Discretisation

Mesh

The computational mesh class that provides the spatial discretisation for finite element calculations.

class underworld3.discretisation.Mesh

Bases: Stateful, uw_object

Unstructured mesh with PETSc DMPlex backend.

The Mesh class provides the spatial discretisation for finite element computations. It wraps PETSc’s DMPlex for unstructured mesh management, supporting various cell types (triangles, quadrilaterals, tetrahedra, hexahedra) and coordinate systems.

Parameters:

• plex_or_meshfile (PETSc.DMPlex or str) – Either a PETSc DMPlex object or path to a mesh file (gmsh, exodus).

• degree (int, optional) – Polynomial degree for the coordinate field (default 1).

• simplex (bool, optional) – True for simplicial elements (triangles/tets), False for quads/hexes.

• coordinate_system_type (CoordinateSystemType, optional) – Coordinate system for vector calculus (Cartesian, cylindrical, etc.).

• qdegree (int, optional) – Quadrature degree for numerical integration (default 2).

• boundaries (list of NamedTuple, optional) – Boundary region definitions with names and values.

• boundary_normals (dict, optional) – Outward normal vectors for each boundary.

• units (str or pint.Unit, optional) – Physical units for mesh coordinates.

• verbose (bool, optional) – Print mesh construction information.

Examples

Meshes are typically created via the meshing module:

>>> mesh = uw.meshing.UnstructuredSimplexBox(
...     minCoords=(0, 0), maxCoords=(1, 1), cellSize=0.1
... )
>>> T = mesh.add_variable("T", vtype=uw.VarType.SCALAR)

See also

underworld3.meshing

Mesh generation utilities.

underworld3.discretisation.MeshVariable

Field variables on meshes.

mesh_instances = 0

__init__(plex_or_meshfile, degree=1, simplex=True, coordinate_system_type=None, qdegree=2, markVertices=None, useRegions=None, useMultipleTags=None, filename=None, refinement=None, refinement_callback=None, coarsening=None, coarsening_callback=None, return_coords_to_bounds=None, boundaries=None, boundary_normals=None, name=None, units=None, verbose=False, *args, **kwargs)

property dim: int

Topological dimension of the mesh.

Returns:

The mesh dimension (2 for 2D, 3 for 3D).

Return type:

int

property cdim: int

Coordinate dimension (embedding space dimension).

For most meshes, cdim == dim. For surface meshes embedded in 3D (e.g., a 2D spherical shell), dim=2 but cdim=3.

Returns:

The coordinate dimension.

Return type:

int

property element: dict

Element type information for the mesh.

Contains details about the finite element discretization including cell type, polynomial degree, and quadrature order.

Returns:

Element information dictionary.

Return type:

dict

Notes

UW3 does not support mixed-element meshes; this applies uniformly to all cells.

property length_scale: float

Length scale for non-dimensionalization.

This property is IMMUTABLE after mesh creation to ensure synchronization with all spatial operators (gradient, divergence, curl, etc.).

The length scale is derived from model reference quantities at mesh creation: - Priority 1: domain_depth from model.set_reference_quantities() - Priority 2: length from model.set_reference_quantities() - Default: 1.0 (no scaling)

Returns:

Length scale value for non-dimensionalization

Return type:

float

Examples

>>> model.set_reference_quantities(domain_depth=uw.quantity(100, "km"))
>>> mesh = uw.meshing.UnstructuredSimplexBox(...)
>>> mesh.length_scale
100000.0  # meters

See also

length_units

Units string for length scale

property length_units: str

Unit string for the length scale.

Returns:

Units for the length scale (e.g., “meter”, “kilometer”)

Return type:

str

Examples

>>> mesh.length_units
'kilometer'

view(level=0)

Displays mesh information at different levels.

Parameters:

level (int (0 default)) – The display level. 0, for basic mesh information (variables and boundaries), while level=1 displays detailed mesh information (including PETSc information)

view_parallel()

returns the break down of boundary labels from each processor

clone_dm_hierarchy()

Clone the dm hierarchy on the mesh

nuke_coords_and_rebuild(verbose=False)

update_lvec()

