Update check mesh function with non-manifold checks

Author: w1th0utnam3

splashsurf/src/reconstruction.rs

@@ -346,7 +346,7 @@ pub struct ReconstructSubcommandArgs {
     /// Optional filename for writing the octree used to partition the particles to disk
     #[arg(help_heading = ARGS_DEBUG, long, value_parser = value_parser!(PathBuf))]
     pub output_octree: Option<PathBuf>,
-    /// Enable checking the final mesh for topological problems such as holes (note that when stitching is disabled this will lead to a lot of reported problems)
+    /// Enable checking the final mesh for holes and non-manifold edges and vertices
     #[arg(
         help_heading = ARGS_DEBUG,
         long,
@@ -356,6 +356,36 @@ pub struct ReconstructSubcommandArgs {
         require_equals = true
     )]
     pub check_mesh: Switch,
+    /// Enable checking the final mesh for holes
+    #[arg(
+        help_heading = ARGS_DEBUG,
+        long,
+        default_value = "off",
+        value_name = "off|on",
+        ignore_case = true,
+        require_equals = true
+    )]
+    pub check_mesh_closed: Switch,
+    /// Enable checking the final mesh for non-manifold edges and vertices
+    #[arg(
+        help_heading = ARGS_DEBUG,
+        long,
+        default_value = "off",
+        value_name = "off|on",
+        ignore_case = true,
+        require_equals = true
+    )]
+    pub check_mesh_manifold: Switch,
+    /// Enable debug output for the check-mesh operations (has no effect if no other check-mesh option is enabled)
+    #[arg(
+        help_heading = ARGS_DEBUG,
+        long,
+        default_value = "off",
+        value_name = "off|on",
+        ignore_case = true,
+        require_equals = true
+    )]
+    pub check_mesh_debug: Switch,
 }
 
 #[derive(Copy, Clone, Debug, PartialEq, Eq, clap::ValueEnum)]
