A

The Astrophysical Journal, 648:435– 451, 2006 September 1# 2006. The American Astronomical Society. All rights reserved. Printed in U.S.A.

DUST DESTRUCTION IN THE HIGH-VELOCITY SHOCKS DRIVEN
BY SUPERNOVAE IN THE EARLY UNIVERSE
Takaya Nozawa1 and Takashi Kozasa1
Division of Earth and Planetary Sciences, Graduate School of Science, Hokkaido University, Sapporo 060-0810, Japan;
nozawa@ep.sci.hokudai.ac.jp, kozasa@ep.sci.hokudai.ac.jp

andAsao Habe1
Division of Physics, Graduate School of Science, Hokkaido University, Sapporo 060-0810, Japan; habe@astro1.sci.hokudai.ac.jp
Received 2006 March 11; accepted 2006 May 5

ABSTRACT
We investigate the destruction of dust grains by sputtering in the high-velocity interstellar shocks driven by
supernovae (SNe) in the early universe to reveal the dependence of the timescale of dust destruction on the gas density
nH;0 in the interstellar medium (ISM ), as well as on the progenitor mass Mpr and explosion energy E51 of SNe. The
sputtering yields for the combinations of dust and ion species of interest to us are evaluated by applying the so-called
universal relation with a slight modification. The dynamics of dust grains and their destruction by sputtering in shocks
are calculated by taking into account the size distribution of each dust species, together with the time evolution of thetemperature and density of the gas in spherically symmetric shocks. The results of the calculations show that theefficiency of dust destruction depends not only on the sputtering yield but also on the initial size distribution of eachgrain species. The efficiency of dust destruction increases with increasing E51 and /or increasing nH;0 but is almost
independent of Mpr as long as E51 is the same. The mass of gas swept up by a shock is an increasing function of E51
and a decreasing function of nH;0 . Combining these results, we present the approximation formula for the timescale ofdestruction for each grain species in the early universe as a function of E51 and nH;0 . This formula is applicable forinvestigating the evolution of dust grains at the early epoch of the universe with the metallicity of Z P 103 Z. The
effects of the cooling processes of gas on the destruction of dust are briefly discussed.
Subject headinggs: dust, extinction — early universe — shock waves — supernova remnants — supernovae: general
Online material: color figures

1. INTRODUCTION

emit it at infrared ( IR) and submillimeter wavelengths. The thermal emission from dust at high redshifts can distort the cosmic microwave background (CMB) radiation (Loeb & Haiman 1997),
and the obscuration and reddening of starlight by dust residing in
the early interstellar and/or intergalactic space lead to a serious
underestimate in evaluating the cosmic SFR from the observations
directed toward higher redshifts (e.g., Hauser & Dwek 2001).
In fact, recent observations have confirmed the presence oflarge amounts of dust grains at high redshifts of z k 5. The submillimeter (Priddey et al. 2003; Robson et al. 2004) and millimeter (Bertoldi et al. 2003) observations have reported continuum
thermal emission of dust grains from some quasars at z > 5 andhave suggested a large amount of cold dust reaching up to 108 –109 M in the host galaxies of quasars. Maiolino et al. (2004a)
discovered the substantial extinction of stellar light by dust in
broad absorption line (BAL) quasars at z > 4:9. They have shownthat their extinction curves could be different from those of lowredshift quasars at z < 4 and have suggested that the difference
may reflect different formation and evolution mechanisms of dust
grains at z k 5.
The main sources of dust grains in the early universe at z > 5
are believed to be supernovae (SNe) evolving from the massive
stars because of their short lifetime. Dust grains formed in theejecta of SNe are injected into the interstellar medium (ISM) andare subjected to destructive processes in the interstellar shocks
induced by the ambient SNe (e.g., Jones et al. 1994; Dwek et al.
1996). Thus, the size distribution and amount of dust in earlyinterstellar space are determined by the balance between their

Population III stars formed from metal-free gas clouds where
the H2 molecule is the main coolant of gas are considered to be
more massive than 102 M (Bromm & Larson 2004 and references therein), although there are still numerous uncertainties
regarding the mass range of the first stars due to a poor understanding of the relevant radiative feedback in the late accretion
phase (see, e.g., Tan & McKee 2004). Once the primordial gas isenriched by dust grains, dust has a great influence on the subsequent formation and evolution history of stars and galaxies; thestar formation rate (SFR) in the metal-poor star-forming clouds isenhanced via the formation of H2 molecules on the surface of dust
grains (Hirashita & Ferrara 2002; Cazaux & Spaans 2004). In addition, the cooling of gas by thermal radiation from dust itself
makes even the gas clouds with metallicity as low as 105 Z fragment into low-mass gas clumps (Schneider et al. 2002, 2003;
Omukai et al. 2005). Thus, dust grains in the early universe may
cause the characteristic mass scale of stars to shift from very highmass to the low mass of 1 M observed at the present time andstrongly affect the evolution of the initial mass function (IMF) atthe early epoch of the universe.
In addition, dust grains in the early universe play crucial roles
in revealing the structure and evolution of the universe from therelevant observations because they absorb stellar light and re1

Current address: Department of Cosmosciences, Graduate School of Science, Hokkaido University, Sapporo 060-0810, Japan.

435

436

NOZAWA, KOZASA, & HABE

production in SNe and their destruction in shocks, and theychange with time according to the star formation activity (e.g.,
Dwek 1998).
The amount of dust that is destroyed in shocks, as well as how
much of stellar light is absorbed and reemitted through thermal
emission by dust grains, heavily depends on their chemical composition, size, and amount. Thus, it is essential to clarify the evolution of dust in the early universe by treating the formation anddestruction processes of dust grains in a consistent manner, inorder to elucidate the SFR and the IMF during the early evolution of the universe from future observations. Nozawa et al.
(2003) investigated the chemical composition, size distribution,
and amount of dust grains formed in the ejecta of Population III
SNe. In this paper, as the second step to revealing the evolutionof dust grains in the early universe, we explore dust destructionin the high-velocity interstellar shock driven by the SN explosion, adopting the dust model by Nozawa et al. (2003) for the
initial dust residing in the early ISM.
We focus on dust destruction in a nonradiative shock where
dust grains are considered to be efficiently destroyed by nonthermal and thermal sputtering (e.g., McKee 1989) because the nonradiative shock has a high shock velocity (>100 km s1) and ahigh gas temperature (>106 K ). However, only a few sets of experimental data on sputtering yield are available for the projectiletarget combinations of astrophysical interest. Hence, we firstevaluate the sputtering yields for the combinations of dust species produced in Population III SNe and ion species of interest
to us, collecting a large amount of sputtering data and applying
the universal relation for sputtering yield proposed by Bohdansky(1984) with a slight modification so as to reproduce the experimental data well. Then we investigate the dependence of the efficiency of dust destruction for each grain species not only on thegas density nH;0 in the ISM but also on the explosion energy E51
and progenitor mass Mpr of the SN driving the interstellar shocks
by applying the hydrodynamic models of Population III SNe by
Umeda & Nomoto (2002). Finally, we evaluate the timescale ofdestruction for each grain species.
The efficiency of dust destruction by sputtering depends on
the temperature and density of the gas and the velocity of dust
relative to gas in the shock, as well as the sputtering yield and thedust size. Thus, in the calculation, we solve the time evolutionof the temperature and density of gas in spherically symmetricshocks, including the cooling of gas by the atomic process, inverse Compton scattering, and the thermal emission of dust. Together with the hydrodynamic calculations for gas, we carefully
treat the erosion of dust by thermal and nonthermal sputtering
caused by the motion of dust relative to gas, taking into accountthe size distribution of each dust species.
This paper is organized as follows: We first define the efficiencyof dust destruction in x 2. In x 3, we describe the properties of the
ISM in the early universe, together with the dust model adopted in
the calculations, and in x 4, we present the basic equations of thehydrodynamic calculations for gas and the method of the simulation for dust destruction in the shock. The sputtering yield ofeach grain species necessary for the calculations of dust destruction is evaluated in x 5, and the grain physics in the shock is described in x 6. The results of calculations of dust destruction are
given in detail for the model of interstellar shock that is driven by a
SN with E51 ¼ 1 and Mpr ¼ 20 M and is propagating into the
ISM with nH;0 ¼ 1 cm3, and the approximation formula for thetimescale of dust destruction as a function of E51 and nH;0 is derived in x 7, where the effects of the cooling of gas on the evolutionof nonradiative shock and the destruction of dust are briefly discussed. The summary is presented in x 8.

Vol. 648

2. THE DEFINITION OF THE EFFICIENCYOF DUST DESTRUCTION
Our main goal is to reveal the timescale of dust destruction byinterstellar shocks in the early universe. In general, the timescaleof destruction SN; j for grain species j is defined as (e.g., McKee
1989)
1
SN;
j ¼ jMswept
SN ;
MISM

ð1Þ

where j is the efficiency of destruction of grain species j, Mswept
is the mass of gas swept up by a shock until the shock velocity
Vshock decelerates below 100 km s1, MISM is the mass of gas
and dust in the ISM, and SN is the effective SN rate. McKee(1989) has estimated j to be nearly constant and found that
Mswept is proportional to the SN explosion energy E51 in units of1051 ergs. However, it is expected that both j and Mswept dependon not only E51 but also the progenitor mass Mpr of the SN andthe number density of gas nH;0 in the ISM. Hence, we aim atderiving the dependences of j and Mswept on E51 , Mpr , and nH;0 .
The truncation time ttr being defined as the time at which Vshock
drops down to 100 km s1, the efficiency of destruction j foreach dust species is calculated byj ¼

Md;destj (t ¼ ttr )Md;sweptj (t ¼ ttr )

;

ð2Þ

swept
where Md;destj (t) and Md; j (t) are, respectively, the mass of dust
destroyed and the mass of dust swept up by the shock, and themass of gas Mswept swept up by a spherically symmetric shock is
given by

Mswept ’40 R3shock (t ¼ ttr );3

ð3Þ

with 0 being the gas density in the ISM and Rshock (t) the travel
distance of the shock front.
3. PROPERTIES OF THE ISM IN THE EARLY UNIVERSE
3.1. Density, Temperature, and Metallicity of the ISM
We assume that the ISM in the early universe is homogeneousand stationary (v ¼ 0), referring to the results of the simulations on the formation of the early H ii regions by the first stars(Kitayama et al. 2004), which have shown that the radiative
feedback from the massive stars causes the ambient ISM to behomogenized with a mean gas density of P1 cm3. In the calculations, we consider the range of hydrogen number density in
the ISM to be 0:03 cm3  nH;0  10 cm3 , in order to reveal
the dependence of the efficiency of dust destruction on the gas
density in the ISM. The ISM surrounding the SN explosion isconsidered to be ionized by the radiation of the massive progenitor stars and heated to T0  104 K ( Kitayama et al. 2004;
Machida et al. 2005), and then T0 is expected to decrease below
104 K, where the recombination rate can be high even for a metallicity as low as Z ¼ 104 Z ( Wolfire et al. 1995). However, theevolution of the gas temperature in the postshock flow is almost
independent of T0 until the end of the Sedov-Taylor phase of thehigh-velocity shock considered in this paper, and the temperature of gas Tshock at the shock front does not decrease below several times 104 K at the truncation time, as is shown in the resultof the calculations (see x 7.1). Thus, we take T0 ¼ 104 K as therepresentative gas temperature in the ISM, independently of thedensity of the gas.

