Hydrogen Isotope Effect on Self-Organized Radial Electric Field Criticality During Electron Internal Transport Barrier Formation in Toroidal Plasmas

Self-organized structure formation in magnetically conﬁned plasmas is one of the most attractive subjects in modern experimental physics. Nonequilibrium media are known to often exhibit phenomena that cannot be predicted by superposition of linear theories. One representative example of such phenomena is the hydrogen isotope effect in fusion plasmas, where the larger the mass of the hydrogen isotope fuel is the better the plasma conﬁnement becomes, contrary to what simple scaling models anticipate. In this article, threshold condition of a plasma structure formation is shown to have a strong hydrogen isotope effect. To investigate the underlying mechanism of this isotope effect, the electrostatic potential is directly measured by a heavy ion beam probe. It is elucidated that the positive radial electric ﬁeld structure can be driven by less input power normalized by plasma density in plasmas with larger isotope mass across the structure formation.


Introduction
A long-standing open question in the magnetically confined plasma study is the hydrogen isotope effect: [1][2][3][4][5] the heavier the fuel particles are the better the confinement becomes, contrary to what simple scaling models anticipate. The hydrogen isotope effect has been classically assessed by means of the global confinement time scaling or the local diffusivity. When evaluating such parameters, it is usually inevitable to rely on specific assumptions in heating profile calculation or transport models, which can be sources of uncertainty. Moreover, the hydrogen isotope effect in global or local transport coefficients is not always distinct but can be marginal particularly in stellarators. 5 Therefore, a large scale statistics is necessary for drawing a conclusive trend, which can hamper interpretation on detailed underlying process.
The hydrogen isotope effect is also known to be evident in the edge transport barriers (ETBs) 4,[6][7][8] or the internal transport barriers (ITBs), 9 where the threshold condition for the transition is eased in plasmas with larger ion mass. The transport barrier criticality can be qualitatively defined without any model assumption, therefore an experiment-driven discussion for unveiling the underlying physics of the hydrogen isotope effect is possible. Furthermore, as the future thermonuclear fusion plant is currently planned to be operated with the tritium-deuterium mixture fuel and with transport barriers, understanding of the background mechanism of the hydrogen isotope effect is important for projecting the reactor plasma performance. When the transport barriers emerge, a localized radial electric field structure is spontaneously formed by the quasi-neutral plasma, 10 which is believed to suppress the anomalous turbulent transport. Direct measurement of the radial electric field in systematic isotope experiments is therefore a key factor for resolving this issue, but has been challenging for a long time.
In this article, the hydrogen isotope effect on the ITB criticality in the Large Helical Device (LHD) stellarator is elaborated, focusing upon the ITB in the electron temperature profile, the so-called electron ITB. The threshold value in the input power normalized by plasma density for triggering the electron ITB is found to be lowered in plasmas with larger isotope mass. Heavy ion beam probe (HIBP) reveals that the radial electric field transition at the ITB formation is induced by a lower input power normalized by density in plasmas with heavier ions. The fact that the self-organized radial electric field criticality has a susceptibility to plasma ion mass can be a missing piece in resolving the hydrogen isotope effect.

Experimental set-up
Hydrogen isotope effect in the electron ITB is dynamically investigated by means of an input power modulation experiment in LHD (see Methods for details). Electron cyclotron resonance heating (ECH) focused at the torus core is modulated with the peak-to-peak amplitude of P ECH = 1 MW and the frequency of f mod = 23 Hz. In order to assess the threshold input power normalized by plasma density for the ITB formation, a shot-to-shot basis density scan is performed in deuterium (D) dominant plasmas, hydrogen (H) dominant plasmas, and D-H mixed plasmas. Here, the D content determined by a passive spectroscopy is maintained approximately constant at 89 %, 4 %, and 60 %, respectively, during the density scan.