This method creates and/or updates the mesh variable local vector. If the local vector is already up to date, this method will do nothing.

property lvec: Vec

Returns a local Petsc vector containing the flattened array of all the mesh variables.

access(*writeable_vars)

Dummy access manager that provides deferred sync for backward compatibility. Uses NDArray_With_Callback.delay_callbacks_global() internally.

This is a compatibility wrapper that allows existing code using the access() context manager to work with the new direct-access variable interfaces. All variable modifications are deferred and synchronized at context exit.

Parameters:

writeable_vars (MeshVariable) – Variables that will be modified (ignored - all variables are writable with the new interface, this parameter is kept for API compatibility)

Return type:

Context manager that defers variable synchronization until exit

Notes

This method is deprecated. New code should access variable.data or variable.array directly without requiring an access context.

property N: CoordSys3D

SymPy coordinate system for symbolic calculus.

The base coordinate system used for gradient, divergence, and curl operations. Access base scalars via mesh.N.x, mesh.N.y, mesh.N.z and base vectors via mesh.N.i, mesh.N.j, mesh.N.k.

Returns:

The SymPy coordinate system object.

Return type:

sympy.vector.CoordSys3D

See also

X

Coordinate system with data access.

r

Tuple of coordinate scalars.

property Gamma_N: CoordSys3D

SymPy coordinate system for boundary/surface coordinates.

Returns:

The boundary coordinate system object.

Return type:

sympy.vector.CoordSys3D

property Gamma: CoordSys3D

Boundary coordinate scalars as a row matrix.

Returns:

Row matrix of boundary coordinate scalars.

Return type:

sympy.Matrix

property X

Coordinate system with symbolic coordinates and data access.

The primary interface for mesh coordinates, providing both symbolic expressions for equations and numerical data for evaluation.

Returns:

Coordinate system object with:

• mesh.X[0], mesh.X[1]: Symbolic coordinate functions

• mesh.X.coords: Coordinate data array (vertex positions)

• mesh.X.units: Coordinate units

• x, y = mesh.X: Unpack symbolic coordinates

Return type:

CoordinateSystem

Examples

>>> x, y = mesh.X  # Symbolic coordinates for equations
>>> coords = mesh.X.coords  # Numerical vertex positions

See also

N

SymPy coordinate system for vector calculus.

property CoordinateSystem: CoordinateSystem

Alias for X (the coordinate system object).

property r: Tuple[BaseScalar]

Tuple of coordinate scalars \((x, y)\) or \((x, y, z)\).

Returns:

Tuple of SymPy base scalars (N.x, N.y[, N.z]).

Return type:

tuple

See also

rvec

Position vector form.

property rvec: Vector

Position vector \(\mathbf{r} = x\hat{i} + y\hat{j} [+ z\hat{k}]\).

Returns:

The position vector in the mesh coordinate system.

Return type:

sympy.vector.Vector

property data: ndarray

The array of mesh element vertex coordinates.

Deprecated since version 0.99.0: Use X.coords instead. mesh.data is deprecated in favor of mesh.X.coords (coordinate-system-aware interface).

This is an alias for mesh.points (which is also deprecated).

property points

Mesh node coordinates in physical units.

Deprecated since version 0.99.0: Use X.coords instead. mesh.points is deprecated in favor of mesh.X.coords (coordinate-system-aware interface).

When the mesh has coordinate scaling applied (via model units), this property automatically converts from internal model coordinates to physical coordinates for user access.

When the mesh has coordinate units specified, returns a unit-aware array.

Returns:

numpy.ndarray or UnitAwareArray: Node coordinates (with units if specified)

property physical_coordinates

Mesh coordinates in physical units.

Returns the mesh coordinate array scaled to physical units using the model’s length scale. Requires the mesh to be associated with a model that has reference quantities set.

Returns:

Coordinates in physical units, or None if no model scaling available

Return type:

UWQuantity or None

Examples

>>> model.set_reference_quantities(domain_length=1000*uw.units.km, ...)
>>> mesh = uw.meshing.StructuredQuadBox(...)
>>> physical_coords = mesh.physical_coordinates  # In kilometers

property physical_bounds

Mesh bounds in physical units.

