Proceedings of ITC/ISHW2007

EMC3-EIRENE implementation on NCSX
and first applications
F. Sardei 1), A. Brooks 2), Y. Feng 1), F. Herrnegger 1),T. Kaiser 3), J. Kisslinger 1),
R. Maingi 4), D. Monticello 2), M.C. Zarnstorff 2)
1)

Max-Planck-Institut für Plasmaphysik, Euratom Association, Greifswald/Garching, Germany
2)
Princeton Plasma Physics Laboratory, Princeton, USA
3)
Lawrence Livermore National Laboratory, Livermore, USA
4)
Oak Ridge National Laboratory, Oak Ridge, USA

The paper presents the status of the implementation of the EMC3-EIRENE code on NCSX.
Particular attention is devoted to the construction of the 3D plasma computational grid in the SOL, which
consists of a 2D radial-poloidal system of finite flux tubes extending half a field period. In a firstapplication, the code simulates the energy transport for electrons and ions in a typical SOL configuration
bounded by model target plates. The distribution of connection lengths defined by the plates and theresulting temperature profiles are shown and discussed.
Keywords: stellarator, NCSX, divertor, transport, 3D modeling, EMC3-EIRENE code

1. Introduction
In helical devices, SOL plasmas and recycling
are intrinsically three-dimensional due to nonaxisymmetric magnetic fields and toroidally
inhomogeneous
divertorplates.
Additional
complexity arises from a large variety of magnetic
topologies encountered in configurations withdifferent magnetic shear and complex magnetic
structures ranging from small or large island chainsto even coexisting slightly or strongly ergodic
regions. For such open topologies standard magnetic
coordinates do not exist. The approach used by theEMC3 code [1,2] is a Monte Carlo (MC) techniqueoptimized for highly anisotropic 3D fluid transport
processes in open magnetic topologies of arbitrarycomplexity. The parallel and cross-field transport
dynamics of mass, momentum and energy are
simulated by advancing MC “particles” along andacross field lines through local orthogonal steps. Upto now, the EMC3 code, coupled self-consistently tothe EIRENE code [3] for neutral transport, has been
applied to W7-AS, W7-X, LHD, TEXTOR-DEDand ITER start-up.
Presently, the NCSX magnetic boundary region[4,5] is being implemented in the code to start theinvestigation of the SOL transport properties for
francesco.sardei@ipp.mpg.de

different configurations on the way to a first divertor
design. A major step of this implementation is theconstruction of a 3D field-aligned grid extending
toroidally half a field period and covering the radial
SOL region of interest. At the boundary of theperiodicity domain, field-line coordinates aretransformed reversibly to the next domain byapplying the RFLM technique [6]. This technique
avoids both a numerical cross-field diffusion of thefast parallel transport and a radial accumulation ofinterpolation errors, which are critical issuesespecially in the high temperature region close to theLCFS.

2. NCSX boundary configuration
NCSX is a quasi-axisymmetric compact
stellarator with three field periods, major radius of
1.4 m and aspect ratio of 4.4 [7]. The plasma
configuration chosen for the present first application

Proceedings of ITC/ISHW2007

Z/cm

Lc /m

Fig. 1 NCSX configuration over half a field period(φ = 0°,30°,60°) with model target plates
extending 1/4 of a field period at stellaratorsymmetric positions.

(Fig. 1) has β = 4.1% and ι = 0.65 at the plasmaedge (Fig. 2). The last closed flux surface (LCFS) issurrounded by a narrow weakly ergodic zone
followed by an island chain at ι = 3/5 (Fig. 1).
Outside this resonance strong radial field-line
excursions to the wall indicate the end of the plasma
region.
In a first approximation, four model target
plates per field period are placed close to the LCFS
at stellarator-symmetric positions. They are locatedat the tips of the bean-shaped cross section, whichhave the largest flux expansion (Fig. 1), and extend
0.8
0.7

6/93/50.6iota

0.5
0.4
0.3
0.2
0.1
00510152520minor radius [ cm]3035Fig. 2 ι profile for the configuration shown in Fig. 1.

1/4 of a field period starting at the bean-shaped
symmetry plane (Fig. 1). The plates cut the LCFS atseveral toroidal positions, defining a limiter SOL
with a highly inhomogeneous distribution ofconnection lengths (Fig. 3). Within a radial distanceof about 3 cm from the LCFS, connection lengths Lc= 15-20 m, i.e. 5-7 field periods, dominate the SOL.
The two relatively large up/down symmetric zones
with Lc≈ 3 m ≈ one field period in the upper andR/cm

