Repair self-intersecting STL surfaces before TetGen meshing

.github/scripts/basic_smoke_test.py

@@ -1,3 +1,4 @@
+import math
 import os
 import sys
 import tempfile
@@ -157,6 +158,7 @@ def main() -> None:
     import pyvista as pv
 
     from svv import Domain, Forest, Simulation, Tree
+    from svv.domain.routines.tetrahedralize import tetrahedralize
     from svv.utils.remeshing.mmg import get_mmg_candidates, get_mmg_exe
     from svv.utils.solvers.solver_0d import get_solver_0d_candidates, get_solver_0d_exe
     from svv.utils.remeshing.remesh import remesh_surface
@@ -185,6 +187,41 @@ def main() -> None:
     cube.solve()
     cube.build()
 
+    # Real geometry-recovery path using a compact, deterministic folded surface.
+    # Future TetGen releases may accept the original directly, so assert portable
+    # mesh invariants rather than a version-specific strategy.
+    _log("SMOKE: tetrahedralize: folded closed surface")
+    folded = pv.Sphere(theta_resolution=20, phi_resolution=20)
+    folded.points[folded.points[:, 2] > 0.2, 2] -= 0.8
+    recovery = tetrahedralize(
+        folded,
+        order=1,
+        nobisect=True,
+        repair_max_distance_ratio=0.4,
+        return_result=True,
+    )
+    assert recovery.report.selected_strategy in {
+        "original",
+        "meshfix",
+        "pyacvd",
+        "pyacvd_meshfix",
+    }
+    assert recovery.elements.shape[0] > 0
+    assert recovery.surface.is_manifold
+    assert recovery.surface.n_open_edges == 0
+    recovery_volumes = recovery.grid.compute_cell_sizes(
+        length=False,
+        area=False,
+        volume=True,
+    ).cell_data["Volume"]
+    assert recovery_volumes.size > 0
+    assert all(
+        math.isfinite(float(value)) and float(value) >= 0
+        for value in recovery_volumes
+    )
+    assert math.isfinite(float(recovery_volumes.sum()))
+    assert float(recovery_volumes.sum()) > 0
+
     # MMG remeshing (validate that packaged/built executables run)
     _log("SMOKE: mmg: remesh_surface(pv.Cube())")
     with _temp_cwd():

.github/workflows/basic-smoke-test.yml

@@ -4,6 +4,9 @@ on:
   push:
     branches:
       - main
+  pull_request:
+    branches:
+      - main
 
 permissions:
   contents: read
@@ -55,7 +58,7 @@ jobs:
         run: |
           python -m pip install --upgrade pip setuptools wheel
           python -m pip install -r requirements.txt
-          python -m pip install cmake
+          python -m pip install cmake pytest
 
       - name: Build MMG (v5.8.0, Release, -O3)
         run: |
@@ -132,6 +135,19 @@ jobs:
           print(f"svzerodsolver -h exit: {proc.returncode}", flush=True)
           PY
 
+      - name: Run domain recovery tests
+        run: |
+          python -m pytest -q \
+            test/test_mesh_diagnostics.py \
+            test/test_tetgen_worker.py \
+            test/test_tetrahedralize_recovery.py \
+            test/test_domain_tetrahedralization_recovery.py \
+            test/test_gui_domain_build_feedback.py
+        env:
+          QT_QPA_PLATFORM: "offscreen"
+          SVV_GUI_DISABLE_VTK: "1"
+          SVV_TELEMETRY_DISABLED: "1"
+
       - name: Run basic smoke test
         run: |
           if [ "${{ runner.os }}" = "Linux" ]; then

docs/api/domain.html

@@ -177,13 +177,23 @@ <h2>Attributes</h2>
                             <tr>
                                 <td><code>boundary</code></td>
                                 <td>pv.PolyData</td>
-                                <td>Boundary mesh of the domain</td>
+                                <td>The exact boundary surface used to create the current interior mesh</td>
+                            </tr>
+                            <tr>
+                                <td><code>original_boundary</code></td>
+                                <td>pv.PolyData</td>
+                                <td>Unmodified imported surface retained as source provenance</td>
                             </tr>
                             <tr>
                                 <td><code>mesh</code></td>
                                 <td>pv.UnstructuredGrid</td>
                                 <td>Interior tetrahedral/triangular mesh</td>
                             </tr>
+                            <tr>
+                                <td><code>mesh_build_report</code></td>
+                                <td>TetrahedralizationReport or None</td>
+                                <td>Ordered strategies, surface checks, versions, and bounded TetGen diagnostics from the latest 3D mesh build</td>
+                            </tr>
                             <tr>
                                 <td><code>patches</code></td>
                                 <td>list</td>
