New type of self-sustained divertor oscillation driven by magnetic island dynamics in Large Helical Device

A new self-sustained divertor oscillation is discovered in magnetic island induced detached plasmas in the Large Helical Device. The divertor oscillation is found to be a self-regulation of the width of an edge magnetic island accompanied by detachment-attachment transitions. The modified Rutherford equation combined with an ad-hoc bootstrap current equation is introduced to describe the divertor oscillation as a predator–prey model between the magnetic island width and a remnant X-point bootstrap current. The model successfully reproduces the experimental observations in terms of the oscillation frequency, the phase relation between variables, and the oscillation amplitude.


Introduction
Divertor detachment operation, in which a radiative boundary plasma in front of divertor target mitigates direct heat flux, is a viable solution for reconciling a high core plasma performance and a tolerable heat load onto divertor targets.In order to acquire a physics-based control manner of the detachment operation, it is essential to investigate the background physics of detachment.Experimental and modeling efforts have been devoted for uncovering key factors that determining detachment conditions [1][2][3][4].As a result of intrinsic nonlinearity, recycling driven self-sustained divertor oscillations (SSDOs) are occasionally observed in tokamaks [5,6], which can provide an opportunity to untangle the background nonlinear physics.Recently, a new type of SSDO was discovered in the Large Helical Device (LHD) stellarator [7], where the peripheral plasma was detached from diverters by means of externally applied perturbation fields (PF) [8].Not only for untangling the nonlinear background physics of the detachment transition but also exploring a possibility of utilizing external PF as a control knob for the detachment operation [3,[9][10][11], in-depth investigation of the PF-driven SSDO in LHD is highly worthwhile.In this paper, detailed experimental observations of the new SSDO are introduced, and a phenomenological modeling accounting penetration and shielding dynamics of the externally applied magnetic island [12,13] is presented.This model potentially explains the SSDO observed in a tokamak island divertor operation [14].

Experimental observation
LHD is a superconducting helical fusion device with toroidal and poloidal magnetic field periods of 10 and 2, respectively.Plasma having an elliptically shaped cross-section is surrounded by a double-nulled scrape-off layer, in which stochastic fields are occasionally enhanced.The present experiment is conducted in the so-called outward shifted magnetic configuration with a vacuum magnetic axis position of R ax = 3.9 m and a toroidal magnetic field strength of B ϕ = 2.54 T (counterclockwise direction viewed from above).A thick magnetic stochastic region at the periphery is a characteristic feature of this configuration, in which m/n = 1/1 magnetic island is driven by external PF, where m and n denote the poloidal and toroidal mode numbers, respectively.The m/n = 1/1 magnetic island leads to enhanced radiation losses at the magnetic island X-point while maintaining the core plasma parameters in the detachment operation [3].
Figure 1 shows the time evolution of the plasma parameters for a discharge having the SSDO (# 163369).The plasma is sustained by three tangentially injected neutral beam (NB) heatings, which predominantly drive the NB current with I p < 0 (counter-clockwise direction, i.e. co-direction to the toroidal magnetic field).A density ramp-up operation up to ne ∼ 6 × 10 19 m −3 is performed by gas puffing to find the detachment condition, which also leads to an adiabatic decay of the electron temperature, (see figures 1(a), (e) and (f )).Note that there is an apparent in-out asymmetry in the density profile, which likely originates from unphysical reasons.The initial detachment transition occurs at t ∼ 3.62 s and ne ∼ 4 × 10 19 m −3 , where the divertor ion saturation current reduction and the X-point radiation enhancement without any deterioration of the plasma stored energy are seen, as summarized in [3]. Figure 1(d) shows the amplitude and phase of the m/n = 1/1 radial magnetic field structure obtained by an on-vessel poloidal saddle loop coil array [13].Amplitude of the poloidal saddle loop coil array and its vacuum value are denoted as Φ and Φ ext , respectively.The square root of the relative amplitude corresponds to the normalized magnetic island width, i.e.W/δ = √ Φ/Φ ext , where δ is the external PF driven magnetic island width in a vacuum.The detachment transition occurs when the magnetic island width exceeds its vacuum value, W/δ > 1, i.e. √ Φ/Φ ext > 1. Immediately after the transition, a transient shrinkage of the island is observed, which can be regarded as a back reaction of the detachment transition.Note that the phase of the m/n = 1/1 radial magnetic field structure is maintained during the entire discharge.
Figure 2 shows the time evolution of the plasma parameters during the SSDO phase.The SSDO of ∼40 Hz appears when the magnetic island width extends above 1.05 times the vacuum PF magnetic island width.The SSDO is observable in multiple quantities, such as X-point radiation loss P X rad , divertor ion saturation current I Div is , normalized magnetic island width W/δ = √ Φ/Φ ext , and plasma current variation Îp .Here, Îp is given by I p − ⟨I p ⟩, where ⟨I p ⟩ is the slowly varying trend obtained by a low-pass filter.No movement of the divertor strike-point is seen, as shown in figure 2(b).Figure 3 shows the evolution of Lissajous diagrams of √ Φ/Φ ext versus P X rad and √ Φ/Φ ext versus Îp in left and right columns, respectively.Lissajous diagrams in the initial detachment transition phase are shown in the top row.Once the magnetic island width reaches a certain value, the X-point radiation loss sharply rises, which corresponds to the detachment transition.The transient shrinkage of the magnetic island follows the detachment transition, while, no back transition is observed in this period.The second and third rows show the Lissajous diagrams during the SSDO phase.After the detachment transition (increasing the radiation), the magnetic island width turns to decay, during which the plasma current is more or less maintained.The shrinkage of the magnetic island leads to a gradual decay of the X-point radiation loss, that is the attachment back transition.Immediately after the back transition, the plasma current surges, which revitalizes the magnetic island growth making a closed loop.The SSDO ceases suddenly with a slight decay of the loop extent only in the last revolution after six cycles from the beginning.
It is worthwhile to examine the possibility of the recyclingdriven SSDO [15], which is considered to be a possible explanation for the SSDO observed in tokamaks [5,6].The key physics of the recycling-driven SSDO is that there is a multi-valued relation between temperature and particle number in the divertor region, on which a repetition of transitions and back-transitions is possible.In this model, it is predicted that the divertor plasma temperature oscillation is the out of phase with respect to the oscillation in neutral pressure and/or plasma density.Time evolutions of those quantities are shown in figure 4. The neutral pressure in the private flux region [16] and the electron temperature and density in the divertor region are all in phase, and therefore the recycling driven bifurcation model [15] may not be relevant for the present observation.Moreover, from an ASDEX-Upgrade observation, the recycling-driven divertor oscillation was pointed out to be electrostatic [17].Here, the oscillation involves a magnetic activity.In the end, the SSDO observed in LHD is considered to originate in upstream plasma.
Another possibility is the SSDO driven by the impurity radiation-condensation instability [18].It is possible that some radiation related mechanisms may play a role in the SSDO in LHD, as there is a phase shift between the radiation intensity and other plasma parameters in figure 4. Detailed assessments with spatially resolved profiles of the radiation and the impurity emission will be presented in future.In this paper, a model of the SSDO related to the magnetic island dynamics will be discussed.
Plasma profile evolution across the SSDO is shown in figure 5. Figure 5(a) is the time evolution of the averaged plasma minor radius a 99 , which is defined by the effective minor radius in which 99% of the electron stored energy is confined and is calculated according to Thomson scattering data.After the initial detachment transition, an abrupt shrinkage of a 99 is observed, which is due to the edge electron temperature decay, as shown in figure 5(b).The repetition rate of the Thomson scattering laser injection is too sparse to resolve the profile evolution in the SSDO.Radial profiles of electron temperature T e , ion collisionality ν * i , and inverse mean free path ν ii /V ti (ν ii and V ti are the ion-ion collision frequency and the ion thermal velocity, respectively) in fully detached and attached phases are shown in figures 5(c)-(e), respectively.The rotational transform profile ι/2π ≡ 1/q, where q is the safety factor, is overlayed in figure 5(c).Around the ι/2π = 1 rational surface, a flattening of the electron temperature profile is visible in the attached phase, which is a signature of the magnetic island externally driven by the PF.After the detachment transition, the electron temperature inside the magnetic island further decreases and becomes marginal to the diagnostics limitation, leading to increases in the ion collisionality and the inverse mean free path.

Phenomenological modeling
The SSDO is regarded as an interplay between the magnetic island and the bootstrap current at the X-point.In LHD, the direction of the bootstrap current is theoretically predicted to depend on plasma parameters.In particular, when the inverse mean free path is not long enough, the bootstrap current is directed in such a way that the rotational transform is reduced [19].In this case, the bootstrap current can be a source of magnetic island expansion [20], analogous to tokamaks with normal shear.The magnetic island expansion by the bootstrap current is illustrated in figure 6.If there is no seed island, the bootstrap current is poloidally homogeneous, since the pressure gradient is a magnetic-flux function (figure 6(a)).Once the seed island forms, the pressure profile is flattened at the O-point region.As a result, the poloidal symmetry in the bootstrap current distribution is violated and a poloidal modulation in the radial magnetic field is induced.If the direction of the bootstrap current is tokamak-like, this radial magnetic field is directed in such a way that the seed island is expanded (figure 6(b)).This is basically what the modified Rutherford equation says [21][22][23][24].When the magnetic island width exceeds a threshold, the detachment transition is triggered by the increased X-point radiation loss, due to the reduced temperature.This gives a rise in plasma collisionarity at the X-point, which results in a decay of the X-point remnant bootstrap current, and the magnetic island starts to shrink.The plasma eventually returns to the attached state with a reduced magnetic island width, which closes the loop.Here, key quantities are the magnetic island width and the X-point remnant bootstrap current.
For testing this concept, a phenomenological model is developed.The modified Rutherford equation [21] is combined with an ad-hoc bootstrap current equation as This set of equations potentially expresses predator-prey dynamics between magnetic island width W and X-point bootstrap current j BS .Here, vacuum PF island width δ, characteristic magnetic diffusion speed V M = |∆ ′ |r 2 s τ −1 R , and coupling coefficient C = 2r 2 s µ 0 |L q |τ −1 R ⟨B p ⟩ −1 are experimentally determined.Parameters in those terms are given as follows: the minor radius of rational surface r s = 0.52 m, resistive magnetic diffusion time τ R = µ 0 r 2 s /η NC where η NC is neoclassical resistivity, the scale length of the safety factor L q = q(dq/dr) −1 , and flux surface averaged poloidal magnetic field ⟨B p ⟩.The tearing stability parameter |∆ ′ | ∼ 2.4 m −1 is calculated from the magnetic flux profile, which is estimated by integrating the MHD equation [22] in two paths, i.e. from the edge and from the core to the rational surface.The edge boundary condition is determined by the magnetic field measurement in a reference discharge without the external PF application.The core boundary condition and the current profile shape are confirmed to have a negligible impact on |∆ ′ | by a wide parameter scan, because the target rational surface is close to the very edge.The choice of bootstrap current equation (equation ( 2)) is rather highly arbitrary.Here, it is given to capture fundamental physics and being as simple as possible, having the linear growth term with growth rate α and a nonlinear saturation term with the saturation coefficient of βW.The shape of equation ( 2) is identical to a constituent of the classical Lotka-Volterra predator-prey model.More realistic modeling is anticipated in the future step.The value of α is given by the inverse of the parallel diffusion time as ti /ν ii being the ambipolar diffusion coefficient, according to a gyrokinetic particle simulation for neoclassical transport [25].The value of β remains undetermined here.
A perturbative expansion is performed with W = W 0 + W 1 and j = j 0 + j 1 (the subscript 'BS' is now dropped for simplicity), where the terms with subscripts 0 and 1 correspond to mean quantities and small oscillatory quantities, respectively.One of the fixed points is obtained as ).From the observation, the equilibrium island width W 0 can be determined as 1.05 × δ (see figure 2(c)) and β = α/W 0 is now given.The first order equation gives the system eigenvalues as, Here, the real-part and sign-inverted imaginary-part of the eigenvalues correspond to the system frequency and damping rate, respectively.For the present experimental parameters, an oscillatory solution with a frequency of 20 Hz and a smaller damping ratio of 8.5 s −1 is found.The predicted system frequency is in the same order as the observed SSDO frequency of 40 Hz, and is sufficiently larger than the damping rate, which can explain a quasi-continuous oscillation.
The behavior of this model is analyzed in figure 7 in detail by changing the parameters δ, W 0 , and α.Open circles in figure 7 show the experimental values.Figures 7(a) and (b) show δ and W 0 dependence of the system frequency and the damping rate, respectively, where α is fixed to the experimental value.When W 0 < δ, a linear growth of the magnetic island is predicted with zero real-frequency.Oscillatory solutions appear when the island width exceeds its vacuum value, i.e. δ < W 0 .When the system frequency is larger than the damping rate, quasi-continuous oscillation is possible.Figure 7(c) shows the slices of figures 7(a) and (b) at δ = 0.07 m of the experimental value.The damping rate decays as the magnetic island expands, in which the system frequency is more or less maintained.Dependence of the system frequency and the damping rate on α is shown in figure 7(d).In the region α > 10 2 s −1 , a finite frequency appears, while in the region α 3 × 10 3 s −1 , a frequency larger than the damping rate is obtained.For a quasi-continuous oscillating solution, large α, i.e, a sufficiently quick response of the bootstrap current, is necessary.
For examining the system response with a finite oscillating amplitude, a set of direct simulations is performed as shown in figure 8. Initial oscillation amplitude for the relative magnetic island width is set as a control parameter, while that for the bootstrap current is given as the zeroth order equilibrium value.Three different cases are shown in figures 8(a)-(c), where the experimental case corresponds to figure 8(b).The simulation result in figure 8(b) is in qualitative agreement with the observation in terms of the oscillation frequency and the phase relation of two quantities.Both the experiment (figures 3(f ) and (g)) and the simulation (figure 8(d)) show clockwise rotation in the normalized island width versus plasma current (density) plot.Moreover, considering the Xpoint geometry, the cross-section surface area in which the modulated bootstrap current flows is estimated to be 0.03 m 2 .The expected amplitude in the total plasma current is given as ∼0.3 kA, which is also in agreement approximately with the experimental value (figure 2(b)).It should be noted that the initial amplitude was given only for the magnetic island width term in this examination, and the oscillation amplitude in the bootstrap current was determined by the system itself.
Examining figures 8(a)-(c), it is found that the system response is nonlinear with respect to the initial amplitude.When the initial amplitude is small, the oscillation is characterized by nearly sinusoidal waves.Distortion in the waveform becomes prominent as the initial amplitude increases, in particular, the bootstrap current is more bursty as shown in figure 8(c).Nevertheless, the phase relation between the magnetic island width and the bootstrap current is maintained, as shown figure 8(d).As the initial amplitude increases, the system frequency decreases, which is clearly shown in figure 8(e).

Spatially localized density perturbations and turbulence dynamics
As the end of this paper, additional observations that might also be important for the SSDO dynamics are presented; spatially localized density perturbations and turbulence dynamics.Although overall properties are successfully captured by the predator-prey model shown above, how the oscillation is triggered and is quenched, as well as how the nonlinear waveform is determined, remain unsolved.In particular, there are time points at which the decreasing trend of the radial magnetic field suddenly changes, as shown by the vertical lines in figure 9(a), which cannot be described by the model.
In LHD, a 13-channel radial array of vertical far-infrared laser interferometers is installed.Figures 9(b) and (c) show the time evolution of the line-averaged density distribution and its fluctuation component, respectively.When the SSDO starts at 3.87 s, an overall increase of the line-averaged density is observed.This increased line-averaged density is maintained until the X-point radiation starts to increase, at which a very localized transient peak forms at the code of R = 3.75 m.This peak moves inward as time goes on, after which an overall Figures 9(d) and (e) show the spatial distribution of turbulence fluctuation power and its horizontal propagation velocity, respectively, measured by a phase-contrast imaging system [26].At the outer minor radii of the present interest, the horizontal propagation velocity corresponds to the poloidal velocity.The positive value of V turb is defined as the velocity in the electron diamagnetic drift direction in r eff /a 99 > 0 (r eff is the effective minor radius).Two different kinds of turbulence seem to coexist: turbulence propagating in the ion diamagnetic drift direction at r eff /a 99 ∼ 1 and that propagating in the electron diamagnetic drift direction towards the inner side.At the points in time when the decreasing trend of the radial magnetic field changes, a suppression of turbulence amplitude and an acceleration of the turbulence propagation velocity in the electron diamagnetic drift direction are found.Incorporating those factors into the model are desirable for deeper understanding of the phenomenon and for better predictability.

Summary
In this paper, we determined a new type of SSDO in LHD.The oscillation of 40 Hz was observed in the divertor ion saturation current, the magnetic island width, the magnetic island X-point radiation, and the plasma current, which was regarded as a sequential repetition of the detachment-attachment transitions and back-transitions.A phenomenological modeling approach was taken to describe the basic property of the oscillation, by combining the modified Rutherford equation that provided the magnetic island dynamics to an ad-hoc model for the bootstrap current at the X-point of the island.Validity of the model was examined in terms of the oscillation frequency and damping rate, the phase relation between two parameters, and the oscillation amplitude.
As the model was designed to capture the overall feature of the oscillation to maintain physical transparency, unresolved issues remained.For example, how the oscillation was triggered and was quenched, as well as how the nonlinear waveform was determined, were unsolved.The model this stage only has the quantitative predictability in the precision of the order of magnitude.Moreover, a critical discrepancy between the observation and the model prediction was that there seemed to be a limit-cycle trajectory in the experiment, while a finite damping term, which can be occasionally much smaller than the real frequency though, was inevitable.In that sense, the present model is not 'self-sustained'.A structural improvement of the model is therefore foreseen.
The present observation is relevant for tokamak island divertor development with PF in terms of the role of the bootstrap current for the island extension.An SSDO similar to LHD was occasionally observed in J-TEXT [14], where nonlinearity in island-plasma interaction was likely essential.Examination of the predator-prey model on the tokamak edge is foreseen.

Figure 1 .
Figure 1.Time evolutions of plasma parameters for a discharge with detachment-attachment transitions and back-transitions: (a) line-averaged density and heating scenarios, (b) divertor ion saturation current and X-point radiation, (c) plasma current and diamagnetic stored energy, (d) normalized amplitude and phase of m/n = 1/1 radial magnetic field structure, and radial profiles of (e) electron density and (f ) electron temperature.

Figure 2 .
Figure 2. Time evolutions of plasma parameters during detachment-attachment transitions and back-transitions: (a) divertor ion saturation current and X-point radiation, (b) profile of divertor ion saturation current across the strike-line (∆x being the distance from the peak position), (c) fluctuation component of plasma current, and (d) normalized amplitude of m/n = 1/1 radial magnetic field structure.

Figure 3 .
Figure 3. Lissajous diagrams of (a)-(d) X-point radiation and (e)-(h) plasma current with respect to normalized amplitude of m/n = 1/1 radial magnetic field.Rows correspond to time periods indicated in figure 2 top by colored rectangles.

Figure 4 .
Figure 4. Time evolutions of (a) neutral pressure in private flux region and X-point radiation and (b) electron temperature and ion saturation current on divertor target plate.

Figure 5 .
Figure 5.Time evolutions of (a) plasma minor radius calculated according to electron kinetic stored energy profile and X-point radiation and (b) electron temperature inside edge m/n = 1/1 magnetic island (R = 4.58 m); and radial profiles of (c) electron temperature and vacuum rotational transform without RMP, (d) ion collisionality, and (e) inverse mean free path for attached and detached phases.Time periods used for obtaining radial profiles are indicated at top by colored rectangles.Dashed curve in (e) corresponds to ϵ −1.5 , where ϵ is inverse aspect ratio.

Figure 6 .
Figure 6.Illustration of magnetic island expansion by bootstrap current.(a) Case without seed island and (b) case with seed island that breaks poloidal symmetry of bootstrap current and induces modulation in radial magnetic field.

Figure 7 .
Figure 7. Eigenvalues of proposed predator-prey model that depend on parameters; vacuum island width δ; equilibrium island width W 0 ; and linear growth rate of bootstrap current evolution α.Here, real-part and sign-inverted imaginary-part of eigenvalues correspond to system frequency and damping rate, respectively.δ-W 0 dependence of (a) real-part and (b) sign-inverted imaginary-part of eigenvalues, (c) their slices at δ = 0.07 m, and (d) α dependence.Experimental parameters are shown by open circles.

Figure 8 .
Figure 8. (a)-(c) Time evolutions of relative island width and bootstrap current calculated with different initial perturbation amplitude, (d) their lissajous diagrams, and initial amplitude dependence of oscillation frequency.Experimental case corresponds to (b).

Figure 9 .
Figure 9.Time evolutions of (a) normalized amplitude of m/n = 1/1 radial magnetic field structure and X-point radiation; (b) radial profile of line-averaged electron density and (c) its fluctuation component; and (d) radial profile of turbulence fluctuation power and (e) its horizontal propagation velocity (nearly poloidal velocity at outer minor radii, where V turb > 0 in r eff /a 99 > 0 corresponds to electron diamagnetic drift velocity direction).Vertical lines in (a) correspond to time points where decreasing trend of radial magnetic field changes.