Fusion Res.
Res. SERIES,SERIES,Vol.
Vol.6 (2004)
6 (2004)314–317
000–000
J. Plasma Fusion

A Numerical Approach to the Localization of Passive Line Integrated
Neutral Particle Measurements on LHD
GONCHAROV Pavel R., LYON James F.1, OZAKI Tetsuo and SUDO Shigeru
National Institute for Fusion Science, Toki 509-5292, Japan
1Oak Ridge National Laboratory, Oak Ridge, TN 37831-8072, USA(Received: 9 December 2003 / Accepted: 24 March 2004)

Abstract
The problem of correct interpretation of passive neutral particle measurements from non-axis-symmetrical plasma
has been investigated on LHD. While in active measurements the charge exchange neutral flux originates from a localized artificially created target, the naturally occurring neutral flux observed by a passive diagnostic is line-integrated. A
general formulation of the measured quantity is given for an arbitrary shape of plasma parameter isolines, including theheliotron magnetic configuration. A mathematical approach is presented based on a combination of the experimental
and computational diagnostic techniques to obtain radial profile data with a passive non-perturbing method.
Keywords:
neutral particle analysis, charge-exchange, locality, ion distribution, ion temperature profile, plasma heating,stellarator/heliotron, Large Helical Device (LHD)of a more ‘native’ coordinate, such as minor radius for acylindrically symmetric plasma column or a generalized radial variable ρ in the non-axis-symmetrical case of a helicalplasma. The integration variable then should be changed andthe formula should contain the module of the corresponding
Jacobian. The detailed formulation of the problem follows in
the next section. Then, possible numerical solution methods
are discussed and finally, the application of these methods tothe silicon detector-based neutral particle analyzer (SDNPA)[1,2] data processing is described.

1. Introduction
Neutral particle analysis currently represents one of thebasic methods of ion parameter measurements in fusion plasma experiments. Due to the variety of particle orbit classes
and complex full 3D field geometry of LHD as well as manyof the modern magnetic plasma confinement devices, spatial
and angular resolution is required. The physical approach to
the localization of such measurements is to introduce an artificial localized target for the charge exchange process in theplasma, such as a diagnostic neutral beam or a solid pellet.
This implies active measurements discussed in [1] for LHD.
However, in toroidal devices with magnetic confinement one
can also retrieve the spatial profiles of ion parameters frompassive neutral particle diagnostic data. This paper describes
a numerical approach to this problem applicable to an arbitrary shape of isolines in the diagnostic cross-section with noassumption of toroidal axial symmetry.
The measured atomic energy distribution dN/dE is aresult of superposition of local ion distributions along thediagnostic sight line, taking into account the spatial distributions of the relevant plasma parameters, such as the densityof charge exchange targets, and the attenuation of the atomic
flux on the way out to the plasma edge. The integral relation
between the plasma ion distribution function and theobserved neutral particle spectrum is given in [1]. This integral has a simple form since the distance along the sight line
is chosen as the integration variable. In practice the plasma
parameter profiles in the integrand are expressed as functions

2. Problem formulation
The measured quantity Γ (E) = dN/dE dt is the energy
resolved flux of atoms to the analyzer’s collimating apertureand the sought quantity proportional to the local ion distribution function g(E, ρ) ∝ ni(ρ)fi(E, ρ) is the local birth rate ofatoms at a given energy [erg-1cm-3s-1]. As mentioned above,
the relationship between the measured and the sought quantities has the form of an integral equation. For the following
mathematical treatment let us introduce a new sought funcΩStion ∼g (E, ρ (X')) = 4π a g(E, ρ (X')), where X' is the distancealong the sight line and ρ (X') is the local value of the effective minor radius defined as a square root of the normalized
magnetic flux (Ψ/ΨLCMS)1/2. The geometrical factor includes
the viewing solid angle Ω and the aperture area Sa of the neutral particle analyzer. The parameter ρ characterizes a system©2004 by The Japan Society of Plasma
Science and Nuclear Fusion Research

Corresponding author’s e-mail: pavel@nifs.ac.jp1314Goncharov Pavel R. et al., A Numerical Approach to the Localization of Passive Line Integrated Neutral Particle Measurements on LHD

The a priori knowledge about the isolines of plasma parameters allows to determine the functions Q+(ρ,ζ ) and Q-(ρ,ζ )and to write the relation (2) in the explicit form.

