 

 

LIBRARY
Michigan State
University

 

 

 

This is to certify that the
thesis entitled

STABILITY OF CARBON FUSION
ON ACCRETING NEUTRON STARS

presented by

PHILIPP GIRICHIDIS

has been accepted towards fulfillment
of the requirements for the

M. S. degree in Physics and Astronomy

 

 

WA.—

Major Professor’s Signature
/ Z , l l . 0 8/

Date

 

MSU is an Afﬁrmative Action/Equal Opportunity Employer

I PLACE IN RETURN Box to remove this checkout from your record.
TO AVOID FINES return on or before date due.

MAY BE RECALLED with earlier due date if requested.

 

DATE DUE DATE DUE DATE DUE

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

5!08 K;lProj/Acc&Pres/ClRC/DateDue indd

STABILITY OF CARBON FUSION
ON ACCRETING NEUTRON STARS

By

Philipp Girichidis

A THESIS

Submitted to
Michigan State University
in partial fulﬁllment of the requirements
for the degree of

MASTER OF SCIENCE
Physics and Astronomy

2008

ABSTRACT

STABILITY OF CARBON FUSION
ON ACCRETING NEUTRON STARS

By

Philipp Girichidis

With observing missions like the Rossi X—my Timing Explorer ( RX TE ), BeppoS'AX,
XMM-Newton and Chandra many thermonuclear activities on neutron stars have
been observed, especially thermonuclear X—ray bursts on accreting neutron stars.
Aside from frequent short type I X-ray bursts, rare and very long enduring high
energetic bursts, the so-called superbursts, have been found. The large total released
energy during a superburst indicates a larger ignition depth and higher ignition
temperatures than it is the case for type I bursts. These ignition conditions lead
to the conclusion, that unstable carbon burning triggers the thermonuclear runaway
for the superburst.

This work focusses on the carbon plasma layer and its nuclear fusion stability.
With numerical simulations a stability analysis of the layer has been performed, in
order to ﬁnd precise conditions for unstable ignitions. The numerical model used
in this thesis combines a full reaction network with a complex number perturbation
stability analysis, in which effects of temperature, energy ﬂux, composition and ac-
cretion rate on the stability were examined. Furthermore, different burning regimes
in the carbon burning process have been investigated in order to determine the na-
ture of the explosion as well as the exact ignition depth. For a few sets of parameters
burning oscillations were investigated. For the neutron star KS 1731-260 the stability

analysis was used to determine the chemical composition of the carbon burning layer.

ACKNOWLEDGMENTS

Many Thanks to Dr. Edward Brown and Dr. Ralf Klessen for supporting my
unconventional combination of an exchange year, a German diploma thesis and the
MSU masters program. Dr. Brown’s patience and constant availablility made it
easy for me to ﬁnd my way through the neutron star theory and the computational
methods. I thank Dr. Klessen for his help with all the fomality concerning my stay
abroad and my external diploma thesis.

Thanks to my committee members Dr. Edward Brown, Dr. Hendrik Schatz and
Dr. Brian O’Shea.

Thanks to all my ofﬁce mates for the very good company and for helping me with
all my language questions and uncertanties.

Thanks to Debbie Barratt for helping me with all the paperwork before and during
my time at the Michigan State University.

Thanks to Dr. Wolfgang Bauer for encouraging me to come to the United States
and the Michigan State University.

Thanks to the Studienstiftung des deutschen Volkes for supporting me throughout

my studies in Heidelberg and especially during my time at MSU.

iii

Contents

1 Introduction 1
1.1 Neutron Stars in Binary Systems .................... 1
1.2 X-ray Bursts ................................ 2
1.3 Superbursts ................................ 3

1.3.1 Observations ........................... 3
1.3.2 Theoretical Conclusions ..................... 4
1.4 Consequences of an Unstable Ignition .................. 5
1.5 Previous Studies on Superbursts ..................... 5
1.6 Goal of this Study ............................ 6

2 Introduction to Neutron Stars 7
2.1 General Properties ............................ 7
2.2 Structure ................................. 8

2.2.1 Atmosphere ............................ 8
2.2.2 Ocean ............................... 9
2.2.3 Crust ................................ 9
2.2.4 Core ................................ 10
2.3 Schematic Thermal Proﬁle ........................ 11

3 Structure of the Outer Layers 13
3.1 Accretion ................................. 13
3.2 Equations of Hydrostatic Equilibrium .................. 14
3.3 Equation of State (EOS) ......................... 17

3.3.1 A Polytropic Gas ......................... 17
3.3.2 Numerical Fit of the EOS .................... 18
3.3.3 Coupling Parameter ....................... 21
3.4 Thermodynamics ............................. 21
3.4.1 Heating .............................. 21
3.4.2 Cooling .............................. 22
3.5 Thermal Structure Equations ...................... 25

4 Nuclear Processes 27
4.1 Triple-a Process .............................. 27
4.2 Rapid-Proton Process (rp—process) .................... 28
4.3 Carbon Burning and Advanced Burning Stages ............ 29
4.4 Released Energy .............................. 32

iv

4.5 Electron Capture .............................
4.6 Reaction Network .............................
4.7 Electron Screening ............................

5 Numerical Methods
5.1 Basic Numerical Concepts ........................

5.2 General Properties of the Code .....................

6 Steady State Solution
6.1 Helium Boundary Conditions ......................
6.1.1 Temperature ............................
6.1.2 Flux ................................
6.1.3 Density ..............................
6.1.4 Sensitivity of the Boundary Conditions .............
6.2 Carbon Boundary Conditions ......................

6.3 Steady State Results ...........................
6.3.1 Composition ............................
6.3.2 Properties of the Material ....................

6.4 Electron Capture .............................

6.5 Validating the Calculations ........................

6.6 Inﬂuence of Heavy Elements .......................

6.7 Inﬂuence of Temperature .........................

7 Linear Stability Analysis
7.1 Linear Perturbation Equations ......................
7.2 Stability and Overstability ........................
7.3 Application to the Stellar Model .....................
7.4 Perturbation Function ..........................
7.5 Calculation Procedure ..........................
7.6 Discretization and Computing Time ...................
7.7 Results of the Linear Stability Analysis .................
7.7.1 Stability Flux ...........................
7.7.2 Inﬂuence of heavy elements ...................
7.7.3 Inﬂuence of the Temperature ..................
7.7.4 Inﬂuence of the Accretion Rate .................
7.8 Burning Regimes .............................
7.8.1 Direct and Delayed Bursts ....................
7.8.2 Overstable Modes .........................

8 Application
8.1 Neutron stars in quiescence .......................
8.2 The neutron star in KS 1731-260 ....................

9 Conclusion and Outlook

A Thermodynamic Equations

33
34
37

41
41
42

45
45
45
46
47
47
48
49
49
50
51
56
57
63

66
66
67
68
69
69
73
75
75
76
77
77
79
79
79

85
85
86

88

90

B Column Depth

C Elements in the Reaction Network

vi

92

94

List of Figures

2.1

2.2
3.1

3.2

3.3
4.1

4.2

4.3
6.1

6.2
6.3
6.4
6.5
6.6

6.7

6.8
6.9
6.10
6.11
6.12
6.13

Schematic structure of a neutron star .................. 11
Schematic proﬁle of an accreting neutron star ............. 12
Fit functions to the EOS ......................... 20
Gravitational ﬂux ............................. 22
Error in equation 3.39 due to an idealized speciﬁc heat ........ 26
Sample rates for some important reactions ............... 31
Rate Ratio of producing 19 + 23Na and 4He + 20Ne ........... 31
Element chart of included isotopes ................... 37
Abundances with pure accreted helium ................. 49
Coupling parameter ............................ 50
Ignition depth for different accretion rates ............... 51
Ignition density for different accretion rates .............. 52
Released energy for different accretion rates .............. 52

Thermal proﬁle for FoLnnéL = 0.076 MeV/nuc, T8 = 4.0, Th 2 0.2 ThEdd
and an initial composition of 50% carbon ................ 53

Thermal proﬁle for F033;? = 1.0 MeV/nuc, T8 = 4.0, m = 0.2 mEdd

and an initial composition of 50% carbon ................ 54
Energy level for electron capture on 23Na ................ 56
Thermal conductivity for different compositions ............ 58
Thermal proﬁle for 12C and 160 ..................... 59
Thermal proﬁle for 12C and 56Fe .................... 60
Thermal proﬁle for 12C and 10’J‘Ru ................... 61
Ignition density for different heavy elements .............. 62

vii

6.14 Released energy for different heavy elements .............. 62

6.15 Ignition depth for different temperatures ................ 64
6.16 Ignition density for different temperatures ............... 64
6.17 Released energy for different temperatures ............... 65
6.18 Released energy for different temperatures ............... 65
7.1 Perturbation function ........................... 70
7.2 Example for a stable and an unstable solution for w .......... 72
7.3 Growth time of the perturbation ..................... 74
7.4 Inﬂuence of heavy elements on the stability ............... 76
7.5 Inﬂuence of the temperature on the stability .............. 77
7.6 Inﬂuence of the accretion rate on the stability ............. 78

7.7 Direct and delayed bursts for an accretion rate of in = 0.2 mEdd and
a composition of 12C : 56Fe = 1 : 1 .................... 80
7.8 Direct and delayed bursts for an a temperature of T8 = 4.0 and a .
composition of 12C : 56Fe = 1 : 1 ..................... 80

7.9 Direct and delayed bursts for an a temperature of T8 = 4.0 and an

accretion rate of m = 0.2 mEdd ..................... 81
7.10 Appearance pattern of the oscillation .................. 82
7.11 Disappearance pattern of the oscillation ................ 83
7.12 Appearance and disappearance of the oscillation ............ 84
8.1 Thermal proﬁle for KS 1731-260 ..................... 86
8.2 Stability ﬂux as a function of composition for KS 1731-260 ...... 87

B.1 Relation between column depth and density .............. 93

viii

List of Tables

1.1 Superburst sources and properties .................... 3
3.1 Parameters for the EOS ﬁt function ................... 19
4.1 Reaction types used in the nuclear network ............... 35
CI List of elements .............................. 95

Chapter 1

Introduction

1.1 Neutron Stars in Binary Systems

Neutron stars have been observed in binary systems with other neutron stars, white
dwarfs or non-degenerate stars. These systems can be divided into two different types:
systems without mass ﬂow and systems with mass ﬂow from the companion star
onto the neutron star. Many of the accreting neutron stars erupt in high luminous
thermonuclear bursts every few hours to days. These events are known as type I
X-ray bursts or simply as X-ray bursts. These ﬂares, which have been discovered
independently by Grindlay et al. [1] and Belian et al. [2], have been observed by
the Rossi X-ray Timing Explorer (RX TE), the BeppoSAX, the XMM-Newton and
the Chandra mission. Since the launch of BeppoSAX and RXTE in 1996, the sky
has been monitored with high sensitivity and frequency in the X-ray band. This
capability has opend the discovery space for bursts with longer recurrence time in
the order of years, which were missed by previous missions. Cornelisse et al. [9]
observed the ﬁrst so—called superburst, a rare long lasting burst event with a total
energy some orders of magnitude higher than in a type I X-ray burst, in the familiar

type I burster 4U 1735-44.

On a neutron star with a typical mass of 1.4M@ and a radius of 10km, where
the gravitaional energy per accreted baryon is ~ 200 MeV, it perhaps does not seem
obvious that nuclear burning with a total energy release of a few MeV per nucleon
plays an important role in the burst evolution and energy release. However, large
amounts of unburned accumulated fuel that undergo explosive burning exceed the
gravitational energy release by a lot and determine the luminosity proﬁle of the
star. What is to be determined is which nuclear reactions under what circumstances
come into consideration as trigger mechanisms for these unstable ignitions. The total
released energy during a burst, the recurrence time and the time evolution of these
explosive events give information about the structure in the neutron stars layers

below the atmoshere.

1.2 X-ray Bursts

The luminosities during a type I X-ray burst are many times larger than the per-

sistent value of 1036 — 1038 ergss_1.

Typical type I bursts exhibit rise times of
several seconds, last from tens to hundreds of seconds and release a total energy of
1039 — 1040 ergs (e.g. Strohmayer & Bildsten [11]). These events are caused by un-
stable ignition of accreted hydrogen and helium on the surface of the neutron star in
constrast to type II bursts that are ignited by sudden accretion events (Lewin et al.
[14]) The accumulated H/ He mixture is compressed and heated in the atmosphere
of the neutron star, forming a layer of layer of several meters thickness. After an
hour of compression the temperture and density are high enough for the H /He to
ignite and burn via the hot carbon-nitrogen-oxygen (CNO) cycle. If the temperature
and density are high enough to trigger unstable H / He burning all the available fuel

is converted to heavier elements within a few seconds which is known as the type I

X-ray burst.

 

 

Source 4U 1820-30 4U1735-44 KS 1731-260

 

 

