This study deals with the equilibrium beta limit for high-beta plasmas of the Large Helical Device using a 3D equilibrium calculation code, HINT, which is an initial value solver based on the relaxation method without an assumption of nested flux surfaces. For the finite beta equilibrium, due to the nonlinear 3D equilibrium response, the flux surface breaking begins, and the magnetic field line becomes stochastic in the peripheral region. However, although the magnetic axis shifts until the conventional theoretical beta limit, the separatrix limiting the equilibrium does not appear in the plasma core. For the high beta equilibrium, breaking of the force balance begins in the strongly stochastized region. To keep the force balance in that region, the pressure gradient decreases and the fixed profile is relaxed. As a result, the volume averaged beta, which is equivalent to the plasma stored energy, saturates, although the peak beta increases. The beta value, where the force balance begins breaking in the inside of the stochastic region, is proposed as a new index of the equilibrium beta limit for the typical heliotron-type stellarator. According to this new index of the equilibrium beta limit, the achieved beta value in the experiment of the Large Helical Device is still lower than the equilibrium beta limit. This suggests that if the heating power can be increased, the achieved beta value may still increase.

## I. INTRODUCTION

Generating and keeping a good magnetohydrodynamic (MHD) equilibrium with the reactor-relevant beta value are aims of the magnetic confinement research. To date, the high-beta equilibrium in the stellarator is studied experimentally^{1–10} and theoretically.^{11–15} In particular, in the Wendelstein 7-AS (W-7AS), the volume averaged beta, $\u27e8\beta \u27e9$, achieved about 3.4% in the medium magnetic field, *B _{t}* = 0.9–1.1 T.

^{1}In those high-beta discharges, the topological change of the magnetic field is expected by the numerical modeling of the 3D equilibrium calculation. Reiman

*et al.*pointed out that the topological change of the magnetic field in the peripheral plasma region may limit the increasing

*β.*

^{12}In the Large Helical Device (LHD), the $\u27e8\beta \u27e9$ achieved up to 4.8%, which is relevant to a parameter of the reactor operation, for the quasi-steady state operation in the low magnetic field,

*B*= 0.41 T.

_{t}^{8}Also, in the medium magnetic field,

*B*= 1 T, the achieved $\u27e8\beta \u27e9$ extended up to 4.1%.

_{t}^{10}Figure 1 shows the Poincaré plot of magnetic field lines for a reconstructed finite-

*β*equilibrium with the kinetic electron profile. Profiles of the connection length of the magnetic field line, $LC$, from the equilibrium calculation, the electron temperature,

*T*, and the electron mean free path (e-mfp),

