Numerical description of mode coupled waves in inhomogeneous magnetized plasmas

. Geometrical optics (GO) ray tracing has been widely used for a description of electron cyclotron resonance waves in inhomogeneous magnetized fusion plasmas. However, this reduced approach is not correct in su ﬃ cient low density plasmas with a sheared magnetic ﬁeld, where mode coupling between two electromagnetic-like cold plasma modes can occur. Here, we extend a ray tracing method based on the ﬁrst-order theory of extended geometrical optics (XGO), which captures mode coupled complex amplitude between O and X mode along the ray trajectory. In our approach, reference ray is calculated with ray equation to satisfy the lowest-order part of XGO theory and an evolution of complex amplitude proﬁle along the reference ray is calculated with partial di ﬀ erential equation derived from ﬁrst-order terms. Calculation results performed by extended ray tracing are in good agreement with 1D full wave analysis. By introducing second-order terms into our numerical approach, di ﬀ raction will be treated.


Introduction
Geometrical optics (GO) ray tracing has been used for a description of coherent wave propagation in inhomogeneous media, and therefore has been incorporated into a number of numerical codes for electron cyclotron resonance heating (ECRH) simulations in fusion plasmas. These codes treat two eigen plasma modes independently and can not manage the coupling between them. Such treatment is valid in sufficiently high density plasmas, where dispersion relations of two cold plasma modes are well separated. However, dispersion relations of two cold plasma modes are not sufficiently separated in low density plasmas, and therefore these modes will be easily mixed by sheared external magnetic field. Unfortunately, because this is not an unique condition in the peripheral region of fusion plasmas, a more general technique is needed for a more realistic EC wave description. Here, we developed a new ray tracing code based on the extended geometrical optics (XGO) theory [1][2][3][4]. Our extended code can calculate the O and X mode coupled vector complex amplitude along the ray. In this paper, a brief outline of a calculation scheme with 1st order XGO theory (Sec.2) and comparison between results performed by our approach and onedimensional full wave (1DFW) analysis [5,6] This basic vector equation is approximated wave equation up to first-order of GO parameter ϵ = λ/L ≪ 1, where λ is a wave length and L is a scale length of the inhomogeneity. The derivation and details of this theory are presented in [1]. Here, we present just an essence of our calculation. D(x, k), Ξ, E, and a are dispersion matrix, arrayed eigenvectors matrix of D H (x, k) for active modes, diagonal matrix of eigenvalues of D H for active modes, and vector complex amplitude for each active mode, correspondingly. Averaged eigenvalue E 0 is lowest-order term and is qualified as Ray-Hamiltonian. This is a very important technique for our approach. Therefore, details are mentioned below. Other terms are first-orders. Γ is dissipation but becomes zero in cold plasmas, v · ∇ + (∇ · v)/2 represents Poynting vector, and U gives polarizational effects including mode coupling. The index H and A denote the Hermitian and anti-Hermitian parts.
In the derivation of XGO, whole wave field Ψ is defined as where s is a number of eigenmodes, ξ s is eigenvector for eigenmode s, and θ(x) is rapidly oscillating phase function. Wave vector k is usually defined as normal vector of phase front. We also define wave vector as k ≡ ∇θ. Then, we introduced an unfamiliar definition, that is, "active." Roughly speaking, when appropriate wave vector k is given, the modes which can be excited in the medium are active modes. If branches of dispersion relations for each eigenmode are well separated, only one mode should be excited for given k. However, if branches are very close, two modes can be excited simultaneously. We can model the former case with XGO equation in one active limit, where complex amplitude a is treated as a scalar amplitude for only one active mode a. However, modeling of two active modes, which treats complex amplitude as a two-dimensional amplitude for convertible two active modes, is suitable for the latter case. More detailed accounts of modeling of two active modes are presented in [1,3].

