Plasma and Fusion Research: Regular Articles

Volume 14, 1401158 (2019)

Finite Larmor Radius Eﬀect on Ion-Temperature-GradientInstability in Cylindrical Plasmas
Naohiro KASUYA1,2) , Tomotsugu OHNO2) , Makoto SASAKI1,2) and Masatoshi YAGI3)
1)
2)

Research Institute for Applied Mechanics, Kyushu University, 6-1 Kasuga-Koen, Kasuga, Fukuoka 816-8580, Japan
Interdisciplinary Graduate School of Engineering Sciences, Kyushu University, 6-1 Kasuga-Koen, Kasuga, Fukuoka
816-8580, Japan
3)
National Institutes for Quantum and Radiological Science and Technology, 2-166 Omotedate, Obuchi,Rokkasho-mura, Aomori 039-3212, Japan(Received 6 June 2019 / Accepted 19 August 2019)
Numerical analyses using a gyro-fluid model have been performed to investigate the finite-Larmor-radius(FLR) eﬀect on ion-temperature-gradient (ITG) instability in cylindrical plasmas. A spectrum code with FourierBessel expansion has been developed for the analysis of global mode structures. The analytical expression ofthe η (ratio between the density and temperature gradient lengths) threshold value for linear ITG instability has
been obtained from the local dispersion relation, whose dependency on the ion temperature comes from the FLR
eﬀect. Dependency of the threshold by the global analysis is reproduced by the local analysis with appropriate
selection of the perpendicular wavenumber. Break of the Boltzmann relation by the FLR eﬀect is not strong as togenerate another unstable branch of the ITG mode.
c 2019 The Japan Society of Plasma Science and Nuclear Fusion Research

Keywords: gyro-fluid model, finite-Larmor-radius eﬀect, linear device, ion-temperature-gradient instability, dispersion relation
DOI: 10.1585/pfr.14.1401158

1. IntroductionNumerical simulations of micro-scale turbulence havebeen carried out to identify the cause for transport in magnetic confined fusion plasmas. Several gyro-kinetic simulation codes have been developed for that purpose, and local linear stability analyses for several kinds of modes, formation mechanism of global turbulent structures and theirdynamics can be studied [1]. On the other hand, gyro-fluid
models have been also used for transport analysis. Computational costs are reduced compared with gyro-kinetic simulations, so are advantage in carrying out parameter scansand global simulations of turbulence. Transport models
GLF23 and TGLF are developed based on the gyro-fluid
model [2–4]. Three-dimensional turbulence simulations
are carried out for toroidal ion-temperature-gradient (ITG)
turbulence [5], and dynamical phenomena, as ELM burst,
are also studied [6]. Spatio-temporal causal relations included in a wide range of plasma regions is one of the keys
for transport mechanism, so a simulation scheme eﬃcientfor global analysis is needed.
Clarifying fundamental mechanisms of turbulence isimportant [7], and researches on structural formation ofplasma turbulence have been progressing by using linear
devices [8]. A linear device has a cylindrical shape without magnetic curvature and can produce high density plasmas, so is suitable for basic experiments to observe plasma
instabilities in details [9,10]. Numerical simulations of lin-ear devices have been also carried out, and a reduced fluid
model for drift wave turbulence is used to describe the3-D turbulence [11]. For understanding structural formation mechanism, identification of the spatial structure ofeigenmodes is important as drive of turbulence. In basic
experimental devices, the typical spatial scales of the device and gyro-motion of ions are close to each other, so
the finite Larmor radius (FLR) eﬀect is not negligible. Amode with the wave length comparable to the eﬀective ion
Larmor radius can become unstable even in low ion temperature plasmas [12]. For the first step, linear characteristic of the ITG instability has been investigated using theparameter set of PANTA linear device [13]. Here we havedeveloped a global simulation code for the ITG instability
in linear devices [14]. The parameter dependency is evaluated with a gyro-fluid model to identify influence of theFLR eﬀect, and is compared with analytical expression obtained from a linear dispersion relation in this article.

2. ModelFluid models are used for analyzing global mode
structures in collisional plasmas. The following set of linearized gyro-fluid equations [15] is used, where the magnetic curvature terms can be eliminated for the analysis of
the cylindrical configuration;author’s e-mail: kasuya@riam.kyushu-u.ac.jp1401158-1⎛⎞∇ˆ 2⊥ ⎟⎟⎟ 1 ∂Ψdn⎜⎜⎜⎜⎟⎠+ ∇// u// + ⎝1 + η⊥
= 0,dt2 Ln ∂y(1)c 2019 The Japan Society of Plasma
Science and Nuclear Fusion Research