Hydrogen isotope effect in ITB threshold
Evolution of the electron temperature (T e ) profile is measured by the Thomson scattering system. 11  , is plotted as a function of the applied ECH power normalized by the line averaged density P ECH /n e in Fig. 1 (c). When P ECH /n e is low, |δ ∇T e | remains low as well. As P ECH /n e is increased, |δ ∇T e | begins to monotonically increase once P ECH /n e exceeds threshold values that depend on the D content.
In Fig. 2, the threshold values are plotted as a function of the D content by taking the value of P ECH /n e at which |δ ∇T e | surpasses 2.5 keV/m. Error bars in the threshold diagram is mainly consequent upon an imperfect reproducibility among  discharges, i.e., vertical scatters of points in Fig. 1 (c). The threshold value in P ECH /n e is found to monotonically decay as the D content is increased. Eased threshold condition for the different types of transport barriers with larger ion mass is ubiquitously observed in various torus plasmas including tokamaks and stellarators, 4, 6-9 implying existence of an analogous underlying physics possibly related to the radial electric field drive. Note that the energy flow from electrons to ions is less than 40 kW, which is negligibly small and plays a minor role for the threshold condition.
Local transport property is examined through the electron thermal diffusivity obtained in a perturbative manner, the socalled incremental diffusivity χ inc e ≡ |δ q e |/n e |δ ∇T e |, where |δ q e | is the heat flux variation having the transport time scale and n e is the electron density. The nonlocal part of the electron heat flux was almost insensitive to the D content. 12 Figure 1 (d) displays P ECH /n e dependence of χ inc e at r eff = 0.15 m. In the lower P ECH /n e regime, χ inc e increases with P ECH /n e . The upper and lower envelopes of the data points are found to be proportional to (P ECH /n e ) 1.9 , as shown by grey dashed curves. Gradual confinement degradation as the density is lowered is a common feature of stellarators, and, for instance, is seen in the energy confinement time scaling. 13 In the larger P ECH /n e regime, χ inc e stays off from the (P ECH /n e ) 1.9 trend, and is maintained low. Those two regimes correspond to non-ITB and ITB regimes, respectively, and the transition threshold value of P ECH /n e clearly depends on the D content, as discussed above.

Electrostatic potential measurement by HIBP
The electron ITB in stellarators is frequently accompanied by a localized positive radial electric field (E r ) structure at the core, [14][15][16][17] which is primarily induced through the neoclassical process. Therefore, the electron ITB in stellarators is sometimes referred to as the core electron root confinement (CERC). In order to investigate the role of E r on the isotope effect in the electron ITB threshold, direct measurement of the electrostatic potential (φ ) at the plasma core is performed by means of HIBP. 18 HIBP in LHD provides a single point measurement. To obtain the radial profile, the measurement position is scanned either in the 10 Hz reciprocation manner or in a shot-to-shot manner, as shown by the thick curve or the open circles in Fig. 3 (a), respectively. Here, a 99 ∼ 0.6 m is the plasma edge in which 99 % of the electron kinetic energy is confined. The local electrostatic potential is translated into the probe beam image position on detector plates [ Fig. 3 (a) insert], which are located in the energy analyzer. Current flowing through the lower and upper detector plates, I L and I U , respectively, are recorded as the raw HIBP data. The electrostatic potential is given as φ = CI diff I −1 HIBP , where I diff = I L − I U and I HIBP = I L + I U are the HIBP difference current and the HIBP total current, respectively, and C is the proportionality factor.
As an example, radial profile of φ measured by HIBP in low density (n e ∼ 0.8 × 10 19 m −3 ) stationary hydrogen discharges is shown. Figure 3 (b) shows φ profiles in a balanced NB discharge and a stationary ECH superposed discharge measured by a reciprocation position scan. Figure 3 (c) exhibits the electron and ion temperature profiles, where the electron ITB is formed in the ECH superposed discharge. It is evident that φ profile with a strong negative gradient, i.e., positive E r , emerges in the core when the ITB is formed. This positive E r is likely ascribed to the electron root solution in the neoclassical transport balance as described in. [14][15][16][17] Meanwhile, the non-ITB plasma has an almost flat φ profile in the core, which corresponds to the ion root solution.

