ordito.boundary

def boundary_edges(
    vertices: wp.array[wp.vec3],
    faces: wp.array[wp.int32],
    *,
    mates: wp.array[wp.int32] | None = None,
) -> odt.Array2dInt32:
    """
    Undirected boundary edges (each row min-first), without preserving orientation.

    An edge belongs to the boundary when it appears exactly once among the sorted edges of
    all triangles.

    Parameters
    ----------
    vertices
        ``(n_vertices,)`` vertex positions; only the device and the length are read. Every face
        index must be below ``n_vertices``.
    faces
        ``(3 * n_faces,)`` face index buffer.
    mates
        ``(3 * n_faces,)`` halfedge mate codes of ``faces``, from
        [`halfedge_mates`][ordito.halfedge.halfedge_mates], to spare the sort that otherwise pairs
        the halfedges. Not checked against ``faces``.
        [`Trimesh.halfedge_mates`][ordito.mesh.Trimesh.halfedge_mates] is the cached one to pass.

    Returns
    -------
    odt.Array2dInt32
        ``(n_boundary, 2)`` sorted boundary edges on ``faces.device``. Empty ``(0, 2)`` when the
        mesh has no boundary.

    Raises
    ------
    RuntimeError
        If ``vertices`` and ``faces`` are not on one device.
    RuntimeError
        If ``vertices``, ``faces`` and ``mates`` are not all on one device.

    Examples
    --------
    A 3 x 3 vertex grid has 8 boundary edges around its rim; a closed cube has none:

    ```python
    import ordito as od

    vertices, faces = od.creation.grid(count=(3, 3))
    print(od.boundary.boundary_edges(vertices, faces).shape)  # (8, 2)
    vertices, faces = od.creation.box()
    print(od.boundary.boundary_edges(vertices, faces).shape)  # (0, 2)
    ```
    """
    require_same_device(vertices=vertices, faces=faces, mates=mates)
    return _boundary_edges_impl(faces, vertices.size, mates, oriented=False)

def oriented_boundary_edges(
    vertices: wp.array[wp.vec3],
    faces: wp.array[wp.int32],
    *,
    mates: wp.array[wp.int32] | None = None,
    key_order: wp.array[wp.int32] | None = None,
) -> odt.Array2dInt32:
    """
    Directed boundary edges, preserving the orientation from the face winding.

    Boundary edges are detected on the sorted edges (appearing exactly once), but the
    directed ``(i, j)`` pairs are returned.

    Parameters
    ----------
    vertices
        ``(n_vertices,)`` vertex positions; only the device and the length are read. Every face
        index must be below ``n_vertices``.
    faces
        ``(3 * n_faces,)`` face index buffer.
    mates
        ``(3 * n_faces,)`` halfedge mate codes of ``faces``, from
        [`halfedge_mates`][ordito.halfedge.halfedge_mates], to spare the sort that otherwise pairs
        the halfedges. Not checked against ``faces``.
        [`Trimesh.halfedge_mates`][ordito.mesh.Trimesh.halfedge_mates] is the cached one to pass.
    key_order
        ``(3 * n_faces,)`` key order the pairing of ``mates`` sorted the halfedges into, from
        [`halfedge_mates`][ordito.halfedge.halfedge_mates] with ``return_key_order``, to spare the
        sort that otherwise orders the boundary rows; read only with ``mates``. ``None`` (what a
        bucketed pairing returns) sorts them.

    Returns
    -------
    odt.Array2dInt32
        ``(n_boundary, 2)`` directed boundary edges on ``faces.device``. Empty
        ``(0, 2)`` when the mesh has no boundary.

    Raises
    ------
    RuntimeError
        If ``vertices``, ``faces``, ``mates`` and ``key_order`` are not all on one device.

    Examples
    --------
    The single triangle ``(0, 1, 2)`` is its own boundary, with its edges in winding order:

    ```python
    import warp as wp

    import ordito as od

    vertices = wp.array([[0.0, 0.0, 0.0], [1.0, 0.0, 0.0], [0.0, 1.0, 0.0]], dtype=wp.vec3)
    faces = wp.array([0, 1, 2], dtype=wp.int32)
    edges = od.boundary.oriented_boundary_edges(vertices, faces)
    print(sorted(edges.numpy().tolist()))  # [[0, 1], [1, 2], [2, 0]]
    ```
    """
    require_same_device(vertices=vertices, faces=faces, mates=mates, key_order=key_order)
    return _boundary_edges_impl(faces, vertices.size, mates, oriented=True, key_order=key_order)

