30. mesh — mesh.py

mesh.py

Definition of the Mesh class for describing discrete geometrical models. And some useful meshing functions to create such models.

class mesh.Mesh(coords=None, elems=None, prop=None, eltype=None)

A mesh is a discrete geometrical model defined by nodes and elements.

In the Mesh geometrical data model, coordinates of all points are gathered in a single twodimensional array ‘coords’ with shape (ncoords,3) and the individual geometrical elements are described by indices into the ‘elems’ array. This model has some advantages over the Formex data model, where all points of all element are stored by their coordinates:

  • compacter storage, because coordinates of coinciding points do not need to be repeated,
  • faster connectivity related algorithms.

The downside is that geometry generating algorithms are far more complex and possibly slower.

In pyFormex we therefore mostly use the Formex data model when creating geometry, but when we come to the point of exporting the geometry to file (and to other programs), a Mesh data model may be more adequate.

The Mesh data model has at least the following attributes:

  • coords: (ncoords,3) shaped Coords array,
  • elems: (nelems,nplex) shaped array of int32 indices into coords. All values should be in the range 0 <= value < ncoords.
  • prop: array of element property numbers, default None.
  • eltype: string designing the element type, default None.

If eltype is None, a default eltype is deived from the plexitude.

A Mesh can be initialized by its attributes (coords,elems,prop,eltype) or by a single geometric object that provides a toMesh() method.

Methods

Initialize a new Mesh.

Methods

Mesh objects have the following methods:

setCoords(coords)

Replace the current coords with new ones.

Returns a Mesh exactly like the current except for the position of the coordinates.

setProp(prop=None)

Create or destroy the property array for the Mesh.

A property array is a rank-1 integer array with dimension equal to the number of elements in the Mesh. You can specify a single value or a list/array of integer values. If the number of passed values is less than the number of elements, they wil be repeated. If you give more, they will be ignored.

If a value None is given, the properties are removed from the Mesh.

getProp()
Return the properties as a numpy array (ndarray)
maxProp()
Return the highest property value used, or None
propSet()
Return a list with unique property values.
copy()
Return a copy using the same data arrays
toFormex()

Convert a Mesh to a Formex.

The Formex inherits the element property numbers and eltype from the Mesh. Node property numbers however can not be translated to the Formex data model.

ndim()
nelems()
nplex()
ncoords()
shape()
nedges()

Return the number of edges.

This returns the number of rows that would be in getEdges(), without actually constructing the edges. The edges are not fused!

centroids()

Return the centroids of all elements of the Mesh.

The centroid of an element is the point whose coordinates are the mean values of all points of the element. The return value is a Coords object with nelems points.

getCoords()
Get the coords data.
getElems()
Get the elems data.
getLowerEntitiesSelector(level=1, unique=False)

Get the entities of a lower dimensionality.

If the element type is defined in the elements module, this returns a Connectivity table with the entities of a lower dimensionality. The full list of entities with increasing dimensionality 0,1,2,3 is:

['points', 'edges', 'faces', 'cells' ]

If level is negative, the dimensionality returned is relative to that of the caller. If it is positive, it is taken absolute. Thus, for a Mesh with a 3D element type, getLowerEntities(-1) returns the faces, while for a 2D element type, it returns the edges. For bothe meshes however, getLowerEntities(+1) returns the edges.

By default, all entities for all elements are returned and common entities will appear multiple times. Specifying unique=True will return only the unique ones.

The return value may be an empty table, if the element type does not have the requested entities (e.g. the ‘point’ type). If the eltype is not defined, or the requested entity level is outside the range 0..3, the return value is None.

getLowerEntities(level=1, unique=False)

Get the entities of a lower dimensionality.

If the element type is defined in the elements module, this returns a Connectivity table with the entities of a lower dimensionality. The full list of entities with increasing dimensionality 0,1,2,3 is:

['points', 'edges', 'faces', 'cells' ]

If level is negative, the dimensionality returned is relative to that of the caller. If it is positive, it is taken absolute. Thus, for a Mesh with a 3D element type, getLowerEntities(-1) returns the faces, while for a 2D element type, it returns the edges. For bothe meshes however, getLowerEntities(+1) returns the edges.

