Active generation and control of radial electric field by local neutral beamlets injection in tokamaks

The radial electric field plays an important role in plasma confinement in tokamaks and can be generated through neutral beam injection. In this study, we propose a model for calculating the radial electric field resulting from tangential local neutral beamlet injection, aiming to externally control and improve plasma confinement. The Neutral beamlet ion and Energetic particles Orbit mover and Electric field solver code has been developed to analyze this issue, and its simulation results have been validated against results from other codes as well as measurements from correlation reflectometers. The charge separation is primarily caused by the redistribution and loss of beam ions due to magnetic gradient and curvature drift as well as collision effects, and it is maintained through continuous beamlet injection. The electric field is calculated using Poisson’s equation, taking into account both classical and neoclassical polarization effects. The results demonstrate that despite the high losses and low heating efficiency associated with localized beamlets, they are capable of generating a significant radial electric field characterized by a steep gradient. This presents opportunities for external control of the electric field, potentially enhancing plasma confinement.


Introduction
The effect of the radial electric field on the plasma confinement in a tokamak is of widespread concern.Prior research has investigated the connection between the L-H transition and the radial electric field [1,2].The E × B flow shear resulting from the radial electric field has the potential to mitigate turbulent transport and enhance plasma confinement [3,4].In the double transport barrier experiment, the internal radial electric field enhances confinement performance within the plasma's central region, leading to internal transport [5].This characteristic of radial electric fields enables the utilization of electric fields to enhance confinement.
Itoh proposed the existence of bifurcated solutions for the boundary non-dipolar particle loss [6].Shaing et al provided an explanation for the phenomenon of spontaneous bifurcation in terms of the poloidal momentum balance equation [1].In 2000, Heikkinen provided insight into the formation process of the negative radial electric field at the boundary by solving the kinetic equation, taking into account ion orbit loss [7].Similarly, Chang et al and Wu et al developed their own ion orbit loss models to explain the generation of the negative radial electric field [8,9].
A common feature among these theories is that the requirement of charge separation within the corresponding region for the generation of a radial electric fields.However, the driving factors for charge separation differ, including anomalous particle transport, ion orbit loss at the boundary, temperature/density gradients, etc.These factors play a crucial role in the formation of charge separation and radial electric field.Nonetheless, this process is typically spontaneous and lacks active control.Furthermore, neutral beam injection (NBI) experiments have also revealed the generation of a radial electric field [10].
This article presents a preliminary investigation into the potential and viability of actively generating and controlling radial electric fields using neutral beams.The neutral beam injects high-energy neutral particles into the plasma.Since the neutral particles are not charged, the high-energy neutral particles can freely traverse the magnetic field and be injected into the plasma.These neutral particles interact with the background plasma, resulting in charge exchange or impact ionization, which then produces fast ions and electrons.As a result of magnetic drift, the fast ions undergo compression or expansion orbits and significantly diverge from the magnetic surface, while electrons remain closer to it.Disparate orbits between fast ions and electrons give rise to charge separation, ultimately establishing an electric field.If instability associated with the rapid rotation of the cold plasma [11] is not taken into consideration over the timescale of beam ion slowing down, the radial transport of background particles is insufficient to offset the charge separation caused by continuous NBI, resulting in the persistence of the radial electric field.
In order to achieve high plasma temperatures, neutral beams used for heating purposes often employ a larger beam cross-sectional area to increase the heating power.For example, in EAST, the beam cross-section measures 12 cm × 48 cm.The beam requires focused injection with a focal length of 10 m and horizontal and vertical divergence angles of 0.6 • and 1.2 • .NBI with a large beam cross-section and a fixed on-axis injection angle generates fast ion distributions spanning the entire radial space, resulting in a relatively fixed distribution of the radial electric field.Additionally, due to the wide spatial dispersion of the beam ions within the plasma, certain small-scale charge separations may be shielded by beam electrons or ions born from other regions, thereby weakening the strength of the electric field.In the case of instabilities occurring within the plasma core, such as neoclassical tearing modes, the suppression of the instabilities requires the generation of a radial electric field within a narrow region, such that the electric field shear exceeds the instability growth rate.In such scenarios, the use of conventional beam injection cannot produce the required electric field to suppress the instabilities effectively.
To achieve targeted deposition of beam ions and induce charge separation at specific radial positions, off-axis injection using beamlets with small cross-sections is required.Furthermore, by exerting control over the beam injection angle, it may be possible to selectively generate the radial electric field in desired regions, consequently providing an effective means of suppressing instabilities.
In this paper, we propose a method for controlling the generation of a radial electric field through local tangential injection of neutral beamlets.These beamlets are smaller in size and do not require focusing.The schematic diagram illustrating tangential local beamlet injection is presented in figure 1.Specifically, the beamlet is injected tangentially on the low field side.It is anticipated that the beam ions will mainly deposit within a specific magnetic surface (the yellow region), as shown in figure 1(c).These particles are expected to induce charge separation through compression or expansion orbits, thereby leading to the formation of radial electric fields in targeted areas.To describe the injection position of the local beamlet, we introduced the parameter r, which represents the minimum distance from the beamline to the magnetic axis.
The radial electric field produced by the local beamlet is closely influenced by the beam's parameters and the equilibrium magnetic field.By systematically optimizing these parameters, the electric field generated by local beamlet injection may achieve or even exceed the strength and shear rate observed in the H-mode electric field.This may lead to the formation of an internal transport barrier and enhances the confinement performance.To investigate this issue, the Neutral beamlet ion and Energetic particles Orbit mover and Electric field solver (NEOE) code has been developed.This paper is structured as follows: section 2 presents the basic model of the NEOE code.In section 3, the verification against the NUBEAM code [12], and validation against experimental data from the EAST 2020 campaign are provided.Section 4 consists of a detailed simulation analysis of the generation of radial electric fields by local beamlets using the code.Finally, in section 5, we summarize our conclusions.

