{'fraction_non_alphanumeric': 0.03262265604931929, 'fraction_numerical': 0.012843565373747753, 'mean_word_length': 4.645591647331787, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Observations of Heteroclinic Bifurcations in Resistive MHD Simulations of the Plasma Response to Resonant Magnetic Perturbations

A new class of static magnetohydrodynamic (MHD) magnetic island bifurcations is identified in rotating spherical tokamak plasmas during single- and two-fluid resistive MHD simulations. As the magnitude of an externally applied non-axisymmetric magnetic field perturbation is increased in these simulations, the internal flux surfaces that make up a sub-set of the resonant helical magnetic islands in the plasma gradually elongate and undergo heteroclinic bifurcations. The bifurcation results in the creation of a new set of hyperbolic-elliptic fixed points as predicted by the Poincar\'e-Birkoff fixed point theorem. Field line calculations without including the resistive MHD plasma response to the applied perturbation field do not undergo this class of bifurcations indicating the importance of plasma self-organization in the bifurcation process.

A new class of static magnetohydrodynamic (MHD) magnetic island bifurcations is identified in rotating spherical tokamak plasmas during single-and two-fluid resistive MHD simulations. As the magnitude of an externally applied non-axisymmetric magnetic field perturbation is increased in these simulations, the internal flux surfaces that make up a sub-set of the resonant helical magnetic islands in the plasma gradually elongate and undergo heteroclinic bifurcations. The bifurcation results in the creation of a new set of hyperbolic-elliptic fixed points as predicted by the Poincaré-Birkoff fixed point theorem. Field line calculations without including the resistive MHD plasma response to the applied perturbation field do not undergo this class of bifurcations indicating the importance of plasma self-organization in the bifurcation process. Resistive magnetohydrodynamic (MHD) theory and simulations are fundamental tools for understanding the physics of topological bifurcations in the trajectories of magnetic field lines that are associated with tearing and reconnection in astrophysical [1], solar [2], space [3] and magnetically confined toroidal plasmas [4]. While magnetic tearing and reconnection theory is conceptually well established and widely accepted, there are open questions concerning the detailed dynamics responsible for these topological changes that are not yet fully understood [5]. In resonantly perturbed toroidal magnetic plasmas, such as those in tokamaks, theory predicts that only three parameters are required to determine the topology of a fully reconnected magnetic state, namely the normalized plasma viscosity, rotation, and resistivity [6]. As discussed in this Letter, linear response studies of resonantly perturbed tokamak plasmas with the m3d-c 1 resistive MHD code [7] have resulted in the discovery of a new class of magnetic field line bifurcations that do not fit into the topological framework of this three-parameter bifurcation theory. These results are particularly important for understanding the underlying physics of selforganization in toroidal plasmas. The topology of the magnetic equilibrium field in toroidal plasmas plays an essential role in determining the properties of the energy, particle and momentum confinement as well as its MHD stability. Therefore, understanding the mechanisms responsible for this new class of topological bifurcations is of particular importance as the plasma parameters in high-power toroidal fusion devices approach those needed to achieve self-sustained burning or ignited states, since the ability of the plasma to self-organize can open access to new types of operating regimes.
Experience has shown that as the magnetic topology changes in tokamaks, unexpected types of plasma dynamics can appear. In particular, bifurcations of the equilibrium magnetic topology are known to be associated with the ability of the plasma to self-organize in unanticipated ways. For example, the spontaneous generation of a bifurcated helical magnetic core equilibria and saturated internal kinks, observed in several conventional aspect ratio tokamaks [8,9] and in a small aspect ratio spherical tokamak [10,11], are found to be associated with peaked pressure profiles that trigger bifurcated MHD equilibrium states reminiscent of saturated internal kink modes found in animec simulations [12]. In the RFX-mod reversed field pinch, it is found that as the plasma current is increased, the hyperbolic (x-) point of the dominant core magnetic island merges with the main magnetic axis to form a new self-organized quasisingle-helicity state [13]. Core magnetic islands have also been found to trigger internal transport barriers [14] and to spontaneously generate cyclical dynamics in the heat transport across low-order rational surfaces near the midradius of the DIII-D tokamak and the large helical device (LHD) [15]. Consequently, understanding the underlying physics responsible for triggering self-organized plasma states is of paramount importance for the development of magnetic fusion energy devices and is of intrinsic interest from a broader scientific perspective.
In magnetically confined toroidal plasmas, such as in tokamaks and stellarators, non-axisymmetric magnetic field perturbations result from vacuum field sources external to the plasma and from internal MHD plasma instabilities. As β increases, changes in the equilibrium magnetic topology, plasma stability and plasma transport become highly sensitive to the properties of the external non-axisymmetric vacuum field perturbations. Here, β is the ratio of the volume averaged plasma pressure p to the magnetic pressure B 2 /2µ 0 and the associated normalized plasma beta is defined as β N = β(aB T /I p ), where a is the minor radius of the plasma, B T is toroidal magnetic field and I P is the plasma current. Changes in the resonant and non-resonant magnetic equilibrium topology, caused by vacuum field perturbations, are referred to as the plasma response. These changes are known to alter neoclassical tearing mode (NTM) and edge localized mode (ELM) stability as well as the particle, energy and momentum confinement of the plasma [16]. Coupling between stable ideal kink modes and resonant modes on rational surfaces with m = nq have been shown to correlate with reductions in the edge pressure gradient resulting in changes to the ELM stability [17,18]. Here, m and n are poloidal and toroidal mode integers respectively, and q(ψ N ), a function of the normalized poloidal magnetic flux (ψ N ), is the safety factor associated with the radial variation in the normalized poloidal flux ψ N of the axisymmetric equilibrium magnetic field. Thus, the coupling between stable ideal kink modes, with (m,n) components that reside at a smaller radius than the resonant q(ψ N ) surface, can significantly modify the resonant plasma response to an external magnetic perturbation field.
In this Letter, the plasma response to n = 3 externally applied non-axisymmetric magnetic perturbation fields (δb ext = δb ext r + δb ext θ + δb ext φ ) is simulated in two high β N spherical tokamak discharges using the m3d-c 1 resistive MHD code. These simulations result in the formation of n = 3 heteroclinic magnetic islands on each low-order rational surfaces from the center of the plasma to the edge. As δb ext r is increased in the simulation, it is found that magnetic islands on resonant surfaces where m/3 is a rational number undergo internal topological bifurcations. In these simulations, the total resonant magnetic field in the plasma is δb total = δb ext + δb plasma . Here, δb plasma is the MHD plasma response field due to the application of δb ext . It is found that the threshold value of δb ext required for m/3 islands to undergo a bifurcation depends on β N . Furthermore, in cases where no plasma response is included, i.e., δb plasma = 0 no bifurcations are observed.
The linear time-independent single-and two-fluid m3d-c 1 simulations in this work were carried out using synthetically generated kinetic Grad-Shafranov [19] plasma equlibria in the NSTX-U spherical tokamak [20]. Two equlibria are simulated, one with β N = 5.5, I p = 1.5 MA, B T = 1.0 T and a q(ψ N ) on the 95% ψ N surface of q 95 = 8.7 and another with β N = 7.6, I p = 2.0, B T = 1.0 and q 95 = 5.5, where q(ψ N ) > 1 everywhere in both cases. Non-axisymmetric vacuum field perturbations are applied using 12 equally spaced toroidal rectangular coils above and below the equatorial plane, referred to as the Non-axisymmetric Control Coil (NCC). This coil set is being considered for instillation in NSTX-U over the next few years. In the simulation discussed here, the normalized currents in this coil have a toroidal distribution pattern in the upper NCC of (+1, 0, -1, 0, +1, -1, 0, +1, 0, -1, 0, 0) with opposite signs in the lower NCC. m3d-c 1 simulations are done with currents in the NCC (I N CC ) ranging from 1 kA to 30 kA using the radial plasma profiles shown in figure 1 as input. FIG. 1: a) ion (black) and electron (red) density, b) ion (black) and electron (red) temperature, c) thermal (black) and total (blue) pressure, and d) poloidal E xB rotation (ω ExB = 2πf ExB ) for β N = 7.6 (solid line) and β N = 5.5 (dashed line) cases.
Due to the relatively high poloidal E xB rotation in these discharges, shown in figure 1(d), the m3d-c 1 simulations show that δb plasma is reduced by at least a factor of 2 on the resonant surfaces of interest compared to that of δb ext . This results in small magnetic islands, with widths ∆ψ island N < 0.02, formed on each rational surface corresponding to the fundamental n = 3 toroidal modes produced by the NCC.
In the unperturbed δb ext = 0 simulations, field lines are confined to the surface of a 3D torus defined by ψ N , the normalized toroidal flux (χ N ), the poloidal angle (θ) and the toroidal angle (φ), where ψ N ,θ,χ N ,φ are referred to as action-angle variables. Using this representation, the action ψ N , when plotted as a function of the poloidal angle θ, describes a set of parallel straight field lines that intersect a 2D Poincaré plane at a φ = constant angle. On rational surfaces with δb ext = 0, field lines form periodic trajectories which result in a finite number of fixed-points when projected on the Poincaré plane. For example, on the q(ψ N ) = 3 surface individual field lines undergo 3 toroidal revolutions (φ = 6π) before returning to their original ψ N ,θ position. Here, each fixed-point undergoes 3 toroidal transits in order to complete 1 poloidal transit resulting in a winding number of m/n = 3/1. Since there are 3 set of these isolated islands on q(ψ N ) = 3, this results in an effective m/n winding number of 9/3. The application of a small δb ext field destroys the toroidal symmetry resulting in the formation of an even multiple of fixed-points as prescribed by the Poincaré-Birkhoff fixed-point theorem [21]. According to the Poincaré-Birkhoff fixed-point theorem, an even number of fixed-points must appear on every rational surface for an arbitrarily small δb ext . In conservative area or magnetic flux preserving (∇ · B = ∇ · ∇xA = 0) Hamiltonian systems, where each rational surface is bounded on both sided by an irrational surface, the KAM theorem [22], in conjunction with the Poincaré-Birkhoff fixed-point theorem mandates that each pair of fixed-points must consist of a both a hyperbolic and an elliptic fixed-points. In general, this allows for the formation of N f p = 2 m fixed-points, where m is the period of each fixed-point, is referred to as the topological index corresponding to the number of isolated fixed-point pairs and 2 is the Poincaré-Birkhoff integer. Since topological island bifurcations must preserve the fixed-point winding number on each resonant surface, = 1 for homoclinic islands and > 1 for heteroclinic islands. An example of a set of heteroclinic islands resulting from the m3d-c 1 β N = 5.5 simulation on the q(ψ N ) = 3 surface in the NSTX-U with I N CC = 2.0 kA is shown in figure 2.
This figure is obtained by integrating the nonlinear magnetic field line differential equations produced by δb total using the TRIP3D code [23]. Here, 3 isolated heteroclinic sets ( = 3) of period 3 (m = 3) islands are formed when the NCC is turned on, each of which is comprised of 3 elliptic fixed-points (o-points) and 3 hyperbolic fixed-points (x-points), resulting in N f p = 2 m = 18. Field lines originating on each fixed-point make m = 3 toroidal revolution before returning to their original ψ N , θ position. As shown by the red, blue and green color coding in figure 2, field line trajectories corresponding to one set of heteroclinic islands do not connect to either of the other two heteroclinic island sets. Figure 3 shows a detailed view of the q(ψ N ) = 3 MHD magnetic island with its elliptic fixed-point at θ = 179.9 • in figure 2 for increasing I N CC levels. The increase in I N CC causes a dramatic change in the topological structure of the inner flux surfaces of the island along with all of its heteroclinic partners and neighbors. Initially, ∆ψ island N increases with I N CC and the elliptic fixed-point moves from θ = 179.9 • to θ = 183.1 • while a sub-set of internal island flux surfaces stretch asymmetrically in the opposite direction toward the screened vacuum island elliptic point located at θ = 163.3 • , as shown in the I N CC = 3.0 kA panel of figure 3. This asymmetric internal flux surface stretching is consistent with a change in δb plasma when linearly superimposed on the increasing δb ext field, and is expected to involve a redistribution of the δb ext screening current. It is noted that the position of MHD island hyperbolic fixed-points remains constant during this process. As shown in the I N CC = 4.0 kA panel, the elongated lobe has bifurcated to form a new set of MHD island elliptic and hyperbolic fixed-points located at θ = 170.0 • and θ = 173.9 • respectively. We refer to this process as a heteroclinic MHD magnetic island bifurcation since it results in the formation of 3 new heteroclinic sets of isolated m = 3 islands, with b = 2 = 6 resulting in N f p = 2 b m = 36, as required to preserve the m/n = 3/1 winding number in each set of islands on this resonant surface. The bifurcation occurs at I N CC = 3.48 kA.
As I N CC continues to increase beyond 3.48 kA the original MHD island elliptic fixed-point moves to θ = 184. These fixed-points continue to move poloidally by small amounts as I N CC is increased to 14.7 kA where a second heteroclinic bifurcation takes place. This second bifurcation results in 3 new heteroclinic sets of isolated m=3 islands, with 2b = 1.5 b = 9 resulting in N f p = 2 2b m = 54. These heteroclinic bifurcations also occurs on each of the period 3 rational surfaces studied in this simulation starting at q(ψ N ) = 5/3 and cascading radially outward, down the pressure profile shown in figure 1c, as I N CC is increased. A similar process is observed during simulations of the β N = 7.6 case but in this case only one bifurcation occurs on the q(ψ N ) = 3 surface when I N CC passes through 4.23 kA as it is increased to 30 kA.
Although not visible in figure 3 due to the generation of a screening current by the MHD plasma response, the vacuum heteroclinic island due to δb ext is located inside the MHD island with its elliptic and hyperbolic fixed points on the q(ψ N ) = 3 surface at θ = 163.3 • and θ = 183.5 • respectively. Here, the hyperbolic fixed-point of the vacuum island is located at approximately the same position as the original elliptic point of the MHD island shown in figure 3 with I N CC = 1 kA. Our hypothesis is that the position of the vacuum island fixed-points along with the non-resonant kink aligned perturbation field influences the bifurcation dynamics of the MHD islands, as the NCC current is increased, by: 1) regulating the screening current needed to prevent the vacuum island from opening and 2) distorting/stretching the flux surfaces. Coupling between the non-resonant kink response, in the vicinity of each rational surface, and the resonant island response is believed to be the mechanism responsible for these bifurcations.
The shift between the vacuum and MHD island elliptic fixed-points is expected to result from a redistribution of the plasma current density inside the MHD islands. This is required to maintain the co-I p vacuum island screening current density filament (j −vac ) while generating a counter-I p current density filament (−j −M HD ) needed to open the MHD island elliptic points. As I N CC increases, j −vac must increase in order to maintain a screening current consistent with ω ExB . A secondary, counter-I p , −j −M HD filament located in the elongated flux surfaces between the screened vacuum and original MHD island elliptic fixed-points is formed as I N CC increases. Eventually, δb plasma due to the secondary −j −M HD filament cancels δb plasma from the original −j −M HD filament at a point on the q(ψ) = 3 surface between the two −j −M HD filaments and a new set of hyperbolic fixed-points with their associated Poincaré-Birkhoff elliptic fixed-points are formed. It is also found that linear two-fluid m3d-c 1 simulations produce qualitatively equivalent results to those of the single-fluid simulations indicating that the Hall term in the m3dc 1 Ohm's law equation does not affect the bifurcation process.
In addition to satisfying the Poincaré-Birkhoff fixedpoint theorem, linear m3d-c 1 simulations demonstrate that the dynamics of these heteroclinic island bifurcations are the manifestation of a smooth/continuous, invertible, process controlled by δb ext that does not result in bistable or hysteresis-like behavior associated with various types of tearing modes [24]. This behavior is consistent with degenerate equilibrium points associated with local saddle-node bifurcations that are found in Hamiltonian systems composed of conservative vector fields [25]. In addition, the process involved in the stretching of the internal MHD island flux surfaces, leading up to the these heteroclinic bifurcations, is consistent with a diffeomorphic type of evolution found in a conservative dynamical system where subsets of differentiable manifolds undergo smooth one-to-one, invertible, topological transformations as a control parameter is varied.
In summary, a new class of heteroclinic magnetic island bifurcations has been identified during variations of an externally applied non-axisymmetric perturbation field in linear single-and two-fluid m3d-c 1 resistive MHD simulations of two spherical tokamak plasmas. These bifurcations preserve the rational flux surface winding numbers by generating new sets of isolated (heteroclinic) magnetic islands inside existing heteroclinic MHD magnetic islands starting on rational surfaces near the center of the plasma and cascading outward as the perturbation field δb ext is increased. Simulations with different normalized plasma pressures (β N ) result in similar bifurcation sequences but with different values of the δb ext bifurcation parameter. As β N is increased, bifurcations occur at larger δb ext values. Based on experimental results during the presence of magnetic islands in tokamaks and stellarators, bifurcations of this type can be expected to result in significant plasma transport and stability changes as well as possible observations of new types of self-organized behaviors.