Plasma and Fusion Research: Regular Articles
du//+ ∇// τn + T // + Ψ = 0,dt
q//
1 ∂Ψ
1 dT //+ ∇// 2u// ++ η//τ dt
τLn ∂y2νii =−T // − T ⊥ ,3τ⎡⎤⎥⎥ 1 ∂Ψq⊥ ⎢⎢⎢ ∇ˆ 2⊥
1 dT ⊥ˆ2ˆ+ ∇//+ ⎣⎢+ η⊥ 1 + ∇ ⎥⎦⎥τ dtτ
2Ln ∂yνii =T // − T ⊥ ,3τ

Volume 14, 1401158 (2019)
(2)

(3)

(4)

where n is the ion density, u// is the ion velocity, T is theion temperature, τ is the ratio between ion and electron
temperature at the plasma center, η is the ratio between
the density and temperature gradient lengths, and the subscripts // and ⊥ represent the quantities in the parallel andperpendicular directions to the magnetic field, respectively.
The other definitions of the parameters are described inRef. [14]. The gyro-averaged potential Ψ is represented
as Ψ ≡ Γ0 1/2 φ with operators Γ0 1/2 = (1 + bτ/2)−1 andb = −∇2⊥ . Operator b gives the square of the perpendicular
wavenumber k⊥2 , which corresponds to the magnitude ofthe FLR eﬀect. Two modified Laplacian operators ∇ˆ 2⊥ and∇ˆˆ 2 are introduced to include the FLR eﬀects. The quasineutrality relation is given to beΨ
b/2T ⊥ − (1 − Γ0 ) = Ψ(5)
Γ0 n −
1 + bτ/2
τequations is the linearized one, but the radial components
with diﬀerent js are not independent to each other, due tothe (1 / r) terms, as in ∂Ψ /∂y = (1/r) ∂Ψ /∂θ.
For simulations, experimental parameters in
PANTA [10] are used; device length Lz = 4.0 m,
plasma radius a = 0.07 m, density n0 = 1.0 × 1019 m−3 ,Ln = 0.07 m, νii = 350 s−1 , magnetic field B = 0.1 T, temperatures T e0 = 3 eV and T i0 = 0.3 eV. The temperaturesand magnetic field give ρs = 1.1 cm, ρi = 3.5 mm andΩci /2π = 3.8 × 104 Hz for argon plasmas. Here Ln (andη) is set to be constant in space. We concentrate on themode structure restricted by the plasma radial boundaries,
so this simplification is used.

3. Finite-Larmor-Radius Eﬀect
To obtain an analytical form of dependency of ITG
instability, the local model with linearization d/dt → −iω,∇// → ikz , ∂/∂y → ikθ gives the following dispersion relation from Eqs. (1 - 5);
ω3 − ω1 ω2 + P1 ω + P2 = 0,
whereto determine the relation between the density and potential.
The Fourier-Bessel expansion is used to calculate 3-Deigenmode structures. The Bessel expansion in the r direction and Fourier expansions in the θ, z directions areapplied on Eqs. (1) - (5) by usingf (r, θ, z) =M
LJ 

fm jl Jm j (r) eimθm=−M j=1 l=1πz(2l − 1)
,2LzM J L
u//m jl Jm j (r) eimθu// (r, θ, z) =× cos(6)m=−M j=1 l=1

× sinπz(2l − 1),2Lz(7)

