Equilibrium Module

The Equilibrium module provides tools to read, construct, and analyze magnetohydrodynamic (MHD) plasma equilibria. It wraps a variety of equilibrium file formats (EFIT, CHEASE, SOL, LAR, and others), runs the appropriate direct or inverse solver, and returns a processed PlasmaEquilibrium object ready for downstream calculations.

Overview

Key responsibilities of the module:

  • Read equilibrium input files and TOML configuration (see EquilibriumConfig).
  • Provide convenient constructors for analytic / model equilibria (Large
Aspect Ratio, Solev'ev).
  • Build spline representations used throughout the code (1D cubic and
2D bicubic splines).
  • Run the direct or inverse equilibrium solver and post-process results
(global parameters, q-profile finding, separatrix finding, GSE checks).

The module exposes a small public API that covers setup, configuration, and common analyses used by other JPEC components (e.g. DCON, vacuum interfaces).

API Reference

JPEC.Equilibrium.DirectBFieldType
DirectBField

Internal mutable struct to hold B-field components and their derivatives at a point. It is used as a temporary workspace to avoid allocations in tight loops.

source
JPEC.Equilibrium.DirectRunInputType
DirectRunInput(...)

A container struct that bundles all necessary inputs for the direct_run function. It is created by the _read_efit function after parsing the raw equilibrium file and preparing the initial splines.

Fields:

  • equil_input: The original EquilInput object.
  • sq_in # x value: psin # Quantity 1: F = R*Bt [m T] # Quantity 2: mu0 * Pressure (non-negative) [nt^2 / m^2 * mu0 = T^2] # Quantity 3: q-profile # Quantity 4: sqrt(psi_norm)
  • psi_in: # x, y value: R, Z [m] # z value : poloidal flux adjusted to be zero at the boundary [Weber/radian] # 1. ψ(R,Z) = ψ_boundary - ψ(R,Z) # 2. if ψ = ψ * sign(ψ(centerR,centerZ))
  • rmin: Minimum R-coordinate of the computational grid [m].
  • rmax: Maximum R-coordinate of the computational grid [m].
  • zmin: Minimum Z-coordinate of the computational grid [m].
  • zmax: Maximum Z-coordinate of the computational grid [m].
  • psio: The total flux difference abs(ψ_axis - ψ_boundary) [Weber / radian].
source
JPEC.Equilibrium.LargeAspectRatioConfigType
LargeAspectRatioConfig(...)

A mutable struct holding parameters for the Large Aspect Ratio (LAR) plasma equilibrium model.

Fields:

  • lar_r0: The major radius of the plasma [m].
  • lar_a: The minor radius of the plasma [m].
  • beta0: The beta value on axis (normalized pressure).
  • q0: The safety factor on axis.
  • p_pres: The exponent for the pressure profile, defined as p00 * (1 - (r / a)^2)^p_pres.
  • p_sig: The exponent that determines the shape of the current-related function profile.
  • sigma_type: The type of sigma profile, can be "default" or "wesson". If "wesson", the sigma profile is defined as sigma0 * (1 - (r / a)^2)^p_sig.
  • mtau: The number of grid points in the poloidal direction.
  • ma: The number of grid points in the radial direction.
  • zeroth: If set to true, it neglects the Shafranov shift
source
JPEC.Equilibrium.PlasmaEquilibriumType
PlasmaEquilibrium(...)

The final, self-contained result of the equilibrium reconstruction. This object provides a complete representation of the processed plasma equilibrium in flux coordinates.

Fields:

  • equil_input: The original EquilInput object used for the reconstruction.
  • sq: The final 1D profile spline (CubicSpline{Float64}). # x value: normalized psi # Quantity 1: Toroidal Field Function * 2π, F * 2π (where F = R * B_toroidal) # Quantity 2: Pressure * μ₀, P * μ₀. # Quantity 3: dVdpsi # Quantity 4: q
  • rzphi: The final 2D flux-coordinate mapping spline (BicubicSpline). # x value: normlized psi # y value: SFL poloidal angle [0,1] # Quantity 1: rcoord² = (R - ro)² + (Z - zo)² # Quantity 2: Offset between the geometric poloidal angle (η) and the new angle (θnew) `η / (2π) - θ_new # Quantity 3: ν in ϕ=2πζ+ν(ψ,θ) # Quantity 4: Jacobian.
  • eqfun: A 2D spline storing local physics and geometric quantities that vary across the flux surfaces. # These are pre-calculated for efficient use in subsequent stability and transport codes. # x value: Normalized poloidal flux, ψnorm ∈ [0, 1]. # y value: SFL poloidal angle, θnew ∈ [0, 1]. # Quantity 1: Total magnetic field strength, B [T] # Quantity 2: (e₁⋅e₂ + q⋅e₃⋅e₁) / (J⋅B²). # Quantity 3: (e₂⋅e₃ + q⋅e₃⋅e₃) / (J⋅B²).
  • ro: R-coordinate of the magnetic axis [m].
  • zo: Z-coordinate of the magnetic axis [m].
  • psio: Total flux difference |Ψ_axis - Ψ_boundary| [Weber / radian].
source
JPEC.Equilibrium.SolevevConfigType
SolevevConfig(...)

A mutable struct holding parameters for the Solev'ev (SOL) plasma equilibrium model.

Fields:

  • mr: number of radial grid zones
  • mz: number of axial grid zones
  • ma: number of flux grid zones
  • e: elongation
  • a: minor radius
  • r0: major radius
  • q0: safety factor at the o-point
  • p0fac: scale on-axis pressure (P-> P+P0*p0fac. beta changes. Phi,q constant)
  • b0fac: scale toroidal field at constant beta (sPhi,sf,s^2*P. bt changes. Shape,beta constant)
  • f0fac: scale toroidal field at constant pressure (s*f. beta,q changes. Phi,p,bp constant)
source
JPEC.Equilibrium._read_1d_gfile_formatMethod

read1dgfileformat(linesblock, numvalues)

Internal helper function to parse Fortran-style fixed-width numerical blocks from a vector of strings.

Arguments:

  • lines_block: A Vector{String} containing the lines to parse.
  • num_values: The total number of Float64 values to read from the block.

Returns:

  • A Vector{Float64} containing the parsed values.
source
JPEC.Equilibrium.direct_fl_der!Method
direct_fl_der!(dy, y, params, eta)

The derivative function for the field-line integration ODE. This is passed to the DifferentialEquations.jl solver.

Arguments:

  • dy: The derivative vector (output, modified in-place).
  • y: The state vector [∫(dl/Bp), rfac, ∫(dl/(R²Bp)), ∫(jac*dl/Bp)].
  • params: A FieldLineDerivParams struct with all necessary parameters.
  • eta: The independent variable (geometric angle η).
source
JPEC.Equilibrium.direct_fl_intMethod
direct_fl_int(ipsi, psifac, raw_profile, ro, zo, rs2)

Performs the field-line integration for a single flux surface.

Arguments:

  • ipsi: The index of the current flux surface.
  • psifac: The normalized psi value for the surface (ψ_norm).
  • raw_profile: The DirectRunInput object containing splines and parameters.
  • ro, zo: Coordinates of the magnetic axis [m].
  • rs2: R-coordinate of the outboard separatrix [m].

Returns:

  • y_out: A matrix containing the integrated quantities vs. the geometric angle η.
  • bf_start: A DirectBField object with values at the integration start point.
source
JPEC.Equilibrium.direct_get_bfield!Method
direct_get_bfield!(bf_out, r, z, psi_in, sq_in, psio; derivs=0)

Calculates the magnetic field and its derivatives at a given (R,Z) point. The results are stored in-place in the bf_out object.

Arguments:

  • bf_out: A mutable DirectBField struct to store the results.
  • r: The R-coordinate [m].
  • z: The Z-coordinate [m].
  • psi_in: The 2D bicubic spline for poloidal flux ψ(R,Z).
  • sq_in: The 1D cubic spline for profiles F(ψ_norm) and P(ψ_norm).
  • psio: The total flux difference |ψ_axis - ψ_boundary|.

Keyword Arguments:

  • derivs: An integer specifying the derivative level to compute.
    • 0: Calculates ψ, F, and P.
    • 1: Also calculates 1st derivatives of ψ and B field components.
    • 2: Also calculates 2nd derivatives of ψ and 1st derivatives of B.

Returns:

  • nothing. The bf_out object is modified in-place.
source
JPEC.Equilibrium.direct_positionMethod
direct_position(psi_in, sq_in, psio, ro_guess, zo_guess, rmin, rmax)

Finds the key geometric locations of the equilibrium: the magnetic axis (O-point) and the inboard/outboard separatrix crossings on the midplane.

Arguments:

  • psi_in: The 2D ψ(R,Z) spline.
  • sq_in: The 1D profile spline.
  • psio: Initial flux difference.
  • ro_guess, zo_guess: Initial guess for the magnetic axis location [m].
  • rmin, rmax: Radial bounds of the computational domain [m].

Returns:

  • ro: R-coordinate of the magnetic axis [m].
  • zo: Z-coordinate of the magnetic axis [m].
  • rs1: R-coordinate of the inboard separatrix crossing [m].
  • rs2: R-coordinate of the outboard separatrix crossing [m].
  • psi_in_new : returns psi_in renormalized by * psio/psi(ro,zo)
source
JPEC.Equilibrium.direct_refineMethod
direct_refine(rfac, eta, psi0, params)

Refines the radial distance rfac at a given angle eta to ensure the point lies exactly on the target flux surface psi0.

Arguments:

  • rfac: The current guess for the radial distance from the magnetic axis.
  • eta: The geometric poloidal angle.
  • psi0: The target ψ value for the flux surface.
  • params: A FieldLineDerivParams struct.

Returns:

  • The refined rfac value.
source
JPEC.Equilibrium.equilibrium_solverMethod
equilibrium_solver(raw_profile)

The main driver for the direct equilibrium reconstruction. It orchestrates the entire process from finding the magnetic axis to integrating along field lines and constructing the final coordinate and physics quantity splines.

Arguments:

  • raw_profile: A DirectRunInput object containing the initial splines (psi_in, sq_in) and run parameters (equil_input).

Returns:

  • A PlasmaEquilibrium object containing the final, processed equilibrium data, including the profile spline (sq), the coordinate mapping spline (rzphi), and the physics quantity spline (eqfun).
source
JPEC.Equilibrium.inverse_extrapMethod
inverse_extrap(xx::Matrix{Float64}, ff::Matrix{Float64}, x::Float64) -> Vector{Float64}

Performs component-wise Lagrange extrapolation for a vector-valued function.

Arguments:

  • xx: A (m × n) matrix where each row contains the x-values for each component.
  • ff: A (m × n) matrix where each row contains function values at the corresponding xx.
  • x: A scalar Float64 value at which to extrapolate.

Returns:

  • A vector of length n representing the extrapolated function values at x.
source
JPEC.Equilibrium.lar_init_conditionsMethod
lar_init_conditions(rmin, sigma_type, params)

Initializes the starting radius and state vector for solving the LAR ODE system. Also evaluates the initial derivative using the analytic model.

Arguments:

  • rmin: Normalized starting radius (as a fraction of lar_a).
  • lar_input: A LargeAspectRatioConfig object containing equilibrium parameters.

Returns:

  • r: Physical radius corresponding to rmin * lar_a.
  • y: Initial state vector of length 5.
source
JPEC.Equilibrium.read_chease2Method
_read_chease2(equil_in)

Parses a chease2 file, creates initial 1D and 2D splines, finds magnetic axis, and bundles them into a InverseRunInput object.

Arguments:

  • equil_in: The EquilInput object containing the filename and parameters.

Returns:

  • A InverseRunInput object ready for the inverse solver.
source
JPEC.Equilibrium.read_efitMethod
_read_efit(equil_in)

Parses an EFIT g-file, creates initial 1D and 2D splines, and bundles them into a DirectRunInput object.

Arguments:

  • equil_in: The EquilInput object containing the filename and parameters.

Returns:

  • A DirectRunInput object ready for the direct solver.
source
JPEC.Equilibrium.setup_equilibriumFunction
setup_equilibrium(equil_input::EquilInput)

The main public API for the Equilibrium module. It orchestrates the entire process of reading an equilibrium file, running the appropriate solver, and returning the final processed PlasmaEquilibrium object.

Arguments:

  • equil_input: An EquilInput object containing all necessary setup parameters.

Returns:

  • A PlasmaEquilibrium object containing the final result.
source
JPEC.Equilibrium.sol_runMethod

This is a Julia version of the Fortran code in sol.f, implementing Soloviev's analytical equilibrium.

Arguments:

  • mr: Number of radial grid zones
  • mz: Number of axial grid zones
  • ma: Number of flux grid zones
  • e: Elongation
  • a: Minor radius
  • r0: Major radius
  • q0: Safety factor at the o-point
  • p0fac: Scales on axis pressure (s*P. beta changes. Phi,q constant)
  • b0fac: Scales on toroidal field (sPhi,sf,s^2*P. bt changes. Shape,beta constant)
  • f0fac: Scales on toroidal field (s*f. bt,q changes. Phi,p,bp constant)

Returns:

  • plasma_eq: PlasmaEquilibrium object
source

Important types

  • EquilibriumConfig — top-level configuration container parsed from a
TOML file (outer constructor `EquilibriumConfig(path::String)` is
provided).
  • EquilibriumControl — low-level control parameters (grid, jacobian
type, tolerances, etc.).
  • EquilibriumOutput — options controlling what output is written.
  • PlasmaEquilibrium — the runtime structure containing spline fields,
geometry, profiles, and computed diagnostics (q-profile, separatrix,
etc.).
  • LargeAspectRatioConfig, SolevevConfig — convenience structs to
construct analytic/model equilibria when using `eq_type = "lar"` or
`eq_type = "sol"`.

Key functions

  • setup_equilibrium(path::String = "equil.toml")
— main entry point that reads configuration, builds the equilibrium,
runs the solver, and returns a `PlasmaEquilibrium` instance.
  • equilibrium_separatrix_find!(pe::PlasmaEquilibrium) — locate
separatrix and related boundary geometry in-place.
  • equilibrium_global_parameters!(pe::PlasmaEquilibrium) — populate
common global parameters (major radius, magnetic axis, volumes, etc.).
  • equilibrium_qfind!(pe::PlasmaEquilibrium) — compute safety factor
(q) information across the grid.
  • equilibrium_gse!(pe::PlasmaEquilibrium) — diagnostics on the
Grad–Shafranov solution.

Example usage

Basic example: read a TOML config and build an equilibrium

using JPEC

# Build from a TOML file (searches relative paths if needed)
pe = JPEC.Equilibrium.setup_equilibrium("docs/examples/dcon.toml")

println("Magnetic axis: ", pe.params.r0, ", ", pe.params.z0)
println("q(0) = ", pe.params.q0)

# Find separatrix (in-place) and inspect results
JPEC.Equilibrium.equilibrium_separatrix_find!(pe)
println("rsep = ", pe.params.rsep)

Analytic / testing example: construct a large-aspect-ratio model

using JPEC

# Create a LAR config from a small TOML fragment or file
larcfg = JPEC.Equilibrium.LargeAspectRatioConfig(lar_r0=10.0, lar_a=1.0, beta0=1e-3)
pe = JPEC.Equilibrium.setup_equilibrium(JPEC.Equilibrium.EquilibriumConfig(control=Dict("eq_filename"=>"unused","eq_type"=>"lar")), larcfg)

println("Built LAR equilibrium with a = ", lorcfg.lar_a)

Notes:

  • The EquilibriumConfig(path::String) constructor parses TOML and
expects `[EQUIL_CONTROL]` to contain at minimum `eq_filename` and
`eq_type` fields. Paths that are not absolute are resolved relative
to the TOML file location.
  • When eq_type == "inverse_testing" a small example inverse run is
constructed (useful in tests and examples).

Notes and Caveats

  • Many routines rely on spline representations; the Splines module is
used heavily and should be initialized where appropriate.
  • The Equilibrium module contains several reader routines for external
formats (EFIT/CHEASE) and also interfaces to older Fortran helpers —
ensure required data files are present for those backends.
  • For programmatic usage, prefer constructing EquilibriumConfig from a
TOML file to ensure all path resolution and defaults are handled.

See also

  • docs/src/splines.md — spline helpers used by equilibrium routines
  • docs/src/vacuum.md — coupling between equilibrium and vacuum solvers

```