Formulation for two active modes
Let us consider the O and X mode mixed waves, which can be modeled by our theory with the assumption that two eigenmodes are active. Such complicated waves can be excited in sufficiently low density plasmas where two eigenvalues of dispersion matrix D H are almost degenerated. (In vacuum, they are degenerated completely.) Also, interaction between these two cold plasma modes can be easily caused by a shear of a magnetic field, because dispersion relations of both modes are close to each other. The evolution of wave field in such a case is calculated as follows.
First, we use conventional ray equation for a conservation of the lowest-order part of XGO as, where r, k, E 0 , and τ are cartesian coordinate, wave vector, averaged active eigenvalue of D H , and small step along a ray trajectory, respectively.
The way to choose initial E 0 is as follows, that provisional wave vectors k 1 which satisfies E 1 = 0 and k 2 for E 2 = 0 are calculated by using well known and corresponding dispersion relations [7]. Then, averaged wavevector k = (k 1 + k 2 )/2 is calculated and by substituting k into Eq. [4], E 0 andẼ are given. Magnitudes of E 1 and E 2 can be very different andẼ, which is deviation from the averaged part of E, can approach the lowest-order magnitude in high density plasmas. In such cases, the assumptionẼ = O(ϵ) is violated. Thus, this model works well in sufficiently low density plasmas, where E 1 and E 2 are almost degenerated. However, when E 1 and E 2 are very different, its branches of dispersion relation are also very different. It means that two active modes cannot excited simultaneously for a given k. Therefore, in high density plasmas, it seems difficult to emerge mode coupling, and you should model such cases with XGO in one active limit. The lowest-order system gives a reference ray trajectory. Also, first-order terms of XGO equation includingẼ, which expresses the deviation from a reference ray, give mode coupled wave field structure along a ray trajectory.
Here, it is convenient to introduce the curvilinear coordinates (ζ, ς) for a calculation of complex amplitudes around the reference ray, where ζ is the longitudinal and ς is the transverse coordinate as, .
Then, by mapping into ray based curvilinear coordinates and a few simplification, v · ∇ in original coordinates is rewritten as v · ∇ ≈ ±∂ τ + (ṽ α ∂ς µ ∂x α ) ∂ ς µ = ±∂ τ + δV µ ∂ ς µ . (6) Also, since we consider waves in cold magnetized plasma where D A is vanishing, dissipation term Γ becomes zero. Eventually, partial differential equation for mode coupled vector complex amplitudes is obtained from higher-order terms of Eq. [1] as, ) a + δV µ ∂a ∂ς µ . (7) The lowest-order balance with Eq. [3] gives a reference ray, and first-order balance with Eq. [7] gives an envelope. This model is valid when the magnitude ofẼ stays in higher-order. In other words, this model works well in sufficiently low density plasmas. (e.g., a peripheral region of fusion plasmas.) We presented equations for a numerical approach up to the first order theory of XGO, which does not include diffraction. (Analytical approach is presented in [1].) However, we can introduce diffraction by using second-order theory. Actually, our numerical approach facilitates to include the diffraction. A more general method which captures mode coupling and diffraction simultaneously will be discussed elsewhere.

Numerical result
In this section, some numerical results simulated by a newly developed ray tracing code based on a model that treats two active modes, which is described in Sec. 2, and by 1DFW analysis code are presented.    Test wave is excited at z = 0.0 [m] and propagates into slab plasma along z axis. Here, wave frequency is 77.0 [GHz] and polarization state of vector complex amplitude is defined initially as follows, where (α, β) = (10.0, −30.0) (deg) in these cases. In ray tracing code, excited vector wave is mapped into eigenmode representation with a = Ξ † ψ and integrated by paraxial equation along a ray. 1DFW code calculates an evolution of vector wave field ψ directly. Profiles of electron density and external magnetic field, depending only on z direction, are n e (z) = n e0 · exp  Fig. 1 -4. Fig. 1 -4 show the evolutions of wave parameters (amplitude, polarization states, and mode contents) of an electromagnetic-like O and X mode coupled wave along its trajectory in a cold slab plasma. Fig. 1 and 3 are results of extended ray tracing code, and Fig. 2 and 4 are results of 1DFW code [5,6]. Fig. 1 and 2 are calculated for the weakly sheared magnetic field (magnetic field shear period τ s = 1.0 × 10 3 ), and Fig. 3 and 4 are calculated for the strong one (τ s = 2.5 × 10 −1 ).
Upper hybrid resonance, right hand cutoff, cyclotron resonance, and plasma cutoff frequencies normalized by test wave frequency (a), each component of the external magnetic field unit vector (b), each component of the wave field intensity (c), local α and β of mode coupled wave (d) and power fraction of two modes (e) as a function of ray path length from excited point are showen in Fig. 1  -4, respectively. (a) shows that waves do not cross any resonance or cutoff layer. Local α and β of the combined modes are plotted in (d), which are given by where γ = tan −1 (|ψ y |/|ψ x |) and δ = arg ψ y − arg ψ x . Power fraction in (e) is given by where q denotes the each cold plasma mode and (α q , β q ), which is polarization state angle for each mode, can be obtained by applying Eq. [11,12] for each mode field ψ q . Here, η 1 and η 2 in Fig. 1 and 3 correspond to η x and η o in Fig. 2 and 4, respectively. Initial mode fraction is η 1,x : η 2,o ∼ 0.2 : 0.8 in Fig. 1 -4. However, The fraction of two modes are reversed in Fig. 3 and 4. Comparison between the calculation for a weakly sheared magnetic field ( Fig. 1 and 2) and for a strong one ( Fig. 3 and 4) demonstrates that mode coupling is caused by a sheared magnetic field, and the coupling depends on the shear. Then, good agreement between extended ray tracing ( Fig. 1 and  3) and 1DFW analysis ( Fig. 2 and 4) demonstrates that extended ray tracing captures mode coupling very well.
Here, our new approach does not reduce refraction compared to 1DFW, and conventional GO ray tracing cannot calculate mode coupled vector waves in anisotropic media. Therefore, extended ray tracing has a great advantage for simulations of waves in peripheral regions of fusion plasmas.

Summary
Extended ray tracing based on first-order XGO theory is a prospective technique for describing the electromagneticlike vector waves in sufficient low density plasmas with a sheared magnetic field. Such a situation is realized in a peripheral region of fusion plasmas. Mode coupling estimation performed by extended ray tracing shows good agreement with the result of 1DFW analysis. Diffraction will be treated, by introducing second-order terms of XGO into our numerical approach. A more general technique which captures refraction, mode coupling, and diffraction simultaneously will be discussed in forthcoming papers [8].