Effect of pressure profile on stochasticity of magnetic field in a conventional stellarator

In this study, the stochastization of the magnetic field in a simplified stellarator, which has a strong magnetic shear in the peripheral region, is studied using a nonlinear 3D equilibrium calculation code, HINT, which calculates the nonlinear resistive steady state using the relaxation method without the assumption of nested flux surfaces. For the finite beta equilibrium, the pressure-induced perturbed field breaks the nested flux surfaces, and the stochastization of the magnetic field appears in the region of the strong magnetic shear. The stochastization penetrates into the plasma core due to the increased beta and the effective plasma confinement degrades; in other words, the effective plasma volume shrinks. The stochastization strongly depends on the profile shape of the plasma pressure. Since the pressure gradient drives the Pfirsch–Schlüter current, the pressure profile affects the stochastization of the magnetic field. In order to understand the equilibrium beta limit due to the stochastization of the magnetic field, the profile shape of the plasma pressure is an important factor.


Introduction
Toroidal magnetic confinement research, and in particular, stellarator research , aims to generate and retain good flux surfaces in the finite high-beta equilibrium. In the tokamak, the perfectly nested flux surface in the finite beta equilibrium can be defined by the axisymmetry, ∂/∂ϕ = 0, but in the stellarator, the perfectly nested flux surface cannot be defined because of the lacking of axisymmetry. In the early history of stellarator research, it was important, but very difficult, to Original Content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. confirmthe generation and maintenance of good flux surfaces in the finite beta equilibrium. Many 3D equilibrium calculation codes were developed without assuming the existence of the nested flux surface, a priori, [1][2][3][4][5][6][7]9] and the high-beta 3D equilibrium is still studied now [10][11][12][13][14][15][16][17][18][19][20][21]. In the conventional theory, the equilibrium beta limit is defined by the Shafranov shift ∆/a, where ∆ = (R ax (β) − R ax (0)) and a is the minor radius. For the equilibrium beta limit of the circular cross section tokamak with large aspect ratio, if the Shafranov shift achieves about 0.5, the equilibrium is limited [22]. In the equilibrium beta limit, the separatrix appears in the plasma core, because the poloidal field is canceled by the external vertical field to keep the MHD equilibrium. The separatrix in the plasma core leads to the flattening of the electron temperature and also the degradation of the confinement. In the stellarator, the equilibrium beta limit can be defined by the Shafranov shift as well as the tokamak [23]. However, Hayashi et al pointed out the equilibrium beta limit is also defined by the stochasticity of magnetic field lines and it is more severe than the limit defined by the Shafranov shift [10]. 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 of the stellarator. Suzuki et al had proposed a new definition that the equilibrium beta limit is defined by the degradation of the effective plasma confinement due to the nonlinear 3D equilibrium response [24]. In that study, the high beta 3D equilibrium of a stellarator was studied and it was found that the strong stochasticity of the magnetic field cannot keep the pressure gradient in the peripheral region. Let us assume the beta sequences of the 3D equilibrium calculation. The relation between the central beta and volume-averaged beta is linear with the fixed pressure profile. However, because of the strong stochasticity of the magnetic field, the magnetic field cannot fix the prescribed pressure profile, and then the linear relation of the central beta and volume-averaged beta is violated. That means, although the central beta increases, the total stored energy equivalent to the volume-averaged beta cannot increase linearly due to the increased central beta.
The nonlinear 3D equilibrium response is a global effect due to the nonlinear coupling of resonant and nonresonant modes on rational surfaces [25]. On a rational surface, the resonant mode of the pressure-induced perturbed field causes the radial deviation of magnetic field lines. If the radial deviation is sufficiently larger and those deviated field lines are overlapping with each other, the magnetic field becomes stochastic. The nonlinear 3D equilibrium response has been studied by nonlinear 3D equilibrium codes not only for the stellarator but also for the tokamak [25].
Usually, the beta limit due to the nonlinear 3D equilibrium response is lower than the conventional beta limit defined by the Shafranov shift [24]. However, the new definition had only been studied for one stellarator, and it should be generalized to all stellarators. The pressure-induced perturbed field is driven by the Pfirsch-Schlüter (P-S) current. Since the P-S current is driven by the pressure gradient, ∇p, the nonlinear 3D equilibrium response may be sensitive to the profile shape of the plasma pressure. In previous studies, the impact of the profile shape of the plasma pressure has not been studied.
In this paper, we study extensively how the stochastization of the magnetic field changes the magnetic topology, in particular, the magnetic field structure in the peripheral region in a conventional stellarator. Special attention is paid to the impact on the profile shape of the plasma pressure on the stochastization of the magnetic field. The magnetic topology is studied by a nonlinear 3D equilibrium code, HINT [7,25], which does not assume nested flux surfaces. Using the HINT code, we can systematically study the stochastization of the magnetic field with different pressure profiles, which are peaked and board profile shapes. In the next section, we discuss a model stellarator configuration used in this study. We then study the nonlinear 3D equilibrium of the model stellarator. We study the impact of the profile shape on the stochastization of the magnetic field structure. In the last section, we summarize this study.

