Effect of injected flux and current temporal phasing on self-organization in the HIT-SI3 experiment

An active topic of research in the development of magnetic confinement fusion devices is the role of non-axisymmetric perturbations in the confinement and stability properties of axisymmetric equilibria.^1,2 For self-organized configurations in particular, the mode spectrum used to inject energy and helicity is expected to play an important role by setting the initial conditions for energy flow into the sustained equilibrium.^3,4

The HIT-SI3 device at the University of Washington forms and sustains spheromak equilibria through the use of Steady Inductive Helicity Injection (SIHI).^5,6 SIHI is performed using a set of three semi-toroidal injectors that drive oscillating helical plasma states that connect through a central flux-conserving confinement volume. Each injector is driven by two orthogonal coils, which drive a locally toroidal magnetic flux and plasma current through its semi-toroid, connecting back through the central volume. The entire device is constructed from chromium copper coated with insulating alumina on the plasma facing surface, to enforce the boundary conditions $B⊥=0$ and $J⊥=0$ for the duration of a discharge. Figure 1 illustrates the interaction among the non-axisymmetric injector fields, the axisymmetric spheromak equilibrium, and the coils used to drive each injector. The confinement region of the device has a diameter of 1.1 meters, with a magnetic axis radius of 34 cm and a minor radius of 23 cm.

One important detail of the SIHI system is that multiple injectors, each driven with an oscillating waveform, can be phased to provide a constant total helicity injection rate in steady-state while retaining purely inductive operation. The helicity injection rate of a single injector, when the two coils operate in temporal phase with one another, is given as $K̇inj=2Vinjψinj$⁠, where $[Vinj,ψinj]=[V0,ψ0] sin(ωinjt)$ are the injector loop voltage and magnetic flux oscillated at a given frequency ω_inj. By driving V_inj and ψ_inj on a single injector in temporal phase, sign-constant helicity injection can be performed. The three injectors, which have been named A, B, and C, can be operated with a phase shift from one another to vary the spatial and temporal properties of the magnetic helicity injection, with an injection rate of

where $ϕi$ is the temporal phase of a given injector.

All HIT-SI3 discharges follow similar time-evolution of bulk plasma parameters largely set by magnetic helicity conservation.⁷ A discharge initially begins with the injectors driving a 3D state, with a symmetry matching the injector drive, that after a short time undergoes a global relaxation event to form the spheromak equilibrium—the minimum energy state of the system.⁸ During this fast relaxation event, magnetic helicity, which decays on a resistive timescale (⁠$K̇loss=K/τK, τK∼τL/R$⁠), is constant, and the initial amplitude of the spheromak state is set by the initial helicity inventory. The helicity injectors continue to operate, with a steady-state being reached where helicity dissipation balances the helicity injection rate,

In the past work,^6,9 discharges were studied in three specific operating cases, corresponding to the two constant helicity injection configurations, with phases spaced 60° and 120° apart, and the all-in-phase case. The recent implementation of a new GPU-based feedback control system¹⁰ has allowed more precise control over the temporal phase of the HIT-SI3 injectors, which enabled detailed study of the available perturbation spectra.

A scan of the parameter space, using both experimental and simulated discharges, has been performed. This paper begins by detailing the process used to both collect and analyze experimental data and simulation results. Next, the results of the scan, characterizing the dependence of perturbation spectra and spheromak evolution with respect to the chosen injector phasing, are presented and discussed.

A series of discharges were taken with the HIT-SI3 device to scan the parameter space of injector relative temporal phasing. Each discharge, which all follow the same general evolution, consists of four distinct phases. First, when plasma breakdown occurs, the helicity injectors drive primarily 3D currents and flux throughout the volume—with a structure dictated by the injector geometry and phasing. Second, on a timescale comparable to the toroidal Alfvén wave transit time, the magnetic structure in the plasma undergoes a relaxation event, where a minimum energy state is formed while conserving magnetic helicity. This relaxation event converts the magnetic field from primarily non-axisymmetric fields to primarily axisymmetric fields. Third, following the relaxation event, a steady-state period begins, where the bulk parameters of the plasma evolve to a steady-state in which the helicity injection is balanced by helicity decay.⁷ Finally, when the injector operation is terminated, the spheromak equilibrium decays on a resistive timescale to complete the discharge.

The study presented here focuses on the steady-state sustainment phase of the discharge. A recently installed GPU-based feedback control system is capable of driving the injectors at a specified temporal phase to within 10° precision. Because all of the injector fields are periodic, and only the steady-state period of the discharge is being considered, the absolute starting phase of the system is arbitrary. Thus, the available parameter space is 2D, defined by the phase of injectors B and C for $ϕB$ and $ϕC$⁠, respectively, relative to a fixed phase of 0° for the A injector. A set of 324 experimental discharges varying both parameters (⁠$ϕB$ and $ϕC$⁠) using a fixed interval of 20°. To prevent long-timescale operating trends, such as wall conditioning, from affecting the data, the order of the discharges taken was randomized.

In order to develop sufficient statistics with comparable plasmas, lower power discharges with poorer wall conditions were used when compared with the high performance discharges studied in past works.^6,11–13 The injectors were operated with an injector power of roughly 7 MW and at a drive frequency of 16.4 kHz. Typical plasma densities for these discharges were $n≈4×1019$ m^–3 measured with the two-color, multi-chord interferometry system.

The mode content of the magnetic equilibrium and perturbations in HIT-SI3 are obtained through a set of surface magnetic probes,¹⁴ as shown in Fig. 2. From this set, 32 probes, indicated as L05 and L06 in the poloidal cross section in Fig. 2, consisting of two sets of 16 probes spaced on the outboard midplane of the device, are used to Fourier decompose the toroidal magnetic field up to n = 7. The toroidal field was selected to study because these probes are located in a recessed gap on the confinement volume, which affects the poloidal field in this region, which is in the cylindrical radial direction at this location. Similarly, the non-circular shape of the poloidal domain complicates interpreting the spectrum of poloidal decompositions of the magnetic structures.^15,16 The measured values for perturbed fields can be compared between different discharges with varying equilibrium amplitudes by normalizing to the equilibrium field, $δB/B0$⁠, with B₀ computed from the average of the four L04 probes in the same method as used in Ref. 6. Typical equilibrium field strengths obtained in the discharges studied were $B0=23±3$ mT.

As an example of the time evolution of the fields in a discharge and the dependence on the injector temporal phases, Fig. 3 illustrates a discharge where the drive parameters were changed halfway through. While the non-axisymmetric mode structure of the injector drive switched from primarily n = 1 to n = 3, the spheromak persisted through the transition with a small change to amplitude and character of the oscillation in the toroidal plasma current.

Composite equilibria formed from combinations of “Taylor states” with vacuum and homogeneous flux boundary conditions have been used in the past as approximations of the magnetic structures in the HIT-SI⁵ and HIT-SI3 devices.¹⁷ While the details are described in depth in Refs. 17 and 18, a short description is provided here. A “Taylor state” refers to the minimum energy equilibria that self-organized plasmas, such as spheromaks, tend to relax toward,^8,19 being described by

with a uniform value of λ being the eigenvalue of the system. The PSI-Tet code, which is capable of solving for these eigenstates on a 3D unstructured tetrahedral mesh for arbitrary geometries, is used to locate the solution for four distinct cases. The first case is the isolated configuration, where the three injector fluxes are held at zero (homogenous boundary conditions), and corresponds to the minimum energy spheromak state of the main confinement volume, which is found by solving for the smallest eigenvalue of Eq. (3)—setting the value of λ used in the remaining solutions. The remaining three cases all solve for the state where one of the individual injectors has unit flux, and the remaining injectors plus the axisymmetric flux in the confinement volume are all zero. The purely non-axisymmetric structure of these injector states ensures that they do not couple to the axisymmetric eigenstate, as with symmetric vacuum fields,^20,21 allowing these states to be combined with λ given by the lowest eigenvalue.

Linear combinations of these four equilibria, possible by the use of the same eigenvalue λ, allow for approximations of the 3D structure of full equilibria, including injector perturbations, to be considered. This takes the form of

where B_i is the corresponding equilibrium state and $Ψi$ are the scaling coefficients, which are equal to the toroidal and injector fluxes. The temporal phase dependence of the perturbations can be computed by scaling these coefficients in a similar way, such as $ΨA(t)=ΨA,0 sin(ωt+ϕA)$⁠. An example of such a combination is seen in Fig. 1(a).

