Effects of magnetic helicity on 3D equilibria and self-organized states in KTX reversed field pinch

The reversed field pinch (RFP) is a toroidal magnetic configuration in which plasmas can spontaneously transform into different self-organized states. Among various states, the ‘quasi-single-helical’ (QSH) state has a dominant component for the magnetic field and significantly improves confinement. Many theoretical and experimental efforts have investigated the transitions among different states. This paper employs the multi-region relaxed magnetohydrodynamic model to study the properties of QSH and other states. The stepped-pressure equilibrium code (SPEC) is used to compute MHD equilibria for the Keda Torus eXperiment (KTX). The toroidal volume of KTX is partitioned into two subvolumes by an internal transport barrier. The geometry of this barrier is adjusted to achieve force balance across the interface, ensuring that the plasma in each subvolume is force-free and that magnetic helicity is conserved. By varying the parameters, we generate distinct self-organized states in KTX. Our findings highlight the crucial role of magnetic helicity in shaping these states. In states with low magnetic helicity in both subvolumes, the plasma exhibits axisymmetric behavior. With increasing core helicity, the plasma gradually transforms from an axisymmetric state to a double-axis helical state and finally to a single-helical-axis state. Elevated core magnetic helicity leads to a more pronounced dominant mode of the boundary magnetic field and a reduced core magnetic shear. This is consistent with previous experimental and numerical results in other RFP devices. We find a linear relationship between the plasma current and helicity in different self-organized states. Our findings suggest that KTX may enter the QSH state when the toroidal current reaches 0.72 MA. This study demonstrates that the stellarator equilibrium code SPEC unveils crucial RFP equilibrium properties, rendering it applicable to a broad range of RFP devices and other toroidal configurations.


