J. Plasma
J.
Plasma Fusion
Fusion Res.
Res. SERIES,SERIES,Vol.
Vol.6 6(2004)(2004)353–356
000–000

Analysis of Pellet Ablation with Atomic Processes
ISHIZAKI Ryuichi, NAKAJIMA Noriyoshi and PARKS Paul B.1
National Institute for Fusion Science, Toki 509-5292, Japan
1
General Atomics, San Diego, CA 92121, USA(Received: 12 December 2003 / Accepted: 22 January 2004)

Abstract
A new type of a magnetohydrodynamics (MHD) code applicable to solid, liquid and gas states, “CAP” has beendeveloped in order to investigate ablation process of a pellet with atomic processes in hot plasmas. One of the mostimportant features of the code is to be able to treat recession of the pellet surface by ablation without any artificial
boundary condition between the pellet and ablation cloud. A region excluding a magnetic field is induced because themagnetic pressure is overcome by the ablation pressure. It is found that a stationary shock wave is driven by ionization.
Keywords:
pellet, ablation, shock, MHD, perfect conductor

1. Introductionpellet consisting of solid hydrogen is heated by an energy
flux, it changes to liquid and gas during phase transition andsubsequently to plasma during atomic process. In the present
work, it is assumed that the ablation cloud consists of four
species particles, namely, molecules, atoms, ions andelectrons, and it behaves as one fluid due to charge neutrality
and strong charge exchange except for incident plasma
electrons encountering the cloud. The incident electron heat
flux is considered by a kinetic treatment with a half-space
Maxwellian distribution [4]. In addition, the equation of state(EOS) [5] is used to treat solid, liquid and gas. Dissociationand ionization are also included to cause deuterium to change
to plasma. Since domain of the numerical simulation includes
solid and gas regions, the Cubic-Interpolated Psudoparticle(CIP) method [6] is used in the code. The most essential pointof the method is that it can calculate simultaneously
incompressible and compressible regions without any artificial
boundary conditions between them.

Refueling is one of essential methods in order to control
plasma density and sustain steady state plasmas. A gas puffing
has been successful for building and sustaining a plasma
density in an experimental system of past generation.
However, in a large scale experimental system, e.g., the Large
Helical Device (LHD), the plasma sources induced by the gas
puffing are strongly localized near the plasma surface. Then,a pellet injection is placed as a fundamental tool and has been
mainly used to obtain a high density plasma and control adensity profile [1].
A theoretical analysis of the pellet injection was started
by Rose [2], and several ablation models for the pellet basedon different physics were developed, e.g., a neutral gas
shielding (NGS) model set up by Parks and Turnbull [3]. The
NGS model is derived by applying one dimensional,
spherically symmetric and quasi-steady gas dynamics to thefluid ablating from a pellet in a thermonuclear plasma. And
the model supposes the ablation cloud to consist of only one
species, such as ground-state deuterium molecules. Though
the model is very simplified, it has been providing many
useful predictions and contributing to experiments of thepellet injection. However, more detailed analysis is requiredin these days when the pellet injection plays an important rolein high-β plasmas in a large scale experimental system. Our
final goal is to investigate the bulk plasma motion inducedby the pellet injection including pellet dynamics itself, bydeveloping a three dimensional fluid code dealing with botha neutral fluid and MHD plasma. As a first step, a new typeof hydrodynamic code “CAP” has been developed where avarious states of pellet are treated without assumptions.
The model of “CAP” is constructed as follows. When a

2. Basic equations
All species in the cloud are assumed to be treated as
fluids, though there is some doubt as to the validity of thisassumption in the bulk plasma. Since the ablation cloud has
a very high density and low temperature, it has a highcollisional rate required for the validity of a fluid treatment.
The hot bulk plasma, which is not accurately treated as a
fluid, has a minimal effect on the ablation process, except foran effect as the source of the heat flux which drives theablation. The ablation cloud consists of molecules, atoms, ionsand electrons. Effects of atomic physics on the ablation©2004 by The Japan Society of Plasma
Science and Nuclear Fusion Research

Corresponding author’s e-mail: ishizaki@nifs.ac.jp
1353

Ishizaki R. et al., Analysis of Pellet Ablation with Atomic Processes

process, in particular dissociation and ionization in theablation cloud, are considered. It is assumed that the species
have same velocity (u) and temperature (T) due to highcollision and charge exchange rate, so that they are regarded
as one fluid with different densities. An electromagnetic force
is ignored for simplicity. The dynamics in not only gas(ablation cloud) region but also solid (pellet) region are givenby the following equations of mass, momentum and energy
conservations:dρ= – ρ∇ ⋅ u,
dt

(1)