Alongside the experimental study, a series of eXtended MagnetoHydroDynamic (xMHD) simulations spanning the parameter space were performed using the NIMROD code.²² NIMROD has been used extensively to model both the HIT-SI^23–26 and HIT-SI3⁹ devices. For simulations in this work, a zero-β Hall MHD model was used, with the application previously described in Ref. 24. This model assumes a uniform and constant temperature and density, solving the following equations for momentum and magnetic field:

where model closures are given by

A deuterium plasma (⁠$mi=3.34×10−27$ kg) is used with resistivity $η=9.4×10−6$ Ω-m, corresponding to an electron temperature T_e = 6 eV and density $n=2×1019$ m^–3 that were selected based on machine operations prior to the simulations being performed. The viscosity ν = 550 m²/s was selected based on values used in prior work with this model, and work in Refs. 24 and 9 found that the bulk dynamics are insensitive to this value. Due to minor changes in the injector power circuits between when the simulations were performed and the experimental campaign, the simulations were performed with an injector frequency of f_inj = 15.6 kHz. While these operating frequencies differ, they are close enough that similar plasma behavior is expected, with major differences occurring in device operation at $finj>40 kHz$⁠.¹¹ Additionally, similar to the work in Refs. 9 and 24, the electron inertial term is included for numerical stability considerations with a factor (⁠$fme$⁠) 100 times the physical value.

The computational geometry used in NIMROD has two significant approximations of the experimental configuration. First, the insulating layer on the inner surface of the experimental flux conserver is approximated by a thin (1 mm) region of high resistivity (⁠$ηwall/ηplasma∼105$⁠). Second, the toroidal dependence of the solution fields is decomposed as a finite Fourier series, and as such is unable to accurately represent the geometry inside of the injectors. Instead, the injectors are implemented as boundary conditions on the axisymmetric, bow-tie confinement volume, and these boundary conditions are equivalent to those described in detail in reference24 for the older HIT-SI configuration. Figure 4 shows the normal magnetic fields applied to approximate the HIT-SI3 injector configuration for a state with all injectors active and in-phase.

A set of 324 simulations, using 20° increments in the temporal phases of the B and C injectors, were completed to create a comparison with the experimental results. The simulations were run to a physical time of 600 μs to allow the system to relax close to the final state (⁠$K̇=0$⁠), and analysis was primarily done in the final injector period of the simulations.

While each experimental discharge lasted approximately 3.2 ms from plasma breakdown to injector shutoff, the focus of this study is the steady-state behavior observed late in time in the discharges, with bulk evolution of the spheromak and injector parameters shown in Fig. 5. The bulk-injector parameters are computed as the injector-period-averaged quadrature sum of the individual components, for example,

where T is the injector period. The two controlled parameters for each injector are the voltage and flux, with the settings for these discharges designed to maintain a constant voltage, while the injector flux is increased throughout the discharge. To extract these results from each diagnostic's temporal waveform, the average value from the final nine injector cycles (⁠$Δt∼548 μ$s) prior to injector shutoff was computed. This period of time in the discharges sees steady behavior, with all signals having primarily a static DC component and constant amplitude oscillations at f_inj. Additionally, it is observed that spheromaks with both orientations of toroidal current are formed; this is typical for operation of the HIT-SI family of devices because there is no externally applied field or other property with a preference for spheromak orientation. In this dataset, it was observed that the spheromak toroidal current tended to be oriented in the “negative” toroidal direction, encompassing approximately 83% of discharges.

Additionally, an extension to the definition of injector current (I_inj) used to calculate current gain (⁠$G≡ItorIinj$⁠) in past analysis of the HIT-SI3 device is required to compare cases where the helicity injection rate is not constant. Previous work used $Iinj(t)=IA2(t)+IB2(t)+IC2(t)$⁠, which, for the constant helicity injection case, produces a waveform with a strictly DC component, but for other cases produces a waveform oscillating at twice the injector frequency. To extend this definition, we instead use $Iinj=16(I0,A+I0,B+I0,C)$⁠, where the values being summed are the amplitudes of the sinusoidal waveforms, which are obtained by fitting the waveform to a sine-wave of the form $IA=I0,A sin(ωt+ϕA)+CA$ in the time window that the data analysis was performed. This form allows direct comparison with prior work by obtaining identical answers for cases studied in those past studies.^6,9