Introduction
In recent decades, numerous magnetic confinement configurations have been proposed for nuclear fusion, with the tokamak and the stellarator emerging as the most widely recognized.The reverse field pinch (RFP), investigated as an alternative concept, features a toroidal magnetic field at the plasma edge that is opposite to the core.It has the advantage of using ohmic heating to achieve fusion conditions [1,2].In contrast to stellarators, the RFP primarily generates its magnetic field through the plasma itself, without heavy reliance on external coils.The poloidal field of RFP has a comparable order of magnitude to the toroidal field, which is different from that of a tokamak.The RFP plasma exhibits diverse states due to the presence of abundant self-organized effects.
In the axisymmetric equilibrium state, the safety factor q of RFP is strongly sheared with q(0) > 0 at the magnetic axis and q(a) ≤ 0 at the plasma edge.There are a large number of m = 1 rational surfaces in the entire region.These rational surfaces are very unstable and it is easy for magnetic islands to emerge.These magnetic islands may overlap and turn the plasma core into chaos [3,4].Destruction of nested magnetic surfaces reduces particle confinement capabilities.Experimentally, this state is characterized by the presence of many m = 1 Fourier components of the magnetic field with comparable amplitude.
A self-organized helical equilibrium could become a new paradigm for ohmically heated fusion plasma [5,6].By increasing the plasma current, a spontaneous transition to a helical equilibrium occurs, with a high-temperature region in the core.A m = 1 and n > 3R/2a Fourier component of the magnetic field grows up and becomes dominant over the other secondary modes, where a and R are the minor and major radius of RFP configuration.This state is called the quasi-single-helical (QSH) state.
Two types of QSH states with distinct topologies have been observed [7,8].When the dominant component of the magnetic field is up to a few percent of the average magnetic field, there are two O-points in the magnetic field and two twisted magnetic axes.This state is known as the double-axis (DAx) state.As the amplitude of the dominant component increases, a new magnetic topology, which is characterized by a single helical magnetic axis has been experimentally observed.This new state is called the single-helical-axis (SHAx) state.
Taylor's relaxation theory can explain many plasma phenomena in RFP, including the reversed toroidal field [9,10].According to Taylor's theory, the final state of relaxation is the state with minimum energy subject only to the single invariant where A is vector potential such that B = ∇×A and H is the total magnetic helicity.For a plasma enclosed by a perfectly conducting toroidal shell, the corresponding equilibrium of this state satisfies where µ is a constant.As a consequence, the plasma current is everywhere parallel to the magnetic field, and the pressure gradient ∇p vanishes.Such a force-free field, meaning no pressure gradient force, can arise in both laboratory and space plasma.
If we ignore the toroidal curvature, the axisymmetric solution of Equation ( 2) can be expressed using the Bessel function, where r, θ, z are cylinder coordinates.Reversal of the toroidal field will occur and B z (a) < 0 when µa > 2.404.While Taylor's theory has achieved significant success, experimental results show that the plasma in the RFP is not a strict global force-free equilibrium.The ratio of parallel current to magnetic field λ = µ 0 j ∥ /B should be a strict constant in a force-free field.However, the measurement of the current density in MST indicates that the λ of the plasma core and edge are different [11].The experiments of RFX-mod show that there is an electron transport barrier, and the gradients of various physical parameters are large near this barrier [12].Classic Taylor's theory cannot describe the 3D phenomena in RFP, so the stellarator codes have become an important tool for studying physics without axial symmetry in RFP.Gobbin et al. used the stellarator code VMEC for the first time to simulate the three-dimensional equilibrium in RFX-mod with the assumption of nested magnetic surfaces [13].There are some investigations of 3D characteristics and dynamo effect based on VMEC equilibrium in KTX [14,15].Furthermore, some studies show that the equilibrium of RFP without nested surfaces can also be reproduced [16,17].
The multi-region relaxed magnetohydrodynamics (MRxMHD) is a new fluid model that combines the ideal MHD and Taylor relaxation [18].MRxMHD considers the helicity constraint of different subvolumes instead of the global constraint.Dennis et al. pointed out that the use of this model requires only a finite number of constraints to describe the self-organized plasma states in RFP [19].However, this has not been fully investigated, and the impact of these constraints on the self-organized state of plasma has not been shown.A multiple domains scheme (MDS) is proposed to investigate the transport problem in the presence of magnetic islands, which also divides the plasma into multiple regions [20].The MDS was applied to RFP plasma, which has been implemented in a new numerical tool [21].
Here, we employ a barrier as the interface to divide the toroidal volume into inner and outer subvolumes.Different self-organized states, including axisymmetric, DAx, and SHAx, can be reproduced numerically with the MRxMHD model.Helicity plays a critical role in this model and magnetic topology changes with helicity.
This paper is structured as follows.We begin with an introduction to the variational principle of MRxMHD in Section 2. Details of the 3D equilibrium simulation of the Keda Torus eXperiment (KTX) are presented in Section 3.There are more magnetic field properties of different self-organized states in Section 4. The topological properties and the decomposition of helicity are shown in Section 5. We summarize in Section 6.

Multi-region relaxed magnetohydrodynamics
The energy functional of the ideal MHD is given by the integral where γ is the adiabatic index and the plasma volume V is bounded by a toroidal surface.The ideal MHD equilibrium equation J × B = ∇p is obtained when Equation ( 5) is extremized by using variational principle, [22,23] δW where ideal displacement ξ vanishes on the plasma boundary.In Taylor's theory, the magnetic energy functional imposing constraint of the total magnetic helicity with a Lagrange multiplier α is expressed as [9,24] The first variation in F Taylor with respect to A is When δF Taylor = 0, the field satisfies Equation ( 2) with µ = 2α/µ 0 .The first step towards constructing the MRxMHD energy functional in RFP is to separate the space into two subvolumes by an interface I.In each subvolume, the magnetic helicity H i is conserved, where i = 1, 2 represents the inner and outer subvolume respectively.The MRxMHD energy principle aims to minimize the plasma energy subject to these constraints, where W i is the energy given in Equation ( 5) and µ i is the Lagrange multiplier.The first variation in the local constrained functional F i due to the arbitrary variations δB = ∇ × δA in the field and the arbitrary variations ξ in the interface I is given by [18] When F is extremized, the magnetic field in each V i follows the Beltrami equation, and the total pressure must be continuous across the interface I, where the bracket [[ * ]] indicate the jump of the scalar across interface I.
The stepped pressure equilibrium code (SPEC) serves as a numerical solver for the MRxMHD model [18,25,26].In SPEC, there could be any number of nested toroidal interfaces and the subvolume V i enclosed by I i−1 and I i .The total pressure across I i satisfies Equation ( 12) and the magnetic field in each V i follows Equation (11).Therefore, the pressure gradient ∇p is zero in each subvolume and the global pressure profile shows a stepped shape.SPEC solves the magnetic field of each subvolume V i by iteratively adjusting the shape of the interface until the total force across the interface I i is balanced.

