From 92e8103fa97d490ddadc58a26794b0b631a219e3 Mon Sep 17 00:00:00 2001 From: d-burg Date: Fri, 23 Oct 2026 18:47:35 -0400 Subject: [PATCH 1/2] Tearing - BUGFIX! - Polish roots with a real 2x2 Jacobian so the anti-holomorphic SLAYER residual converges SLAYER's solve_inner evaluates conj(Delta_hat(g)) by design (Schwarz reflection, to keep the scan-plane contour orientation consistent), so the tearing residual f(Q) = Delta' - S^(1/3) Delta(Q) - Delta_c is a function of conj(Q). _polish_root modelled f with a complex secant slope, which assumes f is holomorphic. On an anti-holomorphic f that step is the mirror image of the true step: an offset along Im Q is doubled, every halving backtrack is also worse, and the polish returns the unpolished point. The validity gate then rejects any root whose raw contour intersection is not already accurate, which hits steep roots first: on an ARC V3A 2/1 surface the least-damped root at the electron diamagnetic frequency (gamma = -0.029 kHz) was discarded and a far root at gamma = -13.5 kHz reported instead. - Model f with the real 2x2 Jacobian of (Re f, Im f) in (Re Q, Im Q) from forward differences along both axes; exact for holomorphic and anti-holomorphic residuals. Trust region, backtracking and the improve-or-no-op contract are unchanged; one extra f evaluation per iteration. - Tests: anti-holomorphic linear and steep nonlinear residuals, a root offset purely along Im Q (all fail with the complex-secant polish), and an extraction-level check that a steep anti-holomorphic root from a coarse scan passes the validity gate at the true zero. Co-Authored-By: Claude Opus 5.5 --- .../Dispersion/GrowthRateExtraction.jl | 40 ++++++++++++------ test/runtests_dispersion_polish.jl | 42 +++++++++++++++++++ 2 files changed, 70 insertions(+), 12 deletions(-) diff --git a/src/Tearing/Dispersion/GrowthRateExtraction.jl b/src/Tearing/Dispersion/GrowthRateExtraction.jl index 45f5bc74dd..254c5a9b16 100644 --- a/src/Tearing/Dispersion/GrowthRateExtraction.jl +++ b/src/Tearing/Dispersion/GrowthRateExtraction.jl @@ -112,7 +112,10 @@ single-surface scans; `mc.surfaces[mc.ref_idx].tauk` for coupled scans). reports the raw marching-squares estimate, which is sensitive to grid and floating-point rounding. - `polish_maxit` -- max polish iterations per root (default 20). Each costs a - handful of `f` evaluations; failures fall back to the unpolished point. + handful of `f` evaluations (two for the real 2×2 Jacobian, one trial step, + at most three backtracks); failures fall back to the unpolished point. The + polish does not assume `f` is holomorphic, so it also converges on the SLAYER + residual, which is a function of conj(Q). - `validity_rtol` -- with `residual` supplied, a polished root is kept only if `|residual| < validity_rtol · median(|Δ|)` (default 1e-3). Intersections that don't polish to a true zero (contour near-misses, or surfaces whose @@ -500,9 +503,19 @@ end # Bounded Gauss-Newton minimization of |f(Q)| within a disk of radius `R` # around `Q0`. Returns (Q_best, n_evals, improved). `improved == false` ⇒ -# Q0 is returned unchanged (hard fallback). Secant derivative + trust-region -# projection + |f|-decrease backtracking keep it strictly "improve-or-no-op": -# it cannot diverge, jump to another root, or be pulled into a pole. +# Q0 is returned unchanged (hard fallback). Trust-region projection + +# |f|-decrease backtracking keep it strictly "improve-or-no-op": it cannot +# diverge, jump to another root, or be pulled into a pole. +# +# The local linear model is the REAL 2×2 Jacobian of (Re f, Im f) with respect +# to (Re Q, Im Q), from forward differences along both axes, not a complex +# derivative. A complex slope assumes `f` is holomorphic, but the SLAYER +# tearing residual is a function of conj(Q) (`solve_inner` evaluates +# conj(Δ̂(ĝ)); see Riccati.jl). On such a residual a complex-secant step is +# the mirror image of the true step, so the polish stalls and the validity +# gate rejects the root, most often the steep, weakly damped root at the +# diamagnetic frequency. The real Jacobian is exact for holomorphic and +# anti-holomorphic residuals alike. function _polish_root(f, Q0::ComplexF64, R::Float64; maxit::Int=8, tol_step::Float64=1e-12, tol_f::Float64=1e-8) (R > 0) || return (Q0, 0, false) @@ -512,17 +525,21 @@ function _polish_root(f, Q0::ComplexF64, R::Float64; maxit::Int=8, a0 = abs(f0) a0 == 0 && return (Q0, nev, false) - # Second point for the secant slope: a small real offset scaled to R. + h = max(R * 1e-3, eps(Float64)) # finite-difference offset, scaled to R Qk, fk = Q0, f0 - Qp = Q0 + complex(max(R * 1e-3, eps(Float64)), 0.0) - fp = ComplexF64(f(Qp)) - nev += 1 Qbest, abest = Q0, a0 for _ in 1:maxit - df = (Qp == Qk) ? ComplexF64(0) : (fp - fk) / (Qp - Qk) - (isfinite(df) && abs(df) > 0) || break - step = -fk / df # zero of local linear model + fx = ComplexF64(f(Qk + h)) + fy = ComplexF64(f(Qk + im * h)) + nev += 2 + (isfinite(fx) && isfinite(fy)) || break + j11, j12 = real(fx - fk) / h, real(fy - fk) / h + j21, j22 = imag(fx - fk) / h, imag(fy - fk) / h + det = j11 * j22 - j12 * j21 + (isfinite(det) && det != 0) || break + # zero of the local linear model: J * [dRe, dIm] = -[Re f, Im f] + step = complex((-real(fk) * j22 + j12 * imag(fk)) / det, (-j11 * imag(fk) + j21 * real(fk)) / det) Qt = Qk + step if abs(Qt - Q0) > R # project into trust region Qt = Q0 + R * (Qt - Q0) / abs(Qt - Q0) @@ -538,7 +555,6 @@ function _polish_root(f, Q0::ComplexF64, R::Float64; maxit::Int=8, bt += 1 end if isfinite(ft) && abs(ft) < abest - Qp, fp = Qk, fk Qk, fk = Qt, ft Qbest, abest = Qt, abs(ft) (abs(step) < tol_step * (abs(Qk) + tol_step) || abest < tol_f * a0) && break diff --git a/test/runtests_dispersion_polish.jl b/test/runtests_dispersion_polish.jl index 8bc4c97025..b97e662944 100644 --- a/test/runtests_dispersion_polish.jl +++ b/test/runtests_dispersion_polish.jl @@ -25,6 +25,35 @@ @test abs(Qp - r) < 1e-7 end + # The SLAYER residual is a function of conj(Q) (solve_inner evaluates + # conj(Δ̂(ĝ))), so the polish must not assume a holomorphic f. + @testset "_polish_root: anti-holomorphic residual converges" begin + r = 0.42 + 0.27im + f(Q) = (2 - 1im) * conj(ComplexF64(Q) - r) + Qp, _, improved = D._polish_root(f, r + (0.013 - 0.009im), 0.1) + @test improved + @test abs(Qp - r) < 1e-10 + @test abs(f(Qp)) < 1e-9 + end + + @testset "_polish_root: anti-holomorphic root offset along Im Q (mirrored complex step)" begin + # A complex-slope step on f = a·conj(Q - r) maps an offset i·y to 2i·y and + # every halving backtrack stays worse, so a complex-secant polish stalls here. + r = -0.75 - 0.009im + f(Q) = 1e3 * conj(ComplexF64(Q) - r) + Qp, _, improved = D._polish_root(f, r + 0.0017im, 0.01) + @test improved + @test abs(Qp - r) < 1e-10 + end + + @testset "_polish_root: steep nonlinear anti-holomorphic residual converges" begin + r = -0.75 - 0.009im + f(Q) = 1e3 * conj((ComplexF64(Q) - r) * (1 + 0.4 * (ComplexF64(Q) - r))) + Qp, _, improved = D._polish_root(f, r + (0.004 + 0.003im), 0.02) + @test improved + @test abs(Qp - r) < 1e-9 + end + @testset "_polish_root: no-op-or-better invariant (never worse)" begin # The polish must never return a point with larger |f| than the start. for r in (0.1 + 0.2im, -0.5 + 0.0im, 0.3 - 0.4im) @@ -181,4 +210,17 @@ @test !(:spurious in gr3.warning_flags) @test !isnan(real(gr3.Q_root)) end + + @testset "find_growth_rates: steep anti-holomorphic root passes the validity gate" begin + # A steep residual leaves the raw contour intersection far above the gate + # threshold, so the root is kept only if the polish reaches the true zero. + r = 0.42 + 0.27im + f(Q) = 1e3 * conj((ComplexF64(Q) - r) * (1 + 0.4 * (ComplexF64(Q) - r))) + scan = brute_force_scan(f, (0.0, 1.0), (0.0, 0.6); nre=9, nim=7, threaded=false) # raw |f| above the gate + # the pole classifier's threshold is absolute: scale it with f, as SLAYER's adaptive one does + gr = find_growth_rates(scan, 1.0; residual=f, validity_rtol=1e-3, pole_threshold=1e3 * 10.0) + @test !(:spurious in gr.warning_flags) + @test !isnan(real(gr.Q_root)) + @test abs(gr.Q_root - r) < 1e-8 + end end From aac35ab68ffc3e6eee4571ee3e8d86c7b0ae146c Mon Sep 17 00:00:00 2001 From: d-burg Date: Wed, 28 Oct 2026 16:23:03 -0400 Subject: [PATCH 2/2] Tearing - TEST - Bound the polish's one-iteration evaluation count and test a residual that is neither holomorphic nor anti-holomorphic Co-Authored-By: Claude Opus 5.5 --- src/Tearing/Dispersion/GrowthRateExtraction.jl | 17 +++++------------ test/runtests_dispersion_polish.jl | 10 +++++++++- 2 files changed, 14 insertions(+), 13 deletions(-) diff --git a/src/Tearing/Dispersion/GrowthRateExtraction.jl b/src/Tearing/Dispersion/GrowthRateExtraction.jl index 254c5a9b16..d16667889c 100644 --- a/src/Tearing/Dispersion/GrowthRateExtraction.jl +++ b/src/Tearing/Dispersion/GrowthRateExtraction.jl @@ -507,15 +507,8 @@ end # |f|-decrease backtracking keep it strictly "improve-or-no-op": it cannot # diverge, jump to another root, or be pulled into a pole. # -# The local linear model is the REAL 2×2 Jacobian of (Re f, Im f) with respect -# to (Re Q, Im Q), from forward differences along both axes, not a complex -# derivative. A complex slope assumes `f` is holomorphic, but the SLAYER -# tearing residual is a function of conj(Q) (`solve_inner` evaluates -# conj(Δ̂(ĝ)); see Riccati.jl). On such a residual a complex-secant step is -# the mirror image of the true step, so the polish stalls and the validity -# gate rejects the root, most often the steep, weakly damped root at the -# diamagnetic frequency. The real Jacobian is exact for holomorphic and -# anti-holomorphic residuals alike. +# The local model is the real 2×2 forward-difference Jacobian of (Re f, Im f), exact for holomorphic +# and anti-holomorphic f alike; SLAYER's residual is a function of conj(Q) (see Riccati.jl). function _polish_root(f, Q0::ComplexF64, R::Float64; maxit::Int=8, tol_step::Float64=1e-12, tol_f::Float64=1e-8) (R > 0) || return (Q0, 0, false) @@ -536,10 +529,10 @@ function _polish_root(f, Q0::ComplexF64, R::Float64; maxit::Int=8, (isfinite(fx) && isfinite(fy)) || break j11, j12 = real(fx - fk) / h, real(fy - fk) / h j21, j22 = imag(fx - fk) / h, imag(fy - fk) / h - det = j11 * j22 - j12 * j21 - (isfinite(det) && det != 0) || break + jdet = j11 * j22 - j12 * j21 + (isfinite(jdet) && jdet != 0) || break # zero of the local linear model: J * [dRe, dIm] = -[Re f, Im f] - step = complex((-real(fk) * j22 + j12 * imag(fk)) / det, (-j11 * imag(fk) + j21 * real(fk)) / det) + step = complex((-real(fk) * j22 + j12 * imag(fk)) / jdet, (-j11 * imag(fk) + j21 * real(fk)) / jdet) Qt = Qk + step if abs(Qt - Q0) > R # project into trust region Qt = Q0 + R * (Qt - Q0) / abs(Qt - Q0) diff --git a/test/runtests_dispersion_polish.jl b/test/runtests_dispersion_polish.jl index b97e662944..43d06363b2 100644 --- a/test/runtests_dispersion_polish.jl +++ b/test/runtests_dispersion_polish.jl @@ -14,7 +14,7 @@ @test improved @test abs(Qp - r) < 1e-10 @test abs(f(Qp)) < 1e-10 - @test nev <= 3 * 8 + 2 # bounded eval budget + @test nev <= 1 + 6 # one iteration: 2 differences, the trial step, ≤ 3 backtracks end @testset "_polish_root: nonlinear residual converges" begin @@ -25,6 +25,14 @@ @test abs(Qp - r) < 1e-7 end + @testset "_polish_root: residual that is neither holomorphic nor anti-holomorphic converges" begin + r = 0.15 - 0.35im + f(Q) = (ComplexF64(Q) - r) + 0.5 * conj(ComplexF64(Q) - r) # zero only at r + Qp, _, improved = D._polish_root(f, r + 0.01 + 0.012im, 0.1) + @test improved + @test abs(Qp - r) < 1e-10 + end + # The SLAYER residual is a function of conj(Q) (solve_inner evaluates # conj(Δ̂(ĝ))), so the polish must not assume a holomorphic f. @testset "_polish_root: anti-holomorphic residual converges" begin