Three dimensional equilibrium solutions for a current-carrying reversed-field pinch plasma with a close-fitting conducting shell

Toroidal magnetically confined plasmas such as the tokamak, spheromak, and reversed-field-pinch (RFP) are designed to be axisymmetric, but toroidally localized device elements as well as plasma dynamics can break this axisymmetry. This makes understanding non-axisymmetric equilibria through three-dimensional equilibrium reconstruction important. Fast ion confinement is particularly sensitive to even small deviations in toroidal symmetry.¹ In tokamaks, Resonant Magnetic Perturbations (RMPs) that deform the plasma boundary three-dimensionally are useful for ELM mitigation² and are planned to be used on ITER.^3,4 ITER itself has planned a helical hybrid-scenario.⁵ Helical “snakes” have been widely observed in tokamak cores,^6–8 and a single-helical-axis (SHAx) state has been observed in RFPs.^9–11 The SHAx state is of particular interest for the RFP community due to its suppression of high order tearing modes and beneficial confinement properties for fusion,¹² and it will be the focus of this paper.

In the axisymmetric state of the RFP, multiple tearing modes arise with comparable amplitude. In the helical state, the resonant mode closest to the core has significantly higher amplitude than all others, a signature of the transition to a helical equilibrium. The spectral purity can be quantified by the spectral index¹³ 

where $b̃n$ is the perturbed magnetic amplitude of a given tearing mode, as measured at the edge of the plasma. Axisymmetric RFPs have high spectral index, N_s > 2, and are referred to as “multi-helicity” (MH) RFPs. Helical RFPs possessing a SHAx with no residual tearing modes have N_s = 1 and are referred to as “single helicity” (SH) RFPs. This is a theoretical state. The intermediate state, where 1 < N_s < 2 and a SHAx may still be present, is referred to as quasi-single helicity (QSH).

The 3D equilibrium reconstruction code V3FIT¹⁴ has been used to generate equilibria for a number of plasma configurations.¹⁵ It has been employed in the past by RFX-mod to reconstruct helical RFP plasmas,¹⁶ but neither polarimetry data nor the eddy current effects due to their thin conducting shell were included. On the Madison Symmetric Torus (MST),¹⁷ external magnetic diagnostic coils have been used in NCT-SHEq to generate equilibria for use with both soft-x-ray and interferometer/polarimeter measurements,^18,19 but each effort focused on a single internal diagnostic, and neither included the internal measurements as constraints to an equilibrium minimization.

This article presents the first full 3D reconstructions of a helical RFP with a close-fitting conducting shell. The adaptations to VMEC-based V3FIT in order to include a perfectly conducting boundary in equilibrium reconstructions are included. This work is also the first to include polarimetry measurements as a constraint on equilibrium reconstruction of the helical RFP. These reconstructions have been performed using V3FIT with experimental data from the MST RFP.

MST is a toroidally symmetric machine with a circular poloidal cross section of minor radius a = 0.52 m and major radius R₀ = 1.5 m. Its 5-cm-thick conducting vessel acts as a toroidal field winding while also supplying vertical field and MHD mode stabilization via eddy currents.¹⁷ The vessel is composed of 6061-T6 aluminum¹⁷ and has a magnetic penetration time of 0.82 s,²⁰ much longer than the QSH duration. There are two toroidally extending, 1.3-cm-thick graphite limiters, at the outboard and inboard midplanes of the device. Typical magnetic equilibria are axisymmetric, albeit with large islands and stochastic field lines from multiple tearing modes.^21,22 The toroidal field on-axis is $Bϕ(0)≤0.5 T$⁠, decreasing to $Bϕ(a)≤0$ near the boundary. This $Bϕ$ reversal is a defining feature of the RFP. The poloidal field increases from B_θ(0) = 0 to $Bθ(a)≈12Bϕ(0)$⁠. A typical safety factor $q=rBϕ/RBθ$ is monotonically decreasing, from $q(0)≈0.2→q(a)≤0$⁠.

At high plasma current, I_p ≃ 500 kA and q(a) ≈ 0, MST plasmas transition to a SHAx equilibrium, often remaining in this state for the duration of the discharge. QSH plasmas are obtainable at lower I_p, but the probability of transition and the duration of QSH are observed to scale with Lundquist number S

