Nmesh section 4 - relaxation mesher internals

PARITY_COMPLETION.md

@@ -0,0 +1,53 @@
+# NMesh Section 4 Parity Completion
+
+**Status:** exact standalone behavioral parity required  
+**Branch:** `nmesh-section4`  
+**Reviewed:** 2026-05-04  
+**Reference sources:** `nmag-src/src/mesh.ml`, `nmag-src/src/snippets.ml`, `nmag-src/interface/nmeshlib/lib1.py`
+
+## Position
+
+The Python 3 `nmesh` port must implement the full behavior of the original OCaml/Python2 mesher without requiring OCaml at runtime. Missing behavior is implementation work, not an accepted end state.
+
+The target is exact behavioral parity, not approximate parity. Topology decisions, region assignment, point-state transitions, add/delete scheduling, boundary recovery, periodic grouping, and file-visible mesh metadata must match the legacy implementation for the same inputs. Numeric tolerances are only allowed for floating-point coordinate comparisons, and each tolerance must be tied to a specific floating-point or triangulator tie-break reason.
+
+The port may use Python, NumPy, SciPy, Numba, Rust, C++, or C as needed. The public Python API should remain stable, and compiled helpers should be isolated behind that API.
+
+## Completed In This Section 4 Pass
+
+- Removed the Python-side 60-step cap; configured `controller_step_limit_max` is honored.
+- Implemented legacy controller cadence: square-number add/delete checks, relaxed add/delete thresholds, topology-threshold retriangulation, minimum-step convergence, force-equilibrium stopping, and 50-step post-change settling.
+- Aligned force metrics with legacy behavior: density-scaled movement, centroid density for simplex forces, tangent-projected boundary effective force, corner suppression, and OCaml-style Voronoi density correction.
+- Classified force-time simplices with the same centroid, near-vertex probe, `Boundary` state, and raw volume-order ratio rules used by the legacy relevant/irrelevant topology split.
+- Restricted irrelevant-element forces to mobile nodes, matching the active legacy path.
+- Replaced deterministic grid seeding with the legacy random density estimator, D-lattice sphere-packing node estimate, rejection sampler, and RNG consumption order.
+- Replaced midpoint point insertion with Gaussian insertion around the source point using local effective rod length.
+- Ported active boundary recovery behavior from `mirror_simplices`: 2D mirrored points and 3D boundary-edge midpoint prevention points.
+- Seeded paired fixed points on periodic outer-box faces, fixed only exact periodic boundary nodes, and canonicalized multi-axis periodic equivalence groups.
+- Added final high-density moving-point cleanup, outside dynamic point cleanup, and final boundary snapping before final assembly.
+- Normalized final simplex orientation to the legacy positive-volume convention.
+- Fixed flat-boundary filtering so it follows the legacy `Boundary` state and raw volume-order ratio rules instead of deleting all geometrically boundary-adjacent multi-region simplices.
+- Changed unsupported density snippets to fail loudly rather than silently falling back to density `1.0`.
+- Added parity coverage for coarse adjacent pieces, concave difference domains, multi-axis periodic groups, final cleanup, and force-time relevant/irrelevant classification.
+- Added canonical mesh signature helpers so parity fixtures compare topology and metadata exactly, with coordinate tolerance explicitly opt-in.
+- Added a parity comparison runner (`tools/nmesh_parity_compare.py`) that generates modern scenario meshes and compares them with legacy `.nmesh` artifacts or a configured legacy runner command.
+- Tightened runtime matching to avoid mesh-size-scaled boundary drift tolerance and rounded non-periodic periodic-group coordinates.
+
+## Validation Gate
+
+These are verification tasks, not accepted limitations:
+
+- Build reference fixtures from legacy examples and compare canonicalized topology, regions, surfaces, links, periodic groups, controller decisions, and mesh metadata exactly using `nmesh.mesher.parity`.
+- Compare coordinates with the narrowest practical floating-point tolerance, documenting each tolerance and why exact binary equality is not the right assertion.
+- Run the modernization scenario matrix through `tools/nmesh_parity_compare.py` and record metric-level acceptance thresholds for physically equivalent meshes when exact topology is not expected.
+- Complete reference validation for periodic outer-box workflows on multiple periodic directions and complex periodic entities.
+- Expand multi-piece/interface parity fixtures to include more interface-prevention cases.
+- Profile large 3D meshes with legacy defaults and move any proven hot path to a compiled Rust, C++, or C helper behind the Python API.
+
+## Sign-Off Criteria
+
+- Legacy reference examples produce exact canonical mesh topology and metadata parity.
+- Periodic, hint-driven, concave, and multi-piece examples match legacy behavior exactly, apart from explicitly documented floating-point coordinate tolerance.
+- Add/delete, retriangulation, boundary recovery, and final cleanup decisions are covered by parity tests.
+- Unsupported density syntax never changes mesh density silently.
+- The full project test suite and parity fixtures pass in the project venv.

README.md

