Global simulation of ion temperature gradient instabilities in a field-reversed configuration
场反转配置中离子温度梯度不稳定性的全局模拟1

Cite as: Phys. Plasmas 26, 042506 (2019); doi: 10.1063/1.5087079
引用： Phys. Plasmas 26， 042506 （2019）;doi： 10.1063/1.5087079Submitted: 27 December 2018. Accepted: 25 March 2019 .
提交日期：2018年12月27日。录用日期： 2019-03-25 .Published Online: 15 April 2019
网络出版日期： 2019-04-15

1. J. Bao, $^{1, 2, 3, a)}$ C. K. Lau, $^{4}$ (D) Z. Lin, $^{7, b)}$ (D) H. Y. Wang, $^{1, 5}$ D. P. Fulton, $^{4}$ S. Dettrick, $^{4}$ and T. Tajima $^{4}$
1. J. Bao， $^{1, 2, 3, a)}$ C. K. Lau， $^{4}$ （D） Z. Lin， $^{7, b)}$ （D） H. Y. Wang， $^{1, 5}$ D. P. Fulton， $^{4}$ S. Dettrick， $^{4}$ and T. Tajima $^{4}$

2. AFFILIATIONS 2. 隶属关系

'University of California, Irvine, California 92697, USA
'加州大学欧文分校，加利福尼亚州 92697，美国

^{2}

Beijing National Laboratory for Condensed Matter Physics and CAS Key Laboratory of Soft Matter Physics, Institute of Physics, Chinese Academy of Sciences, Beijing 100190, China

^{2}

中国科学院物理研究所，北京凝聚态物理国家实验室，中国科学院软物质物理重点实验室，北京 100190

^{3}

University of Chinese Academy of Sciences, Beijing 100049, China

^{3}

中国科学院大学，北京100049，中国

^{4} T A E

Technologies, Inc., 19631 Pauling, Foothill Ranch, California 92610, USA

^{4} T A E

Technologies， Inc.， 19631 Pauling， Foothill Ranch， California 92610，美国

^{5}

Fusion Simulation Center, Peking University, Beijing 100871, China

^{5}

北京大学聚变模拟中心，北京100871，中国

a)jianbao0626@gmail.com

b) zhihongl@uci.edu b） zhihongl@uci.edu

3. ABSTRACT 3. 摘要

We investigate the global properties of drift waves in the beam driven field-reversed configuration (FRC), the C2-U device, in which the central FRC and its scrape-off layer (SOL) plasma are connected with the formation sections and divertors. The ion temperature gradient modes are globally connected and unstable across these regions, while they are linearly stable inside the FRC separatrix. The unstable global drift waves in the SOL show an axially varying structure that is less intense near the central FRC region and the mirror throat areas, while being more robust in the bad curvature formation exit areas.
我们研究了束流驱动场反转配置（FRC）中漂移波的全局特性，即C2-U器件，其中中央FRC及其刮擦层（SOL）等离子体与地层部分和分流器相连。离子温度梯度模式在这些区域之间是全局相连的，并且不稳定，而它们在FRC分离层内是线性稳定的。SOL中不稳定的全局漂移波呈现轴向变化结构，在FRC中心区域和镜像喉咙区域附近强度较低，而在不良曲率形成出口区域则更为稳健。

Published under license by AIP Publishing. https://doi.org/10.1063/1.5087079
经AIP Publishing许可出版。https://doi.org/10.1063/1.5087079

4. INTRODUCTION 4. 引言

A field-reversed configuration (FRC) is an elongated prolate compact toroid (CT) with magnetic fields predominantly along the poloidal direction, which consists of a core with closed field lines and a scrapeoff layer (SOL) with open field lines.

^{1}

As a fusion reactor concept, FRC has many advantages: The average beta (the ratio between plasma kinetic pressure and magnetic energy density) is close to unity, which suggests a much cheaper fusion energy than the tokamak. The compact shape and simple geometry of FRC also lead to construction convenience and high magnetic efficiency. The SOL region extends to the device ends and forms natural divertors, which are far away from the core region and allow extraction of fusion energy without restriction. TAE Technologies, Inc. has launched a series of FRC experiments.

^{2 - 7}

A significant energetic ion population generated from neutral beam injection (NBI) can suppress the macroinstability.

^{8, 9}

Meanwhile, the large orbit size effects of the energetic particle would not destabilize the microturbulence and the ion scale turbulent transport is suppressed.

^{10}

These experimental efforts lead to the sustainment of beam-driven hot FRC plasmas for more than

5 m s

in C-2U experiments, which is in the confinement regime limited by turbulent transport. In FRC, the open field line SOL region is connected to the closed field line core region, and the turbulence in these two regions affects each other. In recent

C 2 - U

experiment, it is found that ion-scale turbulence fluctuation is suppressed in the core, while in the SOL, ionand electron-scale turbulence is observed.

^{10}

Thus, it is important to understand the transport mechanism in FRC for the improvement of plasma confinement. 1D and 2D magnetohydrodynamic (MHD) codes have been built up to model the global FRC transport, which requires turbulence simulation codes to provide the transport coeffcients.

^{11 - 13}

First-principles particle-in-cell simulation is a powerful tool to study the fusion plasmas combining with theory, which has successful applications in understanding the anomalous transport in tokamak plasmas.

^{14 - 24}

For the particle-in-cell study of turbulence in FRC, some pioneer works

^{25 - 27}

have been performed based on the state-of-the-art fusion plasma simulation code: Gyrokinetic Toroidal Code (GTC). GTC has been successfully applied to simulate microturbulence,

^{21, 24}

energetic particle transport,

^{28}

Alfven eigenmodes,

^{29, 30}

and MHD instabilities

^{31, 32}

in toroidal plasmas. An upgrade to the FRC geometry in the Boozer coordinates of GTC has been carried out by Fulton et al.,

^{25, 26}

and local gyrokinetic particle simulation study by Lau et al.

^{10, 27}

shows that the drift wave is stable in the FRC core due to the large orbit size, magnetic well geometry, and short electron transit length. Meanwhile, a new global particle-in-cell FRC code, ANC, has been developed by incorporating core and SOL regions across the separatrix. ANC simulations show that ion scale turbulence can spread from the SOL to the core.

^{33, 34}

磁场反转构型（FRC）是一种细长的匯实紧凑型环形线圈（CT），其磁场主要沿极向方向，由具有闭合磁力线的磁芯和具有开放磁力线的刮擦层（SOL）组成。

^{1}

作为一种聚变反应堆概念，FRC具有许多优点：平均β（等离子体动压和磁能密度之间的比率）接近统一，这表明聚变能量比托卡马克便宜得多。FRC的紧凑形状和简单的几何形状也带来了施工便利性和高磁效率。SOL区域延伸到器件末端并形成天然分流器，远离核心区域，允许不受限制地提取聚变能。TAE Technologies， Inc. 推出了一系列 FRC 实验。

^{2 - 7}

中性束注入（NBI）产生的显著高能离子群可以抑制宏观不稳定性。

^{8, 9}

同时，高能粒子的大轨道尺寸效应不会破坏微湍流的稳定性，抑制离子尺度的湍流输运。

^{10}

这些实验工作导致光束驱动的热FRC等离子体的维持时间超过

5 m s

C-2U实验，C-2U实验处于受湍流传输限制的约束状态。在FRC中，开放场线SOL区域与封闭场线核心区域相连，这两个区域的湍流相互影响。在最近的

C 2 - U

实验中，发现核心中离子尺度的湍流波动受到抑制，而在溶胶中观察到离子和电子尺度的湍流。

^{10}

因此，了解FRC中改善血浆限制的转运机制非常重要。已经建立了一维和二维磁流体动力学（MHD）代码来模拟全球 FRC 传输，这需要湍流模拟代码来提供传输系数。

^{11 - 13}

第一性原理粒子在细胞模拟是研究聚变等离子体与理论相结合的有力工具，在理解托卡马克等离子体的异常输运方面具有成功的应用。

^{14 - 24}

对于FRC中湍流的粒子细胞研究，基于最先进的聚变等离子体模拟代码：回动环形码（GTC）进行了一些开创性的工作

^{25 - 27}

。GTC已成功应用于模拟环形等离子体中的微湍流、

^{21, 24}

高能粒子输运、

^{28}

Alfven特征模态

^{29, 30}

和MHD不稳定性

^{31, 32}

。Fulton等人对GTC的Boozer坐标中的FRC几何进行了升级，

^{25, 26}

Lau等人

^{10, 27}

的局部陀螺动力学粒子模拟研究表明，由于轨道尺寸大、磁井几何形状大、电子传输长度短，漂移波在FRC核心中是稳定的。同时，通过整合整个分离区的核心和 SOL 区域，开发了一种新的全球细胞颗粒 FRC 代码 ANC。ANC模拟表明，离子级湍流可以从SOL扩散到核心。

^{33, 34}

In order to study the turbulent transports globally up to the divertor region for FRC, a new GTC family code, gyrokinetic toroidal code-X (GTC-X), is developed in this work by refactoring the coordinate system and geometry of the original GTC code, i.e., change the Boozer coordinates to cylindrical coordinates and change the geometry from tokamak and stellarator to FRC. GTC-X enables the crossseparatrix simulation with a field aligned mesh covering the whole geometry of the FRC. Compared to the original GTC code, both the particle trajectory and the Poisson solver are newly written as well as the simulation grids in the GTC-X code. The GTC mixed-model OpenMP-MPI parallelization

^{35}

is adopted in GTC-X. This paper mainly presents the numerical developments, code verification, and initial results of ion temperature gradient (ITG) modes in the global FRC geometry. GTC-X global simulations show that ITG is unstable in the SOL and stable in the core, which is consistent with previous local simulations and experimental observations. We find that the ITG mode grows along the field line direction in the SOL and shows an axial variation. The maximum amplitude of the ITG mode is in the formation region with bad curvature, while the mode amplitude is small in the central FRC region. The mode structure in the SOL is sensitive to the parallel domain size, which experiences a transition from even parity to odd parity when increasing the domain size.
为了研究FRC全球至分流器区域的湍流传输，通过重构原始GTC代码的坐标系和几何结构，将Boozer坐标更改为圆柱坐标，将几何形状从托卡马克和恒星器更改为FRC，开发了一种新的GTC系列代码，即陀螺动力学环形码-X（GTC-X）。GTC-X 通过覆盖 FRC 整个几何形状的现场对齐网格实现交叉分离模拟。与原始 GTC 代码相比，粒子轨迹和泊松求解器以及 GTC-X 代码中的仿真网格都是新编写的。GTC-X采用GTC混合模型OpenMP-MPI并行化

^{35}

。本文主要介绍了全球FRC几何中离子温度梯度（ITG）模式的数值发展、代码验证和初步结果。GTC-X全局仿真结果显示，ITG在SOL中不稳定，在核心中稳定，这与之前的局部仿真和实验观测结果一致。我们发现 ITG 模式在 SOL 中沿场线方向增长并显示出轴向变化。ITG模式的最大振幅在曲率不好的地层区域，而模态振幅在FRC中心区域较小。SOL 中的模式结构对并行域大小很敏感，当增加域大小时，并行域大小会经历从偶奇偶校验到奇偶校验的过渡。

This paper is organized as follows. In Sec. II, we introduce the global FRC geometry implementation. The gyrokinetic particle simulation model for FRC is described in Sec. III. The benchmark simulation results are shown in Sec. IV. In Sec. V, the global simulation of ITG modes is described. The conclusion is discussed in Sec. VI.
本文的组织结构如下。在第二部分中，我们介绍了全局FRC几何实现。FRC的陀螺动力学粒子模拟模型在第III节中进行了描述。基准仿真结果显示在第 IV 节中。在第五节中，描述了ITG模式的全局仿真。结论在第六节中讨论。

5. GLOBAL FRC GEOMETRY IMPLEMENTATION
5. 全局FRC几何实现

In order to avoid the singularity of magnetic coordinates at the separatrix,

^{22, 36}

we adapt the cylindrical coordinate system for global FRC simulation with

(R, ζ, Z)

, where the 3 independent unit vectors satisfy the right hand rule:

\hat{R} \times \hat{ζ} \cdot \hat{Z} = 1

. The poloidal magnetic flux

ψ

of FRC equilibrium and cylindrical coordinates used in GTC-X is shown in Fig. 1, which is calculated by an axisymmetric force balance FRC equilibrium solver: LR_eqMI code.

^{37}

The equilibrium box size is normalized by the radial position of magnetic axis:

R = R_{0} = 26.8 c m

, i.e., the distance between the magnetic axis and the cylinder axis. There are several mirror plugs in the SOL region aiming at decreasing the particle end loss, and the expanded divertors are located at the ends of open field lines, where we can apply the edge biasing via the plasma-gun electrodes to improve the confinement.

^{4}

In this section, based on the characteristics of FRC equilibrium, we will introduce the algorithms used in GTC-X for global particle-in-cell modeling of FRC.
为了避免磁坐标在分离点处的奇异性，

^{22, 36}

我们适配圆柱坐标系进行

(R, ζ, Z)

全局FRC模拟，其中3个独立的单位向量满足右手定则：

\hat{R} \times \hat{ζ} \cdot \hat{Z} = 1

。GTC-X中使用

ψ

的FRC平衡和圆柱坐标的极化磁通量如图1所示，该磁通量由轴对称力平衡FRC平衡求解器：LR_eqMI代码计算。

^{37}

平衡箱尺寸由磁轴的径向位置归一化：

R = R_{0} = 26.8 c m

，即磁轴和圆柱轴之间的距离。SOL区域中有几个镜像塞，旨在减少粒子端损耗，扩展的分流器位于开放场线的末端，在那里我们可以通过等离子枪电极应用边缘偏置来改善约束。

^{4}

在本节中，我们将基于FRC平衡的特性，介绍GTC-X中用于FRC全局粒子单元建模的算法。

6. A. Magnetic field representation in cylindrical coordinates
6. A. 圆柱坐标中的磁场表示

The magnetic field and associated derivatives commonly appear in the particle dynamic equations for the simulation of magnetized plasmas; thus, it is important that the magnetic field satisfies

\nabla \cdot B

= 0

numerically. The magnetic coordinates enable the free divergence representation for the magnetic field as

B = \nabla α \times \nabla β

, where

α

and

β

are coordinates which vary along the directions orthogonal to the magnetic field. However, magnetic coordinates fail to address the simulation containing different geometric topologies with a separatrix.

^{36}

Thus, we apply the cylindrical coordinate system as the basic coordinates. In order to guarantee the free divergence property for the magnetic field in cylindrical coordinates, we use the poloidal magnetic flux

ψ

to calculate the magnetic field components and their derivatives and, thus, enforce the consistency between each component. In FRC, the equilibrium magnetic field

B

has no toroidal component and can be expressed as

B = \nabla ψ \times \nabla ζ = B_{R} \hat{R} + B_{Z} \hat{Z}

. The magnetic field strength in radial and axial directions can then be derived as
磁场和相关导数通常出现在磁化等离子体模拟的粒子动力学方程中;因此，磁场在

= 0

数值上满足

\nabla \cdot B

是很重要的。磁坐标使磁场的自由发散表示为

B = \nabla α \times \nabla β

，其中

α

和

β

是沿与磁场正交的方向变化的坐标。然而，磁坐标无法解决包含不同几何拓扑的模拟问题。

^{36}

因此，我们应用圆柱坐标系作为基本坐标。为了保证圆柱坐标中磁场的自由发散性，我们使用极化磁通量

ψ

来计算磁场分量及其导数，从而加强每个分量之间的一致性。在 FRC 中，平衡磁场

B

没有环形分量，可以表示为

B = \nabla ψ \times \nabla ζ = B_{R} \hat{R} + B_{Z} \hat{Z}

。径向和轴向的磁场强度可以推导出为

B_{R} = - \frac{1}{R} \frac{\partial ψ}{\partial Z}

and

B_{Z} = \frac{1}{R} \frac{\partial ψ}{\partial R}

From equilibrium calculated by the LR_eqMI code, we can get the value of poloidal magnetic flux

ψ

over the whole FRC geometry on the equal space grids in the

(R, Z)

plane as shown in Fig. 1 , and

L S R = 150

and

L S Z = 401

are the equilibrium radial and axial grid numbers, respectively. By using the value on coarse equilibrium grids, we can use the quadratic spline function to calculate