_{e}*λ*, from the experiment (shot #69910) are also shown as the reference. The magnetic axis shifts from 3.53 m to about 3.8 m, and the magnetic field line becomes stochastic, in particular, in the peripheral plasma region. Although the magnetic field structure becomes stochastic, it is seemed that

_{e}*T*spreads over in the peripheral plasma region and the gradient of

_{e}*T*exists in that region. Suzuki

_{e}*et al.*pointed out a hypothesis that the stochastic magnetic field can keep

*T*

_{e}^{11}in a following situation. If the $LC$ is sufficiently longer than the

*λ*, the stochastic field line can keep the $\u2207Te$ from the viewpoint of the parallel electron transport. The estimated

_{e}*λ*from the profile in Fig. 1 is less than 10 m in the peripheral plasma region. It is much shorter than the $LC$. Thus, the reconstructed equilibrium may be consistent with the experimental result. An important remark in the high-beta experiment in the LHD is no significant termination or degradation of the confinement by the MHD instability.

_{e}^{8}This means that the beta value, $\u27e8\beta \u27e9$, did not achieve to the beta limit, and the $\u27e8\beta \u27e9$ may increase if the heating power increases. Probably, the beta may be limited by the MHD instability than the equilibrium beta limit. However, it is an important issue regarding how the equilibrium limits the increased

*β*.

Theoretically, the equilibrium beta limit is defined by the Shafranov shift, $\Delta /a$, where $\Delta =(Rax(\beta )\u2212Rax(0))$ and the *a* is the effective minor radius. If the Shafranov shift achieves about 0.5, the equilibrium is limited. For tokamaks, at the equilibrium beta limit, the separatrix appears in the plasma core, because the poloidal field by the plasma current is compensated by the external vertical field to keep the force balance.^{16} The separatrix in the core leads to the degradation of the confinement. In the stellarator, the equilibrium beta limit can be defined by the Shafranov shift, as well as tokamaks. However, Hayashi *et al.* pointed out that the equilibrium beta limit might be defined by the stochasticity of the magnetic field and it is more severe than the limit defined by the Shafranov shift.^{17} Since the pressure-induced perturbed field breaks the nested flux surface, the stochastization of the magnetic field by the nonlinear 3D equilibrium response is an intrinsic property in the stellarator. The problem regarding how to define and represent the equilibrium beta limit is still an open question. Thus, understanding the equilibrium beta limit is a critical and urgent issue from a viewpoint to aim the steady-state reactor on the stellarator line.

Recently, the equilibrium beta limit for the stellarator was studied theoretically and numerically.^{14,15} In those studies, the significant stochastization of the magnetic field in the peripheral plasma region degrades the plasma confinement, and then the plasma beta achieved up to the limit. In particular, Suzuki pointed out that the beta limit in a heliotron-type stellarator is not disruptive, and it appears as the degradation of the increased *β.*^{15} In this study, extending the previous study^{15} for the high-beta LHD, we theoretically deal with the equilibrium beta limit in the LHD plasmas as a typical heliotron-type stellarator. Special notice is how to define and represent the equilibrium beta limit. Analyses were done by a 3D MHD equilibrium code, HINT,^{11,18} which does not assume nested flux surfaces *a priori*. Using HINT, we can study not only the Shafranov shift but also the stochasticity of the magnetic field on the equilibrium beta limit. In Sec. III, we study the beta effects in the 3D equilibrium. Then, we discuss how to define the equilibrium beta limit and summarize.

## II. BRIEF DESCRIPTION OF THE HINT CODE AND NUMERICAL SETUP

The HINT code is a 3D MHD equilibrium calculation code based on the relaxation method solving single-fluid nonlinear MHD equations on the right-handed system of the cylindrical coordinate $(R,\varphi ,Z)$.^{18} This means that the code allows the flux surface breaking by the magnetic island or the stochastization of the magnetic field due to the 3D equilibrium response. The relaxation process consists of two steps: (i) step-A and (ii) step-B.

The first step (step-A) is a relaxation process of the plasma pressure satisfying a condition, $B\xb7\u2207p=0$, with fixed **B**. The condition, $B\xb7\u2207p=0$, is equivalent to the no parallel pressure gradient of *p*. In the HINT code, an “averaged” plasma pressure on a flux tube is calculated instead of directly calculating $B\xb7\u2207p=0$. The relaxed (averaged) pressure on a grid is given by an averaging process:

where *i* means a step number of iterations, *L _{C}* is the connection length of a magnetic field line starting each grid point (

*L*is finite for open magnetic field lines), and $Lin$ is the length along a magnetic field line followed from each grid point and prescribed as an input parameter. That is, the distribution of the relaxed plasma pressure depends on the relative magnitude between the connection length

_{C}*L*and the parameter $Lin$ prescribing the length of the field line trace.

_{C}The second step (B-process) is a relaxation of the magnetic field with fixed *p*, which is calculated in step-A. In the HINT code, the magnetic field is split by $B=B0+B1$. Here, a subscript 0 means the vacuum field fixed in time and another subscript 1 means the perturbed equilibrium component. To speed up the convergence, dissipative parameters, *η* and *ν*, are assumed to be constant in the calculation.convective, and *ρ* = 1 is assumed. Thus, artificial dissipative MHD equations can be derived,

where $v0$ is a given toroidal ow velocity and $v1$ is the MHD ow velocity. Thus, the total plasma velocity is denoted by $v=v0+v1$. The $jnet$ is the net toroidal current like Ohmic or the bootstrap current. In this study, the toroidal rotation does not consider and $v0=0$ assumes. and the zero net toroidal current, $jnet$, assumes in this study. The last term in the inductive equation is a divergence cleaning term and $\kappa divB$ is a small constant. The final equilibrium state is solved by these two steps iteratively until $dv/dt\u21920$ and $dB/dt\u21920$. The 3D MHD equilibrium is obtained as an ohmically steady state.

The pressure profile is approximated by the polynomial function of the normalized toroidal flux, *s*, which is defined by $s=\Phi /\Phi edge$. The normalized minor radius is defined by $\rho =s$. The toroidal flux, $\Phi $, is calculated by the integration along the contour line of the plasma pressure, *p*, not along the flux surface, and $\Phi edge$ is the toroidal flux within the “last closed contour line” of the plasma pressure. In this study, the last closed contour line of the plasma pressure always passes through a fixed point at *R *=* *4.6 m, a torus outboard point on an equator. It is an input parameter.

In this study, a following fifth order polynomial function,

is used to approximate the pressure profile. For the high-beta plasma of the LHD, since the electron temperature, *T _{e}*, is proportional to $\u221d(1\u2212s)$, and the electron density,

*n*, is proportional to $\u221d(1\u2212s4)$, this is a good approximation to prescribe the pressure profile of the high-beta plasma. The

_{e}*p*

_{0}is the plasma pressure on the magnetic axis, and it controls the beta value. Two beta values, the peak beta

*β*

_{0}and the volume averaged beta $\u27e8\beta \u27e9$, define as

Here, the bracket $\u27e8\u27e9$ means the volume averaged quantity defined by

It is noted that in the HINT modeling the plasma volume is defined by an integral of a volume element, $dV=RdRd\varphi dZ$, for $p>p@s=1$ on the cylindrical rectangular grid. Also, the effective minor radius *a* is calculated by $a=V/2\pi 2Rax$.

## III. 3D EQUILIBRIUM PROPERTIES FOR HIGH-*β* PLASMA IN OPTIMIZED CONFIGURATION AND EQUILIBRIUM BETA LIMIT

*β*

In this section, we study the 3D equilibrium properties in an optimized configuration of the LHD. In the LHD experiments, the achieved beta value strongly depends on the preset magnetic axis for the vacuum field $RaxV$ and the plasma aspect ratio $Ap$.^{4} In the standard configuration, $RaxV$ and $Ap$ are set to 3.6 m and 5.8 m, respectively. The achieved beta value in the standard configuration was about 3% and the Shafranov shift, $\Delta /a$, approached about 0.5. Thus, in order to reduce the Shafranov shift, the optimization of $Ap$ was examined because the Shafranov shift is proportional to the inverse of $Ap$.^{6} For $Ap$ = 6.6, a highest beta value was obtained. In this study, we study this optimized configuration optimized for $Ap$ = 6.6.

Figure 2 shows Poincaré plots of magnetic field lines for the 3D equilibrium of the optimized configuration at the horizontally elongated cross section. As a reference, the vacuum field is also plotted in the same cross section. The blue lines in the figures indicate the well-defined last closed flux surface (LCFS) for the vacuum along the *Z *=* *0 plane, and the green lines indicate the LCFS for the finite-*β* equilibrium. For the vacuum field, clear flux surfaces are keeping till the peripheral region. Increasing *β*, the magnetic axis and the inward LCFS shift to the outward of the torus. Also, the magnetic field line becomes stochastic in the peripheral region. However, for $\u27e8\beta \u27e9\u223c$ 5.2%, the outward LCFS is still kept in the vacuum LCFS. For higher $\u27e8\beta \u27e9$ (∼6.7%), the LCFS shrinks in the inside of the vacuum LCFS. Increasing the stochasticity, the plasma volume with clear flux surfaces deceases about 45% for $\u27e8\beta \u27e9\u223c$ 6.7%.

In order to understand the stochasticity of the magnetic field due to increasing *β*, positions of the magnetic axis ($Rax$), inward and outward of the LCFS ($Rin$ and $Rout$) as the function of $\u27e8\beta \u27e9$, are plotted in Fig. 3. The magnetic axis $Rax$ monotonically changes due to increasing *β*. On the other hand, although $Rax$ and $Rin$ shift to the outward due to increasing *β*, $Rout$ is fixed near the vacuum LCFS until $\u27e8\beta \u27e9\u223c$ 5.2%. However, for higher *β*, $Rout$ shrinks inside the vacuum LCFS and the magnetic field becomes strongly stochastic. In more than $\u27e8\beta \u27e9$ = 6.3%, the breaking of the force balance begins on the inside of the vacuum LCFS, and the equilibrium calculation did not converge. To keep the force balance and converge the calculation, the pressure gradient in the stochastic region decreases and the fixed profile is relaxed. For $\u27e8\beta \u27e9<$ 6.3%, the flux surface of *s *=* *1 passes through the point of *R *=* *4.6 m. However, for $\u27e8\beta \u27e9>$ 6.3%, the stochastic field cannot keep the pressure gradient until *R *=* *4.6 m. Thus, the fixed point of *s *=* *1 changes to *R *=* *4.4 m. As the result, the total plasma stored energy *W _{p}* decreases. We summarize that the 3D equilibrium property in the heliotron-type stellarator involves the following three phases in Fig. 3. In phase (I), the plasma shifts due to the finite beta but $Rout$ is fixed near the vacuum LCFS. In another phase (II), increasing

*β*, the stochastic region increases and the LCFS shrinks compared with the vacuum field. In the last phase (III), increasing the stochasticity, the residual force balance breaks. The magnetic configuration cannot keep the steep pressure gradient.

As we discussed, if the magnetic field becomes strongly stochastic, the stochastic field line cannot keep the steep pressure gradient and the fixed pressure profile is relaxed in order to converge the equilibrium calculation. Here, we theoretically discuss why and how the strong stochastization of the magnetic field line leads the equilibrium beta limit. A schematic view of the relaxed pressure profile in the stochastic region is shown in Fig. 4(a), and another schematic view of a relation between the peak beta, *β*_{0}, and the volume averaged beta, $\u27e8\beta \u27e9$, is shown in Fig. 4(b). In Fig. 4(a), a red line indicates a fixed pressure profile as an input parameter. If the magnetic field does not become strongly stochastic, it can keep the pressure gradient along the red line. However, if the stochastization of the magnetic field becomes strong, the pressure gradient decreases to satisfy the force balance, $\u2207p=j\xd7B$. Since the pressure gradient decreases, the total plasma stored energy also decreases. This represents the reduction of the volume averaged beta, $\u27e8\beta \u27e9$, against the peak beta, *β*_{0}. In Fig. 4(b), the schematic view to represent the reduction of the volume averaged beta is shown. Usually, the volume averaged beta is proportional to the peak beta because the pressure profile is prescribed as the function of the normalized flux, the plasma volume, and so on. However, for the higher *β* equilibrium, the magnetic field cannot keep the prescribed pressure profile in order to converge the equilibrium and the linearity is broken. In Fig. 5, the profiles of the plasma pressure obtained from the 3D equilibrium calculation are shown as the function of the normalized minor radius, *ρ*. In Fig. 5, for $\u27e8\beta \u27e9\u2264$ 6.3%, pressure profiles are almost the same. However, for $\u27e8\beta \u27e9$ = 6.7%, the profile is the same as that for $\rho \u223c$ 0.6 but in more than *ρ* = 0.6, the profile is different. This leads to the reduction of the volume averaged beta.

In order to confirm our consideration, the volume averaged beta obtained from the calculation is plotted as the function of the peak beta. Figure 6 shows the relation between *β*_{0} and $\u27e8\beta \u27e9$. Two auxiliary lines are plotted as the reference. For $\u27e8\beta \u27e9\u2264$ 6.3%, the volume averaged beta increases almost linearly. Since the pressure profile is prescribed as the function of the toroidal flux, which means the function of the cross section, it is not completely linear. However, for higher *β* ($\u27e8\beta \u27e9\u2265$ 6.3%), the stochasticity becomes strong and the fixed pressure profile is relaxed to keep the force balance. Thus, an inflection point appeared for $\u27e8\beta \u27e9\u223c$ 6.3%. We define this inflection point as an index of the equilibrium beta limit (between phases II and III in Fig. 3). Thus, the equilibrium beta limit is observed as the degradation of the plasma stored energy. Special notice is this index does not mean the saturation of increasing the volume averaged beta. In the figure, only the ratio between *β*_{0} and $\u27e8\beta \u27e9$ is reduced and the beta does not saturate. In the considered condition, the equilibrium beta limit is more than 6%. This means the experimental result (#69910) for $\u27e8\beta \u27e9\u223c$ 4.8% does not achieve the equilibrium beta limit. Experimentally, the new index of the equilibrium beta limit will be observed by the fixed pressure profile in the peripheral region because the strongly stochastic magnetic field limits the pressure gradient. In the W7-AS experiment, the fixed electron temperature profile was observed in the peripheral region, if the heating power was increased, which is equivalent to increased *β.*^{1} In the LHD experiment, the similar result may be observed in the beta value near the equilibrium limit, but the change of the magnetic configuration, where the equilibrium beta limit is lower than the optimized configuration, is required.

Increasing the stochasticity of the magnetic field, the increase of the volume averaged beta degrades. However, although the gradient of the function of *β*_{0} becomes weak, the volume averaged beta keeps increasing. The relaxation of the fixed pressure profile means the increase of the peaking factor of the profile. Since the pressure profile becomes effectively peaked, the magnetic axis shifts further. In order to study the equilibrium beyond the index of the equilibrium beta limit, we study the 3D equilibrium for higher *β* ($\u27e8\beta \u27e9\u2265$ 6.5%). Figures 7(a) and 7(b) show the Poincaré plot of the magnetic field line for $\beta 0\u223c$ 16% and a detailed plot near the axis, respectively. In Fig. 7(a), magnetic field lines become strongly stochastic. In addition, the plasma shape strongly elongates because of the large Shafranov shift. Since the vertical field driven by the Pfirsh–Schlüter current cancels the poloidal field generated by the external coil system, the separatrix appears near the axis in Fig. 7(b). In order to satisfy the force balance on the separatrix, the plasma pressure flatters on the separatrix. Thus, the confinement degrades because of flattering of the plasma pressure. The green circle in Fig. 6 indicates the generation of the separatrix. Although the peak beta increases from 12% to 16%, the volume averaged beta is almost constant. The generation of the separatrix gives the most severe limit.

In this study, we deal with the theoretical consideration of the equilibrium beta limit from the 3D equilibrium calculation. However, the beta limit may be decided by the MHD instability. In particular, since the typical-heliotron configuration like the LHD is the magnetic hill in the peripheral region, usually the interchange mode is unstable in there. However, in the high-beta experiment of the LHD, the magnetic fluctuations driven by the interchange mode were observed in the peripheral region and it was identified as the resistive interchange mode.^{8} The growth rate of the resistive interchange mode is proportional to $S\u22121/3$, where the *S* is the magnetic Reynolds number. Thus, in the more reactor-relevant condition with the high *S*, the growth rate will be smaller than the rate in the present experimental condition.

## IV. SUMMARY AND OUTLOOK

We deal with the equilibrium beta limit for the optimized LHD configuration as the typical heliotron-type stellarator. Increasing *β*, the magnetic field becomes stochastic by the nonlinear 3D equilibrium response. If the stochastic field line penetrates deeply into the plasma core region, the increase of the volume averaged beta degrades because the stochastic field line cannot keep the steep pressure gradient. We propose this result as the new index of the equilibrium beta limit for high-beta LHD plasmas.

However, some issues are remaining for further understanding of the equilibrium beta limit in the stellarator, in particular, (i) more sophisticated physics of the transport process on the stochastic magnetic field and also (ii) linkages with the experiments. Regarding the first point, Rechester and Rosenbluth discussed the electron heat transport on the stochastic magnetic field.^{19} Suzuki *et al.* studied an anomalous transport by the stochastic magnetic field for a high-*β* plasma of the LHD.^{20} The radial heat conductivity by the strongly stochastic magnetic field in the peripheral region achieved up to nearly 10^{2} [$m2/s$]. Thus, there is a possibility that the strongly stochastic magnetic field limits the pressure gradient for $\u27e8\beta \u27e9\u2264$ 6.3% because of the increase in the transport. Also, this transport process leads to another problem of the plasma pressure in the stochastic magnetic field. In the stochastic magnetic field, since the magnetic field line opens in many cases, the condition, $B\xb7\u2207p=0$, may violate. In such a case, the parallel force balance is not zero; in other words, an additional term must appear in the equation of the motion to balance the parallel pressure gradient, $\u2207\u2225p$. In the HINT model,^{18} see Eq. (2), the viscous term may balance the parallel force balance. However, as pointed out in Refs. 13 and 21, the Pfirsh–Schlüter current on the stochastic field must be considered more physically. One idea is the anisotropic heat transport to model the parallel and perpendicular pressure gradients.^{22,23} Furthermore, as Watanabe *et al.* pointed out,^{7} the pressure anisotropy along the opened field line may be considered. Studies of the equilibrium beta limit including those transport effects are still ongoing issues. We will report these in future papers.

Regarding the second point, the identification of the stochastic plasma boundary was studied in the LHD experiment.^{24–26} Since the radial electric field reflects the electron heat transport on the stochastic field, the strong radial electric field shear well correlates the connection length of the magnetic field line predicted by the HINT code. In this study, we used the fixed point to prescribe the plasma boundary as the theoretical consideration, but we may use the experimentally obtained database of the plasma boundary. That will give more reliable prediction for the equilibrium beta limit.

Finally, in this study, we discuss only the net-toroidal current free equilibrium. Since we observe the net-toroidal current, like the bootstrap current, ECCD current, and so on, in the experiment, it is necessary to discuss the 3D equilibrium including the net-currents. These are also future subjects of research.

## ACKNOWLEDGMENTS

This work was performed with the support and under the auspices of the NIFS Collaborative Research Programs NIFS07KLPH004 and NIFS12KNST045. This work was supported by Grant-in-aid for Scientific Research (C) 25420890 and (B) 18H01202 from the Japan Society for the Promotion of Science (JSPS). Also, this work was supported by “PLADyS,” JSPS Core-to-Core Program A., Advanced Research Networks.

## DATA AVAILABILITY

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