Model of conventional heliotron configuration
In this study, an L = 2/M = 10 simplified stellarator configuration to model the Large Helical Device (LHD) [27] is considered. Here, L is a pole number and M is a toroidal field period. Figure 1 shows a schematic view of the coil geometry of the stellarator configuration. This configuration consists of one pair of helical coils (light blue curves) and three pairs of axisymmetric poloidal coils (light green curves). For a reference, the plasma is also shown in the figure (dark blue surface). The major radius of the helical coil winding is 1.0 m, and the minor radius of the winding is 0.25 m. The helical coil pitch modulation, α, is 0.1. The R 0 B 0 is 1.0 mT, that is, the toroidal magnetic field, B 0 , at the center of the helical coil, R 0 = 1 m, is 1 T. In the real experiment device, the coil has a finite width and the cross section of the coil is a very sophisticated rectangle. Since those rectangular cross sections cannot be approximated by the single current filament, the sophisticated cross section of the coil is approximated by the multicurrent filaments. In such a case, the L = 2 separatrixes of current filaments are overlapping each other, and then the magnetic field structure around the L = 2 separatrix can easily become stochastic. For this reason, the magnetic field structure of the realistic stellarator configuration is very sophisticated in the peripheral region, for example, the LHD [28]. To understand the problem of how the pressure-induced perturbed field in the finite beta equilibrium stochastizes the magnetic field, the vacuum field must have a clear and nested flux surface. Therefore, the single filament approximation without the higher-mode component is used for the coil. Figure 2 shows a Poincaré plot of the vacuum magnetic field at three different poloidal cross sections. The toroidally averaged magnetic axis is 1.0 m and the effective minor radius, r eff , is 0.17 m. The effective minor radius, r eff , is defined by an equation, S = πr eff 2 , where S is the averaged area of the poloidal cross section along the toroidal angle, ϕ. The blue line and red points in the figures indicate the first wall and one pair of helical coils, respectively. The first wall is axisymmetric, and the minor radius is 0.24 m. The colours of the Poincaré plot indicate the connection length of the magnetic field line, L C . The magnetic field line is traced 10 000 m. A magnetic field line longer than 10 000 m can be approximated to the infinite. If a field line reaches to the first wall, that field line opens. Because of the single filament approximation, clear flux surfaces appear toward the last closed flux surface (LCFS), and magnetic islands or stochastic magnetic field lines do not appear in the vacuum field. In the outside of the LCFS, a thin stochastic region appears, but the L C is very long. Figure 3 shows profiles of the rotational transform, ι, and magnetic well/hill as a function of the effective minor radius, r eff . Horizontal dashed lines (blue color) indicate rationals on the n = 10 toroidal mode, which correspond to n/m = 10/20, 10/15, 10/10, 10/9, 10/8, and 10/7. A vertical dashed line (orange color) indicates a boundary of the magnetic well or hill. Rotational transforms at the magnetic axis and LCFS are 0.53 and 1.37, respectively, and the profile monotonically changes toward the LCFS. The magnetic shear becomes strong in the peripheral region, and the rotational transform profile crosses some rational surfaces in the peripheral region. Since distances of rational surfaces in the peripheral region are close to each other, it might be that the radial deviation of rational field lines, due to the pressure-induced perturbed field in the finite beta equilibrium, are overlapping each other [25]. On the other hand, for r eff < 0.094, the inside of the vertical dashed line, there is the magnetic well but for r eff > 0.094, there is the magnetic hill. These are the natural properties of this model configuration.