@@ -282,14 +292,15 @@ <h4>Parameters</h4>
 
                         <div class="api-method">
                             <div class="api-method-signature">
-                                <code>build(resolution=25, skip_boundary=False)</code>
+                                <code>build(resolution=25, skip_boundary=False, **interior_kwargs)</code>
                             </div>
                             <p>Build the implicit function describing the domain and optionally extract boundary/mesh artifacts.</p>
                             <div class="api-method-params">
                                 <h4>Parameters</h4>
                                 <ul>
                                     <li><code>resolution</code> (int, default=25): Grid resolution for boundary extraction.</li>
                                     <li><code>skip_boundary</code> (bool, default=False): If true, only assemble fast-evaluation structures and skip boundary and interior mesh generation.</li>
+                                    <li><code>**interior_kwargs</code>: Options forwarded to <code>get_interior()</code>, including the bounded surface-recovery controls below.</li>
                                 </ul>
                             </div>
                         </div>
@@ -346,7 +357,7 @@ <h3>Mesh Generation</h3>
                                 <h4>Parameters</h4>
                                 <ul>
                                     <li><code>resolution</code> (int): Grid resolution for marching cubes/squares</li>
-                                    <li><code>get_largest</code> (bool, default=True): Reserved (current implementation always returns the largest connected component)</li>
+                                    <li><code>get_largest</code> (bool, default=True): Keep only the largest connected component during initial boundary extraction. Set to false to retain all components.</li>
                                 </ul>
                                 <h4>Returns</h4>
                                 <ul>
@@ -360,17 +371,70 @@ <h4>Returns</h4>
                             <div class="api-method-signature">
                                 <code>get_interior(verbose=False, **kwargs)</code>
                             </div>
-                            <p>Generate tetrahedral (3D) or triangular (2D) mesh of the interior.</p>
+                            <p>Generate a tetrahedral (3D) or triangular (2D) interior mesh. In 3D, geometry-related TetGen failures use a validated, component-preserving recovery pipeline.</p>
                             <div class="api-method-params">
                                 <h4>Parameters</h4>
                                 <ul>
                                     <li><code>verbose</code> (bool): Print mesh generation progress</li>
                                     <li><code>**kwargs</code>: Parameters passed to TetGen/Triangle</li>
+                                    <li>Raw TetGen <code>switches</code> strings are rejected in 3D because Domain enforces linear <code>order=1</code> elements and <code>nobisect=True</code>. Use named options instead.</li>
+                                    <li><code>repair_on_failure</code> (bool, default=True): Retry a geometry rejection after component-preserving MeshFix repair.</li>
+                                    <li><code>repair_max_distance_ratio</code> (float, default=0.01): Maximum symmetric repair displacement and bounds change, measured against the source bounding-box diagonal.</li>
+                                    <li><code>remesh_on_failure</code> (bool, default=True): Enable validated PyACVD as the final fallback after the original and MeshFix candidates.</li>
+                                    <li><code>remesh_subdivisions</code>, <code>remesh_clusters</code>, <code>remesh_clean_tolerance</code>: Configure the optional PyACVD stage.</li>
                                 </ul>
                                 <h4>Returns</h4>
                                 <ul>
                                     <li><code>mesh</code> (pv.UnstructuredGrid): Interior mesh</li>
                                 </ul>
+                                <h4>Recovery and provenance</h4>
+                                <p>
+                                    The original surface is tried first with the caller's TetGen options. A valid
+                                    first attempt returns immediately. After a geometry rejection, MeshFix runs with
+                                    component joining and component removal disabled. Candidates must remain closed,
+                                    manifold, finite, triangular, component-preserving, and within the configured
+                                    displacement bound. PyACVD output is checked by the same policy and is repaired
+                                    before TetGen when necessary; unsafe candidates are recorded and rejected.
+                                </p>
+                                <p>
+                                    Closed, watertight, and manifold are topological properties; they do not prove
+                                    that a surface has no self-intersecting facets. TetGen may therefore reject a
+                                    surface that passes those checks. The structured report records intersection
+                                    diagnostics from a separate verbose pass when the initial failure is opaque.
+                                </p>
+                                <p>
+                                    When recovery succeeds, <code>boundary</code>, boundary sampling arrays, and the
+                                    extracted volume-mesh surface all describe the selected recovery geometry.
+                                    <code>original_boundary</code> remains unchanged. Inspect
+                                    <code>mesh_build_report.selected_strategy</code> and
+                                    <code>mesh_build_report.attempts</code> for provenance and diagnostics.
+                                </p>
+                            </div>
+                        </div>
+
+                        <div class="api-method">
+                            <div class="api-method-signature">
+                                <code>tetrahedralize(surface, *tet_args, return_result=False, **tet_kwargs)</code>
+                            </div>
+                            <p>
+                                Lower-level public surface meshing helper. Unlike <code>Domain.get_interior()</code>,
+                                this compatibility API forwards named TetGen options or a raw
+                                <code>switches</code> string and can return linear or quadratic tetrahedra.
+                            </p>
+                            <div class="api-method-params">
+                                <h4>Returns</h4>
+                                <ul>
+                                    <li>Default: the historical <code>(grid, nodes, elements)</code> tuple.</li>
+                                    <li><code>return_result=True</code>: <code>TetrahedralizationResult</code> with <code>grid</code>, <code>nodes</code>, <code>elements</code>, the exact selected <code>surface</code>, and a structured <code>report</code>.</li>
+                                </ul>
+                                <h4>Recovery controls</h4>
+                                <p>
+                                    The helper tries the original surface, bounded MeshFix repair, and validated
+                                    PyACVD recovery in order. Set <code>repair_on_failure=False</code> or
+                                    <code>remesh_on_failure=False</code> to disable a stage. Raw switches are unchanged
+                                    for real meshing attempts; an isolated diagnostic-only rerun removes quiet mode
+                                    and adds TetGen diagnostic verbosity when needed.
+                                </p>
                             </div>
                         </div>
                     </div>

