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

Đang mở
  • Có gì không ổn?
#481 3 bình luận 0 reaction 1 người được giao Xem trên GitHub

Maintainer thường phản hồi trong vòng 3 ngày

@jhalpern30 đang làm issue này rồi.

Từ ngày 6/10/2026.

Đánh giá

Độ khó
5/5
Thời gian dự kiến
Hơn một tuần
Mức phù hợp với người mới
30/100
Loại issue
Lỗi
Độ rõ ràng
Khá rõ ràng
Mức độ hoạt động
Sôi nổi
Công nghệ
julia
Lĩnh vực
backend

Hướng nghiên cứu

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.

Do mô hình lập chỉ mục viết ra từ nội dung của issue.

Mô tả

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.
Ngôn ngữ chính
Julia
Star
11
Fork
7
Merge trung bình
3 ngày 11 giờ
Pull request đã merge (30 ngày)
68

Chuẩn bị môi trường

  • Không có Dockerfile hay tệp Docker Compose
  • Có mẫu pull request
  • Không có hướng dẫn đóng góp

Bắt đầu từ đâu

  1. Đọc hết issue, rồi đọc hướng dẫn đóng góp của dự án.
  2. Bình luận trên issue rằng bạn sẽ nhận — tránh hai người làm cùng một việc.
  3. Fork repository và làm thay đổi trên một nhánh.
  4. Mở pull request có tham chiếu số hiệu của issue.

Issue khác của OpenFUSIONToolkit/GPEC

Tất cả issue của OpenFUSIONToolkit/GPEC

Issue tương tự

Thêm issue về Julia

Nhận issue mới trong hộp thư của bạn

Bản tóm tắt ngắn những issue GitHub phù hợp với người mới.