so the probability of a saturated QSH state is lower.¹⁰ The signature of the QSH transition in MST is a rapid increase, then saturation, of the amplitude of the n = 5 tearing mode, coupled with cessation of the plasma rotation. Fig. 1(b) shows the saturation of the n = 5 mode at 25 ms, while the other fluctuating components of the field remain at small amplitude in the new equilibrium, which persists until 45 ms. The shift of the magnetic axis during this transition has been observed previously using the far-infrared (FIR) polarimetry diagnostic, which recorded an outward radial shift of 13 cm in addition to the pre-existing axisymmetric Shafranov shift of 5 cm, for a total shift of 18 cm.¹⁹ 

When the plasma rotation halts, the phase of the n = 5 helical axis with respect to machine coordinates is semi-random, having a preferred locking phase aligned with the uncorrected portion of the radial field error.²³ As will be discussed later, this complicates analysis since the diagnostics, whose positions are fixed, vary in effectiveness when viewing the helix at different locking angles. This is a particularly noticeable effect for the FIR interferometer/polarimeter system, whose chords pass vertically through MST to measure the line integrated product of density and vertical magnetic field.²⁴ Line integration causes the measurement to be insensitive to up-down asymmetry. The FIR system is important for characterizing the helical shift of the equilibrium profiles. The uncertainty in the helical shift due to line integration is most effectively mitigated when the helical magnetic axis is locked fully inboard or outboard (θ = 0°, 180°) at the toroidal location of the chords. Similarly, the Thomson scattering system consists of point measurements spaced vertically below the mid-plane at R = 1.5 m²⁵ and gains the most information when the helical axis is also locked downward (θ = 270°) at the location of the diagnostic.

Although the core of the plasma becomes helical, the conducting boundary maintains a circular last closed flux surface (LCFS). Eddy currents in the highly conductive shell are driven by the shifting of magnetic flux within the plasma. These currents exclude field from the shell, effectively closing and circularizing the flux surfaces at the boundary. The resulting equilibrium has a helical axis embedded in an axisymmetrically bounded system, seen in Fig. 2, a potentially favorable configuration for energy and particle confinement.²⁶ Previous V3FIT reconstructions performed without a conducting shell boundary condition showed a discrepancy between the observed and modeled perturbed n = 5 magnetic field measured at the edge of the plasma of ≈60%. This large discrepancy is attributed to the measurement location. MST's magnetic diagnostic coils are mounted on the vessel wall, <1 cm from the shell currents that contribute significantly to the total measured perturbed n = 5 magnetic signal. Because these currents are dynamically generated, are difficult to measure directly, and contribute the dominant portion of the measured magnetic field, their self-consistent calculation within equilibrium reconstruction is critical.

V3FIT is a 3D equilibrium reconstruction code built around the Variational Moments Equilibrium Code (VMEC).²⁷ VMEC solves the J × B = ∇P MHD force balance equation with spectrally decomposed flux surfaces in three dimensions using an energy minimization principle and assuming the existence of closed nested flux surfaces. To solve an equilibrium, a pressure profile and either a current or rotational transform/q profile as a function of flux are specified. The plasma boundary is either determined by specifying a fixed boundary shape or self-consistently determined by matching to vacuum fields.²⁸ In the latter free boundary solutions, vacuum fields are specified on a grid using a coil model. For stellarators and tokamaks, the external currents producing the field are well known, while the plasma shape will vary depending on field coil currents. By contrast, fields in MST are produced by eddy currents induced in the close fitting conducting shell. The induced eddy current distribution will vary depending on the resulting equilibrium. This makes defining an appropriate coil set to model this current distribution challenging. On the other hand, the close fitting conducting shell of MST defines a well characterized boundary making fixed boundary solutions ideal for MST.

Fixed boundary solutions in VMEC do not contain any information about external fields. As a result, VMEC alone cannot solve for the equilibrium that is self-consistent with induced eddy currents. V3FIT reconstructs a 3D equilibrium by adjusting the free parameters such that modeled signals, computed from the equilibrium model, best match with experimentally observed diagnostic signals. V3FIT model signals for magnetic diagnostic coils, line integrated interferometry and polarimetry, line integrated soft-x-ray emission, and point-measured Thomson scattering diagnostics. For some signals, models beyond the VMEC equilibrium are required. In these cases, density, temperature, and soft-x-ray emissivity are defined as flux-surface-constant profiles extending the VMEC solution.

