Common infrastructure

exception meshmode.Error
exception meshmode.DataUnavailable
exception meshmode.FileExistsError

Mesh management

Design of the Data Structure

Why does a Mesh need to be broken into MeshElementGroup instances?

Elements can be of different types (e.g. triangle, quadrilateral, tetrahedron, what have you). In addition, elements may vary in the polynomial degree used to represent them (see also below).

All these bits of information could in principle be stored by element, but having large, internally homogeneous groups is a good thing from an efficiency standpoint. (So that you can, e.g., launch one GPU kernel to deal with all order-3 triangles, instead of maybe having to dispatch based on type and size inside the kernel)

What is the difference between ‘vertices’ and ‘nodes’?

Nodes exist mainly to represent the (potentially non-affine) deformation of each element, by a one-to-one correspondence with MeshElementGroup.unit_nodes. They are unique to each element. Vertices on the other hand exist to clarify whether or not a point shared by two elements is actually identical (or just happens to be “close”). This is done by assigning (single, globally shared) vertex numbers and having elements refer to them.

Consider the following picture:

_images/nodes-vertices.png

Mesh Data Structure

class meshmode.mesh.MeshElementGroup(order, vertex_indices, nodes, element_nr_base=None, node_nr_base=None, unit_nodes=None, dim=None, **kwargs)

A group of elements sharing a common reference element.

order
vertex_indices

An array (nelements, ref_element.nvertices) of (mesh-wide) vertex indices.

nodes

An array of node coordinates with shape (mesh.ambient_dim, nelements, nunit_nodes).

unit_nodes

(dim, nunit_nodes)

element_nr_base

Lowest element number in this element group.

node_nr_base

Lowest node number in this element group.

dim

The number of dimensions spanned by the element. Not the ambient dimension, see Mesh.ambient_dim for that.

is_affine

A bool flag that is True if the local-to-global parametrization of all the elements in the group is affine.

face_vertex_indices()

Return a tuple of tuples indicating which vertices (in mathematically positive ordering) make up each face of an element in this group.

vertex_unit_coordinates()

Return an array of shape (nfaces, dim) with the unit coordinates of each vertex.

nfaces
__eq__(other)

Return self==value.

__ne__(other)

Return self!=value.

class meshmode.mesh.SimplexElementGroup(order, vertex_indices, nodes, element_nr_base=None, node_nr_base=None, unit_nodes=None, dim=None, **kwargs)
class meshmode.mesh.TensorProductElementGroup(order, vertex_indices, nodes, element_nr_base=None, node_nr_base=None, unit_nodes=None, dim=None, **kwargs)
class meshmode.mesh.Mesh(vertices, groups, *, skip_tests=False, node_vertex_consistency_tolerance=None, skip_element_orientation_test=False, nodal_adjacency=None, facial_adjacency_groups=None, boundary_tags=None, vertex_id_dtype=<class 'numpy.int32'>, element_id_dtype=<class 'numpy.int32'>, is_conforming=None)
ambient_dim
dim
vertices

None or an array of vertex coordinates with shape (ambient_dim, nvertices). If None, vertices are not known for this mesh.

nvertices
groups

A list of MeshElementGroup instances.

nelements
nodal_adjacency

An instance of NodalAdjacency.

Referencing this attribute may raise meshmode.DataUnavailable.

facial_adjacency_groups

A list of mappings from neighbor group IDs to instances of FacialAdjacencyGroup.

facial_adjacency_groups[igrp][ineighbor_group] gives the set of facial adjacency relations between group igrp and ineighbor_group. ineighbor_group and igrp may be identical, or ineighbor_group may be None, in which case an InterPartitionAdjacencyGroup group containing boundary faces is returned.

Referencing this attribute may raise meshmode.DataUnavailable.

_images/facial-adjacency-group.png

For example for the mesh in the figure, the following data structure would be present:

[
    {...},  # connectivity for group 0
    {...},  # connectivity for group 1
    {...},  # connectivity for group 2
    {       # connectivity for group 3
        1: FacialAdjacencyGroup(...)  # towards group 1, green
        2: FacialAdjacencyGroup(...)  # towards group 2, pink
        None: FacialAdjacencyGroup(...)  # towards the boundary, orange
    }
]

(Note that element groups are not necessarily geometrically contiguous like the figure may suggest.)

boundary_tags