ψ

at the arbitrary location

(R, Z)

inside the equilibrium domain as
根据LR_eqMI码计算的平衡，我们可以得到

(R, Z)

平面内等空间网格上整个FRC几何形状的极化磁通量

ψ

值，如图1所示，

L S R = 150

分别

L S Z = 401

是平衡径向和轴向网格数。通过使用粗平衡网格上的值，我们可以使用二次样条函数在平衡域内的任意位置

(R, Z)

计算

ψ

为

FIG. 1. Contour plot of

ψ / | ψ_{0} |

for global FRC geometry, where

ψ_{0}

is the poloidal magnetic flux value at the magnetic axis as shown by the green star. The black solid lines represent the different field lines (contour line of

ψ

), and the red line represents the separatrix. The arrows denote the directions of cylindrical coordinates.
图 1.全局FRC几何形状的

ψ / | ψ_{0} |

等值线图，其中

ψ_{0}

是磁轴处的极化磁通量值，如绿色星形所示。黑色实线表示不同的场线（等值线），

ψ

红线表示分离线。箭头表示圆柱坐标的方向。

\begin{aligned} ψ (R, Z) & = ψ (1, i, j) + ψ (2, i, j) Δ R + ψ (3, i, j) Δ R^{2} \\ + ψ (4, i, j) Δ Z + ψ (5, i, j) Δ R Δ Z + ψ (6, i, j) Δ R^{2} Δ Z \\ + ψ (7, i, j) Δ Z^{2} + ψ (8, i, j) Δ R Δ Z^{2} + ψ (9, i, j) Δ R^{2} Δ Z^{2}, \end{aligned}

where

i \in [1, L S R]

and

j \in [1, L S Z]

are the radial and axial indexes of equilibrium grids, and

Δ R = R - R_{i}, Δ Z = Z - Z_{j}, R_{i} \leq R < R_{i + 1}

and

Z_{j} \leq Z < Z_{j + 1}

.
其中

i \in [1, L S R]

和

j \in [1, L S Z]

是平衡网格的径向和轴向索引，和

Δ R = R - R_{i}, Δ Z = Z - Z_{j}, R_{i} \leq R < R_{i + 1}

Z_{j} \leq Z < Z_{j + 1}

。

It is straightforward to show that

\nabla \cdot B = 0

is guaranteed theoretically and numerically
很明显，

\nabla \cdot B = 0

这在理论上和数值上都是有保证的

\nabla \cdot B = \frac{1}{R} \frac{\partial}{\partial R} (R B_{R}) + \frac{\partial B_{Z}}{\partial Z} = 0 .

7. B. Field line coordinates for the perturbed field calculation
7. B. 扰动场计算的场线坐标

In magnetized plasmas, the gyrocenter drift motion across the magnetic field is much slower than the parallel motion along the field line; thus, the wave pattern is always anisotropic in the parallel and perpendicular directions with

k_{∥} ≪ k_{⊥} (k_{∥}

and

k_{⊥}

are the parallel and perpendicular wave vectors). In order to improve the numerical efficiency and accuracy, the field aligned mesh is widely adapted for particle-in-cell simulation of magnetized plasmas, i.e., the grids are aligned along the magnetic field direction with only a small number in the parallel direction, which can dramatically suppress the high

k_{∥}

noise and save computational cost without sacrificing key physics dominated by small

k_{∥} \cdot^{23, 38}

In global FRC simulation, we setup a field aligned mesh in both core and scrape-off layer (SOL) regions across the separatrix in cylindrical coordinates. Due to the fact that the magnetic field is not uniform in FRC, the field aligned mesh is not regular in cylindrical coordinates. For solving perturbed fields as well as particle-grid gather-scatter operation, we create the field line coordinates for core

(ψ, S_{c})

and SOL

(ψ, S_{S})

regions on the poloidal plane, separately, where

S_{c}

and

S_{S}

represent the normalized field line distances along the magnetic field line direction in core and SOL regions, and the mesh is regular in the corresponding field line coordinates. Thus, in GTC-X, we use two different coordinate systems: cylindrical coordinates and field line coordinates to represent the location.
在磁化等离子体中，陀螺中心在磁场上的漂移运动比沿磁场线的平行运动慢得多;因此，波型在平行和垂直方向上