@@ -2,6 +2,8 @@
 
 This is an in-progress port of [nmag](https://github.com/nmag-project/nmag-src) that updates old dependencies and packages nmag as a standalone python 3 module.
 
+The Python 3 port is intended to be standalone and must not require OCaml at runtime. See [PARITY_COMPLETION.md](PARITY_COMPLETION.md) for the Section 4 parity record and validation gate.
+
 ### Installation
 * Create a virtual environment: `python -m venv venv`
 * Activate the virtual environment:
@@ -12,13 +14,7 @@ This is an in-progress port of [nmag](https://github.com/nmag-project/nmag-src)
 ### Testing
 * To install the optional test dependencies run `pip install -e ".[test]"`
 * Once the project has been initialized simply run `pytest` to run your tests
+* To generate modern nmesh parity scenario outputs run `python tools/nmesh_parity_compare.py --output-dir parity/new`
+* To compare against legacy `.nmesh` artifacts run `python tools/nmesh_parity_compare.py --legacy-dir parity/legacy`
 * To view coverage open `htmlcov/index.html` in your browser
 * [Specifying what tests to run](https://docs.pytest.org/en/latest/how-to/usage.html#specifying-which-tests-to-run)
-
-### OCaml
-* Install latest Ocaml/Opam version 4 https://ocaml.org/install
-  * Recommended to run `opam init` and setup a switch with version 4
-* Install ocaml-in-python `opam install ocaml-in-python`
-* Register the package in python ``pip install --editable "`opam var ocaml-in-python:lib`"``
-* Tell python where to look for the OCaml library `export OCAMLPATH=${DUNE_DIR}/_build/install/default/lib` where DUNE_DIR is the directory of the dune project
-* Run this to explictly activate the Opam switch `eval $(opam env)`

pyproject.toml

@@ -12,7 +12,7 @@ requires-python = ">=3.10"
 classifiers = [
   "Programming Language :: Python :: 3",
   "Operating System :: OS Independent" ]
-dependencies = ["pint", "numpy", "tabulate", "meshio", "h5py"]
+dependencies = ["pint", "numpy", "scipy", "numba", "tabulate", "meshio", "h5py"]
 
 [project.optional-dependencies]
 test = ["pytest", "pytest-cov", "pytest-watch"]

src/nmesh/backend.py

@@ -167,7 +167,23 @@ def mesh_bodies_raw(
         cache,
         hints,
     ):
-        return RawMesh(dim=len(bb_min))
+        from .mesher.relaxation import mesh_bodies_raw as python_mesh_bodies_raw
+
+        return python_mesh_bodies_raw(
+            driver,
+            mesher,
+            bb_min,
+            bb_max,
+            mesh_ext,
+            objects,
+            a0,
+            density,
+            fixed,
+            mobile,
+            simply,
+            periodic,
+            hints,
+        )
 
     def mesh_from_points_and_simplices(
         self,

src/nmesh/io/meshio_support.py

@@ -7,12 +7,7 @@
 
 from ..backend import RawMesh
 from .legacy_nmesh_hdf5 import load_raw_mesh_from_legacy_nmesh_hdf5
-
-try:
-    from meshio._exceptions import ReadError as MeshioReadError
-except ImportError:
-    # Fallback for older meshio versions
-    MeshioReadError = Exception
+from meshio._exceptions import ReadError as MeshioReadError
 
 
 _CELL_TYPE_BY_DIM = {

src/nmesh/mesher/__init__.py

@@ -1,2 +1,4 @@
 from .meshing_parameters import *
 from .driver import *
+from .parity import *
+from .relaxation import *

src/nmesh/mesher/meshing_parameters.py

@@ -7,7 +7,7 @@
 
 from mock_features import MockFeatures
 
-from ..utils.constants import MIN_DIVISION_MAGNITUDE
+from ..utils.constants import BOUNDARY_FUZZ, MIN_DIVISION_MAGNITUDE
 
 log = logging.getLogger(__name__)
 
@@ -248,7 +248,7 @@ def _setup_defaults(self):
             # Volume determination
             "nr_probes_for_determining_volume": 100000,
             # Boundary condition parameters
-            "boundary_condition_acceptable_fuzz": 1e-6,
+            "boundary_condition_acceptable_fuzz": BOUNDARY_FUZZ,
             "boundary_condition_max_nr_correction_steps": 200,
             "boundary_condition_debuglevel": 0,
             # Relaxation parameters
@@ -341,6 +341,9 @@ def to_mesher_config(self, dim):
                 continue
             resolved[spec.internal_name] = spec.cast(value)
 
+        for key, value in self._params.items():
+            resolved.setdefault(key, value)
+
         return resolved
 
     def apply_to_mesher(self, mesher, dim):
@@ -354,6 +357,9 @@ def apply_to_mesher(self, mesher, dim):
                 continue
             mesher["parameters"][spec.internal_name] = spec.cast(value)
 
+        for key, value in self._params.items():
+            mesher["parameters"].setdefault(key, value)
+
         return mesher
 
     def _set_parameter(self, name, value):