diff --git a/.github/workflows/codecov-upload.yml b/.github/workflows/codecov-upload.yml index 45d64d811..91cd69cbb 100644 --- a/.github/workflows/codecov-upload.yml +++ b/.github/workflows/codecov-upload.yml @@ -56,7 +56,9 @@ jobs: - name: Test working-directory: build - run: ctest --no-tests=error --output-on-failure -j 2 + # Coverage instrumentation makes the monolithic unit-test binary + # substantially slower than it is in normal CI. + run: ctest --no-tests=error --output-on-failure --timeout 600 -j 2 - name: Submit to codecov.io uses: codecov/codecov-action@fb8b3582c8e4def4969c97caa2f19720cb33a72f # v7.0.0 diff --git a/Justfile b/Justfile index aefe49247..27eb2fbf7 100644 --- a/Justfile +++ b/Justfile @@ -5,7 +5,7 @@ set minimum-version := "1.57.0" set shell := ["bash", "-euo", "pipefail", "-c"] just_version := "1.57.0" -uv_version := "0.11.29" +uv_version := "0.11.30" pinact_version := "4.1.0" pinact_module := "github.com/suzuki-shunsuke/pinact/v4/cmd/pinact@v" + pinact_version llvm_version := "22" diff --git a/README.md b/README.md index 6e55efc67..91264b014 100644 --- a/README.md +++ b/README.md @@ -370,6 +370,9 @@ Online documentation is at . The scientific transition, proposal-ratio, geometry-delta, counter, and precision contracts are recorded in [`docs/metropolis-hastings.md`](docs/metropolis-hastings.md). +The literature-backed contracts, exact deltas, inverse relationships, and +failure-atomicity rules for the complete 2+1D move set are recorded in +[`docs/ergodic-moves.md`](docs/ergodic-moves.md). Seed replay, PCG stream ownership, checkpoint metadata, and the parallel stream policy are recorded in [`docs/reproducibility.md`](docs/reproducibility.md). The repository-wide scientific bibliography is maintained in diff --git a/docs/ergodic-moves.md b/docs/ergodic-moves.md new file mode 100644 index 000000000..2d6917059 --- /dev/null +++ b/docs/ergodic-moves.md @@ -0,0 +1,103 @@ +# 2+1D CDT ergodic move audit + +This note records the scientific and implementation contract for the five +supported 2+1-dimensional causal dynamical triangulation (CDT) moves. It is the +audit record for issue #106 and distinguishes three separate claims: + +1. the local replacement is one of the Lorentzian Monte Carlo moves; +2. CGAL preserves the validity of the combinatorial triangulation data + structure; and +3. CDT++ independently verifies the stronger causal, foliation, metadata, and + failure-atomicity invariants. + +The primary CDT source is Ambjørn, Jurkiewicz, and Loll, *Dynamically +triangulating Lorentzian quantum gravity*, section 7.1, equations (56)-(58), +with the corresponding f-vector changes in equations (71)-(73). See +@ref cdt-framework-2001 "CDT framework (2001)". The topological background is +@ref pachner-moves "Pachner's theorem", and the implementation primitive +contracts come from the @ref cgal-triangulations "CGAL 6.2 manual". + +## Conventions and common contract + +A `(p,q)` tetrahedron has `p` vertices on slice `t` and `q` vertices on slice +`t+1`. Every accepted finite tetrahedron must therefore have exactly two +adjacent time values and type `(3,1)`, `(2,2)`, or `(1,3)`. A spacelike +subsimplex has equal endpoint times; a timelike one spans adjacent slices. + +Every high-level `do_*_move()` and `propose_*_move()` operation mutates an +owning snapshot. An inapplicable site, a checked CGAL rejection, or a failed +postcondition leaves the source manifold's canonical topology, vertex and cell +metadata, geometry counts, time bounds, foliation parameters, and handle caches +unchanged. A successful move constructs a new `FoliatedTriangulation`, which +rebuilds all handle-bearing caches and scalar geometry from the moved canonical +triangulation. + +CGAL 6.2 documents that checked three-dimensional flips preserve vertex handles +and invalidate only affected cell handles. CDT++ never carries affected cell +handles past a successful flip. A copied triangulation has different handles; +the `(4,4)` and `(6,2)` implementations therefore capture stable point values +and re-resolve vertices and edges in their private copies before mutation. + +None of these moves is required to preserve the Euclidean empty-sphere +(Delaunay) property of the representative coordinates. The scientific state is +a valid causal combinatorial triangulation. `tds().is_valid()` establishes the +CGAL incidence and adjacency contract; `Manifold::is_correct()` adds foliation +and causal cell-metadata checks. + +## Move-by-move record + +| Move | Local cavity and time assignment | Independent delta `(N0, N1_SL, N1_TL, N2, N3_31, N3_22, N3_13, N3)` | Implementation and admissibility | +| --- | --- | --- | --- | +| `(2,3)` | Equation (58): `(3,1) 1345 + (2,2) 2345 -> (3,1) 1234 + (2,2) 1235 + (2,2) 1245`, or its time reflection. The shared triangle `345` is timelike and the new edge `12` is timelike. | `(0, 0, +1, +2, 0, +1, 0, +1)` | A checked `Triangulation_3::flip(facet)` is attempted only from correctly labelled `(2,2)` cells whose neighbor is a correctly labelled `(3,1)` or `(1,3)` cell and whose opposite vertices lie on adjacent slices. CGAL additionally rejects infinite, nonflippable, or geometrically inverted cavities. | +| `(3,2)` | The inverse of equation (58): two `(2,2)` cells and exactly one `(3,1)` or `(1,3)` cell meet at the timelike edge removed by the move. | `(0, 0, -1, -2, 0, -1, 0, -1)` | Before `Triangulation_3::flip(edge)`, CDT++ independently requires three finite incident cells spanning adjacent slices with the exact `2 x (2,2) + 1 x ((3,1) or (1,3))` causal composition. CGAL then rejects hull edges and enforces the total degree-three and geometric flippability contracts. Other degree-three timelike cavities can be combinatorially flippable but causally invalid. | +| `(2,6)` | Equation (56): `(1,3) 1345 + (3,1) 2345` share the spacelike triangle `345`. A new vertex `6` is inserted on the same slice and joined to all five old vertices, producing three tetrahedra above and three below. | `(+1, +3, +2, +8, +2, 0, +2, +4)` | `tds().insert_in_facet()` subdivides the common spacelike facet. The input cells and their metadata must be `(1,3)` and `(3,1)`, the shared vertices must have one time value, and the new star must contain six valid cells. The new point is the facet centroid and receives the facet time. | +| `(6,2)` | The inverse of equation (56): a degree-five vertex has six incident cells, exactly three `(3,1)` and three `(1,3)`, with no `(2,2)` cell. | `(-1, -3, -2, -8, -2, 0, -2, -4)` | On a private copy, a checked timelike edge flip reduces the candidate to degree four; `tds().remove_from_maximal_dimension_simplex()` then applies its documented degree-`dimension+1` removal. Exact finite incidence, causal types, metadata, output counts, and output cell types are checked before publication. | +| `(4,4)` | Equation (57): two `(1,3)` and two `(3,1)` tetrahedra form a diamond. The diagonal of the spatial quadrilateral is exchanged; the move is its own inverse. | `(0, 0, 0, 0, 0, 0, 0, 0)` | The pivot must be a spacelike edge with exactly four finite incident cells, two `(3,1)` and two `(1,3)`, all correctly labelled. A checked TDS facet flip creates the new diagonal, then a checked TDS edge flip removes the old one. The composition is performed on a private copy because the transient geometry can fail a `Triangulation_3` geometric flip even when the final abstract diamond is valid. | + +## Independent delta derivations + +For `(2,3)`, two tetrahedra sharing one face initially have five vertices, nine +edges, seven faces, and two cells. Replacing the shared face by the dual edge +adds that one timelike edge. The three resulting tetrahedra have nine distinct +faces, so `(N1, N2, N3)` changes by `(+1, +2, +1)`. Equation (58) changes one +`(3,1)` plus one `(2,2)` into one `(3,1)` plus two `(2,2)`. + +For `(2,6)`, the inserted vertex is joined to all five old vertices. Its three +edges to the shared spatial triangle are spacelike and its two edges to the +opposite vertices are timelike. Thus `N0` changes by `+1` and `N1` by `+5`. +The local complex changes from two to six tetrahedra. With the boundary fixed, +the Euler relation for the three-ball gives +`delta N2 = delta N1 + delta N3 - delta N0 = 5 + 4 - 1 = 8`. Equation (56) +changes one tetrahedron of each time orientation into three of each, giving +`delta N3_31 = delta N3_13 = +2`. The `(6,2)` changes are the negatives. + +For `(4,4)`, the six vertices and the boundary of the four-cell diamond are +unchanged. One spacelike diagonal replaces another, four tetrahedra replace +four tetrahedra, and the two tetrahedra of each time orientation remain two. +The fixed boundary and Euler relation then force every tracked count delta to +zero. + +These derivations agree with the source's f-vectors (71)-(73) but do not call +`ergodic_moves::check_move()` or any production count-delta helper. + +## Verification record + +`tests/Ergodic_moves_3_audit_test.cpp` supplies: + +- deterministic minimal successful fixtures for all five moves; +- direct enumeration of finite vertices, edges, facets, cells, edge types, and + tetrahedron types from vertex time labels; +- reciprocal adjacency, Euler, TDS, foliation, time-bound, metadata, and cache + checks after every successful move; +- canonical point/time/simplex representations independent of CGAL handle + identity; +- exact inverse round trips for `(2,3)/(3,2)`, `(2,6)/(6,2)`, and `(4,4)`; +- malformed-handle, stale-metadata, wrong-cavity, non-applicable, and empty-state + rejection checks with canonical failure-atomicity comparisons; and +- replay of every inverse pair with seeds `0`, `1`, `2`, `92`, `106`, and + `20260721`. + +The wrong-cavity `(3,2)` regression deliberately constructs a geometrically +flippable CGAL edge whose incident CDT types are invalid. It verifies that raw +CGAL accepts the flip while CDT++ rejects before mutation, directly testing the +boundary between combinatorial and scientific validity. diff --git a/include/Ergodic_moves_3.hpp b/include/Ergodic_moves_3.hpp index 771057c99..94cadb469 100644 --- a/include/Ergodic_moves_3.hpp +++ b/include/Ergodic_moves_3.hpp @@ -14,12 +14,15 @@ /// @see [Pachner moves](../REFERENCES.md#pachner-moves) /// @see [Three-dimensional CDT move /// set](../REFERENCES.md#three-dimensional-cdt-2001) +/// @see [2+1D CDT ergodic move audit](../docs/ergodic-moves.md) #ifndef CDT_PLUSPLUS_ERGODIC_MOVES_3_HPP #define CDT_PLUSPLUS_ERGODIC_MOVES_3_HPP -#include - +#include +#include +#include +#include #include #include @@ -40,6 +43,22 @@ namespace ergodic_moves namespace detail { + /// @brief Compare preserved floating-point configuration state exactly. + /// @details Move construction copies these values; arithmetic tolerance is + /// inappropriate because any representation change indicates state drift. + [[nodiscard]] inline auto same_configuration_value( + double const first, double const second) noexcept -> bool + { + using Representation = std::array; + return std::bit_cast(first) == + std::bit_cast(second); + } + + /// @brief Default validator for internal post-mutation test seams. + [[nodiscard]] inline auto accept_post_mutation( + [[maybe_unused]] Delaunay const& triangulation) noexcept -> bool + { return true; } + /// @brief Rebuild all derived topology and geometry state around a value. [[nodiscard]] inline auto make_manifold(Delaunay triangulation, Manifold const& source) -> Manifold @@ -65,6 +84,30 @@ namespace ergodic_moves valid_index(edge.third) && edge.second != edge.third; } + /// @brief Collect the finite cells incident to a checked edge. + /// @details The returned handles borrow from @p triangulation and remain + /// valid only until an affected-cell mutation is performed. + [[nodiscard]] inline auto finite_incident_cells( + Delaunay const& triangulation, Edge_handle const& edge) + -> std::optional + { + if (!is_well_formed_edge(edge) || triangulation.dimension() != 3 || + !triangulation.tds().is_edge(edge.first, edge.second, edge.third)) + { + return std::nullopt; + } + + auto circulator = triangulation.incident_cells(edge, edge.first); + Cell_container incident_cells; + do + { // NOLINT(cppcoreguidelines-avoid-do-while) + if (triangulation.is_infinite(circulator)) { return std::nullopt; } + incident_cells.emplace_back(circulator); + } + while (++circulator != edge.first); + return incident_cells; + } + /// Select exactly one raw proposal site uniformly. template [[nodiscard]] inline auto random_element(Container const& candidates, @@ -85,7 +128,9 @@ namespace ergodic_moves [[nodiscard]] inline auto null_move(Manifold const& t_manifold) -> Expected { return t_manifold; } // null_move - /// @brief Perform a TriangulationDataStructure_3::flip on a facet + /// @brief Perform a checked Triangulation_3 flip on a causal CDT facet + /// @details On success, affected cell handles are invalidated; vertex handles + /// remain valid. Rejection occurs before mutation. /// @param triangulation The triangulation containing the cell to flip /// @param to_be_moved The cell on which to try the move /// @returns True if move succeeded @@ -94,7 +139,14 @@ namespace ergodic_moves [[nodiscard]] inline auto try_23_move(Delaunay& triangulation, Cell_handle const& to_be_moved) -> bool { - if (to_be_moved->info() != 22) { return false; } // NOLINT + if (to_be_moved == nullptr || triangulation.dimension() != 3 || + !triangulation.tds().is_cell(to_be_moved) || + !foliated_triangulations::is_cell_type_correct<3>(to_be_moved) || + foliated_triangulations::expected_cell_type<3>(to_be_moved) != + Cell_type::TWO_TWO) + { + return false; + } auto flipped = false; // Try every facet of the (2,2) cell for (auto i = 0; i < 4; ++i) @@ -102,6 +154,15 @@ namespace ergodic_moves auto const neighbor = to_be_moved->neighbor(i); if (triangulation.is_infinite(neighbor)) { continue; } + auto const neighbor_type = + foliated_triangulations::expected_cell_type<3>(neighbor); + if (!foliated_triangulations::is_cell_type_correct<3>(neighbor) || + (neighbor_type != Cell_type::THREE_ONE && + neighbor_type != Cell_type::ONE_THREE)) + { + continue; + } + // A causal (2,3) move must replace the facet with a timelike edge. // CGAL also permits topological flips that create a spacelike edge. auto const mirror_index = neighbor->index(to_be_moved); @@ -131,13 +192,15 @@ namespace ergodic_moves /// This function calls try_23_move on (2,2) simplices drawn from a /// randomly shuffled container until it succeeds or runs out of simplices. /// - /// If successful, the triangulation is no longer Delaunay. + /// The move guarantees a valid causal combinatorial triangulation. It does + /// not preserve or require the empty-sphere property of the coordinates. /// /// @tparam Generator A uniform random bit generator type /// @param t_manifold The simplicial manifold /// @param generator Caller-owned generator whose state advances during the /// move /// @returns The Expected (2,3) moved manifold or an Unexpected + /// @note The source manifold is unchanged on success and failure. template [[nodiscard]] inline auto do_23_move(Manifold const& t_manifold, Generator& generator) -> Expected @@ -180,15 +243,58 @@ namespace ergodic_moves return std::unexpected("Selected (2,3) proposal site is not movable.\n"); } - /// @brief Perform a TriangulationDataStructure_3::flip on an edge + namespace detail + { + /// @brief Check for the causal cavity inverse to a (2,3) move. + [[nodiscard]] inline auto is_32_movable(Delaunay const& triangulation, + Edge_handle const& candidate) + -> bool + { + auto const incident_cells = + finite_incident_cells(triangulation, candidate); + if (!incident_cells || incident_cells->size() != 3) { return false; } + + auto const first_time = static_cast( + candidate.first->vertex(candidate.second)->info()); + auto const second_time = static_cast( + candidate.first->vertex(candidate.third)->info()); + auto const time_difference = first_time > second_time + ? first_time - second_time + : second_time - first_time; + if (time_difference != 1) { return false; } + + if (!std::ranges::all_of(*incident_cells, [](auto const& cell) { + return foliated_triangulations::is_cell_type_correct<3>(cell); + })) + { + return false; + } + + auto const cell_type_count = [&](Cell_type const type) { + return std::ranges::count_if(*incident_cells, [&](auto const& cell) { + return foliated_triangulations::expected_cell_type<3>(cell) == type; + }); + }; + auto const incident_31 = cell_type_count(Cell_type::THREE_ONE); + auto const incident_22 = cell_type_count(Cell_type::TWO_TWO); + auto const incident_13 = cell_type_count(Cell_type::ONE_THREE); + return incident_22 == 2 && ((incident_31 == 1 && incident_13 == 0) || + (incident_31 == 0 && incident_13 == 1)); + } + } // namespace detail + + /// @brief Perform a checked (3,2) flip on a causal CDT edge cavity. + /// @details On success, affected cell handles are invalidated; vertex handles + /// remain valid. Rejection occurs before mutation. /// @param triangulation The triangulation containing the edge to flip /// @param to_be_moved The edge on which to try the move - /// @returns True if move succeeded + /// @returns True if the CDT cavity is admissible and CGAL performs the flip /// @see /// https://doc.cgal.org/latest/TDS_3/classTriangulationDataStructure__3.html#a5837d666e4198f707f862003c1ffa033 [[nodiscard]] inline auto try_32_move(Delaunay& triangulation, Edge_handle const& to_be_moved) -> bool { + if (!detail::is_32_movable(triangulation, to_be_moved)) { return false; } return triangulation.flip(to_be_moved.first, to_be_moved.second, to_be_moved.third); } // try_32_move @@ -198,12 +304,14 @@ namespace ergodic_moves /// This removes a (2,2) simplex and the timelike edge. /// This function calls try_32_move on timelike edges drawn from a /// randomly shuffled container until it succeeds or runs out of edges. - /// If successful, the triangulation is no longer Delaunay. + /// The move guarantees a valid causal combinatorial triangulation. It does + /// not preserve or require the empty-sphere property of the coordinates. /// @tparam Generator A uniform random bit generator type /// @param t_manifold The simplicial manifold /// @param generator Caller-owned generator whose state advances during the /// move /// @returns The Expected (3,2) moved manifold or an Unexpected + /// @note The source manifold is unchanged on success and failure. template [[nodiscard]] inline auto do_32_move(Manifold const& t_manifold, Generator& generator) -> Expected @@ -251,11 +359,19 @@ namespace ergodic_moves [[nodiscard]] inline auto find_adjacent_31_cell(Cell_handle const& t_cell) -> std::optional { - if (t_cell->info() != 13) { return std::nullopt; } // NOLINT + if (t_cell == nullptr || + !foliated_triangulations::is_cell_type_correct<3>(t_cell) || + foliated_triangulations::expected_cell_type<3>(t_cell) != + Cell_type::ONE_THREE) + { + return std::nullopt; + } for (auto i = 0; i < 4; ++i) { - if (foliated_triangulations::expected_cell_type<3>(t_cell->neighbor(i)) == - Cell_type::THREE_ONE) + auto const neighbor = t_cell->neighbor(i); + if (foliated_triangulations::is_cell_type_correct<3>(neighbor) && + foliated_triangulations::expected_cell_type<3>(neighbor) == + Cell_type::THREE_ONE) { return std::make_optional(i); } @@ -263,29 +379,25 @@ namespace ergodic_moves return std::nullopt; } // find_26_move() - /// @brief Perform a (2,6) move - /// @details A (2,6) move inserts a vertex into the spacelike face between a - /// (1,3) simplex on the bottom connected to a (3,1) simplex on top. - /// This adds 2 (1,3) simplices and 2 (3,1) simplices. - /// It adds 2 spacelike faces and 6 timelike faces. - /// It also adds 2 timelike edges and 3 spacelike edges, as well as the - /// vertex. - /// This function calls find_adjacent_31_cell on (1,3) simplices drawn from a - /// randomly shuffled container until it succeeds or runs out of simplices. - /// If successful, the triangulation is no longer Delaunay. - /// @image html 26.png - /// @image latex 26.eps width=7cm - /// @tparam Generator A uniform random bit generator type - /// @param t_manifold The simplicial manifold - /// @param generator Caller-owned generator whose state advances during the - /// move - /// @param only_first_site Whether to examine only one uniformly selected site - /// @returns The Expected (2,6) moved manifold or an Unexpected - template - [[nodiscard]] inline auto do_26_move_impl(Manifold const& t_manifold, - Generator& generator, - bool const only_first_site) - -> Expected + namespace detail + { + template + requires std::predicate + [[nodiscard]] inline auto do_26_move_impl( + Manifold const& t_manifold, Generator& generator, + bool const only_first_site, + Post_mutation_validator post_mutation_validator) -> Expected; + } // namespace detail + + // Internal validation seam used to test rejection after mutation. + template + requires std::predicate + [[nodiscard]] inline auto detail::do_26_move_impl( + Manifold const& t_manifold, Generator& generator, + bool const only_first_site, + Post_mutation_validator post_mutation_validator) -> Expected { static auto constexpr INCIDENT_CELLS_FOR_6_2_MOVE = 6; Delaunay triangulation{t_manifold.delaunay_snapshot()}; @@ -303,7 +415,8 @@ namespace ergodic_moves } else { - // Shuffle the container to pick a random sequence of (1,3) cells to try. + // Shuffle the container to pick a random sequence of (1,3) cells to + // try. std::ranges::shuffle(one_three, generator); } for (auto const& bottom : one_three) @@ -359,7 +472,8 @@ namespace ergodic_moves return std::unexpected(msg); } - // Each incident cell should be combinatorially and geometrically valid + // Each incident cell should be combinatorially and geometrically + // valid if (auto check_cells = std::ranges::all_of(incident_cells, [&triangulation](auto const& cell) { @@ -382,19 +496,52 @@ namespace ergodic_moves // Final checks // is_valid() checks for combinatorial and geometric validity - if (!triangulation.tds().is_valid(v_center, true, 1)) + if (!post_mutation_validator( + static_cast(triangulation)) || + !triangulation.tds().is_valid(v_center, true, 1)) { std::string const msg = "v_center is invalid.\n"; return std::unexpected(msg); } - return detail::make_manifold(std::move(triangulation), t_manifold); + return make_manifold(std::move(triangulation), t_manifold); } } // We've run out of (1,3) simplices to try std::string const msg = "No (2,6) move possible.\n"; return std::unexpected(msg); - } // do_26_move() + } // do_26_move_impl() + + /// @brief Perform a (2,6) move + /// @details A (2,6) move inserts a vertex into the spacelike face between a + /// (1,3) simplex on the bottom connected to a (3,1) simplex on top. + /// This adds 2 (1,3) simplices and 2 (3,1) simplices. + /// It adds 2 spacelike faces and 6 timelike faces. + /// It also adds 2 timelike edges and 3 spacelike edges, as well as the + /// vertex. + /// This function calls find_adjacent_31_cell on (1,3) simplices drawn from a + /// randomly shuffled container until it succeeds or runs out of simplices. + /// The move guarantees a valid causal combinatorial triangulation. It does + /// not preserve or require the empty-sphere property of the coordinates. + /// @image html 26.png + /// @image latex 26.eps width=7cm + /// @tparam Generator A uniform random bit generator type + /// @param t_manifold The simplicial manifold + /// @param generator Caller-owned generator whose state advances during the + /// move + /// @param only_first_site Whether to examine only one uniformly selected site + /// @returns The Expected (2,6) moved manifold or an Unexpected + /// @note The source manifold is unchanged on success and failure. + + template + [[nodiscard]] inline auto do_26_move_impl(Manifold const& t_manifold, + Generator& generator, + bool const only_first_site) + -> Expected + { + return detail::do_26_move_impl(t_manifold, generator, only_first_site, + detail::accept_post_mutation); + } /// @brief Perform a (2,6) move using a caller-owned random stream. template @@ -460,30 +607,39 @@ namespace ergodic_moves auto const incident_22 = cell_type_count(Cell_type::TWO_TWO); auto const incident_13 = cell_type_count(Cell_type::ONE_THREE); - // All cells should be classified - if (incident_13 + incident_22 + incident_31 != 6) // NOLINT + // All cells should be causally classified and carry matching metadata. + if (incident_13 + incident_22 + incident_31 != 6 || // NOLINT + !std::ranges::all_of(incident_cells, [](auto const& cell) { + return foliated_triangulations::is_cell_type_correct<3>(cell); + })) { - spdlog::warn("Some incident cells on this vertex need to be fixed.\n"); + return false; } return incident_31 == 3 && incident_22 == 0 && incident_13 == 3; } // find_62_moves() - /// @brief Apply the combinatorial (6,2) retriangulation on a private copy. - /// @details A (3,2) flip of either timelike edge incident to the removable - /// vertex leaves that vertex with degree four. Removing it from its maximal - /// simplex then replaces the original six cells with the required two cells. - /// Keeping both operations in a private copy makes rejection failure-atomic. - /// @param source_triangulation The triangulation containing the candidate - /// @param source_candidate The degree-five vertex to remove - /// @param generator The caller-owned random engine used to order flip paths - /// @returns The moved triangulation, or nullopt when the topology is not - /// flippable or the result violates a triangulation or causal-cell invariant - template - [[nodiscard]] inline auto try_62_move(Delaunay const& source_triangulation, - Vertex_handle const source_candidate, - Generator& generator) + namespace detail + { + template + requires std::predicate + [[nodiscard]] inline auto try_62_move_impl( + Delaunay const& source_triangulation, + Vertex_handle const source_candidate, Generator& generator, + Post_mutation_validator post_mutation_validator) + -> std::optional; + } // namespace detail + + // Internal validation seam used to test rejection after mutation. + template + requires std::predicate + [[nodiscard]] inline auto detail::try_62_move_impl( + Delaunay const& source_triangulation, + Vertex_handle const source_candidate, Generator& generator, + Post_mutation_validator post_mutation_validator) -> std::optional { if (!source_triangulation.tds().is_vertex(source_candidate) || @@ -522,7 +678,8 @@ namespace ergodic_moves } tds.remove_from_maximal_dimension_simplex(candidate); - if (!tds.is_valid() || + if (!post_mutation_validator(static_cast(triangulation)) || + !tds.is_valid() || triangulation.number_of_finite_cells() + 4 != old_cells || triangulation.number_of_vertices() + 1 != old_vertices) { @@ -540,6 +697,26 @@ namespace ergodic_moves } return triangulation; + } // try_62_move_impl() + + /// @brief Apply the combinatorial (6,2) retriangulation on a private copy. + /// @details A (3,2) flip of either timelike edge incident to the removable + /// vertex leaves that vertex with degree four. Removing it from its maximal + /// simplex then replaces the original six cells with the required two cells. + /// Keeping both operations in a private copy makes rejection failure-atomic. + /// @param source_triangulation The triangulation containing the candidate + /// @param source_candidate The degree-five vertex to remove + /// @param generator The caller-owned random engine used to order flip paths + /// @returns The moved triangulation, or nullopt when the topology is not + /// flippable or the result violates a triangulation or causal-cell invariant + template + [[nodiscard]] inline auto try_62_move(Delaunay const& source_triangulation, + Vertex_handle const source_candidate, + Generator& generator) + -> std::optional + { + return detail::try_62_move_impl(source_triangulation, source_candidate, + generator, detail::accept_post_mutation); } // try_62_move() /// @brief Perform a (6,2) move @@ -555,14 +732,15 @@ namespace ergodic_moves /// This function calls is_62_movable() on a randomly shuffled container /// of vertices until it succeeds or runs out of vertices. /// - /// If successful, the triangulation remains Delaunay. (Other moves may - /// change this, however.) + /// The move guarantees a valid causal combinatorial triangulation. It does + /// not preserve or require the empty-sphere property of the coordinates. /// /// @tparam Generator A uniform random bit generator type /// @param t_manifold The simplicial manifold /// @param generator Caller-owned generator whose state advances during the /// move /// @returns The Expected (6,2) moved manifold or Unexpected + /// @note The source manifold is unchanged on success and failure. template [[nodiscard]] inline auto do_62_move(Manifold const& t_manifold, Generator& generator) -> Expected @@ -613,25 +791,7 @@ namespace ergodic_moves Delaunay_t<3> const& triangulation, Edge_handle const& edge) -> std::optional { - if (!detail::is_well_formed_edge(edge) || - !triangulation.tds().is_edge(edge.first, edge.second, edge.third)) - { - return std::nullopt; - } - // Create the circulator of cells around the edge, starting with the cell - // containing the edge - auto circulator = triangulation.incident_cells(edge, edge.first); - Cell_container incident_cells; - // Add cells to the container until we get back to the first one in the - // circulator - do - { // NOLINT(cppcoreguidelines-avoid-do-while) - // Ignore cells containing the infinite vertex - if (triangulation.is_infinite(circulator)) { continue; } - incident_cells.emplace_back(circulator); - } - while (++circulator != edge.first); - return incident_cells; + return detail::finite_incident_cells(triangulation, edge); } // incident_cells_from_edge() /// @brief Find a bistellar flip location @@ -646,20 +806,27 @@ namespace ergodic_moves Delaunay const& triangulation, Edge_handle const& t_edge_candidate) -> std::optional { - if (auto incident_cells = - incident_cells_from_edge(triangulation, t_edge_candidate); - incident_cells.has_value() && incident_cells->size() == 4) + if (!detail::is_well_formed_edge(t_edge_candidate)) { return std::nullopt; } + auto incident_cells = + incident_cells_from_edge(triangulation, t_edge_candidate); + if (!incident_cells || incident_cells->size() != 4) { return std::nullopt; } + + auto const first_time = + t_edge_candidate.first->vertex(t_edge_candidate.second)->info(); + auto const second_time = + t_edge_candidate.first->vertex(t_edge_candidate.third)->info(); + if (first_time != second_time) { return std::nullopt; } + + auto const cell_type_count = [&](Cell_type const type) { + return std::ranges::count_if(*incident_cells, [&](auto const cell) { + return foliated_triangulations::is_cell_type_correct<3>(cell) && + foliated_triangulations::expected_cell_type<3>(cell) == type; + }); + }; + if (cell_type_count(Cell_type::THREE_ONE) == 2 && + cell_type_count(Cell_type::ONE_THREE) == 2) { - auto const cell_type_count = [&](Cell_type const type) { - return std::ranges::count_if(*incident_cells, [&](auto const cell) { - return cell->info() == static_cast(type); - }); - }; - if (cell_type_count(Cell_type::THREE_ONE) == 2 && - cell_type_count(Cell_type::ONE_THREE) == 2) - { - return incident_cells.value(); - } + return incident_cells; } return std::nullopt; } // find_bistellar_flip_location() @@ -675,22 +842,27 @@ namespace ergodic_moves return incident_cells_from_edge(triangulation, edge); } // get_incident_cells() - /// @brief Perform a bistellar flip on triangulation via the given edge - /// @details Pass by value to avoid modifying the original triangulation - /// in the event that the flip is unsuccessful. - /// @param source_triangulation The triangulation to flip - /// @param source_edge The edge to pivot on - /// @param source_top Top vertex of the cells being flipped - /// @param source_bottom Bottom vertex of the cells being flipped - /// @returns A flipped triangulation or nullopt - /// @see [Pachner moves](../REFERENCES.md#pachner-moves) - [[nodiscard]] inline auto bistellar_flip(Delaunay const& source_triangulation, - Edge_handle const source_edge, - Vertex_handle const source_top, - Vertex_handle const source_bottom) + namespace detail + { + template + requires std::predicate + [[nodiscard]] inline auto bistellar_flip_impl( + Delaunay const& source_triangulation, Edge_handle const source_edge, + Vertex_handle const source_top, Vertex_handle const source_bottom, + Post_mutation_validator post_mutation_validator) + -> std::optional; + } // namespace detail + + // Internal validation seam used to test rejection after mutation. + template + requires std::predicate + [[nodiscard]] inline auto detail::bistellar_flip_impl( + Delaunay const& source_triangulation, Edge_handle const source_edge, + Vertex_handle const source_top, Vertex_handle const source_bottom, + Post_mutation_validator post_mutation_validator) -> std::optional { - if (!detail::is_well_formed_edge(source_edge) || source_top == nullptr || + if (!is_well_formed_edge(source_edge) || source_top == nullptr || source_bottom == nullptr || !source_triangulation.tds().is_edge( source_edge.first, source_edge.second, source_edge.third) || @@ -803,7 +975,11 @@ namespace ergodic_moves return std::nullopt; } - if (!triangulation.tds().is_valid()) { return std::nullopt; } + if (!post_mutation_validator(static_cast(triangulation)) || + !triangulation.tds().is_valid()) + { + return std::nullopt; + } for (auto const cell : triangulation.finite_cell_handles()) { @@ -812,6 +988,26 @@ namespace ergodic_moves } return triangulation; + } // bistellar_flip_impl() + + /// @brief Perform a bistellar flip on triangulation via the given edge + /// @details Pass by value to avoid modifying the original triangulation + /// in the event that the flip is unsuccessful. + /// @param source_triangulation The triangulation to flip + /// @param source_edge The edge to pivot on + /// @param source_top Top vertex of the cells being flipped + /// @param source_bottom Bottom vertex of the cells being flipped + /// @returns A flipped triangulation or nullopt + /// @see [Pachner moves](../REFERENCES.md#pachner-moves) + [[nodiscard]] inline auto bistellar_flip(Delaunay const& source_triangulation, + Edge_handle const source_edge, + Vertex_handle const source_top, + Vertex_handle const source_bottom) + -> std::optional + { + return detail::bistellar_flip_impl(source_triangulation, source_edge, + source_top, source_bottom, + detail::accept_post_mutation); } // bistellar_flip() /// @return The center edge of a 4-cell complex @@ -866,6 +1062,7 @@ namespace ergodic_moves /// @param generator Caller-owned generator whose state advances during the /// move /// @return The Expected (4,4) moved manifold or Unexpected + /// @note The source manifold is unchanged on success and failure. template [[nodiscard]] inline auto do_44_move(Manifold const& t_manifold, Generator& generator) -> Expected @@ -979,7 +1176,10 @@ namespace ergodic_moves "Selected (4,4) proposal site could not be flipped.\n"); } - /// @brief Check move correctness + /// @brief Check tracked move deltas and the complete CDT manifold invariant + /// @details This verifies causal foliation, cell metadata, TDS validity, + /// geometry counts, time bounds, and preserved foliation parameters. It does + /// not require Euclidean Delaunayhood. /// @param t_before The manifold before the move /// @param t_after The manifold after the move /// @param t_move The type of move @@ -989,10 +1189,19 @@ namespace ergodic_moves move_tracker::move_type const& t_move) -> bool { + if (!t_after.is_correct() || + !detail::same_configuration_value(t_after.initial_radius(), + t_before.initial_radius()) || + !detail::same_configuration_value(t_after.foliation_spacing(), + t_before.foliation_spacing())) + { + return false; + } + switch (t_move) { case move_tracker::move_type::FOUR_FOUR: - return t_after.is_valid() && t_after.N3() == t_before.N3() && + return t_after.N3() == t_before.N3() && t_after.N3_31() == t_before.N3_31() && t_after.N3_22() == t_before.N3_22() && t_after.N3_13() == t_before.N3_13() && @@ -1003,7 +1212,7 @@ namespace ergodic_moves t_after.max_time() == t_before.max_time() && t_after.min_time() == t_before.min_time(); case move_tracker::move_type::TWO_THREE: - return t_after.is_valid() && t_after.N3() == t_before.N3() + 1 && + return t_after.N3() == t_before.N3() + 1 && t_after.N3_31() == t_before.N3_31() && t_after.N3_22() == t_before.N3_22() + 1 && t_after.N3_13() == t_before.N3_13() && @@ -1015,7 +1224,7 @@ namespace ergodic_moves t_after.max_time() == t_before.max_time() && t_after.min_time() == t_before.min_time(); case move_tracker::move_type::THREE_TWO: - return t_after.is_valid() && t_after.N3() == t_before.N3() - 1 && + return t_after.N3() == t_before.N3() - 1 && t_after.N3_31() == t_before.N3_31() && t_after.N3_22() == t_before.N3_22() - 1 && t_after.N3_13() == t_before.N3_13() && @@ -1027,7 +1236,7 @@ namespace ergodic_moves t_after.max_time() == t_before.max_time() && t_after.min_time() == t_before.min_time(); case move_tracker::move_type::TWO_SIX: - return t_after.is_valid() && t_after.N3() == t_before.N3() + 4 && + return t_after.N3() == t_before.N3() + 4 && t_after.N3_31() == t_before.N3_31() + 2 && t_after.N3_22() == t_before.N3_22() && t_after.N3_13() == t_before.N3_13() + 2 && @@ -1039,7 +1248,7 @@ namespace ergodic_moves t_after.max_time() == t_before.max_time() && t_after.min_time() == t_before.min_time(); case move_tracker::move_type::SIX_TWO: - return t_after.is_valid() && t_after.N3() == t_before.N3() - 4 && + return t_after.N3() == t_before.N3() - 4 && t_after.N3_31() == t_before.N3_31() - 2 && t_after.N3_22() == t_before.N3_22() && t_after.N3_13() == t_before.N3_13() - 2 && diff --git a/include/Metropolis.hpp b/include/Metropolis.hpp index cfa9b5ea3..cfa2d6894 100644 --- a/include/Metropolis.hpp +++ b/include/Metropolis.hpp @@ -344,8 +344,7 @@ class MoveStrategy ++m_attempted_moves[move]; auto candidate = propose_candidate(current, move); - if (!candidate || !candidate->is_correct() || - !ergodic_moves::check_move(current, *candidate, move)) + if (!candidate || !ergodic_moves::check_move(current, *candidate, move)) { ++m_failed_moves[move]; ++m_rejected_moves[move]; diff --git a/include/Move_command.hpp b/include/Move_command.hpp index 69dd92b1e..6c211e02f 100644 --- a/include/Move_command.hpp +++ b/include/Move_command.hpp @@ -132,8 +132,7 @@ class MoveCommand generator); result) { - if (result->is_correct() && - ergodic_moves::check_move(m_manifold, *result, move_type)) + if (ergodic_moves::check_move(m_manifold, *result, move_type)) { swap(result.value(), m_manifold); ++m_succeeded[as_integer(move_type)]; diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index d9f64b528..b50eb472c 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -3,6 +3,7 @@ add_executable( ${PROJECT_SOURCE_DIR}/tests/main.cpp Apply_move_test.cpp Bistellar_flip_test.cpp + Ergodic_moves_3_audit_test.cpp Ergodic_moves_3_test.cpp Foliated_triangulation_test.cpp Function_ref_test.cpp diff --git a/tests/Ergodic_moves_3_audit_test.cpp b/tests/Ergodic_moves_3_audit_test.cpp new file mode 100644 index 000000000..26e7a3d81 --- /dev/null +++ b/tests/Ergodic_moves_3_audit_test.cpp @@ -0,0 +1,1027 @@ +/******************************************************************************* + Causal Dynamical Triangulations in C++ using CGAL + + Copyright © 2026 Adam Getchell + ******************************************************************************/ + +/// @file Ergodic_moves_3_audit_test.cpp +/// @brief Independent scientific audit tests for the 2+1D CDT move set + +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include "Ergodic_moves_3.hpp" + +namespace +{ + using Manifold = manifolds::Manifold_3; + using Delaunay = ergodic_moves::Delaunay; + using Cell_handle = ergodic_moves::Cell_handle; + using Edge_handle = ergodic_moves::Edge_handle; + using Vertex_handle = ergodic_moves::Vertex_handle; + + using Vertex_key = std::tuple; + using Double_representation = std::array; + + [[nodiscard]] auto double_representation(double const value) noexcept + -> Double_representation + { return std::bit_cast(value); } + + template + using Simplex_key = std::array; + + using Cell_record = std::pair, Int_precision>; + + struct Direct_counts + { + Int_precision n0{}; + Int_precision n1_sl{}; + Int_precision n1_tl{}; + Int_precision n2{}; + Int_precision n3_31{}; + Int_precision n3_22{}; + Int_precision n3_13{}; + Int_precision n3{}; + + auto operator==(Direct_counts const&) const -> bool = default; + }; + + struct Count_delta + { + Int_precision n0{}; + Int_precision n1_sl{}; + Int_precision n1_tl{}; + Int_precision n2{}; + Int_precision n3_31{}; + Int_precision n3_22{}; + Int_precision n3_13{}; + Int_precision n3{}; + }; + + struct Canonical_state + { + std::vector vertices; + std::vector> edges; + std::vector> facets; + std::vector cells; + Direct_counts counts; + std::vector> spacelike_facets; + Int_precision min_time{}; + Int_precision max_time{}; + Double_representation initial_radius{}; + Double_representation foliation_spacing{}; + + auto operator==(Canonical_state const&) const -> bool = default; + }; + + [[nodiscard]] auto vertex_key(Vertex_handle const vertex) -> Vertex_key + { + return {CGAL::to_double(vertex->point().x()), + CGAL::to_double(vertex->point().y()), + CGAL::to_double(vertex->point().z()), vertex->info()}; + } + + template + [[nodiscard]] auto sorted_simplex( + std::array const& vertices) -> Simplex_key + { + Simplex_key result{}; + std::ranges::transform(vertices, result.begin(), vertex_key); + std::ranges::sort(result); + return result; + } + + [[nodiscard]] auto independent_cell_type(Cell_handle const cell) -> Cell_type + { + std::array times{}; + for (int index = 0; index < 4; ++index) + { + times.at(static_cast(index)) = cell->vertex(index)->info(); + } + auto const [minimum, maximum] = std::ranges::minmax(times); + if (maximum - minimum != 1) { return Cell_type::ACAUSAL; } + + auto const lower_count = std::ranges::count(times, minimum); + if (lower_count == 3) { return Cell_type::THREE_ONE; } + if (lower_count == 2) { return Cell_type::TWO_TWO; } + if (lower_count == 1) { return Cell_type::ONE_THREE; } + return Cell_type::UNCLASSIFIED; + } + + void assign_independent_cell_metadata(Delaunay& triangulation) + { + for (auto cell = triangulation.finite_cells_begin(); + cell != triangulation.finite_cells_end(); ++cell) + { + cell->info() = static_cast(independent_cell_type(cell)); + } + } + + [[nodiscard]] auto direct_counts(Delaunay const& triangulation) + -> Direct_counts + { + Direct_counts result{}; + result.n0 = static_cast(triangulation.number_of_vertices()); + result.n2 = + static_cast(triangulation.number_of_finite_facets()); + result.n3 = + static_cast(triangulation.number_of_finite_cells()); + + for (auto edge = triangulation.finite_edges_begin(); + edge != triangulation.finite_edges_end(); ++edge) + { + auto const first_time = edge->first->vertex(edge->second)->info(); + auto const second_time = edge->first->vertex(edge->third)->info(); + if (first_time == second_time) { ++result.n1_sl; } + else + { + ++result.n1_tl; + } + } + + for (auto cell = triangulation.finite_cells_begin(); + cell != triangulation.finite_cells_end(); ++cell) + { + switch (independent_cell_type(cell)) + { + case Cell_type::THREE_ONE: ++result.n3_31; break; + case Cell_type::TWO_TWO: ++result.n3_22; break; + case Cell_type::ONE_THREE: ++result.n3_13; break; + case Cell_type::ACAUSAL: + case Cell_type::UNCLASSIFIED: break; + } + } + return result; + } + + [[nodiscard]] auto canonical_state(Manifold const& manifold) + -> Canonical_state + { + auto triangulation = manifold.delaunay_snapshot(); + Canonical_state result{}; + result.counts = direct_counts(triangulation); + result.min_time = manifold.min_time(); + result.max_time = manifold.max_time(); + result.initial_radius = double_representation(manifold.initial_radius()); + result.foliation_spacing = + double_representation(manifold.foliation_spacing()); + + for (auto vertex = triangulation.finite_vertices_begin(); + vertex != triangulation.finite_vertices_end(); ++vertex) + { + result.vertices.emplace_back(vertex_key(vertex)); + } + + for (auto edge = triangulation.finite_edges_begin(); + edge != triangulation.finite_edges_end(); ++edge) + { + result.edges.emplace_back(sorted_simplex<2>({ + edge->first->vertex(edge->second), + edge->first->vertex(edge->third), + })); + } + + std::map spacelike_facets; + for (auto facet = triangulation.finite_facets_begin(); + facet != triangulation.finite_facets_end(); ++facet) + { + std::array vertices{}; + auto output = vertices.begin(); + for (int index = 0; index < 4; ++index) + { + if (index != facet->second) { *output++ = facet->first->vertex(index); } + } + auto key = sorted_simplex<3>(vertices); + if (std::get<3>(key[0]) == std::get<3>(key[1]) && + std::get<3>(key[1]) == std::get<3>(key[2])) + { + ++spacelike_facets[std::get<3>(key[0])]; + } + result.facets.emplace_back(std::move(key)); + } + + for (auto cell = triangulation.finite_cells_begin(); + cell != triangulation.finite_cells_end(); ++cell) + { + result.cells.emplace_back( + sorted_simplex<4>({cell->vertex(0), cell->vertex(1), cell->vertex(2), + cell->vertex(3)}), + cell->info()); + } + + std::ranges::sort(result.vertices); + std::ranges::sort(result.edges); + std::ranges::sort(result.facets); + std::ranges::sort(result.cells); + result.spacelike_facets.assign(spacelike_facets.begin(), + spacelike_facets.end()); + return result; + } + + [[nodiscard]] auto canonical_triangulation(Delaunay const& triangulation) + -> Canonical_state + { + Canonical_state result{}; + result.counts = direct_counts(triangulation); + std::map spacelike_facets; + + for (auto vertex = triangulation.finite_vertices_begin(); + vertex != triangulation.finite_vertices_end(); ++vertex) + { + result.vertices.emplace_back(vertex_key(vertex)); + } + if (!result.vertices.empty()) + { + auto const times = + result.vertices | + std::views::transform([](auto const& v) { return std::get<3>(v); }); + auto const [minimum, maximum] = std::ranges::minmax(times); + result.min_time = minimum; + result.max_time = maximum; + } + + for (auto edge = triangulation.finite_edges_begin(); + edge != triangulation.finite_edges_end(); ++edge) + { + result.edges.emplace_back(sorted_simplex<2>({ + edge->first->vertex(edge->second), + edge->first->vertex(edge->third), + })); + } + + for (auto facet = triangulation.finite_facets_begin(); + facet != triangulation.finite_facets_end(); ++facet) + { + std::array vertices{}; + auto output = vertices.begin(); + for (int index = 0; index < 4; ++index) + { + if (index != facet->second) { *output++ = facet->first->vertex(index); } + } + auto key = sorted_simplex<3>(vertices); + if (std::get<3>(key[0]) == std::get<3>(key[1]) && + std::get<3>(key[1]) == std::get<3>(key[2])) + { + ++spacelike_facets[std::get<3>(key[0])]; + } + result.facets.emplace_back(std::move(key)); + } + + for (auto cell = triangulation.finite_cells_begin(); + cell != triangulation.finite_cells_end(); ++cell) + { + result.cells.emplace_back( + sorted_simplex<4>({cell->vertex(0), cell->vertex(1), cell->vertex(2), + cell->vertex(3)}), + cell->info()); + } + + std::ranges::sort(result.vertices); + std::ranges::sort(result.edges); + std::ranges::sort(result.facets); + std::ranges::sort(result.cells); + result.spacelike_facets.assign(spacelike_facets.begin(), + spacelike_facets.end()); + return result; + } + + void check_independent_invariants(Manifold const& manifold) + { + auto const triangulation = manifold.delaunay_snapshot(); + auto const counts = direct_counts(triangulation); + + auto const production_is_foliated = manifold.is_foliated(); + auto const production_checks_vertices = manifold.check_vertices(); + auto const production_checks_simplices = manifold.check_simplices(); + auto const production_is_correct = manifold.is_correct(); + CAPTURE(production_is_foliated); + CAPTURE(production_checks_vertices); + CAPTURE(production_checks_simplices); + CAPTURE(production_is_correct); + + CHECK_EQ(triangulation.dimension(), 3); + CHECK(triangulation.tds().is_valid()); + + std::optional raw_min_time; + std::optional raw_max_time; + std::size_t raw_vertex_count{}; + for (auto vertex = triangulation.finite_vertices_begin(); + vertex != triangulation.finite_vertices_end(); ++vertex) + { + CHECK(triangulation.tds().is_vertex(vertex)); + raw_min_time = raw_min_time ? std::min(*raw_min_time, vertex->info()) + : vertex->info(); + raw_max_time = raw_max_time ? std::max(*raw_max_time, vertex->info()) + : vertex->info(); + ++raw_vertex_count; + } + REQUIRE(raw_min_time.has_value()); + REQUIRE(raw_max_time.has_value()); + CHECK_EQ(raw_vertex_count, static_cast(counts.n0)); + CHECK_EQ(*raw_min_time, manifold.min_time()); + CHECK_EQ(*raw_max_time, manifold.max_time()); + + CHECK_EQ(counts.n0, manifold.N0()); + CHECK_EQ(counts.n0, manifold.vertices()); + CHECK_EQ(counts.n1_sl + counts.n1_tl, manifold.N1()); + CHECK_EQ(counts.n1_sl + counts.n1_tl, manifold.edges()); + CHECK_EQ(counts.n1_sl, manifold.N1_SL()); + CHECK_EQ(counts.n1_tl, manifold.N1_TL()); + CHECK_EQ(counts.n2, manifold.N2()); + CHECK_EQ(counts.n2, manifold.faces()); + CHECK_EQ(counts.n3, manifold.N3()); + CHECK_EQ(counts.n3, manifold.simplices()); + CHECK_EQ(counts.n3_31, manifold.N3_31()); + CHECK_EQ(counts.n3_22, manifold.N3_22()); + CHECK_EQ(counts.n3_13, manifold.N3_13()); + CHECK_EQ(counts.n3_31 + counts.n3_13, manifold.N3_31_13()); + CHECK_EQ(counts.n0 - counts.n1_sl - counts.n1_tl + counts.n2 - counts.n3, + 1); + + std::map spacelike_facets; + for (auto facet = triangulation.finite_facets_begin(); + facet != triangulation.finite_facets_end(); ++facet) + { + std::optional time; + auto spacelike = true; + for (int index = 0; index < 4; ++index) + { + if (index == facet->second) { continue; } + auto const vertex_time = facet->first->vertex(index)->info(); + if (!time) { time = vertex_time; } + else if (*time != vertex_time) { spacelike = false; } + } + if (spacelike) { ++spacelike_facets[*time]; } + } + for (auto time = manifold.min_time(); time <= manifold.max_time(); ++time) + { + CHECK_EQ(manifold.spacelike_face_count(time), spacelike_facets[time]); + } + + std::size_t classified_cell_count{}; + for (auto cell = triangulation.finite_cells_begin(); + cell != triangulation.finite_cells_end(); ++cell) + { + Cell_handle const cell_handle = cell; + auto const type = independent_cell_type(cell); + CHECK(triangulation.tds().is_cell(cell_handle)); + CHECK((type == Cell_type::THREE_ONE || type == Cell_type::TWO_TWO || + type == Cell_type::ONE_THREE)); + CHECK_EQ(cell->info(), static_cast(type)); + ++classified_cell_count; + for (int index = 0; index < 4; ++index) + { + auto const neighbor = cell->neighbor(index); + REQUIRE(neighbor != nullptr); + CHECK(triangulation.tds().is_cell(neighbor)); + REQUIRE(neighbor->has_neighbor(cell_handle)); + CHECK(neighbor->neighbor(neighbor->index(cell_handle)) == cell_handle); + } + } + CHECK_EQ(classified_cell_count, static_cast(counts.n3)); + } + + void check_delta(Manifold const& before, Manifold const& after, + Count_delta const expected) + { + auto const lhs = direct_counts(before.delaunay_snapshot()); + auto const rhs = direct_counts(after.delaunay_snapshot()); + CHECK_EQ(rhs.n0 - lhs.n0, expected.n0); + CHECK_EQ(rhs.n1_sl - lhs.n1_sl, expected.n1_sl); + CHECK_EQ(rhs.n1_tl - lhs.n1_tl, expected.n1_tl); + CHECK_EQ(rhs.n2 - lhs.n2, expected.n2); + CHECK_EQ(rhs.n3_31 - lhs.n3_31, expected.n3_31); + CHECK_EQ(rhs.n3_22 - lhs.n3_22, expected.n3_22); + CHECK_EQ(rhs.n3_13 - lhs.n3_13, expected.n3_13); + CHECK_EQ(rhs.n3 - lhs.n3, expected.n3); + } + + [[nodiscard]] auto make_23_fixture() -> Manifold + { + static auto constexpr radius_2 = 2.0 * std::numbers::inv_sqrt3_v; + static auto constexpr sqrt_2 = std::numbers::sqrt2_v; + std::vector> vertices{ + { 1, 0, 0}, + { 0, 1, 0}, + { 0, 0, 1}, + {radius_2, radius_2, radius_2}, + { sqrt_2, sqrt_2, 0}, + }; + std::vector const times{1, 1, 1, 2, 2}; + return Manifold{manifolds::make_causal_vertices<3>(vertices, times)}; + } + + [[nodiscard]] auto make_26_fixture() -> Manifold + { + static auto constexpr radius_2 = 2.0 * std::numbers::inv_sqrt3_v; + std::vector> vertices{ + { 0, 0, 0}, + { 1, 0, 0}, + { 0, 1, 0}, + { 0, 0, 1}, + {radius_2, radius_2, radius_2}, + }; + std::vector const times{0, 1, 1, 1, 2}; + return Manifold{manifolds::make_causal_vertices<3>(vertices, times)}; + } + + [[nodiscard]] auto make_44_fixture() -> Manifold + { + static auto constexpr inverse_sqrt_2 = 1.0 / std::numbers::sqrt2_v; + std::vector> vertices{ + { 0, 0, 0}, + { inverse_sqrt_2, 0, inverse_sqrt_2}, + { 0, inverse_sqrt_2, inverse_sqrt_2}, + {-inverse_sqrt_2, 0, inverse_sqrt_2}, + { 0, -inverse_sqrt_2, inverse_sqrt_2}, + { 0, 0, 2}, + }; + std::vector const times{0, 1, 1, 1, 1, 2}; + return Manifold{manifolds::make_causal_vertices<3>(vertices, times), 0, 1}; + } + + [[nodiscard]] auto finite_incident_cells(Delaunay const& triangulation, + Edge_handle const edge) + -> std::optional + { + auto circulator = triangulation.incident_cells(edge, edge.first); + auto const done = circulator; + ergodic_moves::Cell_container cells; + do + { // NOLINT(cppcoreguidelines-avoid-do-while) + if (triangulation.is_infinite(circulator)) { return std::nullopt; } + cells.emplace_back(circulator); + } + while (++circulator != done); + return cells; + } + + [[nodiscard]] auto find_inverse_32_edge(Delaunay const& triangulation) + -> std::optional + { + for (auto edge = triangulation.finite_edges_begin(); + edge != triangulation.finite_edges_end(); ++edge) + { + auto const first_time = edge->first->vertex(edge->second)->info(); + auto const second_time = edge->first->vertex(edge->third)->info(); + if (first_time == second_time) { continue; } + auto const cells = finite_incident_cells(triangulation, *edge); + if (!cells || cells->size() != 3) { continue; } + auto const count = [&](Cell_type const type) { + return std::ranges::count_if(*cells, [&](auto const cell) { + return independent_cell_type(cell) == type; + }); + }; + if (count(Cell_type::TWO_TWO) == 2 && + count(Cell_type::THREE_ONE) + count(Cell_type::ONE_THREE) == 1) + { + return *edge; + } + } + return std::nullopt; + } + + [[nodiscard]] auto find_44_pivot(Delaunay const& triangulation) + -> std::optional + { + for (auto edge = triangulation.finite_edges_begin(); + edge != triangulation.finite_edges_end(); ++edge) + { + auto const first_time = edge->first->vertex(edge->second)->info(); + auto const second_time = edge->first->vertex(edge->third)->info(); + if (first_time != second_time) { continue; } + auto const cells = finite_incident_cells(triangulation, *edge); + if (!cells || cells->size() != 4) { continue; } + auto const count = [&](Cell_type const type) { + return std::ranges::count_if(*cells, [&](auto const cell) { + return independent_cell_type(cell) == type; + }); + }; + if (count(Cell_type::THREE_ONE) == 2 && count(Cell_type::ONE_THREE) == 2) + { + return *edge; + } + } + return std::nullopt; + } + + [[nodiscard]] auto find_vertex(Delaunay const& triangulation, + Point_t<3> const& point) + -> std::optional + { + for (auto vertex = triangulation.finite_vertices_begin(); + vertex != triangulation.finite_vertices_end(); ++vertex) + { + if (vertex->point() == point) { return vertex; } + } + return std::nullopt; + } + + void check_23_32_round_trip(std::uint64_t const seed) + { + auto const original = make_23_fixture(); + auto const before = canonical_state(original); + cdt::Random random{seed}; + CAPTURE(seed); + auto forward = ergodic_moves::do_23_move(original, random); + REQUIRE(forward.has_value()); + check_independent_invariants(*forward); + auto inverse = ergodic_moves::do_32_move(*forward, random); + REQUIRE(inverse.has_value()); + check_independent_invariants(*inverse); + CHECK_EQ(canonical_state(*inverse), before); + } + + void check_26_62_round_trip(std::uint64_t const seed) + { + auto const original = make_26_fixture(); + auto const before = canonical_state(original); + cdt::Random random{seed}; + CAPTURE(seed); + auto forward = ergodic_moves::do_26_move(original, random); + REQUIRE(forward.has_value()); + check_independent_invariants(*forward); + auto inverse = ergodic_moves::do_62_move(*forward, random); + REQUIRE(inverse.has_value()); + check_independent_invariants(*inverse); + CHECK_EQ(canonical_state(*inverse), before); + } + + void check_44_round_trip(std::uint64_t const seed) + { + auto const original = make_44_fixture(); + auto const before = canonical_state(original); + cdt::Random random{seed}; + CAPTURE(seed); + auto forward = ergodic_moves::do_44_move(original, random); + REQUIRE(forward.has_value()); + check_independent_invariants(*forward); + auto inverse = ergodic_moves::do_44_move(*forward, random); + REQUIRE(inverse.has_value()); + check_independent_invariants(*inverse); + CHECK_EQ(canonical_state(*inverse), before); + } +} // namespace + +SCENARIO("Every 2+1D CDT move has the literature-derived local delta" * + doctest::test_suite("ergodic-audit")) +{ + GIVEN("the minimal causal cavity for a (2,3) move") + { + auto const before = make_23_fixture(); + cdt::Random random{10623}; + WHEN("the cavity is replaced by three tetrahedra") + { + auto const after = ergodic_moves::do_23_move(before, random); + REQUIRE(after.has_value()); + THEN("the independent f-vector delta is (0,0,1,2,0,1,0,1)") + { + check_delta(before, *after, {0, 0, 1, 2, 0, 1, 0, 1}); + check_independent_invariants(*after); + } + AND_WHEN("the inverse (3,2) move is applied") + { + auto const restored = ergodic_moves::do_32_move(*after, random); + REQUIRE(restored.has_value()); + THEN("the canonical complex is restored without comparing handles") + { + check_delta(*after, *restored, {0, 0, -1, -2, 0, -1, 0, -1}); + CHECK_EQ(canonical_state(*restored), canonical_state(before)); + } + } + } + } + + GIVEN("the minimal causal cavity for a (2,6) move") + { + auto const before = make_26_fixture(); + cdt::Random random{10626}; + WHEN("a vertex is inserted into the shared spacelike triangle") + { + auto const after = ergodic_moves::do_26_move(before, random); + REQUIRE(after.has_value()); + THEN("the independent f-vector delta is (1,3,2,8,2,0,2,4)") + { + check_delta(before, *after, {1, 3, 2, 8, 2, 0, 2, 4}); + check_independent_invariants(*after); + } + AND_WHEN("the inverse (6,2) move removes the inserted vertex") + { + auto const restored = ergodic_moves::do_62_move(*after, random); + REQUIRE(restored.has_value()); + THEN("the canonical complex is restored without comparing handles") + { + check_delta(*after, *restored, {-1, -3, -2, -8, -2, 0, -2, -4}); + CHECK_EQ(canonical_state(*restored), canonical_state(before)); + } + } + } + } + + GIVEN("the minimal four-tetrahedron diamond for a (4,4) move") + { + auto const before = make_44_fixture(); + cdt::Random random{10644}; + WHEN("the spacelike diagonal is exchanged") + { + auto const after = ergodic_moves::do_44_move(before, random); + REQUIRE(after.has_value()); + THEN("every tracked count is unchanged") + { + check_delta(before, *after, {}); + check_independent_invariants(*after); + CHECK_NE(canonical_state(*after), canonical_state(before)); + } + AND_WHEN("the inverse diagonal exchange is applied") + { + auto const restored = ergodic_moves::do_44_move(*after, random); + REQUIRE(restored.has_value()); + THEN("the original canonical diamond is restored") + { CHECK_EQ(canonical_state(*restored), canonical_state(before)); } + } + } + } +} + +SCENARIO("CDT move rejection is causal and failure-atomic" * + doctest::test_suite("ergodic-audit")) +{ + GIVEN("an empty source manifold") + { + Manifold const source; + auto const before = canonical_state(source); + cdt::Random random{106}; + WHEN("each unsupported move is attempted") + { + CHECK_FALSE(ergodic_moves::do_23_move(source, random)); + CHECK_FALSE(ergodic_moves::do_32_move(source, random)); + CHECK_FALSE(ergodic_moves::do_26_move(source, random)); + CHECK_FALSE(ergodic_moves::do_62_move(source, random)); + CHECK_FALSE(ergodic_moves::do_44_move(source, random)); + THEN("the complete canonical state and caches are unchanged") + { CHECK_EQ(canonical_state(source), before); } + } + } + + GIVEN("a geometrically flippable timelike edge with the wrong CDT cavity") + { + auto triangulation = make_23_fixture().delaunay_snapshot(); + auto two_two = foliated_triangulations::filter_cells<3>( + foliated_triangulations::collect_cells<3>(triangulation), + Cell_type::TWO_TWO); + REQUIRE_EQ(two_two.size(), 1); + REQUIRE(ergodic_moves::try_23_move(triangulation, two_two.front())); + assign_independent_cell_metadata(triangulation); + auto const inverse_edge = find_inverse_32_edge(triangulation); + REQUIRE(inverse_edge.has_value()); + + auto const first = inverse_edge->first->vertex(inverse_edge->second); + auto const second = inverse_edge->first->vertex(inverse_edge->third); + auto const lower_time = std::min(first->info(), second->info()); + auto const incident = finite_incident_cells(triangulation, *inverse_edge); + REQUIRE(incident.has_value()); + for (auto const cell : *incident) + { + for (int index = 0; index < 4; ++index) + { + auto const vertex = cell->vertex(index); + if (vertex != first && vertex != second) + { + vertex->info() = lower_time; + } + } + } + assign_independent_cell_metadata(triangulation); + auto const before = canonical_triangulation(triangulation); + + auto control = triangulation; + auto const control_first = find_vertex(control, first->point()); + auto const control_second = find_vertex(control, second->point()); + REQUIRE(control_first.has_value()); + REQUIRE(control_second.has_value()); + Cell_handle control_cell = nullptr; + int control_i{}; + int control_j{}; + REQUIRE(control.is_edge(*control_first, *control_second, control_cell, + control_i, control_j)); + REQUIRE(control.flip(control_cell, control_i, control_j)); + + WHEN("the CDT (3,2) helper checks the site") + { + CHECK_FALSE(ergodic_moves::try_32_move(triangulation, *inverse_edge)); + THEN("it rejects before CGAL can create an acausal tetrahedron") + { CHECK_EQ(canonical_triangulation(triangulation), before); } + } + } + + GIVEN("valid topology with stale or incompatible causal metadata") + { + auto const corrupt_metadata = [](Cell_handle const& cell) { + auto const expected = + foliated_triangulations::expected_cell_type<3>(cell); + cell->info() = static_cast(expected == Cell_type::TWO_TWO + ? Cell_type::THREE_ONE + : Cell_type::TWO_TWO); + }; + + auto move_23 = make_23_fixture().delaunay_snapshot(); + auto two_two = foliated_triangulations::filter_cells<3>( + foliated_triangulations::collect_cells<3>(move_23), Cell_type::TWO_TWO); + REQUIRE_EQ(two_two.size(), 1); + corrupt_metadata(two_two.front()); + auto const before_23 = canonical_triangulation(move_23); + + auto move_32 = make_23_fixture().delaunay_snapshot(); + auto move_32_cells = foliated_triangulations::filter_cells<3>( + foliated_triangulations::collect_cells<3>(move_32), Cell_type::TWO_TWO); + REQUIRE_EQ(move_32_cells.size(), 1); + REQUIRE(ergodic_moves::try_23_move(move_32, move_32_cells.front())); + assign_independent_cell_metadata(move_32); + auto const inverse_edge = find_inverse_32_edge(move_32); + REQUIRE(inverse_edge.has_value()); + auto const inverse_cells = finite_incident_cells(move_32, *inverse_edge); + REQUIRE(inverse_cells.has_value()); + corrupt_metadata(inverse_cells->front()); + auto const before_32 = canonical_triangulation(move_32); + + auto move_26 = make_26_fixture().delaunay_snapshot(); + auto one_three = foliated_triangulations::filter_cells<3>( + foliated_triangulations::collect_cells<3>(move_26), + Cell_type::ONE_THREE); + REQUIRE_EQ(one_three.size(), 1); + corrupt_metadata(one_three.front()); + auto const before_26 = canonical_triangulation(move_26); + + auto const source_62 = make_26_fixture(); + cdt::Random random_62{10662}; + auto const expanded_62 = ergodic_moves::do_26_move(source_62, random_62); + REQUIRE(expanded_62.has_value()); + auto move_62 = expanded_62->delaunay_snapshot(); + auto vertices_62 = foliated_triangulations::collect_vertices<3>(move_62); + auto const candidate_62 = + std::ranges::find_if(vertices_62, [&](auto const& vertex) { + return ergodic_moves::is_62_movable(move_62, vertex); + }); + REQUIRE(candidate_62 != vertices_62.end()); + ergodic_moves::Cell_container incident_62; + move_62.tds().incident_cells(*candidate_62, + std::back_inserter(incident_62)); + REQUIRE_EQ(incident_62.size(), 6); + corrupt_metadata(incident_62.front()); + auto const before_62 = canonical_triangulation(move_62); + + auto move_44 = make_44_fixture().delaunay_snapshot(); + auto const pivot = find_44_pivot(move_44); + REQUIRE(pivot.has_value()); + auto const pivot_cells = finite_incident_cells(move_44, *pivot); + REQUIRE(pivot_cells.has_value()); + corrupt_metadata(pivot_cells->front()); + auto const before_44 = canonical_triangulation(move_44); + + WHEN("the corresponding cavity selectors inspect the sites") + { + CHECK_FALSE(ergodic_moves::try_23_move(move_23, two_two.front())); + CHECK_FALSE(ergodic_moves::try_32_move(move_32, *inverse_edge)); + CHECK_FALSE(ergodic_moves::find_adjacent_31_cell(one_three.front())); + CHECK_FALSE(ergodic_moves::is_62_movable(move_62, *candidate_62)); + CHECK_FALSE(ergodic_moves::find_bistellar_flip_location(move_44, *pivot)); + THEN("no topology is mutated") + { + CHECK_EQ(canonical_triangulation(move_23), before_23); + CHECK_EQ(canonical_triangulation(move_32), before_32); + CHECK_EQ(canonical_triangulation(move_26), before_26); + CHECK_EQ(canonical_triangulation(move_62), before_62); + CHECK_EQ(canonical_triangulation(move_44), before_44); + } + } + } + + GIVEN("a topological cavity with an invalid time assignment") + { + auto triangulation = make_23_fixture().delaunay_snapshot(); + auto two_two = foliated_triangulations::filter_cells<3>( + foliated_triangulations::collect_cells<3>(triangulation), + Cell_type::TWO_TWO); + REQUIRE_EQ(two_two.size(), 1); + two_two.front()->vertex(0)->info() += 3; + assign_independent_cell_metadata(triangulation); + auto const before = canonical_triangulation(triangulation); + + WHEN("the (2,3) selector checks the now-acausal tetrahedron") + { + CHECK_FALSE(ergodic_moves::try_23_move(triangulation, two_two.front())); + THEN("the invalid site is rejected without mutation") + { CHECK_EQ(canonical_triangulation(triangulation), before); } + } + } + + GIVEN("a (2,6) move rejected by validation after vertex insertion") + { + auto const source = make_26_fixture(); + auto const before = canonical_state(source); + auto const source_counts = direct_counts(source.delaunay_snapshot()); + cdt::Random random{106'260}; + bool observed_mutation = false; + + WHEN("the injected count validator rejects the private copy") + { + auto const rejected = ergodic_moves::detail::do_26_move_impl( + source, random, false, [&](Delaunay const& mutated) { + auto const mutated_counts = direct_counts(mutated); + observed_mutation = true; + CHECK_EQ(mutated_counts.n0, source_counts.n0 + 1); + CHECK_EQ(mutated_counts.n3, source_counts.n3 + 4); + return mutated_counts.n0 == source_counts.n0; + }); + + THEN("the post-mutation rejection leaves the source unchanged") + { + REQUIRE(observed_mutation); + CHECK_FALSE(rejected.has_value()); + CHECK_EQ(canonical_state(source), before); + } + } + } + + GIVEN("a (6,2) move rejected by validation after vertex removal") + { + auto const seed = make_26_fixture(); + cdt::Random setup_random{106'620}; + auto const expanded = ergodic_moves::do_26_move(seed, setup_random); + REQUIRE(expanded.has_value()); + + auto source = expanded->delaunay_snapshot(); + auto vertices = foliated_triangulations::collect_vertices<3>(source); + auto const candidate = + std::ranges::find_if(vertices, [&](auto const& vertex) { + return ergodic_moves::is_62_movable(source, vertex); + }); + REQUIRE(candidate != vertices.end()); + auto const source_counts = direct_counts(source); + auto const before = canonical_triangulation(source); + cdt::Random random{106'621}; + bool observed_mutation = false; + + WHEN("the injected count validator rejects the private copy") + { + auto const rejected = ergodic_moves::detail::try_62_move_impl( + source, *candidate, random, [&](Delaunay const& mutated) { + auto const mutated_counts = direct_counts(mutated); + observed_mutation = true; + CHECK_EQ(mutated_counts.n0, source_counts.n0 - 1); + CHECK_EQ(mutated_counts.n3, source_counts.n3 - 4); + return mutated_counts.n0 == source_counts.n0; + }); + + THEN("the post-mutation rejection leaves the source unchanged") + { + REQUIRE(observed_mutation); + CHECK_FALSE(rejected.has_value()); + CHECK_EQ(canonical_triangulation(source), before); + } + } + } + + GIVEN("a bistellar flip rejected by validation after both TDS flips") + { + auto source = make_44_fixture().delaunay_snapshot(); + auto const pivot = find_44_pivot(source); + REQUIRE(pivot.has_value()); + auto const incident = finite_incident_cells(source, *pivot); + REQUIRE(incident.has_value()); + + auto const pivot_first = pivot->first->vertex(pivot->second); + auto const pivot_second = pivot->first->vertex(pivot->third); + Vertex_handle top = nullptr; + Vertex_handle bottom = nullptr; + for (auto const cell : *incident) + { + for (int index = 0; index < 4; ++index) + { + auto const vertex = cell->vertex(index); + if (vertex == pivot_first || vertex == pivot_second) { continue; } + if (top == nullptr || vertex->info() > top->info()) { top = vertex; } + if (bottom == nullptr || vertex->info() < bottom->info()) + { + bottom = vertex; + } + } + } + REQUIRE(top != nullptr); + REQUIRE(bottom != nullptr); + auto const before = canonical_triangulation(source); + bool observed_mutation = false; + + WHEN("the injected TDS validator rejects the private copy") + { + auto const rejected = ergodic_moves::detail::bistellar_flip_impl( + source, *pivot, top, bottom, [&](Delaunay const& mutated) { + observed_mutation = true; + CHECK(mutated.tds().is_valid()); + CHECK_NE(canonical_triangulation(mutated), before); + return !mutated.tds().is_valid(); + }); + + THEN("the post-mutation rejection leaves the source unchanged") + { + REQUIRE(observed_mutation); + CHECK_FALSE(rejected.has_value()); + CHECK_EQ(canonical_triangulation(source), before); + } + } + } + + GIVEN("causal manifolds without the required move incidences") + { + auto const source_23 = make_23_fixture(); + auto const source_26 = make_26_fixture(); + auto const before_23 = canonical_state(source_23); + auto const before_26 = canonical_state(source_26); + cdt::Random random{1060}; + + WHEN("degree, composition, or flippability preconditions are absent") + { + CHECK_FALSE(ergodic_moves::do_32_move(source_23, random)); + CHECK_FALSE(ergodic_moves::do_26_move(source_23, random)); + CHECK_FALSE(ergodic_moves::do_44_move(source_23, random)); + CHECK_FALSE(ergodic_moves::do_62_move(source_26, random)); + THEN("each source remains canonically unchanged") + { + CHECK_EQ(canonical_state(source_23), before_23); + CHECK_EQ(canonical_state(source_26), before_26); + } + } + } + + GIVEN("malformed and nonflippable edge representations") + { + auto triangulation = make_44_fixture().delaunay_snapshot(); + auto foreign = make_44_fixture().delaunay_snapshot(); + auto const foreign_pivot = find_44_pivot(foreign); + REQUIRE(foreign_pivot.has_value()); + + std::optional boundary_edge; + for (auto edge = triangulation.finite_edges_begin(); + edge != triangulation.finite_edges_end(); ++edge) + { + if (!finite_incident_cells(triangulation, *edge)) + { + boundary_edge = *edge; + break; + } + } + REQUIRE(boundary_edge.has_value()); + + Edge_handle const malformed{nullptr, 0, 1}; + auto const before = canonical_triangulation(triangulation); + WHEN("the checked edge helpers receive them") + { + CHECK_FALSE(ergodic_moves::get_incident_cells(triangulation, malformed)); + CHECK_FALSE( + ergodic_moves::get_incident_cells(triangulation, *boundary_edge)); + CHECK_FALSE(ergodic_moves::find_bistellar_flip_location(triangulation, + *foreign_pivot)); + CHECK_FALSE(ergodic_moves::try_32_move(triangulation, malformed)); + THEN("the triangulation remains unchanged") + { CHECK_EQ(canonical_triangulation(triangulation), before); } + } + } +} + +TEST_CASE("Seeded CDT move replay preserves invariants and inverse structure" * + doctest::test_suite("ergodic-audit")) +{ + for (auto const seed : + std::array{0, 1, 2, 92, 106, 20'260'721}) + { + check_23_32_round_trip(seed); + check_26_62_round_trip(seed); + check_44_round_trip(seed); + } +} + +TEST_CASE("Move validation rejects exact configuration-state drift" * + doctest::test_suite("ergodic-audit")) +{ + auto const before = make_44_fixture(); + auto const after = Manifold{ + foliated_triangulations::FoliatedTriangulation_3{ + before.delaunay_snapshot(), -0.0, before.foliation_spacing()} + }; + REQUIRE(after.is_correct()); + CHECK_FALSE(ergodic_moves::check_move(before, after, + move_tracker::move_type::FOUR_FOUR)); +}