-
Notifications
You must be signed in to change notification settings - Fork 8
Description
Summary
RK45 integration in reference coordinates (magfie REFCOORDS, integmode=0) shows a missing class of losses compared to the symplectic runs for the same vmec_to_wall chartmap case.
The missing losses appear as an entire footprint cluster in (theta,zeta) / (R,Z) for chartmap_wall_125mm_rk45_refcoords, while the corresponding symplectic chartmap_wall_125mm run shows that cluster clearly.
Context
In the SQUID chartmap validation dataset (/home/ert/data/SQUID/SIMPLE/rerun_2026_chartmap), we compared:
chartmap_wall_125mm(symplectic,integmode=3, canonical MEISS)chartmap_wall_125mm_rk45_refcoords(RK45,integmode=0,isw_field_type=REFCOORDS)
Both use the same physical wall chartmap:
coord_input = chartmaps/squid_5c_vmec_wall_125mm.chartmap.nc(rho_convention=vmec_to_wall, rho=1 is wall)
Observations
- Loss fraction difference
- Symplectic wall: 190/8192 lost (2.319%)
- RK45 refcoords wall: 79/8192 lost (0.964%)
- Missing footprint cluster
A cluster in the approximate angular window
- theta in [80°, 170°]
- zeta in [0°, 55°] (nfp=4)
accounts for 92 of the symplectic losses, but 0 of the RK45/refcoords losses.
-
Those missing-band particles do not appear to be saved by simply tightening tolerances
A quick RK45/refcoords run with tighter settings (512 particles,relerr=1e-10,ntimstep=4000,npoiper2=128) still showed 0 losses in that band. -
Geometric closeness to the wall surface
For RK45/refcoords wall losses, the lost points are close to the chartmap wall boundary (rho=1), consistent with the loss test being tied to rho.
However, the missing-band particles (lost in symplectic) end RK45/refcoords at rho<1 and typically remain O(10 cm) away from the wall at their final (theta,zeta).
This suggests the discrepancy is not just a plotting artifact; it is a real integration/trajectory difference.
Hypothesis / Suspects
- RK45/refcoords trajectory differs qualitatively because it uses
magfie_refcoords(splined A,h,B derivatives + computed hcurl), while symplectic uses canonical field (MEISS) and a different integration scheme. - There may be a bug or inconsistency in
magfie_refcoordsevaluation for chartmap coordinate systems (metric/basis/Jacobian handling, periodic wrapping, derivative scaling), producing different drift behavior in some regions. - Alternatively, the chamber/loss check path differs between the RK45 and symplectic pipelines.
Note: In this repo, chamb_can is currently a stub (src/chamb_m.f90) that flags outside when y(1) >= 1, which means for RK45/refcoords runs it effectively checks rho>=1 (and ignores toroidal angle).
Reproduction (from dataset repo)
- In
SQUID/SIMPLE/rerun_2026_chartmap/, run the RK45/refcoords wall case:
runs/validation_vmec_theta/chartmap_wall_125mm_rk45_refcoords/simple.inusesisw_field_type=5(REFCOORDS) andintegmode=0.- Execute with
~/code/SIMPLE/build/simple.x simple.inin that directory.
- Compare footprints against the symplectic wall case using
tools/plot_validation.py.
Evidence
- Combined plot set (includes symplectic wall + RK45/refcoords wall):
runs/validation_vmec_theta/fig/triage_125mm_rk45_compare_footprints_theta_zeta.png
Desired next steps
- Add a targeted diagnostic for
magfie_refcoords/ REFCOORDS+RK45 in chartmap coordinates:- Evaluate/plot key quantities along a representative orbit from the missing band (e.g., Bmod, sqrtg, hcurl, and their derivatives)
- Compare against the canonical MEISS field evaluation along the same physical points.
- Verify that REFCOORDS+RK45 path uses consistent wall/loss check logic with other integrators.