Impacts of pressure profile on stochastization due to nonlinear 3D equilibrium response
In the net current-free stellarator, the nonlinear 3D equilibrium response causes the pressure-induced perturbed field, which is driven by the P-S current, and that perturbed field makes an opening of the magnetic island and the stochastization of the magnetic field line. Since the P-S current is driven by the pressure gradient, ∇p, here, the profile shape of the plasma pressure affects the magnetic field structure.
In this study, the pressure profile is defined by the equation, where s is the normalized toroidal flux, p 0 is the pressure on the magnetic axis, and the profile shape is prescribed by two power factors, α and β. The peaking factor (pf) is defined by the following equation, where the ⟨ ⟩ means the volume averaged quantity defined by Note that in the HINT modeling the plasma volume is defined by an integral of a volume element, dV = RdRdϕdZ, for p > 0 on the cylindrical rectangular grid. In this study, three shapes of the plasma pressure are studied for the parabolic profile, (1 − s), the peaked profile, (1 − s) 2 , and the broad profile, (1 − s 2 ). Peaking factors of three profiles are 2, 3, and 1.5, respectively. To prescribe the profile of the plasma pressure accurately in the equilibrium calculation, the equidistant grid of ∆R = ∆Z~0.002 m is used on the poloidal cross section, and the equidistant grid of ∆ϕ~0.005 is also used in the toroidal direction. The central plasma beta, β 0 , and volumeaveraged beta, ⟨β⟩, are defined by, and respectively. Figure 4 shows the beta sequences of nonlinear 3D equilibrium calculations for the parabolic pressure profile, pf = 2. Poincaré plots at the horizontally elongated cross section of ϕ = 18 deg are shown for β 0~2 %, 4%, and 6%, respectively. Values of ⟨β⟩ correspond to 1%, 2%, and 2.8%, respectively. For a comparison, the vacuum field is shown at the same cross section too. The colours of the Poincaré plots show the connection length, L C . Vertical lines indicate the position of the magnetic axis along the R coordinate. Also, radial profiles of the rotational transform and normalized plasma pressure, p/p 0 , are shown at the same poloidal cross section on the Z = 0 plane. The rotational transform and normalized pressure are plotted on only the closed field line (not reaching to the first wall). Vertical dashed lines indicate the position of the LCFS on the Z = 0 plane. According to the increased β 0 , the magnetic axis shifts outward because of the Shafranov shift. However, the nonlinear 3D equilibrium response causes not only the Shafranov shift but also the stochastization and opening of the magnetic field line in the plasma edge. For β 0~2 %, magnetic field lines on m = 8 and 9 rational surfaces begin the stochastization, but the L C of the stochastic field line is long enough, which can be approximated to the almost closed field line. In the HINT modeling [25,26], the plasma pressure in the relaxation process is approximated by the following equation, where i means a step number of the relaxation process, the L C is the connection length of the magnetic field line starting from the computational grid, and L in is an input parameter to control the length of the field line tracing. In this study, L in is 100 m. This equation calculates the averaged plasma pressure on the flux tube. If the flux surface is perfectly closed, the averaged plasma pressure can equilibrate on the closed surface as the flux surface quantity. Also, if the magnetic field line is stochastic but L C is longer than L in , the averaged plasma pressure corresponds to the plasma pressure on the 'averaged' flux surface, where the magnetic field line travels in the high and low plasma pressures. However, if the magnetic field line reaches the first wall or L C is shorter than L in , the averaged plasma pressure is set to zero. This means that in the strongly stochastized magnetic field, where many magnetic field lines are opening and the L C is very short, the plasma pressure is diffused and it is not possible to keep the pressure gradient due to the prescribed profile. Therefore, for β 0~2 %, it can be considered as the stochastic region to effectively confine the plasma. However, for β 0~4 %, magnetic field lines on m = 7-12 rational surfaces are almost opening and the penetration of the stochastic region begins into the plasma core. Thus, the pressure gradient in the outside of the LCFS becomes weak compared with the gradient on the inside of the LCFS. For β 0~6 %, the stochastic region widely penetrates into the plasma, and the region of the closed flux surface shrinks to the small region in the plasma core. In the stochastic region, the mixing of the closed and opening field lines appears, and the plasma pressure diffuses in the stochastic region. The pressure profile cannot be kept by the prescribed pressure profile in the stochastic region. That means that the effective plasma confinement degrades significantly and the total plasma stored energy may be decreased.
To understand how the stochstization and opening of the field line degrade the effective plasma confinement more clearly, L C of the stochastic field line is studied. Figure 5 shows radial profiles of the L C for β 0~2 %, 4% and 6%. Also, the L C of the vacuum field is also shown. Vertical lines indicate positions of the LCFS on the Z = 0 plane. The outside of the LCSF can be considered as the stochastic region. L C is calculated by the field line tracing, which starts from the equatorial plane on Z = 0, at the same poloidal cross section corresponding to figure 4. According to the increased β 0 , the field line in the peripheral region begins opening to the first wall. Increasing the β 0 , the region of closed field lines shrinks, and opening field lines penetrate into the plasma core. In such a sense, the boundary of opening or closed field lines can be considered as an approximation of the effective LCFS. For β 0~2 %, the magnetic field line in the plasma edge begins opening, the L C in the plasma edge, in particular, on n/m = 10/8 and 10/7 rational surfaces, is still closed. On the other hand, for β 0~4 % and 6%, the magnetic field line on those rational surfaces are already opening and the radial deviation of the stochastic field line is also large in both cases. Note that the stochastic region for β 0~4 % is not large but the stochastic region for β 0~6 % is large, where the stochastic region is wider than the region of the closed flux surface. Under the condition of a sufficiently long electron mean free path, λ e (κ ∥ = ∞), the plasma pressure immediately diffuses, and then the parallel pressure gradient, ∇ ∥ p, is zero. However, under another condition of κ ∥ ̸ = ∞ and L C > λ e , the parallel pressure gradient may remain, depending on the ratio of κ ∥ /κ ⊥ , if the radial deviation of the stochastic field line is not large. From this consideration, the stochastic field line of the long L C , where the radial deviation is not large, for cases β 0 ≤ 4% can retain the pressure gradient, but the strongly stochastic field line with the large radial deviation (β 0~6 %) clearly limits the pressure gradient. That is consistent with the result in figure 4.
The P-S current is driven by ∇p. Therefore, the nonlinear 3D equilibrium response due to the P-S current is affected by the profile shape of the plasma pressure. To see impacts of the pressure profile on the nonlinear 3D equilibrium response, 3D equilibria of different pressure profiles are studied. In figure  6, three different pressure profiles are compared for β 0~2 % (left column) and 4% (right column), respectively. Peaking factors (pf) of the pressure profile are 2 (top row), 3 (middle row), and 1.5 (bottom row), respectively. Poincaré plots for pf = 2, 3, and 1.5 are shown at the same poloidal cross section corresponding to figure 4. Vertical dashed lines in the figures indicate the magnetic axis position for pf = 2. For β 02 %, the small difference in the Shafranov shift appears. However, the magnetic field structure in the plasma edge is significantly different, in particular, near two rational surfaces on n/m = 10/9 and 10/8. For pf = 3, the pressure gradient is not steep on those rational surface, and the stochastization is relatively small. However, for pf = 1.5, the number of opening field lines increases compared with the cases of pf = 2 and 3. On the other hand, for β 0~4 %, the Shafranov shift is different according to the different peaking factor. The magnetic field structure in the peripheral region is more obvious. For the steep pressure gradient in the peripheral region, the stochastization of the magnetic field becomes strong, and penetrates into the plasma core. In particular, for pf = 1.5, the number of opening and short field lines in the peripheral region increases compared with the other two cases.
From the above considerations,an important hypothesis can be concluded. That is, if the equilibrium beta limit in the stellarator is defined by the stochastization of the magnetic field, it does not depend on only the Shafranov shift, but also on the profile shape of the plasma pressure. There is a possibility that the equilibrium beta limit is relatively low, although the Shafranov shift and the β 0 are low. A question arises as to how that can be seen. To see that, the relationship between the central beta, which is the beta value on the axis, β 0 , and the ⟨β⟩ is studied. In the low beta equilibrium, the relationship between β 0 and ⟨β⟩ is linear. However, as discussed above, the stochastization of the magnetic field in the plasma due to the increased beta may degrade the confinement, because the effective plasma volume shrinks. Figure 7 shows ⟨β⟩, as a function of β 0 for three different peaking factors of the pressure profile. Dashed lines in the figure indicate linear fittings for the low beta values (β 0 < 2%). For pf = 3, the stochastization of the edge magnetic field is relatively small because the pressure gradient in the edge region is relatively small. Also, the relationship between β 0 and ⟨β⟩ is almost linear up to β 06 %. However, for pf = 2 and 1.5, the linear relationship of β 0 and ⟨β⟩ is violated for the high beta equilibrium (β 0 > 4%). Also, from the high beta (β 0~4 %), the gradient for pf = 1.5 is relatively small compared with the gradient for pf = 2. Since the stochastization for pf = 1.5 is stronger than the stochastization for pf = 2, the degradation of the effective plasma confinement appears strongly for the broad pressure profile, which has a steep ∇p, in the peripheral region. It suggests a possibility that, for the conventional stellarator, the equilibrium beta limit for the broad pressure profile might be lower than the limit for the peaked pressure profile. Since the high central beta value is important to the alpha heating in the fusion reactor, the high equilibrium beta limit of the peaked pressure profile is favorable in the design of the fusion reactor.

