efit_to_boozer: tested, clean version with clear COCOS convention

[krystophny]

Migrated from archived itpplasma/efit_to_boozer#2. The efit_to_boozer code now lives in libneo; this tracks providing a tested, clean EQDSK->Boozer path with an explicit, documented COCOS convention. Original: itpplasma/efit_to_boozer#2

[krystophny]

Automated triage pass (Claude), grounded in the repository at fe3f76e. Verify before acting.

Still open and still accurate. The EQDSK→Boozer path has zero direct test coverage, and there is no COCOS handling anywhere in it — the COCOS machinery that libneo does own is in a different module that this path never calls.

Purpose

efit_to_boozer.x is the entry point by which every tokamak equilibrium enters the ITPcp Boozer-based codes (NEO-2, NEO-RT, MEPHIT consume fromefit*.bc). Its output signs — iota, Jpol/nper, Itor, pprime, sqrt g(0,0), and the sine harmonics rmns/zmns/vmns/bmns — depend on the handedness of the input g-file's coordinate convention, i.e. on its COCOS index (Sauter & Medvedev 2013). Today that dependence is handled by a runtime sign heuristic plus a prose comment from 2023, and nothing checks it. A silently wrong sign here does not crash: it produces a plausible-looking .bc file with a flipped rotational transform or a flipped toroidal current, which then propagates into transport coefficients downstream. Pinning the convention and covering the converter with a real test is the difference between "the answer is right" and "the answer looks right".

Current state

The converter is a single 230-line main program with no module boundary and no test. app/efit_to_boozer.f90 reads six positional integers plus psimax from efit_to_boozer.inp (app/efit_to_boozer.f90:24-32), calls flint_for_Boozer → spline_magdata_in_symfluxcoord → per-surface Fourier sums, and writes fromefit.bc, fromefit_neo_lhs.bc, fromefit_neo_rhs.bc, box_size_axis.dat, flux_functions.dat inline. It is built at CMakeLists.txt:427-431. Searching test/ for efit_to_boozer or fromefit returns nothing: efit_to_boozer.x is not registered as a ctest.

COCOS is nowhere in this code path. rg -i cocos src/magfie/field_divB0.f90 src/efit_to_boozer/*.f90 app/efit_to_boozer.f90 returns zero hits. field_eq reads the g-file raw and only shifts and rescales:

    rad = rad*1.d2 ! cm
    zet = zet*1.d2 ! cm
    rtf = rtf*1.d2 ! cm
    psi = psi*1.d8
    psib= psib*1.d8
    btf = btf*1.d4
    psi=psi+psib

(src/magfie/field_divB0.f90:194-202) — no sign normalisation, no 2π/exp_Bpol handling.

The convention is instead inferred at runtime from one field evaluation:

! Direction of poloidal field:
  sigma=sign(1.d0,Bz*Bp)

(src/efit_to_boozer/field_line_integration_for_Boozer.f90:108, sampled one step outboard of the axis). Every sign in the three output files is then built out of that scalar (app/efit_to_boozer.f90:96-97, 123-125, 183-194), documented only by the free-text block at app/efit_to_boozer.f90:70-95 (!<=25.07.2023), which explains the lhs/rhs/legacy split in prose. That block is the only statement of the output convention anywhere in the repository.

libneo already has a COCOS classifier — unused here. src/magfie/geqdsk_tools.f90 defines cocos_t (line 16), geqdsk_classify (line 264, deriving sgn_Bpol, sgn_pol, index ∈ {1,3,5,7} plus +10 from the Grad–Shafranov prefactor) and geqdsk_standardise (line 307, forcing cocos%index = 3 at line 340). Its only consumer is the standalone tool src/magfie/standardise_geqdsk.f90 (wired at CMakeLists.txt:402-405). Nothing links the converter to it.

The shipped fixture is not COCOS 3. test/resources/input_efit_file.dat has simag = -1.08122148e-01, sibry = +1.36480480e-01, bcentr = -2.50046372e+00, current = +9.9361e+05, so sgn_dpsi = +1, sgn_Itor = +1 ⇒ sgn_Bpol = +1 — the COCOS 1/5 family, and libneo#116 reports geqdsk_classify returning index 15 on exactly this file. The converter consumes it unstandardised.

Input plumbing is inconsistent. Three divergent copies of efit_to_boozer.inp exist: python/tests/efit_to_boozer.inp (psimax = 1e10), test/python/efit_to_boozer.inp (psimax = 2119063900.d0), doc/efit_to_boozer/RUN/efit_to_boozer.inp. And doc/example_namelistfile_efit_to_boozer documents a &efit_to_boozer_nml namelist that the program does not read.

Existing indirect coverage: python/tests/test_eqdsk_to_boozer_chartmap.py drives the f2py module on the same g-file and asserts |I/G| < 0.5, I_edge ≈ 0.2·Ip, and sign(A_phi) = -sign(torflux) — good physics checks, but on the chartmap route, not the .bc writer, and it pytest.skips whenever _efit_to_boozer is not importable.

Reproduction / evidence

make                          # cmake --preset default && cmake --build --preset default
cd build && ctest -N | grep -i efit
# -> only test_geqdsk_tools; no test exercises efit_to_boozer.x

mkdir -p /tmp/e2b && cd /tmp/e2b
cp <repo>/test/python/efit_to_boozer.inp .
sed -e "s|'../../test/resources/|'<repo>/test/resources/|g" \
    <repo>/python/tests/field_divB0.inp > field_divB0.inp
<repo>/build/efit_to_boozer.x
head -6 fromefit_neo_lhs.bc

Observed: three .bc files are produced and no line of output states which COCOS the input was in or which convention the files are in. Now the discriminating experiment:

<repo>/build/standardise_geqdsk.x <repo>/test/resources/input_efit_file.dat ./g_cocos3.dat
# edit field_divB0.inp to point gfile at ./g_cocos3.dat, rerun efit_to_boozer.x

Expected after the fix: fromefit_neo_lhs.bc from the standardised and the unstandardised g-file agree to round-off in iota, Itor, bmnc, because the converter normalises internally. Observed today: nothing guarantees this, and the psi sign flip performed by geqdsk_standardise (src/magfie/geqdsk_tools.f90:309-320) feeds straight into the sigma = sign(Bz*Bp) heuristic, which is the only thing that keeps iota > 0. This equality is the acceptance test.

Behaviour (BDD)

Scenario: converter reports and enforces the input COCOS
  Given field_divB0.inp points gfile at test/resources/input_efit_file.dat
  When efit_to_boozer.x runs
  Then stdout contains a line "input COCOS index = 15" (from geqdsk_classify)
  And stdout contains "standardised to COCOS 3"
  And fromefit_neo_lhs.bc line 2 reads
      "CC Boozer file format: left handed set (s,theta,phi), input COCOS 15 -> 3"

Scenario: output is invariant under prior standardisation of the input
  Given g_cocos3.dat produced by "standardise_geqdsk.x test/resources/input_efit_file.dat g_cocos3.dat"
  When efit_to_boozer.x is run once on input_efit_file.dat and once on g_cocos3.dat
  Then the converter regression parser gives iota, Itor, Jpol_over_nper, rmnc, zmns, bmnc agreeing
      to relative tolerance 1e-10 between the two fromefit_neo_lhs.bc files

Scenario: lhs file is right-handed-consistent with iota > 0
  Given fromefit_neo_lhs.bc from test/resources/input_efit_file.dat
  When the converter regression parser reads it
  Then iota(s) > 0 for every surface
  And abs(Itor(nsurf)) is within 5 percent of 0.2 * abs(Ip) in G*cm, Ip taken from the g-file CURRENT record

Scenario: harmonic reconstruction reproduces the field on a surface
  Given fromefit_neo_lhs.bc with mpol = 12 and nsurf = 64
  When Bmod is reconstructed from bmnc/bmns at s = 0.5 on 128 theta points
  Then it matches magdata_in_symfluxcoord_ext bmod at the same points to relative 1e-3

Scenario: input file is read as a namelist
  Given an efit_to_boozer.inp containing "&efit_to_boozer_nml nstep=3600 nlabel=64 ntheta=64 nsurfmax=1000 nsurf=64 mpol=12 psimax=1e10 /"
  When efit_to_boozer.x runs
  Then it starts without error and the legacy positional file is still accepted with a deprecation warning on stdout

Implementation plan

1. Create src/efit_to_boozer/efit_to_boozer_config.f90 — new module efit_to_boozer_config with a derived type efit_to_boozer_config_t (nstep, nlabel, ntheta, nsurfmax, nsurf, mpol, psimax) and read_config(cfg, filename). Try &efit_to_boozer_nml first (schema exactly as in doc/example_namelistfile_efit_to_boozer, plus psimax); on iostat /= 0 rewind and fall back to the positional format currently at app/efit_to_boozer.f90:24-32, printing a deprecation warning. Use newunit=, status='old', action='read' per CODING_STANDARD.md. Add the file to the add_library(efit_to_boozer STATIC ...) list in src/efit_to_boozer/CMakeLists.txt:64-69.

2. Create src/efit_to_boozer/efit_to_boozer_convert.f90 — move the entire body of app/efit_to_boozer.f90 into module efit_to_boozer_convert with one public entry run_conversion(cfg, gfile, outdir) and private helpers write_bc_header, write_bc_surface, write_flux_functions, write_box_size_axis. No behaviour change in this step. Add to the same add_library list. This is the "clean" half of the issue and should be its own PR.

3. Reduce app/efit_to_boozer.f90 to a driver: read the config, call run_conversion. Keeps the efit_to_boozer.x target at CMakeLists.txt:427-431 unchanged.

4. Add COCOS normalisation at the front of the pipeline. In src/efit_to_boozer/efit_to_boozer_convert.f90, explicitly call read_field_input once to load field_divB0.inp, save the original gfile, classify and standardise that file with geqdsk_tools, only: geqdsk_t, geqdsk_read, geqdsk_classify, geqdsk_standardise, geqdsk_write, geqdsk_deinit, assign the scratch filename back to gfile, and mark the field reader initialized before the first call to field; otherwise field() will reread field_divB0.inp and silently discard the scratch path. Print the original index, store it for the header writer, and standardise to COCOS 3. src/efit_to_boozer/CMakeLists.txt:70-75 already links magfie, so no new link edge is needed. Depends on libneo#116 — geqdsk_classify/geqdsk_standardise must be trusted on the shipped fixture before this step lands; fix test_geqdsk_tools fails: COCOS index mismatch in standardization #116 first.

5. Emit the convention into the files. In write_bc_header, extend the CC lines currently at app/efit_to_boozer.f90:99-107 to state the input COCOS index, the output COCOS index (3), the handedness, and the Fourier kernel sign (m*theta + n*phi, per doc/ExtraDocuments/Conventions.md). Keep the column layout of the m0b n0b nsurf nper ... record byte-identical; the existing Fortran read_bc_file only handles its restricted stellarator-symmetric case, so full ten-column output must be consumed by an extended or dedicated parser.

6. Add test/source/test_efit_to_boozer.f90 — a program test_efit_to_boozer that calls run_conversion twice (raw g-file and standardise_geqdsk-preprocessed g-file) into two scratch directories and compares the generated numeric records. The current bc_file::read_bc_file intentionally rejects non-stellarator-symmetric rows with ten columns, while this converter writes rmnc/rmns/zmnc/zmns/vmnc/vmns/bmnc/bmns; either extend bc_file first and parse the full representation or use a dedicated test parser. Do not wire the current reader directly into this test. Assert scenarios 2 and 3 above with error stop on failure. Follow the fixture-loading idiom of test/source/test_bc_file.f90. Use a reduced grid (nlabel = 64, ntheta = 64, nsurf = 64, mpol = 12) to keep the runtime under ~30 s.

7. Wire the test in test/CMakeLists.txt: add_executable(test_efit_to_boozer.x source/test_efit_to_boozer.f90), target_link_libraries(test_efit_to_boozer.x efit_to_boozer ${MAIN_LIB} ${MAGFIE_LIB}), add_test(NAME test_efit_to_boozer COMMAND test_efit_to_boozer.x), and a set_tests_properties(... WORKING_DIRECTORY ...) block mirroring lines 293-295 so the relative ../../test/resources/... paths in field_divB0.inp resolve.

8. Add the harmonic-reconstruction check (scenario 4) to python/tests/test_boozer.py, reusing the existing _efit_to_boozer f2py module and a parser that supports the full ten-column representation (extend the current reader if appropriate, otherwise use a dedicated parser).

9. Consolidate the input fixtures: make test/python/efit_to_boozer.inp and python/tests/efit_to_boozer.inp identical (namelist form, psimax = 1e10) or replace both with a single file referenced from the two test dirs, and update doc/efit_to_boozer/RUN/efit_to_boozer.inp and doc/example_namelistfile_efit_to_boozer to match what the code actually reads.

10. Document in doc/ExtraDocuments/Conventions.md (which already has a COCOS section) a subsection stating: input is standardised to COCOS 3; fromefit_neo_lhs.bc / fromefit_neo_rhs.bc handedness; kernel sign; that fromefit.bc is retained for AUG backward compatibility only. Migrate the substance of the !<=25.07.2023 block there and shorten the inline comment to a pointer.

11. Downstream: ci/downstreams gates SIMPLE, NEO-2, MEPHIT, KAMEL and rabe on every libneo PR. Step 4 can change output signs for non-COCOS-3 g-files, so run the gate and check whether NEO-2 or MEPHIT golden records shift; if they do, the sign change is the bug being fixed, and the downstream reference data must be updated in the same PR train.

Suggested split: steps 1-3 + 9 (mechanical cleanup, no numerics change) as PR 1; steps 4, 5, 10 + libneo#116 as PR 2; steps 6-8 as PR 3 — though writing the test first (TDD, per CODING_STANDARD.md) means step 6 should be authored before step 4 lands.

Acceptance criteria

• ctest -R test_efit_to_boozer --output-on-failure passes.
• ctest -R test_geqdsk_tools --output-on-failure passes (libneo#116 resolved, no XFAIL).
• iota, Itor, Jpol_over_nper, rmnc, zmns, bmnc from fromefit_neo_lhs.bc agree to relative 1e-10 between raw and standardise_geqdsk.x-preprocessed test/resources/input_efit_file.dat.
• iota(s) > 0 on all surfaces of fromefit_neo_lhs.bc for the shipped fixture.
• |Itor(nsurf)| within 5 % of 0.2 * |CURRENT| (G·cm) from the g-file.
• efit_to_boozer.x accepts &efit_to_boozer_nml and still accepts the legacy positional file, printing a deprecation warning for the latter.
• Each fromefit*.bc header names its input and output COCOS index and its handedness; the converter regression parser accepts the full ten-column non-stellarator-symmetric rows, while the existing test_bc_file remains green for its four-column synthetic format.
• python -m pytest python/tests/test_eqdsk_to_boozer_chartmap.py python/tests/test_boozer.py -q passes with no new skips.
• .github/workflows/downstream-gate.yml green for SIMPLE, NEO-2, MEPHIT, KAMEL, rabe.

References

• O. Sauter and S. Yu. Medvedev, Tokamak coordinate conventions: COCOS, Computer Physics Communications 184 (2013) 293-302, doi:10.1016/j.cpc.2012.09.010 — the definition of the 1-18 index, the sigma_Bp, sigma_RphiZ, sigma_rho_theta_phi triple and the exp_Bp 2π factor. geqdsk_classify in src/magfie/geqdsk_tools.f90:264 implements exactly this classification.
• W. D. D'haeseleer, W. N. G. Hitchon, J. D. Callen, J. L. Shohet, Flux Coordinates and Magnetic Field Structure, Springer 1991 — Boozer coordinates, the straight-field-line construction and the Jacobian conventions used throughout app/efit_to_boozer.f90; already cited from doc/ExtraDocuments/Conventions.md.
• A. H. Boozer, Plasma equilibrium with rational magnetic surfaces, Physics of Fluids 24 (1981) 1999 — the original definition.
• Comparable implementations: booz_xform (Landreman's standalone rewrite of the STELLOPT BOOZ_XFORM, already a libneo [dev] dependency and consumed by src/field/boozmn_file.f90) fixes its own kernel and handedness explicitly in its output metadata; OMFIT and the OMAS data model carry an explicit COCOS attribute on every equilibrium and convert on read, which is the pattern step 4 follows; the IMAS/ITER equilibrium IDS mandates COCOS 11 for the same reason. GPEC and the ASCOT5 Bfield readers likewise require the caller to state the input convention rather than inferring it.

Scope and risk

Three PRs (mechanical extraction; COCOS normalisation + docs; tests), sequenced after libneo#116. The risk is concentrated in step 4: introducing standardisation can legitimately flip signs in .bc files for any g-file that is not already COCOS 3, which changes NEO-2 / MEPHIT results — that is the point of the change, but it must be caught by the downstream gate and reflected in downstream reference data rather than discovered later.

Related open issues: #363 ("Automatically check for COCOS 3 in EQDSK input") is a strict subset of step 4 here — implement it as part of this work and close it against the same PR. #116 (test_geqdsk_tools COCOS index mismatch) is a hard blocker for step 4, since it means the classifier is not currently trusted on the very fixture this issue's tests use.