def boundary_loops(
    vertices: wp.array[wp.vec3],
    faces: wp.array[wp.int32],
    *,
    pinches: Literal["split", "fan", "gap"] = "split",
    copy: bool = False,
    mates: wp.array[wp.int32] | None = None,
) -> list[wp.array[wp.int32]]:
    """
    Ordered vertex-index loops along each mesh boundary.

    Each boundary loop is a simple cycle in the directed-boundary-edge successor graph (on a
    manifold boundary every boundary vertex has exactly one outgoing boundary edge). Loops are
    grouped via connected-component labeling, then each vertex's ordinal position within its
    loop is ranked by following the successor chain from itself to its loop's canonical start
    (the smallest vertex index in the loop) — entirely GPU-parallel, mirroring
    ``igl::boundary_loop`` (its first overload).

    All loops are found in one batched pass
    ([`boundary_loops_with_offsets`][ordito.boundary.boundary_loops_with_offsets]); this is
    [`split`][ordito.array.split] over its packed result.

    !!! note "The returned arrays are views"
        Each loop slices the single packed buffer ``boundary_loops_with_offsets`` produced, which
        costs no device memory and no launches. Two consequences: holding on to a single loop keeps
        the *whole* buffer alive, and writing into one loop writes into the shared allocation. Pass
        ``copy=True`` for independent buffers.

    Parameters
    ----------
    vertices
        ``(n_vertices,)`` vertex positions; only the count is used, as the successor-array
        size and the edge-key radix, so every face index must be below it.
    faces
        ``(3 * n_faces,)`` face index buffer.
    pinches
        How a loop passes a vertex where the rim meets itself; see
        [`boundary_loops_with_offsets`][ordito.boundary.boundary_loops_with_offsets].
    copy
        Return independent buffers instead of views into the packed result.
    mates
        ``(3 * n_faces,)`` halfedge mate codes of ``faces``, from
        [`halfedge_mates`][ordito.halfedge.halfedge_mates], to spare the sort that otherwise pairs
        the halfedges. Not checked against ``faces``.
        [`Trimesh.halfedge_mates`][ordito.mesh.Trimesh.halfedge_mates] is the cached one to pass.

    Returns
    -------
    list[wp.array[wp.int32]]
        One array per boundary loop, each holding the ordered vertex indices around that loop
        (not repeating the start vertex), on ``faces.device``. Empty list when the mesh has no
        boundary.

    Raises
    ------
    ValueError
        If ``pinches`` is not ``"split"``, ``"fan"`` or ``"gap"``.
    RuntimeError
        If ``vertices`` and ``faces`` are not on one device.
    RuntimeError
        If ``vertices``, ``faces`` and ``mates`` are not all on one device.

    Notes
    -----
    A vertex where the rim meets itself (a non-manifold **pinch** point, as deleting two faces that
    share only a corner leaves) is walked by halfedge rather than by vertex, so every boundary edge
    appears exactly once and every consecutive pair is a real boundary edge; ``pinches`` decides
    which edge the walk leaves such a vertex by. That needs an edge-manifold mesh; on one that is
    not, the loops are bounded but unspecified.

    Orientability is *not* assumed. On a non-orientable surface the winding cannot orient the
    boundary globally -- at the seam two boundary edges leave the same vertex -- so the loops are
    recovered from the undirected edges instead, and their direction is then an arbitrary but
    reproducible choice rather than the face winding's. Only the direction is arbitrary; the loops
    themselves are exact. Costs one extra host readback to detect the case.

    See Also
    --------
    [`boundary_loops_with_offsets`][ordito.boundary.boundary_loops_with_offsets]
        The packed form this splits: the same loops in one buffer plus offsets.
    [`longest_boundary_loop`][ordito.boundary.longest_boundary_loop]
    [`oriented_boundary_edges`][ordito.boundary.oriented_boundary_edges]
    ``igl.boundary_loop_all``

    Examples
    --------
    Two separate triangles have two boundary loops; each starts at its smallest vertex index:

    ```python
    import warp as wp

    import ordito as od

    vertices = wp.array(
        [[0.0, 0.0, 0.0], [1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [3.0, 0.0, 0.0], [4.0, 0.0, 0.0],
         [3.0, 1.0, 0.0]],
        dtype=wp.vec3,
    )
    faces = wp.array([0, 1, 2, 5, 3, 4], dtype=wp.int32)
    loops = od.boundary.boundary_loops(vertices, faces)
    print([loop.numpy().tolist() for loop in loops])  # [[0, 1, 2], [3, 4, 5]]
    ```
    """
    require_same_device(vertices=vertices, faces=faces, mates=mates)
    flat_loops, offsets = boundary_loops_with_offsets(vertices, faces, pinches=pinches, mates=mates)
    # The offsets' length is host metadata: with one loop they are ``[0, len(flat_loops)]``, and the
    # split needs no readback of them.
    host_offsets = [0, flat_loops.size] if offsets.size == 2 else offsets
    return od.array.split(flat_loops, host_offsets, copy=copy)

def boundary_loops_with_offsets(
    vertices: wp.array[wp.vec3],
    faces: wp.array[wp.int32],
    *,
    pinches: Literal["split", "fan", "gap"] = "split",
    mates: wp.array[wp.int32] | None = None,
) -> tuple[wp.array[wp.int32], wp.array[wp.int32]]:
    """
    Every boundary loop at once, packed into one buffer plus per-loop offsets.

    Same loops, same order, same contents as
    [`boundary_loops`][ordito.boundary.boundary_loops] — but with no per-loop Python and no
    per-loop allocation, which is the only form whose cost is independent of the loop *count*.
    Prefer it when a mesh has many small holes or when the loops feed straight into another
    batched kernel, as [`ordito.holes`][ordito.holes] does. The loops come out as
    [`successor_cycles`][ordito.graph.successor_cycles] orders a successor graph's cycles --
    each from its smallest vertex, the loops in ascending order of it -- ranked by
    [`successor_cycles_from_table`][ordito.graph.successor_cycles_from_table] straight off the
    boundary's successor table.

    Where two holes touch at one vertex (a *pinch*), the boundary halfedges into that vertex have
    two ways out, and ``pinches`` picks the loops:

    - ``"split"`` (default): every loop is simple. The walk of ``"fan"`` is cut into cycles at the
      vertices a loop passes twice, which are the loops a hole filler closes one by one (a
      triangulation of a loop that repeats a vertex is not a surface). The cut costs one readback
      of the boundary, and only a pinched boundary pays it.
    - ``"fan"``: the walk stays in the fan of faces it arrived along, which is the boundary of the
      surface with the pinch vertex split once per fan; two holes touching at a vertex can come
      back as one loop passing it twice.
    - ``"gap"``: the walk crosses the gap to the other fan: each loop is the rim of one empty
      region, so two holes touching at a vertex stay two loops, and a loop passes a vertex twice
      only where one hole meets itself there. At a vertex with three fans or more the walk stays
      in its fan, as ``"fan"`` does.

    On a boundary with no pinch the three agree.

    Parameters
    ----------
    vertices
        ``(n_vertices,)`` vertex positions; only the count is used, as the successor-array
        size and the edge-key radix, so every face index must be below it.
    faces
        ``(3 * n_faces,)`` face index buffer.
    pinches
        ``"split"`` (default), ``"fan"`` or ``"gap"``: how a loop leaves a pinch vertex, as
        described above.
    mates
        ``(3 * n_faces,)`` halfedge mate codes of ``faces``, from
        [`halfedge_mates`][ordito.halfedge.halfedge_mates], to spare the sort that otherwise pairs
        the halfedges. Not checked against ``faces``.
        [`Trimesh.halfedge_mates`][ordito.mesh.Trimesh.halfedge_mates] is the cached one to pass.

    Returns
    -------
    flat_loops
        ``(n_loop_vertices,)`` concatenated ordered vertex indices of every loop, on
        ``faces.device``.
    offsets
        ``(n_loops + 1,)`` total-terminated offsets: loop ``i`` occupies
        ``flat_loops[offsets[i] : offsets[i + 1]]``, and ``offsets == [0]`` when the mesh has no
        boundary.

    Raises
    ------
    ValueError
        If ``pinches`` is not ``"split"``, ``"fan"`` or ``"gap"``.
    RuntimeError
        If ``vertices`` and ``faces`` are not on one device.
    RuntimeError
        If ``vertices``, ``faces`` and ``mates`` are not all on one device.

    See Also
    --------
    [`boundary_loops`][ordito.boundary.boundary_loops]
        The list form: [`split`][ordito.array.split] over this result.
    [`longest_boundary_loop`][ordito.boundary.longest_boundary_loop]
    [`successor_cycles`][ordito.graph.successor_cycles]

    Examples
    --------
    The same loops as [`boundary_loops`][ordito.boundary.boundary_loops], packed into one buffer:

    ```python
    import warp as wp

    import ordito as od

    vertices = wp.array(
        [[0.0, 0.0, 0.0], [1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [3.0, 0.0, 0.0], [4.0, 0.0, 0.0],
         [3.0, 1.0, 0.0]],
        dtype=wp.vec3,
    )
    faces = wp.array([0, 1, 2, 5, 3, 4], dtype=wp.int32)
    flat_loops, offsets = od.boundary.boundary_loops_with_offsets(vertices, faces)
    print(flat_loops.numpy().tolist())  # [0, 1, 2, 3, 4, 5]
    print(offsets.numpy().tolist())  # [0, 3, 6]
    ```
    """
    require_same_device(vertices=vertices, faces=faces, mates=mates)
    if pinches not in ("split", "fan", "gap"):
        raise ValueError(f"pinches must be 'split', 'fan' or 'gap', got {pinches!r}")
    n_faces = faces.size // 3
    device = faces.device
    if n_faces == 0:
        return _launch.empty_packed(wp.int32, device)

    n_vertices = vertices.size
    # One boundary detection for every edge view below; ``boundary_edges`` /
    # ``oriented_boundary_edges`` / ``boundary_vertex_indices`` would each redo the pairing.
    boundary = _BoundaryHalfedges(faces, n_vertices, mates)
    has_seam, has_pinch = _boundary_defects(boundary, n_vertices)
    if boundary.count() == 0:
        return _launch.empty_packed(wp.int32, device)
    if has_pinch:
        if pinches == "split":
            return _split_at_pinches(*_pinched_boundary_cycles(boundary, n_vertices))
        return _pinched_boundary_cycles(boundary, n_vertices, cross_gaps=pinches == "gap")
    if has_seam:
        return _unoriented_boundary_cycles(boundary.edges(sort_pair=True), n_vertices)
    # With neither defect no boundary vertex has two boundary edges out or three incident, so the
    # directed rows are cycles over the boundary vertices -- plus, on a mesh that is not
    # edge-manifold, chains running into a vertex with no edge out, which the ranking drops.
    tails, next_node = boundary.successors(n_vertices)
    return successor_cycles_from_table(tails, next_node)

