Naming patch on a volume

[albatta674]

Dear all,

I am making a cylinder with a round patch on the side. I would like to have boundary names for the top, bottom, side of the cylinder, and for the patch. However, I do not get the side boundary. Here the example:

from netgen.csg import *
from ngsolve import Mesh
from ngsolve.webgui import Draw
from math import cos, sin

r = 0.5
h = 1.0

# patch parameters:
z0 = 0.6           # height of patch center
theta0 = 0.0       # angle of patch center
rp = 0.18          # patch size parameter

# point on cylinder wall where the patch is centered
cx = r * cos(theta0)
cy = r * sin(theta0)
cz = z0

# infinite cylinder side
side_all = Cylinder(Pnt(0,0,0), Pnt(0,0,h), r)

# top and bottom caps
bottom = Plane(Pnt(0,0,0), Vec(0,0,-1)).bc("bottom")
top    = Plane(Pnt(0,0,h), Vec(0,0, 1)).bc("top")

# selector for the side patch
selector = Sphere(Pnt(cx, cy, cz), rp)

# split the cylinder wall into patch and the rest
side_patch = (side_all * selector).bc("patch")
side_rest  = (side_all - selector).bc("side")

# reconstruct full side wall, then cap with top/bottom planes
solid = (side_patch + side_rest) * bottom * top

geo = CSGeometry()

geo.Add(solid)

mesh = Mesh(geo.GenerateMesh(maxh=0.08))
print(mesh.GetBoundaries())

Draw(mesh)

Best,

Alberto

[christopher]

I would recommend you to switching to the newer opencascade interface. Such things had to be done with boundary modifiers, see for example here:

Assigning boundary name

Leonhard is referencing to the netgen.occ interface. I’d also advise to move to this more powerful interface to the opencascade geometry kernel. However in the old csg kernel there are some caveats with boundary conditions: probably there is the same plane somewhere else defined in your geometry without the bc name. Csg then unifies them. To name the boundary save the plane in Python (for example as pl1) and then when adding to the Geo use bcmod: geo2.Add(test, bcmod=[(pl1, "plane")]) This be…

But in the occ interface it is easier.
Best
Christopher

[albatta674]

Thanks for the hint!
here it is an example made with the opencascade geometry kernel netgen.occ. Just for reference (to use in jupyter notebook).

from netgen.occ import *
from ngsolve import Mesh
from math import cos, sin
from ngsolve.webgui import Draw

r = 1.
h = 1.0

# main cylinder
cyl = Cylinder((0,0,0), Z, r=r, h=h)

# patch center on the cylindrical wall
theta0 = 1.0
z0 = 0.6
rp = 0.18   # selector radius

cx = r*cos(theta0)
cy = r*sin(theta0)
cz = z0

# selector volume around the patch location
tools = []
tools.append(Sphere((cx, cy, cz), rp))

# patch center on the cylindrical wall
theta0 = 2.
z0 = 0.5
rp = 0.15   # selector radius

cx = r*cos(theta0)
cy = r*sin(theta0)
cz = z0

# selector volume around the patch location
tools.append(Sphere((cx, cy, cz), rp))


# make bulk (cylinder) and patches
bulk = cyl
patches = []
for ii, tool in enumerate(tools):
    bulk -= tool
    patch = cyl * tool
    patch.faces.name = f"patch_{ii}"
    patch.faces.col = (1, 0, 0)
    patches.append(patch) 

bulk.faces.name  = "side"


# fuse them back together into one geometry
shape = Fuse([bulk, *patches])

# overwrite cap names
shape.faces.Min(Z).name = "bottom"
shape.faces.Max(Z).name = "top"

geo = OCCGeometry(shape)
mesh = Mesh(geo.GenerateMesh(maxh=0.08)).Curve(3)

print(mesh.GetBoundaries())
print(mesh.GetMaterials())

Draw(mesh)