Returns the mesh bounding box scaled to physical units using the model’s length scale.

Returns:

(min_coords, max_coords) in physical units, or None if no model scaling

Return type:

tuple of UWQuantity or None

Examples

>>> physical_min, physical_max = mesh.physical_bounds
>>> print(f"Domain: {physical_min} to {physical_max}")

property physical_extent

Mesh spatial extent in physical units.

Returns the mesh size (max - min) in each dimension scaled to physical units.

Returns:

Extent in physical units, or None if no model scaling

Return type:

UWQuantity or None

Examples

>>> extent = mesh.physical_extent
>>> print(f"Domain size: {extent}")

write_timestep(filename, index, outputPath='', meshVars=[], swarmVars=[], meshUpdates=False, create_xdmf=True)

Write mesh and selected variables for visualisation output.

This writes: - one mesh HDF5 file (shared/static or per-step, depending on meshUpdates) - one HDF5 file per mesh variable - optional proxy files for swarm variables - optional XDMF file linking all output files

When create_xdmf=True (the default), variable files also include ParaView-compatible groups (/vertex_fields or /cell_fields), and an XDMF file is generated on rank 0.

Parameters:

• filename (str)

• index (int)

• outputPath (str | None)

• meshVars (list | None)

• swarmVars (list | None)

• meshUpdates (bool)

• create_xdmf (bool)

petsc_save_checkpoint(index, meshVars=[], outputPath='')

Use PETSc to save the mesh and mesh vars in a h5 and xdmf file.

Parameters:

• meshVars (list | None) – List of UW mesh variables to save. If left empty then just the mesh is saved.

• index (int) – An index which might correspond to the timestep or output number (for example).

• outputPath (str | None) – Path to save the data. If left empty it will save the data in the current working directory.

write_checkpoint(filename, meshUpdates=True, meshVars=[], swarmVars=[], index=0, unique_id=False)

Write data in a format that can be restored for restarting the simulation.

The difference between this and the visualisation is 1) the parallel section needs to be stored to reload the data correctly, and 2) the visualisation information (vertex form of fields) is not stored. This routine uses dmplex VectorView and VectorLoad functionality.

Parameters:

• filename (str)

• meshUpdates (bool)

• meshVars (list | None)

• swarmVars (list | None)

• index (int | None)

• unique_id (bool | None)

write(filename, index=None)

Save mesh data to the specified hdf5 file.

Parameters:

• filename (str) – The filename for the mesh checkpoint file.

• index (int | None) – Not yet implemented. An optional index which might correspond to the timestep (for example).

vtk(filename)

Save mesh to the specified file

Parameters:

filename (str)

generate_xdmf(filename)

This method generates an xdmf schema for the specified file.

The filename of the generated file will be the same as the hdf5 file but with the xmf extension.

Parameters:

filename (str) – File name of the checkpointed hdf5 file for which the xdmf schema will be written.

property vars

A list of variables recorded on the mesh.

property block_vars

A list of variables recorded on the mesh.

points_in_domain(points, strict_validation=True)

Determine if the given points lie in this domain. Uses a mesh-boundary skeletonization array to determine whether the point is inside the boundary or outside. If close to the boundary, it checks if points are in a cell.

Parameters:

• points (array-like) – Coordinate array in any physical unit system (will be auto-converted). Plain numbers are assumed to be in model coordinates.

• strict_validation (bool) – Whether to perform strict validation near boundaries

get_closest_cells(coords)

This method uses a kd-tree algorithm to find the closest cells to the provided coords. For a regular mesh, this should be exactly the owning cell, but if the mesh is deformed, this is not guaranteed. Note, the nearest point may not be all that close by - use get_closest_local_cells to filter out points that are (probably) not within any local cell.

Parameters:

coords:

An array of the coordinates for which we wish to determine the closest cells. This should be a 2-dimensional array of shape (n_coords,dim) in any physical unit system (will be auto-converted). Plain numbers are assumed to be in model coordinates.

Returns:

closest_cells:

An array of indices representing the cells closest to the provided coordinates. This will be a 1-dimensional array of shape (n_coords).

Parameters:

coords (ndarray)