No. 1, 2006

437

DUST DESTRUCTION IN EARLY UNIVERSE

The truncation time, as well as the time evolution of gas temperature and density structures in the postshock flow, is considered to be affected by the cooling of gas depending on themetallicity of the ISM. In particular, when Tshock P105 K, theefficient line coolings of H and He ions cause the shock velocity
to decelerate rapidly and as a result regulate the truncation time.
However, the gas cooling function is almost independent of themetallicity for Z P103 Z corresponding to ½Fe/H  3, since
the metal lines do not contribute to the cooling of gas at T 
104 K (Sutherland & Dopita 1993), and we adopt the gas cooling
function for the zero metal case given in Sutherland & Dopita(1993) (see x 4.1). Note that the efficiency of dust destruction,
which is defined by the mass ratio of dust destroyed to dust swept
up by the shock, does not depend on the metallicity as long as
Z P 103 Z, despite the fact that the amount of dust in the ISM
is proportional to the metallicity, and in the calculations weassume that the metallicity of the ISM in the early universe is
Z ¼ 104 Z to investigate the destruction of dust in the ISM.
The elemental composition of gas in the ISM is given in Table 1,
where the C, N, and O metals are included so as to adjust themetallicity to Z ¼ 104 Z and their abundances are simply
evaluated by multiplying the metallicity by the primordial ratios
tabulated in Sutherland & Dopita (1993); note that these metals
play no role in the physical processes relevant to dust destruction because of the much lower abundances than those of H and
He. The electron number density ne is fixed to 1.128nH given by
Sutherland & Dopita (1993).
3.2. Model of Dust in the Early Universe
Dust grains in the ISM at the early epoch of the universe
considered in this paper are expected to be predominantly produced in the ejecta of SNe. To reveal what kind of dust grains
form in the early universe, Todini & Ferrara (2001) calculated
dust formation in primordial core-collapse supernovae (CCSNe),assuming uniform elemental composition and gas density within
the He core. They have shown that the newly formed grains areamorphous carbon with size around 300 8 and Al2O3, MgSiO3,
Mg2SiO4, and Fe3O4 with size around 10–20 8. By extending
the model of dust formation by Todini & Ferrara (2001) to pairinstability supernovae (PISNe) whose progenitors evolve frommetal-free stars with Mpr ¼ 140260 M (Umeda & Nomoto
2002; Heger & Woosley 2002), Schneider et al. (2004) found thatthe typical sizes of dust formed in PISNe range from 0.001 to
0.3 m depending on grain species. The dust models by Todini& Ferrara (2001) can successfully reproduce the extinction curve
observed toward the BAL quasar SDSS 1048+46 at z ¼ 6:2( Maiolino et al. 2004b) and the IR spectral energy distribution(SED) of the blue compact dwarf galaxy SBS 0335052( Takeuchi et al. 2003) except for a slight overestimate of the far
IR continuum; SBS 0335052 is a good target in the local universe suitable for scrutinizing the properties of dust formed in
SNe at high redshifts because of its young age (<107 yr) andlow metallicity (1/41 Z ; Vanzi et al. 2000).
However, as discussed by Kozasa et al. (1989), the species ofdust formed in the ejecta largely depend on the elemental composition in the He core, and their sizes are affected by the time
evolution of the gas density and temperature. Thus, Nozawa et al.
(2003) extensively investigated dust production in Population III
CCSNe and PISNe, by considering two extreme cases for theelemental composition in the ejecta: the unmixed case with theoriginal onion-like structure and the uniformly mixed case within
the He core. In the calculations of dust formation, the radial profile of the gas density is properly treated and the time evolution of
the gas temperature is calculated by solving the radiative transfer

TABLE 1Elemental Abundances of Gas in the ISM
Z/Z
104 ..........

AHe

AC

AN

AO

Ae

8.04 ; 102 3.62 ; 108 1.12 ; 108 2.70 ; 107 1.128

Notes.—The number abundance Ai relative to that of the hydrogen atom in
the ISM. Note that the abundances of C, N, and O are calculated by multiplying
the metallicity of Z ¼ 104 Z by the primordial ratios in Table 4 of Sutherland& Dopita (1993).

equation including the energy deposition of radioactive elements,
which are not taken into account by Todini & Ferrara (2001).
Nozawa et al. (2003) found that a variety of grain species (C, Si,
Fe, FeS, Al2O3, MgSiO3, Mg2SiO4, SiO2, and MgO) condense in
the unmixed ejecta, reflecting the difference of the elemental composition at the formation site, and that only oxide and silicate
grains (Al2O3, MgSiO3, Mg2SiO4, SiO2, and Fe3O4) form in themixed ejecta with C/O < 1. The average sizes of these grains span
the range of 0.001–1 m depending on the elemental composition
and gas density at the location of formation in the ejecta. The species of dust formed in Population III SNe and the behaviors of
their size distributions are almost independent of the progenitor
mass, but the relative mass fraction of each grain species is dependent on it.
Based on the dust models by Nozawa et al. (2003), Hirashita
et al. (2005) have shown that the extinction curve of the quasar
SDSS 1048+46 at z ¼ 6:2 can be better explained by the dust
grains produced in the unmixed Population III CCSNe thanthose in the mixed Population III CCSNe. Takeuchi et al. (2005)
have also shown by using the dust models for Mpr ¼ 20 M that
only dust grains formed in the unmixed ejecta but not in themixed ejecta can give an excellent fit to the observed SED of
SBS 0335052 over near- to far-IR wavelengths. These works
strongly suggest that the comparisons with the observations relevant to dust in the early universe prefer dust grains produced in
the unmixed Population III CCSNe rather than those in the mixed
Population III CCSNe.
Therefore, as the representative model of the initial dust grains
residing in the early ISM, we adopt the size distribution andabundance of each grain species obtained by the calculations ofdust formation in the unmixed ejecta of a Population III SN with
Mpr ¼ 20 M by Nozawa et al. (2003). For comparison, we also
consider the model of dust formed in the mixed ejecta of thesame SN. Hereafter, the models of dust formed in the unmixed
and mixed ejecta are referred to as the unmixed grain model andthe mixed grain model, respectively. The initial values of the depletion factor defined by the ratio of the total dust mass to thetotal metal mass produced in the SN are 0.224 for the unmixed
grain model and 0.297 for the mixed grain model. The bulk density d; j , mass fraction Aj , and average size aave; j of each grain
species are summarized in Table 2, and their size spectra are
given in Figure 6a for the unmixed grain model and in Figure 7afor the mixed grain model, where the bin width is 0.1 dex.
It is useful for the following discussion to address here thebehavior of the size distribution of each dust species producedin the unmixed ejecta. C, Fe, and SiO2 grains have a lognormallike size distribution with the relatively large average radius of>0.01 m, while the size distribution functions of FeS, Mg2SiO4,
and MgO grains are approximately power law. The average radiiof Al2O3 and MgSiO3 grains are very small (<5 ; 103 m). Notethat the Si grain shows a bimodal size distribution: the grains withsize less than 0.03 m occupy almost 100% in number but only
0.2% in mass. Thus, as the average radius of the Si grain, we take

438

NOZAWA, KOZASA, & HABE
TABLE 2
Dust Models Used in the Calculations
Unmixed Grain Model

Mixed Grain Model

Dust Species

d; j(g cm3)

Aj

aave; j(m)

Ajaave; j(m)

C....................
Si ...................
Fe...................
FeS ................
Al2O3 .............
MgSiO3 .........
Mg2SiO4 ........
SiO2 ...............
MgO ..............
Fe3O4 .............

2.26
2.32
7.89
4.84
3.98
3.18
3.20
2.64
3.56
5.21

0.103
0.304
0.050
0.156
0.002
0.010
0.219
0.100
0.056
...

0.029
0.250
0.148
0.00880.00070.00450.00660.047
0.0058
...

...
...
...
...
0.001
0.102
0.205
0.588
...
0.103

...
...
...
...
0.0003
0.011
0.010
0.033
...
0.0022

Notes.—The models of dust grains are adopted from the results of the calculations of dust formation in the ejecta of a Population III SN with Mpr ¼ 20 M by
Nozawa et al. (2003). The bulk density of grain species j is denoted by d; j. The
mass fraction Aj and the average size aave; j of each grain species are given for
the unmixed and mixed grain models. Note that the average radius of Si grains in theunmixed grain model is calculated for the grains larger than 0.03 m, which occupy
99.8% of the total mass.

0.25 m calculated for the grains larger than 0.03 m including
99.8% of the mass. The size distribution function summed up overall grain species is well fitted with a broken power-law formula
whose index is 3.5 for the radius larger than 0.06 m and 2.5
for the smaller one.
4. PHYSICS OF INTERSTELLAR SHOCK
The deceleration rate due to drag force and the erosion rate
due to sputtering of dust grains depend on the time evolution ofgas temperature and density in the postshock flow. In this section, we present the basic equations and initial conditions for thehydrodynamic calculation of the time evolution of the shock andthe method of the simulation for dust destruction in the shock.
4.1. The Basic Equations for the Gas Dynamics
The hydrodynamic equations for gas to calculate the time evolution of a spherically symmetric interstellar shock driven by a
SN explosion, using conventional physical variables, are given by
@
1 @ 2
þ
(r v) ¼ 0;
@t r 2 @r@
1 @ 2 2@P(v) þ 2(r v ) ¼ 
;
@t
r @r
@r

@U
1 @  2þ 2
r (U þ P)v ¼ (nH ; T );
@t
r @r

ð4Þ
ð5Þ

1 2
P
v þ
2
1

ic (T ) ¼ 5:41 ; 1032 (1 þ z)4 ne T4 ;

ð9Þ

where T4 is the gas temperature in units of 104 K, and we takez ¼ 20 in the calculation. The cooling function d (nH ; T ) throughthe thermal emission from dust heated by collisions with gas
particles in the shock is calculated by ( Dwek 1987)
ZX
d (nH ; T ) ¼
nd; j Hj (a; T ; nH ) fj (a) da;ð10Þjwhere nd; j and fj (a) are, respectively, the number density andsize distribution function of dust species j at a given t and r. The
collisional heating rate Hj (a; T ; nH ) of grain species j is given inx 6.3.
The time evolution of gas temperature and density in a shock
driven by a SN is characterized by the explosion energy E51 andthe progenitor mass Mpr , which specify the structures of densityand velocity of gas in the ejecta. Thus, as the initial condition forinterstellar shocks, we adopt the hydrodynamic models of theejecta of Population III SNe with various E51 and Mpr by Umeda& Nomoto (2002) to investigate the dependence of the efficiencyof dust destruction on E51 and Mpr . Note that in the early universe with Z ¼ 104 Z, the progenitors of SNe are not only
Population III stars but also massive stars of the subsequent generation. However, the structures of density and velocity in a SN
with a given E51 and Mpr do not depend very much on the initial metallicity of the progenitor star if its metallicity is less than
104 Z (H. Umeda 2002, private communication). In the calculations, we consider ordinary CCSNe with an explosion energyof E51 ¼ 1, hypernovae (HNe) whose explosion energy is morethan 10 times that of CCSNe, and PISNe whose progenitors arevery massive (Mpr ¼ 150, 170, and 200 M ); hereafter CCSNe
and HNe with progenitor mass of Mpr ¼ 13Y30 M are referred toas Type II SNe (SNe II). The ejecta models used in the calculations are given in Table 3, where the labels C, H, and P represent
CCSNe, HNe, and PISNe, respectively, and the numerical value
denotes the mass of the progenitor in units of solar mass.
4.2. The Method of the Simulations

