This page was generated from unit-4.2-csg/csg.ipynb.
4.2 Constructive Solid Geometry (CSG)¶
These geometries are bases on primitives (e.g. sphere, cylinder, plane) which are used to build solids by performing boolean operations. Netgen offers the following primitives
primitive |
csg syntax |
meaning |
---|---|---|
half-space |
Plane(Pnt a,Vec n) |
point p in plane, normal vector |
sphere |
Sphere(Pnt c,float r) |
sphere with center c and radius r |
cylinder |
Cylinder(Pnt a, Pnt b, float r) |
points a and b define the axes of a infinite cylinder of radius r |
brick |
OrthoBrick (Pnt a, Pnt b) |
axes parallel brick with minimal coordinates a and maximal coordinates b |
and the boolean operators
operator |
set operation |
---|---|
\(*\) |
intersection |
\(+\) |
union |
\(-\) |
intersection with complement |
[1]:
# import netgen.gui
from ngsolve import Draw, Redraw # just for visualization
Using these primitives and operations, we can easily construct a cube. First we import the netgen.csg
module, create 6 plane and intersect them to get the solid cube
.
[2]:
from netgen.csg import *
left = Plane (Pnt(0,0,0), Vec(-1,0,0) )
right = Plane (Pnt(1,1,1), Vec( 1,0,0) )
front = Plane (Pnt(0,0,0), Vec(0,-1,0) )
back = Plane (Pnt(1,1,1), Vec(0, 1,0) )
bot = Plane (Pnt(0,0,0), Vec(0,0,-1) )
top = Plane (Pnt(1,1,1), Vec(0,0, 1) )
cube = left * right * front * back * bot * top
Then we create a CSGeometry
object and add the solid.
[3]:
geo = CSGeometry()
geo.Add (cube)
mesh = geo.GenerateMesh(maxh=0.25)
Redraw()
# mesh.Save("cube.vol")
[3]:
True
[4]:
from netgen.csg import *
cube = OrthoBrick( Pnt(0,0,0), Pnt(1,1,1) )
hole = Cylinder ( Pnt(0.5, 0.5, 0), Pnt(0.5, 0.5, 1), 0.2)
geo = CSGeometry()
geo.Add (cube-hole.maxh(0.05))
mesh = geo.GenerateMesh(maxh=0.1)
Redraw()
[4]:
True
Setting properties of solids¶
A solid has members which we can set to define the desired properties.
[5]:
sphere = Sphere(Pnt(0,0,0),1)
Now we can set a boundary name and a maximal mesh size on the surface of this sphere
[6]:
sphere.bc("sphere").maxh(0.25)
[6]:
<netgen.libngpy._csg.Solid at 0x7f80c01b69b0>
and define a material for the volume
[7]:
sphere.mat("iron")
[7]:
<netgen.libngpy._csg.Solid at 0x7f80c01b69b0>
In case we want to visualize the geometry we can define the color (using rgb values) and transparency of the solid.
[8]:
sphere.col([1,0,0])#.transp()
[8]:
<netgen.libngpy._csg.Solid at 0x7f80c01b69b0>
[9]:
geo = CSGeometry()
geo.Add(sphere)
geo.Draw()
[10]:
ngmesh = geo.GenerateMesh()
print(type(ngmesh))
Redraw()
<class 'netgen.libngpy._meshing.Mesh'>
[10]:
True
To improve the approximation of curved geometries it is possible to use curved elements. This can be done within NGSolve
. Thus we have to convert the Netgen
mesh to a NGSolve
mesh before curving it.
[11]:
from ngsolve.comp import Mesh
mesh = Mesh(ngmesh)
print(type(mesh))
Redraw()
<class 'ngsolve.comp.Mesh'>
[11]:
False
[12]:
mesh.Curve(3)
Draw(mesh)
Setting the mesh size¶
There are the following options to set the mesh size:
globally as argument
maxh
ofGenerateMesh
to the surface of one solid (
maxh
property as above mentioned)for the volume of a solid as optional argument when adding it to the geometry
Add(...,bc)
restrict the mesh size for one point using
RestrictH
use
CloseSurfaces
to generate anisotropic meshes
Global mesh size¶
The global mesh size can be set with the named argument maxh
. The following two versions are equivalent since all arguments of the of the GenerateMesh
function are parsed to the MeshingParameters
if no named argument mp
is given.
[13]:
unit_cube.GenerateMesh(maxh=0.4)
[13]:
<netgen.libngpy._meshing.Mesh at 0x7f80c15b4570>
[14]:
from netgen.meshing import MeshingParameters
mp = MeshingParameters(maxh=0.4)
unit_cube.GenerateMesh(mp = mp)
[14]:
<netgen.libngpy._meshing.Mesh at 0x7f80c1525c70>
Mesh size for one solid¶
To set the mesh size for one domain of the mesh we have to add the desired maxh
as argument when adding the solid to the geometry
[15]:
geo = CSGeometry()
brick = OrthoBrick(Pnt(-2,-2,-2),Pnt(2,2,2))
sphere = Sphere(Pnt(0,0,0),1)
geo.Add(brick-sphere)
geo.Add(sphere,maxh=0.1)
ngmesh = geo.GenerateMesh(maxh=0.4)
Mesh size on a surface¶
If we want to refine just on a surface we define it as property of the solid.
[16]:
geo = CSGeometry()
brick = OrthoBrick(Pnt(-2,-2,-2),Pnt(2,2,2))
sphere = Sphere(Pnt(0,0,0),1)
geo.Add(brick-sphere)
geo.Add(sphere.maxh(0.1))
ngmesh = geo.GenerateMesh()
Mesh size in points¶
This can be done with the MeshingParameters
. Using RestrictH
we can define the mesh size in an arbitrary point.
[17]:
geo = CSGeometry()
brick = OrthoBrick(Pnt(-2,-2,-2),Pnt(2,2,2))
sphere = Sphere(Pnt(0,0,0),1)
mp = MeshingParameters(maxh=0.4)
mp.RestrictH (x=0, y=0, z=1, h=0.025)
geo.Add(brick-sphere)
geo.Add(sphere)
ngmesh = geo.GenerateMesh(mp = mp)
Anisotropic meshes¶
If the geometry contains thin layers we can use CloseSurfaces
to avoid elements with small angles.
[18]:
from netgen.csg import *
geo = CSGeometry()
box = OrthoBrick(Pnt(0,0,0),Pnt(1,1,1))
top = Plane(Pnt(0,0,0.52),Vec(0,0,1))
bot = Plane(Pnt(0,0,0.48),Vec(0,0,-1))
plate = box * top * bot
geo.Add((box-top).mat("air"))
geo.Add(plate.mat("plate"))
geo.Add((box-bot).mat("air"))
slices = [2**(-i) for i in reversed(range(1,6))]
# define the close surfaces
geo.CloseSurfaces(bot,top)#,slices)
nmesh = geo.GenerateMesh(maxh=0.3)
# refine the mesh between the close surfaces
# ZRefinement(nmesh,geo)
Setting boundary conditions¶
Boundary condition on the surface of a solid¶
Setting a boundary condition on the whole surface of solid can be achieved by adding it as property to the solid.
[19]:
brick = OrthoBrick(Pnt(-2,-2,-2),Pnt(2,2,2)).bc('outer')
sphere = Sphere(Pnt(0,0,0),1).bc('sphere')
Modify boundary between two solids¶
This can be done by adding the named argument bcmod
when adding the solid to the geometry. Here we change the boundary condition on the surface between the halfsphere
and the already added box
.
[20]:
halfsphere = sphere * Plane(Pnt(0,0,0),Vec(1,0,0)).bc('plane')
box = brick-sphere
geo = CSGeometry()
geo.Add(box.col([1,0,0]).transp())
geo.Add(halfsphere.col([0,0,1]),bcmod=[(box,"halfsphere")])
geo.Draw()
[21]:
ngmesh = geo.GenerateMesh()
mesh = Mesh(ngmesh)
mesh.GetBoundaries()
[21]:
('outer',
'outer',
'outer',
'outer',
'outer',
'outer',
'sphere',
'halfsphere',
'plane')
[ ]: