🔬 physics

Derivation of the Exact Source Functions of the Grad–Shafranov equation

This paper presents a novel method for deriving exact source functions of the Grad–Shafranov equation by transforming it into a fourth-order partial differential equation to obtain independent solutions, which are then validated through numerical computation and applied to the ITER baseline scenario.

Original authors: Esam H. Abdul-Hafidh, Mohamed H Abdelhamed

Published 2026-08-21

📖 4 min read☕ Coffee break read

CC BY 4.0

Original authors: Esam H. Abdul-Hafidh, Mohamed H Abdelhamed

Original paper licensed under CC BY 4.0 (https://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Imagine a future where the energy that powers our cities comes not from burning fossil fuels, but from replicating the process that makes the sun shine. This is the promise of nuclear fusion, a clean and nearly limitless energy source that scientists have been trying to harness for decades. To make this work on Earth, researchers must trap superheated gas, known as plasma, inside a magnetic cage. This gas is so hot that no physical container could hold it; instead, powerful magnetic fields squeeze and shape it into a stable form. The challenge lies in predicting exactly how this plasma will behave, balancing the outward push of its own heat against the inward squeeze of the magnetic fields. If the balance tips even slightly, the plasma escapes, the reaction stops, and the energy is lost.

At the heart of this balancing act is a mathematical description called the Grad–Shafranov equation. Think of this equation as a rulebook that tells scientists how the magnetic fields and the pressure of the plasma must relate to one another to stay in equilibrium. For a long time, solving this rulebook has required making educated guesses about the shape of the pressure and the flow of electric current inside the plasma. These guesses were often chosen simply because they made the math easier to solve, rather than because they perfectly matched the physical reality of the plasma. This left a gap in our understanding: we could simulate the plasma, but we weren't entirely sure if the underlying rules we used were the exact ones nature follows.

In a new study, researchers Esam H. Abdul-Hafidh and Mohamed H Abdelhamed have closed this gap by deriving the exact mathematical expressions for these hidden rules. Instead of guessing the shape of the pressure and current, they worked backward from the fundamental laws of physics to find the precise formulas that must exist for the plasma to remain stable. They started with the standard equation governing the plasma and, through a series of rigorous mathematical steps, stripped away the unknowns until they were left with a single, complex equation that describes the system without needing any assumptions. Solving this new equation revealed four distinct mathematical patterns that the plasma's pressure and current must follow. These patterns include specific logarithmic relationships, a type of curve that had not been previously identified as the exact solution for this problem.

To prove that their new formulas actually work, the team applied them to the International Thermonuclear Experimental Reactor, known as ITER. This massive machine, currently under construction in France, is designed to be the world's largest fusion experiment. The researchers used the known specifications of ITER—such as its size, the strength of its magnetic fields, and the amount of electric current it will carry—to test their new method. They wrote a computer program to calculate the exact shape of the magnetic fields and the distribution of pressure inside the reactor based on their derived formulas. The results were striking. The simulation produced a stable plasma configuration that matched the expected operating conditions of ITER perfectly. The calculated pressure peaked at the center of the reactor at 638 kilopascals, and the magnetic fields held the plasma in a shape that matched the machine's design goals.

The study also explored how stable the plasma remains under different conditions. By adjusting the ratio of pressure to magnetic force, the researchers found that the plasma behaves erratically if this ratio is too low or too high. However, within a specific range that matches the design of ITER, the plasma settles into a steady, predictable state. This confirms that the exact formulas they discovered are not just mathematical curiosities, but practical tools that describe the real physics of fusion reactors. The team also calculated key safety numbers, known as safety factors, which indicate how well the magnetic field lines twist around the plasma to keep it contained. Their calculations showed that the plasma would be stable, with the safety factors falling exactly where engineers need them to be for safe operation.

This work provides a new, more reliable way to model fusion reactors. By replacing old guesses with exact mathematical descriptions, scientists can now simulate plasma behavior with greater confidence. The method is robust enough to be applied to ITER and other similar fusion devices, offering a clearer path toward understanding how to keep the plasma stable long enough to generate electricity. While the path to commercial fusion energy remains long and complex, having the exact rules for how the plasma behaves removes a layer of uncertainty. The researchers have shown that with the right mathematical tools, we can predict the behavior of the most extreme environments in the universe with precision, bringing the dream of fusion power one step closer to reality.

Technical Summary: Derivation of the Exact Source Functions of the Grad–Shafranov Equation

Problem Statement
The Grad–Shafranov equation (GSE) is the fundamental nonlinear second-order partial differential equation used to model magnetohydrodynamic (MHD) equilibrium in axisymmetric fusion devices like tokamaks. A critical challenge in solving the GSE is determining the "source functions"—plasma pressure ( $p$ ), poloidal current ( $f$ ), and toroidal current density ( $J_\phi$ )—which define the right-hand side of the equation. Historically, researchers have resolved this by assuming specific functional dependencies of $p$ and $f$ on the poloidal flux ( $\psi$ ), such as linear or quadratic forms (e.g., Solov'ev or Herrnegger-Maschke models). The authors argue that these assumptions lack firm physical justification and are primarily introduced to render the equation mathematically separable, thereby limiting the generality of the solutions.

Methodology
The paper proposes a rigorous mathematical derivation to determine the exact source functions without arbitrary assumptions about their dependence on $\psi$ . The methodology proceeds as follows:

Transformation to a Fourth-Order PDE: The authors systematically differentiate the GSE with respect to the spatial coordinates ( $r$ and $z$ ) to eliminate the derivatives with respect to the flux ( $\partial/\partial \psi$ ). This process transforms the original second-order equation into a single fourth-order nonlinear partial differential equation governing the source term function $G(\psi, r) = -\mu_0 r J_\phi$ .
Derivation of Exact Solutions: By solving this fourth-order equation, the authors identify four independent solutions involving logarithmic terms. This leads to an exact expression for $G$ that naturally imposes all powers of $\psi$ through its series expansion.
Analytical Expressions: Substituting the derived $G$ $G$ back into the governing relations yields exact analytical expressions for:
- Plasma Pressure ( $p$ ): A function involving $\psi^2 \ln(\beta \psi)$ and linear $\psi$ terms.
- Poloidal Current ( $f$ ): A function involving $\psi^2 \ln(\beta \psi)$ and a constant on-axis term.
- Toroidal Current Density ( $J_\phi$ ): A function combining $r$ , $1/r$ , and the logarithmic flux terms.
Numerical Validation: To determine the specific constants ( $a, b, \lambda, \beta$ $a, b, λ, β$ ) for a real-world scenario, the authors formulated three integral constraints based on:
- The global plasma current ( $I_{pl}$ ).
- The poloidal beta ( $\beta_p$ ).
- The internal inductance ( $l_i$ ).
  These constraints were solved numerically using a customized MATLAB code employing the finite difference method.

Key Results
The method was applied to the International Thermonuclear Experimental Reactor (ITER) baseline inductive scenario ( $I_{pl} = 15$ MA, $B_\phi = 5.3$ T). The study produced the following results:

Equilibrium Parameters: The code successfully determined unique values for the equilibrium constants ( $a, b, \lambda$ ) that satisfy the ITER constraints, with an error margin of less than 1%.
Flux and Profiles: The numerical solution yielded poloidal flux contours reaching a maximum of 11.39 Wb/rad at the magnetic axis. The resulting pressure profile peaked at 638 kPa on-axis, and the toroidal current density peaked at 1.48 MA/m².
Stability Metrics: The calculated safety factors were $q_0 \approx 0.99$ and $q_{95} \approx 3.02$ , consistent with stability requirements. The poloidal beta ( $\beta_p \approx 0.60$ ) and internal inductance ( $l_i \approx 0.83$ ) fell within the required ranges for stable confinement.
Sensitivity Analysis: The study demonstrated that values of $\beta_p$ outside the 0.6–0.65 range (e.g., 0.35 or 0.8) led to non-steady MHD behavior or pressure shifts, confirming the necessity of the derived parameters for stable equilibrium.

Significance and Claims
The authors claim that this approach is novel, as it derives the exact source functions of the Grad–Shafranov equation through systematic mathematical manipulation rather than relying on pre-assumed functional forms. By eliminating the need for arbitrary assumptions about the dependence of pressure and current on flux, the method provides a self-consistent framework for generating physically consistent equilibrium solutions. The paper asserts that this methodology can be readily applied to ITER and other axisymmetric tokamak fusion reactors to improve the accuracy of plasma equilibrium modeling.

More like this