Computational code for neutral beam particles analysis
On the foundation of the GYCAVA [13] trajectory and electric field computation modules, several additions have been added.These include a geometric system for neutral beamlet injection, a neutral particle deposition module, and an angular momentum calculation module.We have integrated these code components into NEOE code to facilitate our research on beamlet injection-driven radial electric fields.

Neutral beam ionization
A neutral beam fast ion tracking code has been developed to investigate the properties of neutral beam fast ions.The code sets the spatial distribution of neutral particles using a double Gaussian distribution at the source.The energies of the injected particles are randomly selected with probability proportional to the measured beam intensity at E, E/2 and E/3.The Monte Carlo method [14] is employed to determine the deposition positions of the neutral particles in the plasma.Each particle is associated with a random number η, uniformly distributed between 0 and 1.Subsequently, the particles are pushed along the direction of beam incidence and the code performs the integral: In the equation, ν represents the particle ionization macroscopic reaction cross section [15]: n i and n e represent the ion and electron densities, respectively.σ ch denotes the cross section for charge exchange reactions, σ ii represents the cross section for ion impact ionization reactions, and σ e denotes the cross section for electron impact ionization reactions.Meanwhile, v e and v b refer to the velocities of electrons and beam particles, respectively.The term ⟨v e σ e ⟩ corresponds to the average rate of electron impact ionization under the Maxwell distribution.The reaction cross-section data were obtained from [16].At the point along the trajectory where I = ln(1/η) the beam neutral is taken to be ionized and is converted to a fast ion ready to be treated by the orbit code.

Calculation of particle guiding-center orbit
The fast ion orbits are computed using gyrokinetic equations based on canonical variables, similar to the approach employed in the GYCAVA code [17], which includes collision effects [18].The ion orbit calculated by the code is illustrated in figure 2, where the particle undergoes motion across the magnetic surface due to magnetic drift.The ion orbits in the co-current (plasma current) direction shift inward from the low-field side, while the ion orbits in the counter-current direction shift outward.For beam fast ions, the orbit direction (co/counter current) can be obtained by altering the NBI injection direction or reversing the plasma current.
The maximum distance for fast ions to migrate from the birth magnetic surface can be estimated by [19]: Here, the maximum crossing distances of passing ions and trapped ions are denoted as d p and d t respectively.The safety factor is represented by q, while the inverse aspect ratio is given by ϵ = a/R, where a and R correspond to the minor and major radii of the device respectively.The Larmor radius is defined as r L = mv ⊥ /ZeB.This indicates that ions with larger Larmor radii are more affected by magnetic drift, resulting in greater track offsets.In contrast, electrons with smaller Larmor radii exhibit negligible deviations from the magnetic surface.The difference in the radial deviation directions between the co-current and counter-current beam ions results in the formation of a positive or negative radial electric field.
The code statistics all the particle orbits, obtains different distributions of electrons and ions, and transfers this information to the electric field calculation module for determining the electric field.

Calculation of the radial electric field produced by NBI
After a certain slowing-down time, the newly injected beam ions and the thermalized beam ions reach an equilibrium state, resulting in a slowed-down distribution of beam ions.Regardless of the anomalous transport of beam electrons, the distribution of the beam electrons is assumed to be the same as the birth distribution of the beam ions.Using the distribution of beam ions and electrons, the electric field can be determined by solving the Poisson equation: Here, we categorize the particles into four species: the ions and electrons of the background plasma, denoted by the subscript j = bulk; the ions and electrons generated by NBI, denoted by the subscript j = beam.In this paper, we neglect the influence of α particles and impurity ions, which are minority particles.Z j represents the charge number of the ion.n i,j and n e,j denote the correspondingion and electron densities, while ε 0 represents the vacuum permittivity.In bounce average gyrokinetics, the particle density in the presence of an electric field can be expressed as [20]: where N gc,j represents the particle density in the absence of electric field disturbance, and n cl,j and n nc,j denote the classical and neoclassical polarization densities, respectively.The classical polarization density accounts for the discrepancy between the actual particle density and the gyro-averaged density under the influence of an electric field.It can be expressed as [21]: The neoclassical polarization density comes from the difference between the gyroaverage density and the bounce average density in the presence of an electric field.Taking into account the effect of the elongation ratio κ, the corrected expression for the neoclassical polarization density is [22]: When κ = 1, the magnetic field configuration corresponds to the circular configuration.It can be observed that as the elongation ratio κ increases, the impact of neoclassical polarization diminishes.Since the polarization density of electrons is negligibly small, it is also disregarded in this paper.According to the above relationship, the Poisson equation ( 5) can be written as: Due to our focus on the electric field produced by the beam ions and electrons and the assumption of quasi-neutrality in the background plasma, we neglect the electric field arising from the charge separation of the background ions and electrons.
The expression for the radial electric field E r induced by NBI can be simplified as follows [23]: ˚−e (Z beam (N gc,beam +n nc,beam ) − n e,beam ) Jdψ dθdζ, (10) Here, J represents the Jacobian matrix, and S denotes the magnetic surface area in flux coordinates.The left-hand side of the equation (10) consists of three terms.The first term represents the relative permittivity in vacuum.The second term ) accounts for the relative permittivity resulting from the classical polarization effect of background ions and beam fast ions.The third term 3.27 ϵ 1/2 represents the neoclassical polarization effect of background ions.
The neoclassical polarization term n nc,beam for the beam ions remains in the density term on the right-hand side of the equation.This is because we utilized the computed radial electric field in to recalculate the particle orbits, and the neoclassical effects arising from magnetic drift are already incorporated in the provided densities.It has been observed that the radial electric field has minimal impact on the orbit of fast ions, and the change of the fast ions distribution is very small.Furthermore, in the background plasma, the ratio ϵ 1/2 ) dominating the shielding effect of the electric field.In the NEOE code, equation ( 10) is utilized to compute the radial electric field induced by local beamlet injection on the slowing-down time scale.