def loop_perimeters(
    vertices: wp.array[wp.vec3], loops: Sequence[wp.array[wp.int32]]
) -> wp.array[wp.float32]:
    """
    Perimeter of every closed loop, in one launch.

    Takes what [`boundary_loops`][ordito.boundary.boundary_loops] returns -- a list of vertex-index
    cycles whose last entry joins back to the first -- and measures them all together, so the cost
    is one launch and one packing pass rather than one call per loop. Equivalent to
    [`polyline_length`][ordito.polyline.polyline_length] with ``closed=True`` on each loop's
    gathered positions, and the segmented form is why the fill's ``preserve_largest_hole`` can rank
    every rim for the price of one readback. Under
    [`is_deterministic`][ordito.config.is_deterministic] each loop is summed in one fixed order,
    so two rims of equal length rank the same way on every run.

    Parameters
    ----------
    vertices
        ``(n_vertices,)`` mesh vertex positions.
    loops
        Closed vertex-index cycles. An empty sequence gives an empty result; a loop of fewer than
        two entries measures ``0``.

    Returns
    -------
    wp.array[wp.float32]
        ``(n_loops,)`` one perimeter per loop, in the order given, on ``vertices.device``.

    Raises
    ------
    TypeError
        If any loop is not a rank-1 ``wp.int32`` array.
    RuntimeError
        If ``vertices`` and ``loops`` are not all on one device.

    See Also
    --------
    [`loop_directed_areas`][ordito.boundary.loop_directed_areas]
        The vector measure of the same loops: norm is the spanned area, direction is its normal.
    [`boundary_loops`][ordito.boundary.boundary_loops]
        Produces the loops.
    [`polyline_length`][ordito.polyline.polyline_length]
        The single-loop form, over positions rather than indices.

    Examples
    --------
    The rim of a unit grid has perimeter 4:

    ```python
    import ordito as od

    vertices, faces = od.creation.grid(count=(3, 3))
    loops = od.boundary.boundary_loops(vertices, faces)
    perimeters = od.boundary.loop_perimeters(vertices, loops)
    print(f"{perimeters.numpy()[0]:.4f}")  # 4.0000
    ```
    """
    require_same_device(vertices=vertices, loops=loops)
    packed = _pack_loop_segments(loops)
    if packed is None:
        return _launch.empty(0, dtype=wp.float32, device=vertices.device)
    flat_loops, offsets = packed
    return _launch_loop_measure(
        kernel_boundary.loop_perimeters, wp.float32, vertices, flat_loops, offsets
    )

