ForceFreeStates - BUGFIX - Decouple the ideal-crossing jump width from singfac_min

Abierto
  • ¿Algo no está bien?
#481 3 comentarios 0 reacciones 1 asignado Ver en GitHub

Los mantenedores suelen responder en 3 días

@jhalpern30 ya está trabajando en esto.

Desde el 6/10/2026.

Evaluación

Dificultad
5/5
Tiempo estimado
Más de una semana
Aptitud para principiantes
30/100
Tipo de issue
Error
Claridad
Bastante claro
Estado de actividad
Activo
Stack tecnológico
julia
Área
backend

Línea de trabajo

Start by reading cross_ideal_singular_surf! and the sing_get_ca / sing_get_ua routines mentioned in the issue, along with the discussion on #93. Compare the candidate approaches for preserving U₁†U₂ − U₂†U₁ = 0 across each crossing. Done means the crossing no longer introduces a persistent non-Hermiticity defect; the issue does not specify a chosen approach or a test file.

Escrito por el modelo de indexación a partir del texto del issue.

Descripción

Rewritten 2026-10-13. The original body blamed the trapezoid step and proposed re-projection, ODE sub-steps or extended precision; none of those addresses the actual source. The original text is in the edit history.

Problem

At each ideal singular surface, both integrators (forward and Riccati) stop at ψ_s − d and jump to ψ_s + d, where d = singfac_min/(n·|q′|). The jump zeroes the resonant component of the regular solutions, but that component only vanishes like x at the surface. Zeroing it at x = d instead of 0 leaves an error in W, δW and the energies that is first order in singfac_min:

Case singfac_min 1e-4 (default) 1e-6 1e-8
DIIID-like n=1 forward et₁ 0.801550 0.798864 0.798834
  • 100× smaller jump → 100× smaller change in et₁ (DIIID 2.69e-3 → 3.0e-5; Solovev 3.40e-3 → 3.4e-5), so first order, as expected.
  • Riccati et₁ drifts by 0.7% (n=1) to 10% (n=6) on DIIID-like at the default.
  • The non-Hermiticity of W after each crossing (the original symptom here) is the same truncation error.

Why lowering singfac_min isn't the fix

singfac_min sets four things at once:

  1. the jump half-width (wants to be small);
  2. the Riccati Δ′ matching distance (wants to be large: precision degrades as tol·x^(−2α); DIIID Riccati breaks at 1e-7);
  3. the chunk ends and kinetic chunking;
  4. the knot-free equilibrium grid bracket around each rational (GridRefinement).

After #520 the forward path stays accurate at small singfac_min (Solovev to 1e-8), so forward-only runs can lower it as a workaround, at the cost of moving the grid, the per-surface ca_left/ca_right stubs and kinetic chunking. Riccati cannot.

Proposed fix

Give the jump its own distance and leave singfac_min for everything else. Two candidate designs; prototype both and pick one:

A. ODE legs (prototyped). Add δ_jump (default 1e-6/(min(n)·|q′|)). Inside both crossings (cross_ideal_singular_surf!, riccati_cross_ideal_singular_surf!), integrate the same ODE from ψ_s − d to ψ_s − δ_jump, jump across ±δ_jump, and integrate back out to ψ_s + d. Every handoff point stays at d.

B. Crossing in asymptotic-coefficient space (not prototyped). Write u = V c_A (Glasser 2016 Eq. 50), set the large resonant coefficient to exactly zero, carry the non-resonant coefficients across (Glasser–Kolemen 2018 Eq. 36) and rebuild with ua on the right. This removes the O(d) error at today's 1e-4 with no legs and no new key, and also the d^(2α) residual that resonant-row elimination leaves for α < ½. Open risk: extracting coefficients through an ill-conditioned ua.

Prototype evidence (design A, δ_jump = 1e-6)