j × B torque injection
Taking the time partial derivative of equation ( 9), integrating it over space, and simplifying the expression yields: (12) The variable j represents radial currents.NBI results in the active generation of fast ion radial flow denoted as j b .This portion of the current begins to form at the moment of NBI injection.Simultaneously, background ions respond with a polarization current jp in the opposite direction to j b .As discussed above, the polarized current from the beam fast ions is negligible and can be disregarded.The term j NA refers to the non-ambipolar diffusion flow of background ions.Assuming that the background plasma remains quasi-neutral on the slowing-down time scale, j NA due to the diffusion of background particles can be neglected.
As the NBI continuously generates j b , a portion of j b is shielded by the opposing polarization current jp.The remaining current generates an electric field.The NEOE code calculates the residual electric field resulting from the removal of the background ion polarization flow by the fast ion radial flow, without considering the impact of non-ambipolar diffusion.This computational approach aligns with Hinton's approach [24], but incorporates a modified polarization coefficient that considers the classical and neoclassical polarization effects of background ions under a non-circular cross-section magnetic field.
The polarization flow (jp) induced by the plasma generates a j p × B θ torque, which contributes to plasma toroidal rotation.It is worth noting that, keeping the geometry of NBI injection unchanged, changing the direction of plasma current will simultaneously change the direction of B θ and jp (that is determined by co-or counter-current NBI), and the direction of the j p × B θ torque would remain unchanged.Additionally, the fast ions transfer momentum to the plasma through collision.These two processes dominate the toroidal torque injection by the neutral beam, which can be expressed as: Here, angle brackets denote magnetic surface averaging.The origin of the polarization current term ⟨RjpB θ ⟩ can be attributed to the charge separation of the beam fast ions, which occurs due to orbit deviations or losses caused by magnetic drift and the pitch angle scattering in the collision term.The slowing down term ⟨Rê ϕ • F col ⟩ arises from the collisional interaction between fast ions and background particles.
In the case that the time scale is significantly smaller than the momentum confinement time and the ion viscosity term is neglected, the toroidal momentum equation is written according to equation (3) of the [24]: According to the equation ( 14), the NEOE code calculates the plasma torque by analyzing the change in the beam ion's toroidal angular momentum over a short period of time (0.12 ms): Within a brief duration (approximately 10 ms) following NBI injection, the torque predominantly arises from the j × B effect.Estimating the rate of momentum change allows us to estimate and validate the calculated results of the radial electric field.

Comparison of torque injection results
A comparison has been performed between the momentum injection results obtained from the NEOE code and NUBEAM code.The comparison uses the NBI parameters [25] and plasma equilibrium parameters from the EAST discharge #77759.The plasma profile for this discharge is shown in figure A1 in the appendix A. The simulations were carried out for a previous EAST setup, predating 2020.EAST was equipped with two sets of neutral beam injectors, with each injector having a left and right source.The co-current injection sources are AL (tangency radius is 1.264 m) and AR (0.731 m), while the countercurrent injection sources are FL (0.606 m) and FR (1.141 m). Figure 3 presents the birth distribution of the beam ion resulting from tangential injection AL, calculated using NEOE and NUBEAM codes.
Within a slowing-down time of approximately 50 ms, the torque generated by NBI is calculated using equation (15) in NEOE, and the comparison with the results obtained from NUBEAM is shown in figure 4. The positive value of torque indicates the direction of the plasma current, which is counterclockwise from a top-down view.
It can be seen from the figure 4 that the distribution and trend of the injected torque density are roughly the same.In figure 4(a), with the co-current injection, the j × B torque is greater than zero in the region of ρ < 0.6, while in the outer region of ρ > 0.6, the torque is negative due to orbit loss.In figure 4(b), with the counter-current injection, the j × B torque is negative, and the distribution is consistent with NUBEAM.The difference in value between the two distributions may result from variations in collision coefficients, orbit following models, and numerical errors.
The results of total injection torque and loss fraction are presented in table 1, showing relatively consistent values.By using formula (12), the time integral of the j × B torque also corresponds to the radial electric field, thereby indirectly validating the accuracy of the radial electric field calculation.

Comparison with experimental results
The radial electric field generated by NBI drives the plasma to rotate with a velocity perpendicular to the magnetic field, as described by the E × B/B 2 rotation.By measuring the vertical rotation velocity of the plasma in experiments, we can estimate the radial electric field: EAST is equipped with a set of X-mode polarized W-band correlation reflectometers with four different probing frequencies (79.2 GHz, 85.2 GHz, 91.8 GHz, 96 GHz).This system utilizes multi-frequency microwave injection, enabling simultaneous measurements of density fluctuations and turbulence vertical velocity at 2 (poloidal) × 4 (radial) positions [26].
In the experiments, the W-band correlation reflectometer has been employed for turbulence measurements in the plasma core region.Figure 5 presents an example of the measurements.In the shot #80757 of the L-mode operation, the plasma current (Ip) was 450 kA and the toroidal magnetic field (B 0 ) was 2.37 T. Co-current NBI and Electron Cyclotron Resonance Heating (ECRH) were utilized as auxiliary heating methods.The basic plasma parameters of this discharge, including the plasma line-averaged density, auxiliary heating power, and plasma energy, are shown in figures 5(a)-(c).
From 3 to 6 s, co-current tangential NBI was continuously applied, pulsed NBIs of 1 MW were performed at 3.2 s, 4.2 s, and 5.2 s in the near-perpendicular direction.ECRH with a power of 1 MW was mainly used to heat the background electrons, having minimal impact on the radial electric field generated by the neutral beams.The core temperature of the ion was approximately 1.04 keV, and the electron density n 0 was around 3.4 × 10 19 m −3 .The density profiles measured by the profile reflectomtery show a relatively small increase during the NBI heating, which is consistent with the trend of the lineaveraged density.
As shown in figures 5(d) and (e), the spectrum of reflected signals for the 91.8 and 96 GHz probing frequencies shows an modulation of the spectral width of high-frequency fluctuation, which indicate the significant modulation of velocity of turbulence in the core due to the NBI.The vertical velocity of turbulence can be calculated using the following equation: Here, d represents the distance between two measurement points with the two different poloidal positions on the cutoff  surface corresponding to the plasma frequency [27], and τ denotes the propagation delay of turbulence between these two measurement points.Typically, the phase velocity extracted from the reflectometer includes the following aspects: In the equation, v phase represents the intrinsic phase velocity of the density fluctuations and the v phase term is usually much smaller than the v E×B term in tokamaks [28].The estimated vertical velocity of the turbulence, obtained through poloidal correlation analysis [27], is shown in figure 6.It should be noted that the value of the error bar only takes into account the error of the correlation analysis.The previous work [29] showed that only when the vertical velocity of turbulence is large enough (>1 km s −1 ), the estimated velocity will be close to the actual plasma velocity.In contrast, when the vertical velocity of turbulence is small (<1 km s −1 ), the estimated velocity will be higher than the actual plasma velocity.
During the ECRH heating periods at 2.8 s and 6.5 s, the vertical velocity of the core turbulence is approximately −2 km s −1 , where the negative sign indicates the rotation in the electron diamagnetic drift (EDD) direction.
After the application of NBI, noticeable modifications are observed in the vertical velocity of turbulent structures.At 4 s, 5 s, and 6 s, the injection of 1 MW tangential NBI generates a radial electric field that induces rotation in the direction of ion diamagnetic drift (IDD).The plasma velocity shifts the direction from EDD to IDD, with a velocity below 1 km s −1 , could potentially cause an overestimation of velocity in the co-tangential case.The acutal rotational velocity of the turbulence might be approximately 0 km s −1 .Additionally, the subsequent injection of a pulsed 1 MW perpendicular NBI causes the radial electric field to drive the turbulent structures to rotate in the IDD direction, resulting in an increase in velocity to around 2 km s −1 .
We assume a constant relative rotation of turbulence with respect to the plasma, denoted as v phase .By calculating the variation of rotational velocity relative to different time points, we obtain the velocity of v E×B (i.e.radial electric field) for different NBI scenarios.Additionally, we employ NEOE to compute the radial electric field distribution for NBI injection under the same experimental parameters.The comparison between the experimental and simulation results is presented in figure 7, with ρ is the square root of the normalized toroidal magnetic flux.
The rotational measurement data from the core region provides radial electric field profiles at three time points with NBI.The radial electric field generated by co-perpendicular NBI was obtained by comparing the velocity variation at the moments when both co-perpendicular and co-tangential neutral beams were injected (at 3.5 s, 4.5 s, and 5.5 s) with the moments when only co-tangential NBI (at 4 s, 5 s, and 6 s).The results are shown in figure 7.
The results of co-tangential NBI are displayed in figure 7(a).The results at R ∼ 2.02 m (ρ ∼ 0.5, circles in figure 7(a)) within the core exhibit larger measurement errors attributable to inaccuracies in the measurement system.There might be several possible reasons for the larger error bars: (1) the system response is different between these two probing frequency points; (2) if the antenna is not exactly perpendicular to the magnetic field, it will introduce additional errors.This same applies to turbulent structures that are not aligned in the direction of the magnetic field [30]; (3) variations in phase velocity between the two probing locations might also contribute to this discrepancy.Additionally, it also brought uncertainty, especially when the velocity was close to 0 km s −1 in the co-tangential case.
Figure 7(b) presents the results of injecting coperpendicular NBI.The measured data at R ∼ 2.08 m (ρ ∼ 0.67, squares in figure 7(b)) under 1 MW NBI conditions agree well with the calculated results.Figure 7(c) presents the results of injecting both co-perpendicular and co-tangential NBIs.This analysis involved comparing velocity variations at 3.5 s,  The measured data at R ∼ 2.08 m (squares in figure 7(c)) align well with the calculated results.The theoretical calculations at R ∼ 2.02 m (circles in figure 7(b)) fall within the experimental error range, and the gradient directions of the profiles are also consistent.Overall, despite the larger errors in the inner location, the gradient direction of the profiles still can be found in the case of co-perpendicular and co-tangential+coperpendicular and the theoretical calculation aligns with the experimental results.
Although the experiment observed the influence of neutral beams on the rotational velocity of turbulence, the simulation results indicate that the shear generated by NBI remains fixed at a radial position.It is not suitable to generate great velocity shears in specific regions or actively control the location of velocity shear within the turbulence structure.Therefore, we propose a new method of localized neutral beamlet injection to selectively induce the desired electric field shear within specific regions.

Simulation setup
In this section, we employ a validated NEOE code to calculate the generation and control of radial electric fields resulting from local neutral beamlet injections.In the simulations presented in this section, we employed parameters resembling those of EAST tokamak.The major radius is 1.84 m, the  minor radius is 0.44 m, and the elongation is 1.6.The equilibrium magnetic field is shown in figure 8.The magnetic axis is located at R = 1.91 m, where the toroidal magnetic field is Bt = 2.44 T and the toroidal current is Ip = 400 kA.The gap between separatrix position and outer mid plane limiter is 7.5 cm.The plasma profiles used in our simulations are shown in figures 9(a)-(c).The central electron temperatures and densities are 4.2 keV and 3.4 × 10 13 cm −3 , respectively.The safety factor profile is shown in figure 9(d), and the boundary safety factor is q 95 = 6.51.
The beamlet section width was set to 3 cm, and the height was set to 12 cm.Deuterium ions were used as the beam ions, with a total energy of 50 keV.Taking into account the neutralization efficiency and transmission losses of approximately 50%, the neutral beam density was set to 0.12 A cm −2 and the corresponding power of the beamlet is 0.22 MW.The energy components of the beam were distributed in ratios of 0.8:0.14:0.06 for full, one half, and one third energy, respectively.Figure 10 illustrates the layout of co-current and counter-current beamlets at injection positions of r = 0.1 m, 0.2 m, and 0.3.In each case, a total of 20 000 particles were tracked to calculate their motion for a duration of 80 ms, which is sufficiently long for all markers to thermalize or be lost.The computation time for each case, utilizing a machine with 80 cores, was approximately 1 h.