k_{∥} ≪ k_{⊥} (k_{∥}

始终是各向异性的，并且是

k_{⊥}

平行和垂直的波矢量）。为了提高数值效率和精度，磁场对齐网格被广泛应用于磁化等离子体的粒子在单元模拟中，即网格沿磁场方向排列，平行方向上只有少量数字，可以显著抑制高

k_{∥}

噪声，节省计算成本，同时又不牺牲以小

k_{∥} \cdot^{23, 38}

为主的关键物理场在全局FRC仿真中，我们在圆柱坐标的分离线的核心和刮擦层（SOL）区域设置了一个场对齐的网格。由于 FRC 中的磁场不均匀，因此圆柱坐标中的磁场对齐网格不规则。为了求解扰动场以及粒子网格聚集散射操作，我们分别在极化平面上创建核心

(ψ, S_{c})

和 SOL

(ψ, S_{S})

区域的场线坐标，其中

S_{c}

和

S_{S}

表示核心和 SOL 区域中沿磁力线方向的归一化场线距离，网格在相应的场线坐标中是规则的。因此，在 GTC-X 中，我们使用两种不同的坐标系：圆柱坐标和场线坐标来表示位置。

The simulation domain is different from equilibrium shown in Fig. 1. Because drift wave instabilities and associated transports are anisotropic in perpendicular and parallel directions, we choose the simulation domain based on perpendicular coordinates

ψ

(we do not need to consider about

ζ

domain because it is toroidally symmetric from

(0, 2 π)

), i.e., the inner boundary in the core and the outer boundary in SOL are labeled by poloidal magnetic flux:

ψ_{0}

and

ψ_{1}

. Furthermore, left and right boundaries in the SOL region are given by

Z_{0}

and

Z_{1}

, where

Z_{0} = - Z_{1}

is symmetric with respect to outer midplane

Z = 0

.
仿真域与图1所示的平衡域不同。由于漂移波的不稳定性和相关的输运在垂直和平行方向上是各向异性的，因此我们根据垂直坐标选择仿真域

ψ

（我们不需要考虑

ζ

域，因为它是环形对称

(0, 2 π)

的），即磁芯中的内边界和SOL中的外边界由极化磁通量标记：

ψ_{0}

和

ψ_{1}

.此外，SOL 区域中的左右边界由

Z_{0}

和

Z_{1}

给出，其中

Z_{0} = - Z_{1}

相对于外中平面

Z = 0

是对称的。

First, we define

S_{c}

and

S_{S}

by tracing each field line on each

ψ

grid in the core and SOL regions as
首先，

S_{S}

我们将核心区域和 SOL 区域中每个

ψ

网格上的每条场线定义为

S_{c}

ψ (i_{c}) = ψ_{0} + δ ψ_{c} \times (i_{c} - 1),

where

δ ψ_{c} = (ψ_{X} - ψ_{0}) / (l s p c - 1), ψ_{X}

is the value of

ψ

at the separatrix,

l s p c

is the spline resolution in the core region, and

1 \leq i_{c} \leq l s p c

. And
其中

δ ψ_{c} = (ψ_{X} - ψ_{0}) / (l s p c - 1), ψ_{X}

是

ψ

的值在分离点上，

l s p c

是核心区域中的样条分辨率，以及

1 \leq i_{c} \leq l s p c

。和

ψ (i_{S}) = ψ_{X} + δ ψ_{S} \times (i_{S} - 1),

where

δ ψ_{S} = (ψ_{1} - ψ_{X}) / (l s p s - 1)

, lsps is the spline resolution in SOL region, and

1 \leq i_{S} \leq l

sps. Both

S_{c}

and

S_{S}

are normalized by the field line length at each

ψ

grid and range from 0 to 1 . In the core region,

S_{c}

starts at the outer midplane with

Z = 0

and increases along the clockwise direction, and grids at

S_{c} = 0

and

S_{c} = 1

are overlapped. In the SOL region, the parallel coordinate

S_{S}

starts on the left boundary

Z_{0}

and ends on the right boundary

Z_{1}

.
其中

δ ψ_{S} = (ψ_{1} - ψ_{X}) / (l s p s - 1)

，lsps 是 SOL 区域中的样条分辨率，sps

1 \leq i_{S} \leq l

.两者都

S_{c}

S_{S}

由每个

ψ

网格处的场线长度归一化，范围从 0 到 1。在核心区域中，

S_{c}

从外中平面开始，沿顺时针方向增加

Z = 0

，网格在

S_{c} = 0

和

S_{c} = 1

处重叠。在 SOL 区域中，平行坐标

S_{S}

从左边界

Z_{0}

开始，到右边界

Z_{1}

结束。

Next, considering the properties of geometry topology of core and SOL regions in FRC, forward spline functions

S_{c} = S_{c} [ψ (R, Z)

θ (R, Z)]

and

S_{S} = S_{S} [ψ (R, Z), Z]

are created for the transformation from cylindrical coordinates to field line coordinates, where

\sin (θ)

= Z / \sqrt{Z^{2} + {(R - R_{0})}^{2}}

is the geometric angle with respect to the magnetic axis position

(R = R_{0}, Z = 0)

. The uniform scale of geometric angle

θ

for creating spline function

S_{c} (ψ, θ)

is
接下来，考虑FRC中磁芯和SOL区域的几何拓扑性质，正向样条函数

S_{c} = S_{c} [ψ (R, Z)

，

θ (R, Z)]

并

S_{S} = S_{S} [ψ (R, Z), Z]

创建了从圆柱坐标到场线坐标的转换，其中

\sin (θ)

= Z / \sqrt{Z^{2} + {(R - R_{0})}^{2}}

是相对于磁轴位置

(R = R_{0}, Z = 0)

的几何角度。用于创建样条函数

S_{c} (ψ, θ)

的几何角度

θ

的均匀比例为

θ (j_{c}) = δ θ \times (j_{c} - 1)

where

δ θ = 2 π / (l s t c - 1)

, lstc is the spline resolution in the

θ

direction, and

1 \leq j_{c} \leq l s t c

. By tracing each field line and calculating the

S_{c}

value on the uniform scales of

ψ

and

θ

, the spline function

S_{c} (ψ, θ)

can be derived as
其中

δ θ = 2 π / (l s t c - 1)

， lstc 是方向上的

θ

样条分辨率，而

1 \leq j_{c} \leq l s t c

.通过跟踪每条场线并计算

ψ

和

θ

的均匀刻度上的

S_{c}

值，样条函数

S_{c} (ψ, θ)

可以导出为

\begin{aligned} S_{c} (ψ, θ) = & S_{c} (1, i_{c}, j_{c}) + S_{c} (2, i_{c}, j_{c}) Δ ψ + S_{c} (3, i_{c}, j_{c}) Δ ψ^{2} \\ + S_{c} (4, i_{c}, j_{c}) Δ θ + S_{c} (5, i_{c}, j_{c}) Δ ψ Δ θ + S_{c} (6, i_{c}, j_{c}) Δ ψ^{2} Δ θ \\ + S_{c} (7, i_{c}, j_{c}) Δ θ^{2} + S_{c} (8, i_{c}, j_{c}) Δ ψ Δ θ^{2} + S_{c} (9, i_{c}, j_{c}) Δ ψ^{2} Δ θ^{2} . \end{aligned}

The uniform scale of

Z

for creating spline function

S_{S} (ψ, Z)

is
用于创建样条函数

S_{S} (ψ, Z)

的

Z

均匀比例为

Z (j_{S}) = δ Z_{S} \times (j_{S} - 1)

where

δ Z_{S} = (Z_{1} - Z_{0}) / (l s z s - 1)

, lszs is the spline resolution in the

Z

direction, and

1 \leq j_{S} \leq l s z s

. By tracing each field line and calculating the

S_{S}

value on the uniform scales of

ψ

and

Z, S_{S} (ψ, Z)

is given as
其中

δ Z_{S} = (Z_{1} - Z_{0}) / (l s z s - 1)

，lszs 是方向上的

Z

样条分辨率，而

1 \leq j_{S} \leq l s z s

。通过跟踪每条场线并计算

ψ

和

Z, S_{S} (ψ, Z)

的均匀刻度上的

S_{S}

值，给出为

\begin{aligned} S_{S} (ψ, Z) = & S_{S} (1, i_{S}, j_{S}) + S_{S} (2, i_{S}, j_{S}) Δ ψ + S_{S} (3, i_{S}, j_{S}) Δ ψ^{2} \\ + S_{S} (4, i_{S}, j_{S}) Δ Z + S_{S} (5, i_{S}, j_{S}) Δ ψ Δ Z + S_{S} (6, i_{S}, j_{S}) Δ ψ^{2} Δ Z \\ + S_{S} (7, i_{S}, j_{S}) Δ Z^{2} + S_{S} (8, i_{S}, j_{S}) Δ ψ Δ Z^{2} + S_{S} (9, i_{S}, j_{S}) Δ ψ^{2} Δ Z^{2} . \end{aligned}

Then, we create the inverse spline functions

R_{c} (ψ, S_{c})

and

Z_{c} (ψ, S_{c})

for the core region, and

R_{S} (ψ, S_{S})

and

Z_{S} (ψ, S_{S})

for the SOL region. The uniform scales of

S_{c}

and

S_{S}

for the spline function are
然后，我们为核心区域

R_{S} (ψ, S_{S})

Z_{S} (ψ, S_{S})

创建反样条函数

R_{c} (ψ, S_{c})

和

Z_{c} (ψ, S_{c})

以及 SOL 区域。样条函数

S_{c}

的均匀尺度为

S_{S}

S_{c} (j_{c}) = δ S_{c} \times (j_{c} - 1)

where

δ S_{c} = 1 / (l s s c - 1), l s s c

is the spline resolution in the parallel direction in the core region, and

1 \leq j_{c} \leq l s s c

其中

δ S_{c} = 1 / (l s s c - 1), l s s c

是核心区域平行方向上的样条分辨率，以及

1 \leq j_{c} \leq l s s c

S_{S} (j_{S}) = δ S_{S} \times (j_{S} - 1),

where

δ S_{S} = 1 / (l s s s - 1)

, lsss is the spline resolution in the parallel direction in the SOL region, and

1 \leq j_{S} \leq l s s s

. It is straightforward to get the values

R_{c} [ψ (i_{c}), S_{c} (j_{c})], Z_{c} [ψ (i_{c}), S_{c} (j_{c})], R_{S} [ψ (i_{S}), S_{S} (j_{S})]

, and

Z_{S} [ψ (i_{S}), S_{S} (j_{S})]

on the uniform scales of

(ψ, S_{c})

and

(ψ, S_{S})

, and the quadratic spline functions

R_{c} (ψ, S_{c}), Z_{c} (ψ, S_{c}), R_{S} (ψ, S_{S})

, and

Z_{S} (ψ, S_{S})

can then be created as
其中

δ S_{S} = 1 / (l s s s - 1)

，lsss 是 SOL 区域中平行方向上的样条分辨率，而

1 \leq j_{S} \leq l s s s

。在和的均匀尺度上，得到值

R_{c} [ψ (i_{c}), S_{c} (j_{c})], Z_{c} [ψ (i_{c}), S_{c} (j_{c})], R_{S} [ψ (i_{S}), S_{S} (j_{S})]

和以及二次样条函数

R_{c} (ψ, S_{c}), Z_{c} (ψ, S_{c}), R_{S} (ψ, S_{S})

是很简单的，然后

Z_{S} (ψ, S_{S})

可以创建为

(ψ, S_{S})

(ψ, S_{c})

Z_{S} [ψ (i_{S}), S_{S} (j_{S})]

\begin{aligned} R_{c} (ψ, S_{c}) = & R_{c} (1, i_{c}, j_{c}) + R_{c} (2, i_{c}, j_{c}) Δ ψ + R_{c} (3, i_{c}, j_{c}) Δ ψ^{2} \\ + R_{c} (4, i_{c}, j_{c}) Δ S + R_{c} (5, i_{c}, j_{c}) Δ ψ Δ S + R_{c} (6, i_{c}, j_{c}) Δ ψ^{2} Δ S \\ + R_{c} (7, i_{c}, j_{c}) Δ S^{2} + R_{c} (8, i_{c}, j_{c}) Δ ψ Δ S^{2} + R_{c} (9, i_{c}, j_{c}) Δ ψ^{2} Δ S^{2}, \end{aligned}

and

\begin{aligned} Z_{c} (ψ, S_{c}) = & Z_{c} (1, i_{c}, j_{c}) + Z_{c} (2, i_{c}, j_{c}) Δ ψ + Z_{c} (3, i_{c}, j_{c}) Δ ψ^{2} \\ + Z_{c} (4, i_{c}, j_{c}) Δ S + Z_{c} (5, i_{c}, j_{c}) Δ ψ Δ S + Z_{c} (6, i_{c}, j_{c}) Δ ψ^{2} Δ S \\ + Z_{c} (7, i_{c}, j_{c}) Δ S^{2} + Z_{c} (8, i_{c}, j_{c}) Δ ψ Δ S^{2} + Z_{c} (9, i_{c}, j_{c}) Δ ψ^{2} Δ S^{2}, \end{aligned}

where

Δ ψ = ψ - ψ (i_{c})

and

Δ S = S_{c} - S (j_{c}) \cdot ψ (i_{c}) \leq ψ < ψ (i_{c} + 1)

, and

S (j_{c}) \leq S_{c} < S (j_{c} + 1)

.
其中

Δ ψ = ψ - ψ (i_{c})

和

Δ S = S_{c} - S (j_{c}) \cdot ψ (i_{c}) \leq ψ < ψ (i_{c} + 1)

和

S (j_{c}) \leq S_{c} < S (j_{c} + 1)

和。

\begin{aligned} R_{S} (ψ, S_{S}) = & R_{S} (1, i_{S}, j_{S}) + R_{S} (2, i_{S}, j_{S}) Δ ψ + R_{S} (3, i_{S}, j_{S}) Δ ψ^{2} \\ + R_{S} (4, i_{S}, j_{S}) Δ S + R_{S} (5, i_{S}, j_{S}) Δ ψ Δ S + R_{S} (6, i_{S}, j_{S}) Δ ψ^{2} Δ S \\ + R_{S} (7, i_{S}, j_{S}) Δ S^{2} + R_{S} (8, i_{S}, j_{S}) Δ ψ Δ S^{2} + R_{S} (9, i_{S}, j_{S}) Δ ψ^{2} Δ S^{2}, \end{aligned}

and

\begin{aligned} Z_{S} (ψ, S_{S}) = & Z_{S} (1, i_{S}, j_{S}) + Z_{S} (2, i_{S}, j_{S}) Δ ψ + Z_{S} (3, i_{S}, j_{S}) Δ ψ^{2} \\ + Z_{S} (4, i_{S}, j_{S}) Δ S + Z_{S} (5, i_{S}, j_{S}) Δ ψ Δ S + Z_{S} (6, i_{S}, j_{S}) Δ ψ^{2} Δ S \\ + Z_{S} (7, i_{S}, j_{S}) Δ S^{2} + Z_{S} (8, i_{S}, j_{S}) Δ ψ Δ S^{2} + Z_{S} (9, i_{S}, j_{S}) Δ ψ^{2} Δ S^{2}, \end{aligned}

where

Δ ψ = ψ - ψ (i_{S})

and

Δ S = S_{S} - S (j_{S}) . ψ (i_{S}) \leq ψ < ψ (i_{S} + 1)

and

S (j_{S}) \leq S_{S} < S (j_{S} + 1)

. By using Eqs. (7)-(10), we could compute the field aligned mesh in cylindrical coordinates with a given regular mesh in the core region

(ψ, S_{c})

and in the SOL region

(ψ, S_{S})

, respectively. An example of global field aligned mesh is given in Fig. 2(a). The field aligned grids are irregular in cylindrical coordinates. For the overlap part of core and SOL regions at the separatrix

(R > 0)

, the grid positions are determined by using

R_{c} (ψ, S_{c})

and

Z_{c} (ψ, S_{c})

, which are shared by both core and SOL regions with the scopes of two coordinates:

(ψ, S_{c})

and

(ψ, S_{S})

as shown by the black stars. In magnetic coordinates, the grids in SOL and core regions are regular inside each domain, respectively. The shared grids at the separatrix are designed regularly with the interior grids in the core region (keep the same parallel coordinate

S_{c}

value), which are shown in Figs. 2(b) and 2(c).
其中

Δ ψ = ψ - ψ (i_{S})

和

Δ S = S_{S} - S (j_{S}) . ψ (i_{S}) \leq ψ < ψ (i_{S} + 1)

和

S (j_{S}) \leq S_{S} < S (j_{S} + 1)

.通过使用方程。（7）-（10）中，我们可以分别计算在核心区域

(ψ, S_{c})

和SOL区域

(ψ, S_{S})

中给定规则网格的圆柱坐标中的场对齐网格。图2（a）给出了全局场对齐网格的示例。场对齐的网格在圆柱坐标中是不规则的。对于分离点

(R > 0)

处核心区域和 SOL 区域的重叠部分，网格位置由和

R_{c} (ψ, S_{c})

确定，它们由核心区域和

Z_{c} (ψ, S_{c})

SOL 区域共享，范围为两个坐标：

(ψ, S_{c})

和

(ψ, S_{S})

，如黑色星号所示。在磁坐标中，SOL 和核心区域中的网格分别在每个域内是规则的。分离点处的共享网格与核心区域的内部网格有规律地设计（保持相同的平行坐标

S_{c}

值），如图2（b）和图2（c）所示。

In GTC-X, particle dynamic equations are evolved in cylindrical coordinates

(R, ζ, Z)

to avoid the singularity at the separatrix. However, the field aligned mesh is irregular in

(R, Z)

space, which is difficult to carry out particle-grid gather-scatter operation for PIC simulation. It is noted that the field aligned grids are labeled by both cylindrical and magnetic coordinates:

(R, Z)

and

(ψ, S_{c})

in the core region and

(R, Z)

and

(ψ, S_{S})

in the SOL region. The mesh is regular in magnetic coordinates in the core and SOL regions, respectively, as shown in Fig. 2. Thus, we can transform the particle coordinate from

(R, Z)

(ψ, S_{c})

and

(ψ, S_{S})

by using Eqs. (3), (5), and (6), and then a simple linear interpolation between the particle location and the regular rectangular mesh can be used for particle-grid gather-scatter operation,
在GTC-X中，粒子动力学方程在圆柱坐标中演化，

(R, ζ, Z)

以避免分离点处的奇点。然而，场对齐网格在空间上

(R, Z)

是不规则的，难以进行粒子-网格聚集-散射操作进行PIC仿真。需要注意的是，磁场对齐的网格由圆柱坐标和磁坐标标记：

(R, Z)

(ψ, S_{c})

在核心区域和

(R, Z)

(ψ, S_{S})

SOL 区域。网格分别在核心区域和 SOL 区域的磁坐标中呈规则，如图 2 所示。因此，我们可以使用方程将粒子坐标从

(R, Z)

到和

(ψ, S_{c})

(ψ, S_{S})

通过进行转换。（3）、（5）和（6），然后使用粒子位置和规则矩形网格之间的简单线性插值进行粒子网格聚集-散射操作，

FIG. 2. (a) Global field aligned mesh in cylindrical coordinates. Field aligned mesh mapping from cylindrical coordinates to magnetic coordinates: (b) core region grids in

(ψ, S_{c})

coordinates and

(c)

SOL region grids in

(ψ, S_{S})

coordinates. The black stars represent the shared grids at the separatrix. The grids shown here are only for illustrating the algorithm, which are much sparser than the ones used in realistic simulation. where we implement

31 D

linear interpolations along

ψ

ζ

, and

S_{c} / S_{S}

for the simulations in 3 dimensional

(ψ, ζ, S_{c} / S_{S})

space.
图 2.（a）圆柱坐标中的全局场对齐网格。从圆柱坐标到磁坐标的场对齐网格映射：（b）坐标中

(ψ, S_{c})

的核心区域网格和

(c)

坐标中的

(ψ, S_{S})

SOL 区域网格。黑色星星代表分离点的共享网格。此处显示的网格仅用于说明算法，这些网格比实际模拟中使用的网格要稀疏得多。我们实现沿

ψ

、

ζ

和

S_{c} / S_{S}

的

31 D

线性插值，用于 3 维

(ψ, ζ, S_{c} / S_{S})

空间的模拟。

8. Laplacian operator 8. 拉普拉斯算子

In the global simulation of FRC geometry with the field aligned mesh, the perpendicular Laplacian operator is discretized in

(ψ, ζ, S_{c})

space for the core region and in

(ψ, ζ, S_{S})

space for the SOL region, respectively, since the irregular mesh in cylindrical coordinates becomes regular in magnetic coordinates.
在具有场对齐网格的 FRC 几何体的全局仿真中，垂直拉普拉斯算子分别在

(ψ, ζ, S_{c})

核心区域的空间和 SOL 区域的空间

(ψ, ζ, S_{S})

中离散化，因为圆柱坐标中的不规则网格在磁坐标中变得规则。

The Laplacian operator can be expanded in a generalized coordinates
拉普拉斯算子可以在广义坐标中展开

\nabla^{2} f = \frac{1}{J} \frac{\partial}{\partial ξ^{α}} (J \nabla ξ^{α} \cdot \nabla ξ^{β} \frac{\partial f}{\partial ξ^{β}})

where

f

represents an arbitrary scalar field,

ξ^{α}, ξ^{β}

refer to the three dimensional coordinates, and

J

is Jacobian and defined as

J^{- 1} ϵ^{α β γ}

= \nabla ξ^{α} \times \nabla ξ^{β} \cdot \nabla ξ^{γ}, \nabla ξ^{α} = \frac{\partial ξ^{α}}{\partial R} \hat{R} + \frac{\partial ξ^{α}}{\partial ζ} \frac{1}{R} \hat{ζ} + \frac{\partial ξ^{α}}{\partial Z} \hat{Z}

, and

ϵ^{α β γ} = 1

, and it becomes antisymmetric when the indices change.

^{39}

其中

f

表示任意标量场，

ξ^{α}, ξ^{β}

指三维坐标，是

J

雅可比坐标，定义为

J^{- 1} ϵ^{α β γ}

= \nabla ξ^{α} \times \nabla ξ^{β} \cdot \nabla ξ^{γ}, \nabla ξ^{α} = \frac{\partial ξ^{α}}{\partial R} \hat{R} + \frac{\partial ξ^{α}}{\partial ζ} \frac{1}{R} \hat{ζ} + \frac{\partial ξ^{α}}{\partial Z} \hat{Z}

和

ϵ^{α β γ} = 1

，当索引发生变化时，它变为反对称。

^{39}

In the core region with

(ψ, ζ, S_{c})

, considering

S_{c} = S_{c} [ψ (R, Z), θ (R, Z)]

, the Laplacian is written as
在核心区域中

(ψ, ζ, S_{c})

，考虑到

S_{c} = S_{c} [ψ (R, Z), θ (R, Z)]

，拉普拉斯被写成

\begin{aligned} \nabla_{c}^{2} f = & g^{ψ ψ} \frac{\partial^{2} f}{\partial ψ^{2}} + 2 g_{c}^{ψ S} \frac{\partial^{2} f}{\partial ψ \partial S_{c}} + g_{c}^{S S} \frac{\partial^{2} f}{\partial S_{c}^{2}} + g^{ζ ζ} \frac{\partial^{2} f}{\partial ζ^{2}} \\ + \frac{1}{J_{c}} [\frac{\partial (J_{c} g^{ψ ψ})}{\partial ψ} + \frac{\partial (J_{c} g_{c} ψ S)}{\partial S_{c}}] \frac{\partial f}{\partial ψ} \\ + \frac{1}{J_{c}} [\frac{\partial (J_{c} g_{c}^{ψ S})}{\partial ψ} + \frac{\partial (J_{c} g_{c}^{S S})}{\partial S_{c}}] \frac{\partial f}{\partial S_{c}}, \end{aligned}

where 哪里

\begin{aligned} J_{c}^{- 1} = \nabla ψ \times \nabla ζ \cdot \nabla S_{c} \\ = (\frac{\partial ψ}{\partial R} \hat{R} + \frac{\partial ψ}{\partial Z} \hat{Z}) \times \frac{1}{R} \hat{ζ} \cdot [\frac{\partial S_{c}}{\partial ψ} (\frac{\partial ψ}{\partial R} \hat{R} + \frac{\partial ψ}{\partial Z} \hat{Z}) \\ + \frac{\partial S_{c}}{\partial θ} (\frac{\partial θ}{\partial R} \hat{R} + \frac{\partial θ}{\partial Z} \hat{Z})] \\ = \frac{1}{R} \frac{\partial ψ}{\partial R} \frac{\partial S_{c}}{\partial θ} \frac{\partial θ}{\partial Z} - \frac{1}{R} \frac{\partial ψ}{\partial Z} \frac{\partial S_{c}}{\partial θ} \frac{\partial θ}{\partial R}, \\ g^{ψ ψ} = \nabla ψ \cdot \nabla ψ = {(\frac{\partial ψ}{\partial R})}^{2} + {(\frac{\partial ψ}{\partial Z})}^{2}, \\ g_{c}^{S S} = \nabla S_{c} \cdot \nabla S_{c} = {(\frac{\partial S_{c}}{\partial ψ} \frac{\partial ψ}{\partial R} + \frac{\partial S_{c}}{\partial θ} \frac{\partial θ}{\partial R})}^{2} + {(\frac{\partial S_{c}}{\partial ψ} \frac{\partial ψ}{\partial Z} + \frac{\partial S_{c}}{\partial θ} \frac{\partial θ}{\partial Z})}^{2}, \\ g_{c}^{ψ S} = \nabla ψ \cdot \nabla S_{c} \\ = \frac{\partial S_{c}}{\partial ψ} [{(\frac{\partial ψ}{\partial R})}^{2} + {(\frac{\partial ψ}{\partial Z})}^{2}] + \frac{\partial S_{c}}{\partial θ} [\frac{\partial ψ}{\partial R} \frac{\partial θ}{\partial R} + \frac{\partial ψ}{\partial Z} \frac{\partial θ}{\partial Z}], \\ g^{ζ ζ} = \nabla ζ \cdot \nabla ζ = \frac{1}{R^{2}} . \end{aligned}

In the SOL region with

(ψ, ζ, S_{S})

, considering

S_{S} = S_{S} [ψ (R

Z), Z]

, the Laplacian is written as
在 SOL 区域中

(ψ, ζ, S_{S})

，考虑到

S_{S} = S_{S} [ψ (R

，

Z), Z]

拉普拉斯被写成

\begin{aligned} \nabla_{S}^{2} f = & g^{ψ ψ} \frac{\partial^{2} f}{\partial ψ^{2}} + 2 g_{S}^{ψ S} \frac{\partial^{2} f}{\partial ψ \partial S_{S}} + g_{S}^{S S} \frac{\partial^{2} f}{\partial S_{S}^{2}} + g^{ζ ζ} \frac{\partial^{2} f}{\partial ζ^{2}} \\ + \frac{1}{J_{S}} [\frac{\partial (J_{S} g^{ψ ψ})}{\partial ψ} + \frac{\partial (J_{S} g_{S}^{ψ S})}{\partial S_{S}}] \frac{\partial f}{\partial ψ} \\ + \frac{1}{J_{S}} [\frac{\partial (J_{S} g_{S}^{ψ S})}{\partial ψ} + \frac{\partial (J_{S} g_{S}^{S S})}{\partial S_{S}}] \frac{\partial f}{\partial S_{S}}, \end{aligned}

where 哪里

\begin{aligned} J_{S}^{- 1} = & \nabla ψ \times \nabla ζ \cdot \nabla S_{S} \\ = & (\frac{\partial ψ}{\partial R} \hat{R} + \frac{\partial ψ}{\partial Z} \hat{Z}) \times \frac{1}{R} \hat{ζ} \cdot [\frac{\partial S_{S}}{\partial ψ} (\frac{\partial ψ}{\partial R} \hat{R} + \frac{\partial ψ}{\partial Z} \hat{Z}) + \frac{\partial S_{S}}{\partial Z} \hat{Z}] \\ = & \frac{1}{R} \frac{\partial ψ}{\partial R} \frac{\partial S_{S}}{\partial Z}, \\ g_{S}^{S S} = \nabla S_{S} \cdot \nabla S_{S} = {(\frac{\partial S_{S}}{\partial ψ} \frac{\partial ψ}{\partial R})}^{2} + {(\frac{\partial S_{S}}{\partial ψ} \frac{\partial ψ}{\partial Z} + \frac{\partial S_{S}}{\partial Z})}^{2}, \\ g_{S}^{ψ S} & = \nabla ψ \cdot \nabla S_{S} = \frac{\partial S_{S}}{\partial ψ} [{(\frac{\partial ψ}{\partial R})}^{2} + {(\frac{\partial ψ}{\partial Z})}^{2}] + \frac{\partial ψ}{\partial Z} \frac{\partial S_{S}}{\partial Z} . \end{aligned}

For the shared grids between core and SOL regions at the separatrix, the perpendicular Laplacian operator is expanded as
对于分离点处核心区域和 SOL 区域之间的共享网格，垂直拉普拉斯算子扩展为

\nabla_{⊥, X}^{2} f = {g^{ψ ψ ψ} \frac{\partial^{2} f}{\partial ψ^{2}} |}_{b_{0} \times ζ} + {\frac{1}{J_{X}} \frac{\partial J_{X} g^{ψ ψ}}{\partial ψ} \frac{\partial f}{\partial ψ} |}_{b_{0} \times ζ} + g^{ζ ζ} \frac{\partial^{2} f}{\partial ζ}

where

|_{b_{0} \times ζ ζ}

represents the partial derivative with respect to

ψ

among the orthogonal grids along the

b_{0} \times ζ

direction, and
其中

|_{b_{0} \times ζ ζ}

表示相对于沿

b_{0} \times ζ

方向的正交网格

ψ

之间的偏导数，以及

\begin{aligned} J_{X}^{- 1} & = \nabla ψ \times \nabla ζ \cdot \nabla S_{X} \\ = (\frac{\partial ψ}{\partial R} \hat{R} + \frac{\partial ψ}{\partial Z} \hat{Z}) \times \frac{1}{R} \hat{ζ} \cdot [\frac{\partial S_{X}}{\partial R} \hat{R} + \frac{\partial S_{X}}{\partial Z} \hat{Z}] \\ = \frac{1}{R} \frac{\partial ψ}{\partial R} \frac{\partial S_{X}}{\partial Z} - \frac{1}{R} \frac{\partial ψ}{\partial Z} \frac{\partial S_{X}}{\partial R} . \end{aligned}

Taking advantage of the toroidal symmetry of FRC, we could transform the Laplacian operator into Fourier space with respect to the toroidal angle

ζ

, which can avoid solving the 3 dimensional matrix. In this paper, we can also simplify Eqs. (12)-(14) assuming

k_{⊥} ≫ k

and

k_{∥} ≪ k_{⊥}

, where

κ = \nabla B_{0} / B_{0}

. Thus, for each toroidal mode

n

(\partial / \partial ζ = i n ζ)

, the Laplacian operators in the core region, SOL region, and at the separatrix can be written as
利用FRC的环形对称性，我们可以将拉普拉斯算子变换为相对于环形角

ζ

的傅里叶空间，从而避免求解三维矩阵。在本文中，我们还可以简化方程。（12）-（14）假设

k_{⊥} ≫ k

和

k_{∥} ≪ k_{⊥}

，其中

κ = \nabla B_{0} / B_{0}

.因此，对于每个环形模

n

(\partial / \partial ζ = i n ζ)

态，核心区域、SOL 区域和分离点处的拉普拉斯算子可以写成

\begin{matrix} \nabla_{⊥ c}^{2} f_{n} \approx g^{ψ ψ} \frac{\partial^{2} f_{n}}{\partial ψ^{2}} + 2 g_{c}^{ψ S} \frac{\partial^{2} f_{n}}{\partial ψ \partial S_{c}} + g_{c}^{S S} \frac{\partial^{2} f_{n}}{\partial S_{c}^{2}} - n^{2} g^{ζ ζ} f_{n}, \\ \nabla_{⊥ S}^{2} f_{n} \approx g^{ψ ψ} \frac{\partial^{2} f_{n}}{\partial ψ^{2}} + 2 g_{S}^{ψ S} \frac{\partial^{2} f_{n}}{\partial ψ \partial S_{S}} + g_{S}^{S S} \frac{\partial^{2} f_{n}}{\partial S_{S}^{2}} - n^{2} g^{ζ ζ} f_{n}, \\ {\nabla_{⊥ X f_{n}}^{2} \approx g^{ψ ψ ψ} \frac{\partial^{2} f_{n}}{\partial ψ^{2}} |}_{b_{0} \times ζ} - n^{2} g^{ζ ζ} f_{n}, \end{matrix}

where

f_{n}

is the

n

toroidal mode component of

f

.
其中

f_{n}

是的

n

f

环形模态分量。

The grids shown in Figs. 2 and 3 are only for illustrating the field aligned grid algorithm in cylindrical coordinates clearly, which are much sparser than the realistic simulation.
图 2 和图 3 中所示的网格仅用于清楚地说明圆柱坐标中的场对齐网格算法，比实际模拟要稀疏得多。

FIG. 3. The green + symbols denote the orthogonal grids in the SOL and core regions with respect to the shared grids at the separatrix, and the dotted line in magenta is along the perpendicular direction

b_{0} \times \hat{ζ}

from the shared grids as shown by the black stars. The grids shown here are only for illustrating the algorithm, which are much sparser than the ones used in realistic simulation.
图 3.绿色 + 符号表示 SOL 和核心区域中相对于分离点处共享网格的正交网格，洋红色的虚线沿共享网格的垂直方向

b_{0} \times \hat{ζ}

，如黑色星形所示。此处显示的网格仅用于说明算法，这些网格比实际模拟中使用的网格要稀疏得多。

9. PHYSICS MODEL 9. 物理模型

10. A. Formulation 10. A. 配方

We use the electrostatic Vlasov-Poisson system for physics simulation in this paper. Particle dynamics is described by the gyrokinetic equation using gyrocenter position

R

, magnetic moment

μ

, and parallel velocity

v_{∥}

as independent variables in five dimensional phase space
本文使用静电Vlasov-Poisson系统进行物理模拟。粒子动力学由陀螺动力学方程描述，使用陀螺中心位置

R

、磁矩

μ

和平行速度

v_{∥}

作为五维相空间中的自变量

\begin{matrix} (\frac{\partial}{\partial t} + \dot{R} \cdot \nabla + {\dot{v}}_{∥} \frac{\partial}{\partial v_{∥}}) f_{α} (R, v_{∥}, μ, t) = 0, \\ \dot{R} = v_{∥} b + v_{E} + v_{d}, \\ \dot{v_{∥}} = - \frac{1}{m_{α}} \frac{B^{*}}{B_{∥}^{*}} \cdot (Z_{α} \nabla ⟨ δ ϕ ⟩ + μ \nabla B), \end{matrix}

where

Z_{α}, m_{α}

, and

f_{α}

are the charge, mass, and distribution function of

α

species.

B

is the equilibrium magnetic field,

b = B / B, B^{*} = B

+ (B v_{∥} / Ω_{c α}) \nabla \times b

, and

B_{∥}^{*} = b \cdot B^{*} . δ ϕ

is the electrostatic potential.

⟨ \dots ⟩ = (1 / 2 π) \int d x d ξ (\dots) δ (R + ρ_{α} - x)

represents the gyrophase average,

ξ

is the gyrophase angle,

x

represents the particle position,

ρ_{α} = b \times v_{⊥} / Ω_{c α}

is the particle gyroradius, and

Ω_{c α}

is the particle cyclotron frequency.

v_{E}

is the

E \times B

velocity, and

v_{d}

is the magnetic drift velocity, which are given as
其中

Z_{α}, m_{α}

，和

f_{α}

是物种的

α

电荷、质量和分布函数。

B

是平衡磁场，

b = B / B, B^{*} = B

+ (B v_{∥} / Ω_{c α}) \nabla \times b

B_{∥}^{*} = b \cdot B^{*} . δ ϕ

是静电势。

⟨ \dots ⟩ = (1 / 2 π) \int d x d ξ (\dots) δ (R + ρ_{α} - x)

表示陀螺相位平均值，

ξ

表示陀螺相位角，

x

表示粒子位置，

ρ_{α} = b \times v_{⊥} / Ω_{c α}

表示粒子回转半径，

Ω_{c α}

表示粒子回旋加速器频率。

v_{E}

是

E \times B

速度，

v_{d}

是磁漂移速度，其表示为

\begin{matrix} v_{E} = \frac{c b \times \nabla ⟨ δ ϕ ⟩}{B_{∥}^{*}}, \\ v_{d} = \frac{c m_{α} v_{∥}^{2}}{Z_{α} B_{∥}^{*}} b \times (b \cdot \nabla b) + \frac{c μ}{Z_{α} B_{∥}^{*}} b \times \nabla B . \end{matrix}

In order to minimize the particle noise, the perturbative

δ f

simulation method

^{16, 17}

is applied. The particle distribution is decomposed into equilibrium and perturbed parts as

f_{α} = f_{α 0} (R, μ, v_{∥}) + δ f_{α} (R, μ, v_{∥}, t)

, and the equilibrium

f_{α 0}

satisfies the following equation:
为了最小化粒子噪声，采用扰动

δ f

模拟方法

^{16, 17}

。粒子分布分解为平衡部分和扰动部分

f_{α} = f_{α 0} (R, μ, v_{∥}) + δ f_{α} (R, μ, v_{∥}, t)

，均衡

f_{α 0}

满足以下方程：

L_{0} f_{α 0} = 0

where

L_{0} = \partial / \partial t + (v_{∥} b + v_{d}) \cdot \nabla - (μ / m_{α}) B^{*} \cdot \nabla B / B_{∥}^{*} (\partial / \partial v_{∥})

. Because

v_{d}

is only in the

\hat{ζ}

direction, the particle drift orbit width is zero in FRC geometry, then

f_{α 0} = n_{α 0} {(\frac{m_{α}}{2 π T_{α 0}})}^{1.5} \exp [- \frac{m_{α} v_{∥}^{2} + 2 μ B}{2 T_{α 0}}]

is the exact solution of Eq. (21), and

n_{α 0} (ψ)

and

T_{α 0} (ψ)

are the

1 D

function of magnetic flux surface. Subtracting Eq. (18) by Eq. (21), we have the equation for perturbed distribution

δ f_{α}

其中

L_{0} = \partial / \partial t + (v_{∥} b + v_{d}) \cdot \nabla - (μ / m_{α}) B^{*} \cdot \nabla B / B_{∥}^{*} (\partial / \partial v_{∥})

.因为

v_{d}

只是在方向上

\hat{ζ}

，粒子漂移轨道宽度在FRC几何中为零，所以

f_{α 0} = n_{α 0} {(\frac{m_{α}}{2 π T_{α 0}})}^{1.5} \exp [- \frac{m_{α} v_{∥}^{2} + 2 μ B}{2 T_{α 0}}]

是方程（21）的精确解，和

n_{α 0} (ψ)

是

T_{α 0} (ψ)

磁通量面的

1 D

函数。用方程（21）减去方程（18），我们得到扰动分布

δ f_{α}

方程

L δ f_{α} = - δ L f_{α 0},

where

δ L = v_{E} \cdot \nabla - (Z_{α} / m_{α} / B_{∥}^{*}) B^{*} \cdot \nabla ⟨ δ ϕ ⟩ (\partial / \partial v_{∥})

and

L = L_{0}

+ δ L

. Defining particle weight as

w_{α} = δ f_{α} / f_{α}

, we can derive the weight equation from Eq. (22) as
其中

δ L = v_{E} \cdot \nabla - (Z_{α} / m_{α} / B_{∥}^{*}) B^{*} \cdot \nabla ⟨ δ ϕ ⟩ (\partial / \partial v_{∥})

和

L = L_{0}

+ δ L

.将粒子重量定义为

w_{α} = δ f_{α} / f_{α}

，我们可以从方程（22）中推导出重量方程为

\begin{aligned} \frac{d w_{α}}{d t} = & (1 - w_{α}) [- {v_{E} \cdot \frac{\nabla f_{α 0}}{f_{α 0}} |}_{v_{⊥}} - \frac{Z_{α} v_{∥} b \cdot \nabla ⟨ δ ϕ ⟩}{T_{α 0}} \\ - \frac{Z_{α}}{T_{α 0}} (\frac{μ b \times \nabla B_{0}}{m_{α} Ω_{c α}} + \frac{v_{∥}^{2}}{Ω_{c α}} b \times (b \cdot \nabla b)) \cdot \nabla ⟨ δ ϕ ⟩], \end{aligned}

where we have used the chain rule

{\nabla f_{α 0} |}_{v_{⊥}} = {\nabla f_{α 0} |}_{μ} + (μ f_{α 0} \nabla B / T_{α 0})

in the derivation of Eq. (23) from Eq. (22).
其中，我们在从方程（22）推导方程（23）时使用了链式法则

{\nabla f_{α 0} |}_{v_{⊥}} = {\nabla f_{α 0} |}_{μ} + (μ f_{α 0} \nabla B / T_{α 0})

。

The gyrokinetic Vlasov equation is used for ion species, and its perturbed density is
陀螺动力学弗拉索夫方程用于离子种类，其扰动密度为

⟨ δ n_{i} (x, t) ⟩ = \int δ f_{i} (R, μ, v_{∥}, t) d v d R d ξ δ (R + ρ_{i} - x) / (2 π)

Electron dynamics is assumed as adiabatic for simplicity, and the elec tron perturbed density is
为简单起见，假设电子动力学为绝热，电子加速器扰动密度为

δ n_{e} (x, t) = \frac{e δ ϕ}{T_{e 0}} n_{e 0} .

The gyrokinetic Poisson's equation is

^{14}

陀螺动力学泊松方程为

^{14}

\frac{Z_{i}^{2} n_{i 0}}{T_{i 0}} (δ ϕ - \tilde{δ ϕ}) = Z_{i} ⟨ δ n_{i} ⟩ - e δ n_{e}

where

n_{i 0}

and

T_{i 0}

are the equilibrium ion density and temperature.

\tilde{δ ϕ}

is the double gyrophase average of electrostatic potential for ion species, which is given as
其中

n_{i 0}

和

T_{i 0}

是平衡离子密度和温度。

\tilde{δ ϕ}

是离子种类静电势的双陀螺相平均值，其表示为

\tilde{δ ϕ} (x, t) = \frac{1}{n_{i 0}} \int f_{i 0} (R, μ, v_{∥}) ⟨ δ ϕ ⟩ d v d R d ξ δ (R + ρ_{i} - x) / (2 π),

where

\int d v = \frac{2 π B}{m_{i}} \int d v_{∥} d μ

. 其中

\int d v = \frac{2 π B}{m_{i}} \int d v_{∥} d μ

11. B. Implementation of dynamic equations in cylindrical coordinates
11. B. 圆柱坐标中动力学方程的实现

The gyrocenter equation of motion is
陀螺中心运动方程为

\dot{R} = \underset{{parallel motion}}{\underset{⏟}{v_{∥} b}} + \underset{{curvaturedrift}}{\underset{⏟}{\frac{v_{∥}^{2}}{Ω_{c α}} \nabla \times b}} + \underset{{grad - B drift}}{\underset{⏟}{\frac{μ b \times \nabla B}{m_{α} Ω_{c α}}}} \underset{{E \times B drift}}{\underset{⏟}{+ \frac{c b \times \nabla ⟨ δ ϕ ⟩}{B}}},

with 跟

\begin{aligned} \underset{{parallel motion}}{\underset{⏟}{v_{∥} b}} = & v_{∥} (b_{R} \hat{R} + b_{Z} \hat{Z}), \\ \underset{{curvature drift}}{\underset{⏟}{\frac{v_{∥}^{2}}{Ω_{c α}} \nabla \times b}} = & [\frac{1}{m_{α} Ω_{c α}} \frac{m_{α} v_{∥}^{2}}{B} (b_{Z} \frac{\partial B}{\partial R} - b_{R} \frac{\partial B}{\partial Z}) \\ + \frac{v_{∥}^{2}}{Ω_{c α} B} (\frac{\partial B_{R}}{\partial Z} - \frac{\partial B_{Z}}{\partial R})] \hat{ζ}, \end{aligned}

\begin{aligned} \underset{{grad - B drift}}{\underset{⏟}{\frac{μ b \times \nabla B}{m_{α} Ω_{c α}}}} = & \frac{μ}{m_{α} Ω_{c α}} (b_{Z} \frac{\partial B}{\partial R} - b_{R} \frac{\partial B}{\partial Z}) \hat{ζ}, \\ \underset{E \times B drift}{\underset{⏟}{\frac{c b \times \nabla ⟨ δ ϕ ⟩}{B}}} = & \frac{c}{B} [- \frac{1}{R} b_{Z} \frac{\partial ⟨ δ ϕ ⟩}{\partial ζ} \hat{R} + (b_{Z} \frac{\partial ⟨ δ ϕ ⟩}{\partial R} - b_{R} \frac{\partial ⟨ δ ϕ ⟩}{\partial Z}) \hat{ζ} \\ + \frac{1}{R} b_{R} \frac{\partial ⟨ δ ϕ ⟩}{\partial ζ} \hat{Z}], \end{aligned}

and

\begin{aligned} {\dot{v}}_{∥} & = \underset{{I}}{\underset{⏟}{- \frac{1}{m_{α}} \frac{1}{B} (B + \frac{B v_{∥}}{Ω_{c α}} \nabla \times b) \cdot (Z_{α} \nabla ⟨ δ ϕ ⟩ + μ \nabla B)}} \\ = \underset{{I}}{\underset{⏟}{- \frac{1}{m_{α}} b \cdot (Z_{α} \nabla ⟨ δ ϕ ⟩ + μ \nabla B)}} \underset{{I I}}{\underset{⏟}{- \frac{1}{m_{α}} \frac{v_{∥}}{Ω_{c α}} \nabla \times b \cdot (Z_{α} \nabla ⟨ δ ϕ ⟩ + μ \nabla B)}} \\ = \underset{{I I}}{\underset{⏟}{- \frac{Z_{α}}{m_{α}} \frac{\partial ⟨ δ ϕ ⟩}{\partial S} - \frac{μ}{m_{α}} (b_{R} \frac{\partial B}{\partial R} + b_{Z} \frac{\partial B}{\partial Z})}} \\ \underset{⏟}{- \frac{Z_{α}}{m_{α}} \frac{v_{∥}}{Ω_{c α}} \frac{1}{R} \frac{\partial ⟨ δ ϕ ⟩}{\partial ζ} [(b_{Z} \frac{\partial B}{\partial R} - b_{R} \frac{\partial B}{\partial Z}) + (\frac{\partial B_{R}}{\partial Z} - \frac{\partial B_{Z}}{\partial R})]} . \end{aligned}

The weight equation is
权重方程为

where 哪里

\begin{aligned} \underset{wdrive}{\underset{⏟}{- {v_{E} \cdot \frac{\nabla f_{i 0}}{f_{i 0}} |}_{v_{⊥}}}} = & \frac{c}{B_{0}} b_{Z} \frac{1}{R} \frac{\partial ⟨ δ ϕ ⟩}{\partial ζ} [\frac{1}{n_{i 0}} \frac{\partial n_{i 0}}{\partial R} + (\frac{ϵ}{T_{i 0}} - 1.5) \frac{1}{T_{i 0}} \frac{\partial T_{i 0}}{\partial R}] \\ - \frac{c}{B_{0}} b_{R} \frac{1}{R} \frac{\partial ⟨ δ ϕ ⟩}{\partial ζ} [\frac{1}{n_{i 0}} \frac{\partial n_{i 0}}{\partial Z} + (\frac{ϵ}{T_{i 0}} - 1.5) \frac{1}{T_{i 0}} \frac{\partial T_{i 0}}{\partial Z}], \end{aligned}

with

ϵ = 0.5 m_{i} v_{∥}^{2} + μ B

, 和

ϵ = 0.5 m_{i} v_{∥}^{2} + μ B

，

\underset{wpara}{\underset{⏟}{\frac{Z_{i}}{T_{i 0}} v_{∥} E_{∥}}} = - \frac{Z_{i}}{T_{i 0}} v_{∥} b \cdot \nabla ⟨ δ ϕ ⟩ = - \frac{Z_{i} v_{∥}}{T_{i 0}} \frac{\partial ⟨ δ ϕ ⟩}{\partial S},

with

S = S_{c} \times L

S = S_{S} \times L

, and

L

is the field line length, and
with

S = S_{c} \times L

或

S = S_{S} \times L

，是

L

场线长度，并且

\begin{aligned} \underset{wdrift}{\underset{⏟}{\frac{Z_{i}}{T_{i 0}} v_{d} \cdot E_{⊥}}} = & - \frac{Z_{i}}{T_{i 0}} (\frac{μ b \times \nabla B}{m_{i} Ω_{c i}} + \frac{v_{∥}^{2}}{Ω_{c α}} \nabla \times b) \cdot \nabla ⟨ δ ϕ ⟩ \\ = & - \frac{Z_{i}}{T_{i 0}} \frac{1}{R} \frac{1}{B} \frac{\partial ⟨ δ ϕ ⟩}{\partial ζ} [\frac{c}{Z_{i}} (μ + \frac{m_{i} v_{∥}^{2}}{B}) (b_{Z} \frac{\partial B}{\partial R} - b_{R} \frac{\partial B}{\partial Z}) \\ + \frac{v_{∥}^{2}}{Ω_{c i}} (\frac{\partial B_{R}}{\partial Z} - \frac{\partial B_{Z}}{\partial R})] . \end{aligned}

In GTC-X, the gyro-average is performed analytically by multiplying the fields with Bessel function as

⟨ δ ϕ ⟩ = δ ϕ J_{0} (k_{ζ} ρ_{i})

and

⟨ δ n_{i} ⟩ = δ n_{i} J_{0} (k_{ζ} ρ_{i})

, where

ρ_{i} = v_{t h, i} / Ω_{c i}

is the ion gyroradius,

v_{t h, i} = \sqrt{T_{i} / m_{i}}, Ω_{c i} = Z_{i} B / c m_{i}, k_{ζ} = n / R

is the toroidal wave vector,

n

is the toroidal mode number, and

R

is the radial position, which is only valid for single toroidal mode simulation. The radial component of perpendicular wave vector

k_{r}

is not considered for the gyroaverage in the current implementation for simplicity. In the future work, a more realistic gyro-average using the 4 point average method

^{21}

will be implemented in GTC-X for multiple mode simulation, which is used in another FRC drift wave code ANC.
在 GTC-X 中，陀螺平均值是通过将贝塞尔函数为

⟨ δ ϕ ⟩ = δ ϕ J_{0} (k_{ζ} ρ_{i})

和

⟨ δ n_{i} ⟩ = δ n_{i} J_{0} (k_{ζ} ρ_{i})

的场相乘来执行的，其中

ρ_{i} = v_{t h, i} / Ω_{c i}

是离子陀螺半径，

v_{t h, i} = \sqrt{T_{i} / m_{i}}, Ω_{c i} = Z_{i} B / c m_{i}, k_{ζ} = n / R

是环形波矢量，

n

是环形模态数，

R

是径向位置，仅对单环形模态模拟有效。为简单起见，在当前实现中，陀螺平均值不考虑垂直波矢量

k_{r}

的径向分量。在未来的工作中，将在 GTC-X 中使用 4 点平均法

^{21}

实现更真实的陀螺仪平均值，用于多模仿真，该仿真用于另一个 FRC 漂移波代码 ANC。

The comparisons between gyrocenter and fully kinetic particle trajectories with different particle energies and locations are shown in Fig. 4. Both gyrocenter and fully kinetic particle pitch angles at the outer midplane are

\tan θ = v_{∥} / v_{⊥} = 1

, which are trapped between two mirror throats. It is seen that gyrokinetic description is suitable for SOL region simulation with low temperature and high magnetic field, and its fidelity decreases with increasing the particle energy and moving toward the core region. Recently, a theoretical work has illustrated that the fidelity of the gyrokinetic equation is well achieved in the nonuniform magnetic field with

ϵ = ρ / L < 1,^{40}

where

ρ

is the gyroradius and

L

is the magnetic field gradient scale length. With regard to FRC geometry, the gyrokinetic description of thermal ion can capture the gyro-orbit but not betatron and figure- 8 orbit.

^{41}

However, a strong finite Larmor radius effect stabilizes the drift wave instability near or inside the core region, which has been observed in experiment

^{10}

and by previous local simulations.

^{25 - 27}

Thus, in this paper, we focus on gyrokinetic simulation of ITG instability in FRC as the first step to demonstrate the code capability for GTC-X as well as reveal the global nature of ITG instability in the SOL region.
具有不同粒子能量和位置的陀螺中心和全动力学粒子轨迹之间的比较如图4所示。外中平面的陀螺中心和全动力学粒子俯仰角都被困在

\tan θ = v_{∥} / v_{⊥} = 1

两个镜喉之间。可以看出，陀螺动力学描述适用于低温高磁场的SOL区域模拟，其保真度随着粒子能量的增加和向核心区域移动而降低。最近，一项理论工作表明，陀螺动力学方程的保真度在非均匀磁场中得到了很好的实现，

ϵ = ρ / L < 1,^{40}

其中

ρ

陀螺半径和

L

磁场梯度尺度长度。关于FRC几何形状，热离子的陀螺动力学描述可以捕获陀螺轨道，但不能捕获βtron和8字形轨道。

^{41}

然而，强烈的有限拉莫尔半径效应稳定了核心区域附近或内部的漂移波不稳定性，这在实验

^{10}

和先前的局部模拟中已经观察到。

^{25 - 27}

因此，本文重点研究了FRC中ITG不稳定性的陀螺动力学仿真，作为展示GTC-X代码能力的第一步，并揭示了SOL区域ITG不稳定性的全局性质。

12. Poisson solver 12. 泊松求解器

Poisson's equation is solved in a semispectral form. Applying Padé approximation,

\hat{δ ϕ} \approx \frac{1}{1 + k_{⊥}^{2} ρ_{i}^{2}} δ ϕ

, Eq. (26) can be written as
泊松方程以半光谱形式求解。应用帕德近似，

\hat{δ ϕ} \approx \frac{1}{1 + k_{⊥}^{2} ρ_{i}^{2}} δ ϕ

方程（26）可以写成

\begin{matrix} \frac{e^{2} n_{i 0}}{T_{e 0}} [- \frac{Z_{i}^{2}}{e^{2}} ρ_{s}^{2} (1 + \frac{e^{2}}{Z_{i}^{2}} \frac{n_{e 0}}{n_{i 0}} \frac{T_{i}}{T_{e}}) \nabla_{⊥}^{2} + \frac{n_{e 0}}{n_{i 0}}] δ ψ \\ = Z_{i} δ n_{i} \frac{1}{1 + \frac{e^{2}}{Z_{i}^{2}} \frac{n_{e 0}}{n_{i 0}} \frac{T_{i 0}}{T_{e 0}}}, \end{matrix}

FIG. 4. Comparison of the trajectory between the fully kinetic particle and the gyrocenter for

44.4 e V

and

1022 e V

deuterium cases. The black solid line represents the separatrix, and green solid lines are the contour of magnetic flux. where

δ ψ = δ ϕ - \frac{T_{i 0}}{Z_{i} n_{i 0}} δ n_{i} / (1 + \frac{e^{2}}{Z_{i}^{2}} \frac{n_{e 0}}{n_{i 0}} \frac{T_{i 0}}{T_{e 0}}), ρ_{s} = C_{s} / Ω_{c i}

, and

C_{s}

= \sqrt{T_{e} / m_{i}}

is the ion sound speed. By applying the Fourier transform in the toroidal direction to Eq. (31), we have
图 4.完全动力学粒子与氘

44.4 e V

1022 e V

形体的陀螺中心之间的轨迹比较。黑色实线表示分离线，绿色实线表示磁通量的轮廓。其中

δ ψ = δ ϕ - \frac{T_{i 0}}{Z_{i} n_{i 0}} δ n_{i} / (1 + \frac{e^{2}}{Z_{i}^{2}} \frac{n_{e 0}}{n_{i 0}} \frac{T_{i 0}}{T_{e 0}}), ρ_{s} = C_{s} / Ω_{c i}

，和

C_{s}

= \sqrt{T_{e} / m_{i}}

是离子声速。通过将环形方向的傅里叶变换应用于方程（31），我们有

\begin{matrix} \frac{e^{2} n_{i 0}}{T_{e 0}} [- \frac{Z_{i}^{2}}{e^{2}} ρ_{s}^{2} (1 + \frac{e^{2}}{Z_{i}^{2}} \frac{n_{e 0}}{n_{i 0}} \frac{T_{i}}{T_{e}}) \nabla_{⊥}^{2} (n) + \frac{n_{e 0}}{n_{i 0}}] δ ψ (n) \\ = Z_{i} δ n_{i} (n) \frac{1}{1 + \frac{e^{2}}{Z_{i}^{2}} \frac{n_{e 0}}{n_{i 0}} \frac{T_{i 0}}{T_{e 0}}} . \end{matrix}

For each toroidally spectral component with mode number

n

, we carry out central finite difference of

\nabla_{⊥}^{2} (n)

as shown by Eqs. (15)-(17) on the

(R, Z)

plane and construct the tridiagonal matrix on field aligned simulation grids. Then, the sparse discrete matrix equation for Eq. (32) is solved using the Krylov method implemented in PETSc software.
对于每个具有模态数

n

的环形谱分量，我们进行中心有限差分，

\nabla_{⊥}^{2} (n)

如方程所示。（15）-（17）在场

(R, Z)

对齐的仿真网格上构造三对角矩阵。然后，使用PETSc软件中实现的Krylov方法求解方程（32）的稀疏离散矩阵方程。

13. CODE VERIFICATION 13. 代码验证

14. A. Slab limit 14. A. 板坯限制

In this section, we show the benchmark results of GTC-X against the analytical dispersion relation in an approximately uniform magnetic field (ignore the ion diamagnetic and curvature drifts). Applying Fourier transform

\partial_{t} = - i ω, b_{0} \cdot \nabla = i k_{∥}

, and

\nabla_{⊥} = i k_{⊥}

for Eq. (22) in linear plasmas in a uniform magnetic field, we can derive the linear perturbed distribution of gyrokinetic ion species
在本节中，我们展示了GTC-X在近似均匀的磁场中与解析色散关系的基准结果（忽略离子抗磁性和曲率漂移）。应用傅里叶变换

\partial_{t} = - i ω, b_{0} \cdot \nabla = i k_{∥}

，

\nabla_{⊥} = i k_{⊥}

对于均匀磁场中线性等离子体中的方程（22），我们可以推导出陀螺动力学离子种类的线性扰动分布

\begin{aligned} δ f_{i} = & \frac{1}{ω - k_{∥} v_{∥}} {- \frac{Z_{i}}{T_{i}} ω_{i}^{*} [1 + (\frac{m_{i} v_{∥}^{2} + 2 μ B_{0}}{2 T_{i}} - 1.5) η_{i}] \\ + \frac{Z_{i} k_{∥} v_{∥}}{T_{i}}} ⟨ δ ϕ ⟩ f_{i 0}, \end{aligned}

where

⟨ δ ϕ ⟩ = \frac{1}{2 π} ∮ δ ϕ (x) d ζ = \sum_{k} δ ϕ (k) \exp (i k \cdot R) J_{0} (\frac{k_{⊥} v_{⊥}}{Ω_{c i}}), x

is the particle position, and

R

is the gyrocenter position,

ω_{i}^{*} = \frac{c T_{i}}{Z_{i} B_{0}}

\frac{1}{n_{i 0}} \frac{\partial n_{i 0}}{\partial r} k_{ζ}

is the ion diamagnetic drift frequency (negative value for the normal profile which decreases along the

\hat{R}

direction),

k_{ζ} = n / R

is the perpendicular wave vector, and

η_{i} = d \ln T_{i} / d l n n_{i 0}

. Integrating Eq. (33) to the 0th order, the ion density perturbation becomes
其中

⟨ δ ϕ ⟩ = \frac{1}{2 π} ∮ δ ϕ (x) d ζ = \sum_{k} δ ϕ (k) \exp (i k \cdot R) J_{0} (\frac{k_{⊥} v_{⊥}}{Ω_{c i}}), x

是粒子位置，

R

是陀螺中心位置，

ω_{i}^{*} = \frac{c T_{i}}{Z_{i} B_{0}}

\frac{1}{n_{i 0}} \frac{\partial n_{i 0}}{\partial r} k_{ζ}

是离子抗磁漂移频率（沿

\hat{R}

方向减小的法线剖面为负值），

k_{ζ} = n / R

是垂直波矢量，和

η_{i} = d \ln T_{i} / d l n n_{i 0}

。将方程（33）积分到0阶，离子密度扰动变为

\begin{aligned} ⟨ δ n_{i} ⟩ = & {- \frac{Z_{i} n_{i 0}}{T_{i}} - \frac{Z_{i} n_{i 0}}{T_{i}} \frac{ω - ω_{i}^{*}}{ω} ξ_{i} Z (ξ_{i}) \\ + \frac{Z_{i} n_{i 0}}{T_{i}} η_{i} \frac{ω_{i}^{*}}{ω} ξ_{i}^{2} [1 + ξ_{i} Z (ξ_{i})] - \frac{3}{2} \frac{Z_{i} n_{i 0}}{T_{i}} η_{i} \frac{ω_{i}^{*}}{ω} ξ_{i} Z (ξ_{i})} \\ \times Γ_{0} (k_{⊥}^{2} ρ_{i}^{2}) δ ϕ + \frac{Z_{i} n_{i 0}}{T_{i}} η_{i} \frac{ω_{i}^{*}}{ω} ξ_{i} Z (ξ_{i}) [(1 - k_{⊥}^{2} ρ_{i}^{2}) \\ \times Γ_{0} (k_{⊥}^{2} ρ_{i}^{2}) + k_{⊥}^{2} ρ_{i}^{2} Γ_{1} (k_{⊥}^{2} ρ_{i}^{2})] δ ϕ \end{aligned}

where

Γ_{0} (k_{⊥}^{2} ρ_{i}^{2}) = I_{0} (k_{⊥}^{2} ρ_{i}^{2}) \exp (- k_{⊥}^{2} ρ_{i}^{2}), Γ_{1} (k_{⊥}^{2} ρ_{i}^{2}) = I_{1} (k_{⊥}^{2} ρ_{i}^{2})

\exp (- k_{⊥}^{2} ρ_{i}^{2})

and

Z (ξ_{i}) = π^{- 1 / 2} \int_{- \infty}^{+ \infty} \frac{dtexp (- t^{2})}{t - ξ_{i}}

is the plasma function.
其中

Γ_{0} (k_{⊥}^{2} ρ_{i}^{2}) = I_{0} (k_{⊥}^{2} ρ_{i}^{2}) \exp (- k_{⊥}^{2} ρ_{i}^{2}), Γ_{1} (k_{⊥}^{2} ρ_{i}^{2}) = I_{1} (k_{⊥}^{2} ρ_{i}^{2})

\exp (- k_{⊥}^{2} ρ_{i}^{2})

和

Z (ξ_{i}) = π^{- 1 / 2} \int_{- \infty}^{+ \infty} \frac{dtexp (- t^{2})}{t - ξ_{i}}

是等离子体函数。

Combing Eqs. (34), (25), and (26), we can derive the linear dispersion relation as

^{15}

梳理方程。（34）、（25）和（26），我们可以推导线性色散关系为：

^{15}

FIG. 5. Frequency and growth rate of the ITG mode dependence on

k_{⊥} ρ_{i}

in panels (a) and (b) and dependence on

τ = T_{e} / T_{i}

in panels (c) and (d).
图 5.ITG模式

k_{⊥} ρ_{i}

的频率和增长率取决于图（a）和（b）以及图（c）和（d）中的依赖性

τ = T_{e} / T_{i}

。

FIG. 6. (a) Magnetic poloidal flux normalized by a separatrix value:

ψ / | ψ_{0} |

, (b)-(d) are total magnetic field

B

, radial magnetic field

B_{R}

, and axial magnetic field

B_{Z}

normalized by

B_{0} = 531 G

. The circles represent the separatrix.
图 6.（a）由分离值归一化的磁极化通量：

ψ / | ψ_{0} |

，（b）-（d）是总磁场

B

，径向磁场

B_{R}

，轴向磁场

B_{Z}

归一化为

B_{0} = 531 G

。圆圈代表分离。

\begin{aligned} 0 = & [1 - Γ_{0} (k_{⊥}^{2} ρ_{i}^{2})] + \frac{T_{i}}{T_{e}} \frac{e^{2} n_{e 0}}{Z_{i}^{2} n_{i 0}} - η_{i} \frac{ω_{i}^{*}}{ω} ξ_{i} Z (ξ_{i}) \\ \times [(1 - k_{⊥}^{2} ρ_{i}^{2}) Γ_{0} (k_{⊥}^{2} ρ_{i}^{2}) + k_{⊥}^{2} ρ_{i}^{2} Γ_{1} (k_{⊥}^{2} ρ_{i}^{2})] \\ - {- 1 - (1 - \frac{ω_{i}^{*}}{ω}) ξ_{i} Z (ξ_{i}) + η_{i} \frac{ω_{i}^{*}}{ω} ξ_{i}^{2} [1 + ξ_{i} Z (ξ_{i})] \\ - \frac{3}{2} η_{i} \frac{ω_{i}^{*}}{ω} ξ_{i} Z (ξ_{i})} Γ_{0} (k_{⊥}^{2} ρ_{i}^{2}) . \end{aligned}

Next, in order to compare with theory Eq. (35) in the slab limit, we carry out simulations in the domain of Fig.

1 : Z / R_{0} \in [- 0.2, 0.2]

and

ψ

ranges from

ψ (R = 1.702 R_{0}, Z = 0)

ψ (R = 1.705 R_{0}, Z = 0)

, where the magnetic field variation is small and magnetic drift effect can be ignored. In the simulation,

κ_{T i} = {\frac{1}{T_{i}} \frac{\partial T_{i}}{\partial R} |}_{Z = 0} = 5.0, κ_{T e} = \frac{1}{T_{e}}

{\frac{\partial T_{e}}{\partial R} |}_{Z = 0} = 0.0, κ_{n i} = κ_{n e} = {\frac{1}{T_{i}} \frac{\partial n_{i}}{\partial R} |}_{Z = 0} = 2.5, η_{i} = κ_{T i} / κ_{n i} = 2.0

, and the parallel vector

k_{∥}

is fixed:

k_{∥} ρ_{i} = 1.08 \times 10^{- 2}

, and the perpendicular vector is determined by the toroidal component:

k_{⊥} = k_{ζ}

. We first fixed

τ = T_{e} / T_{i} = 0.35

and scanned

k_{⊥} ρ_{i}

(we use deuterium as ion species in the simulation) by increasing

k_{ζ}

; and the simulation results agree well with theory as shown in Figs. 5(a) and 5(b). Then,
接下来，为了与理论方程（35）在板极限中进行比较，我们在图域中进行了模拟。

1 : Z / R_{0} \in [- 0.2, 0.2]

范围

ψ

从

ψ (R = 1.702 R_{0}, Z = 0)

到

ψ (R = 1.705 R_{0}, Z = 0)

，其中磁场变化很小，磁漂移效应可以忽略不计。在仿真中，

κ_{T i} = {\frac{1}{T_{i}} \frac{\partial T_{i}}{\partial R} |}_{Z = 0} = 5.0, κ_{T e} = \frac{1}{T_{e}}

{\frac{\partial T_{e}}{\partial R} |}_{Z = 0} = 0.0, κ_{n i} = κ_{n e} = {\frac{1}{T_{i}} \frac{\partial n_{i}}{\partial R} |}_{Z = 0} = 2.5, η_{i} = κ_{T i} / κ_{n i} = 2.0

，，平行向量

k_{∥}

是固定的：

k_{∥} ρ_{i} = 1.08 \times 10^{- 2}

，垂直向量由环形分量确定：

k_{⊥} = k_{ζ}

。我们首先通过

τ = T_{e} / T_{i} = 0.35

k_{⊥} ρ_{i}

增加

k_{ζ}

;仿真结果与理论吻合较好，如图5（a）和图5（b）所示。然后

FIG. 7. (a) Ion temperature and density gradients at outer midplane

Z = 0

, where

T_{i}

is normalized by electron temperature at the axis

T_{e a} = 80.0 e V, n_{i}

is normalized by electron density at the axis

n_{e a} = 2.44 \times 10^{13} {c m}^{- 3}

, and

r

is normalized by

R_{0}

= 26.8 c m

. (b) Magnetic field strength

B, \nabla B

, and

\nabla \times B

scale length along the parallel direction at field line

ψ = ψ (R = 2.02 R_{0}, Z = 0)

.
图 7.（a）外中平面

Z = 0

处的离子温度和密度梯度，其中

T_{i}

由轴

T_{e a} = 80.0 e V, n_{i}

上的电子温度归一化，由轴上的电子密度归一化，并由

r

归

n_{e a} = 2.44 \times 10^{13} {c m}^{- 3}

一化。

R_{0}

= 26.8 c m

（b）磁场强度

B, \nabla B

和

\nabla \times B

沿磁力线

ψ = ψ (R = 2.02 R_{0}, Z = 0)

平行方向的刻度长度。

FIG. 8. (a) Electrostatic potential of the ITG mode in FRC. The circles represent the separatrix. The mode structures of real and imaginary parts (b) along the radial direction at the outer midplane

(Z = 0)

as shown by the black dotted line in

(a)

, and the separatrix radial location is

R / R_{0} = 1.42

, and (c) along the parallel direction [field line with

ψ (R = 2.02 R_{0}, Z = 0)

as shown by the black dashed line in

(a)]

. TABLE I. Parallel domain size effects on the ITG mode.
图 8.（a） FRC中ITG模式的静电势。圆圈代表分离。实部和虚部的模态结构（b）沿外中平面

(Z = 0)

的径向方向，如中的

(a)

黑色虚线所示，分离的径向位置为

R / R_{0} = 1.42

，以及（c）沿平行方向 [场线

ψ (R = 2.02 R_{0}, Z = 0)

如中的

(a)]

黑色虚线所示。表 I. 并行域大小对 ITG 模式的影响。

Domain size 域名大小	Frequency 频率	Growth rate 增长率	Parity of real part 实部奇偶校验
$Z / R_{0} \in [- 9.37, 9.37]$	$- 0.0108 Ω_{c p}$	$0.0022 Ω_{c p}$	Even 甚至
$Z / R_{0} \in [- 13.6, 13.6]$	$- 0.0046 Ω_{c p}$	$0.0046 Ω_{c p}$	Even 甚至
$Z / R_{0} \in [- 16.2, 16.2]$	$- 0.0022 Ω_{c p}$	$0.0030 Ω_{c p}$	Even 甚至
$Z / R_{0} \in [- 19.2, 19.2]$	$- 0.0052 Ω_{c p}$	$0.0035 Ω_{c p}$	Even 甚至
$Z / R_{0} \in [- 21.0, 21.0]$	$- 0.0043 Ω_{c p}$	$0.0037 Ω_{c p}$	Odd

we fixed

k_{⊥} ρ_{i} = 2.0

and scanned

τ

; and GTC-X simulation results are in good agreement with theory as shown in Figs. 5(c) and 5(d). In Fig. 5 , the frequency and growth rate are normalized by

v_{t h, i} / L

, where

v_{t h, i} = \sqrt{T_{i} / m_{i}}

is the ion thermal speed, and

L = 0.4 R_{0}

is the parallel length of simulation domain.
我们修复

k_{⊥} ρ_{i} = 2.0

并扫描

τ

;GTC-X仿真结果与理论吻合较好，如图5（c）和图5（d）所示。在图5中，频率和增长速率归一化为，

v_{t h, i} / L

其中

v_{t h, i} = \sqrt{T_{i} / m_{i}}

是离子热速度，

L = 0.4 R_{0}

是模拟域的平行长度。

15. B. FRC geometry 15. B. FRC几何形状

Next, we carry out simulation in the middle part of a realistic FRC geometry. The deuterium plasma is simulated in this case which is consistent with experiments.

^{10}

The inner boundary and outer boundary of the simulation domain are

ψ_{0} = 0.96 ψ_{O}

and

ψ_{1} = 13.8 | ψ_{O} |

, where

ψ_{O}

is the poloidal magnetic flux value at the magnetic axis,

R / R_{0}

\in [0.0, 2.45]

in the radial direction at the outer midplane,

Z / R_{0}

\in [- 9.37, 9.37]

in the axial direction,

R_{0} = 26.8 c m

, and the region close to the magnetic axis

ψ \in [0, ψ_{0}]

is excluded for the gyrokinetic model validation. Both particle and field boundary conditions at the field line ends

(Z / R_{0} = \pm 9.37)

are periodic since the FRC center is far away from the divertor region, where the sheath effect becomes important. In this paper, we do not take into account the presheath and sheath effects on the ITG in the FRC center and formation region. In the radial boundaries, i.e.,

ψ = ψ_{0}

and

ψ = ψ_{1}

, we apply the reflected boundary condition for particle and zero boundary condition for fields since the fluctuations of ITG are considered to be zero there. The

ψ, B

B_{R}

, and

B_{Z}

value are shown in Fig. 6, which are normalized by the magnetic field strength

B_{0} = 531 G

at the cross point between the separatrix and the outer midplane:

(R = R_{X} = 1.42 R_{0}, Z = 0)

. The ion temperature and density gradients are shown in Fig. 7 . We choose to simulate the ITG mode with the toroidal mode number

n = 20

, and the corresponding

k_{ζ} ρ_{i} = 0.36

. The elongated parallel mode structure and finite radial structure of the ITG mode are found in our simulation as shown in Fig. 8, i.e.,

k_{∥} ≪ k_{r}

and

k_{∥} ≪ k_{ζ}

, and we can estimate

k_{r} ρ_{i} \approx 0.57

. According to the parallel mode structure along the field line

ψ

= ψ (R = 2.02 R_{0}, Z = 0)

as shown in Fig.

8 (c)

, both the real and imaginary parts of electrostatic potential perturbation of

n = 20

mode are even parity. The frequency and growth rates are

ω_{r} = - 5.98 v_{t h, i} /

L_{S} = - 0.0108 Ω_{c p}

and

γ = 1.22 v_{t h, i} / L_{S} = 0.0022 Ω_{c p}

, where

L_{S}

= 18.8 R_{0}

is the field line distance between left and right boundaries and

Ω_{c p} = e B_{0} / ({c m}_{p})

is the proton cyclotron frequency. GTC-X simulation of the ITG mode in the local FRC geometry agrees with another gyrokinetic particle code ANC, which gives the result

ω_{r, A N C} = - 5.30 v_{t h, i} / L_{S}

and

γ_{A N C} = 1.17 v_{t h, i} / L_{S .}^{33, 34}

The location of the unstable ITG mode is determined by a balance of temperature drive,

η_{i}

value, and magnetic drift strength. In FRC, the magnetic gradient drift is opposite to diamagnetic drift, and magnetic curvature drift is the same as the diamagnetic drift in the outer midplane and opposite to diamagnetic drift in mirror throats.
接下来，我们在逼真的FRC几何图形的中间部分进行仿真。在这种情况下模拟了氘等离子体，这与实验一致。

^{10}

仿真域的内边界和外边界为

ψ_{0} = 0.96 ψ_{O}

和

ψ_{1} = 13.8 | ψ_{O} |

，其中

ψ_{O}

是磁轴处的极化磁通量值，

R / R_{0}

\in [0.0, 2.45]

在外中平面处的径向方向上，

Z / R_{0}

\in [- 9.37, 9.37]

轴向上的极化磁通量值，

R_{0} = 26.8 c m

并且在陀螺动力学模型验证中排除了靠近磁轴

ψ \in [0, ψ_{0}]

的区域。场线末端

(Z / R_{0} = \pm 9.37)

的粒子和场边界条件都是周期性的，因为FRC中心远离分流器区域，在分流器区域中，鞘效应变得很重要。本文不考虑FRC中心和地层区域ITG的预鞘和鞘效应。在径向边界中，

ψ = ψ_{0}

即和

ψ = ψ_{1}

，我们对粒子应用反射边界条件，对场应用零边界条件，因为 ITG 的涨落在那里被认为是零。

ψ, B

、

B_{R}

和

B_{Z}

值如图 6 所示，它们由分离线和外中平面之间交叉点处的磁场强度

B_{0} = 531 G

归一化：

(R = R_{X} = 1.42 R_{0}, Z = 0)

。离子温度和密度梯度如图7所示。我们选择使用环形模态数

n = 20

和相应的

k_{ζ} ρ_{i} = 0.36

.ITG模态的细长并行模态结构和有限径向模态结构如图8所示，即

k_{∥} ≪ k_{r}

和

k_{∥} ≪ k_{ζ}

，我们可以估计

k_{r} ρ_{i} \approx 0.57

。根据并联模式结构沿场线

ψ

= ψ (R = 2.02 R_{0}, Z = 0)

如图1所示。

8 (c)

，静电势扰动模式

n = 20

的实部和虚部都是偶偶校验的。频率和增长速率是

ω_{r} = - 5.98 v_{t h, i} /

L_{S} = - 0.0108 Ω_{c p}

和

γ = 1.22 v_{t h, i} / L_{S} = 0.0022 Ω_{c p}

，其中

L_{S}

= 18.8 R_{0}

是左右边界之间的场线距离，

Ω_{c p} = e B_{0} / ({c m}_{p})

是质子回旋加速器频率。GTC-X在局部FRC几何结构中对ITG模式的仿真与另一个陀螺动力学粒子代码ANC一致，该代码给出了

ω_{r, A N C} = - 5.30 v_{t h, i} / L_{S}

结果和

γ_{A N C} = 1.17 v_{t h, i} / L_{S .}^{33, 34}

不稳定ITG模式的位置由温度驱动、

η_{i}

值和磁漂移强度的平衡决定。在FRC中，磁梯度漂移与抗磁漂移相反，磁曲率漂移与外中平面的抗磁漂移相同，与镜喉中的抗磁漂移相反。

16. GLOBAL GYROKINETIC SIMULATION OF ITG INSTABILITY
16. ITG不稳定性的全局陀螺动力学模拟

In this section, we study the global effects of FRC on the ITG. Using FRC equilibrium as shown in Fig. 1, we choose the parallel domain

Z / R_{0}

[- 13.6, 13.6]

[- 16.2, 16.2], [- 19.2, 19.2]

, and

[- 21.0, 21.0]

, respectively, for simulations of deuterium plasmas. The plasma profile and radial domain are the same with Sec. IV B as shown in Fig. 7(a). We still focus on the

n = 20

toroidal mode

(k_{ζ} ρ_{i} = 0.36)

and compare the mode structure, frequency, and growth rate among simulations with different parallel domain lengths, which are shown in Table I. The 2D poloidal mode structure and

1 D

parallel mode
在本节中，我们研究了FRC对ITG的全球影响。使用如图1所示的FRC平衡，我们分别选择并行域

Z / R_{0}

为

[- 13.6, 13.6]

、

[- 16.2, 16.2], [- 19.2, 19.2]

和

[- 21.0, 21.0]

，用于模拟氘等离子体。等离子体轮廓和径向结构域与Sec.IV B相同，如图7（a）所示。我们仍然关注

n = 20

环形模

(k_{ζ} ρ_{i} = 0.36)

态，并比较了不同平行域长度仿真的模态结构、频率和增长速率，如表I所示。二维极化模态结构和

1 D

并行模态

FIG. 9. Comparison of 2D poloidal mode structures of ITG instability between different parallel domain simulations: (a)

Z / R_{0} \in [- 13.6

, 13.6], (b)

Z / R_{0} \in [- 16.2

, 16.2], (c)

Z / R_{0} \in [- 19.2, 19.2]

, and (d)

Z / R_{0} \in [- 21.0, 21.0]

. The dashed lines show the flux surfaces with the maximum mode amplitude. The blue solid line is the separatrix.
图 9.不同并行域模拟中ITG不稳定性的二维极化模结构比较：（a）、

Z / R_{0} \in [- 13.6

13.6]、（b）、

Z / R_{0} \in [- 16.2

16.2]、（c）

Z / R_{0} \in [- 19.2, 19.2]

和（d）。

Z / R_{0} \in [- 21.0, 21.0]

虚线表示具有最大振幅的磁通量表面。蓝色实线是分离线。

FIG. 10. Comparison of

1 D

parallel mode structures on the diagnostic flux surfaces between different parallel domain simulations.

L / ρ_{i}

is the field line distance normalized by

ρ_{i}

at local flux surface. The diagnostic flux surfaces are shown by the dashed lines in Fig. 9.
图 10.不同并行域仿真之间诊断通量面上的

1 D

并行模式结构比较。

L / ρ_{i}

是局部

ρ_{i}

通量表面归一化的场线距离。诊断磁通面如图 9 中的虚线所示。

structure are shown in Figs. 9 and 10, respectively. The global FRC simulation shows that the ITG mode grows along the field line direction with long parallel wavelength mode structure in SOL crossing the regions of the central chamber and formation sections with the transition mirror regions included. The maximum amplitude of this mode is in the formation region with bad curvature. Meanwhile, in central FRC, the amplitude is much lower as shown in Fig. 9. In the core, ITG is stable. This is due to the large electron parallel dynamics, large Larmor radius effect, and large grad-B drift the direction of which is opposite to ion diamagnetic drift.

^{27}

We also find that the ITG parity changes from even to odd by increasing the parallel domain length. For the even parity modes, when increasing the domain size, the frequency shows a down shift, while the growth rate first increases and then decreases, which is caused by the balance between good curvature and bad curvature. The odd parity mode becomes more unstable than the even parity mode by including more formation region

(Z / R_{0} \leq - 10

and

Z / R_{0} \geq 10)

into the parallel domain, and the frequency and growth rate vary with the domain size. Furthermore, by comparing the mode structures from simulations with different domain sizes, the maximum amplitude location of the modes is in the formation region. From our global gyrokinetic particle simulation, it shows that ITG type drift waves are stable in the core region and unstable in the SOL region in FRC.
结构分别如图9和图10所示。全局FRC仿真表明，在SOL中，ITG模式沿场线方向增长，具有长平行波长模式结构，穿过中心腔室和编队部分区域，包括过渡镜区域。该模式的最大振幅在曲率不好的地层区域。同时，在FRC中心，振幅要低得多，如图9所示。在核心方面，ITG是稳定的。这是由于大的电子平行动力学、大的拉莫尔半径效应和大的 B 级漂移，其方向与离子抗磁漂移相反。

^{27}

我们还发现，通过增加并行域长度，ITG奇偶校验从偶数变为奇数。对于偶奇偶校验模式，当增加域大小时，频率呈现下移，而增长速率先增后减，这是由好曲率和坏曲率之间的平衡引起的。奇偶校验模式比偶偶校验模式更不稳定，包括更多的形成区域

(Z / R_{0} \leq - 10

并

Z / R_{0} \geq 10)

进入并行域，并且频率和增长速率随域大小而变化。此外，通过比较不同域大小的模拟模式结构，模式的最大振幅位置位于形成区域。从我们的全局陀螺动力学粒子模拟中可以看出，在FRC中，ITG型漂移波在核心区域是稳定的，在SOL区域是不稳定的。

17. CONCLUSIONS 17. 结论

In this paper, we have studied the global dynamics of drift waves in the FRC plasma that consists of various elements such as the closed field FRC central region with its squeezed SOL, strong mirrors, the formation sections, another set of strong mirrors, and divertors. It is important that we can now look at the globally developed drift waves that are extended in the axial direction, while localized only in the SOL region. This study was enabled by developing the global GTC particle model: a new global particle-in-cell code GTC-X. Two sets of coordinates, field line coordinates and cylindrical coordinates, are used in the simulation, which enable the maximum numerical efficiency for crossseparatrix simulation. A field aligned mesh is applied to suppress the unphysical short wave length noise and dramatically decrease the computational cost. GTC-X is well benchmarked with theory and ANC for ITG instability with adiabatic electrons. Global FRC simulations show that the ITG mode is unstable in the SOL distributed along the field line and stable in the core, which is consistent with local simulations and experiments. By extending the simulation domain along the SOL field lines, the parity of the ITG mode changes from even to odd, which shows that the parallel domain size is important for determining the most unstable eigenmode of the drift wave. The nonlinear simulation of ITG transport with kinetic electrons in global FRC geometry will be reported in a future paper.
在本文中，我们研究了FRC等离子体中漂移波的全局动力学，该等离子体由各种元素组成，例如具有压缩SOL的封闭场FRC中心区域，强镜，形成部分，另一组强镜和分流器。重要的是，我们现在可以查看全球发展的漂移波，这些漂移波沿轴向延伸，而仅局限于 SOL 区域。这项研究是通过开发全局GTC粒子模型来实现的：一种新的全局粒子在细胞中代码GTC-X。仿真中使用了场线坐标和圆柱坐标两组坐标，从而实现了交叉分离模拟的最大数值效率。应用场对齐网格来抑制非物理的短波长噪声，并显著降低计算成本。GTC-X 在绝热电子的 ITG 不稳定性方面进行了理论和 ANC 的良好基准测试。全局FRC模拟表明，ITG模式在沿场线分布的SOL中不稳定，在磁芯中稳定，这与局部模拟和实验一致。通过沿SOL场线扩展仿真域，ITG模式的奇偶性由偶数变为奇数，这表明并行域大小对于确定漂移波最不稳定的特征模态非常重要。ITG在全局FRC几何结构中具有动力学电子的传输的非线性模拟将在未来的论文中报告。

18. ACKNOWLEDGMENTS 18. 致谢

J.B. and Z.L. would like to acknowledge useful discussions with Dr. L. Shi, Dr. W. L. Zhang, Dr. H. S. Xie, and Dr. Y. Xiao. This work was supported by TAE Grant No. TAE-200441, DOE SciDAC ISEP center, and the Strategic Priority Research Program of Chinese Academy of Sciences under Grant No. XDB16010300. GTC-X simulations used resources on the Oak Ridge Leadership Computing Facility at Oak Ridge National Laboratory (DOE Contract No. DEAC05-00OR22725) and the National Energy Research Scientific Computing Center (DOE Contract No. DE-AC02-05CH11231).
J.B. 和 Z.L. 感谢与 L. Shi 博士、W. L. Zhang 博士、H. S. Xie 博士和 Y. Xiao 博士的有益讨论。这项工作得到了TAE资助号的支持。TAE-200441，美国能源部SciDAC ISEP中心，中国科学院战略优先研究计划，资助号。XDB16010300。GTC-X模拟使用了橡树岭国家实验室（DOE合同编号）的橡树岭领导计算设施（DOE合同编号）上的资源。DEAC05-00OR22725）和国家能源研究科学计算中心（DOE合同编号。DE-AC02-05CH11231）。

19. REFERENCES 19. 参考资料

^{1}

M. Tuszewski, "Field reversed configurations," Nucl. Fusion 28(11), 2033 (1988).

^{1}

M. Tuszewski，“场反向配置”，Nucl。融合 28（11），2033 （1988）。

^{2}

M. W. Binderbauer, H. Y. Guo, M. Tuszewski, S. Putvinski, L. Sevier, D. Barnes, N. Rostoker, M. G. Anderson, R. Andow, L. Bonelli et al., "Dynamic formation of a hot field reversed configuration with improved confinement by supersonic merging of two colliding high-compact toroids," Phys. Rev. Lett. 105, 045003 (2010)

^{2}

M. W. Binderbauer， H. Y. Guo， M. Tuszewski， S. Putvinski， L. Sevier， D. Barnes， N. Rostoker， M. G. Anderson， R. Andow， L. Bonelli et al.， “Dynamic formation of a hot field reversed configuration with improved confinement by supersonic mergeing of two colliding high-compact toroids，” Phys. Rev. Lett. 105， 045003 （2010）

^{3}

M. Tuszewski, A. Smirnov, M. C. Thompson, S. Korepanov, and T. Akhmetov, "Field reversed configuration confinement enhancement through edge biasing and neutral beam injection," Phys. Rev. Lett. 108, 255008 (2012).

^{3}

M. Tuszewski、A. Smirnov、M. C. Thompson、S. Korepanov 和 T. Akhmetov，“通过边缘偏置和中性光束注入增强场反向配置约束”，Phys. Rev. Lett. 108,255008 （2012）。

^{4}

M. W. Binderbauer, T. Tajima, L. C. Steinhauer, E. Garate, M. Tuszewski, L. Schmitz, H. Y. Guo, A. Smirnov, H. Gota, D. Barnes et al., "A high performance field-reversed configuration," Phys. Plasmas 22, 056110 (2015).

^{4}

M. W. Binderbauer、T. Tajima、L. C. Steinhauer、E. Garate、M. Tuszewski、L. Schmitz、H. Y. Guo、A. Smirnov、H. Gota、D. Barnes 等人，“高性能场反转配置”，Phys. Plasmas 22， 056110 （2015）。

^{5}

H. Y. Guo, M. W. Binderbauer, T. Tajima, R. D. Milroy, L. C. Steinhauer, X Yang, E. G. Garate, H. Gota, S. Korepanov, A. Necas et al., "Achieving a longlived high-beta plasma state by energetic beam injection," Nat. Commun. 6, 6897 (2015).

^{5}

H. Y. Guo， M. W. Binderbauer， T. Tajima， R. D. Milroy， L. C. Steinhauer， X Yang， E. G. Garate， H. Gota， S. Korepanov， A. Necas et al.， “Achieving a long-life high-beta plasma state by energetic beam injection，” Nat. Commun.6, 6897 (2015).

^{6}

H. Gota, M. W. Binderbauer, T. Tajima, S. Putvinski, M. Tuszewski, S. Dettrick, E. Garate, S. Korepanov, A. Smirnov, M. C. Thompson et al., "Achievement of field-reversed configuration plasma sustainment via 10 MW neutral-beam injection on the C-2U device," Nucl. Fusion 57, 116021 (2017).

^{7}

M. E. Griswold, E. M. Granstedt, M. C. Thompson, K. Knapp, B. Koop, and TAE Team, "Particle and heat flux diagnostics on the C-2W divertor electrodes," Rev. Sci. Instrum. 89, 10J110 (2018).

^{6}

H. Gota、M. W. Binderbauer、T. Tajima、S. Putvinski、M. Tuszewski、S. Dettrick、E. Garate、S. Korepanov、A. Smirnov、M. C. Thompson 等人，“通过 C-2U 器件上的 10 MW 中性光束注入实现场反向配置等离子体维持”，Nucl。Fusion 57,116021（2017 年）。

^{7}

M. E. Griswold、EM Granstedt、MC Thompson、K. Knapp、B. Koop 和 TAE 团队，“C-2W 分流器电极上的粒子和热通量诊断”，Rev. Sci. Instrum。89， 10J110 （2018）.

^{8}

D. C. Barnes, J. L. Schwarzmeier, H. Ralph Lewis, and C. E. Seyler, "Kinetic tilting stability of field-reversed configurations," Phys. Fluids 29, 2616 (1986).

^{8}

D. C. Barnes、J. L. Schwarzmeier、H. Ralph Lewis 和 C. E. Seyler，“场反转构型的动力学倾斜稳定性”，Phys. Fluids 29， 2616 （1986）。

^{9} R

. Horiuchi and T. Sato, "Full magnetohydrodynamic simulation of the tilting instability in a field-reversed configuration," Phys. Fluids B: Plasma Phys. 1, 581 (1989).

^{9} R

。Horiuchi 和 T. Sato，“场反转配置中倾斜不稳定性的全磁流体动力学模拟”，Phys. Fluids B： Plasma Phys. 1， 581 （1989）。

^{10}

L. Schmitz, D. P. Fulton, E. Ruskov, C. Lau, B. H. Deng, T. Tajima, M. W. Binderbauer, I. Holod, Z. Lin, H. Gota, M. Tuszewski, S. A. Dettrick, and L. C. Steinhauer, "Suppressed ion-scale turbulence in a hot high-

β

plasma," Nat. Commun. 7, 13860 (2016).

^{10}

L. Schmitz， D. P. Fulton， E. Ruskov， C. Lau， B. H. 邓， T. Tajima， M. W. Binderbauer， I. Holod， Z. Lin， H. Gota， M. Tuszewski， S. A. Dettrick， and L. C. Steinhauer， “热高

β

等离子体中抑制离子级湍流”， Nat. Commun.7, 13860 (2016).

^{17}

S. Gupta, D. C. Barnes, S. A. Dettrick, E. Trask, M. Tuszewski, B. H. Deng, H. Gota, D. Gupta, K. Hubbard, S. Korepanov, M. C. Thompson, K. Zhai, T. Tajima, and TAE Team, "Transport studies in high-performance field reversed configuration plasmas," Phys. Plasmas 23, 052307 (2016).

^{17}

S. Gupta、DC Barnes、SA Dettrick、E. Trask、M. Tuszewski、BH 邓、H. Gota、D. Gupta、K. Hubbard、S. Korepanov、MC Thompson、K. Zhai、T. Tajima 和 TAE 团队，“高性能场反向配置等离子体的传输研究”，Phys. Plasmas 23,052307 （2016）。

^{12}

M. Onofri, P. Yushmanov, S. Dettrick, D. Barnes, K. Hubbard, and T. Tajima, "Magnetohydrodynamic transport characterization of a field reversed configuration," Phys. Plasmas 24, 092518 (2017).

^{12}

M. Onofri、P. Yushmanov、S. Dettrick、D. Barnes、K. Hubbard 和 T. Tajima，“磁反转构型的磁流体动力学输运表征”，Phys. Plasmas 24,092518 （2017）。

^{13}

L. C. Steinhauer, H. L. Berk, and TAE Team, "Coupled transport in fieldreversed configurations," Phys. Plasmas 25, 022503 (2018).

^{13}

L. C. Steinhauer、H. L. Berk 和 TAE 团队，“场反构型中的耦合传输”，Phys. Plasmas 25,022503 （2018）。

^{14}

W. W. Lee, "Gyrokinetic approach in particle simulation," Phys. Fluids 26, 556 (1983).

^{14}

W. W. Lee，“粒子模拟中的回动力学方法”，Phys. Fluids 26， 556 （1983）。

^{15}

W. W. Lee and W. M. Tang, "Gyrokinetic particle simulation of ion temperature gradient drift instabilities," Phys. Fluids 31(3), 612-624 (1988).

^{15}

W. W. Lee 和 W. M. Tang，“离子温度梯度漂移不稳定性的旋动力学粒子模拟”，Phys. Fluids 31（3），612-624 （1988）。

^{16}

A. M. Dimits and W. W. Lee, "Partially linearized algorithms in gyrokinetic particle simulation," J. Comput. Phys. 107, 309-323 (1993).

^{16}

A. M. Dimits 和 W. W. Lee，“陀螺动力学粒子模拟中的部分线性化算法”，J. Comput。Phys. 107， 309-323 （1993）.

^{17}

S. E. Parker and W. W. Lee, "A fully nonlinear characteristic method for gyrokinetic simulation," Phys. Fluids B: Plasma Phys. 5, 77 (1993).

^{17}

S. E. Parker 和 W. W. Lee，“陀螺动力学模拟的全非线性特性方法”，Phys. Fluids B： Plasma Phys. 5， 77 （1993）。

^{18}

B. I. Cohen, D. C. Barnes, J. M. Dawson, G. W. Hammett, W. W. Lee, G. D. Kerbel, J.-N. Leboeuf, P. C. Liewer, T. Tajima, and R. E. Waltz, "The numerical tokamak project: Simulation of turbulent transport," Comput. Phys. Commun. 87(1-2), 1-15 (1995).

^{18}

B.I.科恩，DC巴恩斯，JM道森，GW哈米特，WW李，GD科贝尔，J.-N.Leboeuf、P. C. Liewer、T. Tajima 和 R. E. Waltz，“数值托卡马克项目：湍流传输的模拟”，Comput。物理。87(1-2), 1-15 (1995).

^{19}

J. Y. Kim, Y. Kishimoto, M. Wakatani, and T. Tajima, "Poloidal shear flow effect on toroidal ion temperature gradient mode: A theory and simulation," Phys. Plasmas 3(10), 3689-3695 (1996).

^{19}

J. Y. Kim、Y. Kishimoto、M. Wakatani 和 T. Tajima，“环形离子温度梯度模式的极化剪切流效应：理论与模拟”，Phys. Plasmas 3（10）， 3689-3695 （1996）。

^{20}

Y. Kishimoto, T. Tajima, W. Horton, M. J. LeBrun, and J. Y. Kim, "Theory of self organized critical transport in tokamak plasmas," Phys. Plasmas 3(4), 1289-1307 (1996)

^{20}

Y. Kishimoto、T. Tajima、W. Horton、M. J. LeBrun 和 J. Y. Kim，“托卡马克等离子体中的自组织临界输运理论”，Phys. Plasmas 3（4），1289-1307 （1996）

^{21}

Z. Lin, T. S. Hahm, W. W. Lee, W. M. Tang, and R. B. White, "Turbulent transport reduction by zonal flows: Massively parallel simulations," Science 281, 1835 (1998).

^{21}

Z. Lin、T. S. Hahm、W. W. Lee、W. M. Tang 和 R. B. White，“纬向流减少湍流传输：大规模并行模拟”，Science 281,1835 （1998）。

^{22}

C. S. Chang, S. Ku, P. H. Diamond, Z. Lin, S. Parker, T. S. Hahm, and N. Samatova, "Compressed ion temperature gradient turbulence in diverted tokamak edge," Phys. Plasmas 16, 056108 (2009).

^{22}

C. S. Chang、S. Ku、P. H. Diamond、Z. Lin、S. Parker、TS Hemm 和 N. Samatova，“转移托卡马克边缘的压缩离子温度梯度湍流”，Phys. Plasmas 16,056108 （2009）。

^{23}

S. Ethier, C.-S. Chang, S.-H. Ku, W. Lee, W. Wang, Z. Lin, and W. Tang, "NERSCs impact on advances of global gyrokinetic PIC codes for fusion energy research," Comput. Sci. Eng. 17, 10-21 (2015).

^{24}

H. S. Xie, Y. Xiao, and Z. Lin, "New paradigm for turbulent transport across a steep gradient in toroidal plasmas," Phys. Rev. Lett. 118, 095001 (2017).

^{23}

S.埃西尔，C.-S.张， S.-H.Ku、W. Lee、W. Wang、Z. Lin 和 W. Tang，“NERSCs 对聚变能研究全球陀螺动力学 PIC 代码进展的影响”，Comput。科学工程 17， 10-21 （2015）.

^{24}

H. S. Xie、Y. Xiao 和 Z. Lin，“环形等离子体陡峭梯度湍流传输的新范式”，Phys. Rev. Lett. 118， 095001 （2017）。

^{25}

D. P. Fulton, C. K. Lau, I. Holod, Z. Lin, and S. Dettrick, "Gyrokinetic particle simulation of a field reversed configuration," Phys. Plasmas 23, 012509 (2016).

^{25}

D. P. Fulton、C. K. Lau、I. Holod、Z. Lin 和 S. Dettrick，“场反构型的回旋动力学粒子模拟”，Phys. Plasmas 23,012509 （2016）。

^{26}

D. P. Fulton, C. K. Lau, L. Schmitz, I. Holod, Z. Lin, T. Tajima, M. W. Binderbauer, and TAE Team, "Gyrokinetic simulation of driftwave instability in field-reversed configuration," Phys. Plasmas 23, 056111 (2016).

^{26}

D. P. Fulton、C. K. Lau、L. Schmitz、I. Holod、Z. Lin、T. Tajima、MW Binderbauer 和 TAE 团队，“场反演配置中漂移波不稳定性的回转动力学模拟”，Phys. Plasmas 23,056111 （2016）。

^{27}

C. K. Lau, D. P. Fulton, I. Holod, Z. Lin, M. Binderbauer, T. Tajima, and L. Schmitz, "Drift-wave stability in the field-reversed configuration," Phys. Plasmas 24, 082512 (2017).

^{27}

C. K. Lau、D. P. Fulton、I. Holod、Z. Lin、M. Binderbauer、T. Tajima 和 L. Schmitz，“场反转配置中的漂移波稳定性”，Phys. Plasmas 24,082512 （2017）。

^{28}

W. Zhang, Z. Lin, and L. Chen, "Transport of energetic particles by microturbulence in magnetized plasmas," Phys. Rev. Lett. 101, 095001 (2008).

^{28}

W. Zhang、Z. Lin 和 L. Chen，“磁化等离子体中微湍流传输高能粒子”，Phys. Rev. Lett. 101， 095001 （2008）。

^{29}

H. S. Zhang, Z. Lin, and I. Holod, "Nonlinear frequency oscillation of Alfven eigenmodes in fusion plasmas," Phys. Rev. Lett. 109, 025001 (2012).

^{29}

H. S. Zhang、Z. Lin 和 I. Holod，“聚变等离子体中 Alfven 特征模态的非线性频率振荡”，Phys. Rev. Lett. 109， 025001 （2012）。

^{30}

Z. Wang, Z. Lin, I. Holod, W. W. Heidbrink, B. Tobias, M. Van Zeeland, and M. E. Austin, "Radial localization of toroidicity-induced Alfven eigenmodes," Phys. Rev. Lett. 111, 145003 (2013).

^{30}

Z. Wang、Z. Lin、I. Holod、W. W. Heidbrink、B. Tobias、M. Van Zeeland 和 M. E. Austin，“toroidicity 诱导的 Alfven 特征模态的径向定位”，Phys. Rev. Lett. 111， 145003 （2013）。

^{31}

J. McClenaghan, Z. Lin, I. Holod, W. Deng, and Z. Wang, "Verification of gyrokinetic particle simulation of current-driven instability in fusion plasmas. I. Internal kink mode," Phys. Plasmas 21, 122519 (2014).

^{31}

J. McClenaghan，Z. Lin，I. Holod，W. 邓和Z. Wang，“聚变等离子体中电流驱动不稳定性的陀螺动力学粒子模拟的验证。I. 内部扭结模式，“Phys. Plasmas 21， 122519 （2014）。

^{32}

D. Liu, W. Zhang, J. McClenaghan, J. Wang, and Z. Lin, "Verification of gyrokinetic particle simulation of current-driven instability in fusion plasmas. II. Resistive tearing mode," Phys. Plasmas 21, 122520 (2014).

^{32}

D. Liu、W. Zhang、J. McClenaghan、J. Wang 和 Z. Lin，“聚变等离子体中电流驱动不稳定性的陀螺动力学粒子模拟的验证。II. 电阻撕裂模式，“Phys. Plasmas 21， 122520 （2014）。

^{33} C

. K. Lau, "Electrostatic turbulence and transport in the field-reversed configuration," Ph.D. dissertation (University of California, Irvine, 2017).

^{33} C

。K. Lau，“场反构中的静电湍流和输运”，博士论文（加州大学欧文分校，2017 年）。

^{34}

C. K. Lau, D. P. Fulton, J. Bao, Z. Lin, T. Tajima, L. Schmitz, and TAE Team, "First simulations of turbulent transport in the field-reversed configuration," in Proceedings of the 26th International Conference on Plasma Physics and Controlled Nuclear Fusion Research, Ahmedabad (International Atomic Energy Agency, Vienna, Austria, 2018), Paper No. IAEA-CN-258/TH/p6-37.

^{34}

C. K. Lau、D. P. Fulton、J. Bao、Z. Lin、T. Tajima、L. Schmitz 和 TAE 团队，“场反转配置中湍流传输的首次模拟”，载于艾哈迈达巴德第 26 届等离子体物理和受控核聚变研究国际会议论文集（国际原子能机构，奥地利维也纳，2018 年），论文编号。原子能机构-CN-258/TH/P6-37。

^{35} S

. Ethier, W. Tang, and Z. Lin, "Gyrokinetic particle-in-cell simulations of plasma microturbulence on advanced computing platforms," J. Phys.: Conf. Ser. 16, 1-15 (2005).

^{35} S

。Ethier、W. Tang 和 Z. Lin，“先进计算平台上等离子体微湍流的回转动力学粒子细胞模拟”，J. Phys.：Conf. Ser. 16， 1-15 （2005）。

^{36}

Y. Nishimura and Z. Lin, "A finite element mesh in a tokamak edge geometry," Plasma Phys. 46(7-9), 551-556 (2006).

^{36}

Y. Nishimura 和 Z. Lin，“托卡马克边缘几何中的有限元网格”，Plasma Phys. 46（7-9）， 551-556 （2006）。

^{37}

L. Galeotti, D. C. Barnes, F. Ceccherini, and F. Pegoraro, "Plasma equilibria with multiple ion species: Equations and algorithm," Phys. Plasmas 18, 082509 (2011).

^{37}

L. Galeotti、DC Barnes、F. Ceccherini 和 F. Pegoraro，“具有多种离子种类的等离子体平衡：方程和算法”，Phys. Plasmas 18,082509 （2011）。

^{38}

Y. Xiao, I. Holod, Z. Wang, Z. Lin, and T. Zhang, "Gyrokinetic particle simulation of microturbulence for general magnetic geometry and experimental profiles," Phys. Plasmas 22, 022516 (2015).

^{38}

Y. Xiao、I. Holod、Z. Wang、Z. Lin 和 T. Zhang，“用于一般磁几何形状和实验剖面的微湍流的回旋动力学粒子模拟”，Phys. Plasmas 22,022516 （2015）。

^{39}

R. B. White and M. S. Chance, "Hamiltonian guiding center drift orbit calculation for plasmas of arbitrary cross section," Phys. Fluids 27, 2455 (1984).

^{39}

R. B. White 和 M. S. Chance，“任意横截面等离子体的哈密顿引导中心漂移轨道计算”，Phys. Fluids 27， 2455 （1984）。

^{40} A

. J. Brizard, "On the validity of the guiding-center approximation in the presence of strong magnetic gradients," Phys. Plasmas 24, 042115 (2017).

^{40} A

。J. Brizard，“在强磁梯度存在的情况下指导中心近似的有效性”，Phys. Plasmas 24,042115（2017）。

^{41}

Y. Hayakawa, T. Takahashi, and Y. Kondoh, "Classification of particle orbits and related stochasticity of plasma ion motion in a field-reversed configuration with D3 He advanced fuel," Nucl. Fusion 42, 1075-1084 (2002).

^{41}

Y. Hayakawa、T. Takahashi 和 Y. Kondoh，“使用 D3 He 先进燃料在场反转配置中粒子轨道的分类和等离子体离子运动的相关随机性”，Nucl。融合 42， 1075-1084 （2002）。

Global simulation of ion temperature gradient instabilities in a field-reversed configuration 场反转配置中离子温度梯度不稳定性的全局模拟1

2. AFFILIATIONS 2. 隶属关系

3. ABSTRACT 3. 摘要

4. INTRODUCTION 4. 引言

5. GLOBAL FRC GEOMETRY IMPLEMENTATION5. 全局FRC几何实现

6. A. Magnetic field representation in cylindrical coordinates6. A. 圆柱坐标中的磁场表示

7. B. Field line coordinates for the perturbed field calculation7. B. 扰动场计算的场线坐标

8. Laplacian operator 8. 拉普拉斯算子

9. PHYSICS MODEL 9. 物理模型

10. A. Formulation 10. A. 配方

11. B. Implementation of dynamic equations in cylindrical coordinates11. B. 圆柱坐标中动力学方程的实现

12. Poisson solver 12. 泊松求解器

13. CODE VERIFICATION 13. 代码验证

14. A. Slab limit 14. A. 板坯限制

15. B. FRC geometry 15. B. FRC几何形状

16. GLOBAL GYROKINETIC SIMULATION OF ITG INSTABILITY16. ITG不稳定性的全局陀螺动力学模拟

17. CONCLUSIONS 17. 结论

18. ACKNOWLEDGMENTS 18. 致谢

19. REFERENCES 19. 参考资料

Global simulation of ion temperature gradient instabilities in a field-reversed configuration
场反转配置中离子温度梯度不稳定性的全局模拟1

5. GLOBAL FRC GEOMETRY IMPLEMENTATION
5. 全局FRC几何实现

6. A. Magnetic field representation in cylindrical coordinates
6. A. 圆柱坐标中的磁场表示

7. B. Field line coordinates for the perturbed field calculation
7. B. 扰动场计算的场线坐标

11. B. Implementation of dynamic equations in cylindrical coordinates
11. B. 圆柱坐标中动力学方程的实现

16. GLOBAL GYROKINETIC SIMULATION OF ITG INSTABILITY
16. ITG不稳定性的全局陀螺动力学模拟