Fig. 3 Distribution of connection lengths as definedby the target plates shown in Fig. 1.

lower inboard corners of the φ = 60° symmetry plane(Fig. 3) are due to a local inward protrusion of theinboard plates into the core at toroidal angles of φ =20°-25°. Deviations from up/down symmetry in theoutermost short-Lc regions of the distribution are due
to technical reasons and can be ignored.
3. Grid construction
In this section we shortly describe the structureof a typical 3D computational grid used by theEMC3-EIRENE code for stellarator-symmetricconfigurations and the basic steps of its construction.
An extension to non-symmetric configurations isstraightforward. Precomputed field lines are started
from a fine 2D up/down symmetric grid of nested
contours (Fig. 4) in the symmetry plane where theplates are located (bean-shaped plane of Fig. 1). The
up/down symmetry of the grid is advisable because
it allows a reduction of the computational domain tohalf a field period by a flux-coordinaterepresentation of the boundary conditions andprevents an accumulation of radial interpolation
errors [6]. The field lines are integrated numericallywith high accuracy over half a field period andstored at fine-spaced toroidal intervals. The resulting
two-parametric set of finite flux tubes covering theradial-poloidal region of interest contains all
relevant information about the field structure
including islands and ergodic regions. The radial andpoloidal grid lines defining the cross-sections of theflux tubes become increasingly distorted in toroidal
direction due to field shear. For high-shear devices,a strong distorsion may require a splitting of theflux-tube domain in two or more toroidal
subdomains [6]. However, the distorsion does notaffect the accuracy of the transport, as the flux-tube
geometry does not enter the MC procedure, which isProceedings of ITC/ISHW2007
φ = 0°

Wall (enlarged)
Z/cm

Plasma grid

R/cm

Fig. 4 Up/down symmetric cross-section of thecomputational grid at the starting plane. The
wall (enlarged) is the radial boundary of thegrid.

φ = 60°Z/cm

Plasma gridWall(enlarged)
R/cm

Fig. 5 Non-up/down symmetric cross section of thecomputational grid resulting from mapping of
the starting grid (Fig. 4) over half a fieldperiod. The wall (enlarged) is the radialboundary of the grid.

based on local orthogonal displacements. In case ofNCSX, the distortion is moderate and a single
toroidal flux-tube domain is sufficient. However, thegrid must be sufficiently fine to resolve small lengthscales of the plasma parameters especially near thetarget plates and to provide sufficient accuracy in thefield-line interpolation.
The finite flux-tube SOL domain is boundedradially inwards by two flux surfaces close to theLCFS and outwards by a non-flux surface ("lastsurface") enclosing the relevant SOL plasma region(Figs. 4,5).
The two innermost flux surfaces are needed tocorrectly apply poloidal-dependent flux boundary
conditions for the radial fluxes. They must be very
close to each other and well defined, i.e. sufficiently
sharp to form continuous, up/down symmetrictoroidal surfaces in both planes of symmetry.
The "last surface" is the outer boundary of thefinite flux-tube system defining the plasmacomputational domain. Together with the innermost
flux surfaces, it should enclose a sufficiently large