Distribution of fast ions from Local Beamlets
The contour plot in figure 11 illustrates the birth distribution of co-current local beamlets injected in different directions.The birth distributions for co-current and counter-current injections are similar.It should be noted that the beamlet's focus and divergence effects are not considered in the analysis.This is due to the fact that beamlets do not require convergence and have relatively low divergence.The small divergence has a negligible effect on the beam ion deposition, and it could be mitigated by using suitable collimators or apertures to achieve a low divergence beam.
Due to the small height (12 cm) of the local beamlet, beam ions tend to deposit in a narrow region along the Z-direction.The deposition distribution of beam ions in the R-direction varies depending on the off-axis injection direction of the local beamlet, with ions depositing in the plasma as they travel from the outer to inner region and then to the outer region again.When the injection position is r = 0.1 m, the plasma experiences the longest deposition region, leading to a broader distribution of beam ion deposition in the R-direction.As the injection radius r increases, the radial deposition region decreases in size.
The deposition density of beam ions is influenced by the plasma density.As the plasma density decreases toward the plasma edge, the deposition density of beam ions in the edge region also decreases, particularly for beams injected at larger values of r.Consequently, the shine-through loss increases.Simulations indicate that the shine-through losses  for co-current beamlet injection at positions 0.1 m, 0.2 m, and 0.3 m are 0.81%, 3.1%, 15.78%, respectively.The shinethrough losses for co-current or counter-current beamlets yield similar results.The estimated thermal load on the wall facing the beam injection section are 0.12 MW m −2 , 0.47 MW m −2 , 2.4 MW m −2 respectively.
The orbit loss fraction induced by local beamlet injection in co/counter-current directions are shown in table 2. The loss fraction for counter-current injection are significantly higher.
Most of the beam ions resulting from local beamlet injection in the counter-current direction deposit in the low-field side loss region.In this region, a majority of trapped ions and some passing ions experience orbit losses and lead to charge separation.Specifically, when the injection position is at r = 0.3 m, the beam ions deposit just inside the last closed flux surface (LCFS), and nearly all of them undergo orbit losses.Although only a few beam ions from the co-current direction beam experience orbit losses on the low-field side, the magnetic drift of beam ions still result in charge separation and the formation of radial electric fields.
The toroidal distribution of wall heat load resulting from local beamlet injection is provided in appendix B. The highest thermal load on the wall caused by the counter-current local beamlet is 0.23 MW m −2 , which is still a relatively small thermal load on the first wall [31].Figure 12 illustrates the slowing-down density distribution of beam electrons and ions resulting from a beamlet injected at different radii.The solid line represents the distribution of beam ions, while the dashed line represents the distribution of beam electrons as well as the birth distribution of beam ions.The magnetic axis ρ = 0 corresponds to the major radius R = 1.91 m.The distribution reveals that tangential injection of the local beamlet leads to density profiles with peaks.Furthermore, the radial position of the beam ions' distribution is determined by the injection radius of the beamlet.
The electron distributions resulting from co-current or counter-current beamlets exhibit similarities.As the injection radius r increases, the deposition of beam ions decreases.In the case of co-current beamlets, the distribution of beam ions shifts inward compared to the beam electrons.However, in the case of counter-current beamlets, the distribution of beam ions shifts outward and then peak density decreases due to orbit loss and larger outer volume.

The radial electric field
The NEOE code utilizes equation (10) to calculate the radial distribution of the electric field resulting from local beamlet injection, and the calculated results are presented in figure 13.
The simulation results demonstrate that tangential injection of a local beam can establish a radial electric field with a sharp gradient.Co-current beamlet injection leads to a positive radial electric field, while counter-current beamlet injection results in a negative radial electric field.It is observed that when the beam energy is just 50 keV and the injected local beam area is 3 cm × 12 cm, the magnitude of the electric field can reach the order of 10 3 V m −1 with a steep gradient (represented by the curve on the left side of the peak), and the shear strength of the electric field is on the order of 10 4 V m −2 .
The injection position r of the beamlet corresponds to the peak position of the radial electric field, the shear position of the electric field can be controlled by the injection radius of the local beam.In the case of beamlet injection in the co-current direction, the peak positions of the electric field are consistently located to the left of the injection position, whereas for counter-current beamlets, the electric field peak values are consistently located to the right of the injection position.It was observed that the width of the electric field of the countercurrent beamlet was broader compared to that of the co-current beamlet.This was due to an additional bump caused by collisions, which enhanced ion orbit losses of counter-current beam ions near the plasma edge.
The strength of the electric field is influenced by the polarization effect.We calculate the polarizability using the magnetic configuration and density profile shown in figures 8 and 9, respectively, with parameters κ = 1.6 and q 95 = 6.51.The profile of polarizability in this case is provided in figure 14.We can see that the polarizability of the plasma χ = ϵ 1/2 is of the order of 10 5 .The polarizability near  the magnetic axis becomes significantly large as ε approaches 0. Near the plasma edge, the safety factor experiences a rapid increase owing to the presence of the X-point in the divertor magnetic configuration, leading to an increase in polarizability.
Although the local beamlet at the injection position r = 0.1 m exhibits a higher peak density, the background polarization effect near the magnetic axis is stronger, resulting in a weaker electric field.Similarly, at the plasma edge, the electric field formed by the counter-current beamlet with an injection position of r = 0.3 m is also shielded by the boundary polarization effect.
In order to enhance the strength of the electric field generated by local beamlet and achieve the desired values of 10 4 V m −1 for electric field strength and 10 5 V m −2 for electric field shear, adjustments can be made to the parameters of the beamlet.These adjustments include increasing the beam energy, density, cross-section, or adjusting the magnitude of the magnetic field.The individual effects of these parameters on the radial electric field are discussed below.