docs/gui.html

@@ -43,6 +43,7 @@ <h3 id="on-this-page">On this page</h3>
       <li><a href="#capabilities">What the GUI can do</a></li>
       <li><a href="#navigation">How to navigate the interface</a></li>
       <li><a href="#workflow">A typical workflow</a></li>
+      <li><a href="#domain-meshing">Domain meshing and recovery</a></li>
     </ul>
 
     <h2 id="capabilities">What The GUI Can Do</h2>
@@ -187,6 +188,39 @@ <h2 id="workflow">Typical Workflow</h2>
       <p style="margin:0">If point-picking does not work, ensure the domain was created/solved/built before loading.</p>
     </div>
 
+    <h2 id="domain-meshing">Domain Meshing And Recovery</h2>
+    <p>
+      Loading a 3D mesh runs the domain create, solve, and build stages. The imported surface is sent to TetGen first.
+      If TetGen detects a geometry problem, the GUI automatically tries a component-preserving MeshFix repair and then
+      a validated PyACVD fallback. A recovered surface is accepted only when it is finite, triangular, closed,
+      manifold, preserves connected components, and stays within the configured geometric-fidelity bound.
+    </p>
+    <p>
+      A surface can be closed, watertight, and manifold while still containing self-intersecting facets. Those
+      topology checks alone therefore do not guarantee that TetGen can construct a volume mesh. When this occurs,
+      the expanded diagnostics identify representative intersecting facets or segments when TetGen reports them.
+    </p>
+    <p>
+      A successful recovery is reported in the status bar as a build completed after surface repair; it does not open
+      a blocking success dialog. The working domain boundary is updated to the exact surface that generated the volume
+      mesh, while the imported surface remains available as <code>original_boundary</code>.
+    </p>
+    <div class="callout warning">
+      <strong>If interior meshing still fails:</strong>
+      <ul>
+        <li>Read the short cause and suggested action in the warning dialog.</li>
+        <li>Expand <em>Show Details</em> to see each attempted strategy, topology counts, TetGen version/options, return code or signal, and representative facet or segment identifiers.</li>
+        <li>For intersecting facets, repair self-intersections in the source model and retry.</li>
+        <li>For an open or non-manifold candidate, close holes and correct non-manifold edges without deleting required components.</li>
+        <li>For a worker or dependency failure, verify the TetGen, PyMeshFix, PyACVD, PyVista, and Python environment.</li>
+      </ul>
+      <p style="margin-bottom:0">
+        Python callers can disable individual recovery stages or tighten the default 1% displacement limit with
+        <code>repair_on_failure</code>, <code>remesh_on_failure</code>, and
+        <code>repair_max_distance_ratio</code>. The GUI retains the structured build report for troubleshooting.
+      </p>
+    </div>
+
     <div class="callout warning">
       <strong>NumPy compatibility:</strong>
       <p style="margin:.5rem 0 0">If importing the GUI fails with <code>_ARRAY_API not found</code> or <code>numpy.core.multiarray failed to import</code>, recreate the environment with a NumPy version supported by your Python version. On Python 3.9–3.12, use <code>python -m pip install --force-reinstall "numpy&lt;2" svv</code>. Python 3.13 requires NumPy 2.1 or newer; use <code>python -m pip install --force-reinstall "numpy&gt;=2.1" svv</code>.</p>

pyproject.toml

@@ -10,3 +10,6 @@ requires = [
 ]
 
 build-backend = "setuptools.build_meta"
+
+[tool.pytest.ini_options]
+testpaths = ["test"]