where f implies {n, T // , T ⊥ , φ}. Bessel functions
Jm j (r) ≡ Jm (λm j r/a)(8)

are used for the expansion to satisfy the radial boundary
conditions {φ, n, u, T // , T ⊥ } = 0 at r/a = 0, 1 for perturbation components, where Jm is the m-th Bessel function,
λm j is the j-th point satisfying Jm (λm j ) = 0, m and l are
the azimuthal and axial mode numbers, respectively, a isthe plasma radius and Lz is the device length. Fourier expansions are used in the azimuthal and axial directions for
simplicity, in spite of existence of the axial boundary determined at the position of the end plate. The set of model3ω∗2τ
− F (1 + ηb1 ),ω1 =−A2 F +32τA1


51 32τ
P1 = − kZ 2+−F ,32 A1 2τω∗
2τ
P2 = − kZ 2 (1 + Fτ) A2 ,3A1b2 τ(τ + 1) + 1,
A1 = b +4
3A2 = −1 + b1 +− b1 + b2 η,
2b/2,
F=
1 + bτ/2bτ/2b1 = −,
1 + bτ/2b2 τ2 /2
b2 = 2b1 +,(1 + bτ/2)2kθω∗ = .
Ln(9)

(10a)
(10b)
(10c)
(10d)
(10e)
(10f)
(10g)
(10h)
(10i)

Here T ⊥ = T // is used for simplicity. In the case withkz 2  1, Eq. (9) is rewritten to be(ω − ω1 − δω) ω2 + δωω + T 1 = 0,
(11)
by use of small δω, whereP1 ω1 + P2,
ω1 2 + P1
P1 2 − P2 ω1
T1 =.
ω1 2 + P1

δω = −(12a)(12b)

An unstable ITG mode (slow wave mode) is given in thecondition;

1401158-2δω2 − 4T 1 ∼ −4T 1 < 0,(13)Plasma and Fusion Research: Regular Articles

Volume 14, 1401158 (2019)and in the case with η ∼ 1 and small τ,
T1 ∼ −P2
ω1(14)

gives the following reduced condition;
A2 > 0.

(15)

From this inequality, the threshold value ηc is given asηc =2 5+ bτ.
3 9(16)This is the threshold value of the unstable ITG mode forthe parameters of the linear device in Sec. 2.
The critical value ηc changes depending on the magnitude of τ. This dependency is related to the FLR eﬀect.
In the case with b2 = 0 the FLR eﬀect is smaller and the τdependency of the threshold becomes weaker;
ηc =2 1+ bτ.
3 9(17)

Without the FLR eﬀect, b = 0, no dependency on τ appearsasηc =2,3(18)radial eigenmode structures of n, T ⊥ and φ. These radially broad modes with azimuthal mode number m = 2and axial mode number l = 1 are most unstable with theplasma parameters in Sec. 2. The values of the perpendicular wavenumber, which is a parameter in local stability
analyses, can be evaluated from the radial structure [14].
The mode includes several Bessel components. Figure 2
shows spectra with the radial mode number j. The component with j = 1 is dominant, and the dependency onη shows that amplitudes of the components with larger j
become smaller with larger η. The dominant mode hask⊥ ρi = 0.6 < 1, so the Padé approximation for the gyroaveraged potential is appropriate. Existence of components with larger j makes the mode structure rather radiallysteeper in the cases with larger η.
The dependency of the growth rate on τ and η areshown in Fig. 3 (a). The ITG instability is unstable when
the ratio of the ion temperature gradient to the density gradient η exceeds the threshold value ηc . The dashed line inFig. 3 (b) represents the threshold in Eq. (16), which welldescribes the numerical result. Figure 4 shows comparison of ηc in Figs. 3 (a) and 3 (b). The local analysis, asin Fig. 3 (b), confirms to reproduce the global analysis re-which is the well-known threshold of the ITG mode.
In the case of the resistive drift wave instability, which
can be unstable in the parameter set of linear devices, breakof the Boltzmann relation n  Ψ is important for the destabilization [16]. In this model, the FLR eﬀect can break theBoltzmann relation, as Eq. (5), which aﬀects the instability.
Equation (5) is transformed toΨ=1(n − FT ⊥ ) .
A1(19)The 2nd term in the right hand side of Eq. (19) corresponds
to the break. Inclusion of higher order δω terms in Eq. (11)gives the relation for δω as(20)2ω1 δω2 + P1 + ω1 2 δω + P1 ω1 + P2 = 0.
IfP1 + ω1 2

2− 8ω1 (P1 ω1 + P2 ) < 0(21)Fig. 1 Radial eigenmode structures with m = 2, l = 1, η = 1.0,
τ = 0.5. Radial profiles of n, T ⊥ and φ are shown, wheresolid and dashed lines represent the real and imaginaryparts, respectively.

is satisfied, a ITG mode (fast wave mode) becomes unstable, but in η = 0 case, this condition cannot be satisfied, asfunction F in Eq. (10f) gives τF < 1. On the other hand,small τ limit gives F ∼ b/2, and using this form, condition(21) can be satisfied with τ ∼ 1, but this is incorrect application of the small τ approximation. The break by theFLR eﬀect is not strong enough, so this eﬀect alone cannot
make the mode unstable.

4. Mode Structure
The numerical analysis for the global model using theset of Eqs. (1 - 5) considers the boundary condition to determines the radial mode structure. Figure 1 shows theFig. 2 Dependency of the radial mode number spectra on η withm = 2, l = 1, τ = 0.5.

1401158-3

Plasma and Fusion Research: Regular Articles

Volume 14, 1401158 (2019)

Fig. 3 Dependency of the growth rate on τ and η with m = 2,l = 1. The results of the (a) global analysis, (b) local,(c) local with b2 = 0 and (d) local without the FLR effect are plotted. The dashed lines represent the analyticalexpressions of the thresholds.

with stronger threshold dependence on τ [19]. In linear devices the magnetic field is homogeneous, so the target inthis article is a slab ITG mode, which is induced by compressibility with parallel ion motion. The FLR terms include τ dependency related to orbit averaging on Ψ , which
gives the threshold dependency. Break of the Boltzmann
relation by the FLR eﬀect is not strong as to generate another unstable branch of the ITG mode, so resistivity israther important with the parameter set of linear devices tobreak the Boltzmann relation, which excites the resistive
drift wave [11].
Linear characteristics have been confirmed with constant background density and temperature gradients. We
are developing the code to include the background profiles, in which couplings between diﬀerent Bessel components must be considered. The routine can be extended toapply to nonlinear calculation, which will be carried out
for investigating dynamics in turbulent plasmas. Torus devices are used for fusion researches, and for investigationof global structures and their dynamics in torus, PLATO
project is now progressing to establish a plasma turbulence
observatory [20]. Construction of the platform for analyzing 3-D dynamics of plasma instabilities with a variety
of simulation codes is necessary for comparison between
experiments and simulations both in cylindrical and torus
plasmas.

Acknowledgments
Fig. 4 Comparison of the thresholds for the ITG instability
from the global and local models with m = 2, l = 1.
The dashed line represents the analytical expression in
Eq. (16).

sult by adjusting the wavenumber to that obtained from theglobal analysis.
The dependency of ηc on τ is confirmed numerically in the smaller or no FLR eﬀect cases with b2 = 0(Fig. 3 (c)) and (Fig. 3 (d)), whose analytical forms with thelocal model are given in Eqs. (17) and (18), respectively.

5. Summary
Numerical analyses using a gyro-fluid model havebeen performed to investigate the FLR eﬀect on ITG instability in cylindrical plasmas. A spectrum code withFourier-Bessel expansion has been developed for the analysis of global mode structures. The analytical expressionof ηc is obtained from the local dispersion relation. Itsτ dependency comes from the FLR eﬀect. The ratio between ion and electron temperatures is an important parameter for plasma performance in torus devices [17], andthe ITG instability threshold is investigated using fluid andkinetic approaches [18]. In torus configurations, resonancewith magnetic drift appears to give a toroidal ITG modeThis work was partly supported by JSPS KAKENHI Grant Number (JP16K06938, JP17H06089), by
the collaboration program of NIFS (NIFS19KNST144,NIFS18KNXN373) and RIAM of Kyushu University.
[1] X. Garbet et al., Nucl. Fusion 50, 043002 (2010).
[2] R.E. Waltz et al., Phys. Plasmas 4, 2482 (1997).
[3] J.E. Kinsey et al., Phys. Plasmas 12, 052503 (2005).
[4] G.M. Staebler et al., Phys. Plasmas 14, 055909 (2007).
[5] R.E. Waltz et al., Phys. Plasmas 1, 2229 (1994).
[6] A. Kendl et al., Phys. Plasmas 17, 072302 (2010).
[7] P.H. Diamond et al., Plasma Phys. Control. Fusion 47, R35(2005).
[8] G.R. Tynan et al., Plasma Phys. Control. Fusion 51, 113001(2009).
[9] A.K. Sen et al., Phys. Rev. Lett. 66, 429 (1991).
[10] S. Inagaki et al., Sci. Reps. 6, 22189 (2016).
[11] N. Kasuya, M. Yagi, K. Itoh and S.-I. Itoh, Phys. Plasmas
15, 052302 (2008).
[12] Y. Miwa et al., Plasma Fusion Res. 8, 2403133 (2013).
[13] G. Hattori et al., Plasma Fusion Res. 10, 3401060 (2015).
[14] T. Ohno et al., Plasma Fusion Res. 13, 1401081 (2018).
[15] W. Dorland and W. Hammet, Phys. Fluids B 5, 812 (1993).
[16] A. Hasegawa and M. Wakatani, Phys. Rev. Lett. 50, 682(1983).
[17] C.C. Petty et al., Phys. Rev. Lett. 83, 3661 (1999).
[18] J. Weiland, Collective Modes in Inhomogeneous Plasma(Institute of Physics, Bristrol, 2000).
[19] A. Casati et al., Phys. Plasmas 15, 042310 (2008).
[20] A. Fujisawa, AIP Conf. Proc. 1993, 020011 (2018).

1401158-4