@@ -450,7 +480,9 @@ mod arguments {
     use walkdir::WalkDir;
 
     pub struct ReconstructionRunnerPostprocessingArgs {
-        pub check_mesh: bool,
+        pub check_mesh_closed: bool,
+        pub check_mesh_manifold: bool,
+        pub check_mesh_debug: bool,
         pub mesh_cleanup: bool,
         pub decimate_barnacles: bool,
         pub keep_vertices: bool,
@@ -601,7 +633,11 @@ mod arguments {
             }
 
             let postprocessing = ReconstructionRunnerPostprocessingArgs {
-                check_mesh: args.check_mesh.into_bool(),
+                check_mesh_closed: args.check_mesh.into_bool()
+                    || args.check_mesh_closed.into_bool(),
+                check_mesh_manifold: args.check_mesh.into_bool()
+                    || args.check_mesh_manifold.into_bool(),
+                check_mesh_debug: args.check_mesh_debug.into_bool(),
                 mesh_cleanup: args.mesh_cleanup.into_bool(),
                 decimate_barnacles: args.decimate_barnacles.into_bool(),
                 keep_vertices: args.keep_verts.into_bool(),
@@ -1455,11 +1491,18 @@ pub(crate) fn reconstruction_pipeline_generic<I: Index, R: Real>(
         info!("Done.");
     }
 
-    if postprocessing.check_mesh {
+    if postprocessing.check_mesh_closed
+        || postprocessing.check_mesh_manifold
+        || postprocessing.check_mesh_debug
+    {
         if let Err(err) = match (&tri_mesh, &tri_quad_mesh) {
-            (Some(mesh), None) => {
-                splashsurf_lib::marching_cubes::check_mesh_consistency(grid, &mesh.mesh)
-            }
+            (Some(mesh), None) => splashsurf_lib::marching_cubes::check_mesh_consistency(
+                grid,
+                &mesh.mesh,
+                postprocessing.check_mesh_closed,
+                postprocessing.check_mesh_manifold,
+                postprocessing.check_mesh_debug,
+            ),
             (None, Some(_mesh)) => {
                 info!("Checking for mesh consistency not implemented for quad mesh at the moment.");
                 return Ok(());
@@ -1468,7 +1511,7 @@ pub(crate) fn reconstruction_pipeline_generic<I: Index, R: Real>(
         } {
             return Err(anyhow!("{}", err));
         } else {
-            info!("Checked mesh for problems (holes, etc.), no problems were found.");
+            info!("Checked mesh for problems (holes: {}, non-manifold edges/vertices: {}), no problems were found.", postprocessing.check_mesh_closed, postprocessing.check_mesh_manifold);
         }
     }
 

splashsurf_lib/src/marching_cubes.rs

@@ -195,21 +195,39 @@ pub(crate) fn triangulate_density_map_to_surface_patch<I: Index, R: Real>(
     })
 }
 
-/// Checks the consistency of the mesh (currently only checks for holes) and returns a string with debug information in case of problems
+/// Checks the consistency of the mesh (currently checks for holes, non-manifold edges and vertices) and returns a string with debug information in case of problems
 pub fn check_mesh_consistency<I: Index, R: Real>(
     grid: &UniformGrid<I, R>,
     mesh: &TriMesh3d<R>,
+    check_closed: bool,
+    check_manifold: bool,
+    debug: bool,
 ) -> Result<(), String> {
     profile!("check_mesh_consistency");
-    let boundary_edges = mesh.find_boundary_edges();
-
-    if boundary_edges.is_empty() {
+    let edge_info = mesh.compute_edge_information();
+
+    let boundary_edges = edge_info
+        .iter()
+        .filter(|ei| ei.incident_faces == 1)
+        .map(|ei| (ei.edge, ei.face, ei.local_edge_index))
+        .collect::<Vec<_>>();
+    let non_manifold_edges = edge_info
+        .iter()
+        .filter(|ei| ei.incident_faces > 2)
+        .map(|ei| (ei.edge, ei.face, ei.local_edge_index))
+        .collect::<Vec<_>>();
+
+    let non_manifold_vertices = mesh.find_non_manifold_vertices();
+
+    if (!check_closed || boundary_edges.is_empty())
+        && (!check_manifold || (non_manifold_edges.is_empty() && non_manifold_vertices.is_empty()))
+    {
         return Ok(());
     }
 
-    let mut error_string = String::new();
-    error_string += &format!("Mesh is not closed. It has {} boundary edges (edges that are connected to only one triangle):", boundary_edges.len());
-    for (edge, tri_idx, _) in boundary_edges {
+    let add_edge_errors = |error_string: &mut String, edge: ([usize; 2], usize, usize)| {
+        let (edge, tri_idx, _) = edge;
+
         let v0 = mesh.vertices[edge[0]];
         let v1 = mesh.vertices[edge[1]];
         let center = (v0 + v1) / (R::one() + R::one());
@@ -221,13 +239,43 @@ pub fn check_mesh_consistency<I: Index, R: Real>(
             let cell_center = grid.point_coordinates(&point_index)
                 + &Vector3::repeat(grid.cell_size().times_f64(0.5));
 
-            error_string += &format!("\n\tTriangle {}, boundary edge {:?} is located in cell with {:?} with center coordinates {:?} and edge length {}.", tri_idx, edge, cell_index, cell_center, grid.cell_size());
+            *error_string += &format!("\n\tTriangle {}, boundary edge {:?} is located in cell with {:?} with center coordinates {:?} and edge length {}.", tri_idx, edge, cell_index, cell_center, grid.cell_size());
         } else {
-            error_string += &format!(
-                "\n\tCannot get cell index for boundary edge {:?} of triangle {}",
+            *error_string += &format!(
+                "\n\tCannot get cell index for edge {:?} of triangle {}",
                 edge, tri_idx
             );
         }
+    };
+
+    let mut error_string = String::new();
+
+    if check_closed && !boundary_edges.is_empty() {
+        error_string += &format!("Mesh is not closed. It has {} boundary edges (edges that are connected to only one triangle).", boundary_edges.len());
+        if debug {
+            for e in boundary_edges {
+                add_edge_errors(&mut error_string, e);
+            }
+        }
+        error_string += &format!("\n");
+    }
+
+    if check_manifold && !non_manifold_edges.is_empty() {
+        error_string += &format!("Mesh is not manifold. It has {} non-manifold edges (edges that are connected to more than twi triangles).", non_manifold_edges.len());
+        if debug {
+            for e in non_manifold_edges {
+                add_edge_errors(&mut error_string, e);
+            }
+        }
+        error_string += &format!("\n");
+    }
+
+    if check_manifold && !non_manifold_vertices.is_empty() {
+        error_string += &format!("Mesh is not manifold. It has {} non-manifold vertices (vertices with more than one triangle fan).", non_manifold_vertices.len());
+        if debug {
+            error_string += &format!("\n\t{:?}", non_manifold_vertices);
+        }
+        error_string += &format!("\n");
     }
 
     Err(error_string)
@@ -240,7 +288,13 @@ fn check_mesh_with_cell_data<I: Index, R: Real>(
     marching_cubes_data: &MarchingCubesInput<I>,
     mesh: &TriMesh3d<R>,
 ) -> Result<(), String> {
-    let boundary_edges = mesh.find_boundary_edges();
+    let edge_info = mesh.compute_edge_information();
+
+    let boundary_edges = edge_info
+        .iter()
+        .filter(|ei| ei.incident_faces == 1)
+        .map(|ei| (ei.edge, ei.face, ei.local_edge_index))
+        .collect::<Vec<_>>();
 
     if boundary_edges.is_empty() {
         return Ok(());
@@ -263,7 +317,7 @@ fn check_mesh_with_cell_data<I: Index, R: Real>(
             let cell_data = marching_cubes_data
                 .cell_data
                 .get(&grid.flatten_cell_index(&cell_index))
-                .expect("Unabel to get cell data of cell");
+                .expect("Unable to get cell data of cell");
 
             error_string += &format!("\n\tTriangle {}, boundary edge {:?} is located in cell with {:?} with center coordinates {:?} and edge length {}. {:?}", tri_idx, edge, cell_index, cell_center, grid.cell_size(), cell_data);
         } else {

splashsurf_lib/src/mesh.rs

@@ -987,27 +987,6 @@ impl<R: Real> TriMesh3d<R> {
         }
     }
 
-    /// Returns all boundary edges of the mesh
-    ///
-    /// Returns edges which are only connected to exactly one triangle, along with the connected triangle
-    /// index and the local index of the edge within that triangle.
-    ///
-    /// Note that the output order is not necessarily deterministic due to the internal use of hashmaps.
-    pub fn find_boundary_edges(&self) -> Vec<([usize; 2], usize, usize)> {
-        let MeshEdgeInformation {
-            edge_counts,
-            edge_info,
-        } = self.compute_edge_information();
-
-        // Take only the faces which have a count of 1, which correspond to boundary faces
-        edge_counts
-            .into_iter()
-            .map(|(_edge, value)| value)
-            .filter(|&(_, count)| count == 1)
-            .map(move |(edge_idx, _)| edge_info[edge_idx].clone())
-            .collect()
-    }
-
     /// Returns all non-manifold vertices of this mesh
     ///
     /// A non-manifold vertex is generated by pinching two surface sheets together at that vertex
@@ -1170,23 +1149,6 @@ mod tri_mesh_tests {
         }
     }
 
-    #[test]
-    fn test_tri_mesh_find_boundary() {
-        let mesh = mesh_one_tri();
-
-        let mut boundary = mesh.find_boundary_edges();
-        boundary.sort_unstable();
-
-        assert_eq!(
-            boundary,
-            vec![
-                ([0usize, 1usize], 0, 0),
-                ([1usize, 2usize], 0, 1),
-                ([2usize, 0usize], 0, 2),
-            ]
-        );
-    }
-
     #[test]
     fn test_tri_mesh_edge_info() {
         let mesh = mesh_non_manifold_edge();

splashsurf_lib/tests/integration_tests/test_full.rs

@@ -133,10 +133,11 @@ macro_rules! generate_test {
             let reconstruction =
                 reconstruct_surface::<i64, _>(particle_positions.as_slice(), &parameters).unwrap();
 
-            write_vtk(reconstruction.mesh(), output_file, "mesh").unwrap();
+            write_vtk(reconstruction.mesh(), &output_file, "mesh").unwrap();
 
             println!(
-                "Reconstructed mesh from particle file '{}' with {} particles has {} triangles.",
+                "Reconstructed mesh '{}' from particle file '{}' with {} particles has {} triangles.",
+                output_file.display(),
                 $input_file,
                 particle_positions.len(),
                 reconstruction.mesh().triangles.len()
@@ -158,10 +159,10 @@ macro_rules! generate_test {
 
             if test_for_boundary(&parameters) {
                 // Ensure that the mesh does not have a boundary
-                if let Err(e) = check_mesh_consistency(reconstruction.grid(), reconstruction.mesh())
+                if let Err(e) = check_mesh_consistency(reconstruction.grid(), reconstruction.mesh(), true, true, true)
                 {
                     eprintln!("{}", e);
-                    panic!("Mesh contains boundary edges");
+                    panic!("Mesh contains topological/manifold errors");
                 }
             }
         }