Creating a mesh with different (sub)-domains

bhaveshshrimali
Topic Author
Offline
Junior Member

4 years 3 weeks ago - 4 years 3 weeks ago #3274 by bhaveshshrimali

Creating a mesh with different (sub)-domains was created by bhaveshshrimali

Hi everyone,
I want to run a run a elasticity problem with two materials -- a matrix and two spherical solids embedded in it. I need a refined mesh between the spheres as I plan to extend the code to do phase-field fracture. For getting an idea of the geometry, consider something like the attached figure:

[attachment=undefined]trialGeometry.png[/attachment] [attachment=undefined]trialGeometry.jpg[/attachment] [attachment=undefined]trialGeometry.jpg[/attachment]

The problem is that I am getting stuck at the mesh generation part itself. I am using NGSolve-6.2.1909

Consider the following snippet:

Code:

import numpy as np
from ngsolve import *
from netgen.csg import *
Do = 2.286e-3
Ho = 0.21 * Do
H1 = Do - Ho/2.
lo, bo, to = 4*Do, 2*Do, 2*Do
eps = 5.e-6 # do two refinements
left = Plane(Pnt(0,0,0), Vec(-1,0,0))
bottom = Plane(Pnt(0,0,0), Vec(0,-1,0))
back = Plane(Pnt(0,0,0), Vec(0,0,-1))
right = Plane(Pnt(lo,0,0), Vec(1,0,0))
top = Plane(Pnt(0,bo,0), Vec(0,1,0))
front = Plane(Pnt(0,0,to), Vec(0,0,1))

box = left * right * bottom * top * back * front

inleft = Plane(Pnt(H1 + Do/2., Do/2., Do/2.), Vec(-1,0,0))
inbottom = Plane(Pnt(H1 + Do/2., Do/2., Do/2.), Vec(0,-1,0))
inback = Plane(Pnt(H1 + Do/2., Do/2., Do/2.), Vec(0,0,-1))
inright = Plane(Pnt(H1 + 1.5 * Do + Ho, Do/2., Do/2.), Vec(1,0,0))
intop = Plane(Pnt(H1 + Do/2., 3.*Do/2., Do/2.), Vec(0,1,0))
infront = Plane(Pnt(H1 + Do/2., Do/2., 3.*Do/2.), Vec(0,0,1))

inbox = inleft * inright * inbottom * intop * inback * infront

spheres = [
    Sphere(Pnt(H1 + Do/2., Do, Do), Do/2.), 
    Sphere(Pnt(H1 +1.5 * Do + Ho, Do, Do), Do/2.),
]
Spheres = spheres[0] + spheres[1]
geo = CSGeometry()
matrix = box - Spheres
innerpdms = inbox - Spheres
outerpdms = matrix - innerpdms
# onlyspheres = box * Spheres
geo.Add(outerpdms)
geo.Add(Spheres)
geo.Add(innerpdms)
outerpdms.maxh(10*eps)
innerpdms.maxh(eps/3.)
Spheres.maxh(10*eps)
msh = geo.GenerateMesh(grading=0.8)

But this simply stalls at:

Code:

importing NGSolve-6.2.1909
 Start Findpoints
 main-solids: 3 
 Found points 24    
 Analyze spec points
 Find edges

Would anyone have a clue as to what could be going wrong?

Thanks
Bhavesh

Last edit: 4 years 3 weeks ago by bhaveshshrimali. Reason: Image not visible

Time to create page: 0.093 seconds