The simulation data are processed in a similar method, instead with only the final injector cycle being used to construct the steady-state values. For comparisons with experimental results, the NIMROD simulations and PSI-Tet composite equilibria are sampled using the experimental locations of the surface magnetic probes and used to calculate the toroidal plasma current with the same method as the experiment. The composite equilibria are unable to produce a prediction of the current gain, since it is an input to the model.

The two computational models (NIMROD and composite equilibria) can also be used to compare the volume-integrated energy spectrum, providing a quantification of how indicative measured surface probes activity is of bulk plasma dynamics for a given model. Individual mode energies are computed as

which can be converted into a comparable form to the surface probe modes with $δBi/B=WB,i/WB,0$⁠, where i denotes the toroidal mode number.

The primary experimental result of this work was the generation of a map detailing both the available perturbation spectra on the HIT-SI3 device and the performance considerations, which is seen in Figs. 6–11. These figures contain the results from the experiment, the composite Taylor states, and both methods of measurement for the non-linear xMHD simulations.

The current gain, which serves as the primary performance metric for device operation, indicates a $≈30%$ in variation between the peak (⁠$≈1.8$⁠) and minimum (⁠$≈1.3$⁠) values observed in the dataset. These values are approximately half the values achieved (⁠$G≈3$⁠) in the high-performance operation mode,⁶ with the main factor driving this difference being reduced values of j/n resulting from poorer wall conditioning.¹³ 

Mode activity for the higher order modes is not shown, as the amplitude drops considerably for n > 3, as shown in Fig. 12. It is clear that the peak amplitude seen for each toroidal Fourier mode for a variety of measurement methods occurs for $n<=3$⁠. This lower mode activity, in particular n = 1 activity, is the result of directly applied higher-order structure coupling to these lower-n structures through dynamic non-linear interactions, in the case of NIMROD, and static linear interactions, in the case of the composite Taylor model. In particular, the boundary conditions applied in the NIMROD simulations contain significant amounts of n > 5 structure, due to the small toroidal footprint (⁠$40°$⁠) of the injector mouths. While the nonlinear model sees the formation of lower-n structures through dynamic processes, we note that the spectrum of the composite equilibrium is produced through linear interactions with eigenmodes of the homogeneous force-free system.

Overall, the simulations do a reasonable job of generating a map of the available perturbations, suggesting that the models are capable of capturing the primary non-linear interactions of the injector perturbations.

From Figs. 7–11, a variety of observations and interpretations of the data can be made.

Beginning with the experimental data, it is clear the applied perturbation spectra of the helicity injection system are both well-described and controllable by the 2D space created by the relative temporal phase of the injectors. Additionally, all possible perturbation configurations on HIT-SI3 lead to spheromak formation and sustainment, consistent with models used in the past that described the bulk evolution of the spheromak in terms of global magnetic helicity conservation.^7,12 Previous studies additionally focused on parameter regions that are places of interest on these maps, with the peak of n = 1 and n = 2 occurring at 120° phase, peak n = 3 occurring when the injectors are in phase, and equal parts of all three modes occurring with 60° phase spacing.

Past studies that examined differences in injector-phasing on the measured results were performed using the assumption that permuting the values of $ϕB$ and $ϕC$ would produce identical results. However, these experimental results suggest that this assumption was invalid, and the resulting perturbation spectrum depends on the rotation direction of the injector fields. While two orientations of spheromaks can be produced by the device, set by the direction of the toroidal plasma current, the resulting perturbation mode structure does not vary between the two. The most likely source of the symmetry break in the system is the vertical asymmetry caused by the injectors all being located on the top of the device.

The rotation asymmetry serves as a good comparison point between the two theoretical models that are used. As seen in Figs. 8 and 9, the composite equilibrium model does not capture this observation, while the dynamic simulations do. In an effort to further isolate the physics in the xMHD simulations that produces this asymmetry, a set of the two 120° phased cases were reproduced while removing the Hall term $(J×B/ne)$⁠, and the effect was still observed. Additionally, it was noted, while performing these simulations, that HIT-SI3 is not nearly as dependent on including the Hall term to reproduce experimentally observed phenomena, as was the case with the HIT-SI configuration.^9,24 While the full explanation is a subject of future work, we note that the reason to include the Hall term found by Akcay et al.²⁴ was to provide a symmetry break for the purely n = 1 perturbation, which is not necessary for HIT-SI3 due to its mixed-mode injection system.