volume to completely cover the expected radialextension of the plasma SOL. If target plates arepresent, like in NCSX, the plasma SOL is boundedby the shadows of target plates catching most of theoutflowing plasma. In the absence of target plates,
like in LHD, the plasma SOL is bounded by flux
tubes "hot-wired" to the wall, i.e. leading to wallcontact within very short connection lengths. To agood starting approximation, the last surface could
be identified with the radial boundary contour of afine Poincaré distribution of particle orbits traced
along field lines with cross-field diffusion. This
could be the first step towards an optimization of thelast surface, which eventually requires an iterative
procedure based on the full plasma transport.
The space between the innermost and outermost
surfaces is filled with nested surfaces obtained byradial interpolation between the bounding surfaces atthe starting plane (Fig. 4).
In the present NCSX grid approximation, thelast surface was obtained by adding a small toroidal
field, thereby smoothing and shifting outwards theLCFS, and then radially extrapolating the basic
island-free geometry of the edge configuration.
Since the 3D grid reflects the basic island-free
topology of the boundary, it is not sensitive to small
variations of rotational transform, ergodicity andedge island structures. Therefore, the same or aslightly modified grid may be used for different
configurations.
Obviously, the 3D table of stored magnetic
field used to integrate the field lines must fully coverthe 3D plasma region. This requirement should bestrictly controlled by the field-line tracing code to
avoid errors which may be hardly detectable after
the grid construction. Furthermore, both the tableand flux-tube grids should be sufficiently fine toresolve magnetic flux conservation within a few %.
This is an essential requirement not only for the
intrinsic quality of the flux-tube system but also forthe correct description of the parallel transport.
The wall represents the outer radial boundary of
the 3D computational grid. No magnetic field isneeded in the space between the wall and the plasmagrid domain, this region being only populated byneutrals. Both the wall and the plasma facing
components (target plates and baffles) are defined bytheir cross sections at φ = const planes.
Overlapping of grid cells arising from strong
radial/poloidal distorsions, which may occur in theProceedings of ITC/ISHW2007
whole computational grid up to the wall, must beavoided. Within the plasma grid it is sufficient tocheck for flux-tube overlapping at the boundary of
the toroidal domain, whereas in the space between
the plasma grid and the wall surface, overlappingneeds to be checked at any toroidal position, as thepoloidal distributions of the wall grid points are notmagnetically interconnected in toroidal direction.
4. First EMC3 application
In the first test application of the EMC3 code toNCSX, only the 3D energy transport for electrons
and ions were simulated, whereas the plasma density
was kept fixed at ne = 1019 m-3 throughout thecomputational domain. A power of PSOL = 1.2 MW
entering the SOL was equally distributed between
ions and electrons, and the anomalous cross-fieldheat diffusivities were set to χ e = χ i = 3 m2/s. These
numbers have to be considered only as an example,
which does not reflect realistic NCSX parametervalues, as the subject of the present study is the code
implementation, not an edge physics study. As
expected from the Lc contour plot (Fig. 3), therelevant radial energy transport takes place in anarrow radial zone dominated by Lc = 15-20 m. The
corresponding collisionality νe* = Lc/λee = 1-2 < 10indicates that the SOL transport is in the sheathlimited convection regime with no significant

φ = 0°

Z/cm

Te /eVR/cm

Fig. 6 Electron temperature distribution at the φ = 0°symmetry plane.

φ = 60°Z/cm

Ti /eVR/cm

Fig. 7 Electron temperature distribution at the φ =60° symmetry plane.

parallel temperature gradients [8]. The resulting
temperature decay length (Figs. 6,7) is proportional
toL c and is estimated as λT = 2χ ⊥τ e // = 1-1.5

cm. A local energy sink in the two up/downsymmetric SOL regions with Lc≈ 3 m (Fig. 3) isvisible in Fig. 7 by comparing the radial Te profiles
at these two positions with that at the inboard
midplane, which is governed by larger λT. The
spatial distribution of Ti is similar as that of Te, with
Ti being about 10% higher throughout the SOL.
5. Next steps
Major technical steps still to be done towards a
realistic description of the SOL transport in NCSX
with the EMC3-EIRENE code are:• improvement of the target-plate geometry• improvement of the outermost plasma grid surface• grid refinement near the plates• activation of the full transport-equation system(particle, momentum, energy)• activation of the self-consistent neutral gastransport• tests of numeric convergence
Then, dedicated transport studies aiming, forexample, at exploring the divertor potential ofNCSX, could be started. Major tasks towards thisgoal would include• characterizing the basic plasma transport properties
for different 3D SOL configurations• estimating the heat and particle deposition profiles
on the target plates for different configurations andplate geometries• optimizing the target-plate and baffle design
with/without inboard plates with respect to wall
protection, configuration flexibility, homogeneousload, neutral compression• exploring the role of ergodicity and low-order edge
islands for a divertor concept• estimating the impurity retention capability in theSOL• checking the sensitivity of the transport to low-Proceedings of ITC/ISHW2007
order resonances and ergodicity.

References[1] Y. Feng, F. Sardei, J. Kisslinger, J. Nucl. Mater. 266269, 812 (1999)[2] Y. Feng, F. Sardei, J. Kisslinger, P. Grigull,
K. McCormick and D. Reiter, Contrib. Plasma Phys.
44 1-3, 57 (2004)[3] D. Reiter, M. Baelmans, and P. Boerner, Fusion Sci.
Technol. 47, 172 (2005)[4] R. Maingi, et. al., Proc. of 33rd EPS on Plasma Physics
and Controlled Fusion, Roma, Italy, June 19-23 (2006)
paper P5.116[5] T. Kaiser, et. al., Bull. Am. Phys. Soc. 51, 39 (2006)[6] Y. Feng, F. Sardei and J. Kisslinger, Phys. Plasmas 12,052505 (2005)[7] M.C. Zarnstorff, et. al., Plasma Phys. Contr. Fusion 43,
A237 (2001)[8] P.C Stangeby, The Plasma Boundary of Magnetic
Fusion Devices, Plasma Physics Series, IoP 2000