ForceFreeStates - REFACTOR - Single-source the auto chunk-count target

src/ForceFreeStates/EulerLagrange.jl

@@ -61,6 +61,29 @@ function ode_itime_cost(psi1::Float64, psi2::Float64, intr::ForceFreeStatesInter
     return cost
 end
 
+"""
+    min_crossing_chunks(msing) -> Int
+
+Floor on the chunk count set by the `msing` rational-surface crossings; explicit `nchunks`
+requests below it are clamped up (with a warning) by `balance_integration_chunks`.
+"""
+min_crossing_chunks(msing::Integer) = 2 * msing + 3
+
+"""
+    auto_chunk_target(msing) -> Int
+
+The chunk-count target `balance_integration_chunks` resolves `nchunks = 0` to. Beyond the
+crossing floor, BVP propagator conditioning needs at least 8 non-crossing sub-chunks per
+segment (axis→surf₁, surfᵢ→surfᵢ₊₁, surfₙ→edge) plus the crossing chunks — without them,
+`assemble_fm_matrix(condition=true)` cannot keep accumulated products well-conditioned
+because single long-span propagators may already have cond ~ 10²⁴ (STRIDE uses 33 intervals
+for comparable problems). Derived from `msing` alone — never from the thread count — which
+is what keeps the chunk list, and hence every Riccati output, identical whatever `julia -t`
+provides. Tests that steer decompositions relative to the auto target call this instead of
+mirroring the formula.
+"""
+auto_chunk_target(msing::Integer) = max(min_crossing_chunks(msing), 8 * (msing + 1) + msing)
+
 """
     balance_integration_chunks(chunks, ctrl, intr) -> Vector{IntegrationChunk}
 
@@ -79,20 +102,14 @@ each original chunk retains `needs_crossing=true` and the original `ising`, so t
 rational surface crossing still fires at the correct ψ in the serial assembly phase.
 """
 function balance_integration_chunks(chunks::Vector{IntegrationChunk}, ctrl::ForceFreeStatesControl, intr::ForceFreeStatesInternal)
-    min_chunks = 2 * intr.msing + 3
-    # Ensure enough sub-chunks for BVP propagator conditioning: at least 5 non-crossing
-    # sub-chunks per segment (axis→surf₁, surfᵢ→surfᵢ₊₁, surfₙ→edge), plus crossing
-    # chunks. STRIDE uses 33 intervals for comparable problems. Without enough sub-chunks,
-    # assemble_fm_matrix(condition=true) can't keep accumulated products well-conditioned
-    # because single long-span propagators may already have cond ~ 10²⁴.
-    min_bvp_intervals = 8 * (intr.msing + 1) + intr.msing
+    min_chunks = min_crossing_chunks(intr.msing)
     if ctrl.nchunks > 0
         if ctrl.nchunks < min_chunks
             @warn "nchunks = $(ctrl.nchunks) is below the $min_chunks chunks required by $(intr.msing) singular surfaces; clamping up."
         end
         target_n = max(ctrl.nchunks, min_chunks)
     else
-        target_n = max(min_chunks, min_bvp_intervals)
+        target_n = auto_chunk_target(intr.msing)
     end
 
     result = collect(chunks)

test/runtests_decomposition_invariance.jl

@@ -49,11 +49,10 @@ end
     @test auto.delta_prime !== nothing
     msing = size(auto.delta_prime.matrix, 1)
 
-    # The nchunks=0 target, mirroring balance_integration_chunks' internal formula (as
-    # runtests_parallel_integration.jl does). Invariance is asserted ABOVE this floor only:
-    # fewer chunks than the msing-derived minimum is not a different decomposition but a
-    # structurally deficient one (the floor gives the rational-surface crossings room).
-    auto_target = max(2 * msing + 3, 8 * (msing + 1) + msing)
+    # The nchunks=0 target, from its single source. Invariance is asserted ABOVE this floor
+    # only: fewer chunks than the msing-derived minimum is not a different decomposition but
+    # a structurally deficient one (the floor gives the rational-surface crossings room).
+    auto_target = GP_TI.ForceFreeStates.auto_chunk_target(msing)
     finer = _solve_at_nchunks(dir, auto_target + 11)
     @test finer.delta_prime !== nothing
 

test/runtests_parallel_integration.jl

@@ -134,11 +134,11 @@ using TOML
         base_chunks = GeneralizedPerturbedEquilibrium.ForceFreeStates.chunk_el_integration_bounds(odet, ctrl, intr)
         balanced = GeneralizedPerturbedEquilibrium.ForceFreeStates.balance_integration_chunks(base_chunks, ctrl, intr)
 
-        # Must mirror balance_integration_chunks' internal target_n formula for nchunks = 0
-        # (src/ForceFreeStates/EulerLagrange.jl). Keep this in sync. The formula reads only
-        # intr.msing — no thread count enters it, which is what makes Riccati outputs
-        # independent of how many threads `julia -t` provides.
-        target_n = max(2 * intr.msing + 3, 8 * (intr.msing + 1) + intr.msing)
+        # The nchunks = 0 target, from its single source. It reads only intr.msing — no
+        # thread count enters it, which is what makes Riccati outputs independent of how
+        # many threads `julia -t` provides.
+        target_n = GeneralizedPerturbedEquilibrium.ForceFreeStates.auto_chunk_target(intr.msing)
+        @test target_n == max(2 * intr.msing + 3, 8 * (intr.msing + 1) + intr.msing)
 
         # After balancing, chunk count equals target_n: the while-loop adds exactly one
         # chunk per iteration (a bisection split) and exits when length(result) >= target_n,
@@ -183,7 +183,7 @@ using TOML
         @test length(balanced_more) == target_n + 7
 
         # An nchunks below the singular-surface floor is clamped up, with a warning.
-        min_chunks = 2 * intr.msing + 3
+        min_chunks = GeneralizedPerturbedEquilibrium.ForceFreeStates.min_crossing_chunks(intr.msing)
         ctrl_few = GeneralizedPerturbedEquilibrium.ForceFreeStates.ForceFreeStatesControl(;
             (Symbol(k) => v for (k, v) in inputs["ForceFreeStates"])..., nchunks=1)
         balanced_few = @test_logs (:warn,) match_mode=:any GeneralizedPerturbedEquilibrium.ForceFreeStates.balance_integration_chunks(base_chunks, ctrl_few, intr)