3. Solution Methods
Kinetic Modeling: The problem of retrieving the ion
parameter spatial profiles from Γ (E, ζ ) can be solved by iterative comparing and joining the results of numerical modeling with the experimental data. The most complete and elaborate approach is the kinetic one. It is based on the known
cross-sections of the elementary processes determining theformation of the atoms within the plasma and their passage tothe periphery. The necessary electron temperature and plasma density distributions, as a rule, are known, e.g. from
Thomson scattering and microwave diagnostic measurements. The unknown distributions of required values can betaken as free parameters. The free parameters, such as the ion
temperature profile, should be variated until the numericalresult matches the experimentally observed dN/dE.
The general scheme of this method is as follows. First,
the spatial distribution of neutral atoms n a(r) should beobtained. Consider the quasistationary kinetic equation for
the neutral atom distribution na(r)fa(r, v) in the plasma:Fig. 1 Integration along a viewing chord crossing the system of nested isolines.

of nested contour lines corresponding to the magnetic surfaces with ρ = 0 referring to the magnetic axis. The structureof magnetic surfaces can be obtained from a numerical solution of Grad-Shafranov equation [3].
Suppose that the contour lines are convex closed curves
without self-intersections and the neutral particle analyzer
performs a scan over the angle ζ. This situation is depicted inFig. 1. Assuming that the geometry of measurements is originally described in the laboratory coordinate system XZ,introduce a new coordinate system X’ Z’ so that the axis X’
is parallel to the sight line L at the given angle ζ. The sightline enters the plasma at X1’ crossing the contour line ρ = 1,
then reaches the deepest region ρ = ρmin for the given ζ andfinally exits the plasma at X2’ crossing the contour line ρ = 1
again. Taking into account the atomic flux attenuation, weobtainΓ (ζ ) =∫X2'
X1'X'~(1)where λmfp is the mean free path of atoms with respect to ionizing collisions. The functions ρ (X') and X'(ρ) along the integration path determined by the angle ζ are known from thestructure of the contour lines ρ = const. Let us denote Q+ (ρ,ζ ) = dX'/dρ on the interval between the point ρ = ρmin and thepoint (2) where dX'/dρ > 0 and Q-(ρ, ζ ) = dX'/dρ on the interval between the point (1) and the point ρ = ρmin where dX'/dρ< 0. Changing the integration variable in (1) from X' to ρyields
1Q (ρ~, ζ )λ (E, ρ~)d ρ~Γ (E, ζ ) = e ∫ ρ−−1mfp

1∫ ~g(E, ρ) ×min

ρmin

⎡ +−
∫ Q (ρ~, ζ )λ (E, ρ~)dρ~⎢ Q (ρ, ζ )e ρ⎣ρ−1
mfp+minQ (ρ , ζ )λ (E, ρ )d ρ⎤− Q (ρ, ζ )e ∫ ρ⎥ dρ.
−−ρmin

−~

−1
mfp~~

⎦

∂ na(r)fa(r, v)= [s(r, v) + g(r, v)] − ν na(r)fa(r, v). (3)∂rThe first two terms in the round brackets on the right handside of this equation represent the sources of atoms and thethird term reflects the losses of atoms. The source functions(r, v) corresponds to the slow atoms penetrating into theplasma from the peripheral regions. The source function g(r,v) already mentioned above is the local differential birth rateof fast atoms within the plasma column. The losses of atoms
due to charge exchange, ion impact ionization and electron
impact ionization depend on the summary total rate of theseprocesses ν = νa /λmfp.
Iterative Join: The integration of eq. (3) over the velocity space allows to obtain the neutral atom density distribution na(r). The next step of the numerical modeling is to use
this distribution to calculate the escaping neutral particle
energy spectrum. Thus, at the certain chosen free parameters,
one should obtain na(r) from the kinetic equation and thenplug it into formula (2) to calculate dN/dE. The actual valuesof the unknown free parameters are determined by iterating
this procedure until the calculated atomic spectrum matches
the measured one.
This method was successfully applied to tokamak dataprocessing on the basis of the kinetic equation in cylindricalcoordinates assuming the toroidal axial symmetry. In thisconnection the early pioneer work [4] should be mentioned.
However, in case of the arbitrary isolines the calculations
would be much more complicated. Alternatively, a Monte
Carlo simulation method such as AURORA numerical code
can be employed to calculate the neutral atom distribution.
Algorithms of this kind and their applications are explained~−λ (E, X' )dX'~
g [E, ρ (X')] × e ∫ X1'dX',−1
mfpv

(2)
2315

Goncharov Pavel R. et al., A Numerical Approach to the Localization of Passive Line Integrated Neutral Particle Measurements on LHDe.g. in [5,6].
Profile Reconstruction: Another method is possiblethat does not require any assumptions about the unknown ion
temperature profile. Neglecting the attenuation of the atomic
flux by setting λ-1
mfp = 0 in (2) yields

Γ (E, ζ ) =1

∫ ~g(E, ρ)[Q (ρ, ζ ) − Q (ρ, ζ )]dρ,+ρmin−(4)which is the generalized Abel equation for arbitrary non-circular isolines. This equation belongs to the class of incorrectly posed problems, i.e. a small perturbation on the experimental data function Γ (E, ζ) can give rise to an arbitrarily
large perturbation in the solution ∼g (E, ρ ). Methods of solving
incorrectly posed problems presently are abundant in literature. An algorithm of numerical solution of such integral
equations by Tikhonov regularization method is described in
the work [7].

Fig. 2 VMEC data points and analytical fit curves ρ(X') at Rax= 3.6 m, β = 0.22% for SDNPA detector #3.

4. Practical application
For the correct interpretation of line-integrated measurements it is essential to know the structure of the magnetic
surfaces for the given experimental conditions. It can be calculated by such a numerical magnetohydrodynamic (MHD)
equilibrium code as VMEC [3]. Thus, the local values of theeffective minor radius ρ (X') along the diagnostic sight line
can be obtained in order to use the formalism of Sec. 2.
From a practical viewpoint the profile reconstruction
method appears to be disadvantageous since it is preferable
not to neglect the attenuation. Moreover, a large data set Γ (ζ,E) is required for this method. It can be measured either during a steady state of a long plasma discharge by continuously
scanning the plasma over the angle ζ or by changing theangle ζ in a sequence of similar plasma discharges.
The integration of the kinetic equation for neutrals is notpractical because of its complexity for helical configurations.
AURORA code has been used instead.
The function ρ (X') for SDNPA sight line geometry has
been calculated from VMEC numerical data. For the practical computational treatment of (2) an analytic fit has beenfound in the form

ρ = Aj(ζ )|Bj(ζ ) − X' | r ( ζ ) + Cj(ζ ).
jFig. 3 Calculated (solid) and experimental (dots) atomic
energy distribution from ne = 4 × 1013 cm-3 NBI plasma.

ured distribution. In the calculation the experimental crosssections were used from [8]. The attenuation was assumed tobe mainly due to the secondary charge exchange and ion
impact ionization. Density data from FIR interferometer wasused. Maxwellian ion distribution was assumed for the bulk
plasma. The agreement between the calculation and theexperiment was reached for the trial Ti profile in the form
Ti(ρ ) = Ti(0)(1 − ρ 2)2 with Ti(0) 2.75 keV as it is shown inFig 3. The lower energies are out of the measurable range.
The central Ti value agrees with the crystal spectrometer data.
The high-energy part of the spectrum corresponds tosuprathermal pitch-angle scattered particles from NBI.

(5)

Here j = 1,2 corresponds to the intervals where dX'/dρ hasdifferent signs. Fig. 2 shows the examples of VMEC data
points and solid curves calculated using (5) for three different
values of the scanning angle ζ for detector #3 of SDNPA. Aquite good fitting (5) can be achieved as shown in Fig. 2, sothat the average |ρ VMEC − ρ fit| along the sight line is about
0.01 or less.
Formula (2) was used to match the experimental dN/dE
measured from NBI#1 and NBI#3 tangentially heated hydrogen plasma. Data from tangential sight line looking against
neutral beam injector (NBI) was used. The model estimation
was iteratively joined to the bulk thermal part of the meas-5. Summary
The formulation of the neutral particle flux from a nonaxisymmetrical plasma has been discussed for an arbitrary
magnetic configuration. Ti profile retrieving from NPA dataon LHD has been demonstrated.

Acknowledgements
The authors would like to express their gratitude to all
the colleagues of LHD Experimental Groups for their cooper3316

Goncharov Pavel R. et al., A Numerical Approach to the Localization of Passive Line Integrated Neutral Particle Measurements on LHDation. We are also thankful to NIFS technical personnel forthe constant assistance.

[4] A.B. Izvozchikov and M.P. Petrov, Sov. J. Plasma Phys.
2, 117 (1976).
[5] M.H. Hughes and D.E. Post, J. Comput. Phys. 28, 43(1978).
[6] S. Tamor, J. Comput. Phys. 40, 104 (1981).
[7] J. Weese, Comput. Phys. Comm. 69, 99 (1992).
[8] C.F. Barnett ed., Atomic Data for Fusion (ORNL-6086,USA, 1990).

References[1] P.R. Goncharov, J.F. Lyon et al., J. Plasma Fusion Res.
SERIES 5, 159 (2002).
[2] J.F. Lyon, P.R. Goncharov et al., Rev. Sci. Instrum. 74,1873 (2003).
[3] L.L. Lao, S.P. Hirshman et al., Phys. Fluids 24, 1431
(1981).

4317