Tearing - MERGE - Merge feature/tearing-growthrates into FKR-width

Brings FKR-width up to date after 362 commits on feature/tearing-growthrates,
including two develop merges and the InnerLayer relocation.

Two conflicts, both from that relocation. src/Tearing/InnerLayer/SLAYER/ now
lives at src/InnerLayer/SLAYER/ (InnerLayer moved to the top level so it loads
before ForceFreeStates, which calls it for the matched-Delta' Galerkin solve):

  - LayerThickness.jl: the LayerWidths constructor call conflicted because this
    branch added the delta_FKR / delta_visco arguments while the other side
    reindented. Kept both this branch's extra arguments and the target branch's
    formatter indentation. The struct itself merged cleanly and carries all nine
    fields.
  - SetResistiveWidthPsihigh.jl was added inside the renamed directory, so git
    flagged a file-location conflict; accepted its relocation to
    src/InnerLayer/SLAYER/.

Also stripped trailing whitespace from the two touched files, which would fail
the pre-commit hygiene hook.

Verified: the package precompiles clean, LayerWidths carries
(ising, m, n, dels_db, delta_s, delta_s_m, d_beta, delta_FKR, delta_visco), and
runtests_slayer_{params,riccati,inputs,runner}.jl and runtests_innerlayer.jl all
pass.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

Author: d-burg

.claude/agent-memory/fortran-physics-reviewer/MEMORY.md

@@ -0,0 +1,7 @@
+# Agent Memory Index
+
+- [Fortran-Julia Correspondence Map](fortran_correspondence_map.md) — Maps Fortran GPEC source files to Julia JPEC modules
+- [Known Physics Issues in PerturbedEquilibrium](known_issues_perteq.md) — xss_mn not implemented, xi_modes placeholders
+- [KineticForces (NTV) Audit Checklist](kinetic_ntv_map.md) — pentrc->KineticForces map, Logan 2015 matrices, what to verify
+- [InnerLayer (Resistive) Audit Checklist](resistive_layer_map.md) — rmatch->InnerLayer map, GGJ Wasow basis / Δ′, what to verify
+- [Galerkin Δ′ Assembly Map](galerkin_assembly_map.md) — gal.f<->GalerkinAssembly.jl; resonant sign chain verified; PASS

.claude/agent-memory/fortran-physics-reviewer/feedback_annotations.md

@@ -0,0 +1,13 @@
+---
+name: annotation_style
+description: How to cite Fortran code in annotations — use papers and submodule names, not line numbers
+type: feedback
+---
+
+Do NOT cite Fortran source line numbers in code annotations (e.g., "gpeq.f:100"). The Fortran codebase is actively developed and line numbers change. Instead, cite:
+- Paper equations (e.g., "Park et al. 2007, Eq. 8")
+- Fortran subroutine/function names (e.g., "matches Fortran `gpeq_contra`")
+
+**Why:** Fortran GPEC is still under active development; line-number references rot quickly.
+
+**How to apply:** When the fortran-physics-reviewer suggests line-number annotations, convert them to subroutine-name or paper-equation references before writing to code.

.claude/agent-memory/fortran-physics-reviewer/fortran_correspondence_map.md

@@ -0,0 +1,40 @@
+---
+name: Fortran-Julia File Correspondence Map
+description: Mapping between Fortran GPEC source files at ~/Code/gpec and Julia JPEC modules for physics review
+type: reference
+---
+
+## Coil/ForcingTerms
+- `~/Code/gpec/coil/coil.F` (coil_read) -> `src/ForcingTerms/CoilGeometry.jl` (apply_transforms)
+- `~/Code/gpec/coil/field.F` (field_bs_psi) -> `src/ForcingTerms/BiotSavart.jl` + `src/ForcingTerms/CoilFourier.jl`
+
+## PerturbedEquilibrium
+- `~/Code/gpec/gpec/gpeq.f` (gpeq_sol, gpeq_contra, gpeq_surface, gpeq_normal) -> `src/PerturbedEquilibrium/FieldReconstruction.jl` + `src/PerturbedEquilibrium/ResponseMatrices.jl`
+- `~/Code/gpec/gpec/gpresp.f` (gpresp_pinduct, gpresp_sinduct, gpresp_permeab) -> `src/PerturbedEquilibrium/ResponseMatrices.jl`
+- `~/Code/gpec/gpec/gpout.f` (gpout_singcoup, gpout_xbnormal) -> `src/PerturbedEquilibrium/SingularCoupling.jl` + `src/PerturbedEquilibrium/FieldReconstruction.jl`
+- `~/Code/gpec/gpec/gpvacuum.f` (gpvacuum_flxsurf) -> `src/PerturbedEquilibrium/SingularCoupling.jl` (compute_surface_inductance_from_greens)
+
+## KineticForces (NTV, Fortran `pentrc/`) — see kinetic_ntv_map.md for the audit checklist
+- `~/Code/gpec/pentrc/torque.F90` -> `src/KineticForces/Torque.jl` (tpsi! single-surface torque)
+- `~/Code/gpec/pentrc/pitch.f90` -> `src/KineticForces/PitchIntegration.jl` (pitch-angle lambda integration, bounce-averaged operator assembly)
+- `~/Code/gpec/pentrc/energy.f90` -> `src/KineticForces/EnergyIntegration.jl` (energy-space quadrature integrand)
+- `~/Code/gpec/pentrc/pentrc.F90` (orchestration) -> `src/KineticForces/{KineticForces.jl,Compute.jl}`
+- `~/Code/gpec/pentrc/inputs.f90` + `params.f90` -> `src/KineticForces/KineticForcesStructs.jl`
+- `~/Code/gpec/pentrc/dcon_interface.f` -> `src/KineticForces/CalculatedKineticMatrices.jl` (bridge to ForceFreeStates kinetic stability) + `src/ForceFreeStates/{Kinetic.jl,FixedKineticMatrices.jl}`
+- bounce averaging of the 6 perturbed-action matrices -> `src/KineticForces/BounceAveraging.jl`
+
+## InnerLayer (resistive matched-asymptotics, Fortran `rmatch/`) — see resistive_layer_map.md for the audit checklist
+- `~/Code/gpec/rmatch/deltac.f` (Galerkin solver) -> `src/InnerLayer/GGJ/Galerkin.jl`
+- `~/Code/gpec/rmatch/deltar.f` (shooting solver) -> `src/InnerLayer/GGJ/Shooting.jl`
+- `~/Code/gpec/rmatch/inps.f` + `inpso.f` (Wasow asymptotic basis: T,J,P,B,Q,C,D,Y,Z,U matrices) -> `src/InnerLayer/GGJ/InnerAsymptotics.jl`
+- `~/Code/gpec/rmatch/{inner.f,match.f,msing.f,gamma.f}` -> `src/InnerLayer/GGJ/` (parameters, matching data, special functions)
+- SLAYER drift-MHD two-fluid solver: `src/InnerLayer/SLAYER/` is a placeholder pending implementation
+
+## Key Fortran Conventions
+- `sq%f(1)` = F (toroidal field function = R*B_tor), `sq%f(4)` = q (safety factor)
+- `sq%f1(4)` = dq/dpsi
+- `chi1 = twopi*psio` (poloidal flux normalization)
+- `singfac = mfac - nn*q` (resonance factor)
+- `ifac = (0,1)` (imaginary unit)
+- `wegt=0`: no weighting; `wegt=1`: J*|grad psi| weighting (default for bmn)
+- Vacuum Green's functions use reversed theta in Fortran (rtheta = mthsurf - itheta)

.claude/agent-memory/fortran-physics-reviewer/galerkin_assembly_map.md

@@ -0,0 +1,83 @@
+---
+name: galerkin-assembly-map
+description: RDCON singular Galerkin Δ′ solver — gal.f ↔ GalerkinAssembly.jl element-level weak-form correspondence and audit notes
+metadata:
+  type: project
+---
+
+# RDCON Galerkin Δ′ assembly: gal.f ↔ GalerkinAssembly.jl
+
+File: `src/ForceFreeStates/Galerkin/GalerkinAssembly.jl` ↔ `GPEC_fortran/rdcon/gal.f`.
+On branch `feature/rdcon-galerkin`. Reviewed 2026-06-16.
+
+## Correspondence (verified, not just docstring-trusted)
+- gal_hermite (222-262) ↔ gal_hermite: pb/qb EXACT.
+- gal_get_fkg (1013-1071) ↔ gal_get_fkg: F=Fbar·sf[i]·sf[j], K=K·sf[i], G=G. singfac=mlow-nn*q+(0..mpert-1) = DIRECT (m-nq), not reciprocal. Correct.
+- gal_gauss_quad (1079-1138) ↔ gal_gauss_quad!: F·qq + K·qp + conj(K[jpert,ipert])·pq + G·pp. swap pb(2)/pb(3) (0-based) = Julia pb[3]/pb[4]. Correct.
+- gal_extension (1146-1322) ↔ gal_extension!: emat/ediag/rhs/erhs all match incl surface terms 1277-1300 and w1/w2.
+- gal_resonant! ↔ gal_lsode_int(780-900)+gal_lsode_der(908-954)+gal_make_arrays signs(1366-1369). QuadGK replaces LSODE.
+- gal_assemble_mat (694-772) ↔ gal_assemble_mat!: LU offset=kl+ku+1, chol offset=1, conj(e) Hermitian placement. Correct.
+- gal_assemble_rhs (642-686) ↔ gal_assemble_rhs!: isol++ at ext2-left / res-right. Correct.
+- gal_set_boundary (1554-1620) ↔ gal_set_boundary!: axis idx, edge 3-way (rpec/vac/fixed). DOF index 3(0-based)=Julia 4.
+
+## Resonant sign chain (the trap — verified correct)
+Fortran: LSODE integrates x0→x_lsode giving raw u; then line 882 `IF right: u=-u`; then make_arrays:
+erhs=-u_res(1), ediag=+u_res(2), rhs=-u_hermite1, emat=+u_hermite2.
+Net with s=+1(left)/-1(right): erhs=-s·res1, ediag=s·res2, rhs=-s·hbig, emat=s·hsmall.
+Julia sgn_u = s (=-1 right/+1 left): erhs=-sgn_u·raw1, ediag=sgn_u·raw2, rhs=-sgn_u·hbig, emat=sgn_u·hsmall. MATCHES.
+QuadGK limits x0→x_lsode preserve LSODE orientation so raw == pre-negation u. Correct.
+du_res conj-weight uses the SMALL solution (isol(2)) in both. Correct.
+
+## Only open item (MINOR, out-of-file)
+sing_get_ua_gal/dua_gal two-sided z<0 → sig=-1 series at -z with derivative ×(-1) is a REIMPLEMENTATION
+of sing.f's internal singular-surface handling, not a line port. The ×(-1) chain is internally consistent
+but worth a dedicated check against sing.f sing_get_ua/dua if Δ′ ever disagrees with Fortran. Lives outside
+GalerkinAssembly.jl.
+
+## RPEC matching (GalerkinMatch.jl <-> rmatch/match.f match_rpec) — reviewed 2026-06-20
+- match block (mat idx1..idx4, -delta1/delta2 signs) + rmat=-TRANSPOSE(delta coil block) + cout/cin split: EXACT port. PASS.
+- matched outer soln Σcout(isol)·sols + sols(:,:,csol): EXACT (match_output_solution L32-41).
+- resist.f E/F/H/M/G/K/taua/taur/v1: EXACT (Resist.jl). eta/rho caller-supplied = deliberate documented deviation.
+- gal quadrature: jacobi_alloc "gll" = Gauss-Lobatto = Julia gausslobatto. nq1 lower=|nn·q1|, upper=|q1|: match.
+
+## SUSPECT (should-fix) — RPEC forced eigenvalue γ
+GalerkinMatch.jl:66 uses γ = 2π·im·n·gal_rotation. Fortran match_rpec passes `initguess` (namelist,
+default 0.0) DIRECTLY to deltac_run — NOT i·n·rotation. The `guess+ifac*ntor*rotation` Doppler form
+lives only in match_delta (Newton eigenvalue search), NOT the RPEC path. Two issues:
+  (1) extra 2π: Fortran eig=ifac*ntor*rotation has NO 2π (rotation is ω[rad/s], eig is 1/time, q0=x0/taua).
+  (2) RPEC path in Fortran uses initguess (the user forced eigenvalue), not a per-surface i·n·rotation shift.
+Confirm intended convention vs Fortran before trusting resistive Δ(Q); ideal_flag path (γ unused) is unaffected.
+RESOLUTION 2026-06-20: author confirms the Hz convention is DELIBERATE (gal_rotation is f[Hz], so γ=2πi·n·f);
+the overstated comments were corrected. Δ(Q) now guarded by test/runtests_innerlayer.jl (Julia self-pin on
+glasser_wang_2020_eq55). Remaining follow-up (mpharr): upgrade that guard to a Fortran deltac_run / paper
+cross-check to fully close the convention question.
+
+## Verdict: PASS WITH REQUIRED ANNOTATIONS for outer solve; RPEC γ convention is SHOULD-FIX (verify vs match.f).
+
+
+# Second pass
+
+
+## PR #266 "Galerkin Integrator" full review (2026-06-25) — PASS
+Exhaustive line-by-line vs gal.f / match.f / deltac.f / resist.f / sing.f. All physics-bearing code faithful.
+- RESOLVED (was SUSPECT): forced eig γ=2πi·n·rotation IS correct. match.f match_rpec line:
+  `rpec_eigenvalues(ising)=2*pi*rotation(ising)*REAL(ntor,r8)*ifac` — the 2π IS in Fortran. Julia matches exactly.
+  My earlier "no 2π in Fortran" note was wrong (confused match_delta Newton path with match_rpec). Closed.
+- Verified EXACT: gal_get_fkg (F=QF̄Q,K=QK̄,G=Ḡ,singfac=m-nq), gauss_quad weak form (F qq+K qp+K† pq+G pp),
+  gal_extension, gal_resonant (mpert-fastest hermite layout; Fortran's apparent (0:np,mpert,2) reshape is a
+  no-op relabel — der writes mpert-fastest, caller reads mpert-fastest), assemble_mat/rhs, set_boundary,
+  gal_get_solution (restore_uh/us/ul), PEST-3 A/B/Γ/Δ ± combos, Δ′ extraction, match_rpec matrix signs,
+  resist_eval E/F/G/H/K/M/taua/taur, mercier_di=E+F+H-1/4, rescale=sfac^(2p1/3)·v1^(2p1) + delta(1)↔(2) swap
+  (swap applied in BOTH solve_inner and solve_inner_profile), sing_matvec, sing_get_dua, ξ_s=-A⁻¹(Bξ'+Cξ).
+- Documented deviations (all fine): η/ρ caller-supplied (no Spitzer); gal_gamma default 5/3 = resist.f;
+  QuadGK replaces LSODE; sing_get_dua omits Fortran `sig` (moved to sing_get_dua_gal left-side ×(-1)).
+- InnerAsymptotics.jl + Shooting.jl diffs = PURE JuliaFormatter whitespace, zero numeric change.
+- Reference.jl G=8.950e1,H=1.292e-2,K=2.332e2 = corrected GW2020 Eq.55 (prior transcription bug fixed).
+
+## OPEN ITEM (latent, n≥2 only): gal_make_grid resonant-cell width asymmetry
+gal.f gal_make_grid genuinely uses nq1=ABS(nn*sing(ising)%q1) on the LOWER bound but nq1=ABS(sing(ising+1)%q1)
+(NO nn) on the UPPER bound. Julia GalerkinGrid.jl faithfully reproduces. Natural layer scale is 1/(n q1) on BOTH
+sides (singfac'≈-n q1), so the upper/left side of each interior surface gets cells ~n× too wide for n≥2.
+No effect at n=1 (RDCON resistive target). galerkin_solve allows single n≥2 (blocks only multi-n npert>1), so it
+IS reachable. Recommend: keep faithful but warn, or fix both sides to nn*q1 (deviation, more correct) gated by
+validation. Suspected upstream Fortran oversight, not a Julia port error.

.claude/agent-memory/fortran-physics-reviewer/kinetic_ntv_map.md

@@ -0,0 +1,33 @@
+---
+name: KineticForces (NTV) Audit Checklist
+description: Curated map and physics-audit checklist for the KineticForces module (NTV, formerly PENTRC); what to read and what to verify
+type: reference
+---
+
+The KineticForces module computes neoclassical toroidal viscosity (NTV) torque and the
+kinetic energy/matrices from trapped-particle nonambipolar transport. Fortran counterpart is
+`~/Code/gpec/pentrc/` (file mapping in fortran_correspondence_map.md).
+
+## Governing theory
+- Logan & Park (2013) PoP 20, 122507: diamagnetic frequency (Eq. 7), energy integrand (Eq. 8), torque normalization Im(T) = 2n·δW_k (Eq. 19).
+- Logan (2015) PhD Thesis, Ch. 7: the six kinetic matrices A_k, B_k, C_k, D_k, E_k, H_k (Eqs 7.30–7.35) as energy-space integrals of perturbed action operators W_X, W_Y, W_Z; resonance splitting/suppression F_h = (Q−P†)F̄(Q−P)+… (Eq. 7.46). App. C: DCON matrix form of the perturbed action. App. D: numerical treatment of integrable singularities in bounce averages.
+- Park et al. (2009) PRL 102, 065002: trapped-particle nonambipolar transport foundation.
+
+## Key Julia files and roles
+- `Torque.jl` — `tpsi!()`: single-surface torque; poloidal grid sampling; diamagnetic frequencies (Eq. 7).
+- `EnergyIntegration.jl` — energy-space quadrature integrand (Eq. 8); Maxwellian vs JKP distribution options.
+- `PitchIntegration.jl` — pitch-angle (λ) integration; bounce-averaged operator assembly.
+- `BounceAveraging.jl` — bounce-averaging of the perturbed-action matrices (Eqs 7.30–7.35); packs 3 Hermitian + 3 full blocks.
+- `Compute.jl` — `compute_torque_all_methods!()` orchestration over all ψ surfaces.
+- `CalculatedKineticMatrices.jl` — bridge to ForceFreeStates kinetic stability (`Kinetic.jl`, `FixedKineticMatrices.jl`).
+- `KineticForcesStructs.jl` — `KineticForcesControl` (TOML params) and `KineticForcesInternal` (state, hints, mode indexing).
+- `Output.jl` — HDF5 writer.
+
+## What to verify in a review
+- **Matrix structure**: A symmetric; the Hermitian blocks actually Hermitian; the 3-Hermitian + 3-full packing in BounceAveraging matches Logan 2015 Eqs 7.30–7.35.
+- **Resonance handling**: Sokhotski–Plemelj pole shift / singular-denominator gating near rational surfaces (App. D). Check the singular-eps threshold is applied where the denominator vanishes, not globally.
+- **Quadrature tolerances**: outer ψ (atol_psi/rtol_psi) and inner pitch/energy (atol_xlmda/rtol_xlmda) are passed through, not hard-coded to lazy defaults.
+- **Distribution/collision options**: energy integrand (Eq. 8) honors the configured collision operator (harmonic/Krook/zero) and distribution (Maxwellian/JKP/CGL) — not silently fixed to one.
+- **Normalization**: torque Im(T)=2n·δW_k (Eq. 19) and diamagnetic-frequency sign/factor conventions match Fortran `torque.F90`.
+- **Mode indexing**: m, n ranges and block packing over n stay consistent with ForceFreeStates.
+- **Method variants**: FGAR/TGAR/PGAR/RLAR/CLAR/FCGL/TMM/WMM each present, not stubbed to a single fallback.

.claude/agent-memory/fortran-physics-reviewer/known_issues_perteq.md

@@ -0,0 +1,17 @@
+---
+name: Known Physics Issues in PerturbedEquilibrium
+description: Documented physics issues found in perturbed equilibrium module during PR #196 review
+type: project
+---
+
+## xss_mn (toroidal displacement) not implemented
+`FieldReconstruction.jl:316` sets `xss = 0`. This makes b_theta_modes and b_zeta_modes incorrect.
+b_psi_modes is unaffected. The Fortran computes xss_mn from idcon_matrix A/B/C matrices.
+**Why:** Not yet ported from Fortran gpeq_sol (lines 86-96).
+**How to apply:** Any downstream use of b_theta or b_zeta (e.g., tangential displacement, parallel field) requires implementing xss_mn first.
+
+## xi_modes.theta and xi_modes.zeta are zero placeholders
+`FieldReconstruction.jl:136-137`. Fortran computes these via gpeq_contra Jacobian convolution.
+
+## gpeq_contra Jacobian convolution bypassed for xi_n
+Julia skips the Jacobian convolution for xi_psi before computing xi_n. This is mathematically equivalent (J cancels in xwp/(J*delpsi) = xsp/delpsi) but the intermediate quantities differ from Fortran.

.claude/agent-memory/fortran-physics-reviewer/resistive_layer_map.md

@@ -0,0 +1,35 @@
+---
+name: InnerLayer (Resistive Matched-Asymptotics) Audit Checklist
+description: Curated map and physics-audit checklist for the InnerLayer module (GGJ/SLAYER resistive Delta-prime matching); what to read and what to verify
+type: reference
+---
+
+The InnerLayer module performs matched-asymptotic analysis of resistive MHD stability at
+rational surfaces — solving the resistive inner-layer equations and computing the tearing
+stability parameter Δ via parity-projected splitting (Δ_odd, Δ_even). Fortran counterpart is
+`~/Code/gpec/rmatch/` (file mapping in fortran_correspondence_map.md).
+
+## Governing theory
+- Glasser (2016) PoP 23, 072505: computation of resistive instabilities by matched asymptotic expansions; shooting method, Frobenius exponent splitting.
+- Glasser, Wang & Park (2016) PoP 23, 112506: matched-asymptotic matching data (Eqs 34–35, parity projection).
+- Glasser (2020) "Asymptotic solutions and convergence studies of the resistive inner region equations": full Wasow asymptotic basis (T, J, P, B, Q, C, D, Y, Z, U matrices); Wang rescaling X₀^(2√(−D_I)).
+- Glasser (2018) PoP 25, 032501; Wang et al. (2020) PoP 27, 122509: robust Δ′ matrix / asymptotic-matching resistive response.
+
+## Key Julia files and roles
+- `InnerLayer.jl` / `InnerLayerInterface.jl` — pluggable `InnerLayerModel` + `solve_inner()`.
+- `GGJ/GGJ.jl` — Glasser–Greene–Johnson solver selector/exports.
+- `GGJ/GGJParameters.jl` — physical/dimensionless params (E, F, G, H, K); Mercier indices D_I, D_R; Lundquist ratio S = τ_R/τ_A.
+- `GGJ/InnerAsymptotics.jl` — Wasow basis construction (T, J matrices; Lyapunov solve; splitting transform; Y-series Frobenius exponents).
+- `GGJ/Shooting.jl` — backward shoot x_max → 0 (Frobenius basis).
+- `GGJ/Galerkin.jl` — Hermite-cubic finite-element solver with parity BCs.
+- `GGJ/Reference.jl` — benchmark `glasser_wang_2020_eq55`.
+- `SLAYER/Slayer.jl` — placeholder (drift-MHD two-fluid, pending).
+
+## What to verify in a review
+- **Mercier gating**: code must require D_I < 0 (stable Mercier) before proceeding, or error — confirm the assertion exists.
+- **Frobenius origin exponents**: p₁ = √(−D_I); check the stability requirement and the Y-series exponents match Glasser 2020.
+- **Wasow basis convergence**: truncation order kmax vs growth-rate scale Q; verify the basis matrices (T,J,P,B,Q,C,D,Y,Z,U) are built per Glasser 2020, not approximated.
+- **Parity projection**: (Δ_odd, Δ_even) split follows Glasser–Wang–Park 2016 Eqs 34–35.
+- **Rescaling**: physical Δ = matching data × X₀^(2p₁) × v₁^(2p₁) — confirm the rescaling factor is applied.
+- **Two solvers agree**: shooting and Galerkin must produce consistent Δ on `glasser_wang_2020_eq55`; a divergence is a red flag.
+- **Stiff ODE control**: timestep/tolerance handles the exponential-decay modes in the asymptotic regime.