Return type:

ndarray

test_if_points_in_cells(points, cells)

Determine if the given points lie in the suggested cells. Uses a mesh skeletonization array to determine whether the point is with the convex polygon / polyhedron defined by a cell.

Exact if applied to a linear mesh, approximate otherwise.

Parameters:

• points (array-like) – Coordinate array in any physical unit system (will be auto-converted)

• cells (array-like) – Cell indices to test

Returns:

Boolean array indicating if points are in cells

Return type:

numpy.ndarray

get_closest_local_cells(coords)

This method uses a kd-tree algorithm to find the closest cells to the provided coords. For a regular mesh, this should be exactly the owning cell, but if the mesh is deformed, this is not guaranteed. Also compares the distance from the cell to the point - if this is larger than the “cell size” then returns -1

Parameters:

coords:

An array of the coordinates for which we wish to determine the closest cells. This should be a 2-dimensional array of shape (n_coords,dim) in any physical unit system (will be auto-converted).

Returns:

closest_cells:

An array of indices representing the cells closest to the provided coordinates. This will be a 1-dimensional array of shape (n_coords).

Parameters:

coords (ndarray)

Return type:

ndarray

get_min_radius_old()

This method returns the global minimum distance from any cell centroid to a face. It wraps to the PETSc DMPlexGetMinRadius routine. The petsc4py equivalent always returns zero.

Return type:

float

get_min_radius()

This method returns the global minimum distance from any cell centroid to a face. It wraps to the PETSc DMPlexGetMinRadius routine. The petsc4py equivalent always returns zero.

Return type:

float

get_max_radius()

This method returns the global maximum distance from any cell centroid to a face.

Return type:

float

stats(uw_function, uw_meshVariable, basis=None)

Returns various norms on the mesh for the provided function.

• size

• mean

• min

• max

• sum

• L2 norm

• rms

NOTE: this currently assumes scalar variables !

meshVariable_mask_from_label(label_name, label_value)

Extract single label value and make a point mask - note: this produces a mask on the mesh points and assumes a 1st order mesh. Cell labels are not respected in this function.

register_swarm(swarm)

Register swarm as dependent on this mesh for coordinate change notifications

unregister_swarm(swarm)

Unregister swarm (called during swarm cleanup)

register_surface(surface)

Register surface as dependent on this mesh for adaptation notifications.

unregister_surface(surface)

Unregister surface (called during surface cleanup).

adapt(metric_field, verbose=False)

Adapt the mesh discretization based on a metric field.

This method refines or coarsens the mesh in place, automatically transferring all attached MeshVariables, updating Surfaces, and marking Solvers for rebuild on their next solve() call.

Parameters:

• metric_field (MeshVariable) – A scalar MeshVariable containing metric values (1/h² where h is target edge length). Larger values mean finer mesh (smaller elements). Use Surface.refinement_metric() to create this field from distance.

• verbose (bool, optional) – If True, print progress and statistics during adaptation.

Notes

The adaptation uses PETSc’s mesh adaptation with MMG/pragmatic backend.

What happens automatically:

• MeshVariables are interpolated to the new mesh

• Surfaces recompute their distance fields

• Swarms are marked as stale (particle-element associations invalidated)

• Solvers are marked for rebuild (happens lazily on next solve())

Examples

>>> # Define metric from fault distance
>>> metric = uw.discretisation.MeshVariable("H", mesh, 1)
>>> with mesh.access(metric):
...     # Smaller H near fault, larger far away
...     metric.data[:, 0] = 0.01 + 0.09 * fault.distance_from(mesh.data)
>>> mesh.adapt(metric, verbose=True)
>>> stokes.solve()  # Solver rebuilds automatically

MeshVariable

The primary class for field data on meshes.

underworld3.discretisation.MeshVariable

alias of EnhancedMeshVariable

Checkpointing

checkpoint_xdmf

underworld3.discretisation.checkpoint_xdmf(filename, meshUpdates=True, meshVars=[], swarmVars=[], index=0)

Parameters:

• filename (str)

• meshUpdates (bool)

• meshVars (list | None)

• swarmVars (list | None)

• index (int | None)