def loop_perimeters_from_offsets(
    vertices: wp.array[wp.vec3],
    flat_loops: wp.array[wp.int32],
    offsets: wp.array[wp.int32],
    *,
    validate: bool = True,
) -> wp.array[wp.float32]:
    """
    Perimeter of every loop, taking the loops in the packed form rather than as a list.

    Same measure as [`loop_perimeters`][ordito.boundary.loop_perimeters] and the same launch; the
    two arguments after ``vertices`` are exactly what
    [`boundary_loops_with_offsets`][ordito.boundary.boundary_loops_with_offsets] returns, in that
    order. Reach for this whenever the loops are already packed, which is the form
    [`ordito.holes`][ordito.holes] carries them in throughout: the list form would have to be
    split back out and repacked to be measured, and the split is a per-loop Python object.

    Parameters
    ----------
    vertices
        ``(n_vertices,)`` mesh vertex positions.
    flat_loops
        ``(n_loop_vertices,)`` concatenated ordered vertex indices of every loop.
    offsets
        ``(n_loops + 1,)`` total-terminated offsets: loop ``i`` is
        ``flat_loops[offsets[i] : offsets[i + 1]]``.
    validate
        When ``True`` (default), check that ``offsets`` is non-decreasing from ``0`` to
        ``flat_loops``' length before launching -- a real cost (a host readback of ``offsets``),
        paid because a hand-built or stale pair otherwise drives the kernel's per-loop index past
        ``flat_loops``'s end with no exception raised. Pass ``False`` only when the pair is known
        correct by construction, as [`ordito.holes`][ordito.holes]'s internal packing already is.

    Returns
    -------
    wp.array[wp.float32]
        ``(n_loops,)`` one perimeter per loop, in the packed order, on ``vertices.device``.

    Raises
    ------
    RuntimeError
        If ``vertices``, ``flat_loops`` and ``offsets`` are not all on one device.
    ValueError
        If ``validate`` and the packed pair is malformed -- see ``validate`` above.

    See Also
    --------
    [`loop_perimeters`][ordito.boundary.loop_perimeters]
        The list form, which packs and then launches the same kernel.
    [`loop_directed_areas_from_offsets`][ordito.boundary.loop_directed_areas_from_offsets]
        The vector measure of the same packed loops.

    Examples
    --------
    ```python
    import ordito as od

    vertices, faces = od.creation.grid(count=(3, 3))
    flat_loops, offsets = od.boundary.boundary_loops_with_offsets(vertices, faces)
    perimeters = od.boundary.loop_perimeters_from_offsets(vertices, flat_loops, offsets)
    print(perimeters.shape)  # (1,)
    print(f"{perimeters.numpy()[0]:.4f}")  # 4.0000
    ```
    """
    require_same_device(vertices=vertices, flat_loops=flat_loops, offsets=offsets)
    if validate:
        _validate_packed_loops(flat_loops, offsets, "loop_perimeters_from_offsets")
    return _launch_loop_measure(
        kernel_boundary.loop_perimeters, wp.float32, vertices, flat_loops, offsets
    )

def loop_directed_areas(
    vertices: wp.array[wp.vec3], loops: Sequence[wp.array[wp.int32]]
) -> wp.array[wp.vec3]:
    """
    Directed area vector of every loop, in one launch.

    Half the sum of ``p_i x p_{i+1}`` around each cycle. Its **norm** is the area of the planar
    polygon the loop spans and its **direction** is that polygon's normal, oriented by the loop's
    own winding -- so it answers "how big is this hole" and "which way does it face" at once, and
    the sign flips if the loop is reversed. A vector rather than a scalar for that reason: the area
    alone loses the orientation, which is what tells an outer boundary from an inner one.

    Origin-independent, because the cross products of a closed ring cancel any shift of the origin,
    so no centroid pass is needed. Exact for a planar loop; for a non-planar one it is the area of
    the loop's projection onto the plane normal to the result, which is the standard convention.

    Parameters
    ----------
    vertices
        ``(n_vertices,)`` mesh vertex positions.
    loops
        Closed vertex-index cycles, as [`boundary_loops`][ordito.boundary.boundary_loops] returns.

    Returns
    -------
    wp.array[wp.vec3]
        ``(n_loops,)`` one directed area per loop, in the order given, on ``vertices.device``.

    Raises
    ------
    TypeError
        If any loop is not a rank-1 ``wp.int32`` array.
    RuntimeError
        If ``vertices`` and ``loops`` are not all on one device.

    See Also
    --------
    [`loop_perimeters`][ordito.boundary.loop_perimeters]
        The scalar measure of the same loops.
    [`polyline_normal`][ordito.polyline.polyline_normal]
        The single-loop direction, normalized and over positions.

    Examples
    --------
    The rim of a unit grid in the ``z = 0`` plane spans area 1; reversing the loop flips the vector:

    ```python
    import warp as wp

    import ordito as od

    vertices, faces = od.creation.grid(count=(3, 3))
    loop = od.boundary.boundary_loops(vertices, faces)[0]
    reversed_loop = wp.array(loop.numpy()[::-1].copy(), dtype=wp.int32)
    areas = od.boundary.loop_directed_areas(vertices, [loop, reversed_loop])
    print(areas.numpy().round(4).tolist())  # [[0.0, 0.0, 1.0], [0.0, 0.0, -1.0]]
    ```
    """
    require_same_device(vertices=vertices, loops=loops)
    packed = _pack_loop_segments(loops)
    if packed is None:
        return _launch.empty(0, dtype=wp.vec3, device=vertices.device)
    flat_loops, offsets = packed
    return _launch_loop_measure(
        kernel_boundary.loop_directed_areas, wp.vec3, vertices, flat_loops, offsets
    )