ð7Þ

with adiabatic index , and we adopt  ¼ 5/3 in the calculation.
The cooling rate of gas (nH ; T ) in units of ergs s1 cm3 isexpressed by
(nH ; T ) ¼ ne nH gas (T ) þ ic (T ) þ d (nH ; T );where each term represents the rate of cooling by the atomic
process, by inverse Compton scattering, and by the thermal
emission of dust. In the calculations, we adopt the gas cooling
function gas (T ) for the zero metal case given in Sutherland &Dopita (1993) and ignore the contribution of cooling by metal
ions ejected from dust by sputtering because, as mentioned in
x 3.1, the metals do not contribute to the cooling of gas for a
metallicity of Z ¼ 104 Z (Sutherland & Dopita 1993). For the
inverse Compton scattering caused by collisions of hot electrons in the shock with CMB photons, the cooling rate ic (T ) is
given by ( Ikeuchi & Ostriker 1986)

ð6Þ

where P ¼ kT /mH is the gas pressure with mean molecularweight  and hydrogen mass mH , and the sum of kinetic energy
and thermal energy per unit volume isU¼

Vol. 648

ð8Þ

Equations (4)–(6) are solved by using the flux-splitting methodwith second-order accuracy in space and first-order accuracy intime (van Albada et al. 1982; Mair et al. 1988). This algorithm is
one of the upwind schemes for the Euler equations and is well
suited to solving problems involving a shock.
The calculations for the time evolution of the shock and dust
destruction are performed via the following two steps. In the firststep, we assume that the SN ejecta collides with the ISM in 10 yr
after the explosion and calculate the time evolution of gas density and temperature without including dust, based on the ejecta
models of Umeda & Nomoto (2002). The spatial mesh number
is 500–5000 depending on the gas density in the ISM, and theNo. 1, 2006

DUST DESTRUCTION IN EARLY UNIVERSETABLE 3Models of Supernova Ejecta

Model

Progenitor Mass(M )Explosion Energy E51(1051 ergs)

C13.................................
C20.................................
C25.................................
C30.................................
H25.................................
H30.................................
P150 ...............................
P170 ...............................
P200 ...............................

13
20
25
30
25
30
150
170
200

1
1
1
1
10
30
1
20
28

Notes.—The models of the SN ejecta are adopted from Umeda & Nomoto(2002). In the model names, the labels C, H, and P represent CCSNe, HNe, andPISNe, respectively, and the numerical value denotes the mass of the progenitor
in units of solar masses.

inner 20 mesh points are assigned to the ejecta. Note that dust
destruction in this first step does not contribute to the evaluationof the destruction efficiency because the mass of gas swept up by
the shock is less than 0.1% of that by t ¼ ttr . In addition, the cooling of gas by dust has no effects on the subsequent structure of
the shock. In the second step, fixing the spatial mesh number to430, we rearrange the density, temperature, and velocity of thepostshock gas obtained from the first step into the inner 30 linear
mesh points. Then, including dust, we calculate the time evolution of the shock, as well as the dynamics and destruction of thedust.
By treating dust as a test particle, the physical processes relevant to dust are calculated as follows. Consider a stationary(vd ¼ 0) dust grain characterized by its composition and number
density with radius between a and a þ da at a given position Rd
in the ISM. Once the grain enters the shock front, the relative velocity to gas wd ¼ v  vd is assigned. Then for the temperature
T and density nH of gas at the shock front, we calculate the deceleration rate dwd /dt due to drag force, the erosion rate da/dt by
sputtering, and the heating rate H(a; T ; nH ) by collisions withgas. By updating the relative velocity by wd0 ¼ wd þ (dwd /dt)t
and the grain radius by a 0 ¼ a þ (da/dt)t for a given small time
step t, the velocity and position of the dust at time t 0 ¼ t þ t
are calculated by vd0 ¼ v þ wd0 and Rd0 ¼ Rd þ vd0 t. The procedure described above is repeated for the gas temperature T anddensity nH at the position Rd0 of the grain. Note that the dust grain
is defined to be completely destroyed if the grain radius is smaller
than the nominal monomer radius of condensate given in Nozawa
et al. (2003).
5. EVALUATION OF SPUTTERING YIELD
Sputtering is divided into two categories reflecting the difference in the underlying mechanisms: physical sputtering andchemical sputtering. Physical sputtering invokes a transfer of kinetic energy from the incident particle to target atoms and subsequent ejection of those atoms that have acquired enough kinetic
energy to overcome the binding forces exerted by the target.
Chemical sputtering invokes a chemical reaction induced by theimpinging particles that produces an unstable compound at thetarget surface (Sigmund 1981). In the high-velocity shocks considered in this paper, dust grains are predominantly destroyedby physical sputtering because the dust temperature is at most
100 K in the postshock region ( Dwek 1987); chemical sputtering
is of no importance, since this process occurs at a temperature as

439

high as 800 K (Roth 1983). Thus, in what follows, we simply
refer to physical sputtering as ‘‘sputtering.’’ The efficiency of destruction by sputtering is characterized by the sputtering yielddefined as the mean number of emitted atoms per incident particle. The sputtering yield depends on the impact energy of a projectile and the incident angle to the target surface, as well as on
the projectile-target combination. However, only a few sets ofexperimental data on sputtering yield are available for the combinations of dust materials and ion species of interest to us.
One method of evaluating the sputtering yields for the projectiletarget combinations of astrophysical interest is to apply the Monte
Carlo code for the transport of ions in matter (TRIM; Ziegler et al.
1985) that simulates the energy loss of a projectile ion with given
impact energy and incident angle in a solid (Field et al. 1997; May
et al. 2000; Bianchi & Ferrara 2005). The TRIM code requires as
input parameters not only the surface binding energy of the target
but also the displacement energy Ed and the lattice binding energy
Eb . Particularly, uncertainties in the value of Ed are a significant
source of error in the calculated sputtering yield (May et al. 2000;
Bianchi & Ferrara 2005). In fact, the sputtering yields computedby this simulation code have been known to significantly differ
from the experimental data at low impact energies (Bianchi &Ferrara 2005). Furthermore, it is time consuming to combine theMonte Carlo simulation with the hydrodynamic calculation for
dust destruction.
Another method is to apply the universal relation, which isthe analytic formula derived by Bohdansky (1984) for the energy
dependence of sputtering yield at normal incidence for a monatomic solid. To investigate the destruction of interstellar dust
such as C, Fe, SiC, SiO2, and H2O ice, Tielens et al. (1994) haveapplied the universal relation by treating the constant K in thisformula as a parameter, and they determined it by comparing theavailable experimental data with the sputtering yield calculatedby the universal relation. They have shown that the universal relation can be applied to not only monatomic solids but also compound solids by taking the appropriate average values for theatomic and mass numbers and the surface binding energy under
the assumption that each element comprising the target is sputtered off at the same rate. However, the agreement of the calculated sputtering yield with the experimental data is not very goodfor the ion-target combinations with a mass ratio of target atom
to incident ion ranging from 0.7 to 2.
Therefore, we apply the universal relation for sputtering yieldwith a slight modification so as to reproduce the experimental
data well, collecting and minutely examining a large amount of
sputtering data, in order to evaluate the sputtering yields for the
combinations of dust and ion species of interest to us. Following
Tielens et al. (1994), we treat K as a free parameter and determine
it by fitting to the experimental data for the grain species with theyield data covering a wide range of impact energy of the projectiles. For the grain species for which little or no sputtering data
are available, we deduce the value of K by fitting to the sputtering yield calculated by means of the Monte Carlo code for the erosion and deposition based on a dynamic model (EDDY; Ohya &Kawata 1997). The EDDY code requires only the surface binding
energy as an input parameter and allows us to predict the sputtering yield for a specific impact energy and incident angle moreeasily than the TRIM code, although the energy range applicableto the calculations is limited to 0.1–10 keV. In the following subsection, we first introduce the universal relation for sputtering
yield and then determine the value of K for each grain species
by fitting to the available experimental data and /or the resultsof the EDDY simulation.

440

NOZAWA, KOZASA, & HABE
5.1. Universal Relation for Sputtering Yield

For backward sputtering and with the neglect of inelastic
energy losses, the sputtering yield at normal incidence Yi 0 (E) by
the projectile i impacting with energy E is given by (Bohdansky
1984)"   #
Eth 2=3Eth 2
014 Si (E ) i ( i )11;Yi (E ) ¼ 4:2 ; 10
U0 Ki þ 1
EEð11Þ
where U0 is the surface binding energy in units of eV, i ¼
Md /Mi is the ratio of the mass number of the target atom to theincident ion, i is the energy-independent function of i , and K
is treated as a free parameter to be adjusted to reproduce theexperimental data ( Tielens et al. 1994). The threshold energies
Eth have been obtained by fitting the yield data for low-energy
sputtering and are approximately given by ( Bohdansky et al.
1980; Andersen & Bay 1981)
8
U0
Mi
>
>
for
 0:3;><
gi (1  gi )
Md
ð12Þ
Eth ¼
 1=3
>
Mi
Mi
>
>
for
> 0:3;: 8U0
Md
Md
where gi ¼ 4Mi Md /(Mi þ Md )2 is the maximum fractional energy
transfer possible in a head-on elastic collision. The function Si (E )is the nuclear stopping cross section in units of ergs cm2 and canbe expressed by the following universal relation (Sigmund 1981):
Si (E ) ¼ 4asc Zi Zd e 2Misi (i );Mi þ Md

ð13Þ

where e is the elementary charge and Zi and Zd are the atomic
numbers of the projectile and the target, respectively. The screening length asc for the interaction potential between the nuclei isgiven by

1=2
2=3
2=3
ð14Þ
asc ¼ 0:885a0 Zi þ Zd
with the Bohr radius a0 ¼ 0:529 8. The function si (i ) depends
on the detailed form adopted for the screened Coulomb interaction and can be approximated by ( Matsunami et al. 1980)
pﬃﬃﬃﬃ
3:441 i ln (i þ 2:718)si (i ) ¼
pﬃﬃﬃﬃ
pﬃﬃﬃﬃ ; ð15Þ1 þ 6:35 i þ i 1:708 þ 6:882 i

TABLE 4
Parameters for Calculations of Sputtering Yield
and References for Yield Data

Dust Species

U0
(eV )

Zd

Md

K

References

C............................
Si ...........................
Fe...........................
FeS ........................
Al2O3 .....................
MgSiO3 .................
Mg2SiO4 ................
SiO2 .......................
MgO ......................
Fe3O4 .....................

