subset

def submesh(mesh, *, faces_index=None, vertex_index=None):
    """Re-imlementation of trimesh.submesh that is faster for our use case.

    Notably we:
     - ignore normals (possibly needed) and visuals (definitely not needed)
     - allow only one set of faces to be passed
     - return vertices and faces instead of a new mesh
     - make as few copies as possible
     - allow passing vertex indices instead of faces

    This function is 5-10x faster than trimesh.submesh for our use case.
    Note that the speed of this function was never the bottleneck though,
    it's about the memory footprint.
    See https://github.com/navis-org/navis/issues/154.

    Parameters
    ----------
    mesh :          trimesh.Trimesh
                    Mesh to submesh.
    faces_index :   array-like
                    Indices of faces to keep.
    vertex_index :  array-like
                    Indices of vertices to keep.

    Returns
    -------
    vertices :  np.ndarray
                Vertices of submesh.
    faces :     np.ndarray
                Faces of submesh.

    """
    if faces_index is None and vertex_index is None:
        raise ValueError("Either `faces_index` or `vertex_index` must be provided.")
    elif faces_index is not None and vertex_index is not None:
        raise ValueError("Only one of `faces_index` or `vertex_index` can be provided.")

    # First check if we can return either an empty mesh or the original mesh right away
    if faces_index is not None:
        if len(faces_index) == 0:
            return np.array([]), np.array([])
        elif len(faces_index) == len(mesh.faces):
            if len(np.unique(faces_index)) == len(mesh.faces):
                return mesh.vertices.copy(), mesh.faces.copy()
    else:
        if len(vertex_index) == 0:
            return np.array([]), np.array([])
        elif len(vertex_index) == len(mesh.vertices):
            if len(np.unique(vertex_index)) == len(mesh.vertices):
                return mesh.vertices.copy(), mesh.faces.copy()

    # Use a view of the original data
    original_faces = mesh.faces.view(np.ndarray)
    original_vertices = mesh.vertices.view(np.ndarray)

    # If we're starting with vertices, find faces that contain at least one of our vertices
    # This way we will also make sure to drop unreferenced vertices
    if vertex_index is not None:
        faces_index = np.arange(len(original_faces))[
            np.isin(original_faces, vertex_index).all(axis=1)
        ]

    # Get unique vertices in the to-be-kept faces
    faces = original_faces[faces_index]
    unique = np.unique(faces.reshape(-1))

    # Generate a mask for the vertices
    # (using int32 here since we're unlikey to have more than 2B vertices)
    mask = np.arange(len(original_vertices), dtype=np.int32)

    # Remap the vertices to the new indices
    mask[unique] = np.arange(len(unique))

    # Grab the vertices in the order they are referenced
    vertices = original_vertices[unique].copy()

    # Remap the faces to the new vertex indices
    # (making a copy to allow `mask` to be garbage collected)
    faces = mask[faces].copy()

    return vertices, faces