Impact of injection radius on electric field
The ion deposition resulting from the beamlet at different radial positions varies, consequently generating different radial electric fields.We performed a radius scan for beamlet injection under the same simulation parameters.The results are depicted in figure 15.The location of the electric field generated by the beamlet is controlled by the injection radius r: the smaller the injection radius, the closer the deposition of fast ions from the local beam is to the core region.
Under the condition of equal injection power, the beamlet induced electric field has a maximum value in the radial direction.Simulation results demonstrate that the strongest radial electric field is formed by the co-current beamlet with an injection radius r in the range of 0.2-0.25 m, reaching a magnitude of up to 8 kV m −1 .Similarly, for the counter-current beamlet, the strongest radial electric field is observed in the interval of 0.15 m to 0.2 m, with a maximum electric field magnitude of 6.4 kV m −1 .Those region corresponds to the intermediate region with weaker polarization effects, according to the polarizability profile provided in figure 14.

Impact of vertical positions of beamlet on electric field
We investigated the electric field resulting from local beamlet injection at various vertical positions (in the Z-direction).The beamlets are injected at heights of 0.234 m, 0 m, and −0.234 m, with an area of 3 cm × 12 cm and a beam energy of 50 keV.We focus on the local beamlet injection results for an injection radius of r = 0.2 m.
As shown in figure 16(a), the deposition position of the beam ions injected at the mid-plane is closer to the core compared to the local beam deposition positions at the upper and lower positions.The variation in deposition positions results in an inward movement of the radial electric field generated by the local beam injection at the mid-plane, accompanied by an increased width of the electric field, as shown in figure 16(b).
The results of counter-current beamlet injection exhibit similarities to those of co-current beamlets.These results indicate that under the condition of off-axis injection, although beamlets injected at different heights exhibit variations in the radial electric field positions, the electric field and its shear remain local and controlled.Even in cases where the injection angle of the beamlet is constrained by device structure, it is still possible to generate a localized and controllable radial electric field through off-axis injection.

Impact of beamlet cross-section size on electric field
The magnitude of the radial electric field generated by the local beamlet exhibits a clear correlation with the cross-sectional area of the beamlet, a larger cross-sectional area corresponds to greater NBI power (number of fast ions).Simulation results demonstrate that increasing the vertical dimension of the beamlet does not significantly affect the shear position of the radial electric field, and the magnitude of the electric field should increase proportionally with the height of the beamlet's cross-section within a certain range.However, the impact of increasing the width of the beamlet is more complex.
The electric field distribution generated by a local beamlet with a beam energy of 50 keV, a height of 0.12 m, and a width ranging from 0.1 m to 0.5 m is illustrated in figure 17.As the beamlet width expands, the proportion of fast ion deposition increases proportionally, accompanied by an increase in the peak value of the radial electric field.The shear rate of the electric field on the outer side exhibits a continuous growth with the widening of the beamlet, while the shear rate on the inner side remains relatively stable.However, the shear region progressively shifts inward as the local beamlet width expands.
Simultaneously, we observed that as the width increases, the magnitude of the electric field increment diminishes.This can be attributed to the widening distribution of deposited  beam ions, which may lead to partial shielding of charge separation by beam ions in other regions, as depicted in figure 18 (red area).The charge separation in this region does not contribute to the overall radial electric field.Consequently, it is preferable to employ a narrow beamlet to establish a localized ion distribution and electric field, with the beam width typically not exceeding the distance of the beam ions across the magnetic surface, which is approximately 14 cm for our field configuration.For an estimate, in the ITER case with a magnetic field of 5.3 T and q = 2.0 [32], the orbit deviation of an 80-400 keV deuterium ion is approximately 5-10 cm.Additionally, a specific beam width corresponds to a saturated beam energy, beyond which the electric field generated by the beamlet ceases to increase.We further examine this relationship in the subsequent section.