Discussion and summary
We have studied how the magnetic field becomes stochastic and how the profile shape of the plasma pressure affects the stochastization of the magnetic field, for a model of a conventional stellarator. In the model configuration, the magnetic shear is strong in the peripheral region, and then the pressureinduced perturbed field leads to the stochastization of the magnetic field line therein. Increasing the stochasticity of magnetic field lines, the linear relation between β 0 and ⟨β⟩ is violated. That means the effective plasma confinement degrades. However, although the effective plasma confinement degrades due to the increased stochasticity, this degradation does not cause a termination event like the disruption of the tokamak.
As mentioned above, the conventional stellarator in this study has a strong magnetic shear in the peripheral region. Thus, if the steep pressure gradient drives the strong P-S current in that region, the magnetic field easily becomes stochastic. This is a reason that, for pf = 1.5, the strong    stochastization in that region, and the degradation of the effective plasma confinement, appear. In the conventional theory, the equilibrium beta limit is defined by the Shafranov shift. However, by nonlinear 3D equilibrium analyses, which allow the appearance of the stochastization of the magnetic field and magnetic island, the equilibrium beta can be limited by the degradation of the plasma confinement due to the stochastization of the magnetic field. In such a case, a key factor is the pressure profile; in other words, ∇p on the strong magnetic shear region. Since the pressure-induced perturbed field due to ∇p drives the radial deviation of the magnetic field line, that radial deviation of neighboring magnetic field lines can be overlapped easily. Therefore, the equilibrium beta limit should be considered by a coupling of the rotational transform and pressure profile. Loizu et al studied the equilibrium beta limit due to the nonlinear 3D equilibrium response using SPEC code [9]. SPEC could reproduce the huge stochastization and the degradation of the flux surface leads to the beta limit for the conventional stellarator [20]. However, since SPEC uses the stepped pressure profile, detailed studies of the stochastization by different pressure profiles are still difficult.
The target configuration of this study is the conventional stellarator. However, the conclusion of this study can extend to the optimized stellarator because the P-S current cannot be suppressed completely. Also, since the optimized stellarator has a low magnetic shear, the pressure-induced perturbed field can cause the large radial deviation of the magnetic field line. It means that, in the optimized stellarator, the steep pressure gradient on the low-order rational surface should be avoided. In previous studies of Wendelstein 7-X (W7-X) [17,19,29], the stochastization due to the nonlinear 3D equilibrium response appears at the high beta. In those studies, only a limited number of pressure profiles is studied. Thus, the impact of the profile shape of the plasma pressure on the stochastization will become an important issue for a future subject.
The stochastization of the magnetic field line in the plasma is an intrinsic property in the stellarator. Since the transport strongly links the MHD equilibrium, it seems that the stochasticity affects the transport processes. Rechester et al suggested that the stochastization of the magnetic field affects the electron heat transport [30]. Thus, there is a possibility that the stochasticity limits ∇p. Suzuki et al studied the heat transport by stochastic field line based on the Rechester-Rosenbluth model [26]. In HINT code, the transport effect of stochastic field lines is modeled by the perpendicular transport of the stochastic field. Thus, the transport of stochastic field lines strongly depends on a parallel correlation length, where the length of the field line tracing is a given parameter a priori. To treat the transport more physically, the anisotropic heat transport must be involved. The anisotropic heat transport has been studied already [31] but the implementation into HINT code has not yet been completed. Also, a study of the stochastization with the anisotropic heat transport using another 3D nonlinear code, NIMROD [32], has been started. The benchmarking of HINT and NIMROD will give a more extensive understanding of the stochastization of the magnetic field and the equilibrium beta limit.
In this study, we discuss only the zero toroidal current equilibrium. The toroidal current may affect strongly the stochastization of the magnetic field in the stellarator. In general, the optimized stellarator minimizes the P-S current, but the optimized stellarator still has significant bootstrap current. The effect of the toroidal current must be discussed. In a previous study [19], it was found that the toroidal current 'density' is more important than the toroidal current. This also suggests the importance of the profile shape of the toroidal current. The self-consistent nonlinear 3D equilibrium calculation including the toroidal current density is another important future subject. Also, in this study, the stochastization driven by MHD instabilities are not discussed. In the realistic situation, the magnetic field can be changed by not only the equilibrium but also the instability. Since a driving term of the MHD instability is ∇p, not β 0 , the degradation due to the stochastization by the MHD instability may begin from the low beta. Therefore, the beta limit due to the confinement degradation due to the MHD instability would be lower than the beta limit defined by the equilibrium. An integration of the equilibrium and instability on the magnetic field structure is an important issue. It will be discussed in different papers.