Despite the limitations of the composite-equilibrium model on reproducing the full experimental results, it still provides a useful starting point for modeling this style of device due to the lower computational cost. While a simple 1D model, such as the boundary conditions, we apply in NIMROD simulations, at the injector-spheromak interface would incorrectly see high-order (n > 5) structure on the applied perturbation, the composite-equilibrium is able to see the volumetric lower-order structure that appears in the force-free time-invariant limit. This illustrates that any model attempting to predict the mapping of applied perturbation spectra on a similar inductive helicity injection system requires the full 3D spatial representation of the system.

The difference in results between the surface probe and energy diagnostics seen on the NIMROD simulations in Figs. 9 and 11 provides some additional details in the decomposition of these modal structures in the poloidal plane. A significant reduction in the amount of rotation asymmetry is seen on the volume-integrated diagnostic, suggesting that the perturbations closer to the core of the device do not see the same asymmetry that the edge perturbations do. This serves as an illustration of the limitations in measurements of the dynamics at play in the HIT-SI3 device, as core measurements of the 3D magnetic perturbations are not available in the diagnostic set. In contrast to the NIMROD results, the composite equilibrium model sees qualitatively similar spectra with the two diagnostic methods, as seen in Figs. 8 and 10, with differences primarily seen in the amplitudes. This provides further evidence that dynamic processes play an important role in the interactions between the helicity injection system and the spheromak equilibrium.

While this dataset was taken with the device in a reduced performance operating mode to obtain the desired quantity of data, the measurements of current gain point to some previously unexplored operating points in past high-performance campaigns. The gain is observed to be minimized when the injectors are operated in the n = 3 configuration. While HIT-SI3's n = 3 configuration is additionally one of the cases where the helicity injection rate (⁠$K̇tot$⁠) varies the most throughout the injector cycle, we see that the peak current amplification appears near the other case but instead has more n = 1 and n = 2 activity, which is the positions [180,180], [180,0], and [0180] on Fig. 6(a). This peak current gain region has not yet been explored in higher performance operation (G > 2.5), and thus, the question if this result scales up to higher performance is currently unanswered. It is noted that the region of elevated performance (G > 1.7) includes the operational configuration (⁠$ϕB=240°,ϕC=120°$⁠) that previous high performance campaigns had identified as peak performance.⁶ 

While the simulations may capture some of the variation in current amplification, the level of disagreement seen is consistent with past results studying the importance of finite-β effects in simulations of HIT-SI and HIT-SI3. Past work^9,25,26 found that toroidal current gain has a strong dependence on the electron temperature, which sets the plasma resistivity and plasma density, which influences the magnitude of the dynamo terms that provide the majority of current drive. The zero-β model used did not allow for either of these terms to change from their chosen values, failing to capture any effect their evolution may have on these simulations. While finite-β models have been successfully applied to simulations of HIT-SI3,^9,25 they tend to require significantly more computational resources making them prohibitive for the hundreds of simulations performed in this study.

This study demonstrates the range of operating space of the HIT-SI3 device, where spheromak equilibria can be generated with a wide range of non-axisymmetric perturbation spectra. In particular, these results show that the operating space is larger than previously thought, providing incentive for future modeling efforts to explore the space. Additionally, we have demonstrated that non-linear extended MHD models are capable of replicating the structure and asymmetry of these perturbations, validating these models toward use in predictions of other injector configurations. Finally, by comparing the results obtained from a range of models, we were able to demonstrate the role of non-linearities in models of the device. While the cascade of energy to lower mode numbers requires only spatial non-linearity, the symmetry break of the injector system requires the full time-dependent system to be evaluated. Future work will involve applying these same techniques to the now operational HIT-SIU device, which uses a fully connected helicity injector manifold system to provide more flexibility in the available perturbation spectra. Broad maps of the available parameter space will inform the regions of the parameter space to identify promising regions for experimental campaigns.

The authors acknowledge the past and present members of the HIT research team for valuable discussion in performing this research and for the creation of diagnostic systems used in this work. The information, data, or work presented herein were funded in part by the Advanced Research Projects Agency-Energy (ARPA-E), the U.S. Department of Energy under Award No. DE-AR0001266, and by CTFusion, Inc., the primary recipient of ARPA-E Award No. DE-AR0001098.

Two of the authors of this work (K.D.M. and A.C.H.) have a financial interest in the prime awardee and funder of this work, CTFusion, Inc.

The data that support the findings of this study are available from the corresponding author upon reasonable request.