du
1
= – ∇p +(∇ × B) × B,
µ0
dt

(2)

de
η
= – p∇ ⋅ u + H + 2 ∇ × B 2
dt
µ0

(3)

discharges are considered so that the distribution of energies
of electrons and ions of the bulk plasma may be taken as
Maxwellian. Owing to the lower thermal velocity and shortermean free path of ions, the thermal ion heat flux is too small
to penetrate the ablation cloud and drive ablation. Since theelectron gyroradius is much smaller than the size of theablation cloud surrounding the pellet, orbits of the electrons
are strongly tied to a single magnetic field line, when theelectrons penetrate the ablation cloud. Then, the heating of
the cloud can be calculated field line by field line. The heat
flux, q, is solved by using a kinetic treatment for thedeposition energy by plasma electrons incident to the cloudwith a half-space Maxwellian distribution as shown in Ref.
[4]. When z is supposed to be the spatial co-ordinate along
the field line directed into the ablation cloud, the heat source
H in Eq. (3) is given by ∂ q/∂ z.

(4)

3. Ablation in 2D cylindrical axisymmetricsystem

ρ

ρ

⎛
⎞
∂B
η
= ∇ × ⎜u × B –∇ × B⎟∂tµ0⎝⎠

Geometrical effects are investigated by multidimensional simulation where the cylindrical geometry (r, θ,z) is used as shown in Fig. 1. A spherical pellet is placed at
the center of the geometry (r = 0, z = 0) and the heat flux
encounters the pellet along the magnetic field (along zdirection), namely, it is anisotropic heating. Since physicalquantities are assumed to be uniform in the θ-direction, thisproblem is reduced to two dimensional (2D) one. Boundary
condition on the base of the cylinder is periodic one without
flow that can go out. Boundary on the side of it is assumed
to be a perfect conductor.
Figure 2 shows temporal evolution of the contour ofmagnetic pressure: B2/(2µ0). Initial conditions are the pellet
radius rp = 2 mm, electron temperature Te∞ = 2 keV andnumber density ne∞ = 1020 m–3 in the bulk plasma. Magnetic
field is Bz = 1.12 T and Br = 0 T. Te∞ and ne∞ are assumed tobe constant through temporal evolution. When a pellet isheated, an ablation cloud with lower temperature and higher
density than the bulk plasma is constructed surrounding thepellet. The magnetic pressure in the region around the pellet
becomes very small because the magnetic pressure isovercome by ablation pressure. Lines of the magnetic force
are bent by the ablation pressure as shown in Fig. 2. Figure 3
shows temporal evolution of the contour of pressure: p. It isfound that expansion in r-direction is inhibited by thewhere ρ is the total density, u the velocity, p the total pressure,e the total specific internal energy, H the heat source, B themagnetic field, µ0 the permeability of free space and η theresistivity. When a solid is heated, it is transformed to gas(molecules) by absorbing the sublimation energy (ε s =0.01 eV). Subsequently, when the molecules are transformed
to atoms by absorbing the dissociation energy ( ε d =4.48 eV), the number density of the molecules decreases but
one of the atoms increases. When the ablation cloud is moreheated, the atoms are transformed to ions and electrons byabsorbing the ionization energy (εi = 13.6 eV) and the number
density of the atoms decreases but ones of the ions andelectrons increase. Since the ablation cloud can be regarded
as one fluid, p and e are respectively given by the following
equations:ρkT
1
1
p = ps + ⎛ fs + fd + fi ⎞
,
⎝2
⎠ m
2
⎡ f –f
f + fi ⎤ kT
e = (1 – fs ) es + ⎢ s d + d⎥⎣ 2(γ m – 1) γ – 1 ⎦ m
kε
1 kε
1 kε
+ fs s + fd d + fi i ,
m 2
m
2
m

(5)

(6)

where ps and es is the pressure and specific internal energy inthe solid, respectively, that should be determined by the EOSof the solid. γ m and γ are ratios of the specific heats formolecules and atoms, respectively. k and m are the Boltzmann
constant and nuclei mass, respectively. Fractions fs, fd and fi
are defined by fs = (2ng + na + ni)/nt, fd = (na + ni)/nt and fi =ni /nt, respectively, where nt = 2ns + 2ng + na + ni . ns, ng, na
and ni are nuclei number densities of solid, gas, atoms andions, respectively. These fractions are given by assuming the
local thermodynamic equiriblium (LTE). The atomic process
induced by the incident electrons encountering the cloud is
not considered.
In the present model, only Ohmic (no runaway)

EquatorrqθzPerfect conductor

Periodic BC w/o flow
Pole

Fig. 1

2354

