I’m trying to mesh unit cubes containing spherically coated, periodically arranged spheres; the geometries are defined via python. I’ve run into a meshing difficulty I do not understand.
In a simple case the unit cube contains just 2 compound spheres, one in its interior and one crossing a surface of the cube. If I’m content with having three parts (“matrix”, “core” and “interphase”, where core and interphase make up the coated sphere), the following Pyton sequence works fine:
geo = CSGeometry()
celsz = 1.000000
generate cell
left = Plane( Pnt(0.00,0.00,0.00), Vec(-1.0,0.0,0.0) )
front = Plane( Pnt(0.00,0.00,0.00), Vec(0.0,-1.0,0.0) )
bot = Plane( Pnt(0.00,0.00,0.00), Vec(0.0,0.0,-1.0) )
right = Plane( Pnt(celsz,celsz,celsz), Vec(1.0,0.0,0.0) )
back = Plane( Pnt(celsz,celsz,celsz), Vec(0.0,1.0,0.0) )
top = Plane( Pnt(celsz,celsz,celsz), Vec(0.0,0.0,1.0) )
VolEl = left * right * bot * top * front * back
center coordinates and radii of spheres
numSph = 2
centers of spheres
spCc = Matrix(numSph,3)
spCc[0,0] = 0.440789
spCc[0,1] = 0.466176
spCc[0,2] = 0.373534
spCc[1,0] = 0.423044
spCc[1,1] = 0.034480
spCc[1,2] = 0.948371
radii of outer spheres
spRdP = Vector(numSph)
spRdP[0] = 0.212157
spRdP[1] = 0.212157
radii of inner spheres
spRdL = Vector(numSph)
spRdL[0] = 0.180000
spRdL[1] = 0.180000
function for duplicating spheres for periodicity
def DuplicateSphere (c,r):
spheres = None
for x in [-1,0,1]:
for y in [-1,0,1]:
for z in [-1,0,1]:
sph = Sphere(c+celsz*Vec(x,y,z),r)
if spheres:
spheres = spheres + sph
else:
spheres = sph
return spheres
set up compound spheres corresponding to particles
sphere00
cp = Pnt(spCc[0,0],spCc[0,1],spCc[0,2])
sphP00 = DuplicateSphere (cp,spRdP[0])
sphL00 = DuplicateSphere (cp,spRdL[0])
sphrsP = sphP00
sphrsL = sphL00
sphere01
cp = Pnt(spCc[1,0],spCc[1,1],spCc[1,2])
sphP01 = DuplicateSphere (cp,spRdP[1])
sphL01 = DuplicateSphere (cp,spRdL[1])
sphrsP = sphrsP + sphP01
sphrsL = sphrsL + sphL01
#==============================================================================
set up geometry - parts correspond to phases
geo.Add (VolEl-sphrsP)
geo.Add (sphrsL*VolEl)
geo.Add ((sphrsP-sphrsL)*VolEl, maxh=0.030000)
#==============================================================================
geo.PeriodicSurfaces(left, right)
geo.PeriodicSurfaces(top, bot)
geo.PeriodicSurfaces(front, back)
geo.SetBoundingBox(Pnt(-1,-1,-1), Pnt(2,2,2))
mesh = geo.GenerateMesh(maxh=0.0500)
mesh.SecondOrder()
mesh.Save (“test.vol”)
If, however, I try to obtain separate parts for the core and the coating of each sphere by using
NetGen returns an error. If I mesh only sphere00, I get the expected mesh, but if I try meshing only sphere01, NetGen throws a
netgen.libngpy._meshing.NgException: Meshing failed!
Can anybody help me understand what I did wrong?