def loop_directed_areas_from_offsets(
    vertices: wp.array[wp.vec3],
    flat_loops: wp.array[wp.int32],
    offsets: wp.array[wp.int32],
    *,
    validate: bool = True,
) -> wp.array[wp.vec3]:
    """
    Directed area vector of every loop, taking the loops in the packed form rather than as a list.

    The packed counterpart of
    [`loop_directed_areas`][ordito.boundary.loop_directed_areas], exactly as
    [`loop_perimeters_from_offsets`][ordito.boundary.loop_perimeters_from_offsets] is of
    [`loop_perimeters`][ordito.boundary.loop_perimeters]; the arguments after ``vertices`` are what
    [`boundary_loops_with_offsets`][ordito.boundary.boundary_loops_with_offsets] returns.

    Parameters
    ----------
    vertices
        ``(n_vertices,)`` mesh vertex positions.
    flat_loops
        ``(n_loop_vertices,)`` concatenated ordered vertex indices of every loop.
    offsets
        ``(n_loops + 1,)`` total-terminated offsets: loop ``i`` is
        ``flat_loops[offsets[i] : offsets[i + 1]]``.
    validate
        When ``True`` (default), check that ``offsets`` is non-decreasing from ``0`` to
        ``flat_loops``' length before launching -- see
        [`loop_perimeters_from_offsets`][ordito.boundary.loop_perimeters_from_offsets]'s docstring
        for the cost and the reasoning.

    Returns
    -------
    wp.array[wp.vec3]
        ``(n_loops,)`` one directed area per loop, in the packed order, on ``vertices.device``.

    Raises
    ------
    RuntimeError
        If ``vertices``, ``flat_loops`` and ``offsets`` are not all on one device.
    ValueError
        If ``validate`` and the packed pair is malformed -- see ``validate`` above.

    See Also
    --------
    [`loop_directed_areas`][ordito.boundary.loop_directed_areas]
        The list form, which packs and then launches the same kernel.
    [`loop_perimeters_from_offsets`][ordito.boundary.loop_perimeters_from_offsets]
        The scalar measure of the same packed loops.

    Examples
    --------
    ```python
    import ordito as od

    vertices, faces = od.creation.grid(count=(3, 3))
    flat_loops, offsets = od.boundary.boundary_loops_with_offsets(vertices, faces)
    areas = od.boundary.loop_directed_areas_from_offsets(vertices, flat_loops, offsets)
    print(areas.numpy().round(4).tolist())  # [[0.0, 0.0, 1.0]]
    ```
    """
    require_same_device(vertices=vertices, flat_loops=flat_loops, offsets=offsets)
    if validate:
        _validate_packed_loops(flat_loops, offsets, "loop_directed_areas_from_offsets")
    return _launch_loop_measure(
        kernel_boundary.loop_directed_areas, wp.vec3, vertices, flat_loops, offsets
    )

def longest_boundary_loop(
    vertices: wp.array[wp.vec3],
    faces: wp.array[wp.int32],
    *,
    mates: wp.array[wp.int32] | None = None,
) -> wp.array[wp.int32]:
    """
    Ordered vertex-index loop along the longest mesh boundary.

    *Longest* counts vertices, not perimeter, as ``igl::boundary_loop`` does; of loops with equal
    vertex counts the one whose smallest vertex index is lowest wins.

    The name carries the *longest*, because the plural
    [`boundary_loops`][ordito.boundary.boundary_loops] returns every one and two public names
    differing by a single character are a defect even when both are correct.

    Parameters
    ----------
    vertices, faces
        Forwarded to [`boundary_loops_with_offsets`][ordito.boundary.boundary_loops_with_offsets].
    mates
        ``(3 * n_faces,)`` halfedge mate codes of ``faces``, from
        [`halfedge_mates`][ordito.halfedge.halfedge_mates], to spare the sort that otherwise pairs
        the halfedges. Not checked against ``faces``.
        [`Trimesh.halfedge_mates`][ordito.mesh.Trimesh.halfedge_mates] is the cached one to pass.

    Returns
    -------
    wp.array[wp.int32]
        ``(m,)`` ordered vertex indices around the longest boundary loop, ``m`` its length, on
        ``faces.device``. An independent buffer, not a view into the packed result. Empty when the
        mesh has no boundary.

    Raises
    ------
    RuntimeError
        If ``vertices`` and ``faces`` are not on one device.
    RuntimeError
        If ``vertices``, ``faces`` and ``mates`` are not all on one device.

    See Also
    --------
    [`boundary_loops`][ordito.boundary.boundary_loops]
        Every loop, not only the longest.
    [`boundary_loops_with_offsets`][ordito.boundary.boundary_loops_with_offsets]
    ``igl.boundary_loop``

    Examples
    --------
    A unit square made of two triangles has a four-vertex rim; a much larger single triangle has
    the longer perimeter but only three rim vertices:

    ```python
    import warp as wp

    import ordito as od

    vertices = wp.array(
        [[0.0, 0.0, 0.0], [1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [1.0, 1.0, 0.0], [5.0, 0.0, 0.0],
         [9.0, 0.0, 0.0], [5.0, 4.0, 0.0]],
        dtype=wp.vec3,
    )
    faces = wp.array([0, 1, 2, 2, 1, 3, 4, 5, 6], dtype=wp.int32)
    print(od.boundary.longest_boundary_loop(vertices, faces).numpy().tolist())  # [0, 1, 3, 2]
    ```
    """
    require_same_device(vertices=vertices, faces=faces, mates=mates)
    device = faces.device
    flat_loops, offsets = boundary_loops_with_offsets(vertices, faces, mates=mates)
    n_loops = offsets.size - 1
    if n_loops == 0:
        return _launch.empty(0, dtype=wp.int32, device=device)
    if n_loops == 1:
        # The one loop is the whole packed buffer, which nothing else holds.
        return flat_loops
    # The sizes are the offsets' differences on the device, so the winner is an argmax there
    # rather than a Python scan: unpacking the loops first would read the offsets back, build one
    # array view per loop, and then recover the lengths from those views -- a cost linear in the
    # rim count for an answer that is one loop. ``-1`` is below every packed key, so the reduction
    # needs no separate seeding pass.
    best = _launch.full(1, -1, dtype=wp.int64, device=device)
    _launch.launch(
        kernel_boundary.longest_loop_key, dim=n_loops, inputs=[offsets, best], device=device
    )
    # One readback, because the key carries the winner's start in its low half as well as its
    # length in its high half -- see the kernel.
    key = int(read_scalar(best, 0))
    start = INT32_MAX - (key & 0xFFFFFFFF)
    size = key >> 32
    # Only the winner is materialized: the rest of the packed buffer is never copied.
    return _launch.clone(odt.as_dense(flat_loops[start : start + size]))