Duration [hr] 3 7 12
Lpers [LEdd] z 0.1 as 0.25 z 0.1
kTmax [keV] as 3 z 2.6 m 2.4
Lpeak [1038 ergs/s] 3.4 1.5 1.4
Ebumt [1038 ergs] > 1.4 > 0.5 rs 1
References S00, SB02 C00 K02
Source 4U 1636-53 Ser X-l GX 3+1
Duration [hr] > 2 — 3 "~V 4 > 3.3
Lpers [LEddl % 0.1 8 0.2 $3 0.2
kTmax [keV] ? z 2.6 ~ 2
Lpeak [1038 ergs/s] 1.2 1.6 0.8
Eb“rst [1038 ergs] 0.5 -— 1 a 0.8 > 0.6
References W01, SM02 C02 K02

 

Table 1.1: Superburst sources and properties (after Kuulkers et al. [6]). Reference
abbreviations: 800 (Strohmayer [4]), SB02 (Strohmayer & Brown [5]), C00 (Cor-
nelisse et al. [3]), K02 (Kuulkers et al. [6]), W01 (Wijnands [7]), SM02 (Strohmayer
& Markwardt [8]), C02 (Cornelisse et al. [9]), K02b (Kuulkers [10])

1.3 Superbursts

1.3. 1 Observations

Shortly after the discovery of superbursts by Cornelisse et al. [3] six more super-
bursts, lasting 3 - 5 hours from previously known X-ray bursters have been reported:
4U 1820-30 (Strohmayer [4], Strohmayer 81. Brown [5]), KS 1731-260 (Kuulkers et al.
[6]), 4U 1636-53 (Wijnands [7], Strohmayer & Markwardt [8]), Ser X-l (Cornelisse
et al. [9]) and GX 3+1 (Kuulkers [10]). From the neutron star 4U 1636-53 two
bursts were observed, separated by a recurrence time of 4.7 years (Wijnands [7]). An
overview of the superburst properties can be found in table 1.1. These ﬂares show
all hallmarks of thermonuclear explosions; they have thermal spectra which soften
with time. Superbursts show a similar shape in rise and decline of the luminosity,
but are roughly 1000 time more energetic and last much longer, typically several

hours. The recurrence times are much longer, in the order of years, which makes

them rare events. The amount of released energies during a superburst is in the
order Of 1042 ergs s‘l.

As neutron stars are dense objects the relatively large gravitational redshift z in
comparison to other stars inﬂuences the Observational measurements of temperature,
luminosity and radius and can therefore not be neglected. Deﬁning the gravitational

radius as rg = 2G’M/c2 the measured quantities at inﬁnity, indicated by the oo

subscript, are related to the actual values on the neutron stars surface as

Loo = L( 11), (1.1)

R
Teff,oo = Teff( _%)1/2’ (12)
R00 = R(1—%)-1/2. (1-3)

1.3.2 Theoretical Conclusions

The similar spectral evolution and the much larger duration time and energy scales
suggest that the origin of the ignition is in much larger depth and that the ther-
monuclear runaway is triggered by carbon burning processes. The frequent type I
X-ray bursts as well as the steady state burning via the rapid proton (rp) process
typically produce a signiﬁcant amount of heavy ashes, so that the accumulated fuel
in this deeper layer only consists of a small fraction of carbon. The ashes from former
rp processes are in the range of the maximum binding energy per nucleon, so they
do not undergo themonuclear reactions and therefore do not contribute to the total
released energy. However, these heavy elements have much lower thermal conductiv-
ities which results in larger temperature gradients and allow the ignition of bursts

with much less accumulated fuel than it would be the case for pure carbon.

ff

1.4 Consequences of an Unstable Ignition

An unstable ignition of a parcel of hot plasma results in an immediate burning
of the fuel to heavier elements, in this ﬁrst step the immediate products of the
reaction. Depending on the total released energy during the process and the resulting
temperture of the plasma right after this nuclear conversion, further reactions of the
freshly created products can be triggert and can start a nuclear runaway. Having
reached the most bound state of baryons in the nuclei, the thermonuclear burning
stops, which corresponds to an Iron group elements. In the case of an unstable carbon
detonation the resulting upward-propagating shockwave is sufﬁciently strong to reach
freshly accreted material and trigger unstable He burning (Weinberg & Bildsten [15]).
Finally, a carbon explosion convertes all the material above the ignition layer to

maximum bound nuclei like iron and even heavier material (Schatz et al. [39]).

1.5 Previous Studies on Superbursts

Many years before the ﬁrst superburst was observed in 2000, unstable carbon ignition
has been proposed as an explanation for the burst activities of 4U 1820-30. Following
the theoretical ideas of thermonuclear instability for H / He burning by Hansen & van
Horn [28], Taam & Picklum [29] extended the calculations to higher temperatures and
densities which led to the conclusion of a instablitity in a carbon rich environment
being a probable cause for the observed bursts. Due to the lack of computing power,
these calculations were done with many simpliﬁcations, particularly with only a few
combined reaction rates and simple instability assumptions. In a very broad study
about the ocean and crust properties of accreting neutron stars Brown & Bildsten
[12] included more reaction rates to their stellar evolution code, giving a detailed
composition of layers below the H/ He burning region. By balancing heating and

cooling and calculating the effects of compression on the He burning ashes, mostly

carbon, they conﬁrmed the previous suggestion of carbon as an igniter of superbursts
and made precise predictions on the circumstances of high energetic ﬂares. The
instability criterion that was used in this and many other studies on unstable nuclear
burning, e.g. nova ignition (Shen & Bildsten [32]), was, that the heating rate is

greater than the cooling rate (thimoto et al. [31])

61%;? > %. (1.4)
In 2003 Narayan & Heyl [25] investigated burning stabilites by studying the evolu-
tion and consequences of small perturbations applied to the burning processes and
the thermal structure equations for H / He burning in the very outer layers. This
perturbation calculations require complex number arithmetics which increases the
computational effort by a lot and restricted their study again to a simpliﬁed reac-
tion network. More detailed work on stability analysis, e.g. Narayan & Cooper [30],

follow the perturbation method, but concentrate on type I X-ray bursts.

1.6 Goal of this Study

In this thesis the structure of the neutron star layer where unstable carbon ignition
occurs as well as a detailed instablilty conditions for carbon burning were investi-
gated. This was done with numerical calulations combining a full reaction network
and the complex number based perturbation instability analysis. The parameter
study covers the inﬂuence of the temperature, the energy ﬂux, the accretion rate and
the composition of the ignited layer on the stability of the burning process. Goal of
this thesis was to ﬁnd precise conditions under which a steady state shell burning
switches from a stable process to an unstable explosion, resulting in a superburst.
In addition to that a more detailed study on where and when in the burning process

the instability occurs is determined.

1:.-

Chapter 2

Introduction to Neutron Stars

2.1 General Properties

Neutron stars (NS) are very compact stars which consist in average of a large fraction
of neutrons in their interior. The density in such a star exceeds several times the nu-
clear density. With typical masses of 1.4MQ and typical radii of about 10 km (Haensel
[17]) neutron stars are the densest stars known. The large mass in comparison to the

very small volume (IO—ISVQ) results in an enormous surface gravity.

2
Egrav % 921- m 5 x 1053 ergs a: 0.2Mc2 (2.1)

GM _
ngZ-z—m2x1014cms 2 (2.2)
The mean density of a typical neutron star is

3M

47rR3 z 7 x 1014 gcm‘3 z 2 — 3p0, ' (2.3)

 

ﬂ:

where G is the gravitational constant, c is the speed of light and p0 is the mass

density in heavy atomic nuclei.

2.2 Structure

The interior structure of a neutron star can be separated into four main layers:
the very thin atmosphere, followed by a liquid ocean, the crust and the core. The
crystalized crust below the ocean can be divided in the outer and inner crust, as well

as the core which is remarkably large in comparison to other stars.

2 .2. 1 Atmosphere

The photosphere of a neutron star is a very thin plasma layer whose thickness varies
from 10 cm in a hot neutron star with a surface temperature of 3 x106 K down to a few
millimeter in a cold neutron star with a surface temperature of 3 x 105 K. Very cold
neutron stars can even have liquid or solid surfaces without a plasma atmosphere.
The emitted radiation of this layer formes a spectrum of thermal electromagnetic
radiation. Current models assume temperatures of T3 5 106K and very strong
magnetic ﬁelds (B 2 1011 G). This temperature is sufﬁciently high in order to start
nuclear burning processes, so most of the accreted hydrogen is burned into helium
in this outer layer of the neutron star. At the surface the ﬂux is therefore large
enough for the radiation force to exceed the very strong gravitational force which
leads to a turbulent and unstable material convection and a strong plasma outﬂow.
In hot unmagnetized stars where the radiative force is evoked by Thomson scattering
the radiation force equals the gravitational force when the luminosity exceeds the

Eddington Iimi t

LEdd = 47TCGMmp/0’T s 1.3 x 1038 (%> ergss‘l, (2.4)

8

where mp is the proton mass and U7- is the Thomson cross section. As the atmosphere
is very unstable and sensitive to the neutron star’s properties the actual models are

not yet completed.

2.2.2 Ocean

The ocean is often deﬁned as the region below the hydrogen and helium burning
layer, where the accreted matter is decelerated from its free-fall velocity, and above
the depth at which the material crystallizes. The ocean is typically ~ 100m thick
(Brown & Bildsten [12]) and ranges from densities of 105 — 109gcm_3. The fully
ionized plasma is liquid under these conditions (see calculations in this thesis). As
the accreted matter is already decelerated, the gravitational settling and the corre-
sponding released gravitational energy can be neglected in comparison to the much

higher nuclear energy release.

2.2.3 Crust

Outer Crust The outer crust which is about several hundred meters thick consists

of ions and free electrons. Its bottom is deﬁned as the layer at which neutrons start

22

to drip out from nuclei forming a free neutron gas at a density of p z PND
4 x 1011gcm"3. The pressure in the outer crust is mainly provided by electrons.
The thin surface layer can be up to a few meters thick in hot stars and contains a non-
degenerate electron gas. In deeper layers the electrons that are highly degenerate and
ultrarelativistic at a denstity p > 106 gcm‘3 form an almost ideal gas. At densities
larger than 104 gcm”3 the atoms are fully ionized by electron pressure. In deep
layers of the crust the ions form a strong coupled Coulomb system, which can be
liquid or solid. The deeper the layer the larger is the fraction of solidiﬁed matter.
As the electron Fermi energy rises with growing density beta captures in the atomic

nuclei occur and enrich nuclei with neutrons.

Inner Crust The inner crust spans the density range PND < p < 0.5p0 correspond-
ing to a thickness of about one kilometer. The matter consists of electrons, neutron

rich nuclei and free neutrons with a rising fraction of free neutrons for higher densi-

IF.‘.—_ 'TZ“

ties. The neutralization processes via ﬂ-capture softens the equation of state. The
deeper the layer the more important are the n — n-interactions. At the bottom layer
where the density reaches values of 1 / 3 — 1/2p0 the nuclei change into a non-sperical
structure and disappear at the crust core interface. In the inner crust free neutrons

and nucleons that are conﬁned in nuclei can become superﬂuid.

2.2.4 Core

Outer Core The outer core is several kilometers thick and extends from a density
of z 0.5p0 up to % 2p0. The matter consists mainly of neutrons with a small fraction
of protons, electrons and muons and is also called the npep. layer. This ripen plasma
is strongly degenerate. While the electrons and muons form an almost ideal Fermi
gas the neutrons and protons interact via the strong nuclear force and form a Fermi

liquid which can be superﬂuid.

Inner Core The inner core at densities from p 2, 2p0 up to a central density of
~ 4p0 has a radius of several kilopmeters and does not occurs in low-mass neutron
stars. The composition and equation of state is very model dependent. Eventually
strange matter condensates can be found in that layer. The main four theories of the

neutron star central structure are:
0 H yperonization of matter: appearance of hyperons, mostly 2‘ and A hyperons
o Pion condensation: boson condensate Od pion-like excitations
o Kaon condensation: BEC Of Kaon-like exctitations, “strange” matter

0 Phase transition to quark matter, mainly u, d and s quarks with only a small

admixture of electrons or even no electrons.

10

outer layers 4 .
0.3 — 0.5 km Inner crust

eZn

 
    
 

outer core
npeu

 
 
   
 

inner core

Hyperons?
7r condensation?
Kaon condensation?
Quarks?
I l s 0590 e 290 ~ 490 l*-—-9 — 10 km---—-'|'-—’I|
1—2km

 
   
     

 

Figure 2.1: Schematic structure of a neutron star
2.3 Schematic Thermal Proﬁle

The thermal proﬁle of an accreting neutron star follows a characteristic scheme.
An illustration of the structure can be seen in ﬁgure 2.2. Starting with a rising
temperature at the outer boundary due to gravitational energy ﬂux and nuclear
burning of the accreted matter and a constant outgoing ﬂux the temperature reaches
a maximum where the energy ﬂux switches from a net outgoing to an inward directed
ﬂux. The temperature below that layer is getting cooler and reaches a fairly cold,
nearly constant core temperature. In the region where the temperature is dropping
the ﬂux value is estimated to be 17% ~ —0.5MeV/nuc (Brown & Cumming [19]),

where mu denotes the nucleon mass and m the accretion rate per unit area.

11

"F—f—_ “T:

temperature

ﬂux

Schematic Proﬁle of an Accreting Neutron Star

 

r 7 l l T I l r

 

 

 

 

 

 

l l l l l l l

 

depth

Figure 2.2: Schematic proﬁle of an accreting neutron star

12

:FT—‘Tr

Chapter 3

Structure Of the Outer Layers

3. 1 Accretion

Neutron stars in binary systems can accrete material from a nearby non-compact
companion star. The composition of the accreted material depends on the type of
the companion, but is typically a mixture of hydrogen and helium. The amount of
accreted material ranges from 10'10 —- IO—SMQ and results in local inﬂux on the

23"1 under the assumption that the material is accreted

surface of 104 — 105 gem—
across the entire surface area. This accretion rate is on the order of 10 — 30% of the
Eddington limit MEdd which denotes the mass accretion rate at which the accretion
luminosity is equal to the Eddington luminosity. The permanent accretion of light
elements like hydrogen and helium is essential to provide a steady state composition
in the outer envelope which is the basis for the investigation of carbon burning in
the way it is done in this thesis. Although the accretion rate is on the order of 1 / 10
of the Eddington limit the total accreted mass in comparison to the neutron star’s
mass is fairly small. The mass impact of the material on the star can therefore be

neglected which means that the assumption of hydrostatic equlibrium is still valid.

As the accretion ﬂow onto the neutron star is in general asymmetric it is more

13

convenient to use the local accretion rate, deﬁned as the mass accreted per unit area.

The local Eddington rate is deﬁned as

'mEdd = #:7712136 = 7.5 x 104gcm_28—1,ue (

 

 

(3.1)

10km
R i

where 0T is the Thomson cross section, 0 is the speed of light, R is the radius of the

star, mu is the atomic mass unit and the number density of electrons is p/(nemu).

3.2 Equations Of Hydrostatic Equilibrium

As neutron stars are very dense and massive one needs to take effects of general rela-
tivity into account. A very important coefﬁcient is the compactness parameter r9 / R
where R is the Schwarzschild radius. Typical values of the compactness parameter
are very high for netron stars: 0.2 — 0.4. All other stars have values much smaller
than one, i.e. white dwarfts (as 10‘4) or main sequance stars (x 10'6). The metric

that is used to describe the space around the neutron star is
ds2 = C2dt262q) — «‘32’\dr2 — 'r2(c102 + sin2 6d¢2), (3.2)

where t and r are the time and space coordinates, 0 and ()5 are the polar and azimutal
angle and (I) and A are the curvature parameter that are described in detail below.
The angular geometry with respect to 0, (b is the same in ﬂat and curved space-time
because of the spherical symmetry. However the space-time is curved along the time
or space dimension. This space-time curvature is evoked by massive stars. For a
ﬂat space-time the curvature parameter are zero (<I>(r) = /\(r) = 0). For a curved

space-time /\ determines the curvature in radial direction which leads to a change in

14

length

T
12/0 eAdr'. (3.3)

For many calculations it is more convenient to use the mass as a function of radius
m(r) instead of the curvature parameter A(r). This transformation can be done with

the following equation

 

e_’\ = \/1 — 2Gm/(rc2), (3.4)

where m represents the gravitational mass inside a sphere of radius r. In neutron
stars the gravitational mass is smaller than the baryon mass (“rest mass”) due to
the gravitational mass defect.

The parameter (I) corresponds to the time curvature and inﬂuences the proper
time dr = e‘p(")dt. It describes the gravitational redshift. Assuming a periodic
signal with a frequancy u, an observer at inﬁnity will measure a modiﬁed freqency

1100 = u(r)eq’(") or a redshift Of

z(r) = 11—0-2)- — = 641)“) — 1. (3.5)

Voo

Making the assumption that a neutron star’s matter is a perfect ﬂuid (all stresses
are zero except for isotropic pressure), there are four functions of 1' that have to be
solved: the two metric functions <I>(r), Mr), the pressure P(r) and the mass density

p(r). With the Einstein equation
1 87rG
Rik - 59M? = 7TH: (3-6)

we get three relativistic equations for hydrostatic equilibrium for a spherically sym-

15

metric neutron star

 

 

 

3 —1
if: = _szm 1+5) (1+47rPr (1_ 20m ,
dr r 5 me2 czr
$1? = 47r7‘2p, (3.7)
dd) _ 1dP 1 P ‘1
dr — 5dr 8 ’

where the ﬁrst one is the Tolman-Oppenheimer-Volkoff equation (TOV) of hydro-
static equilibrium. In order to solve these equations they have to be supplemented
by an equation of state. As the ﬁrst two equations do not depend on the third one,

they can be solved seperately. In quasi-Newtonian form the TOV can be written as

9‘5 -— — (r) — Gm
d7“ _ gp, g _ r2\/1— 2Gm/(r02)’

 

 

(3.8)

where g(r) denotes the local gravitational acceleration. In the stellar interior applies
P > 0, dP/dr < 0 with the bounardy conditions P(r 2 R) = 0, p(r 2 R) = 0 and
m(r 2 R) = M where M is the total gravitational mass.

As it is more convenient the column density is used as the independent parameter

of the differential equations instead of the radius. The column density is deﬁned as

y — fr Pd? (it a: —- —p- (3.9)

The equations of hydrostatic equilibrium (equations. 3.7) can then be transformed

to the following ones

 

 

dP Gm P 47rPr3 2Gm ‘1

— = — 1 — I l— 3.10
dy r2(+€)(+m02)( 027'), ( )
dm 2

—- = — .1
dy 47rr, (3 1)
(id) ldP P ‘1

_ = __— 1— — . 3.12
dy e dy ( e) ( )

16

3.3 Equation of State (EOS)

In order to examine the interior properties of a neutron star one has to solve the
equations of hydrostatic equilibrium. These have to be supported by an equation
of state that relates the pressure of the material to the actual density. The basic
assumption to make in order to ﬁnd an EOS is local thermodynamic equilibrium
(LTE). This assumption holds for a large range of the star because the particle-
particle and particle-photon mean free paths are very short and the collision rates
very rapid in comparison to other stellar lengths or times. Typically the LTE breaks
down in the stellar atmosphere, but this layer is not of interest in this thesis. The
typical length scale in which the surrounding conditions change signiﬁcantly is the

pressure scale height

61nP —1 P

 

which is many orders of magnitude higher than the mean free path. With an equation
of state one normally integrates starting from the center of the star to the surface.
As the interior of a neutron star, espacially the EOS of the neutron star’s core is
not known very well, examinations of the outer envelopes, especially the ocean, start

with intergrations from the surface to deeper layers.

3.3.1 A Polytropic Gas

In the outer envelopes the plasma can be considered as a polytropic gas with index
n. The pressure depends on the density like P oc p1+1/". In hydrostatic balance the
density can be described as a function of radius by the following equation (Brown &

Bildsten [12])

_ -5" =i0_.
p‘p°(1 A)’ A (n+1)P0’ (3'14)

17

ll—

where A ist the total height of the polytropic atmosphere and n is the polytropic
index, which is about 4/3. The constants with subscript 0 are the values for r = 0.

Inserting A in the density equation gives

M)”. (3.15)

p=po 1-
(

For a given minimum and maximum density (Pmim pmax) at the beginning and end

of the ocean region the fraction PO/po can be calculated with

1:9- «Riff/n”)
5m

 

 

__ (3.16)
1/
P0 (n+1) (( min) ”rout—Tin)
and
l/n “'7’
{(5.5) - 1}
P0 = pmin 1 - ) (3°17)

1/n
(2m) Tout "' Tin

pmin
where Tout, Tin are the outer and inner radius of the integration. The transformation
to column depth ﬁnally yields

2 (n+l)PQr n+1
p09 (1’ 9/30 )

y 2 fr p(r’)dr’ = POW + U2 . (3.18)

 

3.3.2 Numerical Fit Of the EOS

For stars EOSs are usually tabulated and subsequently interpolated between mesh
points. But in order to do analytical calculations it is of great interest to have an
analytical function for the EOS. There are several models for the EOS in deeper

layers that differ slightly for high densities above p = 1012 gem-'3. As the densities

18

l'F—‘T‘n T"

(I). I (11'.

 

 

I

1 6.22 10 11.8421
2 6.121 11 -22.003
3 0.006004 12 1.5552
4 0.16345 13 9.3

5 6.50 14 14.19

6 11.8440 15 23.73

7 17.24 16 -1.508
8 1.065 17 1.79

9 6.54 18 15.13

 

Table 3.1: Parameters for the EOS ﬁt function

in the carbon burning layer are much lower the differences in the EOS between these

models can be neglected. For densities above 105 gcm’3

numerical approximation
with a one-parametric function works well, for lower densities that can be found in
the atmosphere the EOS also depends on the temperature and therefore cannot be

described as a one-parametric function any more. The ﬁt function consists of several

fractional-polynomial parts, matched together by using the function

 

fo(-'L‘) = - (3.19)

The EOS as well as the ﬁt function differ between rotating and non-rotating stars.
All the calculations in this thesis were made with the non-rotating equations. For
the EOS it is instructive to describe the pressure as a function of density. Using the
two abbreviations E = logp and C = log P where the density and pressure are given

2

in gcmT3 and dyn cm— , respectively, the parameterized function is

3
a1 +laifa+ a3€_f0 [a5(£ - a5)] + (a7 + agofo [a9(a10 — 6)] (3.20)
45

+(aII + 012€)fo [013(014 - E)l + (015 + a165)f0 [017(018 - 01- (3-21)

C:

The ﬁtting coefﬁcients a, can be found in table 3.1 (Haensel & Potekhin [33]).

31 r I I I I I

30 * FPS ﬁt —— ..........

29 [- SLy ﬁt ””””” ﬂ
linear ﬁt ............

28 — 4

log P [dyne cm’z]

20

Fits for the EOS (FPS, SLy) and Linear Fit

 

 

_

_

 

l l

 

 

8
103 p [s cm“3l

9 10 11 12

Figure 3.1: Fit functions to the EOS. The abbreviations FPS and SLy correspond
to the Friedman-Pandharipande-Skyrme and the Skyrme Lyon effective interaction
model (Pandharipande & Ravenhall [34] and Douchin & Haensel [35]).

As the EOS in the outer layers of a neutron star is primarily determined by degenerate

electrons, the pressure can be described using a simple power law as a function

of density, and the temperature dependence can be neglected. The pressure of a

non-relativistic degenerate electron gas scales as P8 oc p5/3, for the relativistic case

P3 or p4/3. Including coulomb corrections the EOS can be ﬁtted very well to the

complex ﬁt function mentioned above by the following simple power law (ﬁgure 3.1)

g = (1.378 :t 0.004)§ + (14.10 :I: 0.03),
P = 1.259 x 1014 pl'378,

p = 5.858x10“11P0'7257.

20

(3.22)
(3.23)

(3.24)

rrr—*— _ “

3.3.3 Coupling Parameter

In real gases the interactions between the particles have to be taken into account. A
measure of the interaction energy between ions is the Coulomb potential between two
charged particles. An ionic charge of Z for both particles separated by a distance a
results in a Coulomb potential of Z 2e2 / a. Electromagnetic interactions are expected
to become important when the electrostatic energy is comparable to the thermal
energy kT. The ratio of electrostatic and thermal energy

_ (Z)262

_ akT (3.25)

 

is the so—called coupling parameter, where PC = 1 is the rough demarcation at which
Coulomb effects become important. The mean distance a between the ions is chosen

such that
'3'71'0. = — = —, (326)

where (A) denotes the mean ionic mass number and N A is Avogadro’s constant. The
relation between the thermal and electrostatic energy can also be used to descibe the
phase of the material. In the range of 1 5 1’0 ,3 170 the plasma is in a liquid phase.
For very low values I‘c << 1 the gas can be considered as an ideal gas, values above

~ 170 — 180 cause the material to crystalize.

3.4 Thermodynamics

3.4.1 Heating

The two main heating processes in the outer layers are compressional heating due to

gravitational energy of the accreted material and the energy release due to nuclear

21

Gravitational Flux
0000014 I I I IIIII' I I I IIIII' I I I IIjII

0.00012
0.0001
8e-05
6e-05
4e-05
2e-05

0

 

F gm, [keV]

 

 

 

 

'28‘05 1 1 111111] 1 1 1L1111l 1 #411111

1e+10 1e+11 1e+12 1e+13

 

column depth [g cm

Figure 3.2: Gravitational ﬂux

burning. With a typical radius of 10km and a mass of 1.4MQ an accreted proton can

release

M
cgrav z G—rﬂ z 3 x 10’4 ergs z 200 MeV, (3.27)

which is much higher than the net energy relaesed in nuclear fusion from a proton to
a helium nucleus (6.4 MeV per nucleon). The temperature proﬁle of the atmosphere
is therefore strongly inﬂuenced by the gravitational heating of the accreted material.
However, in deeper layers Where carbon ignites the contribution of the gravitational
settling to the energy ﬂux can be neglected in comparison to the nuclear burning

energy of roughly 0.5 MeV/nuc (see ﬁgure 3.2).

3.4.2 Cooling

The cooling is provided by several different processes, namely the heat transport

by radiation, by thermal conduction and by emission of neutrinos. The neutrino

22

energy release at the given temperatures can be neglected and not included in the
calculations.

Considering the thermal transport through the outer envelope the amount of
energy taken away from this layer is determined by the energy ﬂux F which obeys

Fick’s law and can be written as

F = ——KVT, (3.28)

where T denotes the temprature and K the thermal conductivity. Heat transport is
given by both radiative (mainly Thomson scattering and free—free absorption) and
conductive (electron-electron and electron-ion scattering) processes. The thermal

conductivity is then

4ac

K =
3p’crad

 

T3 + I(conduction, (3'29)

where lewd is the radiative opacity and Kconduction the electron conduction contri-

bution.

Radiative Opacities The radiative Opacities are determined by electron scattering

and free-free absorption.

Conductivity The conductivity can be written as (Wiedemann-F‘ranz law)

3m; ye

. (3.30)

K conduction =

where m; = max/1 + :I:2 is the effective mass of an electron and :1: is the relativistic

parameter given by

1/3
a: = if; m 1.009 (psi—i—i) . (3.31)

23

V5 = um- + Veg is the effective electron collision frequency, where the index ei denotes
the electron-ion and es the electron-electron collisions, respectively. The ﬁrst one
can be described (Yakovlev & Urpin [22]) by

_ 4miljZe4

V6,; — 37Th3 (3.32)

 

with L being the Coulomb logarithm, a slowly varying function of the density and
temperature. After (Potekhin et al. [23]) the Coulomb logarithm can be calculated
by

 

_ 2kF 933(9)F2(9)R(9) _Efg
L— [0 dq 946(9) (1 416%), (3.33)

where 6 = :r/ V1 + 172, hp = pp/h, 6(q) is the static longitudinal dielectric function
of electrons, F (q) is the nuclear form factor and R(q) is a non-Born correction factor.

The effective electron-ion collision frequency can be expressed by

_ 3ﬂ2 (19137")2
68 — EWJ(xayl- (3'34)

9 = x/3Tpe/T, where Tpe = hwpe/kB is the plasma temperature determined by the
plasma frequency Lupe = (Macaw/mg. oz is the ﬁne structure constant and again
ﬂ = 50/ v1 + 11:2. J (y) can be calculated using the ﬁt function by Yakovlev & Urpin

 

[21]
6 2
J($,y) ~ (1+ 33:7 + 5734-) X (3.35)
3 2
y 2.81 0.81vF
l 1 —— — —— . .
[3(1+ 0.07414y)3 n( + y ) y .22] (3 36)

24

Ii“

3.5 Thermal Structure Equations

Having discussed the heating and cooling mechanisms the differential equations that
allow to integrate the thermal structure of the neutron star shell can be deﬁned. The

entropy equation is given by

ds 1
_=__ . 3.7
Tdt pV F+e, (3)

where F is the ﬂux and c is the sum of all sources and sinks of energy. Using the
equation of hydrostatic balance, the ﬂow velocity of the accreted material I) = r'n/ p

and the thermodynamic relation

(as—(29(3).

equation 3.37 and the ﬂux equation from Fick’s law (3.28) can be transformed to the

more useful forms

 

6F __ 0T . 0T CpTr'n

3y— — Cp (at + may) y Vad E, (3.39)
6T F

0y _ pK' (3.40)

Here Cp denotes the speciﬁc heat at constant pressure and

 

__ alnT
Vad : (alnp)s (3.41)

is the adiabat. The speciﬁc heat in general has contributions from the ions, the

electrons and the radiation.

 

06 66 ('96
cp = <—) +<—> +(-—)
6T ions 6T electrons 6T radiation
3NAk 7.21.2 NA(Z)T 4aT3
(2<A>)+(xm.c2 <A> l+( p ) (3'43)

25

IE

dF/dy error

1 I rIiIIIIIl I I rIIIIIl I IjIIII"

0.1
0.01
0.001
0.0001
1e-05
1e-06
1e-07
1e-08 . .......I ..... .I .

1e+10 1e+11 1e+12
2]

 

relative dF/dy error

 

1e+13
column depth [g/cm

Figure 3.3: Error in equation 3.39 due to an idealized speciﬁc heat

where a: is the relativity parameter a: = pF/(mc). In the speciﬁc heat above an ideal
gas is assumed. For the polytropic gas there are corrections that are neglected here.

The internal energy of a Coulomb liquid can be ﬁtted by the function (Hansen [37])

 

, = 3/2 AI As
u, r [anrrr+1+r]’ @AO

where A3 = —\/3/2 — Al/m, A2 = 0.62954, A1 = —0.9070 (Chabrier & Potekhin
[36]) and I‘ is the coupling parameter. The speciﬁc heat derived from that ﬁt function
differs by a factor of a few in comparison to the ideal gas speciﬁc heat. As the
contribution of the terms with speciﬁc in equation 3.39 are very small relative to the
released energy 6 the relative error to dF/dy is small (see ﬁgure 3.3). The feature
at a column depth of z 2 x 1012 g cm‘2 is due to an electron capture on 23Na (see

section 6.4).

26

Chapter 4

Nuclear Processes

4.1 Triple-oz Process

In contrast to the hot CNO cycle in the triple-a process
4He + 4He + 4He _. 120 (4.1)

there is no lack of neutrons, so weak interactions are not needed, which accelerates
the process, because the weak interaction coupling strength (01 ~ 10’s) is much
smaller than the coupling strength for strong interactions (as ~ 1). On the other
hand helium is inhibited by the fact that there are no exothermic two-body reactions
in which only 4He is involved. The simplest reaction would be 4He + 4He H 8Be but
the mass of 8Be is 92 keV higher than the mass of the 4He nucliei, so the 8Be decays
with a lifetime of r = 1 x 10‘ 16 s back to 4He+4He. With a density of p ~ 105 g cm—3

and a thermal energy of kT ~ 10 keV the abundance of 8Be is very tiny

3332 ~ 10’9. (4.2)
n4He

27

Using this small abundance it is possible to produce 12C through the reaction
4He + 8Be —-> 12C + ’7, (4-3)

which is an exothermal reaction with a released energy of Q = 7.366 MeV. Because
of the small 8Be abundance the 12C production is very small. However this rate can
be increased if 12C has an excited state 12C* at a similar energy level. This excited
state has an energy of 7654 keV above the 12C ground state, 283 keV above 4He+8Be
and decays mostly via an a decay returning the original 8Be nucleus. Nevertheless
a fraction of 10‘3 of the 12C* decay results in stable 12C in the ground state. The

irreversible production of 12C ﬁnally proceeds through

34He ——> 4He + 8Be —>12C* —-> 12C + ’7 + '7. (4.4)

4.2 Rapid-Proton Process (rp-process)

For high temperatures and accretion rates at the surface of the neutron star the
accreted H / He mixture does not burn in separate layers from H to He and then via
the triple alpha process to carbon, but via the rapid proton process from H / He to
very heavy elements beyond the iron group, A ~ 60 - 100 (Schatz et al. [39]). In
this fusion chain the rate of capturing protons is larger than the ﬂ—decay rate, which
leads to neutron poor isotopes as fusion products, located above the valey of stability
in the chart of nuclei. The larger the accretion rate and the proton capture rate to
the nuclei, the fewer neutrons the isotopes have in comparison to their stability valey
isotope.
At high accretion rates the temperature at the ignition depth of the fuel exceeds
T Z, 8 x 107 K. The accreted material burns via the hot CNO cycle and reaches

H e ignition conditions before the H is fully consumed, which leads to helium fusion

28

in an hydrogen rich environment (Lamb & Lamb [38], Taam & Picklum [29]). For
subsolar metallicities and accretion rates of 10—10M9yr_1 f, M S, 10’8 M9 yr”1
the burning is unstable and results in a type I burst. If the accretion rate is high
enough for stable burning the helium is fused within 20 min. The enhancement of
CNO seed nuclei combined with a so-called breakout from the CNO cycle leads to a
rapid hydrogen burning via the rp-process. In that case the helium ignition acts as
a trigger for rapid hydrogen burning. At an accretion rate of m z mEdd the ﬁnal

carbon abundance is z 4.1% by mass (Schatz et al. [39]).

4.3 Carbon Burning and Advanced Burning Stages

After the irreversible ignition of helium to carbon the carbon burning processes can

occur. The exothermic reaction
4He + 120 _. 160 + 7, Q = 7.162 MeV (4.5)

competes with the triple-oz process because its reaction rate is linear in the 4He
concentration while the triple-oz rate is proportional to the third power of the 4He
concentration. The helium burning stage thus generates a mixture of 12C and 16O.

In addition to that two more exothermic carbon burning reactions, namely

12C + p —> 13N + 7 Q = 1.943 MeV (4.6)

12C + a ——> 150 + n Q = 7.162 MeV (4.7)

take place. The ﬁrst of these two reaction rates is many orders of magnitude higher
than than any helium reaction. However this reaction is strongly suppressed by the
low number of available protons in the composition at that stage. The second one is

Very low and only comes into play at temperatures around 109 K.

29

Although the reaction rates of the above processes cannot be neglected, in the
stellar interior the triple-oz can convert most of the helium to carbon before further
burning branches can fuse 12C+4He. That results in high carbon abundances before
the temperature and density reaches high enough values to ignite carbon and burn it
to heavier elements. The amount of produced oxygen, however is a priori not small.
But the fact that at temperatures above 3.5 x 108 K the destruction rate for oxygen

via
4He + 16O —» 2ONe + 7 Q = 4.73 MeV. (4.8)

is higher than the production rate (equation 4.7) burnes most of the produced oxygen
to neon which results in a remarkable neon abundance.

A comparison of the reaction rates of the main burning processes can be found in
ﬁgure (4.1). Note that these are only the reaction rates. These have to be multiplied
by the abundances of the input isotopes which affects the branching ratio of the
products.

As the composition resulting from nuclear burning is dominated by carbon the
reactions in which only carbon is involved play an important role. The two main

reactions in a carbon rich environment are

12C + 120 —> p + 23Na Q = 2.242 MeV (4.9)

120+120 —. 4He+20Ne Q=4.621MeV (4.10)

Comparing the reaction rates for 23Na and 20Ne shows that the ratio of the two
rates does not vary very strongly in the temperature range of 108 - 109K and is
roughly one. That means that roughly for every four consumed carbon one of each
p,23Na,4He and 20Ne is produced. The ratio of both reaction rates can be seen in

ﬁgure 4.2

30

$21.]

reaction rate N A (av)

Sample Burning Rates

 

 

 

 

‘ p+b12rqINi2'_-I_;_ I r I l I Iﬁ r:

100000 p + 012 —> F17 ------ 1

He4 + C12 —-> 016 ............ _

1 He4 + He4 + He4 —> C12 ---—- _._1

H64 + 016 —-> N820 -----~---- .

C12 + C12 ——> He4 + Ne20 ------------ .

1e-05 ,3

1

1e-10 _:

1e-15 d:

1e-20 _E

1

1e-25 a

18'30 4 l d... 1 I 1 1 1 I 1 .33: 1 1 1 1 1 1 -
1e+07 le+08 1e+09

temperature [K]

Figure 4.1: Sample rates for some important reactions

Carbon Depletion Rate Ratio

1'3 I I I I I I I I

 

 

 

 

1.23 I I I l I I l 4
0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Temperature [GK]

 

Figure 4.2: Rate Ratio of producing p + 23Na and 4He + 20Ne

31

4.4 Released Energy

The energy that is released in nuclear fusion is because of the slightly lower mass of

the products relatively to the reactants. The Q value is given by
Q = AE = AmCZ. (4.11)

This released energy, which is due to the different binding energy of the nucleus, can

be described by the Weizsiicker mass formula

M(Z, A) = N Mn + Z Mp + Zme mass of the constituents (4.12)
—avA volume term
+asA2/3 surface term
22
+a°:4_1/_3- Coulomb term
(N - 2)2
+aaT asymmetry term
6 . .
+2”? pairing term

where the coefficients av, as, ac and aa depend on the range of masses for which they

are optimised. One possible set of parameters is (Povh [18])

av = 15.67 MeV (4.13)
as = 17.23 MeV
Co = 0.714 MeV

aa = 93.15 MeV

—11.2 MeV for even Z and N (even-even nuclei)
5 = 0 MeV for odd A (odd-even nuclei)
+112 MeV for odd Z and N (odd-odd nuclei).

32

IF“

4.5 Electron Capture

The basic process that occurs during an electron capture to the nucleus is the inverse

neutron decay
p+e_ —> n+ V6. (4.14)

The reaction in an isolated environment is endothermal and therefore does not occur
in a system with low density. However if the reacting proton is bound in a nucleus
the reaction and the corresponding switch from a proton to a neutron inﬂuences the

binding energy of the nucleus. The reaction
e_+(Z,A) —> (Z—1,A)+z/t3 (4.15)

can be endothermal or exothermal depending on the change in binding energy. The
nuclear Q-value is the mass energy difference of the product and the reactant nucelus

plus the captured electron’s rest mass
0 = (lM(Z — 1. A) — M(Z, A)l — me) oz. (416)

The density at which an endothermic electron capture can occur can be obtained from
the Q-value by expressing the electron chemical potential in terms of the relativistic

parameter r F

#6 = mec2 [(1+ 33%)1/2 — I] . (4.17)

33

4.6 Reaction Network

In order to calculate the abundance changes and the released energy during the
burning the set of differential rate equations have to be derived. The continuity

equation

5m + V ° ("2'”) = -dz' + 102'. (4-18)

is determined by the destruction rates d,- and production rates p,- for each isotope 2'

due to thermonuclear reactions. n, is the particle density of isotope i which can be

written as
X-
nz‘ = jN/Ip E YINAP (4-19)
2

with X,- being the mass fraction and A,- the atomic mass number. Because the number
of nucleons and the density in one parcel in the star are conserved, the total particle

ﬂux in a volume dV must be 0
Otp + V(pv) = 0 4:» Z X,- = 1 (4.20)

with p = Z,- niAi/NA. Ignoring particle diffusion, all isotopes comove vi = v and

equation (4.20) can be subtracted fromrequation (4.18) which yields

8.1/.- + v - W,- = (4.21)

P P
— 2(pp‘1NA r(T) H Y1.) + Z(pp_1NAT(T)l—[ Y1.) (4-22)
’D k=1 ’P k=1

For each isotope one has to sum over all the reactions that produce or destroy this
isotope. The set of reactions is indicated by the index ’D for the destruction processes

and ’P for the production processes, respectively. In each process the abundance

34

 

type reactants products

1 el —» 62
2 e1 —> 62 + 83
3 81 -—* 82 + 83 + 64
4 61+ 82 ——> 83
5 el + 62 —-> 63 + e4
6 61+62 ——+ e3+e4+e5
7 el+e2 —-> e3+e4+e5+e6
8 61 + 62 + 63 —-> 64
9 el+e2+e3+e4 —+ e5

10 el+ez+e3+e4+e5 ——> 65

11 61 —-> €2+e3+64+e5

 

Table 4.1: Reaction types used in the nuclear network

product of the reactants has to be computed. The density appears one power less
than this number. The reaction rate r(T) is only a function of temperature. Care
must be taken if reactants or products appear multiple times in a reaction, i.e. 4He+
4He+4He —-> 12C, which will be discussed later in this section. The network includes
the reaction types mentioned in table 4.1

In order to solve equation (4.21) it has to be linearized in time. As the system
of diffenrential equations is very stiff, an implicit integration methods has to be
used, which requires the calculation of the Jacobian matrix J , where Jij = (972/619.
The reaction rates are tabulated by reaction and so it makes sense to go over all
the reactions and add the contribution of this reaction to the reactant and product
entries in J. As an example, consider a reaction of type 5: a + b —> c + d. This

reaction makes a contribution in four abundance equations

Y. = —Yaprr(T) (4.23)
Y, = —Yaprr(T)
Y0 = +Yaprr(T)
I’d = +Yaprr(T)

35

and there are two factors of Y appearing on the right-hand side, so there are 8 terms

inJ

Jaa = — prr(T) (4.24)
Jba = — prr(T)
Jab = — Yapr(T)
Jbb = — Yapr(T)
Jca = — prr(T)
J61, = — apr(T)
Jda. = — prr(T)
Jdb = - aPT(T)

in both sets of equations the ellipsis (...) indicate the contributions from other reac-
tions.

The case of multiple reactants or products of the same isotope requires a special
treatment of the reaction rate. Supposed there are N powers of one abundance on
the right-hand side of the abundance time derivative. Then the reaction rate has to
be devided by the statistical weight of N I. In addition to that this reaction destroys
N particles of this isotope, so the rate has to be multiplied by N. This algorithm
can be done easily by the following calculation procedure:

Let each rection have a number of N reactants and M poducts. Replace each
reaction rate r(T) by wr(T) where w is the statistical weight. Then compute the
derivatives and the Jacobian matrix without taking care of the multiplicity of the
isotopes. This procedure sums up all the quantities correctly.

The reaction network includes the elements p, n, 4He, the isotopes shown in ﬁg-
ure 4.3, 56Fe and 104 Ru. More detailed information of the relevant element properties

is listed in appendix C. The elements are connected through 68 weak and 470 strong

36

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

22AI 23.41 24A1 25A1 26AI 27AI
20Mg 21Mg 3171;23Mg 24Mg 25Mg 26Mg
20Na 21Na 22Na 23Na
l7Ne “Na ”Na 20Ne 21Ne 22Ne 23Ne
17F 18F 19F
130 140 150 160 170 180
12N 13N 14N 15N
120 130 14C

 

 

 

 

 

 

 

Figure 4.3: Element chart of included isotopes without p, n, 4He, 56Fe and 104Ru

reactions, taken from the J INA reaction library (http: //www . nscl .msu . edu/"nero/db/).

4.7 Electron Screening

At sufﬁciently high temperatures all (or practically all) of the atoms are ionized. The
bare nuclei surrounded by the free electrons collide with the kinetic energy arising
from the thermal motion. For the two (in some reactions even three) nuclei to
undergo a nuclear transformation they have to approach each other up to a distance
at which a reaction in the nucleus can take place. This distance is of the order
of 10—13cm, the nuclear radius. Approaching each other the repulsive Coulomb
force increases strongly so that the potential energy is much larger than the thermal
energy kT when the distance of the nuclei is in the order of the nuclear diameter.

The important factor in the reaction is the barrier penetration factor, the probability

37

of the nulcei approaching sufﬁciently close for nuclear forces to become important.

The reaction rate r is proportional to (Salpeter [24])

r oc /000 131/2 exp (TEE?) P(E)anuc(E)dE. (4.25)

The ﬁrst term is the Maxwell-Boltzmann distribution factor, the second (P) is the
barrier penetration factor which depends strongly on the kinetic energy E and the
charge of the approaching ions Zi. June is a nuclear factor which depends on the
details of interaction after barrier penetration and usually varies slowly with energy.

Normally the penetration factor is calculated only with the Coulomb potential
of the bare nucei with unscreened charges Zi. In a neutron star the density is very
high, the mean distance between the ions and electrons is pretty small. Even though
completely ionized each ion polarizes the surrounding gas by attracting electrons
and repelling neighbouring ions. The ion is then completely screened by a cloud of
electrons. The radius of this cloud is in the order of the mean distance between the
ions and differs depending on the ratio of the thermal energy [CT and the repulsive
Coulomb force. The ion takes its electron cloud with it while moving through the
plasma which changes the forces and, of course, the interaction energy between two

approaching ions. The total interaction energy can be written as

Z1Z2e2

+ U (7‘12). (4.26)
7‘12

 

Utot =

where ZlZ2e2/r12 is the unscreened Coulomb energy and U (r12) is the screening
term.
Let deﬁne Emax to be the energy at the classical turning point where two positive

charges with kinetic energies smaller than the penetration energy have the smallest

38

distance:

2 z 2
Emu: 1 26. (4.27)

To

 

Let furthermore a be the distance deﬁned by the relation
47ra3pNA = 1, (4.28)

where p is the gas density and N A is Avogadro’s number. The distance a is a measure
of the interparticle distance, containing an average mass of 1/3 atomic mass units
inside a sphere of radius a. Let R be the radius of the screening cloud.

Let us consider a case where the classical impact parameter rc is much less than

the the screening radius R. The barrier penetration factor mainly depends on

Z z e2
E — U(r12)— In: . (4.29)

 

Now U (r12) must be a function that is small for large distances (r12 > R) and
becomes constant (U0) for distances r12 << R. The value U0 will be of order of

Z1 de2/ R. Horn that it follows that
— ~ — << 1. (4.30)

If this inequality is satisﬁed the screening potential U (r12) can be replaced by U0
which is independent of both E and r12. Therefore the penetration factor P and
the nuclear factor June for energy E without screening are the same as for energy

(E + U0) with screening. The reaction integral can then be replaced by

rscreen 0c [000 [(E + U0)1/2 exp (WET) exp 0%)] P(E)anuc(E)dE. (4.31)

39

As U0 is much smaller than Emax the term (E + U0)”2 can be replaced by E 1/ 2 and
the whole screening effect reduces to a factor exp(-—U0/kT) which the unscreened

reaction rate has to be multiplied with;

U
Tscreen = 7' EXP (—%) - (4.32)

40

Chapter 5

Numerical Methods

5.1 Basic Numerical Concepts

Concerning the stability of nuclear burning, there are two different approaches of
calculation.

In one approach, one simulates the physics of the accreted material in a fully
time-dependent reaction network, including a large number of different nuclear re-
actions and detailed thermodynamics. This kind of calculation is important for the
understanding of the details of these thermonuclear explosions. A very popular ex-
ample for codes like this is the FLASH program, written at Chicago. However, this
method is not very convenient for parameter surveys and the comparison of theoret-
ical predictions of bursts.

In a second approach one starts with a stellar composition in quasi-equlibrium
and focuses on the thermonuclear instability that triggers the explosive burning. As
the general structure of the accreted layer varies only slowly with time one ﬁrst calcu-
lates a steady state solution for the composition before applying small perturbations
on the system. If these perturbations grow with time, the burning process in this

environment is considered to be unstable and leads in the temperature regime of a

41

neutron star’s outer layers to a X-ray burst.

5.2 General Properties Of the Code

The program that was used for this work is written in the object oriented program-
ming language C++ and divided into two main parts. The ﬁrst part is the solver for
the nuclear reactions, the second part consists of the complex number based stability
calculation.

The nuclear network solver is called afterburner and was written by Dr. Edward
Brown. It integrates the stiff set of differential equations with an implicit Bulirsch-
Stoer integrator and includes strong and weak reaction rates, as well as electron
screening effects. The integrator calculates the change of the elemental abundance

vector dY for a given density, temperture and burning time.
dY = afterburner.Burn(time,rho,Y,T)

The main integration variable in all the calculation is the column depth y which is

related to the burning time like

% = m (5.1)
for an inward moving parcel of gas in the Lagrangian coordinate system.

The second part of the program consists of two parts, the calculation procedures
for the steady state solution and the complex stability analysis. Both parts integrate
from outer to inner layers using the column depth as the main integration variable.
The stepsize was choosen constant in logarithmic scale. The thermal structure equa-
tions do not form a set of stiff ODEs and therefore were integrated with a ﬁxed step
fourth order Runge-Kutta integrator.

The ﬁrst and second part of the program are combined as follows: First the

42

afterburner routines determine the change of the element abundances. From that
the total released energy was found by calculating total difference in mass excess of
the system. With the gained energy the thermal structure including the outgoing
energy ﬂux, the temperature and the density were computed. The reason why the
program was split up in these two parts is that the thermal structure equations can be
changed for tests and in order to check the inﬂuence of the single terms without the
need of changing the total Jacobian matrix for the implicit integrator. The following

sequence shows the calculation in pseudo code:

for(coldepth=c0; coldepth=c1; dy)
I

// abundance changes
dtime = dy/mdot
dY = afterburner.Burn(dtime,rho,Y,T)

// thermal structure
d1? =1 ...

dF = ...

drho = ...
rungekutta(dT,dF,drho)

For the stability analysis the thermal structure equations are complex functions. As
the perturbations applyed to the system are small the magnitude of the imaginary
parts of the equations can be neglected in the abundance change calculations. The
afterburner routines therefore were used with the real parts of the thermal equa-
tions.

The program ﬁrst calculates the steady state solution before systematically com-
puting the complex thermal structure for a range of complex perturbation parame-
ters.

As input variables the program uses the temperature, the density, the inital abun-
dances and the column depth range. The output for the steady state solution is a

thermal proﬁle as a function of column depth. The stability analysis gives the com-

43

plex values for temperature, density and ﬂux at the inner boundary as output.

44

Chapter 6

Steady State Solution

6.1 Helium Boundary Conditions

6.1.1 Temperature

In order to ﬁnd the temperature below the photospere one has to balance the net
energy deposition in this layer. The outward heat ﬂux is set by the released energy
from nuclear burning below the layer and from gravitational settling. Assuming that

all helium is burned to carbon the outward ﬂux due to nuclear burning is given by

 

Pm,C = 162230 m = 0.61 MeV (11) , (6.1)

mu mu

where Q30 is the released energy per triple-oz process and m the local accretion rate

per unit area

With typical accretion rates Of 0.1 — 0.3 mm, the sum of the ﬂuxes gives

F ~ 1024 ergs cm"2 s_1. (6.3)

45

The temperature can be calculated with the Stefan-Boltzmann law

4 _ l
T — 0' (Fgrav + Fnuc) . (6.4)
6.1.2 Flux

As the ﬂux value strongly inﬂuences the thermal structure of the layers below the
carbon burning region, the actual value in previous studies is determined by setting
inner boundary ﬂux values and using shooting methods in order to reach the given
inner constrains. Several approaches for inner boundary values have been used:
Brown & Bildsten [12] set the inner ﬂux value at the bottom of the helium burning
to a range of0.1— 0.3 MeV/nuc. Narayan & Heyl [25] integrated the stellar equations
further into the stellar core, which is assumed to have a constant temperature in their
calaulations, and shoot for the ﬂux value that fulﬁlls the desired conditions.

As the energy ﬂux value at outer shells is set by the thermal properties like speciﬁc
heat, thermal conductivity, temperature and energy release of all layers below, the
exact determination of the ﬂux requires exact models for deeper layers and precise
knowledge of the central values. The energy ﬂux therefore varies signiﬁcantly with
the model and the central boundaries, especially because the sum of contributions
from deeper layers is large in comparison to the contribution from outer layers. The
released energy during carbon burning for example is ~ (0.3 — 0.4) MeV/nuc, deep
crustal heating in contrast releases 1.4 - 1.7 MeV/ nuc (Haensel & Zdunik [40]). In
order to avoid including all the different crust and core physics into the calculations,
the outer ﬂux value is used as a free parameter. All the quantities are therefore given
for a wide ﬂux range, which allows to switch between crust and core models and

boundaries without doing the calculations for each speciﬁc constellation again.

46

6.1.3 Density

From an assumed coupling parameter

 

_ (3)262 ~ 4 3_ (A)
— akT NI’ 37m — (6-5)

I"
C PNA

the following density temperature relation can be derived

pl/3 3(A) k
T = (4pm.) (2)2e2' (6'6)

With pure helium and a temperature of ~ 3 x. 107 K the outer density is

 

p0 z 105 gcm73. (6.7)

This surely is just a rough estimate, but as can be seen in appendix B the density and

column density in the carbon burning region is very insensitive to the outer density.

6.1.4 Sensitivity of the Boundary Conditions

Before varying the outer boundary conditions it is interesting to know how sensitive
the stellar structure is to these numbers. Some tests with the integrator clearly
show that the system adjusts inexact values of the outer temperature and the outer
density and column depth, respectively, very quickly and in any case long before the
carbon burning processes come into play. In stellar physics this is also known as
the radiative-zero boundary. A detailed plot of the column depth can be found in
appendix B. The energy ﬂux in contrast strongly determines the peak temperature
of the stellar structure and the column depth where carbon ignites. This sensitivity

supports the idea of using the ﬂux value as a free parameter.

47

6.2 Carbon Boundary Conditions

FI‘OIII the integration through the outer layers of the star to a depth where carbon
ignites one obtains a range of temperature, ﬂux and density values that can be used
for the outer carbon boundary conditions. Instead of starting the integration every
time from the very outer boundary, varying these outer carbon boundary conditions
and starting the integration at the maximum carbon abundance before carbon starts
to burn is less time consuming and allows variations of the composition in that layer.

The temperature range varies depending on the ﬂux between
3 S, T 8 S, 6, (6.8)

where T 8 is denotes the common abbreviation T8 = T/108K. The relevant integration

range for the column depth was set to

1010 gcm‘2 ,S y S, 1013 gem—2. (6.9)
For the outer ﬂux boundary, values in the range of

1keV/nuc f, EDI-:1E $1MeV/nuc (6.10)

seem to be reasonable. The outer flua: value F0 as well as the outer temperature T

in all the steady state and stability calculations is set at the outer boundary for the

010

column depth (yo = 1 gem—2).

48

Pure Structure

 

l I I I l l l I l

 

 

I

1—1
I I 1111"] ‘r I Illllll

D
(D
[\3
o
I
i
I
I
I
1

1 :5
1 1 111111

-2

log abundances by mass

  
   

I IIIII'I'II
B
0‘?
N)
a}.

 

 

he4 ......... . I K..-“—T.-.-....-._..i-.-.-...-..-

5 6 7 8 9 10 1 1
log column depth [g cmT2]

 

 

Figure 6.1: Abundances with pure accreted helium ignoring mixing with heavy ashes

6.3 Steady State Results

6.3. 1 Composition

Starting with pure helium as accreted material the dominant nuclear burning process
is the triple-a process which results in a very high carbon abundance of over 90% (see
ﬁgure 6.1). The 16O abundance is fairly small and vaies with the initial temperature
at the outer boundary. The lower the outer temperature and the resulting temperture
gradient the larger is the gap between the triple-a reaction and the reaction 12C +
4H6 —» 16O. In this case all the available helium is already burned before the the
oxygen producing reaction can ignite. In any case the fraction by mass of oxygen is
smaller than 10% and can be neglected.

However, turbulent mixing in deeper layers and other heavier ashes from previous
burning (bursts, rp-process) can change the composition and the maximum abun-
dance of carbon. In this analysis the carbon abundance is also varied in order to

investigate the inﬂuence of heavier elements. The ashes of bursts that are mixed

49

 

Coupling Parameter
70

 

 
   

 

 

 

C/Fe=3'/7 —' | l I
60 C/Fe=1/1---——-
C.) hC/Fe=7/3 ............
‘2 50 _
G)
E’
a, 40
8.
w 30
.e
E 20 a
8 .............
10 —— ------- ”7:: ..................................
0 ................... I l . l l
10 10.5 11 11.5 12 12.5 13

log column depth [g cm_2]

Figure 6.2: Coupling parameter for an initial temperature of 4 x 108K and outer

ﬂux of F091;? 2 2.5 MeV/nuc. The three curves correspond to different carbon abun-
dances in the initial carbon iron mixture.

with the freshly produced 12C are elements in the iron group with A z 60 (Woolsey
et al. [26]) and elements with A z 104 (Schatz et al. [27]). The major part of the

calculations is done wit a mixture of iron and carbon in varying fractions.

6.3.2 Properties of the Material

The phase of the burning material can be determined with the coupling parameter

_(Z)2e2 47rpNA ”3
Pc— kT ( 3(A) ) , (6.11)

which can be derived from equations 3.25 and 3.26. The coupling parameter for an

 

 

initial mixture of iron and carbon, an outer temperature of T = 4 x 108 K and ﬂux
value of F011"? 2 0.25 MeV/nuc can be seen in ﬁgure 6.2. The values for PC indicate
that the material is in the liquid phase over the entire burning range. The rise at

higher column depth is due to the freshly produced heavier elements with Z w 10

5O

Ignition Depth
4058+12 I I I IIIII' I I rIIIII' I I IIIIII

oooooooooooo
.....
................
'0-0
.-
...
C.-
'-
'0
'0

 

4e+12

 

3.5e+12

 

3e+12 T
2.5e+12 -
2e+12 —
1.5e+12 L-

ignition depth [g/cmz]

1e+12 [-
5e+11 -

 

 

 

 

0 1 1111111' 1 Llllllll 1 1 111111

0.001 0.01 0.1 1
outer ﬂux FOE-777',“ [MeV/nuc]

Figure 6.3: Ignition depth for different accretion rates

which increses the electrostatic energy.

6.4 Electron Capture

The ignition depth, the corresponding density and the released energy per nucleon
for carbon fusion as a function of outer ﬂux for an initial composition of each 50%
12C and 56Fe and an outer temperature of T3 = 4.0 can be seen in ﬁgures 6.3, 6.4
and 6.5. The released energy clearly shows a strong increase for ﬂux values above
~ 0.2 MeV/nuc which corresponds to an ignition density of roughly 1.7 x 109 gem—3.
A closer look at the composition as function of density shows, that at this density the
main product from carbon burning switches from 23Ne at higher ignition densities

to 20Ne at lower densities (see ﬁgures 6.6 and 6.7). This change is due to an electron

capture to 23Na.

51

ignition density [g/cm3]

released energy [MeV/nuc]

Ignition Density
2-5e+09 I I I IIITIIT I I 'ITIIII' I r I IIIII

 

 

 

 

 

2e+09 —

1.5e+09 —

1e+09 -

5e+08 --

 

 

0 l l lljllJL l l IILIILL l l llllll
0.001 0.01 0.1

outer ﬂux F032“ [MeV/nuc]

Figure 6.4: Ignition density for different accretion rates

Released Energy

 

 

0.4 I I I IIIII' I I I IIIIII I I I II'IIE
o

0.39 -
0.38
0.37 -
0.36 e
0.35 —
0.34 —
0.33 _
0.32
0.31 — ................................. ..

03 1 1 l_LJlllI 1 #1411111 1 1 111111
C

0.001 0.01 0.1
outer ﬂux F097;]lL [MeV/nuc]

 

0.1 mEdd
0. 2 mEdd -------
O 3 mEdd ............

 

 

 

 

 

 

Figure 6.5: Released energy for different accretion rates

52

 

Thermal Proﬁle for an Outer Flux F013? = 0.076 MeV/nuc
0.424 . . . . . . .

0.423
0.422
0.421

0.42
0.419
0.418

 

temperature [C K]

 

 

 

0.08 - ' 1 I r I ' ' I —
0.06 - e
0.04 a ..
0.02 - _
0
-0.02 r
-004 L —
-0.06 - 1
-0.08 —
1 ._ u

 

 

ﬂux me“ [MeV/nuc]

 

 

 

 

C13 ‘—
c __________
0.1 : 016 ............
: ne20 “.....-
I- ne22 -------:-"
F 11623 - ...........
- na23 ------'°'

f 6 .-..-..-..-.-
0.01 65 -

1e+12 1e+13
2l

abundance

 

 

 

 

column depth [g cm—

Figure 6.6: Thermal proﬁle for F071;? = 0.076 MeV/nuc, T3 = 4.0, m = 0.2 mEdd
and an initial composition of 50% carbon

53

Thermal Proﬁle for an Outer Flux For—”"34 = 1.0 MeV/nuc

 

 

 

 

 

 

 

     

 

 

 

 

 

 

 

0.65 ....l . .......l . .
5.? 0.6 - ._
Q
g 0.55 e _
g 05
o. ' ” ‘
8
8
0.45 — -
l l l l I I l l l I I I I I l l
1 I I I II I I I I I I I II I I —
F
:1
& 0.95 l— a
%
E 0.9 L _
as
it. 0.85 — —
5
c: 0.8 - —
l l l l I l I I I I I I l I I l
I I I II I I I I I I I II I I
f c12 —— .
~ 013 .......... electron capture -
016 ------------
g ne20 —~——~——
0.1 E- ne21 — ---------- _
'8 : ne22 ............ . . ... .._.'-_- """""" ..
E : ne23 ---------- -- '1 :
or! - na23 ------------- -
- mg24 - -
L mg26 — .
fe56 ----------
0.01 ‘ ‘ ‘ ‘1
1e+11 Ie+12

column depth [g cm‘z]

Figure 6.7: Thermal proﬁle for F037;,“ = 1.0 MeV/nuc, T3 = 4.0, in. = 0.2 mEdd and
an initial composition of 50% carbon

In a carbon rich environment the two main processes are

120+120 —> p+23Na (6.12)

120+12C —» 4He+20Ne (6.13)

As mentioned in section 4.3 the ratio of these two reactions does not vary strongly
in the temperature range between 108 — 109 K and has a value of (4He + 20Ne) / (p +
23Na) z 1.2, which results in almost equal parts of produced p,23Na,4He and 20Ne.

For lower ignition densities the most important following reaction is
23Na + p -+ 4He + 20Ne. (6.14)

20Ne is pretty stable at the given environment and has insigniﬁcant burning rates
which explaines the high abundance at the end of the carbon burning process (see

ﬁgure 6.7). At higher ignition densities the reaction rate for an electron capture on

23Na
23Na + e ——) 23Ne + I76 (6.15)

is much larger than the rate for 23Na+ p —> 4He + 2ONe which turnes a high fraction
of 23Na into 23Ne, which itself is also pretty stable at this temperture range. This
ﬁnally results in an 23Ne rich environment (see ﬁgure 6.6).

The threshold density at which the electron capture becomes important can be
calculated as follows: The mass difference between 23Ne and 23Na is 4.887 MeV.
Subtracting the captured electron’s restmass of 511 keV yields an energy difference
of up = 4.376 MeV (see ﬁgure 6.8). From the the electron chemical potential and
the relativistic parameter the density in an environment with (Z) / (A) = 0.5 ﬁnally

yields p = 1.675 x 109 gcm‘3 for the electron capture on 23Na, which corresponds to

55

23 Ne

 

““3 AB = 4.887 MeV

0.511 MeV

 

 

 

23Na

Figure 6.8: Energy level for electron capture on 23Na

a column depth of y = 2.75 x 1012 gcm‘2 and ﬁts exactly with the burning results
(see ﬁgure 6.7).

6.5 Validating the Calculations

In order to validate the computed structure, especially the column depth and the
total released energy during a superburst, the important quantity is the ignition
depth. In a superburst all the material above the ignited layer is burned to heavy
elements (Weinberg & Bildsten [15]), that do not contribute to further nuclear fusion.

Assuming that most of the accreted hydrogen is burned quickly to helium, the
layers above the carbon burning consist of a mixture of carbon and helium. De-
pending on the fracation of heavy ashes the amount of burnable fuel shrinks. For a
rough estimate equal parts of burnable fuel and inactive ashes are chosen, the fuel
consisting mostly of carbon. The released energy per nucleon is then given by the

difference in the binding energy:

1 3(56Fe) 3(120)
= $- {8.790 MeV/nuc — 7.680 MeV/nuc} (6.17)
= 0.555 MeV/nuc. (6.18)

56

The total mass above the ignition layer is simply given by the column depth. In
the nuclear calculations the layer above the ignited carbon has a mass of a few
1012gcm—2, which gives a number of 1036 nucleons per square centimeter surface
area. With a radius of 10 km the total released energy is ~ 1042 ergs, which matches
the values of former theoretical work and observations.

Estimating the recurrence time by assuming that at constant accretion rate the

012

neutron star has to accumulate the same amount of material (1 g cm”), this time

can be written as
tree = yignition . (6.19)
m

With an accretion rate of m/mEdd z 0.2 the recurrence time is in the order of

tree m 2Yr-

6.6 Influence Of Heavy Elements

As mentioned before mixing between different layers and heavy ashes from previous
explosions can enrich the carbon burning layer with A Z, 60 elements that do not
take part in nuclear fusion any more. However these elements inﬂuence the thermal
quantities like the speciﬁc heat and the thermal conductivity in the layer. Especially
the thermal conductivity, which is much lower for heavier elements causes the ignition
front to move farther out due to the more insulating character of the plasma. The
thermal conductivities for layers composed of 50% carbon and 50% heavier elements
(160, 56Fe, 104Ru) are plotted in ﬁgure 6.9. At a column depth of y ~ 1012 gem-2
the carbon ignites. The resulting inﬂuence of the heavier reaction products on the
conductivity is signiﬁcant for the carbon-oxygen composition. Mixtures of carbon
and iron or ruthenium do not show a remarkable change in the conductivity below

the carbon burning depth. FAIII thermal proﬁles for the same compositions are shown

57

Thermal Conductivity

 

      

 

le+19 : I I I IIIII' I I IiIIIIIl
5 (3-0 __
' C-Fe ..........
[ C-Ru

1e+18

I IIITIIT]

 
 

18+17 E"

l l 1 Llllll

 

thermal conductivity [ergs cm—1 K’l]

 

1€+16 I I I'll." ' ' "lllll l I 111141
1e+10 1e+11 1e+12 1e+13

_2]

 

column depth [g cm

Figure 6.9: Thermal conductivity for an initial abundance of each 50% of 12C+160,
12C+56Fe and 12C+104Ru

in ﬁgures 6.10, 6.11 and 6.12. All three proﬁles were computed with an outer tem-
perature of T3 = 4.0, an outer ﬂux of Fog?!“ = 0.15 MeV/nuc and an accretion rate
of m = 0.2 mm... Figures 6.13 and 6.14 show the ignition density and the released
energy for diﬂerent compositions. The temperature at the outer boundary was set

to T3 = 4.0, the accretion rate again to m = 0.2 mEdd.

58

Thermal Proﬁle for GO Plasma
0.417 I I I I I I I I

 

0.416
0.415
0.414
0.413

temperature [GK]

0.412

 

 

0.16
0.14
0.12

0.1
0.08
0.06
0.04
0.02

 

ﬂux F??? [MeV/nuc]

 

 

 

 

0.1

abundance

 

 

 

21

column depth [g cm—

Figure 6.10: Thermal proﬁle for an initial abundance of each 50% of 12C and 160,
an outer ﬂux 17091;:‘ = 0.15 MeV/nuc and an accretion rate of m = 0.2 mEdd

59

Thermal Proﬁle for C-Fe Plasma

 

 

 

 

 

 

 

 

 

 

  

 

 

 

 

 

0.446 T ' I ' ﬂ I u y 1
E ‘
52, 0.444
3 0.442 —
3
as
3 0.44 _
o.
E, 0.438 _
0.436 _
_ 0.15 —
U
:3
a
\
> 0.1 - _
Q)
E.
as 0.05 F _
EL.
§
:1: 0 \
1 :- * - 1 . . _:
l m ':
‘D r 012 —-
g ’ c13 ----------
"D 01 _ 016 ............
g ' 5 ne20 -——--—-~ :
'8 : ne21 -..“.-- :
- ne22 - -----------
' ne23 ---------. .
"$32 l
0.01 e . .. . . . .
1e+12 1e+13

column depth [g cm‘z]

Figure 6.11: Thermal proﬁle for an initial abundance of each 50% of 120 and 56Fe,
an outer flux Fawn-{i = 0.15 MeV/nuc and an accretion rate of in = 0.2 ThEdd

60

Thermal Proﬁle for C-Ru Plasma

 

 

 

 

 

 

 

 

 

 
 
   
  

 

 

 

 

 

0.474 r I I I I I I I
._ 0.472
:5:
9.. 0.47
0)
5 0.468
E3
33 0.466
a.
E. 0.464
0.462
__ 0.15 _.
U
S
a
\
> 0.1 — _
Q)
E.
as 0.05 - _
in.
§
CI: 0 ‘
1 _ - v . - - - .5
t ~— 2
a) c ----------
a ' 016 ------------ a
33 01 __ ne20 —..-.._.... ,. :..._._.,,___...........-........;
g ' E neg; — ........... f _,..:..'-',--..-r-'--" ------------------------------------------ a
JD - ne ------------ _,..j
‘6 ne23 ---------- ’/ .. .
' na23 ,5 .......
r- mg26 .- :1 . .
0 01 1.11104 if" I 1," .- _,_...-.......-.....--.-...-.-......-.:..
' 1e+12 1e+13

column depth [g cm”2]

Figure 6.12: Thermal proﬁle for an initial abundance of each 50% of 12C and 104Ru,
an outer ﬂux For—"77% = 0.15 MeV/nuc and an accretion rate of m = 0.2 ThEdd

61

2e+09
5""
E
if?
é
"o 1e+09
c:
.9.
'5)
"‘ 5e+08

 

Ignition Density

 

 

I I I L I

 

 

0.1

outer ﬂux Fawn-ii [MeV/nuc]

Figure 6.13: Ignition density for different heavy elements with an initial temperature
of T3 = 4.0 and an accretion rate of in = 0.2 ThEdd

0.42

0.4

0.38

0.36

Q [MeV/nuc]

0.34

0.32

0.3

Released Energy

 

 

 

 

 

 

I T I Ij I I I I I I I ..I."
C- O —
- C- Fe ------ -
C- Ru ............
l.
I I I I l I I l I I I I I
0.1

outer ﬂux For—"n3“ [MeV/nuc]

Figure 6.14: Released energy for different heavy elements with an initial temperature
of T8 = 4.0 and an accretion rate of in = 0.2 T'nEdd

62

6.7 Inﬂuence of Temperature

The ignition results for a ﬁxed accretion rate of Th = 0.2mEdd and a composition
of 50% of each 12C and 56Fe show a strong ignition depth and density dependence
on the temperature for low ﬂux values (see ﬁgures 6.15 and 6.16). For higher ﬂux
values the curves approach each other and drop below 1012 gcm’2 and 109 gcm‘3,
respectively. The released energy (ﬁgure 6.17) changes strongly for the three lower
temperatures as a function of outgoing ﬂux due to crossing the threshold value for
electron capture (see section 6.4). For a temperature of T 8 = 5.0 the ignition depth
is low enough for all ﬂux values in the chosen range to burn all the material before
the critical density for electron capture is reached. Figure 6.18 shows the density as
a function of ignition density. In the overlapping ranges for different temperatures
the released energy clearly shows a very small direct dependence of the Q value on
the temperature.

Note that for all calculations the givien temperatures are the outer boundary

values at y = 1010gcm’2. The temperature at the ignition points is remarkably

higher for high flux values.

63

Ignition Depth

 

 

 

 

 

 

 

 

 

 

7e+12= I IrITIIl rrIIII] I I I IIIII
6e+12 -
0qu 5e+12 —
5i
39.
a 4e+12
*5.
Q)
'2 3e+12 —
.9
a 2e+12 —
le+12 if
0 l I IIIIILI I I IJIIJII #LIIIIII
0.001 0.01 0.1 1

outer flux Foam“ [MeV/nuc]

Figure 6.15: Ignition depth for different temperatures at an accretion rate of m =
0.2 T'nEdd and a composition of 50% of each 120 and 56Fe.

Ignition Density
3.5e+09 .

 

I ITIIII' I rIIIIIT' I I IIIIII

 

T8=3.5 —
3e+09 T8=4.0 ..... -
T8=4.5 ............
2.5e+09 L “3:50

  

 

- -... - .-

2€+09 - “I“...N‘N.

 
 

    
  

1.5e+09 —

1e+09 -

ignition density [g/cm3]

5e+08 —

 

 

0 I I IIIIIII I I LIIIIII

0.001 0.01 0.1 1
outer flux F091;? [MeV/nuc]

Figure 6.16: Ignition density for different temperatures at an accretion rate of m =
0.2 ThEdd and a composition of 50% of each 12C and 56Fe.

64

Released Energy

 

    
  

0:4 I I TIIIII] I I IIIIII' I I IIIII
out.

 

 

 

 

 

 

 

”....

0.38 — °

3'

:‘o

i 0.36 —

>

G)

E- 0.34

is

8

C.‘ 0.32 - ‘_'_”.__‘,,...-"'

G)

"O

O.)

3 0.3 —

ass—-
0’28 " T8;4.5 ............ "

T8=5.0 ------
0.26=——' l lllllll I l 111111L 4 . Ill—111
0.001 0.01 0.1 1

outer flux Foam“ [MeV/nuc]

Figure 6.17: Released energy for different temperatures at an accretion rate of m =
0.2 mEdd and a composition of 50% of each 12C and 56Fe.

Released Energy
0-42 l T l I l l

 

0.4 -
0.38 r-
0.36 —
0.34 h
0.32 L-

0.3 -

released energy [MeV/nuc]

0.28 -

 

 

 

0.26 l I l L J L

0 0.5 1 1.5 2 2.5 3 3.5
3]

 

ignition density [109g/cm

Figure 6.18: Released energy for different temperatures at an accretion rate of fn =
0.2 mEdd and a composition of 50% of each 120 and 56Fe.

65

Chapter 7

Linear Stability Analysis

7.1 Linear Perturbation Equations

In order to be able to do a linear stability analysis a steady state solution is needed.
Although the star’s outer layers are not perfectly in steady state, the time scales at
which the structure of these layers change are much larger than the burning time
scale, so that the assumption of a time independent thermal structure holds for this
model. For an inward moving parcel of burning plasma the Lagrangian Perturbation

in its linearized form can be applied

(My, t) = <1>0(y) + (My, 0 (71)

where (1)0 is the solution obtained by steady state integration and <I>' is an inﬁnitesimal

perturbation, which satisﬁes the condition

I

(I)
3 << 1 (7-2)

This inﬁnitesimal character allows to neglect terms of higher order. The important

property of the Lagrangian method is that a change of the quantity (I) is computed in

66

a mass shell not at a given spatial position (Eulerian perturbation). The perturbation

is arbitrary and can be written generally as a complex variation
(My, t) = (my) exp(wt) (7-3)

with the complex quantity 0) which plays the role of a growth factor of the pertur-

bation. Normalizing the amplitude to the steady state quantitiy (DO yields

<I>(y, t) = <1>o(y) [1 + (My) exp(wt)] (7-4)

with a perturbation function ¢(y).

7.2 Stability and Overstability

The solution is ﬁnally a complex quantity with the complex growth frequency w =
w,» + iwi. The stability of the system can be obtained from the value of this complex
0.1. A closer look at the exponential function shows that for positive values for wr
the perturbation ¢(y) increases with time whereas negative values cause the external
impact to die out with time. The imaginary part of the frequency does not give
any information about the amplitude of the perturbation. The stability criterion can

therefore be summarize as follows: A found solution of the quantity is

stable if wr < 0 (7.5)
indifferent if wr = 0 (7.6)
unstable if our > 0. (7.7)

The real part of 0) does not only give information about the stability, but also about

the time scale of the rise and damping of the perturbation. This time is called

67

e-folding time, is connetcted to w as

_ 27r
lwrl

 

(7.8)

7"e

and denotes the time for the perturbation to increase or decrease about a factor of
e.
Although the imaginary part 01,- does not determine the stability, it indicates

oscillations in the perturbation during the burning
<I>’(y,t) = 6<I>(y) exp(wrt) exp(iwz-t). (7.9)

In case of a positive real part w and an imaginary part unequal zero, the unstable

mode is called overstable. The oscillation period II is given by

n (7.10)

— Iwz‘l'

7.3 Application to the Stellar Model

The important quantities that determine the burning processes are the temperature,
the density and the outcoming energy flux. The stability analysis therefore concen-

trates on perturbations to these quantities, namely

T(.% t) = Tss(y) (1 + 9(1) exp(wt)) , (7-11)
p(y, t) = pss(y) (1 + ((0 exp(wt)), (7-12)
F (y, t) = Fss(y) (1 + f (t) exp(wt)), (713)

where the subscript ss denotes the steady state proﬁle and the time t in the perturbing

functions is related to the column depth like dy/dt = mace. The new stellar equations

68

are now in general complex quantities. The solution of the system of differential
equations therefore yields for every complex growth frequency w = cur + 202,- a set of

three complex functions.

7 .4 Perturbation Function

The perturbation functions (0,(, f) are dimensionless functions with an amplitude
much smaller than one. The duration of the perturbation (tp) should be a non
negligible fraction of the total burning time. In these calculations a fraction of 25%
of the total burning time was chosen and the following shape of the function was

used (see ﬁgure 7.1)

_t_ 1
1 for tp < 3

9(t),€(t).f(t) = A0 ices [7r (g; -710] +§ for ,1, g ,t; 3.3. (7.14)
0 for it; > 3,

where A0 denotes the amplitude of the perturbation and was set to 10‘3. The exact
shape of the function was chosen arbitrarily. Some tests with different functions
show, that smooth functions lead to less numerical noise in the solution for w in
comparison to a simple step function. The shape of the function as well as the exact
duration of the perturbation do affect the values for 6min, (see section 7.5), but do
not affect the value of 02, except for very short (tp S, 0.05 t) or very long (tp 2, 0.8 t)

durations of the perturbation.

7.5 Calculation Procedure

For the complex stability analysis a reasonable range of the complex growth frequency

02 has to be scanned. For each value of w the thermal structure equations have to be

69

Perturbation Function

 

amplitude A(t) /A0

 

 

 

 

l l l I l l l l l

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
time ltC,totl

 

Figure 7.1: Perturbation function

solved. Among the ﬁnal set of solutions 3(0)) the ones have to selected, that match

the inner boundary conditions of the steady state solution. These conditions are:
0 The solutions for p, T and F must be real.

0 The real parts of the complex structure equations have to match the steady

state solution.

As ad has to be discretized, the solutions generally can only ﬁt within a threshold

value ere for the real and 6im for the imaginary part. The numerical conditions can

70

be written as follows

 

 

 

Im(p)
_ S e- ,
lpl ”n
Im(T)
_ < .
ITI — ElI‘n)
Im(F)
__ < . .1
lFl _ 61m: (7 5)
lRe(p)-pss| S Ere,
Ipl
lRe(T)—Tssl
<
lTl _ Ere)
lRe(F) " Fssl
<

where the subscript ss denotes the steady state result. The values for are and Eim
vary with the region in the complex plane for w. For large positive wr min(ere,im) 5
10’2, at the transition from unstable to stable burning min(ere,im) S, 10‘4. The
large values for very unstable burning are caused by the large grid size, that had
to be chosen in order to reduce the computing time (see section 7.6). Interpolating
between the grid points in the complex plane reduces the threshold values min(ere’im)
signiﬁcantly, but does not change the result. Figure 7.2 shows two examples for a
solution for the complex 0). Dark areas indicate the regions where the solution
conditions (equations 7.15) are fulﬁlled. For the growth rate of the perturbation it is
instructive to make comparisons to the total carbon burning time. A perturbation
that grows on much larger time scales than the burning time can not affect the the
burning process signiﬁcantly. Growth times much smaller than the burning time
are considered to cause instabilities. The total carbon burning time therefore sets a
minimum for the modulus of lwl. As the initial perturbation amplitude is small (10—3)

it is reasonable to start with a modulus of |w| that can increase the perturbation to

at least a few percent. In the following calculations the minimum value for was set

71

Stable c0

 

 

 

 

 

 

 

 

 

 

10 I I I
5 — —.

TU)

‘9

e o - -

3

E
-5 )— ...
_10 l J I

-10 -5 0 5 10
Re(u)[10'6 s"])
Unstableco
20 I I I I I I I
15 l- -
10 -— -

Tm 5 ’- ﬂ '-

‘9

2 o - -

'5'

g -5 — 9 ._
-10 — -I
-15 - —I
_20 l J I J l I l

-20 -15 -10 -5 0 5 10 15 20

Re(co[10’5 s"])

Figure 7.2: Example for a stable and an unstable solution for 02. Light gray areas
correspond to values of w, in which equations 7.15 are not fulﬁlled (relative differences
> emim). Dark areas indicate solutions for all equations 7.15.

72

to

3
lwminl = 2": (7°16)
b

where tb is the total burning time and the arbitrary factor of 3 corresponds to an

increase of 2%,
f(t = 0) exp(3) z 10‘3 - 20 = 0.02. (7.17)

The upper value was set to

30
lwmaxl = 10leinl = E» (7.18)

which corresponds to a maximum possible growth factor of z 1010.
In ﬁgure 7.3 the growth time in units of the total carbon burning time is plotted.
The transition between unstable and stable burning is limited to a very narrow flux

range, which gives a distinct number for the threshold value.

7.6 Discretization and Computing Time

The most intuitive way of dicretizing the 0) grid in the complex plane is transforming

w to polar coordinates
a) = wr + iwz- = |w|(cos¢ +i sin <15) (7.19)

and discretize the radius |w| and the angle 03 in the complex plane. The stepsize was

determined after some tests after having found the size of the solution pattern. The

73

 

 

Perturbation Growth Time

 

003 I I l IIIjI' I I rIIIII'! I I I ITIII
1

0.25 b
0.2 — -

0.15 - unstable stable

I

 

Te /tburn

0.1 — _

0.05 L _

 

 

0 I I I IIIJJI I I I IIIIII I I I IIJII

0.001 0.01 0.1 1
outer ﬁux F097;?- [MeV/nuc]

 

Figure 7.3: Growth time of the perturbation in units of the total carbon burning time
for T = 4 x 108 K, a composition of equal parts of carbon and iron and an accretion
rate of in = 0.2mEdd.

radius was divided into 30 steps

71 = Iwminl + kﬂ‘géP—l, k = 0.29, (7.20)

for the angle (15 90 steps were choosen

2 _
(1),, = kgg, k = 0.89. (7.21)

The grid in the complex plane therefore consists of 2700 parcels. The integration time
for one wk was roughly one second. As the Bulirsch-Stoer integrator in the reaction
network uses an adaptive stepsize and the time of the fourth order Runge—Kutta
integrator can be neglected, this integration time was pretty much independent of

the given stepsize. The errors in the reaction network concerning the abundances 61/

74

with
EY = 1 - Z Xi (7.22)
i

were less than 10‘15 for all test calulations. The error for the thermal structure
equations in the Runge-Kutta integrator were not measurable for stepsizes above
N z 1000. The total computational time in order to determine the stability for
one set of parameters was therefore 45-50 minutes. In order to ﬁnd the threshold
ﬂux between unstable and stable burning for a ﬁxed accretion rate, temperature and
initial composition a set of z 30 stability calculations was practicable, which took

roughly one day.

7.7 Results of the Linear Stability Analysis

7.7 .1 Stability Flux

As mentioned in section 6.1.2 the ﬂux value is used as a free parameter. For a ﬁxed
set of outer boundary values the flux, at which the fusion switches from an unstable
to a stable burning, was determined. This ﬂux value is called stability flux in this
thesis. For all ﬂuxes lower than this threshold the burning is unstable, for all higher
values it is stable. The big advantage of determining the stability ﬂux is, that for any
given ﬂux, determined by various models and boundary values, one can directly ﬁnd
out whether the carbon burning is stable or not, without doing the stability analysis
again for every set of parameters. In addition to that, the computed ﬂux value from
other calculations does not only show whether stable or not, but also how close to

the stability this flux is.

75

Ignition Density

 

2e+09

1.5e+09

1e+09

ignition density [g cm—3]

5e+08

 

0..
I I I I I I ’L I I I I I I

 

 

 

 

 

 

0.1 1
outer ﬂux Fog-5n“ [MeV/nuc]

Figure 7.4: Influence of heavy elements on the stability. The vertical lines indicate

the threshold ﬂux value for stable burning. Higher ﬂux values are stable, lower are
unstable.

7 .7.2 Inﬂuence of heavy elements

In ﬁgure 7.4 the ignition densities can be seen for three initial compositions, each
consisting to equal parts of 12C and the heavy element (160, 56Fe, 104Ru). The
stability threshold flux is indicated by the vertical line, giving the following values

for a temperature of T3 = 4.0 and an accretion rate of Th = 0.2 ThEdd

120 — 16o : E; = 0.12 MeV/nuc ya = 3.6 x 1012gcm'2 (7.23)
12C — 56Fe: PC = 0.13 MeV/nuc yo = 2.65 x 1012 gem—2

12C — 104Ru : C = 0.16 MeV/nuc yo = 1.75 x 1012 gcm‘2

with strongly decreasing ignition column depth for heavier elements.

76

Temperature Dependence
0'35 I I I

 

0.3

0.25

0.2

0.15

outer ﬂux Fomnvii [MeV/nuc]

0.1

 

 

 

0.05 I l
0.34 0.36 0.38 0.4 0.42 0.44 0.46 0.48 0.5 0.52

temperature [GK]

 

Figure 7.5: Inﬂuence of the temperature on the stability. For ﬂux values higher than
the shown curves, the burning is stable, for lower values it is unstable

7 .7 .3 Inﬂuence of the Temperature

The calculations for different temperatures are all set up with an initial composition
of 12C and 56Fe at different ratios. The results (ﬁgure 7.5) clearly show a strong
temperature dependence for all compositions. Even more sensitive to the stability
flux is the composition of the material. Higher carbon abundances require much

higher flux values for a stable burning.

7.7 .4 Inﬂuence of the Accretion Rate

The dependence of the accretion rate on the stability ﬂux (ﬁgure 7.6) is fairly small.

In the range of 10 — 30% the Eddington limit the'threshold values barely increase.

77

Dependence on Accretion Rate
0.3

 

 

 

 

 

 

 

 

 

C l l l l l I I
0.275 — (328 — _
”g 0.25 — 838 ____'__'_ _
; 0.225 _
(D
a 02 _ ............................................ .1
515 0.175 : ..
In.
0.15 '— —I
s
“5 0.125 - _
s
g 0.1 — _
0.075 1.
0.05 I 4 J I l I I

 

 

 

0.1 0.125 0.15 0.175 0.2 0.225 0.25 0.275 0.3
accretion rate lmEddl

Figure 7.6: Inﬂuence of the accretion rate on the stability. For ﬂux values higher
than the shown curves, the burning is stable, for lower values it is unstable

78

7.8 Burning Regimes

7 .8.1 Direct and Delayed Bursts

For an unstable burning regime the question arises, at what point in the carbon
burning process the perturbation starts growing signiﬁcantly and causes the ignition
of the burst. In order to determine the ignition region, the [burning process the
is devided into two parts. The ﬁrst part covers the burning regime from 90% to
50% of the initial carbon abundance (hereafter region I), the second one refers to
the burning region from 50% to 10% (hereafter region 11). If an instability already
occurs in region I the burst is called an direct burst, if region I is stable but region II
is unstable the burst is called a delayed burst. In order to determine the ignition
time, a complex stability analysis is done separartetly for each region. It turns
out, that burning region 11 shows an unstable character for all ﬂux values below
the stability threshold ﬂux, whereas burning region I is only unstable for low ﬂux
values, corresponding to deeper ignition depths, and switches to a stable mode for
intermediate ﬂuxes. That means that that for increasing ﬂuxes the unstable ignition
front moves lower carbon abundances until all the carbon is burned before reaching
the unstable depth. Figure 7.7 shows the direct and delayed bursts for a ﬁxed initial
composition of equal parts of carbon and iron and an accretion rate of m = 0.2 T'nEdd.
In ﬁgure 7.9 the fraction of carbon in the initial composition is varied. The switching
line from direct to delayed bursts roughly lies parallel to the threshold stability line.
Figure 7.8 shows the dependence on the accretion rate for the same composition and

a temperature of T8 = 4.0.

7.8.2 Overstable Modes

As mentioned in section 7.2 the ﬂame is overstable if the imaginary part of the growth

frequency is unequal zero. This leads to oscillations in the increasing perturbation.

79

Direct and delayed bursts

 

 

 

 

 

 

 

 

L I d 9; dye d I a I I I I I
—~ ‘ irect I 7
a 0'32 .- threshold ﬂux 1
\ .- I
T - stable ~
2.. 0.16 L E, a -
l' E] 2 £1 a ‘
€15 -—r"‘"':" E: g a g 3 1
LE ' g 3 3 3 “’1' 3 E' 3 l
g 0.08 - g g g estableg E El ‘-
u: n a g 3 . as :
H I g 3 E1 " ' : I
8 I I z : : I l
5 0-04 - 3 S : : : : : -
§ I I I K I I
"‘ l " l ‘I‘ I "' I "‘ I " I I
0.34 0.36 0.38 0.4 0.42 0.44 0.46 0.48 0.5 0.52

temperature [GK]

Figure 7.7: Direct and delayed bursts for an accretion rate of in = 02de and a

composition of 12C : 56Fe = 1 : 1

Direct and delayed bursts

 

 

 

 

 

 

I I I I I
0'32 _ delayed a 1
? _ direct a: -
a threshold ﬂux .
\ _ .
a; 0.16 _ stable d
a —5——'ﬂ—_—E———
El"? r E E E un- 3 3 1
LS 0.08 — g g g stableE s
.2 a E a 1 a .
.32 5’ 3 : : :
8 0.04 [- a S : :1 : _
'I‘ ll 1' 1' ’l'
0.05 0.1 0.15 0.2 0.25 0.3 0.35

accretion rate [mEddl

Figure 7.8: Direct and delayed bursts for an a temperature of T8 = 4.0 and a com-

position of 12C : 56Fe = 1 : 1

80

Direct and delayed bursts

 

  
 
      
     
 

 

    
 
   

       

 

 

 

: I I I I I I I I

_ 0384 f“ delayed a

g . direct a!
& t threshold ﬂux —— E, 3
> 0.192 L- a g g j
g . stable 3 g g g .
‘_‘ : a a 3 S S :
as 0.096 r a g g a 3 i: 1
LE 1 D 3 3 3 5 : z: 1
g I E g a : : l I :
:1: (1048 f g E g : I : : T
5 ' E : : i 2 i 2 1
g F = ’- " ' : : I; :
O 0.024 a; : unstable ,7, ,. . . Ii
1 I T l T I I I 1
0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

carbon fraction

Figure 7.9: Direct and delayed bursts for an a temperature of T8 = 4.0 and an
accretion rate of m = 0.2 ThEdd

As at a given temperature, density, accretion rate and composition the ﬂux range,
within which these oscillations appear and disappear is very small, it is not possible
to make a general statement about the occurence of these oscillations based on the
relatively large ﬂux steps, that were used in the stability calculations. A wider ﬂux
range, however, could be determined for high temperatures above T 3 > 4.5 and high
carbon abundances (70%, 90%), where they were easy to detect. Focussing on this
boundary values, the evolution of the oscialltions follows an interesting pattern. For
very low ﬂux values the unstable ﬂame does not show overstable modes. At the lower
threshold ﬂux value the solution for w undergoes a bifurcation and branches out to
two overstable modes (see ﬁgure 7.10). The imaginary part increases with increasing
ﬂux up to a maximum of 60 cycles during the total C burning time. Increasing the
ﬂux a real unstable mode appears and replaces the overstable modes, see ﬁgure 7.11.
In ﬁgure 7.12 the frequency in units of total carbon burning time is shown. As

expected the oscillations appear in pairs with positive and negative imaginary part

81

 

 

Im(m)
I
.-
I
Im(co)
I
«In
I

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

I I I I I I I I I I L I I l I I I L
Fle(03) Rem)
r I I I I I I I I I I I I I I I I I
._ 3 1 ._ 4 _
_ 1 .. _
[_ _ .. _
s: 1 I s ' 1" ‘
E , E ‘ . ‘
_ r I _ _ _ I _
_ _ _ J
I l l I I l I l l I l l l I l l l l
Re(oa) Re(co)
I I I I I I I I I I I I I I I I I I
E _ é _ _
E r E .
_ _ - I _
I I I I I I I I I I I I I I I I I I
Re(co) Re(m)

Figure 7.10: Appearance pattern of the oscillation for Foam“ m 0.107 MeV/nuc, T8 =
5.0, 70% carbon. Light gray areas correspond to values of w, in which equations 7.15
are not fulﬁlled (relative differences > swim). Dark areas indicate solutions for all
equations 7.15.

82

 

 

 

 

 

 

 

 

 

 

 

 

I I I I I I I I I l
.—1 — l— 2 1
" 0 _ 7 3

'5 ’55

E ' I E t ‘
_ I — — 3 —
_ 1 _ _

I I I I I I I I I I
Re(m) Re(0))
I I I I I I I I I I
_ 3 _ _ 4 ._
_ ‘ _ .. a; _

’5 , '87

E ' ! * E - i ‘
I ' - - A

I I I I I I I I I I
Re(co) R900)

 

 

Im(co)
I
I
Im(m)
I
I

 

 

 

 

 

 

I I
Re(co) Re(a))
Figure 7.11: Disappearance pattern of the oscillation for F097;? m 0.153 MeV/nuc,
T3 = 5.0, 70% carbon. Light gray areas correspond to values of w, in which equa-
tions 7.15 are not fulﬁlled (relative differences > €re,im)- Dark areas indicate solutions
for all equations 7.15.

83

Osciallation frequency

 

60 I I I I I

A
O

[\D
O

 

['0
o

frequency [total burning time]
is.
o o

 

 

 

 

I I I
0.08 0.09 0.1 0.11 0.12 0.13 0.14 0.15 0.16 0.17

outer ﬂux For—rﬁn‘ [MeV/nuc]

 

Figure 7.12: Appearance and disappearance of the oscillation for T8 = 5.0 and 70%
carbon.

w and the same absolute value: to: = —wz- .

A closer look at the region where these oscillations are generated shows that
burning region I does not have overstable modes at all. The patterns in region II
exactly match the pattern of the total burning process which leads to the conclusion,
that the oscillations are driven by reactions in the region, where most of the carbon
is already burned and a signiﬁcant amount of the nuclear energy is already released.
In all the tested cases, the ﬂame switches from a direct to a delayed burst roughly

at the ﬂux value at which the oscillation appears. This transition is indicated by the

vertical line in ﬁgure 7.12.

84

Chapter 8

Application

8.1 Neutron stars in quiescence

Aside from direct burst oberservation and analysis of the ﬂare itself, the further time
evolution of the neutron star during cooling gives information about the burning
circumstances in the carbon burning layer. Brown & Cumming [19] found that
the temperature proﬁle of the steady state burning during accretion and the proﬁle
at the end of the outburst differ by less than 4% over the entire column depth
range. Mapping the crustal heating with the observed cooling lightcurves allows to
determine the temperature and ﬂux right behind the ignition front of the steady state
burning. Running a stability analysis with the given temperature and accretion rate
for different compositions of the ocean before the unstable ignition gives a relation
between the composition and the threshold ﬂux values for a superburst. The mapped
flux value in the ignition front sets a limit on the composition for unstable carbon

burning.

85

Temperature proﬁle

 

ET
9.. 4
(D
H
:3
4—7
G —
I...
8.
E
3 -I

 

 

 

 

9 10 11 12 13 14 15 16 17 18 19
column depth [g/cmz]

Figure 8.1: Thermal proﬁle for KS 1731-260

8.2 The neutron star in KS 1731-260

The neutron star has been observed to cool following a 2.5 yr-long outburst (Wijnands
et al. [20]). The star has a mass of 1.6 MG and a radius of 11.6 km, which leads to
a redshift factor of 1 + z = 1.32 and a surface gravity of g = 2.3 x 1014 cm 3'2.
The outburst accretion rate is observed to be M m 1017 gs‘l, yielding an accretion
rate per unit area of in z 0.095 T'nEdd. Fitting a cooling curve to the observed data
the temperature behind the ignition front was determined to be Tb = 3.8 x 108K
(Brown & Cumming [19]), corresponding to a peak temperature at the ignition of
Tp = 6 x 108 K (see ﬁgure 8.1). The ﬂux value at the outer boudary was determined
to be 0.22 MeV/nuc. The result of the stability analysis can be seen in ﬁgure 8.2.
The difference between T = 5.5 x 108 K and T = 6 x 108 K can be neglected. Possible
compositions are indicated by the thick line, which yields a fraction of carbon of over

z 65%.

86

Stability Flux

 

 

 

outer ﬂux 170%“ [MeV/nuc]

 

 

 

T8=5.5 El
0-3 - T8=6.0 x '-
observed flux value

0.25 — -
0.2 _ -
0.15 - -
0.1 - E/ +
0.05 b -
0 _.

-0.05 l ' ' '
0 0.2 0.4 0.6 0.8 1

carbon fraction

Figure 8.2: Stability flux as a function of composition for KS 1731-260

87

Chapter 9

Conclusion and Outlook

The work presented in this thesis is a broad theoretical parameter study of the carbon
burning processes in the ocean of accreting neutron stars. The main interest focused
on the conditions where the nuclear burning switches from a steady state burning
to an unstable ignition, leading to a thermonuclear superburst. The study was per-
formed with numerical simulations under the assumption of a sperically symmetric
star in one radial dimension. The burning processes were calculated by a full reaction
network, including the 45 most relevant isotopes, connected by over 500 weak and
strong reactions. The stability analysis was based on a complex number perturbation
evolution, giving precise stability threshold values for temperature, ﬂux and accretion
rate. For a more precise determination of the point in time, at which the unstable
ﬂame ignites, the burning process was not only examined as a whole, but also devided
into two sub—regions, associated with direct and delayed bursts. In addition to that
oscillations of the perturbation were examined for a set of parameters.

In future work the effects of rotation have to be included in the model as well
as consequences of turbulent mixing between the layers. Especially if the mixing
timescale is in the order or even signiﬁcantly shorter than the burning timescale, the

resulting heat and material transport will certainly affect the burning and therefore

88

the stability conditions. Strong mixing can also lead to locally very high accretion
rates which are not covered in the present work. Furthermore the observed features
during the rise of the ﬂare such as millisecond oscillations are still unconnected to
the nuclear burning before the unstable ignition. Concerning the ignition mechanism
the inﬂuence of shocks and a possible connection to the often observed preceding
helium ﬂash before the superburst are not covered in this thesis. With a systematic
inclusion of density ﬂuctuations, the sensitivity of the ﬂame to this perturbation has
to be examined. In addition to that precise conditions for the appearance and disap-
pearance of burst oscillations are needed as well as a detailed oscillation structure for
a larger range of parameters. Finally a transformation of the stability problem with
a full reaction network and a complex perturbation analysis to a three dimensional

model would certainly approach the model further to the real star.

89

Appendix A

Thermodynamic Equations

Here a detailed derivation of the thermodynamic equations for the temperature and
the ﬂux can be found. Starting with the thermodynamic potential

dE = TdS + pdN — PdV (A.1)

the time derivative of the energy yields

dt dt dt 7t" (A2)

Not changing the number of particles nor the volume, the last two terms vanish.
Transforming to quantities per unit mass gives

E=T§=TdsdT

dt dt as" (A3)

The change in energy is given by the total rate de/dt = e, which can be separated
into heating and cooling éheat?£C001' Heating is provided by nuclear reactions, the
cooling is given by thermal diffusion and neutrino emission. The contribution from
thermal diffusion can be derived from the diffusion equation:

8 K
——T = ——v2:r (A.4)
at pCp
with the total thermal conductivity K. The diffusion equation can be split up into
the continuity equation

0 1 6 1

aT — —;—C—P-V ° F (1} CPET — —'p-V ' F (A.5)
and Fick’s law

F = ——KVT (A5)

90

With the deﬁnition of the speciﬁc heat

T BS (93
C — m7 “ Ta W)
equation A.3 can be written as
d3 ds dT (9 1 ,
_ = __ = _ =_— __ . A.
Tdt Tdet CpatT pV F+e ( 8)

where 5’ denotes the other sources and sinks of energy (nuclear heating and neutrino
cooling)

91

Appendix B
Column Depth

The column depth is deﬁned as the integral over the density and can be written in
spherical symmetry as follows

W) = /°°p(r>dr. (8.1)

In the integration inward the star the accumulated mass above the radius 1‘ therefore
varies with the initial density at the outer boundary and depends on the equation
of state of the material. However, in the neutron star’s ocean in the carbon burning
region, where the column depths reaches values on the order of 1010 — 1013g cm”2,
the column depth can be described by an exponential function which barely depend
on the initial conditions. The following diagram shows the column depth-density
relation for different initial conditions.

92

Density
10 I I I I I

 

 

 

 

 

 

 

logpo = 4 —
10g p0 = 6 ............

 

log denSIty [g / cm3]
q
I

03/)
0.744 . log(y) ".--.--

 

 

 

4 i i i j i
2]

 

log column depth [g/cm

Figure B.1: Relation between column depth and density for various initial conditions

93

Appendix C

Elements in the Reaction Network

In this thesis a full reaction network is used to determine the elemental abundances.
As the temperature and density range is limited and the main focus is on the carbon
burning processes only the relevant elements and their reactions are included in
the network in order to optimize the computing time. Starting with a very large
number of isotopes, this number was reduced by removing elements that change
the abundances by less than 10'4. As the released energy is calculated from the
abundances, the impact on the energy balance is also very small. The isotopes that
are used in the reaction network as well as their mass excess and the total binding
energy are tabulated in table C.1.

 

 

 

 

 

Eleﬁent Mass excess [keV] Binding energy |keV]=

 

n 8071

p 7289
4He 25930 28296
12c: 0 92162
13C 3125 97108
14C 3020 105285
12N 17338 74041
13N 5346 94105
14N 2863 104659
15N 101 115492
130 23111 75558
140 8006 98733
150 2855 111956
160 —4737 127619
170 -809 131763
180 —782 139807
17F 1952 128220
18F 873 137369

 

94

Mass excess [keV] Binding energy [keV]

 

Element

19F —1487 147801
17Ne 16490 112900
”Me 5307 132154
”Me 1751 143781
”Ne —7042 160645
21Ne —5732 167406
”Ne —-8024 177770
”Ne —5154 182971
”Na. 6845 145976
”Na. —2184 163076
”Na —5182 174145
”Na. —9529 186564
”Mg 17570 134470
”Mg 10912 149198
”Mg -397 168578
”Mg -—5473 181725
24Mg —13933 198257
25Mg —13193 202535
26Mg —16214 216681
”AI 18180 149220
23A1 6767 168703
24AI —55 183596
25A1 —8916 200528
26AI —12210 211894
27AI —17197 224952
56Fe —60601 492254
104Ru 893085

—§8091

 

 

 

 

Table C.1: List of elements

95

L

Bibliography

[1] Grindlay, J. E. et al. 1976, ApJ, 205, L127
[2] Belian, R. D., Conner, J. P. & Evans, W. D. 1976, ApJ, 206, L135

[3] Cornelisse, R., Heise, J ., Kuulkers, E., Verbunt, F. & in’t Zant, J. J. M. 2000,
A&A, 357, L21

[4] Strohmayer, T. E. 2000, AAS/ High Energy Astrophysics Division, 32
[5] Strohmayer, T. E. & Brown, E. F. 2002, ApJ, 566, 1045

[6] Kuulkers, E. et al. 2002, A&A, 382, 503

[7] Wijnands, R. 2001, ApJ, 554, L59

[8] Strohmayer, T. E. & Markwardt, C. B. 2002, ApJ, 577, 337

[9] Cornelisse, R., Kuulkers, E., int Zand, J. J. M., Verbunt, F., & Heise, J. 2002,
A&A, 382, 174

[10] Kuulkers, E. 2002, A&A, 383, L5

[11] Strohmayer, T. & Bildsten, L. 2003 arXiv:astro-ph/0301544v2
[12] Brown, E. F., Bildsten, L. 1998, ApJ, 496, 915

[13] Narayan, R. & Heyl, J. S. 2003, ApJ, 599, 419N

[14] Lewin, W. H. G., van Paradijs, J., & Taam, R. E. 1993,

[15] Weinberg, N. N. & Bildsten, L. 2007, ApJ, 670, 1291

[16] Heger, A., Fryer, C. L., Woosley, S. E., Langer, N., Hartmann, D. H. 2003,
ApJ, 591, 288

[17] Haensel, P., Potekhin, A. Y., Yakovlev, D. G. Neutron Stars 1, Equation of
State ans Structure, Springer 2006

[18] Povh, B., Rith, K., Scholz, C., & Zetsche, F., Particles and Nuclei, Springer
2006

96

[19] Brown, E. F., & Cumming, A. 2008 in prep.

[20] Wijnands, R., Nowak, M., Miller, J. M., Homan, J., Wachter, S., & Lewin,
w. H. G. 2003, ApJ, 594, 952

[21] Yakovlev, D., & Urpin, V. 1980, Sov. Astron., 24:3

[22] Yakovlev, D., & Urpin, V. 1980, Sov. Astron., 242303

[23] Potekhin, A. Y., Chabrier, G., Yakovlev, D. G. 1997, A&A, 323, 415
[24] Salpeter, E. E. 1954, AuJPh, 7, 373

[25] Narayan, R, & Heyl, J. S. 2003, ApJ, 599, 419

[26] Woosley, S. et al. 2004, ApJS, 151, 75

[27] Schatz, H. et al. 2001, Phys. Rev. Lett., 86, 3471

[28] Hansen, C. J., & van Horn, H. M. 1975, ApJ, 195, 735

[29] Taam, R. E., & Picklum, R., E. 1978, ApJ, 224, 210

[30] Narayan, R., & Cooper, R. L. 2007, ApJ, 665, 628

[31] Fujimoto, M. Y., Hanawa, T., & Miyaji, S. 1981, ApJ, 247, 267
[32] Shen, K. J., Bildsten, L. 2008, arXiv:0805.2160

[33] Haensel, 13., & Potekhin, A. Y. 2004, A&A, 428, 191 '

[34] Pandharipande, V. R.; Ravenhall, D. G. 1989, Nuclear Matter and Heavy
Ion Collisions. NATO Advanced Study Institutes Series. Volume B205, 1989,
p.103 ’

[35] Douchin, F., & Haensel, P. 2001, A&A, 380, 151

[36] Chabrier, G., & Potekhin, A. Y. 1998, Phys. Rev. E, 58, 4941
[37] Hansen, J. P. 1973, Phys. Rev. A, 8, 3097

[38] Lamb, D. Q, & Lamb, F. K. 1978, ApJ, 220, 291

[39] Schatz et al. 1999, ApJ, 524, 1014

[40] Haensel, P. & Zdunik, J. L. 1990, A&A, 227, 431H

97

[11111111111]

030

VI
U.]Il.
Tl
A“
I“
S“

3 129

 

[Hill]