4.0
4.66
4.31
4.12
8.5
6.0
5.8
6.42
5.17
4.98

6
14
26
21
10
10
10
10
10
15.7

12
28
56
44
20.4
20
20
20
20
33.1

0.61
0.43
0.23
0.18
0.08
0.1
0.1
0.1
0.06
0.15

1, 2, 3, 4, 5
3, 4, 6, 7, 8, 9, 10, 11
3, 4, 5, 7, 8, 12, 13
EDDY
3, 14, 1516
16
3, 14, 15, 17, 18
14, 15, EDDY
EDDY

Notes.—The surface binding energy is represented by U0 . The mean atomic
number and mean mass number of the target material are Zd and Md , respectively, and the free parameter K is determined by fitting to the available experimental data and/or the results of the EDDY simulation. See text for details.
References.—(1) Bohdansky et al. 1978; (2) Bohdansky et al. 1976; (3) Rothet al. 1979; (4) Rosenberg & Wehner 1962; (5) Hechtl et al. 1981; (6) Laegreid &Wehner 1961; (7) Southern et al. 1963; (8) Eernisse 1971; (9) Sommerfeldt et al.
1972; (10) Coburn et al. 1977; (11) Blank & Wittmaack 1979; (12) von Seefeld
et al. 1976; (13) Bohdansky et al. 1977; (14) Bach 1970; (15) Nenadovic et al. 1990;(16) Tielens et al. 1994; (17) Bach et al. 1974; (18) Edwin 1973.

However, equation (17) for i slightly overestimates the sputtering yields for the projectile-target combinations with a mass
ratio of 0:7 P i P 2. Hence, from a comparison with a large
amount of sputtering data, we apply the following formula for i :8for i  0:5;
>
< 0:2
2
1
0:1
þ
0:25
ð


0:5
Þ
for 0:5 < i  1;
i ¼
i
i
>
:
2=3
0:3ði  0:6Þfor 1 < i :ð18ÞAs we show below, this modified function i can produce excellent agreement with the sputtering data for any values of i
considered here (0:3  i  56).
By treating K as a parameter, the surface binding energy U0must be specified to evaluate the normal-incidence sputtering
yield Yi 0 (E ) from equation (11). We assume U0 to be equal to thesublimation energy and evaluate it from the JANAF Thermochemical Tables (Chase et al. 1985) except for C and Al2O3, whose
values of U0 are discussed below. The values of U0 for each grain
species are summarized in Table 4 along with the values of Zd
and Md .

where the reduced energy i is given by
i ¼

Md
asc
E:
Mi þ Md Zi Zd e 2

5.2. Sputtering Yield of Each Grain Species
ð16Þ