Reconstructing the VMEC equilibrium with signals sensitive to the induced eddy currents allows V3FIT to find the equilibrium self-consistent with its own induced currents. In general, magnetic diagnostic signals are implemented by computing the response, R, between a diagnostic coil and a source of current I

where the response function R depends on the vector potential A and a factor f to convert from flux to magnetic field. Response functions for unit current are precomputed using a Biot-Savart integration of each current source using the V3RFUN code. Magnetic signals from plasma current and external coils are treated separately.²⁹ 

For the plasma current, response functions are precomputed on a 3D grid encapsulating the plasma volume. The total plasma signal response is computed from a volume integration of the plasma fields

where the current density J is determined from the curl of the VMEC fields.

When external sources of current are well defined, the response function is recast as a mutual inductance matrix M_ij between each external field coil i and the diagnostic coil j. The induced signal is determined from the external coil current I_i and mutual inductance matrix

If the current distribution in the shell were well known, eddy currents could be modeled as an additional set of external field coils. For MST, however, the shell current distribution will change depending on the equilibrium. To model the eddy currents as external coils, a set of coils modeling every possible equilibrium state would be required. In 3D, this quickly makes the problem intractable. Instead, a new method of calculating magnetic signals for a conducting shell is implemented that parallels the virtual casing method.^30,31

Knowing that the magnetic penetration time of the shell is much longer than the duration of the QSH state, the simplifying assumption can be made that the vacuum vessel is a perfect conductor. The resulting eddy currents in the shell will shield out the plasma fields. Thus the surface current K in the conducting shell is determined from the jump condition in the parallel magnetic field at the plasma edge H₁ and the shell H₂

In the MST device for the q(a) = 0 case, the only source of magnetic field is the plasma current with no fields outside the shell, making H₂ = 0. The signal contribution due to the conducting shell can be computed from this surface current density and a response function of the shell geometry

Shell response functions R_shell are precomputed on a 2D grid representing a conducting surface in the shell. For the MST response function, the conducting surface is placed half way within the 5-cm-thick shell to avoid singularities near the magnetic diagnostic coils.

To compute the conducting shell signal contribution, only the fields parallel to the shell (B_θ and $Bϕ$⁠) are needed. Since the MST vacuum vessel is close fitting to the plasma, the fields at the last closed flux surface are used, and the vacuum region is ignored for simplicity. In the VMEC coordinate system, the eddy current signal from Equation (7) becomes

where $⁣$ is the Jacobian of the VMEC coordinate system composed of radial s, poloidal u, and toroidal v coordinates. The e^s is the contravariant basis vector and points normal to a flux surface at any u,v location. Note that e^s is not a unit vector and must be normalized when used as the $n̂$ direction. Equation (8) reduces to

where B_u and B_v are the covariant components of the magnetic fields for the u and v coordinates. The covariant basis vectors e_u and e_v are orthogonal to the flux surface normal. The total signal for a MST magnetic diagnostic is the sum of the plasma and shell contributions, seen in Fig. 3. The contribution to the signal from the shell reduces the total signal from the plasma-only signal. Figure 4 shows the resulting eddy current distribution of the shell-currents along the entire torus for a helical equilibrium.

For the MST RFP, V3FIT was benchmarked by a series of lower current, axisymmetric cases against MSTFit, a 2D equilibrium reconstruction code developed for MST.³² With a conducting boundary, V3FIT successfully reproduces axisymmetric equilibrium solutions obtained using MSTFit. While VMEC solves the 3D MHD force balance equation using spectrally decomposed flux surfaces, the equilibrium solution portion of MSTFit solves the Grad-Shafranov equation assuming shifted circular flux surfaces. The conducting shell in MSTFit is modeled with 72 toroidal current filaments. Pre-calculated unit current vector potentials are used to determine the current distribution in the shell that satisfies b_r = 0 at the plasma edge, a constraint that follows from a perfectly conducting shell assumption. Similar solutions from the two codes provide confidence in the capability of V3FIT to reconstruct RFP plasmas with eddy currents included.

