{'fraction_non_alphanumeric': 0.034843661419039725, 'fraction_numerical': 0.01255941738258253, 'mean_word_length': 4.602736379183972, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Variational Symplectic Particle-in-cell Simulation of Nonlinear Mode Conversion from Extraordinary waves to Bernstein Waves

In this paper, the nonlinear mode conversion of extraordinary waves in nonuniform magnetized plasmas is studied using the variational symplectic particle-in-cell simulation. The accuracy of the nonlinear simulation is guaranteed by the long-term accuracy and conservativeness of the symplectic algorithm. The spectra of the electromagnetic wave, the evolution of the wave reflectivity, the energy deposition profile, and the parameter-dependent properties of radio-frequency waves during the nonlinear mode conversion are investigated. It is illustrated that nonlinear effects significantly modify the physics of the radio-frequency injection in magnetized plasmas. The evolutions of the radio-frequency wave reflectivity and the energy deposition are observed, as well as the self-interaction of the Bernstein waves and mode excitations. Even for waves with small magnitude, nonlinear effects can also become important after continuous wave injections, which are common in the realistic radio-frequency wave heating and current drive experiments.


I. INTRODUCTION
Radio-frequency (rf) wave injection plays an important role in auxiliary heating and current drive for modern magnetic confinement fusion devices [1][2][3][4]. A typical method for rf wave injection into overdense plasmas, where the plasma frequency ω pe is much larger than the electron cyclotron frequency ω ce , is to excite the electron Bernstein waves (EBWs). Two different approaches are commonly used for coupling rf power to EBWs, i.e., the ordinaryextraordinary-Bernstein (O-X-B) and extraordinary-Bernstein (X-B) mode conversion. The O-X-B scheme is usually used for the electron Bernstein wave current drive (EBCD) in spherical tokamaks (STs), where the electron cyclotron resonance current drive (ECCD) method is ineffective [4,5]. The X-B scheme requires steeper plasma density profile, which is also available in some STs [6][7][8][9]. During the X-B mode conversion, the electromagnetic wave is initially injected into the plasma perpendicularly to the background magnetic field, and propagates as an extraordinary (X) wave within the plasma. When the X wave reaches the cutoff location, where the local right-handed cutoff frequency ω R is equal to the frequency of the incident wave ω, part of the energy of the incident wave penetrates the cutoff layer and continues propagating as EBWs. The X-B conversion will be significantly enhanced if the electron density gradient becomes large enough.
The linear X-B mode conversion can be theoretically investigated using the cold plasma model [10,11]. In this model, the kinetic effects of the electrons are neglected. At the same time, the ions offer a positive charged background. The brief description of the linear theory is presented in Appendix A. The linear theory predicts that the reflectivity of the injected wave is a constant that does not change with time.
The cold plasma assumption is only accurate when the wave length is much larger than the electron gyroradius. As the electron temperature increases, the kinetic effect becomes significant and should be taken into account. The kinetic theory of mode conversion is so complicated that the analytical solution is unaccessible except for some special cases [6]. In addition, the linear condition no longer holds for large-amplitude or long-term wave evolution. In modern magnetic confinement devices, the power of the rf wave injection is so large that the nonlinear effects can not be neglected. The nonlinear effects during the propagation of Bernstein waves have been directly observed in some experiments [4,12,13]. Those nonlinear effects can modify the reflectivity during the wave injection, which cannot be treated by the linear theory. The simulation study on the nonlinear physics of rf waves in plasmas has become an active research field [14][15][16][17]. However, the analytical and numerical studies on nonlinear rf physics are severely constrained by existing mathematical and algorithmic tools.
Containing both kinetic and nonlinear effects, first-principle particle simulation is regarded as a necessary and powerful tool for dealing with the nonlinear rf physics. Different aspects of the X-B mode conversion have lately been investigated using particle-in-cell (PIC) codes [14][15][16][17]. For nonlinear mode conversions, where the timescale of the mode conversion is much larger than the periods of the rf waves, long-term simulations are inevitably required.
Since the number of time steps in a long-term simulation is large, the coherent accumulation of the numerical error from each time step turns out to be the most serious difficulty. To ensure the reliability of the simulation results, long-term numerical accuracy and stability should be guaranteed. To overcome the difficulty, a newly developed variational symplectic PIC method [18] is employed in this paper to study the nonlinear X-B mode conversion. The symplectic PIC algorithm bounds the global numerical error for an arbitrary long time by preserving the discrete symplectic structure of the original Lagrangian system. Its long-term accuracy, stability, and conservativeness have been theoretically discussed and practically verified [18,19].
Suffering no accumulation of the global errors, the symplectic PIC method shows its unparalleled advantage in the long-term simulations and provides accurate nonlinear results.
In the simulation of nonlinear rf mode conversion, it is illustrated that the nonlinear effects significantly modify the physics of the rf wave injection in magnetized plasmas. The wave reflectivity, which is regarded as a constant in linear theory, varies in the nonlinear mode conversion process. The evolutions of the rf wave reflectivity and the energy deposition are observed, as well as the self-interaction of the Bernstein waves and mode excitations.
The spectra of the electromagnetic wave, the energy deposition profile, and the parameterdependent properties of radio-frequency waves during the nonlinear mode conversion are also investigated. It is found that even for waves weak in magnitude, nonlinear effects can become significant after continuous wave injections. When simulating rf waves with small amplitudes, nonlinear phenomena become prominent after several hundred wave periods. New modes are excited over a sufficient long time even if the amplitude of incident wave is small. To study the parameter dependency of the nonlinear X-B conversion, we compare the rf wave injection processes with different electron temperatures, amplitudes of injection waves, and plasma density gradients. It is discovered that a lower electron temperature can enhance the reflectivity, and the reflectivity depends on the plasma density gradient in a sensitive way.
The paper is organized as follows. In Sec. II, we introduce the physical model and the problem of the nonlinear X-B mode conversion. The numerical method used in the longterm nonlinear simulation, i.e., the variational symplectic PIC algorithm, is also introduced.
In Sec. III, the symplectic PIC simulation of the nonlinear X-B conversion is carried out.
Different factors, which effect the nonlinear rf wave processes, are analyzed in detail. Section IV is a the brief summary and future research plan.

II. THE MODEL OF X-B CONVERSION AND THE SYMPLECTIC PIC METHOD
In the X-B conversion process, the incident wave is perpendicular to the background magnetic field. This problem can thus be described by a planar model. Consider a planar plasma with the density nonuniformity in the x-direction and a uniform background magnetic field along the z-direction. The wave source is located outside the plasma. The source of the electromagnetic wave, with the electric field component along the y-direction, excites an X wave in the plasma boundary. The X wave propagates in the plasma until it reaches a cut-off location, and part of the wave converts to the Bernstein wave. There are two main restrictions in the PIC simulation of the X-B mode conversion. Firstly, the wave length of the Bernstein wave is much shorter than the electromagnetic wave with the same frequency in the vacuum. So the grid length in the simulation should be small enough to resolve it, which means the time step is restricted to a small value correspondingly. A large time step will lead to numerical instability because the particles can pass through multiple grids within a single time step. Another issue is that the group velocity of the Bernstein wave is much slower than the speed of light in the vacuum, so the mode conversion process lasts for a long time. It takes hundreds of periods before the wave reaches the core region of plasma. Then at least ten thousand time steps are required to represent a full X-B mode conversion process. Traditional numerical methods may fail to obtain useful results due to the accumulation of numerical errors. In order to study the mode conversion process with more accuracy, a variational symplectic method, which has excellent long-term conservation properties [19], is adopted in the particle-in-cell simulation.
The basic idea of constructing a variational symplectic scheme is to discretize the Lagrangian and perform the discrete variational principle [20]. The starting point is to consider the self-consistent system consisting of charged particles and electromagnetic fields. The Lagrangian is the integral where is the Lagrangian density. Here m s , q s , x s , andẋ s denote the mass, electric charge, velocity, and position of the s-th particle respectively, ǫ 0 and µ 0 are the permittivity and permeability of free space, and A and φ are vector potential and scalar potential of the electromagnetic fields. The time derivative of vector potential is denoted asȦ.
We follow the method used in Ref. [18] to discretize the Lagrangian L presented in Eq. (1) in cubic meshes and perform the discrete variational principle to obtain the numeric scheme.
The discrete Lagrangian is where A Here ∆x and ∆t are the grid size and time step respectively. The subscript l stands for variables in the lth time interval, and i, j, k are integer indices for three directions of the Cartesian coordinate frame. The electric potential φ is set to be zero because of the adoption of temporal gauge. The corresponding discrete action integral A d is and the motion equations are obtained from the discrete variational principle The discretization of Lagrangian and the discrete variational principle lead to a PIC iteration scheme. Theoretically, the variational symplectic PIC scheme automatically conserves a discrete symplectic two-form [20,21]. For any discrete Lagrangian L d (q l , q l+1 ), where q l denotes the generalized coordinates discretized at the lth time step, there exists an discrete symplectic two-form defined as where the operator d denotes the exterior derivative. The discrete Lagrangian of Vlasov-Maxwell system is defined in Eq. (??) with q l = (A i,j,k,l , x s,l ). The two-form in Eq. (15) is a geometric structure, called the discrete symplectic structure. Because the original selfconsistent Vlasov-Maxwell system is Hamiltonian, it is naturally equipped with a continuous symplectic structure, which is conserved exactly during the time evolution. The conserved discrete symplectic structure in the variational symplectic PIC method is an approximation of the original continuous geometric structure. The good long-term properties of the symplectic algorithm are guaranteed through the conservation of the discrete symplectic two-form, that is Ω l−1,l = Ω l,l+1 for any l, during time advance in the numerical simulation, which can be proved exactly [20]. Therefore, the numerical method has long-term conservativeness and accuracy, which are necessary when simulating long-term nonlinear rf wave problems and have been verified numerically [18].
To decrease the numerical noise caused by the lack of sampling points, the shape function is applied. In this paper it is chosen as where and The factor a is set to be 0.25∆x.
With the variational symplectic PIC method, we are able to simulate the nonlinear process of X-B mode conversion. The plasma is inhomogeneous in x direction, and as an example we assume the electron density profile is where n r ∆x is the thickness of the plasma boundary, ∆x = 2.773 × 10 −5 m, and n 0 = 2.0 × The density profile is plotted in Fig. 1. The physical domain is 0 < x < 1300∆x.
At x = 0 we use the Mur's absorbing boundary condition (MABC) [22] to eliminate the unwanted wave reflections. However at x = 1300∆x the density of plasma is not zero and the MABC fails. To construct an absorbing boundary, we apply an auxiliary calculation domain 1300∆x < x < 1792∆x with the decreasing electron density During the simulation, the Bernstein wave propagates in the plasma but not fast enough where E 1 is the amplitude of the oscillating electric field, and ω is the frequency of incident wave. The time step is set to ∆t = ∆x/c ≈ 9.24 × 10 −14 s. According parameters, the corresponding electron cyclotron frequency is ω ce = 9.67 × 10 10 rad/s, and the plasma frequency is ω pe 0 = 25.2 × 10 10 rad/s in the core plasma. The frequency of the incidence wave ω is chosen to be ω = 0.0145/∆t ≈ 1.62ω ce . According to the linear mode conversion  forms in a short time, and the reflection is stable (Fig. 2). After about 100 periods, the wave structure to the right of UHR becomes unstable, as marked by the dashed ellipse in Fig. 3.
The wave structure varies with time irregularly, which is an implication of the nonlinear effect. There is no direct evidence of nonlinear generation of other frequencies in Fig. 3.
Further results of nonlinear generation of EBWs with other frequencies will be discussed in Sec. III B.

B. New mode generation during X-B mode conversion
The primary signature of nonlinear processes is the generation of additional modes with frequencies other than the incident frequency. In this subsection, the new mode generation The generation of the mode at frequency 2ω can be explained as the self-interaction of the initial Bernstein wave. This second harmonic generation (SHG) of Bernstein wave has been studied extensively [12,17,25,26]. Then the wave with frequency ω ′ = 0.77ω is generated near the location where the SHG of Bernstein wave happens, which implies that the generation of the wave with frequency ω ′ is closely related to the SHG of Bernstein wave.
According to the theory of parametric decay instability (PDI), there should be a nother wave with frequency ω ′′ , which satisfies either 2ω + ω ′ = ω ′′ or 2ω − ω ′ = ω ′′ . To verify this PDI process, we plot in Fig. 5 the dispersion relation of electromagnetic EBWs in hot plasma, where the electron density of the plasma is n e (608∆x), and the electron temperature is in the X-B mode conversion simulation at x = 608∆x where the wave with frequency ω ′ occurs. We can see that these modes satisfy the relation k ′′ + k ′ = 2k and ω ′′ + ω ′ = 2ω.
In Fig. 5, it is evident that the frequency and wave number of the third mode satisfy ω ′′ = 2ω − ω ′ and k ′′ = 2k − k ′ , respectively. This mode can be visualized through the spectra of E y plotted in Fig. 6, where the mode at frequency ω ′′ = 1.23ω manifest itself clearly.
The evolution of E x spectrum in plasmas with higher electron temperature (123.9 eV) or weaker incident wave (E 1 =0.25 MV/m) are also investigated.
When temperature is higher, the location where the wave with frequency ω ′ = 0.77ω appears shifts to x ≈ 680∆x, which is much deeper in the plasma. This means that the wave has a greater tendency to penetrate the plasma with higher temperature. At the same time, the reflectivity is lower.
If we use a smaller source amplitude E 1 , the waves at frequencies ω ′ and ω ′′ generated during the mode conversion are much weaker than the wave with frequency ω. This is of course due to the nonlinear nature of the PDI. The scenario is the same as Fig. 4.

C. Dependence of reflectivity on wave amplitude
In wave heating and current drive the energy conversion efficiency C is a key index. Given the energy reflectivity R, the conversion efficiency can be calculated according to C = 1 − R.
During the nonlinear X-B mode conversion, the modes with different frequencies coexist.
The reflectivity defined for linear case, as in Ref. [14], is no longer valid. In this work, we calculate R based on the Poynting vectors (see Appendix B for detail).
We now consider the reflectivity of two incident electromagnetic waves with different amplitudes, i.e.,

D. The energy deposition profile
Since the location of wave energy deposition is crucial to rf heating, it is important to study the energy deposition profile. Generally speaking, the energy depositing inside the plasma core region will be helpful for heating and current driving. In this subsection, the energy deposition is obtained by comparing the kinetic energy distribution of particles at different times. In plasmas, the wave energy is stored in both electromagnetic fields and oscillating particles during each period of the wave. To calculate the energy deposition rate, we average the energies of particles over 4000 time steps and record the energy evolution in The result in Fig. 8(b), which shows no evident electron heating, is consistent with the prediction of the linear theory. According to the linear theory, in the cold plasma approximation the wave energy should accumulate at x = x U HR exactly. But in hot plasma, it is expected that the wave energy will be converted into the Bernstein mode and be carried into On the other hand, for waves with large amplitude in Fig. 8(a), obvious energy deposition occurs to the right of the UHR position, which is caused by the nonlinear effect. Simulation results show that a large amount of wave energy is accumulated at x ∼ 608∆x, where the SHG of EBW appears. From Fig. 4(c) and Fig. 6(a) it can be observed that the waves with frequencies ω ′ = 0.77ω and ω ′′ = 1.23ω are generated at x ∼ 608∆x. Notice that ω ′′ = 1.23ω ∼ 2ω ce (see Fig. 5), and the wave with frequency ω ′′ will lead to electron heating through the second-harmonic resonance. The complete process of the nonlinear energy accumulation for large amplitude incident waves is provided as follows.
1. The X wave converts to Bernstein wave with the same frequency ω. Because the nonlinear effect strongly depends on the parameters of the plasma and incident wave, different nonlinear deposition processes will appear in different situations. We will investigate the large variety of nonlinear energy deposition phenomena in future studies.

E. Reflectivity influenced by electron temperature
To find out the dependence of the energy conversion efficiency on the background electron temperature, waves propagating in plasmas with different electron temperatures are investigated. The reflectivity evolution in plasmas with electron temperatures T e = 123.9 eV and T e = 57.6 eV is plotted in Fig. 9. We can observe that a long time (1200/ω) after the wave is injected, the energy reflectivity in the lower T e is about twice as large as the higher T e . This suggests that in wave heating experiments a higher electron temperature should be kept in order to maintain a larger mode conversion efficiency.
To obtain more information of the relation between electron temperature and reflectivity, the case of T e ≈0 eV is also tested. Two different wave amplitudes E 1 = 0.066kV/m and wave is injected, E y (x) of the wave with small amplitude agrees well with the linear theory.
On the other hand, Fig. 11 implies that even for waves with small amplitude, the mode conversion structure is not stable and the reflectivity varies significantly after a long time (t ≈ 400/ω). Figure 11 also suggests that X-B mode conversion in cold plasmas for wave with large amplitude (E 1 =1.0 MV/m) is unstable, and the reflectivity increases quickly after t > 400/ω. The variation of reflectivity is very sensitive to the value of density gradient. Therefore, when applying the X-B rf heating in the magnetic confinement devices, it is preferable to maintain a large edge electron density gradient in order to maximize the conversion efficiency.

IV. SUMMARY AND DISCUSSION
The variational symplectic PIC scheme has been applied to simulate the nonlinear X-B mode conversion. The preservation of symplectic structure guarantees the accuracy of the long-term simulation results. Simulation results show that when the amplitude of incident wave is sufficiently large, the reflectivity will increase during the mode conversion, and modes with frequencies other than ω are exited. Even for waves with small amplitudes, the reflectivity changes significantly after a long time. In the meantime, we find that the energy deposits within a wide spatial range in the plasma. The simulation study also suggests that a higher electron temperature can improve the efficiency of energy conversion. Moreover, a In this paper, the X-B mode conversion is discussed, and a series of interesting and useful nonlinear phenomena are investigated. The conservation property of the variational symplectic PIC plays a crucial role in achieving the new results during the long-term simulations.
In future studies, the long-term simulations of mode conversion of other waves, such as lower hybrid wave (LHW) and ion cyclotron wave (ICW), will be carried out. We expect that, the variational symplectic PIC simulation is effective for these nonlinear rf processes and is able to suggest optimal parameters for rf heating and current drive schemes.
In fluid theory, the evolution of plasma is determined by the following equations, ∂ ∂t n e = −∇ · (n e u) , 1 c 1 c For the propagation of X-waves, we only need to consider the extraordinary component (E 1 · B 0 = 0). For 1D inhomogeneous, field variable can be expressed as After linearizing Eqs. (A1)-(A4), we obtain a second-order ODE where Here n e (x) is decided by Eq. (18), n r = 380, and the external magnetic field is B 0 = 0.55 Te z .
Note that this ODE may have a singular point ξ U HR where According to previous results of Refs. [10] and [11], the connection condition at ξ = ξ U HR is E y | ξ U HR −0 = E y | ξ U HR +0 = E y (ξ U HR ) , ate ω for the incident wave, we calculated the reflectivity for ω ∈ [1.2ω ce , 2.0ω ce ] and plotted in Fig. 14, which shows that the lowest reflectivity is at ω ≈ 1.37ω ce . For ω ∈ [1.2ω ce , 2.0ω ce ], the reflectivity will not be too high to make the mode conversion insignificant. For all our examples, ω is set to 1.62ω ce , and the corresponding reflectivity is R ≈ 0.15 when n r = 380.
Therefore, the energy reflectivity is All Poynting vectors here are time averaged in 20π/ω, which is significantly larger than the period of waves under investigation.