Impact of beamlet energy on electric field
The influence of beam energy on the generated radial electric field can be analyzed from the following aspects: (1) the neutralization efficiency of positive deuterium ions decreases continuously as the beam energy increases.When the ion energy exceeds 100 keV, the neutralization efficiency drops below 20%.(2) The increase in beam energy directly affects the ionization cross-section ν, consequently influencing the deposition position of the fast ions, i.e. the location of charge separation.(3) The cross-magnetic surface width of the beam ions is determined by the particle's pitch and energy, which also affects the distribution of the radial electric field.
In practical devices, a beamline consists of multiple ion sources, and the decrease in beam current caused by the reduction in neutralization efficiency can be compensated by injecting multiple sources simultaneously.Hence, we assume a constant beam current density (0.12 A cm −2 ) injected into the plasma and investigate the impact of the latter two factors on the radial electric field.These factors have already been includede in the neutral particle deposition and orbit calculations within the NEOE code.
The electric field resulting from local beamlet injection at different beam energies is shown in figure 19.The beam has an area of 3 cm × 12 cm, and is injected at a radius of r = 0.2 m.At a beam energy of 50 keV, the maximum electric field generated by the beamlet in the co-current direction is 5 kV m −1 .As the beam energy increases to 80 keV, the maximum electric field reaches 9 kV m −1 .While the peak value increases, the width of the electric field remains  Figure 20(a) illustrates the variation of the maximum electric field value with beam energy.We determine that the electric field has reached saturation when the increase in the electric field is less than 2%.The results in the figure indicate that the electric field gradually saturates as the beam energy increases.When the deposition width is held constant, there exists an energy threshold where the ions' orbits across the magnetic surface become sufficiently large, leading to a separation between the birth deposition area and the redistribution area of beam ions, as shown in figure 18. Continuing to increase the ion energy only widens the gap between these regions, without affecting the quantity of charge separation or the maximum electric field, which reaches saturation.However, due to the increased separation distance, the electric field shear weakens.
The beam energy corresponding to the local beam saturation electric field with different beam widths is shown in figure 20(b), while the height is kept at 12 cm.When the beam width is less than 8 cm, the wider the local beam, the higher the beam energy required to reach the saturation electric field.At the same time, it should be noted that when the beam width is greater than 8 cm, as the beam energy increases, the beam ion shine-through loss will also increase, and the beam ions deposited in the plasma will become less, and the electric field will no longer increase.
In the future, with the larger size and higher plasma density of tokamak devices, greater energy is required for local beamlet penetration into the plasma.This allows us to use larger  cross-sectional beamlets to enhance the strength of the radial electric field.However, careful parameter selection is necessary to avoid excessive high-energy ion losses to the first wall.