Hydrogen isotope effect in radial electric field criticality
HIBP measurement is applied to the ECH modulation discharge to study the critical behavior in E r across the ITB formation. Aiming at capturing a marginal response in φ particularly in higher density discharges discussed below, the measurement position is fixed during the discharge. The radial profile of the φ modulation amplitude is obtained in a shot-to-shot HIBP position scan. Figure 4 (a) illustrates the time evolution of a low density hydrogen target discharge. The ECH modulation is applied on the nearly stationary plasma. Figure 4  In order to extract the repetitive nature of φ variation, the conditional averaging is performed asφ mod is the period of the ECH modulation. The value t i indicates the i-th turn-on time of ECH and N is the total number of the modulation. Figures 4 (d) and (e) show the conditionally averaged φ signals at r eff = 0.08 m and 0.26 m, i.e., inside and outside the ITB, respectively. In addition, reconstructed T e variation measured by the Thomson scattering system 19 are overlaid. Inside the ITB, φ and T e simultaneously increase when the ECH is applied, while outside their modulation amplitudes are much smaller. This oscillation amplitude difference in φ is indicative of the existence of E r oscillation that plays a role for the electron ITB formation. It is worthwhile to note that there is an asymmetry in time evolutions of φ and T e , i.e., a stepwise increase and a nearly linear decay. The ITB begins to form at the moment of the ECH application, and quickly saturates making the step-like waveform. When the ECH is terminated, the ITB slowly weakens and eventually returns to the non-ITB plasma. Figures 5 (a) and (b) compare the conditionally averaged φ at different radial positions in D and H low density plasmas of n e ∼ 0.8 × 10 19 m −3 . The oscillation amplitude gradually decays in the radially outward direction in H plasmas, while it stays nearly unchanged in D plasmas. In order to more distinctly quantify this difference in D and H plasmas, the radial profile of the peak-to-peak φ oscillation amplitude for the fundamental modulation frequency component (23 Hz) is shown in Fig. 5 (c). Although the oscillation amplitude is approximately equivalent at the core, its radial extent in D plasmas is wider than that in H plasmas. The E r peak locations in H and D plasmas are approximately at r eff ∼ 0.2 m and r eff ∼ 0.3 m, respectively. A stronger but more localized E r structure at the ITB region is seemingly excited in H plasmas in contrast to a moderate and wider structure peaked at slightly outer radius in D plasmas. Figure 5 (d) shows T e variation by the ITB formation. The location where the electron temperature gradient modulation is largest is r eff ∼ 0.1 m, which is clearly inner than the E r peak location. D plasmas have a stronger ITB with a slightly wider radial extent. Since E r is zero in the core by definition, the E × B flow shear at the ITB location is estimated as E larger E peak r and smaller r peak eff . Nevertheless, the resultant ITB is stronger in D plasmas, indicating that the E × B flow shear paradigm is insufficient for describing the observation. Phenomenologically, a wide radial extent of the positive E r structure and its shear in D plasmas seem to be more effective for transport suppression and significant ITB formation, rather than the large but localized structure in H plasmas. Importance of positive E r in diminishing the zonal flow damping is argued in, 20 which can be an interpretation of this observation. According to this model, the zonal flow is a hidden key factor, although its detection in the core is still challenging. Enhancement of the zonal flow activity in D plasmas was experimentally addressed in. 21 Regardless of the radial gradient scale, the φ modulation amplitude at the plasma core is a good proxy for the integrated magnitude of the positive radial electric field modulation. Therefore, response in φ to the ECH modulation at the plasma core is explored in the density scan experiment to discuss the underlying mechanism of the isotope effect in the ITB threshold. A diagnostic difficulty for this experiment is that the HIBP signal intensity (total current I HIBP ) exponentially decays as the density is increased. In such a low signal intensity case, the value of φ ∝ I diff I −1 HIBP can be uncertain because of small value in the denominator. In order to settle this problem, the HIBP difference current I diff is used as a proxy of local φ . Moreover, level of φ response to the ECH modulation pulse is statistically characterized as γ 2 sin α, where γ 2 and α are the squared cross coherence and the cross phase between I diff and P ECH , respectively. If φ modulation is induced by the ECH modulation, α must be π/2 as shown in Fig. 4 and Fig. 5. Figures 6 (a) and (d) show γ 2 sin α at the modulation ECH frequency of 23 Hz plotted against P ECH /n e for D and H plasmas. An example of the frequency spectrum is given in the insert in Fig. 6 (d). In both cases γ 2 sin α decays as P ECH /n e is decreased (n e is increased).
The decay in γ 2 can be caused either due to disappearance of φ modulation in the core or due to deteriorated signal-tonoise ratio. The latter possibility is examined through a synthetic data approach. Decay characteristic of γ 2 is estimated using a mimic HIBP difference current I HIBP ≡Î HIBP sin 2π f mod t + WGN(σ 0 ), where WGN(σ 0 ) is the white Gaussian noise with a standard deviation of σ 0 . The modulation amplitudeÎ HIBP is decayed following the known density dependence,Î HIBP = I 0 exp[−n e /0.55 × 10 19 ], obtained by a dedicated test. The ratio between I 0 and σ 0 is determined in such a way that the γ 2 curve passes through the experimental data points in the low density domain, where φ oscillation is unambiguously observed. The synthetically derived squared cross coherence, referred to as γ 2 max , represents the case where the coherence decays only due to the descent in the HIBP signal intensity.
As shown in Fig. 6 (a), the measured value follows γ 2 max in the entire P ECH /n e range in D plasmas. This indicates that φ modulation occurs in all the cases in D plasmas and the variation in the coherence is solely due to the declined HIBP signal intensity in higher density range. In contrast, the measured value stays off from the γ 2 max curve at P ECH /n e ∼ 0.8 MW/10 19 m −3 in H plasmas as displayed in Fig. 6 (d). It turns out that φ at the core scarcely responds to the ECH pulse in P ECH /n e < 0.7 MW/10 19 m −3 . This difference in D and H plasmas is more evident in the relative plot to the synthetic data, γ 2 sin α/γ 2 max , as shown in Figs. 6 (b) and (e). Figures 6 (c) and (f) indicate the peak-to-peak modulation amplitude in T e gradient. H plasmas begin to form the ITB structure at P ECH /n e ∼ 0.8 MW/10 19 m −3 , where φ oscillation becomes substantial. In contrast, the ITB in D plasmas is maintained in the entire range of P ECH /n e with the meaningful φ oscillation. This observation indicates that the critical behavior in the E r structure formation that depends on the isotope mass plays a dominant role for the isotope  Figure 6. (a,d) Squared coherence of the HIBP difference current with respect to the ECH modulation pulse, (b,e) that normalized by its maximum possible value, and (c,f) the peak-to-peak modulation amplitude of the electron temperature gradient as a function of the applied ECH power normalized by the line averaged density for deuterium plasmas and hydrogen plasmas, respectively; and (insert) squared coherence spectrum for a low density hydrogen discharge. dependence of the ITB threshold. Note that the threshold value in P ECH /n e is not perfectly identical to the case shown in Fig. 1 due to different experimental conditions.