Case Today (et₁ / flagged steps) With fix
Forward DIIID n=1 0.801550 / 14 0.798864 / 0
Forward multi-n DIIID −0.575053 / 64 −0.575094 / 0
Riccati DIIID n=1 0.804133 0.798892
  • DIIID-like n=1–6: ≤ 1.1e-3 error at δ_jump = 1e-6; forward and Riccati extrapolations in δ_jump agree within 2e-5.
  • Riccati Δ′ matrix is bit-identical in every case (the STRIDE BVP reads only ua at the matching distance and the chunk propagators). Kinetic is identical.
  • Cost: +5–10% runtime at 1e-6.
  • Developer page with figures and the draft diff: https://claude.ai/artifact/MaMkwYpbRxALD99cT6GAWB

Caveat: the prototype predates #520's final trigger (cond(U₁) on the ucrit key, U₂ / lu(U₁)); re-run the table on #520's head. The earlier "DIIID n=3 forward 2.7× wrong" was mostly forward basis collapse, which #520 fixes, not this jump error (Riccati at the default is within ~1%).

Required corrections before a PR (design A)

  • Capture ca_left/ca_right and the Δca stub at d, outside the legs. They go to HDF5 and SLAYER's fallback, and ua is far worse conditioned at δ_jump.
  • Riccati: the column zeroing assumes (S, I) form, and the left leg ends without a renormalization, so the large solution is effectively removed at d. Renormalize at ψ_s − δ_jump before zeroing, or reuse #489's resonant-row elimination. This is the likely cause of the remaining 7e-5 forward/Riccati gap at DIIID n=3.
  • Scale δ_jump with min(n) in multi-n.
  • Attribute non-Hermiticity flags with a fixed axis start (frobenius_psi_max = 0); the default Frobenius start adds its own m=0 floor (#464).
  • The kinetic crossing (cross_kinetic_singular_surf!) has the same structure at singfac_min; its convergence is untested. Leave it unchanged and say so, or open a follow-up issue.
  • Check reverse-shear double rationals (#492) where two surfaces sit within d.

Impact

  • Every free-boundary energy moves toward convergence by roughly singfac_min relative (DIIID et₁ −0.34%), so the PR takes !. Shipped Fortran inputs use 1e-4, so results move away from Fortran benchmarks by about that amount.
  • Re-pin the harness and add a small-jump Solovev case.

Dependencies

  • #489 (crossing pivots on the resonant row) → #520 (cond(U₁) reduction trigger, resizable reduction storage, stable W). Both open; XB stacks on #520. Without #520, a small jump collapses the forward basis after high-α surfaces.
  • Independent of the Δ′-input work (#506, #517, field-line and grid PRs), since Δ′ is unchanged by construction. Those PRs do move the equilibrium, so take the final before/after numbers on the base this lands on. If design B is chosen, land #517 (trace-aware exponents) first, since B relies on ua.
  • The grid-bracket PR must keep its bracket ≥ d so the matching points and legs stay in one cubic cell.

Related

  • #93: the original report, closed by #480.
  • #490: forward basis collapse undetected by ucrit; fixed by #520 (open).
  • #464: Frobenius axis-start floor; separate from this.
Lenguaje dominante
Julia
Estrellas
11
Forks
7
Merge medio
3 d 14 h
PR fusionados (30 d)
69

Preparar el entorno

  • Sin Dockerfile ni archivo de Docker Compose
  • Tiene una plantilla de pull request
  • Sin guía de contribución

Primeros pasos

  1. Lee el issue completo y luego la guía de contribución del proyecto.
  2. Comenta en el issue que vas a ocuparte — evita que dos personas hagan lo mismo.
  3. Haz un fork del repositorio y trabaja en una rama.
  4. Abre un pull request que haga referencia al número del issue.

Más de OpenFUSIONToolkit/GPEC

Todos los issues de OpenFUSIONToolkit/GPEC

Issues similares

Más issues de Julia

Recibe los nuevos issues en tu correo

Un resumen breve de issues de GitHub para principiantes.