def boundary_vertex_indices(
    vertices: wp.array[wp.vec3],
    faces: wp.array[wp.int32],
    *,
    mates: wp.array[wp.int32] | None = None,
) -> wp.array[wp.int32]:
    """
    Sorted unique vertex indices lying on the mesh boundary.

    A vertex belongs to the boundary when it belongs to at least one boundary edge.

    Parameters
    ----------
    vertices
        ``(n_vertices,)`` vertex positions; only the count is used (every boundary vertex index is
        below it).
    faces
        ``(3 * n_faces,)`` face index buffer.
    mates
        ``(3 * n_faces,)`` halfedge mate codes of ``faces``, from
        [`halfedge_mates`][ordito.halfedge.halfedge_mates], to spare the sort that otherwise pairs
        the halfedges. Not checked against ``faces``.
        [`Trimesh.halfedge_mates`][ordito.mesh.Trimesh.halfedge_mates] is the cached one to pass.

    Returns
    -------
    wp.array[wp.int32]
        ``(n_boundary_vertices,)`` sorted unique boundary vertex indices on ``faces.device``.
        Empty when the mesh has no boundary.

    Raises
    ------
    RuntimeError
        If ``vertices`` and ``faces`` are not on one device.
    RuntimeError
        If ``vertices``, ``faces`` and ``mates`` are not all on one device.

    Examples
    --------
    Every vertex of a 3 x 3 grid except the centre one (index 4) is on the rim:

    ```python
    import ordito as od

    vertices, faces = od.creation.grid(count=(3, 3))
    rim = od.boundary.boundary_vertex_indices(vertices, faces)
    print(rim.numpy().tolist())  # [0, 1, 2, 3, 5, 6, 7, 8]
    ```
    """
    require_same_device(vertices=vertices, faces=faces, mates=mates)
    if faces.size // 3 == 0:
        return _launch.empty(0, dtype=wp.int32, device=faces.device)
    # The endpoints are vertex indices below ``len(vertices)``, so a per-vertex flag array and its
    # scan give the sorted unique set directly -- no edge list and no ``unique_1d``.
    n_vertices = vertices.size
    return _BoundaryHalfedges(faces, n_vertices, mates).vertex_indices(n_vertices)

