ForceFreeStates - BUGFIX - Decouple the ideal-crossing jump width from singfac_min
メンテナーはふだん 3 日以内に返信
@jhalpern30 がすでに取り組んでいます。
2026年10月6日 から。
評価
調査の方向性
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.
索引モデルが issue の本文から書いたものです。
説明
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:
- the jump half-width (wants to be small);
- the Riccati Δ′ matching distance (wants to be large: precision degrades as tol·x^(−2α); DIIID Riccati breaks at 1e-7);
- the chunk ends and kinetic chunking;
- 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
uaat 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_rightand the Δca stub at d, outside the legs. They go to HDF5 and SLAYER's fallback, anduais 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 atsingfac_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_minrelative (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.
- 主要言語
- Julia
- スター
- 11
- フォーク
- 7
- 平均マージ
- 4日 1時間
- マージ済み PR(30日)
- 68
環境構築
- Dockerfile・Docker Compose ファイルなし
- プルリクエストのテンプレートあり
- コントリビューションガイドなし
はじめの一歩
- issue を最後まで読み、次にプロジェクトのコントリビューションガイドを読みます。
- 着手することを issue にコメントします — 二人が同じ作業をするのを防げます。
- リポジトリをフォークし、ブランチを切って変更します。
- issue 番号を参照したプルリクエストを送ります。
OpenFUSIONToolkit/GPEC のほかの issue
-
難易度 1/5 1時間未満 初心者へのやさしさ 82/100
OpenFUSIONToolkit/GPEC#644 ·
メンテナーはふだん 3 日以内に返信
-
難易度 4/5 3〜5日 初心者へのやさしさ 35/100
OpenFUSIONToolkit/GPEC#643 ·
メンテナーはふだん 3 日以内に返信
-
PerturbedEquilibrium - FEATURE - Complete the resistive perturbed-equilibrium workflow from Galerkin and matched results対応中かも @matt-pharr が 2 日前に担当しました。 オープンfeature
OpenFUSIONToolkit/GPEC#635 · コメント 1 件 · 担当者 1 名 ·
メンテナーはふだん 3 日以内に返信
-
test
難易度 4/5 3〜5日 初心者へのやさしさ 25/100
OpenFUSIONToolkit/GPEC#633 · コメント 2 件 ·
メンテナーはふだん 3 日以内に返信
-
Tearing - FEATURE - Add a general modified Rutherford equation for island evolution対応中かも @StuartBenjamin が 2 日前に担当しました。 オープンfeature
OpenFUSIONToolkit/GPEC#630 · 担当者 1 名 ·
メンテナーはふだん 3 日以内に返信
OpenFUSIONToolkit/GPEC の issue をすべて見る
似ている issue
-
authority-de-dommel ready-for-agent
難易度 1/5 1時間未満 初心者へのやさしさ 78/100
Deltares/Ribasim-NL#997 ·
メンテナーはふだん 1 日以内に返信
-
難易度 2/5 1〜3時間 初心者へのやさしさ 76/100
Sienna-Platform/InfrastructureSystems.jl#675 ·
メンテナーはふだん 2 日以内に返信
-
TargetEncoder reuses the first level's shrinkage factor for all levels when m > 0対応中かも @parth-patel010 が今日担当しました。 オープン
難易度 2/5 1〜3時間 初心者へのやさしさ 78/100
JuliaAI/MLJTransforms.jl#95 · コメント 1 件 ·
-
難易度 1/5 1時間未満 初心者へのやさしさ 78/100
epiforecasts/BVDOutbreakSize#1232 · コメント 1 件 ·
メンテナーはふだん 1 日以内に返信
-
難易度 2/5 1〜3時間 初心者へのやさしさ 72/100