Self-organized equilibrium simulation of KTX
The KTX with a measurement system of helicity flux density located at the University of Science and Technology of China is a large RFP device [27,28].The KTX vacuum chamber has a circular cross-section with major radius R = 1.4 m and minor radius a = 0.4 m.
We use SPEC to simulate the three-dimensional equilibrium in KTX.There are two subvolumes and an interface in the simulation.The state of the magnetic field in each subvolume is the state of minimum energy with invariant helicity, respectively.The magnetic field is obtained by solving Equation (11) in SPEC and choosing an appropriate interface shape to balance the total pressure across I.According to Taylor's theory, the entire plasma should be in force-free equilibrium and λ is a constant globally.This is unattainable at the plasma edge.Although a polynomial function model (PFM) with λ = λ 0 (1 − (r/a) α ) is proposed [29], our simulation based on the MRxMHD uses two radially discontinuous λ profiles to solve this problem.Figure 1 illustrates the radial λ profiles in different models.
Without the assumption of nested surfaces, there could be magnetic islands and chaotic fields.To enforce various boundary conditions, it is convenient to use toroidal coordinates (s, θ, ζ), which are adapted to the plasma boundary and the interfaces.For stellarator symmetric equilibrium, the toroidal surface x s can be described by where e R and e Z are the basis vectors in cylindrical coordinates.Similarly, the vector potential on a toroidal surface,     As mentioned above, we use trigonometric function basis in the poloidal and toroidal directions to express physical quantities, and we use polynomial function basis in the radial direction to discretize these quantities.The input data in these simulations are toroidal flux Ψ i , stepped pressure p i and helicity H i , where i = 1, 2 represent the inner and outer subvolume, respectively.The difference between the two stepped-pressures is set to p 1 − p 2 = 560 Pa and the total toroidal flux is Ψ = Ψ 1 + Ψ 2 = 0.1298 Wb according to the discharge and design parameters of KTX.The step of the pressure profile is determined by the magnetic flux of the subvolume.Based on previous results [4,12], the transport barrier locates at a relatively outer position in the radial direction.In our simulaton, the flux of inner subvolume is Ψ 1 = 0.86Ψ, which keeps the position of the interface relatively out.The SPEC outputs are the shape of the interface and the vector potential A. The magnetic field can be obtained by B = ∇ × A.
In our simulations, we found different core magnetic helicity density ⟨H⟩ results in magnetic fields with different topologies.Here, core magnetic helicity density ⟨H⟩ is the average of the helicity in the inner subvolume, The system is axisymmetric when the helicity is at a lower density.As ⟨H⟩ increases, the plasma spontaneously generates three-dimensional structures.The magnetic field gradually transforms from an axisymmetric state to DAx state and then to SHAx state.Figure 2 shows the flux surfaces of these different self-organized states.There is a new magnetic axis generated and the old magnetic axis disappears in the process of changing the magnetic topology.Figure 3 shows more details of the magnetic field through the Poincaré cross-section of six cases with different helicity densities.It should be noted that there is no Poincaré plot in the outer subvolume because the toroidal field is very small at the edge.

