freegsnke.GSstaticsolver module

Module for solving the static forward and inverse GS problems.

Copyright 2025 UKAEA, UKRI-STFC, and The Authors, as per the COPYRIGHT and README files.

This file is part of FreeGSNKE.

FreeGSNKE is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.

FreeGSNKE is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

You should have received a copy of the GNU Lesser General Public License along with FreeGSNKE. If not, see <http://www.gnu.org/licenses/>.

class freegsnke.GSstaticsolver.NKGSsolver(eq, l2_reg=1e-06, collinearity_reg=1e-06, seed=42)

Bases: object

Newton–Krylov solver for the nonlinear static Grad–Shafranov equilibrium problem.

This class solves the nonlinear free-boundary Grad–Shafranov equation:

Δψ = − μ₀ R Jtor(ψ)

written as a nonlinear root-finding problem:

F(ψ) = 0

where:

ψ : plasma poloidal flux Jtor(ψ) : toroidal plasma current density profile μ₀ R : geometric operator coefficients

The solver supports:

• Forward equilibrium solving (fixed coil currents)

• Inverse magnetic configuration optimisation

Solver structure:

The nonlinear system is solved using:

• Newton–Krylov nonlinear iteration

• Finite difference Jacobian approximations

Solution domain:

The spatial grid is fixed at instantiation using the provided equilibrium object.

Notes

The solver assumes:

• Fixed computational domain

• Fixed grid spacing

• Compatible equilibrium geometry

F_function(plasma_psi, tokamak_psi, profiles)

Compute the nonlinear Grad–Shafranov residual written as a root-finding problem.

The Grad–Shafranov equation is solved in residual form:

F(ψ_plasma) = ψ_plasma − G(ψ_boundary, Jtor(ψ_plasma))

where:

G(…) represents the solution of the linearised GS operator with given boundary conditions and toroidal current source.

The equilibrium solution satisfies:

F(ψ_plasma) = 0

Physically, this measures the mismatch between:

• The current plasma flux field

• The flux predicted by solving the linear GS equation with the computed toroidal current density.

Parameters:

• plasma_psi (ndarray) – Flattened plasma poloidal flux vector. Shape = (nx * ny,).

• tokamak_psi (ndarray) –

  Vacuum flux contribution from:

  • Active coils

  • Passive conducting structures

  Same flattened shape as plasma_psi.

• profiles (freegsnke profile object) –

  Plasma profile model used to compute toroidal current density:

  Jtor(ψ).

Returns:

residual – Nonlinear GS residual vector. Root of this function corresponds to GS equilibrium.

Return type:

ndarray

Notes

This formulation is used for:

• Newton–Krylov nonlinear solves

• Picard fixed-point iterations

• Residual diagnostics during equilibrium solving

__init__(eq, l2_reg=1e-06, collinearity_reg=1e-06, seed=42)

Initialise the Grad–Shafranov nonlinear solver.

The constructor prepares all numerical operators required for nonlinear GS solving, including:

• Linear GS multigrid solver

• Green’s function boundary response operator

• Newton–Krylov nonlinear solver backend

• Random generator for Krylov direction perturbations

Parameters:

• eq (FreeGSNKE equilibrium object) –

  Defines the computational domain and geometry.

  Required attributes:

  eq.R : 2D radial grid eq.Z : 2D vertical grid

  The solver grid is fixed to this domain.

• l2_reg (float, optional (default=1e-6)) – Tikhonov regularisation coefficient applied to nonlinear least-squares systems.

• collinearity_reg (float, optional (default=1e-6)) – Additional regularisation penalising collinear search directions in Krylov space. Improves numerical stability of NK iteration.

• seed (int, optional (default=42)) –

  Random number generator seed used for:

  • Krylov perturbation generation

  • Directional exploration in nonlinear solve

self.R, self.Z

Computational grid coordinates.

Type:

ndarray

self.nx, self.ny

Grid dimensions.

Type:

int

self.dRdZ

Differential area element used for integration.

Type:

float

self.linear_GS_solver

Multigrid solver for linearised GS equation.

self.greenfunc

Boundary response Green’s function matrix.

self.nksolver

Newton–Krylov nonlinear solver backend.

Notes

The solver domain cannot be changed after construction. A new solver instance must be created for different grids.

forward_solve(eq, profiles, target_relative_tolerance, max_solving_iterations=100, Picard_handover=0.11, step_size=2.5, scaling_with_n=-1.0, target_relative_unexplained_residual=0.2, max_n_directions=16, max_rel_update_size=0.2, clip=10, force_up_down_symmetric=False, verbose=False, suppress=False)

Solve the forward static Grad–Shafranov (GS) equilibrium problem.

This method computes the plasma poloidal flux ψ_plasma that satisfies the nonlinear static Grad–Shafranov equation:

F(ψ_plasma; ψ_tokamak, profiles) = 0

where:

ψ_tokamak : vacuum flux generated by metal currents profiles : plasma current profile specification (Jtor(ψ))

The total flux is:

ψ_total = ψ_tokamak + ψ_plasma

The nonlinear system is solved using a hybrid scheme:

• Picard iterations when far from convergence

• Newton–Krylov (NK) iterations when close to solution

The final solution is written to:

eq.plasma_psi

Parameters:

• eq (FreeGSNKE equilibrium object) –

  Provides:

  • metal currents (via eq.tokamak)

  • initial guess eq.plasma_psi

  • Green’s functions eq._vgreen

• profiles (FreeGSNKE profile object) – Defines plasma current model Jtor(ψ).

• target_relative_tolerance (float) – Iterations stop when the relative nonlinear residual falls below this value.

• max_solving_iterations (int, optional) – Maximum allowed nonlinear iterations.

• Picard_handover (float) – Relative residual threshold above which Picard iteration is used instead of Newton–Krylov.

• step_size (float) – Proposed NK step magnitude scaling factor.

• scaling_with_n (float) – Additional scaling applied to NK step depending on iteration count: (1 + n_it)**scaling_with_n.

• target_relative_unexplained_residual (float in (0,1)) – NK internal Arnoldi iteration stops once this fraction of the initial residual is explained.

• max_n_directions (int) – Maximum number of Krylov directions explored.

• max_rel_update_size (float) – Maximum allowed relative update to ψ_plasma. Larger updates are rescaled.

• clip (float) – Maximum magnitude of update contribution from each explored Krylov direction.

• force_up_down_symmetric (bool) – If True, enforces up-down symmetry at each step.

• verbose (bool) – Enables detailed iteration logging.

• suppress (bool) – Suppresses all print output.

Returns:

• None

• Side Effects

• ————

• • Updates eq.plasma_psi

• • Updates profiles via Jtor(…)

• • Updates internal convergence diagnostics

Notes

This solver includes:

• Adaptive initialisation

• Residual-based step control

• Trust-region-style update limiting

• Automatic fallback between Picard and NK

• Residual collinearity detection

• Symmetry enforcement (optional)

The algorithm is designed for robustness in highly nonlinear free-boundary equilibrium problems.

freeboundary(plasma_psi, tokamak_psi, profiles)

Apply free-boundary Grad–Shafranov boundary conditions and compute plasma current source terms.

This routine constructs the nonlinear source term and boundary flux conditions required to solve the Grad–Shafranov equation:

Δψ = − μ₀ R Jtor(ψ)

in free-boundary form.

The algorithm performs three main tasks:

1. Compute toroidal current density:

   Jtor(ψ_total)

2. Compute RHS source term for linear GS solve

3. Compute boundary flux contributions using Green’s functions

The total flux used is:

ψ_total = ψ_tokamak + ψ_plasma

Parameters:

• plasma_psi (ndarray) – Flattened plasma poloidal flux vector. Shape = (nx * ny,).

• tokamak_psi (ndarray) –

  Vacuum flux contribution generated by:

  • Active coils

  • Passive structures

• profiles (FreeGSNKE profile object) – Provides plasma current model Jtor(ψ).

Returns:

Updates internal solver state:

self.jtor self.rhs self.psi_boundary

Return type:

None

Notes

• This method is called before solving the linearised GS equation.

• Boundary flux is computed using Green’s function convolution.

get_rel_delta_psi(new_psi, previous_psi, profiles)

Compute the relative change between two flux states in the plasma core.

This function measures the relative difference between two poloidal flux fields, restricted to the plasma core region. It is used as a convergence and stability metric in both forward and inverse solves.

The relative change is defined as:

rel_delta = || (ψ_new − ψ_old) * M ||₂

|| (ψ_new + ψ_old) * M ||₂

where:

ψ_new, ψ_old : flattened flux vectors M : binary core-region mask ||·||₂ : Euclidean (L2) norm

The mask ensures that only the plasma core region contributes to the norm. If no core mask is available, the full domain is used.

This symmetric normalisation (using ψ_new + ψ_old in the denominator) prevents bias when one field is small and provides a scale-invariant measure of change.

Parameters:

• new_psi (ndarray) – Flattened array of the updated poloidal flux values.

• previous_psi (ndarray) – Flattened array of the previous poloidal flux values. Must be the same shape as new_psi.

• profiles (freegsnke profile object) –

  Used to obtain the plasma core mask:

  • profiles.diverted_core_mask if available

  • otherwise the entire computational domain is used.

Returns:

rel_delta_psit – Relative L2 change between the two flux fields restricted to the plasma core region.

Return type:

float

Notes

• If profiles.diverted_core_mask is None or not present, the relative change is computed over the full domain.

• This quantity is dimensionless.

• Used to:

  • Detect convergence of forward solves

  • Limit coil-induced flux changes in inverse solves

  • Control damping and adaptive tolerances

get_rel_delta_psit(delta_current, profiles, vgreen)

Estimate the relative core-region flux perturbation induced by a requested coil current change.

This function computes the magnitude of the poloidal flux perturbation Δψ generated by a proposed change in coil currents, and expresses it relative to the magnitude of the existing tokamak flux.

The induced flux perturbation is computed using the coil Green’s functions:

Δψ(R, Z) = Σ_i ΔI_i G_i(R, Z)

where:

ΔI_i : requested current change in coil i G_i(R, Z) : Green’s function of coil i (R, Z) : spatial grid coordinates

The relative change metric is then defined as:

rel_delta = || Δψ_core ||₂

|| ψ_tokamak ||₂ + ε

where:

ψ_tokamak : current total flux field Δψ_core : Δψ restricted to the plasma core ||·||₂ : Euclidean (L2) norm ε : small regularisation constant (1e-6)

If a diverted core mask is available, only the plasma core region contributes to the numerator. Otherwise, the full domain is used.

Parameters:

• delta_current (ndarray) – One-dimensional array of requested coil current changes. Shape (n_coils,).

• profiles (freegsnke profile object) – Provides the plasma core mask via profiles.diverted_core_mask, if available.

• vgreen (ndarray) – Array of coil Green’s functions. Shape (n_coils, Nx, Ny), where each slice vgreen[i] corresponds to the flux response of coil i per unit current.

Returns:

rel_delta_psit – Dimensionless estimate of the relative magnitude of the induced flux perturbation.

Return type:

float

Notes

• This quantity measures how large the requested current update is

in terms of its effect on the plasma core flux.

• It is typically used to:

  • Control step sizes in inverse solves,

  • Enforce trust-region constraints,

  • Prevent excessively large coil updates,

  • Provide stabilisation in iterative current optimisation.

• The denominator is computed over the full domain

(not masked), matching the current implementation.

inverse_solve(eq, profiles, constrain, target_relative_tolerance, target_relative_psit_update=0.001, max_solving_iterations=100, max_iter_per_update=5, Picard_handover=0.15, step_size=2.5, scaling_with_n=-1.0, target_relative_unexplained_residual=0.3, max_n_directions=16, clip=10, max_rel_update_size=0.15, threshold_val=0.18, l2_reg=1e-09, forward_tolerance_increase=100, max_rel_psit=0.02, damping_factor=0.995, use_full_Jacobian=True, full_jacobian_handover=[1e-05, 0.007], l2_reg_fj=1e-08, force_up_down_symmetric=False, verbose=False, suppress=False)

Inverse solver for static free-boundary Grad–Shafranov equilibria.

This routine solves a coupled inverse problem:

• Unknowns:

  • Control coil currents

  • Plasma poloidal flux (plasma_psi)

• Constraints:

  • Magnetic constraints provided via constrain

  • Grad–Shafranov equation consistency

  • Target plasma profiles (via profiles)

The algorithm alternates between:

1. Updating coil currents to better satisfy magnetic constraints.

2. Performing a forward GS solve to update plasma_psi for the new currents.

Two Jacobian models are available for the current update:

• Green’s function Jacobian (plasma fixed)

• Full Jacobian (plasma response included)

Convergence requires BOTH:

• GS residual < target_relative_tolerance

• Relative tokamak flux update < target_relative_psit_update

The solve is therefore a nonlinear block-iterative scheme coupling:

• Coil optimisation

• Forward equilibrium solve

Upon success:

• eq.plasma_psi contains the converged plasma flux

• eq.tokamak currents contain the converged coil currents

Parameters:

• eq (FreeGSNKE equilibrium object) – Contains the tokamak geometry and coil system. - eq.tokamak holds coil currents. - eq.plasma_psi holds the plasma poloidal flux. This object is updated in-place.

• profiles (FreeGSNKE profile object) –

  Specifies the desired plasma profiles, including:

  • Total plasma current

  • Pressure profile p(psi)

  • F(psi) profile

  These define the toroidal current density Jtor(psi).

• constrain (freegsnke inverse_optimizer object) –

  Defines:

  • Which coils are controllable

  • The magnetic constraints to be satisfied

  • The loss function measuring constraint violation

• target_relative_tolerance (float) – Target relative tolerance for the GS residual. The forward GS problem must satisfy this tolerance for the inverse problem to be considered solved.

• target_relative_psit_update (float, optional (default=1e-3)) – Maximum allowed relative change in tokamak flux (in the plasma core) between successive inverse iterations. Acts as a stability criterion.

• max_solving_iterations (int, optional (default=100)) – Maximum number of outer inverse iterations.

• max_iter_per_update (int, optional (default=5)) – Maximum number of forward GS iterations performed after each coil current update.

• Picard_handover (float, optional (default=0.15)) – Relative residual threshold above which Picard iteration is used instead of Newton–Krylov in the forward solve.

• step_size (float, optional (default=2.5)) – Scaling of exploratory steps in the forward nonlinear solve. Expressed in units of the residual norm.

• scaling_with_n (float, optional (default=-1.0)) –

  Additional scaling factor applied to exploratory steps:

  (1 + iteration_number)**scaling_with_n

• target_relative_unexplained_residual (float, optional (default=0.3)) – Forward solver termination criterion based on how much of the residual can be explained by current search directions.

• max_n_directions (int, optional (default=16)) – Maximum number of search directions used in forward solve.

• clip (float, optional (default=10)) – Maximum update magnitude allowed per direction (in units of exploratory finite-difference step).

• max_rel_update_size (float, optional (default=0.15)) – Maximum allowed relative update to plasma_psi in forward solve. If exceeded, the update is rescaled.

• threshold_val (float, optional (default=0.18)) – Threshold relative GS residual below which the solver switches to more conservative update logic.

• l2_reg (float or 1D numpy array, optional (default=1e-9)) – L2 regularisation factor applied when using the Green’s function Jacobian (simplified current optimisation). If array, must match number of control coils.

• forward_tolerance_increase (float, optional (default=100)) – Controls how tightly the forward problem is solved relative to the size of the coil-induced flux change. Smaller flux changes trigger tighter forward tolerances.

• max_rel_psit (float, optional (default=0.02)) – Maximum relative change in tokamak flux allowed due to coil current updates. Used as a step limiter.

• damping_factor (float, optional (default=0.995)) –

  Exponential damping applied to allowed tokamak flux change across iterations:

  allowed_change *= damping_factor**iteration

• use_full_Jacobian (bool, optional (default=True)) – If True, enables switching to full Jacobian optimisation when sufficiently close to convergence.

• full_jacobian_handover (list[float], optional (default=[1e-5, 7e-3])) –

  Two thresholds:

  [GS_residual_threshold, tokamak_flux_update_threshold]

  Below these values, the solver switches to full Jacobian mode.

• l2_reg_fj (float, optional (default=1e-8)) – L2 regularisation factor applied when using the full Jacobian.

• force_up_down_symmetric (bool, optional (default=False)) – If True, enforces up–down symmetry during forward solve.

• verbose (bool, optional (default=False)) – If True, prints iteration progress and diagnostic information.

• suppress (bool, optional (default=False)) – If True, suppresses all output (overrides verbose).

Returns:

The solution is written in-place into:

eq.plasma_psi eq.tokamak coil currents

Return type:

None

Notes

• The method is a nonlinear block-iterative inverse solver.

• Early iterations use a frozen-plasma approximation for robustness and efficiency.

• Close to convergence, the full Jacobian is used for quadratic-like convergence.

• For diverted equilibria, current updates are automatically resized to prevent loss of X-points or O-points.

optimize_currents(eq, profiles, constrain, target_relative_tolerance, relative_psit_size=0.001, l2_reg=1e-12, verbose=False)

Compute coil current updates using the full (plasma-aware) Jacobian.

This routine builds the Jacobian of the magnetic constraint residual vector with respect to the control coil currents:

A = d b / d I

where:

b = constraint residual vector I = control coil currents

Unlike the simplified Green’s-function approach (which assumes the plasma response is frozen), this method recomputes the full plasma equilibrium for each finite-difference perturbation of the control currents. Therefore the Jacobian includes:

• Direct magnetic effect of coil current changes

• Indirect effect via plasma response (GS equation)

The Jacobian is constructed using forward finite differences:

A[:, i] ≈ (b(I + δI_i) − b(I)) / δI_i

Once A is constructed, the Newton step is computed by solving the Tikhonov-regularised least-squares problem:

min || A ΔI + b0 ||² + ||R ΔI||²

where:

b0 = current constraint residual R = regularisation matrix

If current or flux limits are active, a quadratic optimisation routine is used instead of the closed-form normal equations.

Parameters:

• eq (freegsnke equilibrium object) –

  Contains:

  • Current coil currents

  • Current plasma_psi

  Modified only through auxiliary copies inside this routine.

• profiles (freegsnke profile object) – Provides plasma properties required to solve the forward Grad–Shafranov problem.

• constrain (freegsnke inverse_optimizer object) –

  Defines:

  • Control coils

  • Constraint residual vector b

  • Loss function

  • Optional current or flux limits

• target_relative_tolerance (float) – Forward GS tolerance used when recomputing equilibria for each Jacobian column.

• relative_psit_size (float, optional (default=1e-3)) – Target relative change in tokamak flux used to scale the finite-difference perturbations δI. Ensures perturbations are neither too small (noise-dominated) nor too large (nonlinear).

• l2_reg (float or 1D array, optional (default=1e-12)) – Tikhonov regularisation factor. If scalar: isotropic regularisation. If array: per-coil regularisation weights.

• verbose (bool, optional) – If True, prints progress during Jacobian construction.

Returns:

• Newton_delta_current (ndarray) – Optimal Newton update for control coil currents.

• loss (float) – Norm of the constraint residual after applying the computed Newton step (predicted reduction).

Notes

• Computational cost scales with number of control coils because one full GS solve is required per column.

• Provides quadratic convergence when close to solution.

• Intended to be used only when already near equilibrium.

port_critical(eq, profiles)

Transfer critical equilibrium and topology information from the plasma profile solver to the equilibrium object after solving the Grad–Shafranov equation.

This function synchronises diagnostic and structural information between the profile model and equilibrium container. It is typically called after a successful GS solve.

The following information is transferred:

• Magnetic axis (O-point)

• X-point locations (if present)

• Plasma boundary and limiter flags

• Total plasma current estimate

• Plasma profile state snapshot

• Tokamak vacuum flux (if available)

Parameters:

• eq (FreeGSNKE equilibrium object) – Equilibrium object that will be updated in-place.

• profiles (FreeGSNKE profile object) –

  Profile object used during GS solution. Must contain:

  • profiles.xpt : X-point locations (if diverted plasma)

  • profiles.opt : O-point / magnetic axis data

  • profiles.psi_bndry : Boundary flux value

  • profiles.flag_limiter : Limiter configuration flag

  • profiles.jtor : Toroidal current density profile

Returns:

Updates are performed in-place on the equilibrium object.

Return type:

None

Notes

This routine ensures consistency between:

• Plasma topology diagnostics

• Flux solution storage

• Current profile bookkeeping

It is primarily used internally by forward and inverse solvers.

relative_del_residual(res, psi)

Compute a relative residual measure based on the range (max − min) of the residual and flux field.

This metric measures relative variation using amplitude range rather than L2 norms, and is useful for detecting large-scale flux structure errors.

The relative residual is defined as:

relative_residual =

( max(ψ) − min(ψ) )

where:

ψ = plasma poloidal flux field res = nonlinear GS residual vector

This quantity is dimensionless and provides a scale-invariant measure of solution error magnitude.

Parameters:

• res (ndarray) – Flattened nonlinear Grad–Shafranov residual vector.

• psi (ndarray) – Flattened plasma poloidal flux field.

Returns:

• rel_residual (float) – Relative residual defined using amplitude range normalisation.

• del_psi (float) –

  Flux amplitude range:

  max(ψ) − min(ψ)

Notes

• Range-based metrics are less sensitive to localised spikes than pointwise residuals.

• If ψ has very small dynamic range, this metric may become ill-conditioned.

• Used internally as a convergence and stability diagnostic.

relative_norm_residual(res, psi)

Compute a relative residual using Euclidean (L2) norm normalisation.

This function measures the magnitude of the nonlinear GS residual relative to the magnitude of the plasma flux field using:

relative_residual =

||ψ||₂

where:

res = nonlinear Grad–Shafranov residual vector ψ = plasma poloidal flux field

This provides a scale-invariant convergence metric for nonlinear equilibrium solves.

Parameters:

• res (ndarray) – Flattened nonlinear Grad–Shafranov residual vector.

• psi (ndarray) – Flattened plasma poloidal flux field.

Returns:

Dimensionless relative L2 residual.

Return type:

float

Notes

• This metric is sensitive to small ||ψ|| values.

• Used as a primary convergence indicator during nonlinear solves.

solve(eq, profiles, constrain=None, target_relative_tolerance=1e-05, target_relative_psit_update=0.001, max_solving_iterations=100, max_iter_per_update=5, Picard_handover=0.1, step_size=2.5, scaling_with_n=-1.0, target_relative_unexplained_residual=0.3, max_n_directions=16, clip=10, max_rel_update_size=0.15, l2_reg=1e-09, forward_tolerance_increase=100, max_rel_psit=0.01, damping_factor=0.98, use_full_Jacobian=True, full_jacobian_handover=[1e-05, 0.007], l2_reg_fj=1e-08, force_up_down_symmetric=False, verbose=False, suppress=False)

Unified entry point for solving Grad–Shafranov problems (forward or inverse).

This method dispatches to either:

• forward_solve — fixed coil currents, solve for plasma equilibrium

• inverse_solve — adjust control coil currents to satisfy magnetic constraints

The mode is determined by the constrain argument:

constrain is None → Forward equilibrium solve constrain is provided → Inverse free-boundary optimisation

FORWARD MODE (constrain=None)

Solves the nonlinear free-boundary Grad–Shafranov equation:

F(ψ, I) = 0

for the plasma poloidal flux ψ, with fixed coil currents I.

The solve terminates when the relative GS residual falls below target_relative_tolerance or when max_solving_iterations is reached.

The result is written in-place into:

eq.plasma_psi

INVERSE MODE (constrain provided)

Solves the PDE-constrained optimisation problem:

min_I ½ || b(ψ(I), I) ||² subject to F(ψ(I), I) = 0

where:

ψ = plasma flux I = control coil currents b = magnetic constraint residual vector

The algorithm alternates between:

1. Optimising control coil currents

2. Performing forward GS solves to update plasma response

Convergence requires BOTH:

• GS residual < target_relative_tolerance

• Relative tokamak flux update < target_relative_psit_update

Upon convergence:

eq.plasma_psi contains the final plasma solution eq.tokamak contains the optimised coil currents

param eq:

Contains:

• Tokamak geometry and coils

• Current coil currents

• Current plasma flux (initial guess)

Modified in-place.

type eq:

FreeGSNKE equilibrium object

param profiles:

Defines plasma profiles used to compute toroidal current density:

• Pressure profile p(ψ)

• F(ψ) profile

• Total plasma current

type profiles:

FreeGSNKE profile object

param constrain:

If None:

Forward equilibrium solve is performed.

If provided:

Specifies control coils and magnetic constraints for inverse optimisation.

type constrain:

freegsnke inverse_optimizer object or None, optional

param target_relative_tolerance:

Target relative GS residual tolerance for forward solves. In inverse mode, this is the required final equilibrium accuracy.

type target_relative_tolerance:

float, optional (default=1e-5)

param target_relative_psit_update:

(Inverse mode only) Maximum allowed relative change in tokamak flux between successive current updates.

type target_relative_psit_update:

float, optional (default=1e-3)

param max_solving_iterations:

Maximum number of outer iterations. In forward mode: maximum GS iterations. In inverse mode: maximum inverse iterations.

type max_solving_iterations:

int, optional (default=100)

param max_iter_per_update:

(Inverse mode only) Maximum forward GS iterations performed after each coil current update.

type max_iter_per_update:

int, optional (default=5)

param Picard_handover:

Relative GS residual threshold above which Picard iteration is used instead of Newton–Krylov in forward solves.

type Picard_handover:

float, optional (default=0.1)

param step_size:

Scaling factor for nonlinear update steps in forward solve. Expressed in units of residual norm.

type step_size:

float, optional (default=2.5)

param scaling_with_n:

Additional scaling of nonlinear steps as:

(1 + iteration_number)**scaling_with_n

type scaling_with_n:

float, optional (default=-1.0)

param target_relative_unexplained_residual:

Forward solver termination criterion based on how much of the residual is explained by accumulated search directions.

type target_relative_unexplained_residual:

float, optional (default=0.3)

param max_n_directions:

Maximum number of search directions in forward solve.

type max_n_directions:

int, optional (default=16)

param clip:

Maximum update magnitude per search direction in forward solve.

type clip:

float, optional (default=10)

param max_rel_update_size:

Maximum allowed relative change in plasma_psi per forward iteration. Updates exceeding this are rescaled.

type max_rel_update_size:

float, optional (default=0.15)

param l2_reg:

(Inverse mode only) Tikhonov regularisation applied when using Green’s-function Jacobians for coil optimisation.

type l2_reg:

float or 1D ndarray, optional (default=1e-9)

param forward_tolerance_increase:

(Inverse mode only) Controls how tightly forward problems are solved relative to the magnitude of coil-induced flux changes.

type forward_tolerance_increase:

float, optional (default=100)

param max_rel_psit:

(Inverse mode only) Maximum allowed relative tokamak flux change due to a single coil current update.

type max_rel_psit:

float, optional (default=0.01)

param damping_factor:

(Inverse mode only) Exponential damping applied to the allowed tokamak flux change across iterations to encourage convergence.

type damping_factor:

float, optional (default=0.98)

param use_full_Jacobian:

(Inverse mode only) If True, switches to full plasma-aware Jacobian when sufficiently close to convergence.

type use_full_Jacobian:

bool, optional (default=True)

param full_jacobian_handover:

(Inverse mode only) Thresholds [GS_residual, tokamak_flux_update] below which full Jacobian optimisation is activated.

type full_jacobian_handover:

list[float], optional (default=[1e-5, 7e-3])

param l2_reg_fj:

(Inverse mode only) Regularisation factor applied when full Jacobian is used.

type l2_reg_fj:

float, optional (default=1e-8)

param force_up_down_symmetric:

If True, enforces up–down symmetry during forward solves.

type force_up_down_symmetric:

bool, optional (default=False)

param verbose:

If True, prints iteration progress information.

type verbose:

bool, optional (default=False)

param suppress:

If True, suppresses all printed output.

type suppress:

bool, optional (default=False)

returns:

The solution is written in-place into:

eq.plasma_psi eq.tokamak coil currents (inverse mode only)

rtype:

None

Notes

• This is the recommended high-level API for solving equilibria.

• Forward mode solves a nonlinear free-boundary PDE.

• Inverse mode solves a nonlinear PDE-constrained least-squares problem.

• Internally dispatches to forward_solve or inverse_solve.