A list of boundary tag identifiers. BTAG_ALL and BTAG_REALLY_ALL are guaranteed to exist.

btag_to_index

A mapping that maps boundary tag identifiers to their corresponding index.

Note

Elements of boundary_tags that do not cover any part of the boundary will not be keys in this dictionary.

vertex_id_dtype
element_id_dtype
is_conforming

True if it is known that all element interfaces are conforming. False if it is known that some element interfaces are non-conforming. None if neither of the two is known.

copy(**kwargs)
__eq__(other)

Return self==value.

__ne__(other)

Return self!=value.

class meshmode.mesh.NodalAdjacency(valuedict=None, exclude=None, **kwargs)

Describes nodal element adjacency information, i.e. information about elements that touch in at least one point.

neighbors_starts

element_id_t [nelements+1]

Use together with neighbors. neighbors_starts[iel] and neighbors_starts[iel+1] together indicate a ranges of element indices neighbors which are adjacent to iel.

neighbors

element_id_t []

See neighbors_starts.

__eq__(other)

Return self==value.

__ne__(other)

Return self!=value.

class meshmode.mesh.FacialAdjacencyGroup(valuedict=None, exclude=None, **kwargs)

Describes facial element adjacency information for one MeshElementGroup, i.e. information about elements that share (part of) a face.

_images/facial-adjacency-group.png

Represents (for example) one of the (colored) interfaces between MeshElementGroup instances, or an interface between MeshElementGroup and a boundary. (Note that element groups are not necessarily contiguous like the figure may suggest.)

igroup
ineighbor_group

ID of neighboring group, or None for boundary faces. If identical to igroup, then this contains the self-connectivity in this group.

elements

element_id_t [nfagrp_elements]. elements[i] Group-local element numbers.

element_faces

face_id_t [nfagrp_elements]. element_faces[i] indicate what face index of the opposite element indicated in neighbors[iface] touches face number iface of element number iel_grp in this element group.

neighbors

element_id_t [nfagrp_elements]. neighbors[i] gives the element number within ineighbor_group of the element opposite elements[i].

If this number is negative, then this indicates that this is a boundary face, and the bits set in -neighbors[i] should be interpreted according to Mesh.boundary_tags.

neighbor_faces

face_id_t [nfagrp_elements]. neighbor_faces[i] indicate what face index of the opposite element indicated in neighbors[i] has facial contact with face number element_faces[i] of element number elements[i] in element group igroup.

Zero if neighbors[i] is negative.

__eq__(other)

Return self==value.

__ne__(other)

Return self!=value.

class meshmode.mesh.InterPartitionAdjacencyGroup(valuedict=None, exclude=None, **kwargs)

Describes boundary adjacency information of elements in one MeshElementGroup.

igroup

The group number of this group.

ineighbor_group

None for boundary faces.

elements

Group-local element numbers. Element element_id_dtype elements[i] and face face_id_dtype element_faces[i] is connected to neighbor element element_id_dtype partition_neighbors[i] with face face_id_dtype global_neighbor_faces[i]. The partition number it connects to is neighbor_partitions[i].

element_faces

face_id_dtype element_faces[i] gives the face of element_id_dtype elements[i] that is connected to partition_neighbors[i].

neighbors

Since this is a boundary, element_id_dtype neighbors[i] is interpreted as a boundary tag. -neighbors[i] should be interpreted according to :class:Mesh.boundary_tags.

partition_neighbors

element_id_dtype partition_neighbors[i] gives the volume element number within the neighboring partition of the element connected to element_id_dtype elements[i] (which is a boundary element index). Use ~meshmode.mesh.processing.find_group_indices to find the group that the element belongs to, then subtract element_nr_base to find the element of the group.

If partition_neighbors[i] is negative, elements[i] is on a true boundary and is not connected to any other :class:Mesh.

neighbor_faces

face_id_dtype global_neighbor_faces[i] gives face index within the neighboring partition of the face connected to element_id_dtype elements[i]

If neighbor_partitions[i] is negative, elements[i] is on a true boundary and is not connected to any other :class:Mesh.

neighbor_partitions

neighbor_partitions[i] gives the partition number that elements[i] is connected to.

If neighbor_partitions[i] is negative, elements[i] is on a true boundary and is not connected to any other :class:Mesh.

New in version 2017.1.

meshmode.mesh.as_python(mesh, function_name='make_mesh')

Return a snippet of Python code (as a string) that will recreate the mesh given as an input parameter.

meshmode.mesh.check_bc_coverage(mesh, boundary_tags, incomplete_ok=False, true_boundary_only=True)

Verify boundary condition coverage.

Given a list of boundary tags as boundary_tags, this function verifies that

1. the union of all these boundaries gives the complete boundary, 1. all these boundaries are disjoint.

Parameters
  • incomplete_ok – Do not report an error if some faces are not covered by the boundary conditions.

  • true_boundary_only – only verify for faces tagged with BTAG_ALL.

meshmode.mesh.is_boundary_tag_empty(mesh, boundary_tag)

Return True if the corresponding boundary tag does not occur as part of mesh.

Predefined Boundary tags

class meshmode.mesh.BTAG_NONE

A boundary tag representing an empty boundary or volume.

class meshmode.mesh.BTAG_ALL

A boundary tag representing the entire boundary or volume.

In the case of the boundary, BTAG_ALL does not include rank boundaries, or, more generally, anything tagged with BTAG_NO_BOUNDARY.

In the case of a mesh representing an element-wise subset of another, BTAG_ALL does not include boundaries induced by taking the subset. Instead, these boundaries will be tagged with BTAG_INDUCED_BOUNDARY.

class meshmode.mesh.BTAG_REALLY_ALL

A boundary tag representing the entire boundary.

Unlike BTAG_ALL, this includes rank boundaries, or, more generally, everything tagged with BTAG_NO_BOUNDARY.

In the case of a mesh representing an element-wise subset of another, this tag includes boundaries induced by taking the subset, or, more generally, everything tagged with BTAG_INDUCED_BOUNDARY

class meshmode.mesh.BTAG_NO_BOUNDARY

A boundary tag indicating that this edge should not fall under BTAG_ALL. Among other things, this is used to keep rank boundaries out of BTAG_ALL.

class meshmode.mesh.BTAG_PARTITION(part_nr)

A boundary tag indicating that this edge is adjacent to an element of another Mesh. The partition number of the adjacent mesh is given by part_nr.

part_nr

New in version 2017.1.

class meshmode.mesh.BTAG_INDUCED_BOUNDARY

When a Mesh is created as an element-by-element subset of another (as, for example, when using the Firedrake interop features while passing restrict_to_boundary), boundaries may arise where there were none in the original mesh. This boundary tag is used to indicate such boundaries.

Mesh generation

Curves

meshmode.mesh.generation.make_curve_mesh(curve_f: Callable[[numpy.ndarray], numpy.ndarray], element_boundaries: numpy.ndarray, order: int, *, unit_nodes: Optional[numpy.ndarray] = None, node_vertex_consistency_tolerance: Optional[Union[float, bool]] = None, closed: bool = True, return_parametrization_points: bool = False)
Parameters
  • curve_f – parametrization for a curve, accepting a vector of point locations and returning an array of shape (2, npoints).

  • element_boundaries – a vector of element boundary locations in \([0, 1]\), in order. \(0\) must be the first entry, \(1\) the last one.

  • order – order of the (simplex) elements. If unit_nodes is also provided, the orders should match.

  • unit_nodes – if given, the unit nodes to use. Must have shape (2, nnodes).

  • node_vertex_consistency_tolerance – passed to the Mesh constructor. If False, no checks are performed.

  • closed – if True, the curve is assumed closed and the first and last of the element_boundaries must match.

  • return_parametrization_points – if True, the parametrization points at which all the nodes in the mesh were evaluated are also returned.

Returns

a Mesh, or if return_parametrization_points is True, a tuple (mesh, par_points), where par_points is an array of parametrization points.

Curve parametrizations

meshmode.mesh.generation.circle(t: numpy.ndarray)
Parameters

t – the parametrization, runs from \([0, 1)\).

Returns

an array of shape (2, t.size).

meshmode.mesh.generation.ellipse(aspect_ratio: float, t: numpy.ndarray)
Parameters

t – the parametrization, runs from \([0, 1)\).

Returns

an array of shape (2, t.size).

meshmode.mesh.generation.cloverleaf(t: numpy.ndarray)
Parameters

t – the parametrization, runs from \([0, 1)\).

Returns

an array of shape (2, t.size).

meshmode.mesh.generation.starfish
meshmode.mesh.generation.drop(t: numpy.ndarray)
Parameters

t – the parametrization, runs from \([0, 1)\).

Returns

an array of shape (2, t.size).

meshmode.mesh.generation.n_gon(n_corners, t)
Parameters

t – the parametrization, runs from \([0, 1)\).

Returns

an array of shape (2, t.size).

meshmode.mesh.generation.qbx_peanut(t: numpy.ndarray)
Parameters

t – the parametrization, runs from \([0, 1)\).

Returns

an array of shape (2, t.size).

meshmode.mesh.generation.apple(a: float, t: numpy.ndarray)
Parameters
  • a – roundness parameter in \([0, 1/2]\), where \(0\) gives a circle and \(1/2\) gives a cardioid.

  • t – the parametrization, runs from \([0, 1)\).

Returns

an array of shape (2, t.size).

class meshmode.mesh.generation.WobblyCircle(coeffs: numpy.ndarray)
static random(ncoeffs: int, seed: int)
__call__(t: numpy.ndarray)
Parameters

t – the parametrization, runs from \([0, 1)\).

Returns

an array of shape (2, t.size).

class meshmode.mesh.generation.NArmedStarfish(n_arms: int, amplitude: float)

Inherits from WobblyCircle.

__call__(t: numpy.ndarray)
Parameters

t – the parametrization, runs from \([0, 1)\).

Returns

an array of shape (2, t.size).

Surfaces

meshmode.mesh.generation.generate_icosahedron(r: float, order: int, *, node_vertex_consistency_tolerance: Optional[Union[float, bool]] = None, unit_nodes: Optional[numpy.ndarray] = None)
meshmode.mesh.generation.generate_icosphere(r: float, order: int, *, uniform_refinement_rounds: int = 0, node_vertex_consistency_tolerance: Optional[Union[float, bool]] = None, unit_nodes: Optional[numpy.ndarray] = None)
Parameters
  • r – radius of the sphere.

  • order – order of the (simplex) elements. If unit_nodes is also provided, the orders should match.

  • uniform_refinement_rounds – number of uniform refinement rounds to perform after the initial mesh was created.

  • node_vertex_consistency_tolerance – passed to the Mesh constructor. If False, no checks are performed.

  • unit_nodes – if given, the unit nodes to use. Must have shape (3, nnodes).

meshmode.mesh.generation.generate_torus(r_major: float, r_minor: float, n_major: int = 20, n_minor: int = 10, order: int = 1, node_vertex_consistency_tolerance: Optional[Union[float, bool]] = None, unit_nodes: Optional[numpy.ndarray] = None)

Generate a torus.

_images/torus.png

Shown: A torus with major circle (magenta) and minor circle (red). Source: https://commons.wikimedia.org/wiki/File:Torus_cycles.svg (public domain image by Krishnavedala).

The torus is obtained as the image of the parameter domain \((u, v) \in [0, 2\pi) \times [0, 2 \pi)\) under the map

\[\begin{split}\begin{align} x &= \cos(u) (r_\text{major} + r_\text{minor} \cos(v)) \\ y &= \sin(u) (r_\text{major} + r_\text{minor} \sin(v)) \\ z &= r_\text{minor} \sin(v) \end{align}\end{split}\]

where \(r_\text{major}\) and \(r_\text{minor}\) are the radii of the major and minor circles, respectively. The parameter domain is tiled with \(n_\text{major} \times n_\text{minor}\) contiguous rectangles, and then each rectangle is subdivided into two triangles.

Parameters
  • r_major – radius of the major circle.

  • r_minor – radius of the minor circle.

  • n_major – number of rectangles along major circle.

  • n_minor – number of rectangles along minor circle.

  • order – order of the (simplex) elements. If unit_nodes is also provided, the orders should match.

  • node_vertex_consistency_tolerance – passed to the Mesh constructor. If False, no checks are performed.

  • unit_nodes – if given, the unit nodes to use. Must have shape (3, nnodes).

Returns

a Mesh of a torus.

meshmode.mesh.generation.refine_mesh_and_get_urchin_warper(order: int, m: int, n: int, est_rel_interp_tolerance: float, min_rad: float = 0.2, uniform_refinement_rounds: int = 0)
Parameters
  • order – order of the (simplex) elements.

  • m – order of the spherical harmonic \(Y^m_n\).

  • n – order of the spherical harmonic \(Y^m_n\).

  • est_rel_interp_tolerance – a tolerance for the relative interpolation error estimates on the warped version of the mesh.

Returns

a tuple (refiner, warp_mesh), where refiner is a Refiner (from which the unwarped mesh may be obtained), and whose get_current_mesh() returns a locally-refined Mesh of a sphere and warp_mesh is a callable taking and returning a mesh that warps the unwarped mesh into a smooth shape covered by a spherical harmonic of order \((m, n)\).

meshmode.mesh.generation.generate_urchin(order: int, m: int, n: int, est_rel_interp_tolerance: float, min_rad: float = 0.2)
Parameters
  • order – order of the (simplex) elements. If unit_nodes is also provided, the orders should match.

  • m – order of the spherical harmonic \(Y^m_n\).

  • n – order of the spherical harmonic \(Y^m_n\).

  • est_rel_interp_tolerance – a tolerance for the relative interpolation error estimates on the warped version of the mesh.

Returns

a refined Mesh of a smooth shape covered by a spherical harmonic of order \((m, n)\).

meshmode.mesh.generation.generate_surface_of_revolution(get_radius: Callable[[numpy.ndarray, numpy.ndarray], numpy.ndarray], height_discr: numpy.ndarray, angle_discr: numpy.ndarray, order: int, *, node_vertex_consistency_tolerance: Optional[Union[float, bool]] = None, unit_nodes: Optional[numpy.ndarray] = None)

Return a cylinder aligned with the “height” axis aligned with the Z axis.

Parameters
  • get_radius – A callable function that takes in a 1D array of heights and a 1D array of angles and returns a 1D array of radii.

  • height_discr – A discretization of [0, 2*pi).

  • angle_discr – A discretization of [0, 2*pi).

  • order – order of the (simplex) elements. If unit_nodes is also provided, the orders should match.

  • node_vertex_consistency_tolerance – passed to the Mesh constructor. If False, no checks are performed.

  • unit_nodes – if given, the unit nodes to use. Must have shape (3, nnodes).

Volumes

meshmode.mesh.generation.generate_box_mesh(axis_coords, order=1, coord_dtype=<class 'numpy.float64'>, group_cls=None, boundary_tag_to_face=None, mesh_type=None)

Create a semi-structured mesh.

Parameters
  • axis_coords – a tuple with a number of entries corresponding to the number of dimensions, with each entry a numpy array specifying the coordinates to be used along that axis.

  • group_cls – One of meshmode.mesh.SimplexElementGroup or meshmode.mesh.TensorProductElementGroup.

  • boundary_tag_to_face

    an optional dictionary for tagging boundaries. The keys correspond to custom boundary tags, with the values giving a list of the faces on which they should be applied in terms of coordinate directions (+x, -x, +y, -y, +z, -z, +w, -w).

    For example:

    boundary_tag_to_face={"bdry_1": ["+x", "+y"], "bdry_2": ["-x"]}
    

  • mesh_type

    In two dimensions with non-tensor-product elements, mesh_type may be set to "X" to generate this type of mesh:

    _______
    |\   /|
    | \ / |
    |  X  |
    | / \ |
    |/   \|
    ^^^^^^^
    

    instead of the default:

    _______
    |\    |
    | \   |
    |  \  |
    |   \ |
    |    \|
    ^^^^^^^
    

    Specifying a value other than None for all other mesh dimensionalities and element types is an error.

Changed in version 2017.1: group_factory parameter added.

Changed in version 2020.1: boundary_tag_to_face parameter added.

Changed in version 2020.3: group_factory deprecated and renamed to group_cls.

meshmode.mesh.generation.generate_regular_rect_mesh(a=(0, 0), b=(1, 1), *, nelements_per_axis=None, npoints_per_axis=None, order=1, boundary_tag_to_face=None, group_cls=None, mesh_type=None, n=None)

Create a semi-structured rectangular mesh with equispaced elements.

Parameters
  • a – the lower left hand point of the rectangle.

  • b – the upper right hand point of the rectangle.

  • nelements_per_axis – an optional tuple of integers indicating the number of elements along each axis.

  • npoints_per_axis – an optional tuple of integers indicating the number of points along each axis.

  • order – the mesh element order.

  • boundary_tag_to_face – an optional dictionary for tagging boundaries. See generate_box_mesh().

  • group_cls – see generate_box_mesh().

  • mesh_type – see generate_box_mesh().

Note

Specify only one of nelements_per_axis and npoints_per_axis.

meshmode.mesh.generation.generate_warped_rect_mesh(dim, order, *, nelements_side=None, npoints_side=None, group_cls=None, n=None)

Generate a mesh of a warped square/cube. Mainly useful for testing functionality with curvilinear meshes.

Tools for Iterative Refinement

meshmode.mesh.generation.warp_and_refine_until_resolved(unwarped_mesh_or_refiner, warp_callable, est_rel_interp_tolerance)

Given an original (“unwarped”) meshmode.mesh.Mesh and a warping function warp_callable that takes and returns a mesh and a tolerance to which the mesh should be resolved by the mapping polynomials, this function will iteratively refine the unwarped_mesh until relative interpolation error estimates on the warped version are smaller than est_rel_interp_tolerance on each element.

Returns

The refined, unwarped mesh.

New in version 2018.1.

Mesh input/output

class meshmode.mesh.io.ScriptSource(source, extension)

New in version 2016.1.

class meshmode.mesh.io.FileSource(filename)

New in version 2014.1.

class meshmode.mesh.io.ScriptWithFilesSource(source, filenames, source_name='temp.geo')

New in version 2016.1.

source

The script code to be fed to gmsh.

filenames

The names of files to be copied to the temporary directory where gmsh is run.

meshmode.mesh.io.read_gmsh(filename, force_ambient_dim=None, mesh_construction_kwargs=None)

Read a gmsh mesh file from filename and return a meshmode.mesh.Mesh.

Parameters
  • force_ambient_dim – if not None, truncate point coordinates to this many dimensions.

  • mesh_construction_kwargsNone or a dictionary of keyword arguments passed to the meshmode.mesh.Mesh constructor.

meshmode.mesh.io.generate_gmsh(source, dimensions=None, order=None, other_options=None, extension='geo', gmsh_executable='gmsh', force_ambient_dim=None, output_file_name='output.msh', mesh_construction_kwargs=None, target_unit=None)

Run gmsh on the input given by source, and return a meshmode.mesh.Mesh based on the result.

Parameters
  • source – an instance of either gmsh_interop.reader.FileSource or gmsh_interop.reader.ScriptSource

  • force_ambient_dim – if not None, truncate point coordinates to this many dimensions.

  • mesh_construction_kwargsNone or a dictionary of keyword arguments passed to the meshmode.mesh.Mesh constructor.

  • target_unit – Value of the option Geometry.OCCTargetUnit. Supported values are the strings ‘M’ or ‘MM’.

meshmode.mesh.io.from_meshpy(mesh_info, order=1)

Imports a mesh from a meshpy mesh_info data structure, which may be generated by either meshpy.triangle or meshpy.tet.

meshmode.mesh.io.from_vertices_and_simplices(vertices, simplices, order=1, fix_orientation=False)

Imports a mesh from a numpy array of vertices and an array of simplices.

Parameters
  • vertices – An array of vertex coordinates with shape (ambient_dim, nvertices)

  • simplices – An array (nelements, nvertices) of (mesh-wide) vertex indices.

meshmode.mesh.io.to_json(mesh)

Return a JSON-able Python data structure for mesh. The structure directly reflects the meshmode.mesh.Mesh data structure.

Mesh processing

meshmode.mesh.processing.find_group_indices(groups, meshwide_elems)
Parameters
  • groups – A list of MeshElementGroup instances that contain meshwide_elems.

  • meshwide_elems – A numpy.ndarray of mesh-wide element numbers. Usually computed by elem + element_nr_base.

Returns

A numpy.ndarray of group numbers that meshwide_elem belongs to.

meshmode.mesh.processing.partition_mesh(mesh, part_per_element, part_num)
Parameters
  • mesh – A Mesh to be partitioned.

  • part_per_element – A numpy.ndarray containing one integer per element of mesh indicating which part of the partitioned mesh the element is to become a part of.

  • part_num – The part number of the mesh to return.

Returns

A tuple (part_mesh, part_to_global), where part_mesh is a Mesh that is a partition of mesh, and part_to_global is a numpy.ndarray mapping element numbers on part_mesh to ones in mesh.

New in version 2017.1.

meshmode.mesh.processing.find_volume_mesh_element_orientations(mesh, tolerate_unimplemented_checks=False)

Return a positive floating point number for each positively oriented element, and a negative floating point number for each negatively oriented element.

Parameters

tolerate_unimplemented_checks – If True, elements for which no check is available will return NaN.

meshmode.mesh.processing.flip_simplex_element_group(vertices, grp, grp_flip_flags)
meshmode.mesh.processing.perform_flips(mesh, flip_flags, skip_tests=False)
Parameters

flip_flags – A numpy.ndarray with meshmode.mesh.Mesh.nelements entries indicating by their Boolean value whether the element is to be flipped.

meshmode.mesh.processing.find_bounding_box(mesh)
Returns

a tuple (min, max), each consisting of a numpy.ndarray indicating the minimal and maximal extent of the geometry along each axis.

meshmode.mesh.processing.merge_disjoint_meshes(meshes, skip_tests=False, single_group=False)
meshmode.mesh.processing.split_mesh_groups(mesh, element_flags, return_subgroup_mapping=False)

Split all the groups in mesh according to the values of element_flags. The element flags are expected to be integers defining, for each group, how the elements are to be split into subgroups. For example, a single-group mesh with flags:

element_flags = [0, 0, 0, 42, 42, 42, 0, 0, 0, 41, 41, 41]

will create three subgroups. The integer flags need not be increasing or contiguous and can repeat across different groups (i.e. they are group-local).

Parameters

element_flags – a numpy.ndarray with nelements entries indicating how the elements in a group are to be split.

Returns

a Mesh where each group has been split according to flags in element_flags. If return_subgroup_mapping is True, it also returns a mapping of (group_index, subgroup) -> new_group_index.

meshmode.mesh.processing.map_mesh(mesh, f)

Apply the map f to the mesh. f needs to accept and return arrays of shape (ambient_dim, npoints).

meshmode.mesh.processing.affine_map(mesh, A: Optional[Union[numbers.Real, numpy.ndarray]] = None, b: Optional[Union[numbers.Real, numpy.ndarray]] = None)

Apply the affine map \(f(x) = A x + b\) to the geometry of mesh.

meshmode.mesh.processing.rotate_mesh_around_axis(mesh, *, theta: numbers.Real, axis: Optional[numpy.ndarray] = None)

Rotate the mesh by theta radians around the axis axis.

Parameters

axis – a (not necessarily unit) vector. By default, the rotation is performed around the \(z\) axis.

Mesh refinement

class meshmode.mesh.refinement.Refiner(mesh)

An older that mostly succeeds at preserving adjacency across non-conformal refinement.

Note

This refiner is currently kind of slow, and not always correct. See RefinerWithoutAdjacency for a less capable but much faster refiner.

__init__(mesh)

Initialize self. See help(type(self)) for accurate signature.

refine(refine_flags)
Parameters

refine_flags – an ndarray of dtype bool and length meshmode.mesh.Mesh.nelements, indicating which elements should be split.

refine_uniformly()
get_current_mesh()
get_previous_mesh()
class meshmode.mesh.refinement.RefinerWithoutAdjacency(mesh)

A refiner that may be applied to non-conforming meshmode.mesh.Mesh instances. It does not generate adjacency information, and it is typically faster than meshmode.mesh.refinement.Refiner.

Note

If the input meshes to this refiner are not conforming, then the resulting meshes may contain duplicated vertices. (I.e. two different numbers referring to the same geometric vertex.)

__init__(mesh)

Initialize self. See help(type(self)) for accurate signature.

refine(refine_flags)
Parameters

refine_flags – an ndarray of dtype bool and length meshmode.mesh.Mesh.nelements indicating which elements should be split.

refine_uniformly()
get_current_mesh()
get_previous_mesh()
meshmode.mesh.refinement.refine_uniformly(mesh, iterations, with_adjacency=False)

Mesh visualization

meshmode.mesh.visualization.draw_2d_mesh(mesh, draw_vertex_numbers=True, draw_element_numbers=True, draw_nodal_adjacency=False, draw_face_numbers=False, set_bounding_box=False, **kwargs)
meshmode.mesh.visualization.draw_curve(mesh, el_bdry_style='o', el_bdry_kwargs=None, node_style='x-', node_kwargs=None)
meshmode.mesh.visualization.write_vertex_vtk_file(mesh, file_name, compressor=None, overwrite=False)
meshmode.mesh.visualization.mesh_to_tikz(mesh)
meshmode.mesh.visualization.vtk_visualize_mesh(actx, mesh, filename, vtk_high_order=True)