def boundary_vertices(
    vertices: wp.array[wp.vec3],
    faces: wp.array[wp.int32],
    *,
    mates: wp.array[wp.int32] | None = None,
) -> wp.array[wp.vec3]:
    """
    Coordinates of the vertices lying on the mesh boundary.

    Parameters
    ----------
    vertices
        ``(n_vertices,)`` vertex positions.
    faces
        ``(3 * n_faces,)`` face index buffer.
    mates
        ``(3 * n_faces,)`` halfedge mate codes of ``faces``, from
        [`halfedge_mates`][ordito.halfedge.halfedge_mates], to spare the sort that otherwise pairs
        the halfedges. Not checked against ``faces``.
        [`Trimesh.halfedge_mates`][ordito.mesh.Trimesh.halfedge_mates] is the cached one to pass.

    Returns
    -------
    wp.array[wp.vec3]
        ``(n_boundary_vertices,)`` boundary vertex positions on ``vertices.device``,
        ordered by ascending vertex index. Empty when the mesh has no boundary.

    Raises
    ------
    RuntimeError
        If ``vertices`` and ``faces`` are not on one device.
    RuntimeError
        If ``vertices``, ``faces`` and ``mates`` are not all on one device.

    Examples
    --------
    ```python
    import ordito as od

    vertices, faces = od.creation.grid(count=(3, 3))
    rim = od.boundary.boundary_vertices(vertices, faces)
    print(rim.shape)  # (8,)
    print(rim.numpy()[0])  # [-0.5 -0.5  0. ]
    ```
    """
    require_same_device(vertices=vertices, faces=faces, mates=mates)
    if faces.size // 3 == 0:
        return _launch.empty(0, dtype=wp.vec3, device=faces.device)
    return _BoundaryHalfedges(faces, vertices.size, mates).vertex_positions(vertices)

def ears(
    faces: wp.array[wp.int32],
    *,
    n_vertices: int | None = None,
    mates: wp.array[wp.int32] | None = None,
) -> tuple[wp.array[wp.int32], wp.array[wp.int32]]:
    """
    Find ear faces (triangles with exactly two boundary edges).

    For each ear face, ``ear_opp`` is the local index of the corner opposite the non-boundary
    edge -- the ear's tip -- so the interior edge is
    ``(faces[f, (ear_opp + 1) % 3], faces[f, (ear_opp + 2) % 3])``. This is ``igl.ears``'
    numbering, edge ``i`` being the one facing corner ``i``.

    Parameters
    ----------
    faces
        ``(3 * n_faces,)`` face index buffer.
    n_vertices
        Optional vertex count, forwarded to
        [`halfedge_mates`][ordito.halfedge.halfedge_mates]: it narrows the sort that pairs the
        halfedges. Trusted, not checked: it must exceed every index in ``faces``.
    mates
        ``(3 * n_faces,)`` halfedge mate codes of ``faces``, from
        [`halfedge_mates`][ordito.halfedge.halfedge_mates], to spare the sort that otherwise pairs
        the halfedges. Not checked against ``faces``.
        [`Trimesh.halfedge_mates`][ordito.mesh.Trimesh.halfedge_mates] is the cached one to pass.

    Returns
    -------
    ear : wp.array[wp.int32]
        ``(n_ears,)`` face indices of ear triangles on ``faces.device``, ascending. Empty when no
        ears exist.
    ear_opp : wp.array[wp.int32]
        ``(n_ears,)`` local index of the corner opposite each ear face's interior edge, same length
        as ``ear``.

    Raises
    ------
    RuntimeError
        If ``faces`` and ``mates`` are not all on one device.

    See Also
    --------
    ``igl.ears``

    Examples
    --------
    A 2 x 2 vertex grid is two triangles, and each has two edges on the boundary:

    ```python
    import ordito as od

    vertices, faces = od.creation.grid(count=(2, 2))
    ear, ear_opp = od.boundary.ears(faces)
    print(ear.numpy().tolist())  # [0, 1]
    print(ear_opp.numpy().tolist())  # [1, 2]
    print(faces.numpy().tolist())  # [0, 2, 3, 0, 3, 1]
    ```
    """
    require_same_device(faces=faces, mates=mates)
    n_faces = faces.size // 3
    device = faces.device
    empty = _launch.empty(0, dtype=wp.int32, device=device)
    if n_faces == 0:
        return empty, empty

    # A boundary edge is a halfedge with no mate, read per corner: no mask, scan or readback.
    mates = mates if mates is not None else od.halfedge.halfedge_mates(faces, n_vertices)

    # Flag, scan in place, and emit at the scan's steps: ascending face order, one readback.
    inclusive = _launch.empty(n_faces, dtype=wp.int32, device=device)
    _launch.launch(kernel_boundary.mark_ears, dim=n_faces, inputs=[mates, inclusive], device=device)
    _launch.array_scan(inclusive, out_array=inclusive, inclusive=True)
    # Sizes the output: the one host readback.
    n_ears = int(read_scalar(inclusive))
    if n_ears == 0:
        return empty, empty
    # Both outputs in one allocation, every entry written by the emit.
    pair = _launch.empty(2 * n_ears, dtype=wp.int32, device=device)
    ear = _launch.view(pair, 0, n_ears)
    ear_opp = _launch.view(pair, n_ears, 2 * n_ears)
    _launch.launch(
        kernel_boundary.emit_ears,
        dim=n_faces,
        inputs=[mates, inclusive, ear, ear_opp],
        device=device,
    )
    return ear, ear_opp