Modelling Spherical Inclusions in 3D
Hi,
I'm trying to generate a mesh representing multiple spherical regions in 3D using gmsh Python api. I have a working example in 2D. It generates 3 concentric disks, with 3 separate regions, corresponding to "source", "vacuum" and "wall". Each region is labelled for later use in FEniCS. Finally, the mesh is refined, such that it is more fine near the boundaries of each region. Here is the code in 2D:
import gmsh
gmsh.initialize()
gmsh.model.add("2D_mesh")
geom = gmsh.model.occ
'Geometry settings------------------------------------------------------------'
wt = 0.1 # The width of the wall
r = 0.3 # The radius of the spherical source
'Defining the source, vacuum and the wall regions:----------------------------'
# Define a surface for the source. (Also creates a boundary labelled 1.)
source_boundary_line = geom.addCircle(x=0, y=0, z=0, r=r)
# Define chamber walls and vacuum.
inner_wall = geom.addCircle(x=0, y=0, z=0, r=1)
outer_wall = geom.addCircle(x=0, y=0, z=0, r=1+wt)
source_boundary = geom.addCurveLoop([source_boundary_line])
inner_wall_boundary = geom.addCurveLoop([inner_wall])
outer_wall_boundary = geom.addCurveLoop([outer_wall])
source = geom.addPlaneSurface([source_boundary])
vacuum = geom.addPlaneSurface([source_boundary, inner_wall_boundary])
wall = geom.addPlaneSurface([inner_wall_boundary, outer_wall_boundary])
geom.synchronize()
'Refining the mesh along the different region boundaries: --------------------'
# Create "field" to control the mesh size in different regions of the mesh
gmsh.model.mesh.field.add("Distance", 1)
gmsh.model.mesh.field.setNumbers(1, "CurvesList", [source_boundary, inner_wall_boundary])
gmsh.model.mesh.field.setNumber(1, "NumPointsPerCurve", 100)
gmsh.model.mesh.field.add("Threshold", 2)
gmsh.model.mesh.field.setNumber(2, "InField", 1)
gmsh.model.mesh.field.setNumber(2, "SizeMin", 0.001)
gmsh.model.mesh.field.setNumber(2, "SizeMax", 0.1)
gmsh.model.mesh.field.setNumber(2, "DistMin", 0.0)
gmsh.model.mesh.field.setNumber(2, "DistMax", 0.1)
gmsh.model.mesh.field.add("MathEval", 3)
gmsh.model.mesh.field.setString(3, "F", "0.1*F1 + 0.01")
gmsh.model.mesh.field.add("Min", 4)
gmsh.model.mesh.field.setNumbers(4, "FieldsList", [2]) # to apply field[3] put 3 in list. e.g. [2,3]
gmsh.model.mesh.field.setAsBackgroundMesh(4)
gmsh.option.setNumber("Mesh.MeshSizeExtendFromBoundary", 0)
gmsh.option.setNumber("Mesh.MeshSizeFromPoints", 0)
gmsh.option.setNumber("Mesh.MeshSizeFromCurvature", 0)
'Marking the different regions. This will be later used to define density----'
# Mark the physical domains. This will be later used to
gmsh.model.addPhysicalGroup(dim=2, tags=[source], tag=1)
gmsh.model.addPhysicalGroup(dim=2, tags=[vacuum], tag=2)
gmsh.model.addPhysicalGroup(dim=2, tags=[wall], tag=3)
gmsh.model.addPhysicalGroup(dim=1, tags=[source_boundary], tag=1)
gmsh.model.addPhysicalGroup(dim=1, tags=[inner_wall_boundary], tag=2)
gmsh.model.addPhysicalGroup(dim=1, tags=[outer_wall_boundary], tag=3)
'Generating the mesh:--------------------------------------------------------'
# Generate and save mesh.
gmsh.model.mesh.generate(dim=2)
gmsh.write("mesh.msh")
# Preview mesh.
gmsh.fltk.run()
# Clear mesh and close gmsh API.
gmsh.clear()
gmsh.finalize()
The picture of the generated mesh is attached. My question is how to create an analogous mesh but in 3D, which contains 3 cocentric spheres, with the mesh boundaries refined. Here is my attempt so far. The result (attached) looks pretty close, but I don't know how to correctly do the second boolean cut in order to create the wall region. Somehow the mesh is duplicated between the wall and the source region. And if I attempt the second boolean cut, I get an error saying no elements corresponding to tag 6. Here is the code:
import gmsh
gmsh.initialize()
gmsh.model.add("3D_mesh")
geom = gmsh.model.occ
wt = 0.1 # Wall thickness
r = 0.2 # Source readius
source = geom.addSphere(xc= 0, yc=0, zc=0, radius = r, tag = 1)
inner_wall = geom.addSphere(xc= 0, yc=0, zc=0, radius = 1, tag = 2)
outer_wall = geom.addSphere(xc= 0, yc=0, zc=0, radius = 1+wt, tag = 3)
## Boolean difference (defining the different regions):
geom.cut(objectDimTags = [(3,2)], toolDimTags= [(3,1)], tag = 4 ,removeObject = True, removeTool=False)
#geom.cut(objectDimTags = [(3,3)], toolDimTags= [(3,4)], tag = -1 ,removeObject = True, removeTool=False)
geom.synchronize()
gmsh.model.addPhysicalGroup(dim=3, tags=[1], tag=1)
gmsh.model.addPhysicalGroup(dim=3, tags=[2], tag=2)
#gmsh.model.addPhysicalGroup(dim=3, tags=[3], tag=3)
gmsh.model.addPhysicalGroup(dim=3, tags=[4], tag=4)
gmsh.model.addPhysicalGroup(dim=2, tags=[1], tag=1)
gmsh.model.addPhysicalGroup(dim=2, tags=[2], tag=2)
#gmsh.model.addPhysicalGroup(dim=2, tags=[3], tag=3)
gmsh.model.addPhysicalGroup(dim=2, tags=[4], tag=4)
gmsh.model.addPhysicalGroup(dim=1, tags=[1], tag=1)
gmsh.model.addPhysicalGroup(dim=1, tags=[2], tag=2)
#gmsh.model.addPhysicalGroup(dim=1, tags=[3], tag=3)
gmsh.model.addPhysicalGroup(dim=1, tags=[4], tag=4)
#gmsh.model.addPhysicalGroup(dim=3, tags=[6], tag=4)
gmsh.model.mesh.setTransfiniteCurve(tag = 1, numNodes=40)
gmsh.model.mesh.setTransfiniteCurve(tag = 2, numNodes=100)
gmsh.model.mesh.setTransfiniteCurve(tag = 3, numNodes=70)
#gmsh.model.mesh.setTransfiniteCurve(tag = 6, numNodes=70)
gmsh.model.mesh.generate(dim=3)
gmsh.write("mesh3D_test.msh")
# Preview mesh.
gmsh.fltk.run()
# Clear mesh and close gmsh API.
gmsh.clear()
gmsh.finalize()
Does anyone know what is wrong with my approach? Thank you kindly for your help.