# 3D mesh with labeled boundary components

I would like to solve a problem with mixed boundary conditions in 3D. I want to label different subsets of the boundary to have different boundary conditions on each component. In the following code I tried to label the boundaries of the upper and lower hemispheres of a sphere:

from ngsolve import *
from netgen.csg import *
R=1
sph=Sphere(Pnt(0,0,0),R)
top=Plane(Pnt(0,0,0),Vec(0,0,1))
bot=Plane(Pnt(0,0,0),Vec(0,0,-1))
geo=CSGeometry()
mesh=Mesh(geo.GenerateMesh(maxh=.2))
Draw(mesh)
print("Boundaries: ",mesh.GetBoundaries())
print("Materials: ",mesh.GetMaterials())

The output is
Boundaries: (‘screen’, ‘default’, ‘screen’)
Materials: (‘s1’, ‘s2’)

As expected there are two volume subdomains, but the entire boundary is labeled “screen” rather than part being “screen” and part “artificial”. How can I label the upper and lower parts of the boundary?

Also I would like to do the same thing for an ellipsoid. Is there a way to create an ellipsoid with CSG or must I use OpenCascade (a new venture for me!).

from ngsolve import *
from netgen.csg import *

R=1
top=Plane(Pnt(0,0,0),Vec(0,0,1))
bot=Plane(Pnt(0,0,0),Vec(0,0,-1))
sph=Sphere(Pnt(0,0,0),R)
sph1=sph*top
sph1.mat(“s1”).bc(“screen”)
sph2=sph*bot
sph2.mat(“s2”).bc(“artificial”)
geo=CSGeometry()
mesh=Mesh(geo.GenerateMesh(maxh=.2))
Draw(mesh)
print("Boundaries: ",mesh.GetBoundaries())
print("Materials: ",mesh.GetMaterials())

Results
Boundaries: (‘screen’, ‘artificial’, ‘screen’)
Materials: (‘s2’, ‘s1’)

Thank you for this suggestion. But your code also labels the surface of the entire sphere as “screen”, and in addition labels the internal surface as “artificial”. In my application I want the upper surface labeled “screen” and lower “artificial” and the internal surface left as “default”.

The trick is to use bcmod: 4.2 Constructive Solid Geometry (CSG) — NGS-Py 6.2.2302 documentation

``````from ngsolve import *
from netgen.csg import *
R=1
sph=Sphere(Pnt(0,0,0),R).bc("screen")
top=Plane(Pnt(0,0,0),Vec(0,0,1))
bot=Plane(Pnt(0,0,0),Vec(0,0,-1))
geo=CSGeometry()
mesh=Mesh(geo.GenerateMesh(maxh=.2))
Draw(mesh)
print("Boundaries: ",mesh.GetBoundaries())
print("Materials: ",mesh.GetMaterials())
``````

But I would rather change to the new opencascade interface for the future, makes defining such geometries easier and more stable, your code with occ geometries:

``````rom ngsolve import *
from netgen.occ import *
R=1
sph=Sphere((0,0,0),R)
sph.faces.name = "sphere"
cut=HalfSpace((0,0,0),(0,0,1))
top = sph-cut
bot = sph*cut
top.faces["sphere"].name = "artificial"
bot.faces["sphere"].name = "screen"
top.solids.name = "s1"
bot.solids.name = "s2"
shape = Glue([top, bot])
geo = OCCGeometry(shape)
mesh=Mesh(geo.GenerateMesh(maxh=.2))
Draw(mesh)
print("Boundaries: ",mesh.GetBoundaries())
print("Materials: ",mesh.GetMaterials())
``````

Best
Christopher

Hi Peter,
with the latest nightly release we have also ellipsoids in netgen-occ: Either with Axis (i.e. point + direction) and two radii, or with Axes (=coordinate system) + 3 radii:

``````ell1 = Ellipsoid (Axis ((0,0,0), Z), 3, 1)
ell2 = Ellipsoid (Axes ((3,3,0), X, Y+Z), 3, 2, 1)
DrawGeo (ell1+ell2);
``````

best, Joachim