Impact of toroidal magnetic field on electric field
In future tokamak devices, the central magnetic field intensity can exceed 6 T. We performed a magnetic field scan ranging from 2 T to 6 T to investigate the electric field distribution generated by beamlet injection, while maintaining the qprofile depicted in figure 9(d).The beamlet injection was performed at the midplane position (r = 0.2 m) with a beam energy of 60 keV and a beam cross-sectional size of 3 cm × 12 cm.The scan results are illustrated in figures 21 and 22.
According to equation (3), the distance of beam ions acrossing the magnetic surface d, is inversely proportional to the magnetic field strength B. As the magnetic Figure 23.Variation of the boundary q 95 with plasma current (cross labels) and the maximum radial electric field with respect to plasma current (square labels).These results correspond to beamlet injection at the midplane position (r = 0.2 m, with a beam energy of 60 keV and a beam cross-sectional size of 3 cm × 12 cm.field strength increases, the deviation of the beam ions becomes smaller, resulting in reduced charge separation.However, on the other hand, the polarizability of the plasma ϵ 1/2 in equation (10), is inversely proportional to B 2 and this relationship is demonstrated in figure 21.Therefore, although the charge separation decreases, the plasma background polarization effect decays more rapidly, leading to an increasing trend of the electric field with the enhancement of the magnetic field.In magnetic field B 0 > 6 T, a radial electric field of 40 kV m −1 may be generated, corresponding to a beamlet power of 0.26 MW.
Simultaneously, the improved confinement due to the high magnetic field reduces the orbit losses of countercurrent beamlets from 49% to 13%.Consequently, the impact of the beam ions on wall heat flux is reduced.This implies that under future high magnetic field conditions, high-energy localized beams can be utilized to ensure charge separation and generate localized radial electric fields without causing significant effects on the wall.

Impact of plasma current on electric field
We conducted a parameter scan of the plasma current, which affects the generation of the radial electric field by influencing the safety factor.On one hand, it can be observed from the equation (10) that the neoclassical polarizability is proportional to q 2 .As q decreases, the shielding effect of background ions becomes weaker, resulting in a larger radial electric field.
On the other hand, according to equation (3), a smaller q results in smaller ion orbit deviations, which can potentially lead to a reduced radial electric field.These two mechanisms compete with each other.
We performed the scan with a beamlet injection at the midplane position (r = 0.2 m), beam energy of 60 keV, and beam cross-sectional size of 3 cm × 12 cm, while keeping the toroidal magnetic field (Bt) constant at 2.44 T. The plasma current was varied from 200 to 800 kA.
The scan results, as shown in figure 23, indicate that as the plasma current increases and the safety factor decreases, the maximum radial electric field generally exhibits an upward trend with increasing current.This suggests that, under our parameter conditions, the reduction in q leading to a decrease in polarization effects prevails, resulting in an increase in the electric field.However, in practical scenarios, maintaining such a low q profile is not feasible.It would be more appropriate to retain a safe q distribution while maximizing the toroidal magnetic field.

Summary and discussions
This paper presents the NEOE code, a numerical calculation code for the radial electric field generated by NBI.The code has been verified against the well-known NUBEAM code, and validated against experimental data from the EAST 2020 campaign.Based on this foundation, we propose a novel method for generating radial electric fields through localized beamlet injection.
The simulation results, based on the current level of neutral beam parameters, demonstrate that off-axis tangential local beamlet injection can generate a radial electric field in the localized region within the plasma.The electric field strength resulting from a 0.26 MW beamlet injection can reach 10 4 V m −1 , with a shear rate exceeding 10 5 kV m −2 .These values of radial electric field strength and electric field shear can match the the edge radial electric field levels of the H-mode observed experimentally.By controlling the angle and width of the local beamlet, we have the ability to manipulate the radial position and width of the electric field, enabling the formation and maintenance of a more significant electric field and electric field shear in a specific region.
In engineering applications, the beamlets may exhibit high losses and consequently low heating efficiency.Certain beamlet injections near the edge could result in significant shinethrough losses, necessitating avoidance or reinforcement of localized wall cooling capabilities.Additionally, the superconducting coils surrounding the vacuum chamber impose constraints on port accessibility.Therefore, achieving a large angle of tangential injection for the horizontal beam may be challenging.As an alternative, vertical injection from the bottom to the top of the device, specifically at the low field side, can still deposit beam ions into the plasma.Since the localized beam's required cross section is small, the necessary port size should not compromise the integrity of the magnet structure, the proposed approach demonstrates feasibility.
In future fusion devices, the low polarization effect under high magnetic fields enables the generation of a controllable radial electric field by a positive ion source beamlet.By increasing the cross-sectional area and beam energy of the beamlet, the strength of the radial electric field can be enhanced.
Simultaneously injecting multiple co-current or countercurrent beamlets at different toroidal positions allows us to modulate the shape of the electric field.The superposition of electric fields generated by these beamlets with different injection directions may lead to substantial changes in the electric field, potentially exceeding 10 6 V m −1 .If an electric field shear can be deliberately formed at the inner region of a large magnetic island, it becomes possible to suppress the growth of the magnetic island and enhance internal confinement.

Figure 1 .
Figure 1.Schematic diagram of tangential local neutral beamlet injection: (a) small cross-section beamlet, (b) tangential injection of local beam on the low-field side, (c) deposition of beam ions in the specific region.

Figure 2 .
Figure 2. Orbits of passing particles (a) and trapped particles (b) under magnetic drift, for particles with an energy of 50 keV and Bt = 2.4 T. The cross labels (red) indicate the birth positions of the particles.

Figure 3 .
Figure 3. Birth distribution of the beam ions from AL tangential co-current NBI calculated by NEOE (red) and NUBEAM (black) Codes.

Figure 5 .
Figure 5. (a) External heating powers, (b) plasma stored energy, (c) line averaged density and Dα signal of EAST shot #80757.(d), (e) Frequency spectra of the reflected signals measured by the correlation reflectometry at two different locations.

Figure 6 .
Figure 6.Vertical velocity profiles of turbulence measured by the reflectometer in the shot #80757.

4. 5 s
, and 5.5 s with the velocity at 2.8 s during ECRH-only.

Figure 7 .
Figure 7.Comparison of simulated radial electric field with experimental measurements: (a) co-current tangential NBI, (b) co-current perpendicular NBI, (c) combined co-perpendicular and co-tangential NBI.The corresponding positions are R ∼ 2.08 m with ρ ∼ 0.67 and R ∼ 2.02 m with ρ ∼ 0.5.

Figure 8 .
Figure 8.The magnetic equilibrium configuration (depicted in red) and the shape of the simulation wall (represented in black) utilized in the simulations.

Figure 9 .
Figure 9.The plasma profiles used in the simulations: (a) the electron density, (b) the electron temperature, (c) the ion temperature, (d) the safety factor q profile.ρ is the square root of the normalized toroidal magnetic flux.

Figure 10 .
Figure 10.Local beamlet layout, r represents the distance from the beamline to the magnetic axis.

Figure 11 .
Figure 11.Contours of the birth deposition of local beamlets in the plasma (blue region) and the distribution of wall losses (red region) vary with the injection radius r for co-current injection (a)-(c), considering a beam injection energy of 50 keV and a beam cross-section of 3 cm × 12 cm.

Figure 12 .
Figure 12.Distribution of ions (solid line) and electrons (dashed line, also representing the beam ion birth distribution) resulting from local beamlet injection: (a) co-current injection, (b) counter-current injection, with arrows indicating the injection radii r of the local beamlets.These simulations consider a beam injection energy of 50 keV and a beam cross-section of 3 cm × 12 cm.

Figure 13 .
Figure 13.The radial electric field generated by local beamlet injections: (a) co-current beamlets and (b) counter-current beamlets, with arrows indicating the injection radii r of the local beamlets.These simulations consider a beam injection energy of 50 keV and a beam cross-section of 3 cm × 12 cm.

Figure 14 .
Figure 14.The distribution of polarizability of the plasma: χ = ∑ j

Figure 15 .
Figure 15.The radial electric field distributions of beamlets at different injection radii r: (a) co-current beamlets and (b) counter-current beamlets, with a beam injection energy of 50 keV and a beam cross-section of 3 cm × 12 cm.

Figure 16 .
Figure 16.(a) The birth distribution of beam ions and (b) the radial electric field distribution resulting from co-current local beamlet injection at different height positions.The injection radius is r = 0.2 m, and a 50 keV beam with a cross-section of 3 cm × 12 cm is employed.

Figure 17 .
Figure 17.The electric field distributions generated by local beamlets of different widths: (a) co-current beamlets and (b) counter-current beamlets.The injection radius is r = 0.2 m, and a 50 keV beam with a height of 12 cm is utilized.

Figure 18 .
Figure 18.(a) The schematic diagram of the electric field distribution generated by local beams of different widths.(b) Illustrated the widening distribution of deposited beam ions lead to partial shielding of charge separation by beam ions.

Figure 19 .
Figure 19.The radial electric fields generated by local beam injections at different energies: (a) co-current beamlets and (b) counter-current beamlet.The midplane injection radius is r = 0.2 m, with a cross-section of 3 cm × 12 cm is employed.

Figure 20 .
Figure 20.(a) Variation of the maximum radial electric field generated by a local beam a cross-section of 3 cm × 12 cm as a function of beam energy, and (b) the corresponding beam energy required to reach the saturated electric field for different widths of local beams with a height of 12 cm.

Figure 22 .
Figure 22.Electric field distribution generated by beamlets at different toroidal magnetic field strengths: (a) co-current beamlets and (b) counter-current beamlet.The injection radius is r = 0.2 m, and a 60 keV beam with a cross-section of 3 cm × 12 cm is employed.

Figure B1 .
Figure B1.Variation of the heat load on the first wall with injection radius for co-current beamlets (a), (c), (e), and counter-current beamlets (b), (d), (f ).These simulations consider a beam injection energy of 50 keV and a beam cross-section of 3 cm × 12 cm.The boxes represent the guard limiters and the units of the heat loads are MW m −2 MW −1 .

Table 1 .
Neutral beam orbit losses and total torque.

Table 2 .
Loss fraction of co/counter current beamlets.