The first comparison is between model signals for the poloidal array of 16 magnetic diagnostic coils measuring B_θ, located at $ϕ=0°,θ=0°−360°$⁠. This diagnostic is the ideal test because it is sensitive to both the eddy currents in the shell and the poloidal asymmetry of magnetic field in the axisymmetric RFP equilibrium. A scan in I_p, and thus a scan of B_θ, was performed in $Bϕ(a)=0$ discharges, seen in Fig. 5. MSTFit and V3FIT results show qualitative agreement with each other and quantitative agreement with the goodness-of-fit measurements, $⟨χMSTFit2⟩=11.2$ and $⟨χV3FIT2⟩=15.8$⁠. The differences are likely due to different implementations of the perfectly conducting shell and different methods for parameterizing plasma characteristics. The consistent overestimation of both V3FIT and MSTFit reconstructed signals compared to the observed signals is attributed to the location of the measurements, which occur at a poloidal gap in the conducting shell. Both reconstruction codes assume a smooth conducting shell with no gap and, therefore, do not account for any flux expansion that might arise at this toroidal location.

The second comparison is between profiles reconstructed at different $Bϕ(a)$⁠. MSTFit and V3FIT were run for improved confinement discharges, where pressure is peaked on-axis, and the magnetic field is deeply reversed such that q(a) ∼ −q(0).^33,34 These equilibria are compared to non-reversed (⁠$Bϕ(a)=0,q(a)=0$⁠) equilibria in Fig. 6, where pressure is flatter in the core. Safety factor q and pressure profiles with uncertainties are shown for example, V3FIT cases, with the MSTFit result included. For RFPs, it is necessary to use flux in the helical-poloidal direction as a flux surface label since the toroidal field is characteristically reversed in the plasma periphery. In the axisymmetric limit, the helical-poloidal flux used in V3FIT is comparable to the poloidal flux surfaces MSTFit uses, making pressure and q-profiles comparable between the two codes. While the general trend and some characteristics agree between the V3FIT and MSTFit equilibria, V3FIT consistently returns a reverse-sheared q-profile for deeply reversed plasmas, a feature that is not present in MSTFit equilibria. The spline parameterizations employed within VMEC for V3FIT are more flexible than the stiffer functional form assumed in the MSTFit run, generating wide deviations between the core and mid-radius. The spline parameterizations also have larger uncertainty in areas where they are poorly constrained by data.

For MST plasmas with N_s < 1.5, V3FIT reliably converges to a helical equilibrium solution. The equilibrium has a shifted, bean-shaped core that precesses along with the helical magnetic axis, with circular flux surfaces at the edge, as in Fig. 2. The V3FIT flux surfaces are roughly similar in shape to previous NCT-SHEq reconstructions.¹⁹ There is also quantitative agreement on the observed helical shift of the magnetic axis of 15 cm and an axisymmetric Shafranov shift of 3 cm. Reconstructions consistently achieve a $χr2≈1.25$⁠. $χr2$ is the reduced χ² measure of goodness-of-fit

where $Si0$ is the ith observed signal, $Sim$ is the ith modeled signal, σ_i is the uncertainty in the ith observed signal, and ν is the degrees of freedom equal to the number of signals input into the reconstruction minus the number of parameters used to create the fit.

Because of the orientation of the fixed diagnostics in MST, there are two locked orientations, separated poloidally by 180°, that allow maximum combined effectiveness of the measured signals. The Thomson scattering and finite impulse response diagnostics are separated toroidally by 30°, corresponding to a 150° poloidal rotation of the m = 1, n = 5 helical axis between the two diagnostics, as shown in Fig. 7. However, as discussed in Section II, obtaining maximum effectiveness from each diagnostic individually would require a poloidal rotation of only 90° between the two diagnostics. Therefore, the first ideal locking phase for measurement is a compromise between these two orientations, with the mode locked at θ = 30° at the FIR diagnostic and θ = 240° at the Thomson scattering diagnostic. This orientation is referred to as the “outboard” orientation for the position of the phase at the FIR system. The second orientation, with the mode ideally at θ = 210° at the FIR diagnostic and θ = 60° at the Thomson scattering diagnostic, is less optimal because the magnetic axis is no longer in view of the Thomson scattering diagnostic. However, the density and magnetic axis shift are still distinguishable by the interferometer-polarimeter, and the helical axis opposite to the temperature diagnostic allows a unique set of measurements across the effective plasma mid-radius. This orientation is referred to as the “inboard” orientation. In practice, it is unlikely for the plasma to naturally lock at one of these two locations.²³ The two plasmas under consideration that represent the inboard and outboard orientations have a phase at the FIR system of θ = 173° and θ = 344°, respectively, as shown in Figs. 8(a) and 8(b).

