meshplex — Simplex meshes for Python

meshplex computes all sorts of interesting points, areas, and volumes in triangular and tetrahedral meshes, with a focus on efficiency. Useful in many contexts, e.g., finite-element and finite-volume computations.

For a quickstart, checkout meshplex’s GitHubPage.

Overview of classes and functions

class meshplex.MeshLine(points, cells)

Class for handling line segment “meshes”.

class meshplex.MeshTetra(points, cells, sort_cells=False)

Class for handling tetrahedral meshes.

plot_edge(edge_id)

Displays edge with ce_ratio.

Parameters

edge_id (int) – Edge ID for which to show the ce_ratio.

property q_min_sin_dihedral_angles

Get the smallest of the sines of the 6 angles between the faces of each tetrahedron, times a scaling factor that makes sure the value is 1 for the equilateral tetrahedron.

property q_vol_rms_edgelength3

For each cell, return the ratio of the volume and the cube of the root-mean-square edge length. (This is cell quality measure used by Stellar <https://people.eecs.berkeley.edu/~jrs/stellar>.)

class meshplex.MeshTri(points, cells, sort_cells=False)

Class for handling triangular meshes.

property angles

All angles in the triangle.

compute_ncurl(vector_field)

Computes the n.dot.curl of a vector field over the mesh. While the vector field is point-based, the curl will be cell-based. The approximation is based on

\[n\cdot curl(F) = \lim_{A\to 0} |A|^{-1} \rangle\int_{dGamma}, F\rangle dr;\]

see https://en.wikipedia.org/wiki/Curl_(mathematics). Actually, to approximate the integral, one would only need the projection of the vector field onto the edges at the midpoint of the edges.

flip_until_delaunay(tol=0.0, max_steps=100)

Flip edges until the mesh is fully Delaunay (up to tol).

plot(show_coedges=True, control_volume_centroid_color=None, mesh_color='k', nondelaunay_edge_color=None, boundary_edge_color=None, comesh_color=(0.8, 0.8, 0.8), show_axes=True, cell_quality_coloring=None, show_point_numbers=False, show_edge_numbers=False, show_cell_numbers=False, cell_mask=None, mark_points=None, mark_edges=None, mark_cells=None)

Show the mesh using matplotlib.

plot_vertex(point_id, show_ce_ratio=True)

Plot the vicinity of a point and its covolume/edgelength ratio.

Parameters
  • point_id (int) – Node ID of the point to be shown.

  • show_ce_ratio (bool, optional) – If true, shows the ce_ratio of the point, too.

save(filename, *args, **kwargs)

Save the mesh to a file.

show(*args, fullscreen=False, **kwargs)

Show the mesh (see plot()).

show_vertex(*args, **kwargs)

Show the mesh around a vertex (see plot_vertex()).

meshplex.from_meshio(mesh)

Transform from meshio to meshplex format.

Parameters

mesh (meshio.Mesh) – The meshio mesh object.

Returns mesh{2,3}d

The mesh data.

meshplex.get_signed_simplex_volumes(cells, pts)

Signed volume of a simplex in nD. Note that signing only makes sense for n-simplices in R^n.

meshplex.read(filename)

Reads an unstructured mesh into meshplex format.

Parameters

filenames (str) – The files to read from.

Returns mesh{2,3}d

The mesh data.