Page 4 - Demo

042504-4
Rath et al.
d U~ w d ~I w d U~ p
Phys. Plasmas 24, 042504 (2017) force exerted by the perturbation on the plasma. The resulting
plasma displacement ~n is given by the constraint that the force
from the perturbed external field is balanced by the force
exerted by the equilibrium field on the displaced plasma cur-
rents, ~n 1⁄4 F 1:F~ . If the plasma is displaced in the ~n direc- vw
tion, this goes along with a change dB~ 1⁄4 ~n r B~ of the pp
plasma-generated magnetic field. The normal component of this field, which is identical to the current potential on the con-
~
The linear model presented so far has a free parameter that may not be immediately obvious. The degree of free- dom arises because there is no clear, unique separation between the plasma and the surrounding vacuum region. Intuitively, we do not expect that rigid motion of “the FRC” will result in, e.g., a displacement of the density pro- file in the divertor region—which is separated from the sep- aratrix by both distance and the magnetic mirror. When calculating the flux that is induced on the wall by the mov- ing plasma currents and the force that is exerted on the plasma currents, we therefore have to define a boundary for the volume integrals. This boundary will be called the “rigid volume” because it defines what part of the plasma we consider as rigidly moving and what part we assume will remain stationary. Theoretically, we expect that the rigid volume encloses at least the separatrix and that it excludes at least part of the SOL (both radially and axially). Numerically, we also require the rigid volume to be con- tained in the control surface.
IV. NON-LINEAR SIMULATIONS A. Outline of the Q2D code
Q2D25 is a hybrid, FRC simulation code that combines a 2-D, extended MHD code (“LamyRidge”) with a 3-D, Monte-Carlo “fast-particle” tracker (“MC”). LamyRidge and MC time steps are executed interleaved. MC follows individ- ual ion and neutral trajectories, while the field components produced by the fluid plasma are held constant. LamyRidge calculates the coupled time-evolution of ion-, electron- and neutral-fluids which exchange energy and momentum with the fast particles. In tokamaks, modeling neutral beam injec- tion using the Monte Carlo method is the standard prac- tice,26,27 and the method has also been previously used to model beam injection into FRCs.28–30 However, while in tokamaks, it is sufficient to use the drift orbit equations (averaged over the Larmor orbit), in FRCs, it is necessary to resolve the full ion orbit using the Lorentz equation.
The total number of particles, energy, and momentum among particles and fluids is conserved. Modeled interac- tions between the fluid and particles include collisions, ioni- zation, and charge exchange. Particles whose energy falls under a threshold become part of the respective fluid. The computational domain is bounded by an arbitrary shaped, axisymmetric resistive wall with Dirichlet boundary condi- tions for the temperature. Impacting particles can be absorbed and recycled.
dt 1⁄4Lw: dt þMwp:R: dt ;
(14)
1⁄4 ðLw þ Mwp:R:MpwÞ
we thus obtain a closed system of linear ODEs describing
ðLw þ Mwp:R:MpwÞ dt 1⁄4 Rw:Iw: (16) C. Numerical implementation
d~I w :
the evolution of the coupled plasma-wall system d ~I w ~
trol surface that would produce the same field, is dBp n^. D. Rigid Volume—A free parameter
dt
(15)
Equations (16) and (11) form the mathematical model of rigid FRC evolution in the presence of an arbitrary shaped, resistive wall. To solve these equations, we have adapted the VALEN code,22 which was developed to simulate the evolu- tion of resistive wall modes in Tokamaks.
VALEN represents the wall as a collection of finite tri- angular and quadrilateral elements with circulating branch currents and varying resistivities. To analyze the positional stability of an FRC, we retain the VALEN functions for cal- culating surface inductances, mutuals, and resistances but replace the (Tokamak specific) use of the DCON code. The calculation of Lw ; Rw ; Mpw , and Mwp is performed by VALEN without special adaptations. The FEniCS software suite23,24 is used to compute the reluctance matrix R. We import the equilibrium plasma current profiles and interpo- late it on a volumetric mesh inside the control surface. The basis functions defined on the control surface are then inter- preted as boundary conditions for the normal magnetic field, and we calculate the resulting ~j B~ force on the plasma for each basis function. Since the force must be a linear operator on the normal field at the boundary, this results in a system of linear equations for a force matrix Fw
ð
FðwÞ1⁄4 ~j B~dV; (17) ij pxi
plasma
where B~ is the (unique) curl-free field that satisfies B~ n^ xx
1⁄4 fjð~xÞ on the boundary. The force matrix can be used to cal- culate the force for any external field B~ (defined through its
x expansion coefficients U~x on the control surface)
ð
~j B~dV1⁄4F:U~: (18) pxwx
plasma
The vacuum force matrix Fv is calculated from Equation (1). Since the basis functions for the control surface correspond to the canonical basis of the finite element representation, the integral in Equation (11) vanishes, and the reluctance matrix R can be identified with the integrand
R1⁄4 1 F 1:F r B~ð~xÞ n^: lvw p
0
(19)
The meaning of this equation becomes clearer if we consider how it acts on a specific external field perturbation U~x. To evaluate R:U~x, we first calculate F~w 1⁄4 Fw:U~x, which is the