By default, all entities for all elements are returned and common entities will appear multiple times. Specifying unique=True will return only the unique ones.

The return value may be an empty table, if the element type does not have the requested entities (e.g. the ‘point’ type). If the eltype is not defined, or the requested entity level is outside the range 0..3, the return value is None.

getEdges(unique=False)

Return the edges of the elements.

This is a convenient function to create a table with the element edges. It is equivalent to `self.getLowerEntities(1,unique)`

getFaces(unique=False)

Return the faces of the elements.

This is a convenient function to create a table with the element faces. It is equivalent to `self.getLowerEntities(2,unique)`

getAngles(angle_spec=Deg)

Returns the angles in Deg or Rad between the edges of a mesh.

The returned angles are shaped as (nelems, n1faces, n1vertices), where n1faces are the number of faces in 1 element and the number of vertices in 1 face.

getBorder()

Return the border of the Mesh.

This returns a Connectivity table with the border of the Mesh. The border entities are of a lower jierarchical level than the mesh itself. This Connectivity can be used together with the Mesh coords to construct a Mesh of the border geometry. See also getBorderMesh.

getBorderMesh()

Return a Mesh with the border elements.

Returns a Mesh representing the border of the Mesh. The new Mesh is of the next lower hierarchical level.

report()
fuse()

Fuse the nodes of a Meshes.

All nodes that are within the tolerance limits of each other are merged into a single node.

The merging operation can be tuned by specifying extra arguments that will be passed to Coords:fuse().

compact()

Remove unconnected nodes and renumber the mesh.

Beware! This function changes the object in place.

select(selected)

Return a mesh with selected elements from the original.

  • selected: an object that can be used as an index in the elems array, e.g. a list of element numbers.

Returns a Mesh with only the selected elements. The returned mesh is not compacted.

meanNodes(nodsel)

Create nodes from the existing nodes of a mesh.

nodsel is a local node selector as in selectNodes() Returns the mean coordinates of the points in the selector.

addNodes(newcoords, eltype=None)

Add new nodes to elements.

newcoords is an (nelems,nnod,3) array of coordinates. Each element thus gets exactly nnod extra points and the result is a Mesh with plexitude self.nplex() + nnod.

addMeanNodes(nodsel, eltype=None)

Add new nodes to elements by averaging existing ones.

nodsel is a local node selector as in selectNodes() Returns a Mesh where the mean coordinates of the points in the selector are added to each element, thus increasing the plexitude by the length of the items in the selector. The new element type should be set to correct value.

selectNodes(nodsel, eltype)

Return a mesh with subsets of the original nodes.

nodsel is an object that can be converted to a 1-dim or 2-dim array. Examples are a tuple of local node numbers, or a list of such tuples all having the same length. Each row of nodsel holds a list of local node numbers that should be retained in the new connectivity table.

withProp(val)

Return a Mesh which holds only the elements with property val.

val is either a single integer, or a list/array of integers. The return value is a Mesh holding all the elements that have the property val, resp. one of the values in val. The returned Mesh inherits the matching properties.

If the Mesh has no properties, a copy with all elements is returned.

splitProp()

Partition aMesh according to its prop values.

Returns a dict with the prop values as keys and the corresponding partitions as values. Each value is a Mesh instance. It the Mesh has no props, an empty dict is returned.

convert(totype)
splitRandom(n)
Split a mesh in n parts, distributing the elements randomly.
convertRandom(choices)
Convert choosing randomly between choices
reduceDegenerate(eltype=None)

Reduce degenerate elements to lower plexitude elements.

This will try to reduce the degenerate elements of the mesh to elements of a lower plexitude. If a target element type is given, only the matching recuce scheme is tried. Else, all the target element types for which a reduce scheme from the Mesh eltype is available, will be tried.

The result is a list of Meshes of which the last one contains the elements that could not be reduced and may be empty. Property numbers propagate to the children.

splitDegenerate(autofix=True)

Split a Mesh in degenerate and non-degenerate elements.

If autofix is True, the degenerate elements will be tested against known degeneration patterns, and the matching elements will be transformed to non-degenerate elements of a lower plexitude.

The return value is a list of Meshes. The first holds the non-degenerate elements of the original Mesh. The last holds the remaining degenerate elements. The intermediate Meshes, if any, hold elements of a lower plexitude than the original. These may still contain degenerate elements.

renumber(order='elems')

Renumber the nodes of a Mesh in the specified order.

order is an index with length equal to the number of nodes. The index specifies the node number that should come at this position. Thus, the order values are the old node numbers on the new node number positions.

order can also be a predefined value that will generate the node index automatically: - ‘elems’: the nodes are number in order of their appearance in the

Mesh connectivity.
extrude(n, step=1., dir=0, autofix=True)

Extrude a Mesh in one of the axes directions.

Returns a new Mesh obtained by extruding the given Mesh over n steps of length step in direction of axis dir. The returned Mesh has double plexitude of the original.

This function is usually used to extrude points into lines, lines into surfaces and surfaces into volumes. By default it will try to fix the connectivity ordering where appropriate. If autofix is switched off, the connectivities are merely stacked, and the user may have to fix it himself.

Currently, this function correctly transforms: point1 to line2, line2 to quad4, tri3 to wedge6, quad4 to hex8.

revolve(n, axis=0, angle=360., around=None, autofix=True)

Revolve a Mesh around an axis.

Returns a new Mesh obtained by revolving the given Mesh over an angle around an axis in n steps, while extruding the mesh from one step to the next.

This function is usually used to extrude points into lines, lines into surfaces and surfaces into volumes. By default it will try to fix the connectivity ordering where appropriate. If autofix is switched off, the connectivities are merely stacked, and the user may have to fix it himself.

Currently, this function correctly transforms: point1 to line2, line2 to quad4, tri3 to wedge6, quad4 to hex8.

sweep(path, autofix=True)

Sweep a mesh along a path, creating an extrusion

Returns a new Mesh obtained by sweeping the given Mesh over a path. The returned Mesh has double plexitude of the original. The operation is similar to the extrude() method, but the path can be any 3D curve.

This function is usually used to extrude points into lines, lines into surfaces and surfaces into volumes. By default it will try to fix the connectivity ordering where appropriate. If autofix is switched off, the connectivities are merely stacked, and the user may have to fix it himself.

Currently, this function correctly transforms: point1 to line2, line2 to quad4, tri3 to wedge6, quad4 to hex8.

classmethod concatenate(clas, meshes)

Concatenate a list of meshes of the same plexitude and eltype

Merging of the nodes can be tuned by specifying extra arguments that will be passed to Coords:fuse().

test(nodes='all', dir=0, min=None, max=None, atol=0.)

Flag elements having nodal coordinates between min and max.

This function is very convenient in clipping a TriSurface in a specified direction. It returns a 1D integer array flagging (with a value 1 or True) the elements having nodal coordinates in the required range. Use where(result) to get a list of element numbers passing the test. Or directly use clip() or cclip() to create the clipped TriSurface

The test plane can be defined in two ways, depending on the value of dir. If dir == 0, 1 or 2, it specifies a global axis and min and max are the minimum and maximum values for the coordinates along that axis. Default is the 0 (or x) direction.

Else, dir should be compaitble with a (3,) shaped array and specifies the direction of the normal on the planes. In this case, min and max are points and should also evaluate to (3,) shaped arrays.

nodes specifies which nodes are taken into account in the comparisons. It should be one of the following:

  • a single (integer) point number (< the number of points in the Formex)
  • a list of point numbers
  • one of the special strings: ‘all’, ‘any’, ‘none’

The default (‘all’) will flag all the elements that have all their nodes between the planes x=min and x=max, i.e. the elements that fall completely between these planes. One of the two clipping planes may be left unspecified.

clip(t)

Return a TriSurface with all the elements where t>0.

t should be a 1-D integer array with length equal to the number of elements of the TriSurface. The resulting TriSurface will contain all elements where t > 0.

cclip(t)
This is the complement of clip, returning a TriSurface where t<=0.
clipAtPlane(p, n, nodes='any', side='+')

Return the Mesh clipped at plane (p,n).

This is a convenience function returning the part of the Mesh at one side of the plane (p,n)

equiAngleSkew()

Returns the equiAngleSkew of the elements, a mesh quality parameter .

It quantifies the skewness of the elements: normalize difference between the worst angle in each element and the ideal angle (angle in the face of an equiangular element, qe).

Functions defined in the module mesh

mesh.vectorRotation(vec1, vec2, upvec=[, 0., 0., 1.])

Return a rotation matrix for rotating vector vec1 to vec2

The rotation matrix will be such that the plane of vec2 and the rotated upvec will be parallel to the original upvec.

This function is like arraytools.rotMatrix(), but allows the specification of vec1. The returned matrix should be used in postmultiplication to the Coords.

mesh.sweepCoords(path, origin=[, 0., 0., 0.], normal=0, upvector=2, avgdir=False, enddir=None, scalex=None, scaley=None)

Sweep a Coords object along a path, returning a series of copies.

origin and normal define the local path position and direction on the mesh.

At each point of the curve, a copy of the Coords object is created, with its origin in the curve’s point, and its normal along the curve’s direction. In case of a PolyLine, directions are pointing to the next point by default. If avgdir==True, average directions are taken at the intermediate points. Missing end directions can explicitely be set by enddir, and are by default taken along the last segment. If the curve is closed, endpoints are treated as any intermediate point, and the user should normally not specify enddir.

At each point of the curve, the original Coords object can be scaled in x and y direction by specifying scalex and scaley. The number of values specified in scalex and scaly should be equal to the number of points on the curve.

The return value is a sequence of the transformed Coords objects.

mesh.defaultEltype(nplex)
Default element type for a mesh with given plexitude.
mesh.mergeNodes(nodes)

Merge all the nodes of a list of node sets.

Each item in nodes is a Coords array. The return value is a tuple with:

  • the coordinates of all unique nodes,
  • a list of indices translating the old node numbers to the new.

The merging operation can be tuned by specifying extra arguments that will be passed to Coords:fuse().

mesh.mergeMeshes(meshes)

Merge all the nodes of a list of Meshes.

Each item in meshes is a Mesh instance. The return value is a tuple with:

  • the coordinates of all unique nodes,
  • a list of elems corresponding to the input list, but with numbers referring to the new coordinates.

The merging operation can be tuned by specifying extra arguments that will be passed to Coords:fuse().

mesh.connectMesh(mesh1, mesh2, n=1, n1=None, n2=None, eltype=None)

Connect two meshes to form a hypermesh.

mesh1 and mesh2 are two meshes with same topology (shape). The two meshes are connected by a higher order mesh with n elements in the direction between the two meshes. n1 and n2 are node selection indices permitting a permutation of the nodes of the base sets in their appearance in the hypermesh. This can e.g. be used to achieve circular numbering of the hypermesh.

mesh.connectMeshSequence(ML, loop=False)
mesh.structuredHexGrid(dx, dy, dz, isophex='hex64')
it builds a structured hexahedral grid with nodes and elements both numbered in a structured way: first along z, then along y,and then along x. The resulting hex cells are oriented along z. This function is the equivalent of simple.rectangularGrid but for a mesh. Additionally, dx,dy,dz can be either integers or div (1D list or array). In case of list/array, first and last numbers should be 0.0 and 1.0 if the desired grid has to be inside the region 0.,0.,0. to 1.,1.,1. If isopHex is specified, a convenient set of control points for the isoparametric transformation hex64 is also returned. TODO: include other options to get the control points for other isoparametric transformation for hex.
mesh.correctHexMeshOrientation(hm)
hexahedral elements have an orientation. Some geometrical transformation (e.g. reflect) may produce inconsistent orientation, which results in negative (signed) volume of the hexahedral (triple product). This function fixes the hexahedrals without orientation.

Documentation

Previous topic

29. curve — Definition of curves in pyFormex.

Next topic

31. surface — Operations on triangulated surfaces.

This Page