This page was generated from unit-4.4-occ/occ.ipynb.

4.4 Open Cascade Technology Geometry (NEW)

Netgen provides a Python wrapper around the Open Cascade Technology (OCCT) geometry kernel. It allows to model complex geometric objects. It also allows to import models via STEP format, explore, and modify the geometry. For that geometry type you have to install Netgen with the cmake-flag -DUSE_OCC=ON.

In the wrapper we aim in bringing most of the structure of OCCT to Python. If you are familiar with the C++ interface of Open Cascade you will recognize many classes. Here you can find the OCCT-Netgen dialect to define the famous OCCT-bottle tutorial.

Construction of 3D objects

[1]:
from netgen.occ import *
from netgen.webgui import Draw as DrawGeo

We define a box aligned with the Cartesian coordinates given by two points with minimal and maximal x/y/z coordinates. A cylinder is given by a point on the axis, a direction vector, the radius \(r\), and the height \(h\). The symbols X, Y, and Z are predefined basis vectors for the Cartesian coordinates.

[2]:
box = Box(Pnt(0,0,0), Pnt(1,1,1))
cyl = Cylinder(Pnt(1,0.5,0.5), X, r=0.3, h=0.5)

The Boolean operations fuse, common, and cut provided by OCCT are made available by the operators +, * and -. More can be found in the OCCT-documentation. Note that the fuse generates one new solid, there is no interface face where the cylinder is touching the box.

[3]:
fused = box+cyl
DrawGeo (fused);

The generated objects are Py-wrapped OCCT objects derived from TopoDS_Shape.

[4]:
print ("object type is:", type(box))
print ("type query tells:", box.type)
object type is: <class 'netgen.libngpy._NgOCC.Solid'>
type query tells: TopAbs_ShapeEnum.SOLID

For mesh generation, we make a Netgen OCC-Geometry from the OCC-shape. Then we can call the Netgen mesh generation as usual:

[5]:
from ngsolve import Mesh
from ngsolve.webgui import Draw

geo = OCCGeometry(fused)
mesh = Mesh(geo.GenerateMesh(maxh=0.2))
mesh.Curve(3)
Draw (mesh, clipping=True);

Instead of fusing, we can glue together shapes. Then, the resulting composite solid contains the interface face between the solids. This is important when dealing with separate material regions:

[6]:
geo = Glue( [box, cyl])
DrawGeo (geo)

mesh = Mesh(OCCGeometry(geo).GenerateMesh(maxh=0.2)).Curve(3)
Draw (mesh, clipping=True);

A third option is to form a compound of shapes, then the component shapes are meshed independently:

[7]:
geo = Compound( [box, cyl])
DrawGeo (geo)

mesh = Mesh(OCCGeometry(geo).GenerateMesh(maxh=0.2)).Curve(3)
Draw (mesh, clipping=True);

Transformation of shapes

We can translate, rotate and mirror shapes. The transformation preserves the original shape, and returns a transformed copy. The translation Move takes a vector as argument, the rotation Rotation needs an axis given by a point and a direction, and an angle. The sign of the angle reflects the right hand screw rule.

[8]:
solid = Box((0,0,0), (5,3,1)) + Sphere((0,0,0), 0.3)
solid2 = solid.Move((5,0,2))
solid3 = solid.Move((0,0,4)).Rotate( Axis((0,0,4), X), 45)
DrawGeo (solid + solid2 + solid3);

Hierarchy of shapes

OCCT keeps track of the hierarchy of shapes. A solid knows its faces, a face its edges, and a edge its vertices. We can ask a shape for a list of sub-shapes.

[9]:
solid = Box((0,0,0), (1,1,1))
explodedF = sum( [f.Move(0.2 * (f.center-Pnt(0.5, 0.5, 0.5))) for f in solid.faces] )
explodedE = sum( [e.Move(0.2 * (e.center-Pnt(0.5, 0.5, 0.5))) for e in solid.edges] )

DrawGeo ( Compound( [explodedF, explodedE] ) );
[10]:
for f in solid.faces:
    print ("center of gravity", f.center)
center of gravity (0, 0.5, 0.5)
center of gravity (1, 0.5, 0.5)
center of gravity (0.5, 0, 0.5)
center of gravity (0.5, 1, 0.5)
center of gravity (0.5, 0.5, 0)
center of gravity (0.5, 0.5, 1)
[11]:
e = solid.edges[0]
print ("I am a", e.type)
print ("with startpoint", e.start, "and endpoint", e.end)
I am a TopAbs_ShapeEnum.EDGE
with startpoint (0, 0, 0) and endpoint (0, 0, 1)
[12]:
for f in solid.faces:
    print ("face")
    for e in f.edges:
        print ("edge:", e.start, "-", e.end)
face
edge: (0, 0, 0) - (0, 0, 1)
edge: (0, 0, 1) - (0, 1, 1)
edge: (0, 1, 0) - (0, 1, 1)
edge: (0, 0, 0) - (0, 1, 0)
face
edge: (1, 0, 0) - (1, 0, 1)
edge: (1, 0, 1) - (1, 1, 1)
edge: (1, 1, 0) - (1, 1, 1)
edge: (1, 0, 0) - (1, 1, 0)
face
edge: (0, 0, 0) - (1, 0, 0)
edge: (1, 0, 0) - (1, 0, 1)
edge: (0, 0, 1) - (1, 0, 1)
edge: (0, 0, 0) - (0, 0, 1)
face
edge: (0, 1, 0) - (1, 1, 0)
edge: (1, 1, 0) - (1, 1, 1)
edge: (0, 1, 1) - (1, 1, 1)
edge: (0, 1, 0) - (0, 1, 1)
face
edge: (0, 0, 0) - (0, 1, 0)
edge: (0, 1, 0) - (1, 1, 0)
edge: (1, 0, 0) - (1, 1, 0)
edge: (0, 0, 0) - (1, 0, 0)
face
edge: (0, 0, 1) - (0, 1, 1)
edge: (0, 1, 1) - (1, 1, 1)
edge: (1, 0, 1) - (1, 1, 1)
edge: (0, 0, 1) - (1, 0, 1)

Selection of shapes

Sometimes we want to set properties of shapes, boundary conditions or mesh-size. For that we can use shape selectors.

  • The Max and Min selectors finds the sub-shapes where the center of gravity has maximal or minimal coordinates in a given direction

More of them will come soon …

[13]:
box = Box((0,0,0), (1,1,1))
box.faces.Max(Z).col = (1,0,0)
box.faces.Min(Y).col = (0,0,1)
DrawGeo (box);

Faces where center of gravity has \(x\) coordinate less than 0.8:

[14]:
DrawGeo (Compound(box.faces[X<0.8]));
[ ]: