Skip to content

Graph & Connectivity#

Derived connectivity tables, graph traversal, and binary morphology on a mesh.

Everything here is pure: each function reads a Mesh2d and returns new arrays or a new Connectivity, never writing a derived table back onto the mesh it was handed. The edge tables are derived from the faces; when a file declares its own edge_node_connectivity those edges keep their positions (edge-located variables index into them) and any face-implied edge the file omits is appended after them, so a partial edge list still yields a complete face topology.

sparse_adjacency is the shared primitive: connected_components, reverse_cuthill_mckee and the binary morphology all read it.

One morphology pass is simultaneous — the neighbour test is taken from the mask as it stands at the start of the pass, and the border_value imputation joins that same update — so one iteration moves exactly one ring, matching scipy.ndimage.

pyramids.netcdf.ugrid.graph #

Mesh graph and connectivity operations for UGRID 2D meshes.

A collection of pure functions that each take a :class:~pyramids.netcdf.ugrid.mesh.Mesh2d (and/or plain numpy arrays) and return plain results — numpy arrays, :class:Connectivity wrappers, scipy sparse matrices, shapely geometries, or a rebuilt Mesh2d for the one topology-transforming helper. None of them mutate their input.

The module groups the graph-flavoured mesh operations behind a single import:

  • Connectivity builders that pyramids does not derive today: node_node_connectivity, node_face_connectivity, node_edge_connectivity, edge_node_connectivity, edge_face_connectivity, face_edge_connectivity and edge_edge_connectivity.
  • sparse_adjacency — a scipy sparse adjacency matrix over faces or nodes, the shared primitive the graph algorithms below build on.
  • Facet remapping (to_node / to_edge / to_face) that averages a per-element field onto another mesh facet via connectivity.
  • Graph algorithms: connected_components, reverse_cuthill_mckee (+ permute_faces), binary_dilation / binary_erosion, voronoi / circumcenters and the boundary helpers exterior_edges / exterior_faces / bounding_polygon.

scipy (scipy.sparse, scipy.sparse.csgraph, scipy.spatial) and shapely are hard dependencies of the UGRID subsystem, so they are imported at module top.

node_face_connectivity(mesh) #

Faces incident to each node (inverse of face_node_connectivity).

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required

Returns:

Name Type Description
Connectivity Connectivity

Ragged node -> faces mapping, -1 padded.

Source code in src/pyramids/netcdf/ugrid/graph.py
def node_face_connectivity(mesh: Mesh2d) -> Connectivity:
    """Faces incident to each node (inverse of ``face_node_connectivity``).

    Args:
        mesh: The mesh to read.

    Returns:
        Connectivity: Ragged ``node -> faces`` mapping, ``-1`` padded.
    """
    return _invert_connectivity(
        mesh.face_node_connectivity, mesh.n_node, "node_face_connectivity"
    )

edge_node_connectivity(mesh) #

The two nodes of each edge, derived from the faces when absent.

A row the file declared but never defined -- one padded with the fill value, or a repeat of an edge already listed -- comes back as a (-1, -1) placeholder. It keeps its position, so an edge-located variable stays aligned with the table, but it describes no edge: :func:_real_edges is the mask of the rows that do, and :func:edge_face_connectivity / :func:edge_edge_connectivity answer an all-fill row at the same positions.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required

Returns:

Name Type Description
Connectivity Connectivity

edge -> (node, node) mapping, with a (-1, -1) row wherever the

Connectivity

file's edge list declared an edge it did not define.

Source code in src/pyramids/netcdf/ugrid/graph.py
def edge_node_connectivity(mesh: Mesh2d) -> Connectivity:
    """The two nodes of each edge, derived from the faces when absent.

    A row the file declared but never defined -- one padded with the fill value, or a repeat of
    an edge already listed -- comes back as a ``(-1, -1)`` **placeholder**. It keeps its
    position, so an edge-located variable stays aligned with the table, but it describes no
    edge: :func:`_real_edges` is the mask of the rows that do, and
    :func:`edge_face_connectivity` / :func:`edge_edge_connectivity` answer an all-fill row at
    the same positions.

    Args:
        mesh: The mesh to read.

    Returns:
        Connectivity: ``edge -> (node, node)`` mapping, with a ``(-1, -1)`` row wherever the
        file's edge list declared an edge it did not define.
    """
    edge_node, _, _ = _edge_tables(mesh)
    return Connectivity(
        data=edge_node,
        fill_value=-1,
        cf_role="edge_node_connectivity",
        original_start_index=0,
    )

node_edge_connectivity(mesh) #

Edges incident to each node (inverse of edge_node_connectivity).

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required

Returns:

Name Type Description
Connectivity Connectivity

Ragged node -> edges mapping, -1 padded.

Source code in src/pyramids/netcdf/ugrid/graph.py
def node_edge_connectivity(mesh: Mesh2d) -> Connectivity:
    """Edges incident to each node (inverse of ``edge_node_connectivity``).

    Args:
        mesh: The mesh to read.

    Returns:
        Connectivity: Ragged ``node -> edges`` mapping, ``-1`` padded.
    """
    edge_node, _, _ = _edge_tables(mesh)
    # Keep the edge ids of the real rows only: a placeholder row is not an edge of any node, and
    # its -1 node index would reach np.bincount, which refuses negatives.
    real = _real_edges(edge_node)
    e_ids = np.arange(edge_node.shape[0], dtype=np.intp)[real]
    pairs = edge_node[real]
    keys = np.concatenate([pairs[:, 0], pairs[:, 1]])
    values = np.concatenate([e_ids, e_ids])
    return _group_to_connectivity(keys, values, mesh.n_node, "node_edge_connectivity")

node_node_connectivity(mesh) #

Neighbouring nodes of each node, derived from the mesh edges.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required

Returns:

Name Type Description
Connectivity Connectivity

Ragged node -> neighbour nodes mapping, -1 padded.

Source code in src/pyramids/netcdf/ugrid/graph.py
def node_node_connectivity(mesh: Mesh2d) -> Connectivity:
    """Neighbouring nodes of each node, derived from the mesh edges.

    Args:
        mesh: The mesh to read.

    Returns:
        Connectivity: Ragged ``node -> neighbour nodes`` mapping, ``-1`` padded.
    """
    edge_node, _, _ = _edge_tables(mesh)
    pairs = edge_node[_real_edges(edge_node)]
    keys = np.concatenate([pairs[:, 0], pairs[:, 1]])
    values = np.concatenate([pairs[:, 1], pairs[:, 0]])
    return _group_to_connectivity(keys, values, mesh.n_node, "node_node_connectivity")

edge_face_connectivity(mesh) #

The (up to two) faces on either side of each edge.

Edges on the mesh boundary carry a -1 fill on their missing side.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required

Returns:

Name Type Description
Connectivity Connectivity

edge -> faces mapping, -1 for an exterior side.

Source code in src/pyramids/netcdf/ugrid/graph.py
def edge_face_connectivity(mesh: Mesh2d) -> Connectivity:
    """The (up to two) faces on either side of each edge.

    Edges on the mesh boundary carry a ``-1`` fill on their missing side.

    Args:
        mesh: The mesh to read.

    Returns:
        Connectivity: ``edge -> faces`` mapping, ``-1`` for an exterior side.
    """
    _, edge_face, _ = _edge_tables(mesh)
    return Connectivity(
        data=edge_face,
        fill_value=-1,
        cf_role="edge_face_connectivity",
        original_start_index=0,
    )

face_edge_connectivity(mesh) #

The boundary edges of each face.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required

Returns:

Name Type Description
Connectivity Connectivity

Ragged face -> edges mapping, -1 padded.

Source code in src/pyramids/netcdf/ugrid/graph.py
def face_edge_connectivity(mesh: Mesh2d) -> Connectivity:
    """The boundary edges of each face.

    Args:
        mesh: The mesh to read.

    Returns:
        Connectivity: Ragged ``face -> edges`` mapping, ``-1`` padded.
    """
    _, _, face_edge = _edge_tables(mesh)
    return face_edge

edge_edge_connectivity(mesh) #

Edges that share a node with each edge.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required

Returns:

Name Type Description
Connectivity Connectivity

Ragged edge -> neighbouring edges mapping, -1 padded.

Source code in src/pyramids/netcdf/ugrid/graph.py
def edge_edge_connectivity(mesh: Mesh2d) -> Connectivity:
    """Edges that share a node with each edge.

    Args:
        mesh: The mesh to read.

    Returns:
        Connectivity: Ragged ``edge -> neighbouring edges`` mapping, ``-1`` padded.
    """
    node_edge = node_edge_connectivity(mesh)
    edge_node, _, _ = _edge_tables(mesh)
    n_edge = int(edge_node.shape[0])
    keys_list: list[np.typing.NDArray] = []
    values_list: list[np.typing.NDArray] = []
    for r in range(node_edge.n_elements):
        edges = node_edge.get_element(r)
        if edges.size < 2:
            continue
        left, right = np.meshgrid(edges, edges)
        off_diag = left != right
        keys_list.append(left[off_diag])
        values_list.append(right[off_diag])
    if keys_list:
        keys = np.concatenate(keys_list)
        values = np.concatenate(values_list)
    else:
        keys = np.empty(0, dtype=np.intp)
        values = np.empty(0, dtype=np.intp)
    return _group_to_connectivity(keys, values, n_edge, "edge_edge_connectivity")

sparse_adjacency(mesh, kind='face') #

Build a symmetric scipy sparse adjacency matrix for the mesh.

For kind="face" two faces are adjacent when they share an edge (an edge whose edge_face_connectivity has two valid faces). For kind="node" two nodes are adjacent when an edge joins them. The result is the shared primitive consumed by :func:connected_components, :func:reverse_cuthill_mckee and the binary morphology.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required
kind str

"face" (default) or "node".

'face'

Returns:

Type Description
Any

scipy.sparse.csr_matrix: The symmetric (n, n) adjacency matrix.

Raises:

Type Description
ValueError

If kind is neither "face" nor "node".

Source code in src/pyramids/netcdf/ugrid/graph.py
def sparse_adjacency(mesh: Mesh2d, kind: str = "face") -> Any:
    """Build a symmetric scipy sparse adjacency matrix for the mesh.

    For ``kind="face"`` two faces are adjacent when they share an edge (an edge whose
    ``edge_face_connectivity`` has two valid faces). For ``kind="node"`` two nodes are
    adjacent when an edge joins them. The result is the shared primitive consumed by
    :func:`connected_components`, :func:`reverse_cuthill_mckee` and the binary morphology.

    Args:
        mesh: The mesh to read.
        kind: ``"face"`` (default) or ``"node"``.

    Returns:
        scipy.sparse.csr_matrix: The symmetric ``(n, n)`` adjacency matrix.

    Raises:
        ValueError: If ``kind`` is neither ``"face"`` nor ``"node"``.
    """
    if kind == "face":
        _, edge_face, _ = _edge_tables(mesh)
        both_valid = (edge_face[:, 0] != -1) & (edge_face[:, 1] != -1)
        rows = edge_face[both_valid, 0]
        cols = edge_face[both_valid, 1]
        n = mesh.n_face
    elif kind == "node":
        edge_node, _, _ = _edge_tables(mesh)
        pairs = edge_node[_real_edges(edge_node)]
        rows = pairs[:, 0]
        cols = pairs[:, 1]
        n = mesh.n_node
    else:
        raise ValueError(f"kind must be 'face' or 'node', got {kind!r}.")

    sym_rows = np.concatenate([rows, cols])
    sym_cols = np.concatenate([cols, rows])
    data = np.ones(sym_rows.size, dtype=np.int8)
    matrix = sparse.coo_matrix((data, (sym_rows, sym_cols)), shape=(n, n)).tocsr()
    matrix.data[:] = 1
    return matrix

to_node(mesh, values, source='face') #

Remap a per-face or per-edge field onto the nodes by averaging.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required
values NDArray

Source field, element axis last.

required
source str

"face" (default) or "edge".

'face'

Returns:

Type Description
NDArray

np.ndarray: The node-located field.

Raises:

Type Description
ValueError

If source is not "face" or "edge".

Source code in src/pyramids/netcdf/ugrid/graph.py
def to_node(
    mesh: Mesh2d, values: np.typing.NDArray, source: str = "face"
) -> np.typing.NDArray:
    """Remap a per-face or per-edge field onto the nodes by averaging.

    Args:
        mesh: The mesh to read.
        values: Source field, element axis last.
        source: ``"face"`` (default) or ``"edge"``.

    Returns:
        np.ndarray: The node-located field.

    Raises:
        ValueError: If ``source`` is not ``"face"`` or ``"edge"``.
    """
    if source == "face":
        conn = node_face_connectivity(mesh)
    elif source == "edge":
        conn = node_edge_connectivity(mesh)
    else:
        raise ValueError(f"source must be 'face' or 'edge', got {source!r}.")
    return _average_over(conn, values)

to_face(mesh, values, source='node') #

Remap a per-node or per-edge field onto the faces by averaging.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required
values NDArray

Source field, element axis last.

required
source str

"node" (default) or "edge".

'node'

Returns:

Type Description
NDArray

np.ndarray: The face-located field.

Raises:

Type Description
ValueError

If source is not "node" or "edge".

Source code in src/pyramids/netcdf/ugrid/graph.py
def to_face(
    mesh: Mesh2d, values: np.typing.NDArray, source: str = "node"
) -> np.typing.NDArray:
    """Remap a per-node or per-edge field onto the faces by averaging.

    Args:
        mesh: The mesh to read.
        values: Source field, element axis last.
        source: ``"node"`` (default) or ``"edge"``.

    Returns:
        np.ndarray: The face-located field.

    Raises:
        ValueError: If ``source`` is not ``"node"`` or ``"edge"``.
    """
    if source == "node":
        conn = mesh.face_node_connectivity
    elif source == "edge":
        conn = face_edge_connectivity(mesh)
    else:
        raise ValueError(f"source must be 'node' or 'edge', got {source!r}.")
    return _average_over(conn, values)

to_edge(mesh, values, source='node') #

Remap a per-node or per-face field onto the edges by averaging.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required
values NDArray

Source field, element axis last.

required
source str

"node" (default) or "face".

'node'

Returns:

Type Description
NDArray

np.ndarray: The edge-located field.

Raises:

Type Description
ValueError

If source is not "node" or "face".

Source code in src/pyramids/netcdf/ugrid/graph.py
def to_edge(
    mesh: Mesh2d, values: np.typing.NDArray, source: str = "node"
) -> np.typing.NDArray:
    """Remap a per-node or per-face field onto the edges by averaging.

    Args:
        mesh: The mesh to read.
        values: Source field, element axis last.
        source: ``"node"`` (default) or ``"face"``.

    Returns:
        np.ndarray: The edge-located field.

    Raises:
        ValueError: If ``source`` is not ``"node"`` or ``"face"``.
    """
    if source == "node":
        conn = edge_node_connectivity(mesh)
    elif source == "face":
        conn = edge_face_connectivity(mesh)
    else:
        raise ValueError(f"source must be 'node' or 'face', got {source!r}.")
    return _average_over(conn, values)

connected_components(mesh) #

Label each face by the connected component it belongs to.

Two faces are in the same component when a path of shared edges joins them — the mesh analogue of a raster cluster labelling.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required

Returns:

Type Description
NDArray

np.ndarray: Per-face component label of length n_face.

Source code in src/pyramids/netcdf/ugrid/graph.py
def connected_components(mesh: Mesh2d) -> np.typing.NDArray:
    """Label each face by the connected component it belongs to.

    Two faces are in the same component when a path of shared edges joins them — the mesh
    analogue of a raster cluster labelling.

    Args:
        mesh: The mesh to read.

    Returns:
        np.ndarray: Per-face component label of length ``n_face``.
    """
    adjacency = sparse_adjacency(mesh, kind="face")
    _, labels = csgraph.connected_components(adjacency, directed=False)
    return np.asarray(labels, dtype=np.intp)

reverse_cuthill_mckee(mesh) #

Return a face permutation that reduces the face-adjacency bandwidth.

Computes the reverse Cuthill-McKee ordering of the face adjacency. Apply it to the mesh with :func:permute_faces and to every face-located variable with arr[..., perm].

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required

Returns:

Type Description
NDArray

np.ndarray: A permutation of range(n_face).

Source code in src/pyramids/netcdf/ugrid/graph.py
def reverse_cuthill_mckee(mesh: Mesh2d) -> np.typing.NDArray:
    """Return a face permutation that reduces the face-adjacency bandwidth.

    Computes the reverse Cuthill-McKee ordering of the face adjacency. Apply it to the mesh
    with :func:`permute_faces` and to every face-located variable with ``arr[..., perm]``.

    Args:
        mesh: The mesh to read.

    Returns:
        np.ndarray: A permutation of ``range(n_face)``.
    """
    if mesh.n_face == 0:
        # scipy refuses a 0x0 matrix. An empty mesh has an empty ordering, which is the
        # answer rather than an error.
        return np.empty(0, dtype=np.intp)
    adjacency = sparse.csr_matrix(sparse_adjacency(mesh, kind="face"))
    perm = csgraph.reverse_cuthill_mckee(adjacency, symmetric_mode=True)
    return np.asarray(perm, dtype=np.intp)

permute_faces(mesh, perm) #

Rebuild a mesh with its faces reordered by perm.

Node coordinates and numbering are unchanged — only the face rows (and any face centre coordinates) are reordered — so the result is congruent to the input. New face k is old face perm[k]; the input mesh is not mutated, and no array is shared with it.

The three face-indexed tables -- face_edge_connectivity, face_face_connectivity and edge_face_connectivity -- are dropped rather than permuted, since their contents index the old face numbering. Derive them again from the result if you need them.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to reorder.

required
perm NDArray

A permutation of range(n_face).

required

Returns:

Name Type Description
Mesh2d Mesh2d

The reordered mesh.

Source code in src/pyramids/netcdf/ugrid/graph.py
def permute_faces(mesh: Mesh2d, perm: np.typing.NDArray) -> Mesh2d:
    """Rebuild a mesh with its faces reordered by ``perm``.

    Node coordinates and numbering are unchanged — only the face rows (and any face centre
    coordinates) are reordered — so the result is congruent to the input. New face ``k`` is
    old face ``perm[k]``; the input mesh is not mutated, and no array is shared with it.

    The three face-indexed tables -- ``face_edge_connectivity``, ``face_face_connectivity``
    and ``edge_face_connectivity`` -- are **dropped** rather than permuted, since their
    contents index the old face numbering. Derive them again from the result if you need
    them.

    Args:
        mesh: The mesh to reorder.
        perm: A permutation of ``range(n_face)``.

    Returns:
        Mesh2d: The reordered mesh.
    """
    perm = np.asarray(perm, dtype=np.intp)
    old_fnc = mesh.face_node_connectivity
    new_fnc = Connectivity(
        data=old_fnc.data[perm].copy(),
        fill_value=old_fnc.fill_value,
        cf_role=old_fnc.cf_role,
        original_start_index=old_fnc.original_start_index,
    )
    face_x = mesh.face_x[perm] if mesh.face_x is not None else None
    face_y = mesh.face_y[perm] if mesh.face_y is not None else None
    return Mesh2d(
        node_x=mesh.node_x.copy(),
        node_y=mesh.node_y.copy(),
        face_node_connectivity=new_fnc,
        # Copy what is carried over, as the permuted arrays above are. Aliasing the
        # untouched tables meant a write through the result reached the input mesh.
        edge_node_connectivity=_copied_connectivity(mesh.edge_node_connectivity),
        face_x=face_x,
        face_y=face_y,
        edge_x=None if mesh.edge_x is None else np.array(mesh.edge_x, copy=True),
        edge_y=None if mesh.edge_y is None else np.array(mesh.edge_y, copy=True),
        crs=mesh.crs,
    )

binary_dilation(mesh, mask, iterations=1, border_value=False, fixed=None) #

Grow a True region of a boolean face mask across shared edges.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required
mask NDArray

Input boolean face mask (not mutated).

required
iterations int

Number of rings to grow.

1
border_value bool

Value imputed to the mesh exterior (grows inward from the boundary when True).

False
fixed NDArray | None

Optional boolean faces held at their input value.

None

Returns:

Type Description
NDArray

np.ndarray: The dilated boolean face mask.

Source code in src/pyramids/netcdf/ugrid/graph.py
def binary_dilation(
    mesh: Mesh2d,
    mask: np.typing.NDArray,
    iterations: int = 1,
    border_value: bool = False,
    fixed: np.typing.NDArray | None = None,
) -> np.typing.NDArray:
    """Grow a ``True`` region of a boolean face mask across shared edges.

    Args:
        mesh: The mesh to read.
        mask: Input boolean face mask (not mutated).
        iterations: Number of rings to grow.
        border_value: Value imputed to the mesh exterior (grows inward from the boundary
            when ``True``).
        fixed: Optional boolean faces held at their input value.

    Returns:
        np.ndarray: The dilated boolean face mask.
    """
    return _binary_iterate(mesh, mask, True, iterations, border_value, fixed)

binary_erosion(mesh, mask, iterations=1, border_value=False, fixed=None) #

Shrink a True region of a boolean face mask across shared edges.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required
mask NDArray

Input boolean face mask (not mutated).

required
iterations int

Number of rings to shrink.

1
border_value bool

Value imputed to the mesh exterior (erodes from the boundary when False).

False
fixed NDArray | None

Optional boolean faces held at their input value.

None

Returns:

Type Description
NDArray

np.ndarray: The eroded boolean face mask.

Source code in src/pyramids/netcdf/ugrid/graph.py
def binary_erosion(
    mesh: Mesh2d,
    mask: np.typing.NDArray,
    iterations: int = 1,
    border_value: bool = False,
    fixed: np.typing.NDArray | None = None,
) -> np.typing.NDArray:
    """Shrink a ``True`` region of a boolean face mask across shared edges.

    Args:
        mesh: The mesh to read.
        mask: Input boolean face mask (not mutated).
        iterations: Number of rings to shrink.
        border_value: Value imputed to the mesh exterior (erodes from the boundary when
            ``False``).
        fixed: Optional boolean faces held at their input value.

    Returns:
        np.ndarray: The eroded boolean face mask.
    """
    return _binary_iterate(mesh, mask, False, iterations, border_value, fixed)

circumcenters(mesh) #

Per-face circumcenter coordinates for a triangular mesh.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read (every face must be a triangle).

required

Returns:

Type Description
NDArray

np.ndarray: (n_face, 2) array of circumcenter (x, y) coordinates.

Raises:

Type Description
ValueError

If the mesh is not purely triangular.

Source code in src/pyramids/netcdf/ugrid/graph.py
def circumcenters(mesh: Mesh2d) -> np.typing.NDArray:
    """Per-face circumcenter coordinates for a triangular mesh.

    Args:
        mesh: The mesh to read (every face must be a triangle).

    Returns:
        np.ndarray: ``(n_face, 2)`` array of circumcenter ``(x, y)`` coordinates.

    Raises:
        ValueError: If the mesh is not purely triangular.
    """
    fnc = mesh.face_node_connectivity
    if not fnc.is_triangular():
        raise ValueError("circumcenters is only defined for triangular meshes.")
    tri = fnc.data[:, :3]
    ax, ay = mesh.node_x[tri[:, 0]], mesh.node_y[tri[:, 0]]
    bx, by = mesh.node_x[tri[:, 1]], mesh.node_y[tri[:, 1]]
    cx, cy = mesh.node_x[tri[:, 2]], mesh.node_y[tri[:, 2]]
    d = 2.0 * (ax * (by - cy) + bx * (cy - ay) + cx * (ay - by))
    a_sq = ax * ax + ay * ay
    b_sq = bx * bx + by * by
    c_sq = cx * cx + cy * cy
    # A degenerate triangle has d == 0 and therefore no circumcentre. Dividing under
    # errstate keeps the inf/nan answer the caller can test for without numpy printing a
    # divide-by-zero or invalid-value warning from inside a library call.
    with np.errstate(divide="ignore", invalid="ignore"):
        ux = (a_sq * (by - cy) + b_sq * (cy - ay) + c_sq * (ay - by)) / d
        uy = (a_sq * (cx - bx) + b_sq * (ax - cx) + c_sq * (bx - ax)) / d
    return np.column_stack([ux, uy])

voronoi(mesh) #

Voronoi tessellation of the mesh nodes (the mesh dual).

Uses :class:scipy.spatial.Voronoi over the node coordinates. For a triangular mesh that is Delaunay the Voronoi vertices coincide with the face :func:circumcenters, and interior nodes receive bounded dual cells while boundary nodes receive unbounded (-1) regions.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required

Returns:

Type Description
Any

scipy.spatial.Voronoi: The Voronoi object; .points are the node coordinates.

Source code in src/pyramids/netcdf/ugrid/graph.py
def voronoi(mesh: Mesh2d) -> Any:
    """Voronoi tessellation of the mesh nodes (the mesh dual).

    Uses :class:`scipy.spatial.Voronoi` over the node coordinates. For a triangular mesh that
    is Delaunay the Voronoi vertices coincide with the face :func:`circumcenters`, and interior
    nodes receive bounded dual cells while boundary nodes receive unbounded (``-1``) regions.

    Args:
        mesh: The mesh to read.

    Returns:
        scipy.spatial.Voronoi: The Voronoi object; ``.points`` are the node coordinates.
    """
    points = np.column_stack([mesh.node_x, mesh.node_y])
    return Voronoi(points)

exterior_edges(mesh) #

Indices of the edges on the mesh boundary.

A boundary edge is one whose edge_face_connectivity has a fill on one side (only one adjacent face).

A dangling edge -- one the file lists that no face implies -- has a fill on both sides, so it is reported here too. That is deliberate: such an edge bounds no face, which makes it part of the mesh's boundary in the only sense this function can measure. A placeholder row (a padded edge the file declared but never defined, or a repeat of an edge already listed) also has two fills, and that one is not reported, because it describes no edge at all.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required

Returns:

Type Description
NDArray

np.ndarray: Boundary edge indices (into the mesh's edge ordering).

Source code in src/pyramids/netcdf/ugrid/graph.py
def exterior_edges(mesh: Mesh2d) -> np.typing.NDArray:
    """Indices of the edges on the mesh boundary.

    A boundary edge is one whose ``edge_face_connectivity`` has a fill on one side (only one
    adjacent face).

    A **dangling** edge -- one the file lists that no face implies -- has a fill on *both*
    sides, so it is reported here too. That is deliberate: such an edge bounds no face, which
    makes it part of the mesh's boundary in the only sense this function can measure. A
    placeholder row (a padded edge the file declared but never defined, or a repeat of an edge
    already listed) also has two fills, and that one is **not** reported, because it describes
    no edge at all.

    Args:
        mesh: The mesh to read.

    Returns:
        np.ndarray: Boundary edge indices (into the mesh's edge ordering).
    """
    edge_node, edge_face, _ = _edge_tables(mesh)
    # A placeholder row (a padded edge the file declared but did not define) touches no
    # face, so it looks exterior on both sides. It is not a boundary edge of the mesh.
    exterior = np.any(edge_face == -1, axis=1) & _real_edges(edge_node)
    return np.flatnonzero(exterior).astype(np.intp)

exterior_faces(mesh) #

Indices of the faces that touch the mesh boundary.

A boundary face is one attached to at least one exterior edge.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required

Returns:

Type Description
NDArray

np.ndarray: Sorted unique boundary face indices.

Source code in src/pyramids/netcdf/ugrid/graph.py
def exterior_faces(mesh: Mesh2d) -> np.typing.NDArray:
    """Indices of the faces that touch the mesh boundary.

    A boundary face is one attached to at least one exterior edge.

    Args:
        mesh: The mesh to read.

    Returns:
        np.ndarray: Sorted unique boundary face indices.
    """
    _, edge_face, _ = _edge_tables(mesh)
    exterior = np.any(edge_face == -1, axis=1)
    faces = edge_face[exterior]
    valid = faces[faces != -1]
    return np.unique(valid).astype(np.intp)

bounding_polygon(mesh) #

Dissolved outline of the mesh as a shapely (Multi)Polygon, holes preserved.

Unions the face polygons so shared interior edges dissolve away; the result's exterior ring is the mesh outline and any mesh hole appears as an interior ring. Its area equals the mesh footprint. The returned geometry carries no CRS — the caller attaches one.

Parameters:

Name Type Description Default
mesh Mesh2d

The mesh to read.

required

Returns:

Type Description
Any

shapely.Polygon | shapely.MultiPolygon: The dissolved outline.

Source code in src/pyramids/netcdf/ugrid/graph.py
def bounding_polygon(mesh: Mesh2d) -> Any:
    """Dissolved outline of the mesh as a shapely (Multi)Polygon, holes preserved.

    Unions the face polygons so shared interior edges dissolve away; the result's exterior
    ring is the mesh outline and any mesh hole appears as an interior ring. Its area equals the
    mesh footprint. The returned geometry carries no CRS — the caller attaches one.

    Args:
        mesh: The mesh to read.

    Returns:
        shapely.Polygon | shapely.MultiPolygon: The dissolved outline.
    """
    return unary_union(list(_face_polygons(mesh)))