Discussion
Since the neoclassical E r transition plays a role for the ITB criticality, [14][15][16][17] whether neoclassical E r has a susceptibility to the plasma ion mass is numerically examined by a local neoclassical transport code DKES/PENTA. 22 Figure 7 shows the P ECH /n e dependence of calculated neoclassical E r at r eff ∼ 0.2 m in the ECH applied phase. Although neoclassical E r is systematically larger in D plasmas in the electron-root solution, critical behavior is visible neither in D nor H plasmas. This examination implies that other mechanisms rather than the simple local neoclassical theory are necessary for explaining the isotope effect in the ITB threshold. Zonal flows 20 having an isotopic mass dependence 23 would be a candidate as the underlying mechanism. Global neoclassical effects are also important particularly close to the magnetic axis, which is dependent on the ion gyroradius.
An extra input power necessary for the E r structure formation at the high-confinement mode transition in hydrogen plasmas is reported in the ASDEX-Upgrade tokamak, 8 which possibly has a link to the present case. In addition, an important role of the fast ions produced by NBs on core transport suppression is recently shown to have an isotope effect, 24 which can also contribute to the isotope effect in the ITB criticality.

6/8 Summary
In conclusion, the isotope effect in the electron ITB threshold favorable to plasmas with heavier ions was found. The radial electric field structure likely responsible for the ITB formation was easily formed in plasmas with heavier ions. A better core plasma performance with ITB in future deuterium-tritium plasmas is therefore foreseen. Transport barrier structures are regarded as a dissipative structure formed spontaneously in nonequilibrium open system by an intensive power input. 25 The present work also contributes to understanding of the dissipative structure affected by basic properties of the media.

Large Helical Device (LHD)
LHD is a plasma current-free magnetically confined torus plasma device. High-temperature plasmas are confined in the externally excited helical magnetic field structure with the toroidal magnetic field of B t = 2.75 T and the vacuum torus axis of R = 3.6 m. In the vacuum magnetic field configuration, the rotational transform, ι/2π = 1/q, where q is the safety factor, monotonically increases with the radius. The rational surfaces of ι/2π = 0.5 and 1 typically exist at the core and the edge, respectively. In the present experiments, the plasma is sustained by balanced tangential neutral beams (NBs) with P NB = 3.6 MW. The deuterium content is determined as I D α /(I D α + I H α ), where I D α and I H α are the intensities of the D α emission and the H α emission measured by a passive spectroscopy, respectively.

Heavy Ion Beam Probe (HIBP)
HIBP provides a direct measurement of the plasma electrostatic potential. 18 Heavy ion beam injected to the plasma is secondary ionized, where the beam trajectory diverges from the original beam trajectory. A part of the secondary ionized beam is selected by a slit, and the beam energy is analyzed. The energy difference between the original beam and the secondary ionized beam is equivalent to the electrostatic potential at which the secondary ionization occurs. To manipulate the measurement location, operational parameters for the HIBP, such as the injection beam energy, the incident angle, and the confinement magnetic field, are chosen according to the result of the beam trajectory calculation.