The value of i (i ) in equation (11) depends on how the distribution of energy deposited in the target is approximated. Sigmund(1969) suggested the three types of approximations according to
the adopted distribution function (Gaussian, corrected Gaussian,and non-Gaussian). The values of i (i ) obtained by these approximations are very similar (Sigmund 1969) and are approximately expressed for i  10 as follows (e.g., Bohdansky 1984):(
0:2
for i  0:5;
i ¼
ð17Þ
2=3
0:3i
for 0:5 < i  10:

In order to evaluate the sputtering yield for the projectiletarget combinations of interest to us, we determine the value of K
by fitting to the experimental data for C, Si, Fe, SiO2, and Al2O3
grains for which a sufficient amount of sputtering data are available. For MgO, FeS, and Fe3O4 grains for which little or no yielddata are available, we evaluate K by fitting to the results of the
EDDY simulation. The experimental and simulated sputtering
data at normal incidence, for which the references are summarizedin Table 4, are given in Figure 1 along with the best-fitting theoretical yield calculated by the universal relation (solid lines). In
Figures 1a–1e, the experimental data are shown for the following
ion species: H+ (open circles), D+ (asterisks), He+ (open squares),
Ne+ (crosses), and Ar+ (open triangles). Even in the case that theFig. 1.—Sputtering yields of each grain species vs. incident energy of projectiles. (a) C, (b) Si, (c) Fe, (d ) SiO2, (e) Al2O3, ( f ) MgO, ( g) FeS, and (h) Fe3O4. The
incident ion species are H+, D+, He+, Ne+, and Ar+. The experimental data on sputtering yield are represented by open circles ( H+), asterisks ( D+), open squares( He+), crosses ( Ne+), and open triangles (Ar+), and the results of the sputtering yield calculated by the EDDY code are denoted by filled circles ( H+), filled squares( He+), and filled triangles (Ar+). The solid curves show the best-fitting theoretical yields calculated by the universal relation. In (b) and (c), the dotted curves show
the theoretical sputtering yields by Ar+ projectile calculated by adopting i given by eq. (17). [See the electronic edition of the Journal for a color version of thisfigure.]442

NOZAWA, KOZASA, & HABE

experimental data are sufficiently available, for comparison, we
also plot the results of the EDDY calculations for the projectiles
H+ ( filled circles), He+ ( filled squares), and Ar+ ( filled triangles).
As for C, the sublimation energy estimated to be U0 ¼ 7:43 eVfrom the thermodynamical data cannot reproduce the experimental data. Tielens et al. (1994) have shown that the ion fluences
efficiently amorphize the surface layer of carbon material and reduce the surface binding energy. Thus, we take U0 ¼ 4:0 eV following Tielens et al. (1994), and this value gives a good fit to themeasured yields with K ¼ 0:61 (Fig. 1a) being different from
K ¼ 0:65 by Tielens et al. (1994); the difference stems from theadopted formula for i. The sputtering yield calculated by the
EDDY code is enhanced for H+ by a factor of 2–3 compared tothe experimental data but gives a good agreement for He+.
The sputtering data for the Si and Fe targets are sufficiently
available. The modified universal relation employing equation (18)for i shows excellent agreement with the data points for all
projectiles considered here (Figs. 1b and 1c), by adopting K ¼
0:43 and 0.23, respectively. The EDDY calculations for Si andFe targets also present good agreement with their measured data.
The dotted lines indicate the sputtering yields by Ar+ projectile
calculated by using the values of K derived above but adopting
i given by equation (17). It can be seen from these figures that
equation (17) overestimates the sputtering yields for the projectiletarget combinations with a mass ratio of 0:7 P i P 2 as mentioned in x 5.1. It should be emphasized that the previous formulacannot reproduce the experimental data better than the modified
universal relation whatever values of K are selected.
A very good fit of the universal relation to the experimental
data is obtained for the SiO2 target by taking K ¼ 0:1 ( Fig. 1d ).
This justifies the application of the universal relation to compound materials as shown by Tielens et al. (1994). The results ofthe EDDY simulations are consistent with the experimental data,
but a comparison with the theoretical curves of sputtering yieldshows a little enhancement for the H+ and Ar+ projectile at E 
200 eV. For the MgSiO3 and Mg2SiO4 materials for which sputtering data cannot be found in the literature, as in Tielens et al.
(1994), we assume K to be 0.1, considering SiO2 to be representative of silicates. In fact, this could be reasonable, since K ’ 0:1
for the compounds including oxygen atoms (see Table 4).
Given the surface binding energy of U0 ¼ 6:37 eV evaluatedfrom the thermodynamical data, the yield data of Al2O3 cannot
be reproduced by the universal relation for any values of K. The
disagreement of the theoretical prediction with the yield data is
expected to be due to an unsuitable estimate of U0 , as is the casefor C. Indeed, Roth et al. (1979) have suggested U0 ¼ 8:5 eV for
Al2O3 from the behavior of the sputtering yield at impact energies below the maximum yield by He+. Adopting this surface
binding energy, we can get a better fit with K ¼ 0:08 (Fig. 1e),although the agreement of the universal relation with the yielddata is the worst among grain materials considered in this paper.
To our knowledge, little or no experimental data on sputtering
yield exist for MgO, FeS, and Fe3O4. Thus, we employ the results by the EDDY calculations to extract the values of K forthese materials. In Figures 1f–1h (1f for MgO, 1g for FeS, and1h for Fe3O4), we show the EDDY results and the best-fitting
theoretical curves. The agreement is not always good for H+ andAr+ at low energies, which is also true for SiO2, but the universal
relation reproduces the yield data of SiO2 well. Thus, we believe
the fits to be reasonable and satisfactory for these materials. The
values of K determined are 0.06, 0.18, and 0.15 for MgO, FeS,
and Fe3O4, respectively.
Table 4 summarizes the value of K for each grain species evaluated by fitting to the available experimental data and/or the re-Vol. 648

sults of the EDDY simulations. We confirm that the derived valueof K can well reproduce the experimental data of compounds, aswell as single-element materials, by adopting appropriate averagevalues for the atomic and mass numbers as proposed by Tielens
et al. (1994). By slightly modifying the function i (i ) in the universal relation, we realized that the calculated sputtering yields
agree better with the available experimental data than before. The
modified universal relation derived here could be applicable to evaluate the sputtering yield for any combinations of targets and projectiles of astrophysical interest as long as the experimental data
are available.
6. PHYSICS OF DUST GRAINS IN SHOCK
The collisional interactions between dust and a hot gas in shock
efficiently transfer charge, momentum, and energy between thetwo phases. As a result, dust grains acquire electric charge, undergo resistance to their motion through the gas, are eroded by
sputtering, and are heated by collisions with gas. In addition,
heated dust grains cool the gas through their thermal emission.
The electric charge on a grain plays an important role in grain
physics in hot plasmas through the Coulomb interaction and theLorentz force in the presence of magnetic field. Here we concentrate on the destruction of dust grains in fast (100 km s1)and hot (105 K ) nonradiative shocks. For T k 106 K, the effectof the electric charge of dust grains can be neglected because thedimensionless potential parameter  defined by the ratio of theelectric potential acquired by dust to the kinetic energy of gas isapproximated by   105/T and is significantly smaller than unity(Draine & Salpeter 1979; McKee et al. 1987). In addition, theerosion rate of refractory grains by thermal sputtering quickly decreases for T < 106 K (Draine & Salpeter 1979). Correspondingly, the erosion rate of dust grains by nonthermal sputtering
decreases quickly for the relative velocity of dust to gas wd <200 km s1. Therefore, it can be expected that the electric charge
does not significantly affect the destruction of dust grains in thenonradiative shocks considered in this paper. Also, it is assumed
that the early ISM is not pervaded by a magnetic field because themagnetic field in the early universe is considered to be very weak
(Gnedin et al. 2000). Thus, we formulate the basic equations describing the motion of dust grains and their erosion by sputteringin postshock flow, neglecting the effect of charge on dust grains.
6.1. Dynamics of Dust Grains
Once dust grains that are at rest in the preshock ISM encounterthe interstellar shock, they move ballistically behind the shock
front with the initial velocity relative to gas wd ’ 34 Vshock , where
Vshock is the shock velocity, and are eroded by nonthermal sputtering. A grain streaming through the ionized gas is decelerated
by a drag force due to the direct collisions with gas particles, andthus the velocity of dust relative to gas decreases. The drag force
acting on a grain results from the momentum transfer from gas todust. Under the assumption that a gas particle reflects specularly
from the grain surface and that the kinetic temperature of gas is
the same for all gas species, the deceleration rate of a sphericalgrain of radius a with velocity relative to gas wd is given by (e.g.,
Draine & Salpeter 1979)
dwd
3nH kT X¼
Ai Gi (si );2ad idt

ð19Þ

where nH is the number density of hydrogen atoms, k is theBoltzmann constant, T is the gas temperature, d is the bulk density of the grains, and si2 ¼ mi wd2 /2kT , with mi being the mass of

No. 1, 2006

443

DUST DESTRUCTION IN EARLY UNIVERSE

gas species i. The summation is taken over all gas species, andAi is the number abundance of gas species i relative to that ofhydrogen atoms. For the function Gi (si ), whose exact formula isgiven by Baines et al. (1965), the following analytical approximation has been proposed by Draine & Salpeter (1979):8si
9 2 1=2si;
ð20Þ
Gi (si )  pﬃﬃﬃ 1 þ
64
3 with the accuracy within 1% for 0 < si < 1. In the calculations, we adopt this formula. It can be seen from equation (19)that grains with smaller size and/or higher relative velocity canbe effectively decelerated by drag force and come to comove
with gas in a shock.
6.2. Grain Destruction by Sputtering
The destruction of dust grains in high-velocity and hightemperature nonradiative shocks is dominated by sputtering. Actually, we can neglect the thermal evaporation of dust because
the temperature of dust cannot rise as high as its sublimation temperature (1000–2000 K ) for the preshock gas density considered in this paper (e.g., Dwek 1987). Furthermore, we neglect
the shattering process and partial vaporization caused by graingrain collisions, since the number density of dust particles is verysmall in postshock flow and collisions between dust grains are
extremely rare without a magnetic field.
The erosion rate of grains by sputtering is calculated by taking
the angle-averaged sputtering yield hYi (Ei )i ¼ 2Yi 0 (Ei ) (Draine& Salpeter 1979) and is given by (e.g., Dwek et al. 1996)
X  8kT 1=2 esi2
msp
da
¼
nH
Ai
dt
mi
2d
2si
i
Z
pﬃﬃﬃﬃ ipﬃﬃﬃﬃ;
i e sinh (2si i )Yi 0 (i ) di ;ð21Þwhere i ¼ Ei /kT and msp is the average mass of the sputteredatoms. Equation (21) can be reduced to the following equations
for the two extreme cases of si ( Dwek & Arendt 1992; Tielens
et al. 1994): One is for si T1, where a stationary grain suffers
the thermal sputtering caused by collisions originating from thethermal motion of gas, and the erosion rate is written as

 Z
msp X
1 da
8kT 1=2¼i ei Yi 0 (i ) di :
Ai
ð22Þ
nH dt
mi
2d i
The other is in the limit of si ! 1, where the nonthermal sputtering erodes a hypersonic grain with the rate given by
X
msp
1 da
¼
wd
Ai Yi 0 (E ¼ 0:5mi wd2 ):
nH dt
2d
i

ð23Þ

Figure 2 shows the erosion rate of each dust species by sputtering in a gas with the primordial elemental composition of
metallicity Z ¼ 104 Z given in Table 1; Figure 2a shows thermal sputtering versus the gas temperature T from equation (22),and Figure 2b shows nonthermal sputtering versus the velocityof dust relative to gas wd from equation (23). The erosion rateby thermal (nonthermal) sputtering steeply increases from T 
105 K (wd  30 km s1), reaches a peak at T ¼ (4 Y20) ; 107 K(wd ¼ 500Y1300 km s1), and then slowly decreases with increasing T (wd ). This behavior of the erosion rate as a functionof T and wd reflects the dependence of sputtering yield on im-

Fig. 2.—Erosion rate of each dust species by sputtering in units of m yr1 cm3
calculated for the elemental composition of gas with Z ¼ 104 Z given in Table 1.
(a) Thermal sputtering calculated by eq. (22) as a function of gas temperature T;(b) nonthermal sputtering calculated by eq. (23) as a function of velocity of dust
relative to gas wd . [See the electronic edition of the Journal for a color version ofthis figure.]

pact energy (see x 5.2). Among dust species considered in thispaper, C grains have the lowest erosion rate at T  2 ; 106 K(wd  200 km s1), which is about 1 order of magnitude lower
than that of FeS grains with the highest rate at T  107 K (wd 
400 km s1). For thermal sputtering, the erosion rate da/dt forthe other dust species is 1:2 ; 106 nH m yr1 cm3 within afactor of 3 at T  2 ; 106 K. For nonthermal sputtering at wd 
200 km s1, the destruction by He+ is dominant, and the erosion
rate of Al2O3 is the lowest.
6.3. Grain Heating
The collisions with gas can also heat dust grains in postshock
flow. The collisional heating rate H(a; T; nH ) of a dust grain withradius a is presented by Dwek & Arendt (1992) as
H(a; T ; nH ) ¼ nH a 2 kT
;

Z

X
i


Ai


2
8kT 1=2 esi
mi
2si

pﬃﬃﬃﬃ
3=2
i ei sinh (2si i ) i (a; i ) di ;

ð24Þ

444

NOZAWA, KOZASA, & HABE

Vol. 648

Fig. 4.—Time evolution of shock velocity Vshock and gas temperature Tshockat the shock front for model C20 with nH;0 ¼ 1 cm3. [See the electronic editionof the Journal for a color version of this figure.]

7.1. Dust Destruction in Interstellar Shock for C20 Model
7.1.1. Time Evolution of Temperature and Density of Gas in ShockFig. 3.—Structures of density (top) and temperature of gas (bottom) at given
times in the interstellar shock that is driven by the SN model of C20 with
E51 ¼ 1 and Mpr ¼ 20 M and is propagating into the ISM with nH;0 ¼ 1 cm3.
The arrow in the bottom panel indicates the position of the shock front at a given
time. [See the electronic edition of the Journal for a color version of this figure.]

where i (a; i ) is the fraction of kinetic energy of gas species i
deposited into the dust grain. The values of i (a; i ) are calculated by comparing the effective grain thickness aeA (4a/3 for aspherical grain) with the stopping range of the incident particle
ls . When aeA  ls , the incident particle deposits almost all of theenergy into dust grains; otherwise, the difference between aeA
and ls is exploited to determine the fraction of deposited energy.
The approximations for the stopping range of ion species considered here are given in Dwek (1987). For the electron stopping range, we adopt the approximation derived by Tabata et al.
(1972) and Iskef et al. (1983). This electron stopping range reduces the heating rate by 30% at 106 K P T P 109 K compared with that given by Dwek (1987). The values and detailed
derivation of i (a; E ) will be presented elsewhere ( T. Nozawa
et al. 2006, in preparation).
7. RESULTS OF CALCULATIONOF DUST DESTRUCTION
In this section, we show the results of numerical simulations
for dust destruction in a nonradiative shock, focusing on the results calculated for the ISM dust specified by the unmixed grainmodel described in x 3. In x 7.1, we demonstrate the time evolution of the temperature and density of the gas and the destructionof dust grains in the interstellar shock that is driven by the SN
model of C20 with E51 ¼ 1 and Mpr ¼ 20 M and is propagating
into the ISM with nH;0 ¼ 1 cm3. In x 7.2, we investigate thedependences of the destruction efficiency and the mass of the gas
swept up by the shock on E51 , Mpr , and nH;0 and derive an approximation formula for the timescale of destruction for each dust
species in the early universe. The effects of the cooling processesof gas on the evolution of the nonradiative shock and the destruction of dust are discussed briefly in x 7.3.

Figure 3 shows the structures of density (top) and temperatureof gas (bottom) in a nonradiative shock at given times, and thesolid lines in Figure 4 give the time evolution of the shock velocity Vshock and the gas temperature Tshock at the shock front. As
shown in Figure 3, the gas density increases to 4 times that in the
ISM and the gas temperature rises steeply at the shock front,
which is indicated by the arrow for a given time step. With theinitial shock velocity ’6000 km s1, the gas temperature at theshock front remains above 108 K until t ’ 200 yr, and both Vshock
and Tshock decrease with time. At t > 5 ; 104 yr after the explosion, the density of the gas in the inner postshock region is1 order of magnitude lower than that in the unshocked region,
and a low-density hot bubble with a temperature of several times
106 K is formed (Fig. 3). In this model, Vshock decelerates below
100 km s1 at ttr ’ 105 yr when the shock front travels the distance of ’30 pc and Tshock drops to several times 104 K. After that,a dense thin shell forms at the shock front because of the effective
line cooling, and the supernova remnant (SNR) enters the radiative phase.
7.1.2. Dynamics and Destruction of Dust Grains in ShockThe motion of dust grains (Fig. 5, top) and the evolution of
their sizes (Fig. 5, bottom) in postshock flow as a function oftime are depicted in Figure 5, for example, for C grains with theinitial size of 0.01 m (dashed line), 0.1 m (dotted line), and1 m (thin solid line). The thick solid curve (Fig. 5, top) indicatesthe position of the shock front. As dust grains initially at rest in
the ISM intrude into the blast wave, they display different trajectories depending on their initial sizes, which clearly demonstrates
that the dust grains are segregated and subjected to different sputtering processes in the postshock flow; small grains with initial
radius of 0.01 m are quickly decelerated by drag force, trappedin gas near the shock front, and completely destroyed by thermal
sputtering. The 0.1 m–sized grains are gradually decelerated to
comove with gas in 7 ; 103 yr after entering the shock. Largesized grains with radius of 1 m do not undergo as much deceleration by drag force and continue to keep a high velocity relative togas. As a result, they are subjected to nonthermal sputtering, but
their sizes are reduced very little, partly because they stay in theNo. 1, 2006

445

DUST DESTRUCTION IN EARLY UNIVERSE

Fig. 5.—Time evolution of position (top) and size (bottom) of dust grains inpostshock flow for model C20 with nH;0 ¼ 1 cm3. The dust species consideredis a C grain with size of 0.01 m (dashed line), 0.1 m (dotted line), and 1 m(thin solid line). The thick solid curve in the top panel denotes the position of theshock front Rshock . [See the electronic edition of the Journal for a color versionof this figure.]

inner region of the postshock flow where the gas density is lower
than that near the shock front (see Fig. 3).
The modification of the size distribution of each dust species
due to destruction by sputtering is illustrated in Figure 6: Figure 6a for the initial size distribution before destruction and Figure 6b for the size distribution after destruction. Since the erosion
rate by sputtering does not strongly depend on grain size (seex 6.2), small grains are predominantly destroyed regardless ofgrain species; the number of small-sized grains such as Al2O3 isgreatly reduced. C, SiO2, and Fe grains whose initial size distributions are lognormal with a relatively large average size are
eroded but not completely destroyed, and the numbers of smaller
ones increase. Grains of size larger than 0.1 m are little affectedby erosion. Note that the size distribution summed up over all
grain species gets flatter for the smaller radius with time comparedwith the corresponding initial size distribution approximated bya power-law formula with index of 2.5, while that for the larger
radius remains almost unchanged.

Fig. 6.—Size distribution of each dust species for the unmixed grain modelfor (a) the initial size distribution before destruction and (b) the size distribution obtained by the calculation of dust destruction for model C20 with nH;0 ¼
1 cm3. [See the electronic edition of the Journal for a color version of this figure.]

than that of SiO2 grains because of the lowest erosion rate
among all grain species at T k 2 ; 106 K (see Fig. 2). In addition, FeS grains have a higher destruction efficiency (0.578)
than those of Mg2SiO4 (0.451) and MgO grains (0.505) despite
the larger initial average radius, which reflects not only the higher
TABLE 5The Efficiency of Dust Destruction for C20 Model with nH;0 ¼ 1
Unmixed Grain Model

Mixed Grain Model

7.1.3. Efficiency of Dust Destruction

Dust Speciesjaave; j( m)jaave; j( m)

The destruction efficiency j of each grain species is given in
Table 5 (see unmixed grain model) along with the initial averageradius. The efficiency of dust destruction is expected to be higher
for grain species with a smaller average size because the smaller
grains are predominantly destroyed by sputtering. In fact, theefficiency of destruction of Al2O3 grains with the smallest average radius is 0.667 and is the highest among the dust species
considered here. Si and Fe grains with initial average radii larger
than 0.1 m have destruction efficiencies of less than 0.2; 0.13
for Si grains and 0.15 for Fe grains. However, the destruction
efficiency of C grains with the smaller average radius is smaller

C....................................
Si ...................................
Fe...................................
FeS ................................
Al2O3 .............................
MgSiO3 .........................
Mg2SiO4 ........................
SiO2 ...............................
MgO ..............................
Fe3O4 .............................
Overall...........................

0.247
0.134
0.154
0.578
0.667
0.637
0.451
0.411
0.505
...
0.340

0.029
0.250
0.148
0.00880.00070.00450.00660.047
0.0058
...
...

...
...
...
...
0.794
0.631
0.586
0.399
...
0.741
0.497

...
...
...
...
0.0003
0.011
0.010
0.033
...
0.0022...

446

NOZAWA, KOZASA, & HABE

Vol. 648

of Mg2SiO4 grains in the mixed case is significantly larger thanthat (0.45) in the unmixed case. The reason is that Mg2SiO4 grains
formed in the unmixed ejecta have a power-law–like size distribution and include the large-sized grains of >0.1 m that are barely
destroyed in the shock; the average size is not always suitable forassessing the feasibility of dust destruction. Thus, we conclude
that the dust destruction efficiency is very sensitive to the initialsize distribution. The mass fraction of dust destroyed reaches up
to 34% for the unmixed grain model and 50% for the mixed grainmodel.
7.2. Timescale of Dust Destruction in the Early Universe

Fig. 7.—Same as Fig. 6, but for the mixed grain model. [See the electronic
edition of the Journal for a color version of this figure.]erosion rate at T k 107 K but also a mass distribution that is much
more weighted toward the smaller grains. Likewise, the destruction efficiency of MgSiO3 grains with average size comparable tothat of Mg2SiO4 and MgO is high (0.637) because of the lack oflarge grains. These facts indicate that the efficiency of dust destruction depends on the initial size distribution of dust grains, aswell as on the sputtering yield.
In order to clarify the effect of the initial size distribution
on the efficiency of dust destruction, we present the results ofdust destruction for the mixed grain model calculated with thesame SN model and value of nH;0 as those for the unmixed grainmodel. Figure 7 shows the size distribution of each grain species(Al2O3, MgSiO3, Mg2SiO4, SiO2, and Fe3O4) before destruction (Fig. 7a) and after destruction (Fig. 7b), and the destruction
efficiencies of these grains are tabulated in Table 5 (see mixed
grain model). As is also the case for the unmixed grain model,
the numbers of Al2O3 and Fe3O4 grains with the small average
size (less than a few tens of angstroms) are considerably reduced,
which results in high efficiencies of destruction (>0.7). For
MgSiO3, Mg2SiO4, and SiO2 grains with lognormal size distributions, the erosion of large grains leads to an increase in thenumber of smaller ones, and their destruction efficiencies span
the range of 0:4 P j P 0:64. Note that although the average size
is twice that in the unmixed case, the destruction efficiency (0.59)

In this subsection, we investigate the dependences of the efficiency of dust destruction j and the mass of the gas swept up by
the shock Mswept on E51 , Mpr , and nH;0 and derive an analytic
formula describing the timescale of dust destruction for the unmixed grain model. Figure 8 shows the destruction efficiency ofeach grain species versus SN explosion energy calculated for
nH;0 ¼ 1 cm3; Figure 8a is for Al2O3, FeS, Mg2SiO4, andFe grains, and Figure 8b is for MgSiO3, MgO, SiO2, C, andSi grains. In addition, the overall efficiency of dust destruction,
which is defined as the ratio of the total mass of dust destroyed
to the total mass of dust swept up by the shock, is plotted in
Figure 8a. The SN models used for the calculations are distinguished by open circles (CCSNe), open squares (HNe), and filled
triangles (PISNe).
The destruction efficiencies for each grain species are almost
the same for CCSNe and a PISN with an explosion energy of
E51 ¼ 1 irrespective of the progenitor mass Mpr , and they increase
with increasing E51 . Note that a high explosion energy with E51 
10 causes the temperature of the gas in the postshock flow to riseas high as 109 K, but this does not directly influence the efficiencyof dust destruction because the increase of gas temperature does
not always lead to the enhancement of the erosion rate by sputtering (see Fig. 2a). The reason for the increased efficiency withincreasing E51 is considered as follows. The high-velocity shock(initial shock velocity k104 km s1) generated by the energetic
SN explosion induces a high velocity of dust relative to gas. Then
dust grains are efficiently decelerated by drag force (see eq. [19]),
trapped in the high-density region near the shock front, and significantly eroded by thermal sputtering. Furthermore, since thehigh shock velocity takes a much longer time to drop down to
100 km s1 (see Fig. 9), dust grains are immersed in a hot plasma
for a long time, which causes even larger grains to be more eroded.
It should be pointed out here that the efficiencies of destruction for the models of P170 with E51 ¼ 20 and P200 with E51 ¼
28 are a little higher than that for H30 with E51 ¼ 30. The reason
is that the ejecta mass of PISNe is more than 6 times larger thanthat of HNe, and the resulting longer duration of the free expansion phase causes the longer truncation time than that for
HNe ( Fig. 9). On the other hand, the efficiency of destruction for
P150 with Mpr ¼ 150 M and E51 ¼ 1 is almost the same asthat for CCSNe with E51 ¼ 1 because the initial shock velocity(3000 km s1) much lower than that of CCSNe makes the truncation time comparable to that of CCSNe. Therefore, although
the destruction efficiency j is almost independent of Mpr as longas E51 ¼ 1, the j for PISNe higher than that for HNe at E51  10
reflects the difference in the explosion mechanism depending on
the progenitor mass.
To examine the dependence of j on E51 for each type of SN,
we calculate the coefficients a1; j and b1; j in equation (A1) givenin the Appendix for SNe II (C13, C20, C25, C30, H25, and H30)
and PISNe (P150, P170, and P200) separately, and the calculated overall efficiencies of dust destruction are indicated by dottedNo. 1, 2006

447

DUST DESTRUCTION IN EARLY UNIVERSE

Fig. 9.—Mass of gas swept up by shock, Mswept , and truncation time ttr vs.
explosion energy of SNe for nH;0 ¼ 1 cm3. The results for CCSNe, HNe, andPISNe are represented by open circles, open squares, and filled triangles, respectively. The linear lines are the power-law formula approximated by eq. (A2)for all SNe (solid lines), SNe II (dotted lines), and PISNe (dashed lines). [See theelectronic edition of the Journal for a color version of this figure.]

dust destruction versus nH;0 for C20 (crosses), H25 (open circles),and H30 models ( filled triangles) with E51 ¼ 1, 10, and 30,
respectively. The efficiency of dust destruction increases withincreasing nH;0 , since higher gas density causes more frequentcollisions between dust and gas to efficiently erode the surface ofdust grains by sputtering; for example, the mass fraction of dust
destroyed for the H30 model is 78% for nH;0 ¼ 10 cm3 but only
23% for nH;0 ¼ 0:1 cm3. In Figure 10b, we present the mass ofthe gas swept up by the shock Mswept as a function of nH;0 for theC20, H25, and H30 models. Note that Mswept decreases withincreasing nH;0 because the shock wave more quickly decelerates and travels only a small distance.
The approximation formulae presented in the Appendix beingcombined, the timescale of destruction for each dust species by
the interstellar shock in the early universe is presented by
Fig. 8.—Efficiency of destruction of each grain species in the unmixed grainmodel calculated for nH;0 ¼ 1 cm3 as a function of SN explosion energy for(a) Al2O3, FeS, Mg2SiO4, and Fe grains and (b) MgSiO3, MgO, SiO2, C, andSi grains. The SN models used for the calculations are distinguished by open
circles (CCSNe), open squares ( HNe), and filled triangles ( PISNe). In (a), theoverall efficiency of dust destruction is also plotted. The linear solid lines foreach grain species are the results calculated by the power-law formula given byeq. (A1) for all SNe. The dotted and dashed lines for the overall efficiencies of dust
destruction are the results of calculations for SNe II and for PISNe, respectively.
[See the electronic edition of the Journal for a color version of this figure.](SNe II) and dashed lines (PISNe) in Figure 8a. Although thedifference in the efficiency of dust destruction between SNe II andPISNe increases with increasing E51 , the deviations from the values calculated by a1; j and b1; j for all SN models are at most about
10% at E51 ¼ 30. Thus, we consider j to be almost independentof Mpr as long as E51 is the same. For reference, we tabulate thevalues of a1; j and b1; j for SNe II and PISNe in Table 6. The mass ofgas swept up by shock Mswept and the truncation time ttr calculated
for nH;0 ¼ 1 cm3 are presented in Figure 9 as a function of SN
explosion energy. As is also the case for j, Mswept and ttr are almost the same for E51 ¼ 1 regardless of Mpr , and they increase withincreasing E51 .
Next, we show the dependences of j and Mswept on the preshock gas density nH;0. Figure 10a plots the overall efficiency of

0:063
0:142 E51

1
SN;
j ¼ j (E51 ; nH;0 )

0:84144 E51n H; 0
MISM

M

SN

ð25Þ

for all SN models, where the dependences of j on E51 and nH;0
are given by equations (A1) and (A3), respectively. We derived
the timescale of dust destruction as a function of not only the
explosion energy of SNe but also the gas density in the ISM.
0:8,
The swept-up gas mass is proportional not to E51 but to E51
being different from the formula proposed by McKee (1989).
The timescale of dust destruction derived here could be applicable to investigating the time evolution of dust grains in theearly universe. In particular, the difference of the efficiency of
dust destruction for each grain species may have a great influence on the amount and size distribution of dust grains residing
in early interstellar space.
7.3. The Effects of the Cooling Processesof Gas on Dust Destruction
Figure 11 shows the cumulative energy lost by the atomic process (dashed line), the inverse Compton scattering (dotted line),
and the thermal emission from dust (thin solid line) calculated
for the SN model of C20 and for the ISM with parameters of448

NOZAWA, KOZASA, & HABE

Vol. 648

TABLE 6
Coefficients for j of the Approximation Formula (A1)
For All SNe

For SNe II

For PISNe

Dust Species

a1; j

b1; j

a1; j

b1; j

a1; j

b1; j

C.................................
Si ................................
Fe................................
FeS .............................
Al2O3 ..........................
MgSiO3 ......................
Mg2SiO4 .....................
SiO2 ............................
MgO ...........................
Overall........................

2.532E01
1.380E01
1.583E01
5.842E01
6.724E01
6.429E01
4.585E01
4.183E01
5.117E01
3.458E01

1.660E01
2.371E01
2.402E01
7.897E025.649E02
6.358E02
1.146E01
1.238E01
1.046E01
1.307E01

2.526E01
1.377E01
1.586E015.850E01
6.733E01
6.438E01
4.590E01
4.189E01
5.122E01
3.460E01

1.407E01
1.979E01
2.005E01
6.758E02
4.803E02
5.421E02
9.827E02
1.054E01
9.020E02
1.099E01

2.580E01
1.414E01
1.594E01
5.834E01
6.709E01
6.417E01
4.595E01
4.190E01
5.124E01
3.476E01

1.820E01
2.631E01
2.707E01
8.857E02
6.393E02
7.165E021.275E01
1.383E011.160E01
1.463E01nH;0 ¼ 1 and Z ¼ 104 Z at a redshift of z ¼ 20. The total energy lost by these cooling processes is drawn by the thick solid
curve. The inverse Compton cooling is comparable to that by theatomic process in the early phase of the SNR (t < 104 yr). As thegas temperature decreases, the atomic process by H and He line
cooling becomes dominant. Compared with the above two cooling processes, the thermal emission from dust is extremely low
and contributes only less than 0.1% to the total energy loss. The
transition of nonradiative shock to radiative shock occurs when
the cumulative energy loss reaches 0.01% of the explosion energy. Thus, only the H and He line cooling processes affect thetime evolution of the nonradiative shock and the dust destruction
efficiencies in the ISM with metallicity less than Z ¼ 104 Z
corresponding to the dust-to-gas mass ratio of 4:5 ; 107 in thispaper.
However, the effects of cooling of gas by dust are expected tocome to be important for the ISM with higher dust-to-gas mass
ratio. The result of the calculation for dust destruction in the ISM
with Z ¼ 104 Z and the dust-to-gas mass ratio 4:5 ; 103
shows that the truncation time and the overall efficiency of dust
destruction decrease by 4.4% and 2.6%, respectively, comparedwith the results for a dust-to-gas mass ratio of 4:5 ; 107 ,Fig. 10.—(a) Overall efficiency of dust destruction vs. nH;0 for models C20(crosses), H25 (open circles), and H30 ( filled triangles). The solid curves arecalculated by the approximation formula of eq. (A3). (b) Mass of gas swept up
by shock for models C20 (crosses), H25 (open circles), and H30 ( filled triangles)
as a function of nH;0 . The linear lines are the power-law approximation formula of
0:063. [See the electronic edition of the Journal
Mswept / ngH;0 with g ¼ 0:142E51
for a color version of this figure.]

Fig. 11.—Time evolution of cumulative energy lost by the atomic process(dashed line), inverse Compton cooling (dotted line), and thermal emission fromall dust grains (thin solid line) in units of E51 for the SN model C20, nH;0 ¼ 1, and
Z ¼ 104 Z at a redshift of z ¼ 20. The total energy loss is indicated by the thicksolid curve. [See the electronic edition of the Journal for a color version of thisfigure.]No. 1, 2006

DUST DESTRUCTION IN EARLY UNIVERSE

although the gas temperature in postshock gas decreases by
about 20%. Even if the dust-to-gas mass ratio is raised by up to
4 orders of magnitude, the cooling of gas by dust in a nonradiative shock has significant effects on neither the efficiency of
dust destruction nor the evolution of the nonradiative shock.
Furthermore, the inverse Compton cooling influences the evolution of nonradiative shocks in the early universe. Because thecooling rate by the inverse Compton scattering is proportional
to ne T (1 þ z)4 and that by the atomic process is proportional tone nH  n2H , the contribution of the inverse Compton cooling is enhanced in postshock flow with a lower gas density and/or higher
gas temperature at higher redshift. At redshift z ¼ 40, for the SN
model C20 with preshock gas density of nH;0 ¼ 1 cm3, the truncation time is only 2% shorter than that calculated at z ¼ 20. Thus,in this case, the inverse Compton cooling seems not to affectthe destruction efficiency of dust grains. However, more systematic studies covering wide ranges of SN explosion energy and gas
density in the ISM are necessary to reveal the effect of the inverse
Compton cooling.
Note that the cooling function of gas for Z ¼ 103 Z is almost
the same as that for the zero metal case (Sutherland & Dopita
1993). Thus, the parameters of the ISM considered in this paper(Z ¼ 104 Z and z ¼ 20) are well fitted for investigating dust
destruction in the early universe, and the timescale of dust destruction derived in x 7.2 is applicable over the wide ranges ofmetallicity of gas in the ISM (Z P 103 Z) and redshift (z  40).
8. SUMMARY
We investigate the destruction of dust grains in the interstellar
shocks driven by SNe as the second step to revealing the evolution of dust in the early universe, based on the dust models obtained by Nozawa et al. (2003). We focus on dust destruction innonradiative shocks where dust grains are predominantly destroyed by nonthermal and thermal sputtering because of hightemperature and high velocity of the gas. The sputtering yields
for the combinations of dust and ion species of interest to us areevaluated by applying the universal relation for sputtering yieldwith a slight modification and by determining the value of K byfitting to the available experimental data and /or the results of the
EDDY simulations. The modified universal relation derived here
can present better fits to the available experimental data than theprevious relation.
In the calculations of dust destruction, we solve the time evolution of gas temperature and density in spherically symmetricshocks, adopting the hydrodynamic models by Umeda & Nomoto(2002) as the initial conditions for interstellar shocks and including the cooling of gas by the atomic process, inverse Compton
scattering, and the thermal emission of dust. The erosion of dustby thermal and nonthermal sputtering caused by motion of dust
relative to gas is carefully treated by taking into account the size
distribution of each dust species. The results of the calculations are
summarized as follows.
1. Because the sputtering predominantly destroys the small
grains, the number of small-sized grains such as Al2O3 andMgSiO3 is greatly reduced. The erosion of C, SiO2, and Fe grains
whose size distributions are lognormal-like, with the average size
larger than 0.01 m, increases the numbers of smaller grains. The
size distribution summed up over all grain species becomes flatterfor small radius compared with the initial size distribution, whilethat for radius larger than 0.2 m remains almost unchanged.
2. The efficiency of dust destruction is higher for the grains
with a small average size such as Al2O3 and MgSiO3 and witha power-law–like size distribution such as FeS, Mg2SiO4, and449

MgO. On the other hand, Si, Fe, C, and SiO2 grains, which havea lognormal-like size distribution and a relatively large averageradius, have a lower efficiency of destruction. A detailed analysisof the behavior of the dust destruction efficiency for each grain
species indicates that not only sputtering yields but also the initial size distribution plays a crucial role in the efficiency of dust
destruction by sputtering.
3. The efficiency of destruction j for each dust species increases with increasing explosion energy E51 and/or increasingpreshock gas density nH;0 but is almost independent of the progenitor mass Mpr of the SN as long as E51 remains the same. The
destruction efficiency j as a function of E51 is reproduced wellby a power-law formula given by equation (A1). The dependenceof j on nH;0 is expressed well by the quadratic equation (A3) interms of log nH;0 .
4. As is also the case for j, the mass of gas swept up by theshock wave, Mswept , is an increasing function of E51 and is approximated by a power-law formula given by equation (A2). However, Mswept decreases with increasing gas density in the ISM, andits dependence is also reproduced by a power-law formula whose
index is given by equation (A4) as a function of E51 . Finally, by
combining these results, we present an analytic formula for thetimescale of destruction for each grain species in the early universe as a function of E51 and nH;0 . The derived timescale of dust
destruction can be used to investigate the time evolution of theamount of dust grains in the early universe.
5. In the early universe, only the H and He line cooling processes affect the time evolution of nonradiative shocks and theefficiency of dust destruction. The thermal emission from dust
grains is not very important even for the ISM with a dust-to-gas
mass ratio of 103 as long as the metallicity in the ISM is Z P
103 Z. In addition, inverse Compton cooling is expected not tohave great effects on the evolution of shocks at high redshifts.
Thus, the timescale of dust destruction derived in this paper isapplicable for ranges of metallicity of gas in the ISM of Z P
103 Z and redshift of z  40.
In this study, we have focused on dust destruction in nonradiative shocks where dust grains in the ISM are considered to be
predominantly destroyed. We should mention here that a part ofthe dust grains formed in the ejecta of SNe are expected to be
destroyed by reverse shocks penetrating the ejecta. Because sputtering is the most dominant destruction process of dust grains inreverse shocks as well, the simulation code of dust destruction
constructed in this paper can be applied to explore the destructionof dust grains by reverse shocks. This work is now in progress.
Furthermore, in radiative shocks with shock velocities below
100 km s1, shattering and partial vaporization by grain-grain
collisions become important in the presence of a magnetic field.
A magnetic field frozen into the gas causes charged grains togyrate and accelerate behind the shock front and results in theenhancement of not only the sputtering rate but also the graingrain collision frequency. The shattering by grain-grain collision mainly leads to the redistribution of grain sizes ( Tielens
et al. 1994; Jones et al. 1996). The sophisticated study of thissubject is left for future work.

The authors are grateful to the anonymous referee for criticalcomments, which have improved the manuscript. The authors
thank H. Umeda and K. Nomoto for making the ejecta model ofPopulation III supernovae available. This work has been supported in part by a Grant-in-Aid for Scientific Research from theJapan Society for the Promotion of Sciences (16340051).

450

NOZAWA, KOZASA, & HABE

Vol. 648

TABLE 7
Coefficients for Mswept of the Approximation Formula (A2)
Coefficient

For all SNe

For SNe II

For PISNe

a2 ............................................
b2 ............................................

4.144E+03
7.967E01

4.118E+03
6.627E01

4.481E+03
8.864E01

APPENDIX
THE APPROXIMATION FORMULAE FOR j AND Mswept
In this appendix, we present the approximate formulae for the efficiency of dust destruction j and the mass of the gas swept up by
the shock, Mswept , as a function of the SN explosion energy E51 and the gas density nH;0 in the ISM. The energy dependence of thedestruction efficiency for each grain species can be reproduced well by a power-law formula given bybj ¼ a1; j E511; j :ðA1ÞThe solid lines in Figure 8 indicate the least-squares fits by this power-law formula for all SNe whose numerical coefficients a1; j andb1; j are tabulated in Table 6. The dependence of Mswept on E51 is well fitted by the power-law formulab2
Mswept ¼ a2 E51

ðA2Þusing the coefficients a2 and b2 given in Table 7 and is depicted in Figure 9 for all SNe (solid line), SNe II (dotted line), and PISNe(dashed line). The destruction efficiency j for each grain species j as a function of nH;0 , as well as the overall destruction efficiency, isexpressed by2log j ¼ cj log nH;0 þ dj log nH;0 þ ej ;ðA3Þand the coefficients cj , dj , and ej derived by the least-squares fits are given in Table 8 for models C20, H25, and H30, respectively,
along with the coefficients for overall destruction efficiency. The solid curves in Figure 10a depict the overall destruction efficiencies
calculated by equation (A3) for models C20, H25, and H30. This formula reproduces the destruction efficiencies calculated by thesimulations with an accuracy of within less than 10% for nH;0  10 cm3. The calculated Mswept can be reproduced by the power-law
formula Mswept / ngH;0 , where the index g is weakly dependent on E51 and is approximated by
0:063
g ¼ 0:142E51:

ðA4Þ

This formula produces good agreement with the results of the simulations, as represented by the solid lines in Figure 10b.

TABLE 8
Coefficients for the Approximation (A3)
C20

H25

H30

Dust Species

cj

dj

ej

cj

dj

ej

cj

dj

ej

C.................
Si ................
Fe................
FeS .............
Al2O3 ..........
MgSiO3 ......
Mg2SiO4 .....
SiO2 ............
MgO ...........
Overall........

0.062
0.044
0.057
0.088
0.031
0.0600.087
0.086
0.079
0.067

0.41
0.54
0.55
0.23
0.14
0.18
0.31
0.33
0.29
0.33

0.61
0.87
0.81
0.23
0.17
0.19
0.35
0.39
0.30
0.47

0.07
0.06
0.078
0.077
0.0300.057
0.092
0.093
0.086
0.067

0.36
0.50
0.49
0.18
0.12
0.14
0.25
0.27
0.23
0.29

0.46
0.67
0.61
0.16
0.13
0.13
0.24
0.27
0.20
0.35

0.08
0.072
0.096
0.073
0.032
0.055
0.097
0.097
0.089
0.072

0.33
0.47
0.46
0.16
0.11
0.12
0.22
0.24
0.20
0.26

0.39
0.58
0.51
0.13
0.11
0.11
0.19
0.22
0.16
0.30

No. 1, 2006

DUST DESTRUCTION IN EARLY UNIVERSE

451

REFERENCES
Andersen, H. H., & Bay, H. L. 1981, in Sputtering by Particle Bombardment I,
Mair, G., Müllar, E., Hillebrandt, W., & Arnold, C. N. 1988, A&A, 199, 114
ed. R. Behrisch ( Berlin: Springer), 145
Matsunami, N., Yamamura, Y., Itakawa, W., Itoh, N., Kazumata, Y., Myagawa, S.,
Bach, H. 1970, Nucl. Instrum. Methods, 84, 4Morita, K., & Shimizu, R. 1980, Radiat. Effects Lett., 57, 15
Bach, H., Kitzmann, I., & Schroder, H. 1974, Radiat. Effects, 21, 31
May, P. W., Pineau des Forets, G., Flower, D. R., Field, D., Allan, N. L., &Baines, M. J., Williams, I. P., & Asebiomo, A. S. 1965, MNRAS, 130, 63
Purton, J. A. 2000, MNRAS, 318, 809
Bertoldi, F., Carilli, C. L., Cox, P., Fan, X., Strauss, M. A., Beelen, A., Omont, A.,
McKee, C. F. 1989, in Interstellar Dust, ed. L. J. Allamandora & A. G. G. M.
& Zylka, R. 2003, A&A, 406, L55
Tielens ( Dordrechent: Kluwer), 431
Bianchi, S., & Ferrara, A. 2005, MNRAS, 358, 379
McKee, C. F., Hollenbach, D. J., Seab, C. G., & Tielens, A. G. G. M. 1987,
Blank, P., & Wittmaack, K. 1979, J. Appl. Phys., 50, 1519
ApJ, 318, 674
Bohdansky, J. 1984, Nucl. Instrum. Methods Phys. Res. B, 2, 587
Nenadovic, T., Perraillon, B., Bogdanov, Z., Djordjevic, Z., & Milic, M. 1990,
Bohdansky, J., Bay, H. L., & Ottenberger, W. 1978, J. Nucl. Math., 76–77, 163
Nucl. Instrum. Methods Phys. Res. B, 48, 538
Bohdansky, J., Bay, H. L., & Roth, J. 1977, In Proc. 7th Int. Vacuum Cong.
Nozawa, T., Kozasa, T., Umeda, H., Maeda, K., & Nomoto, K. 2003, ApJ, 598,
and 3rd Int. Conf. on Solid Surfaces, ed. R. Dobrozemsky et al. ( Vienna:
785
R. Dobrozemsky et al.), 1509
Ohya, K., & Kawata, J. 1997, Japanese J. Appl. Phys., 36, L298
Bohdansky, J., Roth, J., & Bay, H. L. 1980, J. Appl. Phys., 51, 2861
Omukai, K., Tsuribe, T., Schneider, R., & Ferrara, A. 2005, ApJ, 626, 627
Bohdansky, J., Roth, J., & Sinha, M. K. 1976, in Proc. 9th Symp. on Fusion
Priddey, R. S., Isaak, K. G., McMahon, R. G., Robson, E. I., & Pearson, C. P.
Technology ( London: Pergamon), 541
2003, MNRAS, 344, L74
Bromm, V., & Larson, R. B. 2004, ARA&A, 42, 79
Robson, I., Priddey, R. S., Isaak, K. G., & McMahon, R. G. 2004, MNRAS,
Cazaux, S., & Spaans, M. 2004, ApJ, 611, 40
351, L29
Chase, M. W., Davies, C. A., Downey, J. R., Frurip, D. J., McDonald, R. A., &Rosenberg, D., & Wehner, G. K. 1962, J. Appl. Phys., 33, 1842
Syverud, A. N. 1985, J. Phys. Chem. Ref. Data, 14, Suppl. 1
Roth, J. 1983, in Sputtering by Particle Bombardment II, ed. R. BehrischCoburn, J. W., Winters, H. F., & Chuang, T. J. 1977, J. Appl. Phys., 48, 3532
( Berlin: Springer), 91
Draine, B. T., & Salpeter, E. E. 1979, ApJ, 231, 77
Roth, J., Bohdansky, J., & Ottenberger, W. 1979, Rep. IPP 9/26 (Garching:Dwek, E. 1987, ApJ, 322, 812
MPI )———. 1998, ApJ, 501, 643
Schneider, R., Ferrara, A., Natarajan, P., & Omukai, K. 2002, ApJ, 571, 30
Dwek, E., & Arendt, R. G. 1992, ARA&A, 30, 11
Schneider, R., Ferrara, A., & Salvaterra, R. 2004, MNRAS, 351, 1379
Dwek, E., Foster, S. M., & Vancura, O. 1996, ApJ, 457, 244
Schneider, R., Ferrara, A., Salvaterra, R., Omukai, K., & Bromm, V. 2003,
Edwin, R. P. 1973, J. Phys. D, 6, 833
Nature, 422, 869
Eernisse, E. P. 1971, J. Appl. Phys., 42, 480
Sigmund, P. 1969, Phys. Rev., 184, 383
Field, D., May, P. W., Pineau des Forêts, G., & Flower, D. R. 1997, MNRAS,———. 1981, in Sputtering by Particle Bombardment I, ed. R. Behrisch
285, 839( Berlin; Springer), 9
Gnedin, N. Y., Ferrara, A., & Zweibel, E. G. 2000, ApJ, 539, 505
Sommerfeldt, H., Mashkova, E. S., & Molchanov, V. A. 1972, Phys. Lett. A,
Hauser, M. G., & Dwek, E. 2001, ARA&A, 39, 249
38, 237
Hechtl, E., Bohdansky, J., & Roth, J. 1981, J. Nucl. Math., 103–104, 333
Southern, A. L., Willis, W. R., & Robinson, M. T. 1963, J. Appl. Phys., 34, 153
Heger, A., & Woosley, S. E. 2002, ApJ, 567, 532
Sutherland, R. S., & Dopita, M. A. 1993, ApJS, 88, 253
Hirashita, H., & Ferrara, A. 2002, MNRAS, 337, 921
Tabata, T., Ito, R., & Okabe, S. 1972, Nucl. Instrum. Methods, 103, 85
Hirashita, H., Nozawa, T., Kozasa, T., Ishii, T. T., & Takeuchi, T. T. 2005,
Takeuchi, T. T., Hirashita, H., Ishii, T. T., Hunt, L. K., & Ferrara, A. 2003,
MNRAS, 357, 1077
MNRAS, 343, 839
Ikeuchi, S., & Ostriker, J. P. 1986, ApJ, 301, 522
Takeuchi, T. T., Ishii, T. T., Nozawa, T., Kozasa, T., & Hirashita, H. 2005,
Iskef, H., Cunningham, J. W., & Watt, D. E. 1983, Phys. Med. Biol., 28, 535
MNRAS, 362, 592
Jones, A. P., Tielens, A. G. G. M., & Hollenbach, D. J. 1996, ApJ, 469, 740
Tan, J. C., & McKee, C. F. 2004, ApJ, 603, 383
Jones, A. P., Tielens, A. G. G. M., Hollenbach, D. J., & McKee, C. F. 1994,
Tielens, A. G. G. M., McKee, C. F., Seab, C. G., & Hollenbach, D. J. 1994,
ApJ, 433, 797
ApJ, 431, 321
Kitayama, T., Yoshida, N., Susa, H., & Umemura, M. 2004, ApJ, 613, 631
Todini, P., & Ferrara, A. 2001, MNRAS, 325, 726
Kozasa, T., Hasegawa, H., & Nomoto, K. 1989, ApJ, 344, 325
Umeda, H., & Nomoto, K. 2002, ApJ, 565, 385
Laegreid, N., & Wehner, G. K. 1961, J. Appl. Phys., 32, 365
van Albada, G. D., van Leer, B., & Roberts, W. W., Jr. 1982, A&A, 108, 76
Loeb, A., & Haiman, Z. 1997, ApJ, 490, 571
Vanzi, L., Hunt, L. K., Thuan, T. X., & Izotov, Y. I. 2000, A&A, 363, 493
Machida, M. N., Tomisaka, K., Nakamura, F., & Fujimoto, M. Y. 2005, ApJ,
von Seefeld, H., Schmidl, H., Behrisch, R., & Scherzer, B. M. U. 1976, J. Nucl.
622, 39
Math., 63, 215
Maiolino, R., Oliva, E., Ghinassi, F., Pedani, M., Mannucci, F., Mujica, R., &Wolfire, M. G., Hollenbach, D., McKee, C. F., Tielens, A. G. G. M., & Bakes,
Juarez, Y. 2004a, A&A, 420, 889
E. L. O. 1995, ApJ, 443, 152
Maiolino, R., Schneider, R., Oliva, E., Bianchi, S., Ferrara, A., Mannucci, F.,
Ziegler, J. F., Biersack, J. P., & Littmark, U. 1985, The Stopping and Range of
Pedani, M., & Roca Sogorb, M. 2004b, Nature, 431, 533
Ions in Solids ( New York: Pergamon)