Reconstructions for the inboard and outboard orientations were compared in order to verify the ability of diagnostics and subsequent fits to correctly identify the helical structure. For these reconstructions, all available diagnostics were incorporated, comprising 200 unique measurements. V3FIT successfully reconstructed equilibria for both the inboard and outboard cases, as shown in Fig. 9. Most critically, the implementation of the conducting boundary allows the edge magnetic data to be fit successfully, as seen in Fig. 9(a). From the V3FIT reconstruction, 37% of the measured fluctuating field at the plasma edge was due to plasma currents and 63% was due to eddy currents in the shell.

The difference between inboard and outboard cases is evident for each of the internal diagnostics, as shown in Figs. 9(b)–9(e). On the interferometer and soft-x-ray cameras, there is a clear shift due to the θ = 180° difference in magnetic axis location. On polarimetry, it is clear that the inboard case enhances the product of magnetic field and density on the inboard side. The zero-crossing in the Faraday rotation angle, Ψ = 0°, is shifted inboard for the inboard case, and outboard for the outboard case. Lastly, Thomson scattering measurements for the outboard case indicate a plateau of temperature from the geometric center to Z ≈ –0.2 m of the machine because it is sampling near the helical core, over relatively few flux surfaces, seen in Fig. 8(d). In the inboard case, the magnetic axis is further away from the measurement locations (see Fig. 8(c)), so a monotonically decreasing temperature with radius is observed.

The profiles in flux coordinates were compared in order to gauge the consistency of V3FIT reconstructions at multiple orientations, as shown in Fig. 10. As these were different discharges, there is no guarantee that the dynamically achieved equilibria of the plasma are perfectly comparable. Additionally, due to a lack of Thomson scattering coverage near the location of the axis, the initial guess for the pressure in the core for the inboard case was not well informed by measurements. Nonetheless, the magnetic equilibria are similar, with a reversed-shear q-profile that peaks below q = 0.2. The shapes of the pressure profiles are not as well matched as the q-profiles, especially in the core, but this may be due to poor guessing of the core pressure value and poor constraint of the pressure profile by the Thomson scattering diagnostic's measurement uncertainty.

With the addition of a conducting boundary, V3FIT successfully produces equilibrium reconstructions for RFP plasmas. There are two important consequences of these reconstructions. First, the extensible V3FIT model can now be used to produce helical reconstructions of MST plasmas, allowing for clear interpretation of internal diagnostic measurements. Second, the magnitude of the measured magnetic field due to the shell currents can now be quantified and separated from the magnetic field due to currents in the plasma. This separation may allow for better comparisons between RFX-mod helical states and MST helical states, since eddy current contributions to edge measurements of the two devices are very different.

Both the V3FIT eddy current model and the underlying VMEC equilibrium parameterization are different than those in MSTFit. While both models converge to similar χ² for external plasma diagnostics, the resulting pressure profiles diverge, particularly in the core and extreme edge. This behavior is expected, particularly for RFP plasmas, where the pressure profile has a weak effect on the overall equilibrium.

Reconstructions of the inboard and outboard locking positions for internal measurement converged to similar helical equilibria, albeit with pressure profile differences in the core. This area was outside the field of view for Thomson scattering in the inboard case, and therefore underconstrained. With greater diagnostic coverage, it is expected that equilibria would converge with similar pressure profiles. Notably, the eddy current model successfully reproduces the measured magnetic fields at the plasma edge, where the helical field is dominated by contributions from eddy currents. These results also confirm the ability of polarimetry and interferometry to constrain helical reconstructions and distinguish between inboard and outboard mode locking.¹⁹ 

Future experiments could focus on larger datasets with the plasma locked at optimal positions and, if possible, greater internal diagnostic coverage. In particular, a larger Thomson scattering dataset will reduce errors, increasing the signal effectiveness of the diagnostic. Inclusion of the Motional Stark Effect diagnostic to measure magnetic field strength would also help to constrain reconstructions.³⁶ 

We would like to thank K. J. McCollam and M. D. Nornberg for their careful review of the manuscript. MST data shown in this paper can be obtained in digital format (supplementary material).³⁷ This material is based upon work supported by the U.S. Department of Energy Office of Science, Office of Fusion Energy Sciences program under Award No. DE-FC02-05ER54814.