From 6c0a6c069cefef3234f8781ae2bf7d67762f00d5 Mon Sep 17 00:00:00 2001
From: Andreas Kloeckner <inform@tiker.net>
Date: Mon, 29 Apr 2024 12:26:48 -0500
Subject: [PATCH] Add, test remove_unused_vertices, use in Firedrake tests

---
 meshmode/mesh/processing.py    | 34 ++++++++++++++++++++++++++++++++++
 test/test_firedrake_interop.py |  3 ++-
 test/test_mesh.py              | 24 ++++++++++++++++++++++++
 3 files changed, 60 insertions(+), 1 deletion(-)

diff --git a/meshmode/mesh/processing.py b/meshmode/mesh/processing.py
index 90fa2e8b..a61349da 100644
--- a/meshmode/mesh/processing.py
+++ b/meshmode/mesh/processing.py
@@ -56,6 +56,8 @@ __doc__ = """
 .. autofunction:: map_mesh
 .. autofunction:: affine_map
 .. autofunction:: rotate_mesh_around_axis
+
+.. autofunction:: remove_unused_vertices
 """
 
 
@@ -1594,4 +1596,36 @@ def make_mesh_grid(
 
 # }}}
 
+
+# {{{ remove_unused_vertices
+
+def remove_unused_vertices(mesh: Mesh) -> Mesh:
+    if mesh.vertices is None:
+        raise ValueError("mesh must have vertices")
+
+    def not_none(vi: Optional[np.ndarray]) -> np.ndarray:
+        if vi is None:
+            raise ValueError("mesh element groups must have vertex indices")
+        return vi
+
+    used_vertices = np.unique(np.sort(np.concatenate([
+        not_none(grp.vertex_indices).reshape(-1)
+        for grp in mesh.groups
+    ])))
+
+    used_flags: np.ndarray = np.zeros(mesh.nvertices, dtype=np.bool_)
+    used_flags[used_vertices] = 1
+    new_vertex_indices = np.cumsum(used_flags, dtype=mesh.vertex_id_dtype) - 1
+    new_vertex_indices[~used_flags] = -1
+
+    return replace(
+        mesh,
+        vertices=mesh.vertices[:, used_flags],
+        groups=tuple(
+            replace(grp, vertex_indices=new_vertex_indices[grp.vertex_indices])
+            for grp in mesh.groups
+        ))
+
+# }}}
+
 # vim: foldmethod=marker
diff --git a/test/test_firedrake_interop.py b/test/test_firedrake_interop.py
index b9458d80..d32419b1 100644
--- a/test/test_firedrake_interop.py
+++ b/test/test_firedrake_interop.py
@@ -78,7 +78,8 @@ def fspace_degree(request):
 
 def make_mm_mesh(name: str) -> Mesh:
     from meshmode.mesh.io import read_gmsh
-    return read_gmsh(name)
+    from meshmode.mesh.processing import remove_unused_vertices
+    return remove_unused_vertices(read_gmsh(name))
 
 
 def make_firedrake_mesh(name: str):
diff --git a/test/test_mesh.py b/test/test_mesh.py
index 2d54fafa..22cd615e 100644
--- a/test/test_mesh.py
+++ b/test/test_mesh.py
@@ -256,6 +256,30 @@ def test_mesh_copy():
     mesh.copy()
 
 
+def test_remove_unused_vertices():
+    mesh = mgen.generate_box_mesh(3*(np.linspace(0, 1, 5),))
+
+    assert mesh.vertices is not None
+
+    mesh2 = replace(
+        mesh,
+        vertices=np.concatenate([np.zeros((3, 1)), mesh.vertices], axis=1),
+        groups=tuple(
+            replace(
+                grp,
+                vertex_indices=grp.vertex_indices + 1
+                )
+            for grp in mesh.groups
+        ))
+
+    mesh3 = mproc.remove_unused_vertices(mesh2)
+
+    assert np.array_equal(mesh3.vertices, mesh.vertices)
+    assert np.array_equal(
+                  mesh3.groups[0].vertex_indices,
+                  mesh.groups[0].vertex_indices)
+
+
 # {{{ as_python stringification
 
 def test_mesh_as_python():
-- 
GitLab