Magnetic field properties of different self-organized states
We simulated a large number of self-organized states with ⟨H⟩ ranging from 5.2×10 −3 to 9.2 × 10 −3 T 2 • m and analyzed the details of the magnetic field in different states.When ⟨H⟩ is less than 5.84 × 10 −3 T 2 • m, the system is always in the axisymmetrical state.If the core helicity density exceeds this threshold, there will be a new magnetic axis and the magnetic field will be in the DAx state where two flux tubes are twisted with each other.New flux tubes become larger as helicity increases.When ⟨H⟩ > 8.61 × 10 −3 T 2 • m, the old flux tube disappears, there is only a large flux magnetic tube and the plasma is in the SHAx state.Figure 3 illustrates this regulation by six Poincaré plots with different core helicity density.A m = 1, n = 6 Fourier component of the magnetic field at the plasma edge becomes more and more evident as the core magnetic helicity density increases.Figure 4 shows   When the plasma spontaneously generates a three-dimensional structure, B 1,6 as the dominant mode is more obvious than other secondary modes.Figure 5 shows a Poincaré plot and n spectrum of m = 1 magnetic modes in a typical SHAx state and  there is a m = 1, n = 6 dominant mode at the edge in the three-dimensional selforganized state of KTX.The same phenomenon has been experimentally discovered in other RFPs and there are different toroidal mode numbers of dominant mode because of the difference in minor and major radius.RFX-mod has a m = 1, n = 7 dominant mode and MST has a m = 1, n = 5 dominant mode [6,30,31].Although there is a dominant mode in our simulation, the other secondary modes will not disappear.Figure 6 shows that the amplitude of other m = 1 secondary helical modes that still exist here is about 10% of the dominant mode.The magnetic shear of the new flux tube in the DAx state and the core in the SHAx state will vanish.Figure 7 shows the inner subvolume safety factor profiles of three selforganized states with different ⟨H⟩.Both the DAx and SHAx states have lower safety factors than the axisymmetric state.The safety factor of the magnetic island in the DAx state and the core in the SHAx state is a constant of 1/6.This shows that the helical magnetic field in the core of the SHAx state is a large magnetic island.
There is an apparent positive linear correlation between ⟨H⟩ and the plasma current.In different self-organized states, the current is enhanced with increasing magnetic helicity.Figure 8 shows this relationship.The helicity flux density measurement system is implemented in KTX, and a similar linear relationship was obtained [28,32].For KTX, when the average toroidal magnetic field ⟨B t ⟩ ≈ 0.26 T, we estimate that a plasma current of about 0.6 MA is needed to spontaneously transform into a threedimensional state and when the current reaches 0.72 MA, the plasma can self-organize to form the SHAx state.There is also a higher probability of obtaining QSH spectra   with higher plasma current experimentally [1,33,34,35].

Decomposition of the magnetic helicity
Mathematically, the sum of the Gauss linking number over every pair of field lines within a volume is the magnetic helicity [36].The linking number L xy of the two curves x and y is If there are N tubes and the flux carried by each tube is Ψ i , i = 1, . . ., N , we can calculate the sum of linking numbers with flux weighting.If we let N → ∞ with Ψ i → 0, the tubes degenerate into field lines and we have where the Coulomb gauge is used in the vector potential.Helicity measures the link and twist between magnetic field lines together with the writhing of the field lines themselves [37,38,39].We can analyze the components of the helicity at different states in RFP.As shown in Figure 9, the transition of these self-organized states is accompanied by the generation and growth of the new flux tube Ψ n , and the reduction and disappearance of the old flux tube Ψ o .We can use the mutual helicity to describe the link between these two tubes The self helicity within the flux tube can be calculated where ι is the rotational transform and the average Gauss linking number of the field line and the axis in one toroidal period.According to the Cǎlugǎreanu invariant [40,41], the linking and self helicity can be produced by the writhing of the axis itself and the twisting of the field lines around the axis.The writhing helicity H writhing = WrΨ 2 measures the three-dimensional structure of the magnetic axis with the writhing number The twisting helicity H twisting = H self − H writhing reflects the contribution from the twisting of field lines around the axis.Figure 10 illustrates the decomposition of helicity in different DAx states.Note that when calculating magnetic helicity, we ignore chaotic regions.With the changes    10a, the H self of the old flux tube gradually decreases and the H self of the new flux tube gradually increases, both of which are significantly larger than the mutual helicity H mutual .Figure 10b describes more details of the twisting and writhing components in self helicity.The magnetic axis of the new flux tube has a significant three-dimensional structure, thus its writhing helicity H writhing is much larger than that of the old flux tube.The H self of the old flux tube mainly consists of the twisting helicity H twisting , while the writhing helicity H writhing is small enough to be ignored.

Summary and discussion
Different self-organized states on KTX were simulated using the SPEC code.We use the barrier as an interface to divide the plasma into two subvolumes, and the toroidal flux in these two subvolumes remains unchanged, ensuring that the barrier is always in a relatively out position.The relaxed state in each region minimizes the energy functional under the local helicity constraint.We reproduced different self-organized states by changing the core helicity density, which plays an important role in the plasma relaxed state.The system is axisymmetric with low helicity and as the core helicity density ⟨H⟩ increases, the self-organized state of the plasma gradually transforms into the DAx state and then into the SHAx state.
Beyond its three-dimensional geometric structure, the high ⟨H⟩ self-organized state exhibits distinct magnetic field properties.The high helicity makes the dominant mode at the plasma edge more pronounced, characterized by a toroidal number of 6.The emergence of a three-dimensional structure eliminates magnetic shear in the core, corroborating experimental observations on other RFP devices.Additionally, we found an approximate linear relationship between the plasma toroidal current and the core helicity density.With this information, we predict that the KTX plasma will enter the DAx state and have three-dimensional structures when the plasma current is higher than 0.6 MA.Furthermore, when the plasma current reaches 0.72 MA, the plasma will spontaneously transform into the SHAx state.
We further analyzed the components of magnetic helicity.Changes in magnetic topology are accompanied by the generation of a new flux tube and the disappearance of the old flux tube.The total helicity is primarily composed of the self-helicity within these two flux tubes, rather than the mutual helicity.The helicity of the old flux tube was predominantly characterized by twisting, while the new flux tube has a certain writhing helicity with an evident three-dimensional structure.
The MRxMHD model and the stellarator equilibrium code SPEC reveal essential equilibrium properties of RFP.SPEC demonstrates its capability for RFP devices and can be used in other RFP applications, such as equilibrium reconstructions.In the future, numerical results can be further validated through experiments.

Figure 1 :
Figure 1: Radial profiles of λ in different models.

Figure 2 :
Figure 2: Flux surfaces of different self-organized states.

Figure 3 :
Figure 3: Poincaré plot at ζ = 0 in the inner subvolume with different ⟨H⟩.The outmost red curve is the cross section of the interface.

Figure 4 :
Figure 4: The relationship between the amplitude of the m = 1, n = 6 mode and the core helicity density ⟨H⟩.

Figure 5 :
Figure 5: The Poincaré plot and n spectrum of m = 1 modes in a typical SHAx state.

Figure 6 :
Figure 6: The amplitude of the dominant mode and the secondary mode as a function of ⟨H⟩.

Figure 8 :
Figure 8: Plasma current as a function of ⟨H⟩.

Figure 9 :
Figure 9: The growth of the new flux tube (orange) and the reduction of the old flux tube (light blue).

7 HH
twisting of the old tube H writhing of the old tube H twisting of the new tube H writhing of the new tube (b) twisting helicity and writhing helicity

Figure 10 :
Figure 10: Helicity decomposition in different DAx states.Ψ n /Ψ o is the magnetic flux ratio of the new flux tube to the old tube