2D axisymmetric cylindrical geometry. The bases of it
are periodic boundary conditions without flow that cango out. The sides of it are perfect conductor.

Ishizaki R. et al., Analysis of Pellet Ablation with Atomic Processes

µs
0.1µ

0.5µs

ram_1150_bsq_002_text

-40

µs
1.0µ

ram_1150_bsq_003_text

µs
1.5µ

ram_1150_bsq_008_text

ram_1150_bsq_013_textr (mm)

-20
0
20
40
-40

-20

0
20
z (mm)

40

0.0

0.5

1.0
2

1.5

2

B (T )

Fig. 2 Contour of B 2/(2µ0) and vector of B.

-40

µs
0.1µ

µs
0.5µ

ram_1150_pre_002_textr (mm)-40

r (mm)

-20

20

20

40
-40

40
-20

0
20
z (mm)µs
1.5µram_1150_pre_008_text

ram_1150_pre_013_text

0

0

-40

µs
1.0µ

ram_1150_pre_003_text
-20

-20

0
20
z (mm)

40

40
0.0

0.5
1.0
1.5
ram_1150_pre_002_text

P (MPa)

0.0

0.5
1.0
1.5
ram_1150_pre_002_text

Fig. 3 Contour of p.

1 02

2.0

5

ram_1150_gnu4_y3_013.text

ram_1150_gnu2_y3_013.text
1

10

0

ρ

- 1

10
1.0

ρ

p, B 2 / 2

Jθ

1 00

B2/ 2

1 0- 2

-5

Jθ

p
1.5

1 0- 3
1 0- 40.5-10

1 0- 5
0.0
-120

-80

-40

0

40

80

1 0- 6
-120

120

-80

-40

40

-15
120

80

r (mm)r (mm)
Spatial profiles of p (dashed line) and B 2/(2µ0) (solidline) in r-direction.

Spatial profiles of ρ (solid line) and Jθ (dashed line) inr-direction. respectively.

Fig. 5magnetic pressure in comparison with one in z-direction.
Figure 4 shows p and B2/(2µ0) in r-direction at 1.5 µs. The
magnetic pressure is comparable to the ablation pressure atr ~ 30 mm where force balance is satisfied. The ablation
cloud expands slowly as keeping the force balance between
the magnetic pressure and ablation pressure. Figure 5 showsdensity ρ and current Jθ. It is found that there is a surface
current at r ~ 30 mm and no current in the region. In otherwords, surface of the ablation cloud becomes perfectconductor, so that there is neither magnetic field nor current.
Figure 6 shows Mach number and fractional ionization
by dashed and solid lines, respectively. Ablated particles
around the pellet is accelerated by the ablation pressure

1.5

1.5
ram_1150_gnu3_y3_013.text

M

fi
1.0

0.5

0.5

M

1

fi

Fig. 4

0

0
-120

-80

-40

0

40

80

0.0
120

r (mm)
Fig. 6

3355

Spatial profiles of M (solid line) and fi (dashed line) inr-direction. respectively.

Ishizaki R. et al., Analysis of Pellet Ablation with Atomic Processes

without the effect of the magnetic field. On the other hand,the flow is disturbed by the magnetic field out of the ablation
cloud. Then, there are acceleration and deceleration phases in
the ablation cloud as shown in Fig. 6. A jump structure isfound to be driven in the acceleration phase. That is a shock
wave because the flow is supersonic ahead of that structureand subsonic behind it. The shock is found to be driven in
the region where the ablation particles ionize (0 < fi < 1).
Since the neutral particles loose their energy by ionization inthat region, their kinetic and internal energies are hard to theablation cloud and drives a shock.

A region excluding a magnetic field is induced by ablationpressure because it overcomes the magnetic pressure. In otherwords, there is neither magnetic field nor current in theablation cloud because surface of the ablation cloud behaves
as a perfect conductor. It is shown that a stationary shock
wave is induced by ionization.

References[1] L.R. Baylor et al., Nucl. Fusion 32, 2177 (1992).
[2] D.J. Rose, Culham Laboratory Technology Division,
Memorandum, No. 82 (1968).
[3] P.B. Parks and R.J. Turnbull, Phys. Fluids 21, 1735
(1978).
[4] P.B. Parks et al., Phys. Plasmas 7, 1968 (2000).
[5] R.L. More, Phys. Fluids 31, 3059 (1988).
[6] T. Yabe and P.Y. Wang, J. Phys. Soc. Jpn. 60, 2105
(1991).

4. Conclusions
Two dimensional MHD code “CAP”, treating with aneutral in various states of matter and plasma simultaneously,
has been developed to investigate the pellet ablation withatomic processes and anisotropic heating from bulk plasmas.

4356