J. Scheffel, P. Brunsell - Fusion Physics
J. Scheffel, P. Brunsell - Fusion Physics
J. Scheffel, P. Brunsell - Fusion Physics
- introduction to
the physics behind fusion energy
Jan Scheffel
Per Brunsell
I
II
Contents
III
2 THE FUSION PLASMA 43
2.1 Parameters of a fusion plasma . . . . . . . . . . . . . . . . . . . . . . . . . 43
2.2 Plasma models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
2.3 Particle dynamics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
2.3.1 Homogeneous magnetic field . . . . . . . . . . . . . . . . . . . . . . 46
2.3.2 Guiding center drift . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
2.3.3 Vertical drift in a toroidal magnetic field . . . . . . . . . . . . . . . 51
2.4 The extension of the particle fields . . . . . . . . . . . . . . . . . . . . . . 52
2.5 Kinetic plasma model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53
2.5.1 The Vlasov equation∗ . . . . . . . . . . . . . . . . . . . . . . . . . . 53
2.5.2 The Fokker-Planck equation∗ . . . . . . . . . . . . . . . . . . . . . 54
2.6 Coulomb collisions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
2.6.1 Mechanics of the Coulomb interaction . . . . . . . . . . . . . . . . . 55
2.6.2 Coulomb collision times . . . . . . . . . . . . . . . . . . . . . . . . 57
2.6.3 Plasma resistivity . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58
2.7 Fluid models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
2.7.1 Derivation of the fluid model equations∗ . . . . . . . . . . . . . . . 60
2.7.2 Magnetohydrodynamic model (MHD) . . . . . . . . . . . . . . . . . 61
2.7.3 Advanced fluid models∗ . . . . . . . . . . . . . . . . . . . . . . . . . 63
3 EQUILIBRIUM 64
3.1 Cylindrical equilibrium . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64
3.2 Beta value . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66
3.3 Toroidal equilibrium . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66
3.3.1 Equilibrium vertical field . . . . . . . . . . . . . . . . . . . . . . . . 67
3.4 Magnetic flux surfaces . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71
3.4.1 Flux functions∗ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72
3.4.2 Grad-Shafranov equation∗ . . . . . . . . . . . . . . . . . . . . . . . 73
IV
4.6.9 Further instabilities influencing confinement∗ . . . . . . . . . . . . . 102
5 TRANSPORT 105
5.1 Classical transport theory . . . . . . . . . . . . . . . . . . . . . . . . . . . 105
5.1.1 The random walk process . . . . . . . . . . . . . . . . . . . . . . . 106
5.1.2 Particle model for classical diffusion . . . . . . . . . . . . . . . . . . 107
5.1.3 Fluid model for classical diffusion . . . . . . . . . . . . . . . . . . . 108
5.1.4 Relation between confinement time and diffusion coefficient . . . . . 109
5.1.5 Classical heat conduction . . . . . . . . . . . . . . . . . . . . . . . . 110
5.1.6 Energy confinement time for ohmically heated plasma∗ . . . . . . . 111
5.2 Neoclassical transport theory . . . . . . . . . . . . . . . . . . . . . . . . . 112
5.2.1 Guiding-centre orbits in a tokamak . . . . . . . . . . . . . . . . . . 112
5.2.2 Particle model for Pfirsch-Schlüter diffusion . . . . . . . . . . . . . 114
5.2.3 Trapped particles in a tokamak . . . . . . . . . . . . . . . . . . . . 115
5.2.4 Banana orbits . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117
5.2.5 Banana orbit diffusion . . . . . . . . . . . . . . . . . . . . . . . . . 119
5.2.6 Summary of neoclassical transport . . . . . . . . . . . . . . . . . . 120
5.2.7 Fluid model for Pfirsch-Schlüter diffusion∗ . . . . . . . . . . . . . . 121
5.3 Experimentally measured transport . . . . . . . . . . . . . . . . . . . . . . 121
5.3.1 Confinement scaling . . . . . . . . . . . . . . . . . . . . . . . . . . 122
5.4 Models for anomalous transport . . . . . . . . . . . . . . . . . . . . . . . . 123
5.4.1 Particle flux from electrostatic fluctuations . . . . . . . . . . . . . . 124
5.4.2 Transport caused by drift waves . . . . . . . . . . . . . . . . . . . . 125
5.4.3 Stochastic magnetic field∗ . . . . . . . . . . . . . . . . . . . . . . . 127
5.4.4 Transport in a stochastic magnetic field∗ . . . . . . . . . . . . . . . 128
V
7 PLASMA DIAGNOSTICS 156
7.1 Rogowski coil . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 157
7.2 Interferometric measurement of electron density . . . . . . . . . . . . . . . 158
7.2.1 Principle of interferometry measurement . . . . . . . . . . . . . . . 159
7.2.2 A laser interferometer system . . . . . . . . . . . . . . . . . . . . . 160
7.3 Electron temperature measurement by Thomson scattering . . . . . . . . . 161
7.3.1 Principle of Thomson scattering . . . . . . . . . . . . . . . . . . . . 161
7.3.2 Thomson scattering diagnostic system . . . . . . . . . . . . . . . . 164
7.4 Measurement of electron cyclotron radiation . . . . . . . . . . . . . . . . . 165
7.5 Measurement of bremsstrahlung spectrum . . . . . . . . . . . . . . . . . . 167
7.6 Neutral particle analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . 168
7.6.1 Neutral particle analyser system . . . . . . . . . . . . . . . . . . . . 170
1
Figure 1: Professor Lyman Spitzer next to an early stellarator. (Figure credit: FUSION EXPO-
EUROfusion)
Prologue
The idea to copy the energy producing process of the sun arose already in the 1920s.
Researchers found that at very high temperatures, hydrogen nuclei in the solar core are
joined together to form heavier ions, releasing energy as they subsequently split into
particles with lower total mass. This process has become known as thermonuclear fusion.
In 1937, the exile German physicists Hans Bethe and Leo Szilard, discussed the possi-
bility of using deuterium instead of hydrogen. The German Houtermans, the Australian
Thonemann and the Englishman Thomson were thinking along similar lines. Shortly after
World War II, secret experimental programs were initiated both in England and in the
Soviet Union. The military potential of the fusion process was clear. However, the real
starting point for fusion energy research took place in 1951 when the dictator of Argentina,
Juan Peron, declared that a German-Argentinian physicist had demonstrated thermonu-
clear reactions under controlled conditions. A year later the claim was disproved but in
the USA, Lyman Spitzer of Princeton University was selected to lead a project to inves-
tigate the possibility of thermonuclear fusion and in the Soviet Union Lev Artsimovich
became leader of another team with the same mission. In figure 1, Spitzer is posing next
to the magnetic confinement device he invented and named stellarator.
When fusion research was finally declassified at the Geneva conference in 1958, it
turned out that the UK, the Soviet Union and the USA had been focusing their efforts
on different configurations. The UK had been experimenting with so called pinches,
which unfortunately proved unstable. The Russians had realized the strong stabilizing
effect of axial magnetic fields, a discovery that led them to the successful tokamak. The
American investigation of the stellarator did not lead to any significant achievements;
plasma diffusion losses exceeded classical predictions by far. When open conferences on
fusion now started, the researcher’s hopes of quickly overcoming the problems related
to controlled fusion were however dashed. Instabilities and energy losses were tough
problems, coming to dominate research for decades to come.
2
Chapter 1
In this chapter, we give the background and motivation to fusion research; the closely
related problems of global energy supply and climate change. Then follows a brief overview
of currently available energy sources and related technologies. We end the chapter with
an overview of fusion energy and current fusion research. This overview serves as an
introduction to the more in-depth presentation of the physics behind fusion energy, which
follows in the remaining chapters of the book.
3
Figure 1.1: Global primary energy consumption by source calculated with the substitution
method, re-calculating the non-fossil fuel inputs to take into account the fossil fuel energy
conversion inefficiency. Data provided by the BP Statistical Review of World Energy (Figure
credit: Our World in Data)
in which primary energy can be calculated, the ”direct” method and the ”substitution”
method. The direct method does not take into account the energy lost in the conversion
of fossil fuels to usable energy. The substitution method attempts to correct for this loss.
When we burn fossil fuel in a thermal power plant most of the energy we put in is lost in
the form of heat. Most fossil fuel power plants have an efficiency of 33 % to 40 %. When
we measure electricity generation from renewables, we are measuring the direct output,
with no losses.
In order to compare the different energy sources in a way that take into account the
low efficiency of fossil fuel energy conversion, the primary energy consumption in the BP
Statistical Review is calculated using the substitution method, in which the renewable,
non-fossil fuel primary energy consumption values are converted to primary energy values
by dividing them with an efficiency factor. The efficiency factor used is increasing from
36 % to 40 % over time, reflecting the gradual improvement of the fossil fuel power plant
efficiency. The direct method would clearly show a larger share of the fossil fuels in
the world consumption. However, the substitution method is considered to give a more
accurate understanding of how low-carbon energy sources are competing with fossil fuels.
The primary world energy consumption has increased roughly linearly from the 1950
level of 28 000 TWh to reach the current level of 176 000 TWh, according to BP Statistical
Review of World Energy data shown in figure 1.1. During the period from 1950 to 2021
the global energy consumption has grown by a factor of 6.3. How can we explain that large
increase? First we consider the world’s population growth. The world population in 1950
was 2.5 billion, which has grown to 8.0 billion in 2022, an increase by a factor of 3.2. This
population growth factor is about half of that of the energy consumption growth, which
4
indicates that the average energy consumption per capita has approximately doubled over
this period.
The trends of the energy consumption in various regions in the world are quite different,
as can be seen in figure 1.2. The energy consumption in the Asia-Pacific region, has had
a strong increase over the last 50 years, and the pace of the increase is even faster in the
last 20 years. The trend in the western, developed parts of the world such as Europe and
North America, on the other hand is flat, or slightly decreasing over the same period.
From the figure, it is also clear that the energy consumption in the Asia-Pacific region
surpassed that in Europe in the 1990’s, and that in North America a few years later,
and now in 2021 is larger that the two regions (Europe and North America) combined.
Clearly, the Asia-Pacific region is largely responsible for the overall increasing trend of
the world’s energy consumption. This is of course closely linked to the population growth
in the region, but also to the rapid industrial and economical development, which results
in an increase of the energy consumption per capita. The current main sources of energy
Figure 1.2: Trends of energy consumption in different regions of the world. Data from BP
Statistical Review of World Energy. (Figure credit: Our World in Data)
5
• geothermal energy
• nuclear power
These are the major alternatives available with the technological knowledge of today. It
is among these limited options we must seek the solution to the energy problems for the
world’s increasing population with higher and higher demand for energy. We must find a
long term solution for energy production, which provides the possibility of a society with
clean environment, stable climate and without the risk for energy related disasters. On
the long term, in the second half of this century, new nuclear technologies such as fusion
energy, is expected to contribute to the energy production and be a part of the mix of
sources in the energy system.
6
Figure 1.3: Global CO2 emissions by fuel or industry type. Data are shown for emissions due
to the fossil fuels coal, oil, and natural gas, the controlled flaring process of natural gas in oil
extraction, and industry processes such as cement production. Data from Our World in Data
Global Carbon Project (Figure credit: Our World in Data)
In the IEA scenario ”Net Zero by 2050”, the full range of energy issues, including
oil, gas and coal supply and demand, renewable energy technologies, electricity markets,
energy efficiency are examined, with the aim to provide an approach to the goal of net
zero carbon dioxide emission by 2050. According to this scenario, in the year 2050
• the world economy is more than twice as big today,
• the population has increased with 2 billion more people,
• but the global energy demand is 8 % smaller than today, and
• instead of fossil fuels, the energy sector is based largely on renewable energy,
• two-thirds of total energy supply is from wind power, solar energy, bioenergy,
geothermal energy and hydro power,
• solar becomes the largest source, accounting for 20 % of energy supplies,
• solar photo-voltaic capacity increases by a factor of 20,
• wind power increases by a factor of 10 from today,
• use of fossil fuels declines from almost 80 % of total energy supply today to slightly
over 20 %, and
• fossil fuels that remain are used in facilities fitted with carbon capture,
• electricity accounts for almost 50 % of total energy consumption,
• electricity plays a key role in all sectors - from transport and buildings to industry
- and is essential to produce low-emission fuels such as hydrogen,
• total electricity generation increases by a factor of 2.5, and
• almost 90 % of electricity generation comes from renewable sources, with wind and
solar photo-voltaic together accounting for nearly 70 %, and
• the remaining 10 % of electricity comes largely from nuclear power.
7
Table 1.1: Fossil fuel consumption and reserves in the year 2020. Data from BP Statistical
Review of World Energy.
Fuel Yearly consumption (TWh) Fuel reserves left (years)
Coal 44 000 139
Oil 51 000 54
Gas 39 000 49
8
be compressed and transported to storage sites such as deep geological formations. These
formations can be petroleum or natural gas fields or deep underground water-filled porous
rock reservoirs known as saline aquifers. A version of carbon capture that is economically
attractive is CO2 Enhanced Oil Recovery (CO2 -EOR). In this method, CO2 is injected
into petroleum fields that are in their final stage of production, in order to extend the
life-time of the field by raising the underground reservoir pressure, which then enables
more petroleum to be extracted from the field.
There are currently around 35 commercial facilities applying carbon capture to in-
dustrial processes, fuel transformation, and power generation with a total annual capture
capacity of 44 million tonnes of CO2 . One example is the Norwegian Sleipner natural
gas field in the North Sea5 . Since the mid 1990s, one million tonnes of carbon dioxide is
annually stored in an underground reservoir. The carbon dioxide captured in the Sleipner
field is the result of purification of the natural gas extracted from the field. The captured
carbon dioxide is compressed to 65 bar and injected into the reservoir at some 1000 meters
depth below the ocean floor in the Utsira aquifer. The storage is considered safe and the
potential is huge. It is estimated that the Utsira aquifer, which is 200-300 meters thick
and covers an area of 26 000 square kilometers, has storage potential corresponding to
the carbon dioxide emission in waste gas from all European power plants for hundreds of
years.
The captured CO2 can also be used in industrial processes. World wide, around 230
million tonnes of CO2 are currently used each year, of which 130 million tonnes is used
in fertilizer production and 90 million tonnes is used in enhanced oil recovery.
In the IEA Net Zero Emission by 2050 scenario, the world-wide CCS capacity is
envisioned to increase from today’s capacity of 44 million tonnes to 8 billion tonnes in
2050. In the NZE scenario, the use of fossil fuels is dramatically reduced by 2050 so that
this capture capacity would be enough to capture largely all emission from burning of
fossil fuel and thus achieve the goal of net zero CO2 emission.
9
Figure 1.4: Levelized cost-of-electricity (LCOE) for different renewable energy technolo-
gies. Levelized cost estimates the average cost per unit of energy generated across the
life-time of a new power plant. Source: International Renewable Energy Agency (IRENA),
Figure credit: Our world in data
range 0.05 - 0.18 USD/kWh, a price range characteristic of fossil fuels. In most places,
power from new renewable power plants is now cheaper than from new fossil fuel plants.
A striking example is the price per Watt for solar modules, which has declined by
99.6 % from 106 USD/Watt in the mid 1970’s to 0.38 USD/Watt in 2020. The fundamental
driver of this change is that renewable energy technologies follow learning curves, which
means that with each doubling of the cumulative installed capacity, their price declines
by some fraction. For example, the price of solar PV modules has dropped by 20 % for
each doubling of the installed capacity. This price reduction factor is called the learning
rate, which then for solar PV modules is 20 %. Wind power also follows a learning curve,
which for on-shore wind power is 23 %. Off-shore wind power has a learning rate of 10 %,
and is still relatively expensive.
10
lations can be integrated into buildings in urban areas. Solar energy is a variable energy
source, which preferably is combined with energy storage systems or hydrogen production.
Concentrating solar power (CSP) plants, using mirrors to track the path of the sun,
concentrate the solar rays to a solar collector where a heat transfer fluid circulates. The
heat can be used for a variety of applications: electricity generation, water desalination,
or chemical production. Heat can be stored in molten salt during day for electricity
production at night-time. The Noor power station in Morocco is in 2021 the world’s
largest concentrated solar power plant with more than half a million mirrors distributed
over the 25 km2 site area.7 The plant uses molten salt storage to produce electricity
at night. Total capacity is 510 MW electricity generation, and annual production is 1.5
TWh.
Solar photovoltaic (PV) cells are made from semiconductor materials that eject elec-
trons when light strikes the surface, producing an electric current. PV cells are small,
and produce a few watts of direct current (DC) electricity. The cells are combined into
PV modules that typically comprise a rectangular grid of 60 to 72 cells.
World-wide the solar PV generation increased by a record 179 TWh ( 22 %) in 2021, to
reach a total of 1000 TWh per year. It is the second largest growing renewable technology
after wind power in 2021. In order to follow the IEA Net Zero Emission by 2050 Scenario,
the annual growth of solar PV generation by 25 % is needed in the period 2022-2030,
and the total solar PV generation in 2030 should reach 7400 TWh. Utility-scale solar
PV is the least costly option for new electricity generation in a majority of countries.
Distributed, small scale solar PV, such as rooftop solar on buildings is expected to grow
faster as a result of higher electricity prices.
11
Table 1.2: Size evolution of wind turbines over time. Source: Environmental Geology by
Steve Earle, 2021, Thompson Rivers University
Year Maximum power (MW) Rotor diameter (m)
1990 0.5 40
2000 2 80
2010 7 126
2020 14 200
Note the cubic dependence on wind speed, which shows the importance of locating wind
power plants to areas with favourable wind conditions. Wind power technology has devel-
oped significantly during recent years, with bigger and bigger wind power turbines being
constructed, see table 1.2.
In 2021, wind electricity generation increased by a record 273 TWh (17 %), to reach
a total production of almost 1 900 TWh per year. This is the highest annual growth
among all renewable power technologies. In order to get on track with the Net Zero
Emissions by 2050 scenario, which has 7 900 TWh annual wind power generation in 2030,
it is necessary to raise average annual capacity additions to 250 GW, more than double
the 2020 growth. In the IEA analysis of the technical potential of offshore wind, the
study showed that close-to-shore wind sites globally could provide almost 36 000 TWh of
electricity per year, which is close to the global electricity demand projected for 2040.
The cost-of-electricity (COE) for on-shore wind power is at reasonable levels, with low
investment costs as compared to several other energy sources.
Currently, Europe’s largest wind power installation is the Markbygden wind farm
located in northern Sweden 9 . The wind farm is comprised of 329 wind turbines in
operation with a total capacity of 1.12 GW (year 2021), and 307 wind turbines under
construction to be completed in 2022-23 with a total capacity of 1.44 GW. Once completed
the total planned capacity is 4 GW and production is estimated to 8-12 TWh per year.
12
requirement.11 CorPower’s WEC devices generate up to 300 kW each, and are planned to
be combined into clusters, with 10-20 MW capacity. These clusters are in turn planned to
be laid out side-by-side to form larger wave farms that can have a utility-scale production
capacity. Prototype tests of the WEC devices will start in 2022 off the coast of Portugal.
1.3.4 Hydropower
Hydropower, or hydroelectric power, is one of the oldest and largest resources of renewable
energy. It uses the natural flow of moving water to generate electricity. Water reservoirs
next to the plant stores the energy. Hydropower technologies generate power by using
the elevation difference, created by the reservoir dam. It generates electricity in a turbine
driven by water flowing in on one side of the dam, and flowing out on the other. The
power generated from a hydropower plant is proportional to the water flow and the height,
and can be calculated by a simple formula 12 . The potential energy of the water at the
top of the dam is
W = mgh (1.3)
where m is the mass of water in kg, g = 9.81 m/s2 is the acceleration due to gravity,
and h is the height difference between the two sides of the dam in m. The mass of the
water is obtained as m = ρV , where ρ = 1000 kg/m3 is the density of water, and V is
the volume of the water passing the turbine. Consider a small time interval ∆t in which
the water flows though the turbine. The volume of water passing the turbine is given by
the volumetric flow rate V̇ as V = V̇ ∆t. The maximum turbine power is obtained as
the kinetic energy of the water passing the turbine divided by the time interval ∆t. This
kinetic energy is equivalent to the corresponding change in the potential energy of the
water, and can be written as
W
P = = ρV̇ gh (1.4)
∆t
where P is the power in Watt. For example, a reservoir dam height of 100 meter and a
volumetric flow rate of 100 m3 /s produces 98 MW, which must be multiplied with the
efficiency of the plant, generally 80–90 %, when determining the electrical power.
The benefits of hydropower have been recognized and harnessed for thousands of years,
and it is the leading renewable energy source world wide. There are virtually no emissions
during the production phase, and it can be used both as baseload power and as balancing
power, compensating for variability in other renewable energy sources.
Currently (2022), the world’s largest hydropower plant is the Three-Gorges Dam in
China, which has 32 turbines and a total capacity of 22.5 GW of electricity power. The
annual production is 112 TWh.13
Hydropower is still the largest renewable energy source - excluding traditional biomass
- and accounts for more than 60 % of renewable electricity generation. It is expected to
remain the world’s largest source of renewable electricity generation. The global annual
hydropower generation in 2021 is 4 300 TWh. Global hydropower capacity (2021) is
1 360 GW. In the IEA Net Zero Emission by 2050 scenario, hydropower maintains an
average annual growth rate of about 3 % in 2022-2030, to provide approximately 5 700
TWh of electricity in 2030. Pumped storage hydropower plants represents 30 % of the
net hydropower additions in this forecast.
11
https://corpowerocean.com/wave-energy-technology/
12
https://en.wikipedia.org/wiki/Hydropower
13
https://www.power-technology.com/projects/gorges/
13
In summary, hydropower is an economical, continuous and reliable source of energy,
but in many countries, it has already been developed as far as concerns for environment
and population allows. The largest remaining potentials for hydropower expansion are
found in Africa, South America and Asia.
1.3.5 Bioenergy
Bioenergy is created in plants when the photosynthesis transforms water, carbon dioxide
and minerals to biomass using energy from sunlight. Bioenergy does not contribute much
to the greenhouse effect, since the carbon dioxide released during combustion was earlier
captured by the plants themselves. Inevitably, transports for harvesting and fertilizing as
well as processing require energy, which in turn should be sustainably produced.
Traditional use of biomass - burning of charcoal, organic wastes, and crop residues
- was an important energy source for a long period in human history. It remains an
important source in the developing countries today. However, high-quality estimates of
energy consumption from these sources is difficult to find. The BP Statistical Review of
World Energy only publishes data on commercially-traded energy, so traditional biomass
is not included. However, modern biofuels are included in the data. Bioethanol and
biodiesel - fuel made from crops such as corn, sugarcane, hemp and cassava, are now a
key transport fuel in many countries.
Bioenergy is a leading energy source in Sweden, and biomass has a dominant position
in the Swedish heat market, to a large part as fuel in district heating. 14
During recent years, the global biofuel energy consumption has increased from around
100 TWh in 1990, to 1000 Twh in 2021, a ten-fold growth. In the IEA Net Zero Emission
by 2050 scenario, an increase in the global bioenergy supply by 60 % is required from
2020 to 2030. However, the increase must not lead to negative effects, such as expan-
sion of cropland for bioenergy or conversion of existing forested land into bioenergy crop
production.
14
energy, of which 1600 litre/second hot water for district heating and 120 MW of electricity.
Ground source heat pumps are widely used in many countries to provide renewable
local heating and cooling of buildings. In these systems fluid circulates in a pipe in the
ground a modest depths below the Earth’s surface. The fluid absorbs heat from the
warmer soil, rock or ground water around it. A heat exchanger transfers the heat to
the building heating/cooling system, and the fluid then returns to the ground at a lower
temperature.
In 2021, the global geothermal power production was nearly 100 TWh.
Grid-scale storage plays an important role in the IEA Net Zero Emission by 2050 Sce-
nario, providing system services for grid stability and long term energy storage. Pumped
hydroelectric storage is the most widely used storage technology and it has significant
potential in several regions. Batteries are the most scalable type of grid-scale storage and
has seen a strong growth in recent years. Hydrogen is an emerging technology that has
potential for long term seasonal storage of renewable energy. Thermal and mechanical
energy storage play a relatively small role in current power systems.
1.4.1 Batteries
Battery energy storage systems (BESS) are one of the fastest growing technologies in the
sustainable energy industry. On its most basic level, a battery is a device consisting of one
or more electrochemical cells that convert stored chemical energy into electrical energy.
Each cell contains a positive terminal, or cathode, and a negative terminal, or anode.
Electrolytes allow ions to move between electrodes and terminals, which allows current to
flow out of the battery. Advances in technology and materials have greatly increased the
reliability, output, and density of modern battery systems, and economics of scale have
dramatically reduced the cost.
• Lithium-ion (Li-Ion) batteries. Many companies are developing larger-format
cells for use in energy-storage applications. Many also expect there be synergies
with the electric vehicles (EV) powered by Li-ion batteries. Li-ion batteries have
been deployed in multi-megawatt batteries for the grid services. Lithium-ion refers
to electrochemical cell which is characterized by the transfer of lithium ions between
the electrodes during the charge and discharge reactions.
17
https://energystorage.org/why-energy-storage/technologies/
15
• Lead-acid batteries. Lead batteries are the most extensively used rechargeable
battery technology in the world. They have an unrivalled track record for reliability,
which makes them the dominant battery in terms of MWh of production. Lead
batteries are widely deployed to support PV installations both in commercial and
domestic premises. Lead batteries have lower cost than other chemistries.
• Redox flow batteries (RFB). Redox flow batteries uses chemical reduction and
oxidation reactions to store energy in liquid electrolyte solutions which flow through
the battery during charge and discharge. The RFB provides a separation of power
and energy. The system energy is stored in the volume of electrolyte, which can
be in the range of kWh to MWh. The power of the system is determined by the
stack of electrochemical cells. The amount of electrolyte flowing in the stack is
rarely more than a few percent of the total, and flow can easily be stopped during a
fault condition. This is in contrast with packaged batteries such as Li-ion and lead
batteries.
• Sodium-sulphur (NaS) batteries. The system operates at a high temperature
300-350 ◦ C. The active materials in a NaS battery is molten sulphur as the positive
electrode, and molten sodium as the negative. The electrodes are separated by a
solid ceramic, sodium alumina, which also serves as the electrolyte. NaS battery
technology has been demonstrated at over 190 sites in Japan. More than 270 MW
capacity for 6 hours of daily peak shaving have been installed. In Abu Dhabi, 15
NaS systems provide 108 MW/648 MWh. The efficiency is high, typically 89 %.
• Zinc-bromine (Zn/Br) flow batteries. The zinc-bromine battery is a hybrid
redox battery, because much of the total energy is stored by plating zinc metal as
a solid onto the anode plates in the electrochemical stack during charge. During
discharge, the zinc metal is oxidized and dissolved in the aqueous anolyte. The net
DC-DC efficiency of this battery is in the range 65-75 %. Integrated Zn/Br energy
storage systems have been tested up to 1 MW/3 MWh for utility-scale applications.
The Pillswood BESS in UK with power capacity of 98 MW and capacity to store 196 MWh
of energy is currently (2022) Europe’s biggest battery storage system. 18 It is based on
Li-ion battery units developed by Tesla. The storage system will be connected to Dogger
Bank, the world’s largest offshore wind farm in 2023, under development outside the
North East coast of England.
The capacity of grid-scale battery energy storage systems are currently far smaller
than pumped hydroelectric storage, but grid-scale batteries are projected to account for
the majority of storage growth worldwide. Batteries are typically employed for sub-hour,
hourly and daily balancing. Total installed gird scale battery storage capacity is 16 GW
in 2021. In the grid-scale battery technology mix, lithium-ion battery storage is the most
widely used, making up the majority of all new capacity installed. Besides lithium-ion
batteries, flow batteries emerge as a breakthrough technology for stationary storage. In
July 2022, the world’s largest vanadium redox flow battery was commissioned in China,
with a power capacity of 100 MW and a storage volume of 400 MWh.
16
produce all of their energy when the sun is shining during the day. The excess energy
produced during peak sunlight is often stored in these facilities - in the form of molten salt
or other materials - and can be used into the evening to generate steam to drive a turbine
to produce electricity. Alternatively, a facility can use off-peak electricity rates which are
lower at night to produce ice, which can be incorporated into a building’s cooling system
to lower demand for energy during the day.
• Pumped Heat Electrical Storage (PHES). In Pumped Heat Electrical Storage,
electricity is used to drive a storage engine connected to two large thermal stores.
To store electricity, the electrical energy drives a heat pump, which pumps heat from
the ”cold store” to the ”hot store”, similar to a refrigerator. To recover energy, the
heat pump is reversed to become a heat engine. The engine takes heat from the hot
store, delivers waste heat to the cold store, and produces mechanical work. When
recovering electricity, the heat engine drives a generator. The expected ”round-trip”
efficiency is 75-80 %.
• Liquid Air Energy Storage (LAES). Liquid Air Energy Storage uses electricity
to cool air until it liquefies, stores the liquid ari in a tank, bring the liquid air back
to a gaseous state by exposure to ambient air or with waste heat from an industrial
process, and uses that gas to a turbine and generate electricity. Liquid Air Energy
Storage is sometimes referred to as Cryogenic Energy Storage (CES). LAES plants
can provide large-scale, long duration energy storage with 100 MW output.
1.4.4 Hydrogen
Hydrogen is a versatile energy carrier, which can help tackle various critical energy chal-
lenges. Hydrogen can be produced from almost all energy resources, though today’s use
of hydrogen in oil refining and chemical production is covered from fossil fuels. Clean
hydrogen, being produced from renewables, nuclear or fossil fuels with carbon capture,
can help to decarbonise a range of sectors, including transport, chemicals, iron and steel
production, where it is difficult to reduce emissions.
Hydrogen can also support the integration of variable renewables in the energy system,
being one of the very few options of long term storage of electricity over days, weeks or
months. Clean hydrogen can be produced from water in electrolysers, stored, and later
generate electricity in fuel cells. The hydrogen storage system is depicted in figure 1.5.
17
Figure 1.5: Schematic of the hydrogen storage system. Electricity generated by wind or solar
is user to power the electrolyser which produces hydrogen gas. The hydrogen gas is compressed
and stored in a tank. The hydrogen can later be used in a fuel cell to produce electricity and
heat. The heat from the fuel cell is stored in a hot water tank, and used for heating of tap water.
Use and storage. Hydrogen has a high energy content by weight, but not by volume.
The hydrogen density is 0.089 kg/m3 at the atmospheric pressure of 1 bar and room
temperature. The energy content of hydrogen is 33.6 kWh/kg, or around 3 kWh/m3 .
This is a rather low energy density, and hydrogen must be stored physically as a high
19
https://en.wikipedia.org/wiki/Electrolysis of water
18
pressure gas or as a liquid. When the hydrogen is compressed to 700 bar corresponding
to a density of 42 kg/m3 , the energy content is 1 400 kWh/m3 . This energy density is
compatible with fuel tanks of motor vehicles. Hydrogen has a wide flammability range
compared to other fuels, corresponding to hydrogen-to-air ratios between 4 % and 75 %.
Hydrogen is the smallest known molecule, and is prone to leakage. In a confined space,
leaking hydrogen can accumulate and reach flammable concentrations. Proper ventilation
and the use of detection sensors can mitigate this hazard. Several car manufacturers
today use compressed hydrogen tanks which are capable of 350 or 700 bars, depending
on the automotive type. Compressed hydrogen fuel tanks are now made of carbon fibre
composites or carbon fiber and metal alloys and composites. The inner liner of the tank is
a high-molecular weight polymer that serves as a hydrogen gas permeation barrier. The
application of hydrogen in the transport sector is
• fuel cell electrical vehicles
• hydrogen combustion engines
• hydrogen-methane mixtures for combustion engines
Fuel cells. A fuel cell is an electrochemical cell that converts chemical energy of the
hydrogen fuel and an oxidizing agent, such as oxygen in air, into electricity. The fuel cell
is composed of an anode, a cathode and an electrolyte membrane. Hydrogen is supplied
at the anode and air is supplied at the cathode. The hydrogen atoms are stripped of their
electron at the anode, and positively charged protons pass the membrane to the cathode,
and the electron flows as electric current in the circuit. The electrons combine with the
protons and the oxygen from the air at the cathode, to generate water and heat. The
overall chemical reaction is the reverse of that in the electrolyser:
2H2 + O2 → 2H2 O (1.6)
A typical fuel cell produces a voltage of 0.6-0.7 V. To deliver the desired power, the fuel
cells can be combined in series to yield higher voltage, and in parallel to allow higher
current in a fuel cell stack. Typical fuel cell efficiency is 40-60 % 20 Assuming 60 %
efficiency in the electrolyser electricity-to-hydrogen production and 50 % efficiency in the
fuel cell hydrogen-to-electricity conversion, the ”round-trip” efficiency of the hydrogen
storage is of the order of 30 %.
19
Binding energy per nucleon (MeV)
8 84
4 56 Kr
He Ni U
54
7
Li Fe
6 6
Li
4
3
He
2 2
D
0
0 40 80 120 160 200 240
Atomic mass number (A)
Figure 1.6: The average nuclear binding energy per nucleon ε = Eb /A as a function of
atomic mass number.
of the plants in operation today are used to provide daily balancing. The world’s most
significant investment in new pumped hydroelectric storage capacity is in China where
currently (2022) 50 GW of capacity is under construction.
20
Figure 1.7: Illustration of the concept of nuclear cross section
It equates energy and mass, and states that the energy E of a nucleon is equal to m,
the mass of the nucleon, times the square of c, the velocity of light. The energy of the
nucleus is thus the sum of the energies of the individual constituent nucleons calculated
with Einstein’s equation minus the binding energy Eb of the nucleus. Since the energy of
the nucleus is less than the sum of the constituent nucleon energies, it is equivalent to a
reduction in mass of the nucleus relative to the sum of the masses of the nucleons. The
difference is referred to as the mass defect.
Nuclear reactions
A nuclear reaction is a process in which two nuclei collide and produce one or more new
nuclides. The process can be viewed as occurring in two stages: In the first stage, energy
equivalent to the sum of the binding energies of the reacting nuclides is supplied, and the
nuclei are disassembled into free nucleons. In the second stage, the free nucleons assemble
into the final nuclides, and energy equivalent to the sum of the binding energies of the
product nuclides is recovered. The reaction must obey the law of energy conservation,
thus, the difference between the recovered and supplied energies the so called reaction
energy, must be released/supplied as kinetic energy. In the case reaction energy is released,
the reaction products born in the nuclear reaction receive high initial velocities, and this
kinetic energy is then transformed to heat in a nuclear reactor as nuclear energy.
From the nuclear binding energy curve in figure 1.6, it can be seen that the elements
with very low, or very high atomic mass number have the smallest binding energies, and
elements with medium atomic mass number (e.g. Fe, Ni, Cr) have the largest. Nuclear
reactions which transforms nuclides with small binding energies to nuclides with large
binding energies release energy. Clearly, nuclear energy can be released in nuclear reactions
of two types:
• splitting (fission) of a heavy nuclide
• merging (fusion) of two light nuclides
Reaction probability
The nuclear cross section σ is used to describe the probability of a nuclear reaction. The
cross section concept originates from the idea of a projectile hitting a target with a cross
section area σ, as depicted in figure 1.7. The larger the cross section area is, the higher
is the probability that the reaction will occur. Consider a thin layer of the target volume
with area A and thickness dx. The nuclei in the volume has density n, and thus the
number of nuclei in the layer is N = nAdx. Assume that the probability for a reaction
between an incident particle and a target nucleus is equal to the probability of the incident
21
particle colliding with a physical target having the cross section area σ. The probability
P of a reaction in the layer with thickness dx is then equal to the fraction of the layer
that is blocked for the incident particle:
It is here assumed that the layer is very thin so that P is less than unity. An incident par-
ticle moves an average distance λ between collisions in the target and the layer thickness
is a fraction of this distance, which is given by the probability as dx = P λ. The mean
free path is then expressed as
λ = dx/P = 1/nσ (1.9)
An incident particle with velocity v traverses the layer in a time interval dt = dx/v. The
probability for a collision in the time interval dt is then P = nσdt. The average collision
time interval in the target of an incident particle, the collision time is τ , and the layer
traversing time dt is a fraction of this time, which is given by the probability as dt = P τ .
The collision time is then expressed as
The nuclear energy released in the reaction is measured in unit ”MeV” or million
electron volts. One electron volt is equivalent to the work done on an electron with the
charge −q = e = 1.6 × 10−19 C as it is moved across a potential difference of one Volt.
The conversion from eV to Joule is simply 1 eV = 1.6 × 10−19 J. The energy 200 MeV
released in one fission reaction is thus equal to 200 × 1.6 × 10−13 = 3.2 × 10−11 J.
We can compare the energy released in nuclear reactions with that released in chemical
reactions as follows: The weight of one uranium atom is 235 × 1.67 × 10−27 = 3.92 × 10−25
kg. One kilogram of pure uranium-235 contains 2.54 × 1024 atoms, and if all undergo
fission reactions, the nuclear energy released is 2.54 × 1024 × 3.2 × 10−11 = 8.13 × 1013 J =
22.6 GWh. The chemical energy released in combustion of 1 kg of methane is about 15
kWh. The nuclear energy density is thus of the order of one million times higher than
that of chemical.
22
Figure 1.8: Fission cross section for some nuclides in the uranium-plutonium fuel cycle in
unit barn as function of the incident neutron energy in MeV. 1 barn is equal to 10−28 m2
highly energetic particles, receiving on average an energy of 1.5 MeV when they are born.
Some of these secondary neutrons split other fissile nuclei, while others are captured in
collisions with non-fissile neutron-absorbing elements in the fuel, such as U-238. When the
reaction becomes critical, the number of neutrons produced from fission is larger than the
number lost to capture, and a chain reaction is sustained. In a power producing reactor,
neutron absorbing components are inserted, in order to keep the reactor in a marginally
critical state. By varying the amount of neutron absorbing material in the reactor core,
the neutron multiplication factor is adjusted to be close to one, and the reaction can be
maintained self-sustained in a steady state.
23
Figure 1.9: Cross sections for neutron capture in U-238 (blue) and neutron fission of Pu-
239 (red). Cross sections are shown in unit barn as function of the incident neutron energy
in eV. Cross sections are comparable for higher neutrons energies, which is essential for
raising the plutonium breeding ratio to the order of unity in the fast reactor.
the hydrogen atom is almost identical to the neutron mass it works well as a moderator.
Lighter elements are in general more efficient as moderators. Typical moderator materials
are light (regular) water (H2 O), heavy water (D2 O), and graphite (C). A reactor that is
using a moderator to slow down neutrons is referred to as a thermal, or slow-neutron
reactor.
From the graph of the fission cross sections, figure 1.8, it is clear that Pu-239 is also
a fissile nuclide, with a fission probability similar to U-235. Several fission reactions are
possible, such as the one resulting in the nuclides xenon and zirconium, also producing
three secondary neutrons, as shown below
239 134
94 Pu +n→ 54 Xe +103
40 Zr + 3n + 200 MeV (1.15)
The fission of Pu-239 produces nuclear energy, which is added to the energy produced
by fission of U-235. In fact, about one-third of the total energy produced in the common
24
commercial reactors today is due to fission of plutonium. 21
The level of Pu-239 in the reactor is determined by the balance of Pu-239 production
and loss processes. Several reactions are possible, such as the production of Pu-239 due
to neutron capture in U-238, and the loss due to the neutron induced fission of Pu-239.
The U-238 isotope is referred to as a fertile isotope due to it’s ability to convert into
another fissile isotope. The cross-sections for the two processes are shown in figure 1.9. In
a conventional thermal fission reactor operating with slow neutrons, the cross-section for
Pu-239 fission is much higher than that for U-238 capture. This means that the Pu-239 is
lost almost instantly after it is produced, and as a result, the fraction of Pu-239 in the fuel
remains small. The sustainment of the chain reaction would in this case still be mainly
driven by the fission of the U-235 isotope. When the fissile U-235 in the fuel has been
consumed, the reactor fuel has to be replaced. In this time period, much of the U-238 is
still present in the spent fuel, and remains in the fuel that is removed from the reactor.
However, it can be seen in figure 1.9 that, for higher neutron energies, the U-238
capture and the Pu-239 fission cross-sections become comparable, with the consequence
that the production and loss of Pu-239 can be balanced. When the so called breeding
ratio is close to one, the fissile Pu-239 is replenished from U-238 at the same pace as it is
consumed. In this case, there can be a significant amount of Pu-239 present in the fuel,
contributing to the sustainment of the chain reaction. The level of Pu-239 in the fuel
remain in this case constant, while the U-238 is constantly consumed in order to replenish
Pu-239. As a result, the uranium fuel can be used for longer time, and the U-238 isotope
is utilized in the reactor, instead of being left unused in the spent fuel. A reactor that
operates without a moderator, in order to use fast neutrons in the fission process, is
referred to as a fast reactor. Another possibility, to utilize the plutonium fission reaction,
is to use reprocessed reactor fuel, so called mixed-oxide (MOX) fuel, containing a mix of
uranium and plutonium.
Thermal reactors
The main fission reactor technology used today is the thermal Water Cooled Reactor
(WCR). The WCR is referred to as a thermal reactor, since it incorporates a moderator,
which functions to slow down the neutrons produced in the fission reactions, reducing
their initial energies in the MeV range to low, thermal energies below 1 eV. The uranium
fuel pellets are stacked in thin, long cylinders placed vertically in the reactor core, so
called fuel rods. The water circulates around the fuel rods, and functions both as coolant,
removing the heat produced in the fission process, and as moderator, slowing down the
neutrons. The reactor is kept at a marginally critical state by adjusting the neutron
multiplication factor to a value close to one by inserting neutron absorbing material.
Typical materials used are boron, cadmium, silver, hafnium or indium. The material is
placed in long cylinders, so called control rods placed vertically, mixed with the fuel rods.
By moving the control rods in (out) of the core, the amount of absorbing material in the
reactor is increased (decreased), adjusting the neutron multiplication factor down( up).
The water can be regular water (H2 O), or heavy water (D2 O), where the regular hydrogen
isotope 11 H is replaced by a heavier hydrogen isotope deuterium 21 H, often more simply
denoted by the letter ”D”. The reason for using heavy water, which has to be produced
21
https://world-nuclear.org/information-library/current-and-future-generation/fast-neutron-
reactors.aspx
25
from regular water at some cost, is that is a better material as a moderator. Heavy water
is also effective in slowing down neutrons in elastic scattering collisions, while it absorbs
significantly less neutrons in inelastic collisions. This characteristic of the moderator is
important, in order for the reactor to sustain the fission chain reaction, and it enables
the use of natural uranium with a low percentage (0.7 %) of fissile U-235 as fuel, a main
advantage of the Heavy Water Reactor (HWR).
The water cooled reactor technology, dominating today, is the Light Water Reactor
(LWR), which uses regular water. It requires that the uranium fuel is enriched to 3-5 %
of fissile U-235, in order to sustain chain reaction. Plants for enrichment of uranium are
capital intensive, and there are presently only a relatively operating facilities worldwide,
in France, Germany, Netherlands, UK, USA and China. 22 . The commercial process
employed for the enrichment is based on the mass difference between the U-235 and
U-238 isotopes, and involves spinning gaseous uranium fuel in fast rotating centrifuge
devices. Since the enrichment processes requires uranium in a gaseous form, uranium
oxide from the mine is converted to uranium hexafluoride in a preliminary process.
The water cooled reactors are of two types. In the Boiling Water Reactor (BWR), the
water circulates in a single circuit. The water is heated to the boiling temperature, and
the steam produced drives a turbine, which is connected to a generator that produces
electricity. In the Pressurized Water Reactor (PWR), the water circulates in two separate
circuits, which are connected by a heat exchanger. In the primary circuit, the water is
kept at a high pressure of 150 bar, allowing the water to reach temperatures of over 300 ◦ C
without boiling. An advantage of the PWR is that the irradiated water flowing through
the reactor remains separate from the water in the secondary, outer circuit.
Fast reactors
Fast reactors are characterized by operation without a moderator, and as a consequence
the fission process is relying on highly energetic, fast neutrons. At higher neutron energies,
the production of the Pu-239 plutonium isotope due to neutron capture in the uranium
U-238 isotope can balance the loss of Pu-239 due to neutron induced fission, at a PU-239
breeding ratio close to one.
In a thermal reactor, the chain reaction is sustained by fission of the U-235 isotope,
and the reactor fuel has to be replaced when the U-235 content in the fuel is consumed.
At this time, most of the U-238 isotope in the fuel remains, and is left in the spent fuel
removed from the reactor. In a fast reactor, the fraction of plutonium in the reactor
fuel can be at a significant level, and the fission chain reaction can be sustained by the
plutonium reactions. The plutonium is continuously replenished from U-238 isotope,
which is consumed in the process. As a result, the uranium fuel rods can be used for
much longer time before replacement, and the U-238 isotope is utilized during the reactor
operation, instead of being left unused in the spent fuel. A major advantage of the fast
reactor is that the reactor, by using the main constituent of the uranium fuel in the
breeding process, uses the uranium fuel much more efficiently than a thermal reactor.
This is important for the long-term usage of fission reactors, since it enables them to
utilize the world’s limited uranium reserves more efficiently.
The fast reactor can not use water a coolant, since it would act as a moderator. Fast
reactors typically use liquid metals as coolants, such as sodium or lead at elevated tem-
peratures. Due to their high thermal conductivity, metal coolants remove heat efficiently,
enabling high power density in the reactor. This makes them attractive in applications
22
https://world-nuclear.org/information-library/nuclear-fuel-cycle/conversion-enrichment-and-
fabrication/uranium-enrichment.aspx
26
where size and weight must be limited, such nuclear reactors for powering naval ships
or submarines. A problem with sodium is it’s chemical reactivity, which requires special
precautions to prevent fires. If sodium comes into contact with water, it reacts to produce
sodium hydro-oxide and hydrogen, and the hydrogen burns in contact with air. Lead has
a relatively high melting point of 328 ◦ C, which makes refuelling difficult. The melting
point can be lowered to 123 ◦ C by forming an alloy with bismuth, but a problem is then
that the lead-bismuth alloy is highly corrosive to most reactor structural materials.
27
• Nuclear proliferation. There is a concern that the development of nuclear energy
increases the likelihood of proliferation of nuclear weapons.
Accidents. There have been only two major accidents at nuclear power plants, and
their impacts have been less severe than was feared.26 The most serious nuclear accident
took place in 1986 at the Chernobyl nuclear power plant in Ukraine. A failed safety test
caused two explosions, a fire that lasted over a week, and the release of large amounts
of radioactive material. In 2011, an earthquake outside the coast of Japan triggered a
massive tsunami, and three of the reactor cores at the Fukushima Daiichi nuclear power
plant melted, and there were several hydrogen explosions and release of nuclear material
into the environment.
The radioactive materials and ionizing radiation released in a nuclear accident is a
potential health risk for people. Much of what is known about cancer caused by radi-
ation exposure from nuclear accidents comes from research on the Chernobyl accident.
27
Approximately 600 persons who worked at the plant during the emergency received
very high doses of radiation and suffered from Acute Radiation Syndome (ARS), also
known as radiation sickness or radiation poisoning. Among those, the 28 workers who
received the highest doses died. In the hundred of thousands of people worked as part of
the clean-up crews in the years after the accident an increased risk of leukaemia has been
found. Among 5 million residents near Chernobyl it is found that children have got an
increased risk of developing thyroid cancer.
Nuclear waste. When the nucleus splits, most of the energy is released immediately and
carried off by the coolant. However, energy continues to be released for thousands of years
after the atom splits. The afterglow heat is what makes nuclear waste hazardous. This
delayed energetic emission means that the nuclear waste is highly radioactive. In practice
the spent fuel is kept underwater for 5-8 years in spent fuel pools until the radiation has
decayed to a level where active cooling is no longer needed. Then, in today’s nuclear
power plants, the nuclear waste is either recycled in a fuel recycling plant, or kept on site
as a temporary solution.
Currently, the main option for the long-term solution to the waste problem is per-
manent storage in stable underground geological formations. The Onkalo deep geological
repository in Finland will be the world’s first commercial site to come into operation, serv-
ing the Olkilouto power plant. 28 Deep in the bedrock in southwest Finland, 450 meters
below the Earth’s surface, fuel rods containing spent uranium pellets will be inserted into
copper canisters that rest in water-absorbing bentonite clay. This disposal concept has
been developed by the Swedish nuclear fuel and waste management company SKB, and
a similar storage site is being planned in Sweden at the Forsmark nuclear power plant.
Another approach to the nuclear waste problem that has been suggested is trans-
mutation, in which the waste is transformed to short-lived, or even stable isotopes by
irradiation with neutrons. This can be done in two ways, either by using the neutron flux
in a fission reactor, or with an accelerator-powered system.
Nuclear proliferation. Nuclear weapons are the most dangerous weapons on Earth.
One weapon can destroy a whole city, potentially killing millions of people, and jeopar-
26
https://world-nuclear.org/nuclear-essentials/what-are-the-effects-of-nuclear-accidents.aspx
27
https://www.cancer.gov/about-cancer/causes-prevention/risk/radiation/nuclear-accidents-fact-
sheet
28
https://psmag.com/ideas/the-hiding-place-inside-the-worlds-first-long-term-storage-facility-for-
highly-radioactive-nuclear-waste
28
dizing the natural environment and lives for future generations. The dangers from such
weapons arise from their very existence. Although nuclear weapons have only been used
twice in warface - in the bombings of Hiroshima and Nagasaki in 1945 - many still exist
in our world today.
Nuclear proliferation is the spread of nuclear weapons, nuclear technology, fissionable
material or nuclear weapons-making information to countries that do not posses these.
The United Nations has already in 1946, in the first resolution by the UN General As-
sembly established a commission to ensure that atomic energy is only used for peaceful
purposes. The International Atomic Energy Agency (IAEA) is known as the world’s
”Atoms for Peace and Development” organization within the United Nations.
The IAEA is the international centre for cooperations in the nuclear field, and it has an
important role to verify the fulfilment of obligations assumed under the Non-proliferation
of Nuclear Weapons Treaty (NPT) that came into force in 1970. 29 The treaty has
recognized the countries that have nuclear weapons as ”Nuclear Weapon States” NWS)
and those that do not possess nuclear weapons as ”Non-nuclear Weapon States” (NNWS).
Nuclear proliferation takes place when the NWS pass their nuclear weapons or technology
to the NNWS. Currently, there are 5 NWS (US, Russia, UK, China and France), and 190
NNWS. Four other states are believed to possess nuclear weapons, but are not recognized
as NWS (India, Pakistan, North Korea, and Israel). In 1993, North Korea became the
only state to withdraw from the treaty.
1. p + p + e → D +ν + γ+ 1.44 MeV
2. D + p → 3 He+γ+ 5.49 MeV
3. 3 He + 3 He → 4 He + p + p + γ + 12.86 MeV.
Here, e denotes electron, ν neutrino and γ gamma rays. The reaction velocity is very
low at the 15 million degrees, which is the temperature of the Sun. But, since the volume
is huge, 4 million tons of material is transformed into energy every second. This energy,
29
https://www.iaea.org/topics/non-proliferation-treaty
29
Figure 1.10: The fusion processes in the sun. See also section 1.2.1. (Figure credit:
http://FusEdWeb.llnl.gov/CPEP/)
mainly in the form of gamma radiation, is transformed into visible light during a long
(millions of years!) and complicated route through the core of the Sun.
30
design a fusion reactor. Fusion research is indeed a long-term effort, which has engaged
a large number of scientists all around the world for many years, and continue to do so
today.
The enormous potential advantages of fusion energy, which has motivated the resources
and human efforts devoted to fusion research world wide over many years, are mainly the
following:
• Abundant energy. In a fusion reactor, energy would be produced in nuclear fusion
reactions involving the two hydrogen isotopes deuterium (21 H) and tritium (31 H).
Deuterium can be extracted from sea-water, while tritium will have to be generated
in the fusion reactor from nuclear reactions in lithium. Terrestrial reserves of lithium
would permit the operation of fusion power plants for more than 1000 years, while
sea-based reserves of lithium would fulfil needs for millions of years.
• No carbon dioxide emission. As with other kinds of nuclear energy, fusion reactions
do not produce greenhouse gas emissions. The by-product of fusion is helium, an
inert, non-toxic gas.
• No long-lived radioactive waste. Nuclear fusion reactors will not produce high ac-
tivity, long-lived nuclear waste. The activation of components in a fusion reactor is
anticipated to be low enough for the materials to be recycled or reused within 100
years.
• Limited risk of nuclear proliferation. Fusion reactors do not employ fissile materials
like uranium and plutonium. There are no materials in a fusion reactor that could
be used to make nuclear weapons.
31
POTENTIAL
ENERGY
Z1 Z2 e 2
EMAX =
EMAX 4π ε 0 R 0
NUCLEI REPEL
(COULOMB REPULSION)
r
R0 DISTANCE BETWEEN NUCLEI
≈ ≈
NUCLEI ATTRACT
(NUCLEAR FORCE)
Figure 1.11: Potential energy as a function of the distance from the nucleus. Another
incident nucleus requires sufficiently high kinetic energy so that it can pass the positive
Coulomb potential barrier in order come close enough to be pulled in by the negative
potential of strong nuclear force. The peak of the barrier is at a distance of order R0 ∼
10−15 m
potential energy barrier of the Coulomb interaction, see figure ??. It should be mentioned
that quantum mechanical tunnelling allows fusion reactions to take place at lower kinetic
energies.
How do we make the particles come close enough to each other? After some consider-
ation, one might come up with at least three ideas:
1. the particles are confined and heated together (thus forming a plasma) so their
thermal velocity may overcome the necessary Coulomb potential,
2. the particles are brought together at high speed by an accelerator,
3. the positive charge of the nuclei is shielded by negative charges, and the nuclei are
“tricked” to approach each other.
The first idea is thermonuclear fusion, where the fuel atoms are heated to a tempera-
ture of several hundred million degrees. It has so far been dominating, and it has evolved
along two main paths. The first path is some form of magnetic confinement. Most of
this text book will focus on this successful method. The second path is based on the
possibility to ignore any confinement, except that of fuel inertia. If the fuel is compressed
and heated for a few nanoseconds, conditions required for a positive energy exchange may
be achieved. This is called inertial confinement and has, so far, seen a similar progress
as magnetic confinement. The terms laser fusion and ion beam fusion are also used,
depending on which method is used for heating of the fuel.
Possibility 2) might seem tempting. Why not simply use two ion beams of opposite
direction, with a relative energy which optimizes the cross section for fusion? Unfortu-
nately, a simple calculation shows that the cross section for Coulomb collisions is several
orders of magnitude greater, which means the particle beams will simply deflect.
Muons, particles with the same charge as the electron, but with a mass 207 times
greater than the electron’s, have mainly been discussed in connection with possibility 3),
so called muon-catalyzed fusion. By replacing the electrons in heavy hydrogen atoms with
the heavier muons, the nuclei will come sufficiently close to each other so that quantum
mechanical tunnelling will become frequent. Unfortunately, muons have a short lifetime,
which limits the number of fusion reactions each muon can catalyze, and since they are
energy-expensive to produce with accelerators, it is difficult to reach energy break-even.
32
Figure 1.12: Cross sections for different fusion reactions. (Figure credit: By permission of
Oxford University Press)
The released energy in the DT reaction is in the form of the kinetic energy of the α-particle
(used to denote to the 4 He nucleus ) and the neutron. No gamma radiation is released.
The first generation of fusion reactors are expected to use a DT-mixture as fuel. To
get an idea of the energy density of the fusion fuel, the amount of fuel needed every
second in a fusion reactor delivering 1000 MW electric energy can be estimated to be 0.01
g/s, or a little under a kilogram during 24 hours. Deuterium is abundant in nature; one
33
Figure 1.13: Fusion of hydrogen isotopes deuterium and tritium. (Figure credit: EUROfusion)
out of 6700 hydrogen atoms is the the form of deuterium. It exists in heavy water D2 O
molecules which is present in ordinary sea water. Tritium is radioactive, with a half-life of
12.3 years, thus existing only in very small amounts in nature. Therefore, the reactor core
in a DT reactor will be surrounded by a blanket containing lithium, which will produce
tritium by the nuclear fission reactions
6
Li + n → T + 4 He + 4.8 MeV
7
Li + n → T + 4 He + n - 2.5 MeV
The neutrons created in the DT reactions will generate new fuel! To enhance this effect,
some neutron multiplicating material will be included in the blanket. If the produced
deuterium and tritium are combusted at the same rate as it is produced, the reaction is
called fully catalysed, and there is no tritium to take care of.
Figure 1.12 shows the a few other interesting reactions. DD reactions (with the same
probability)
D + D → 3 He + n + 3.27 MeV
D + D → T + p + 4.03 MeV
have the advantage that no tritium is needed. The problem, as we will discover shortly,
is that both higher temperature and higher plasma density is needed in this “advanced
reaction” in order to have the same energy exchange as the DT reaction.
It is also possible to go one step further, to the so called exotic reactions which do not
generate neutrons. Thus, the activation of the reactor walls caused by the fast neutrons
is avoided. As seen in figure 1.12, the reaction
D + 3 He → 4 He + p + 18.3 MeV
has this advantage, and in addition it has a better reaction cross section than the DD
reaction at energies above 100 keV. It is, however, impossible to completely avoid neutron
production from DD reactions in the plasma. Further, 3 He is rare on Earth (but not on
the Moon).
There are other reactions which are completely neutron free. Figure 1.14 shows the
reactivity for some of these reactions. One of the most advantageous seems to be
11
5 B + p → 3 42 He + 8.7 MeV
which, however demands very high temperatures and plasma densities.
It is also important to take in consideration that neutron-free reactions transfer the
fusion energy to kinetic energy in charged particles. Methods for direct transfer of kinetic
energy to electric energy must thus be used. This can for example be achieved by locally
opening the magnetic field lines so that the fast particles can escape into a generator.
34
10-20
LEGEND
= p - 9Be
= 3He - 6Li
10-21 = p - 6Li
= 3He - 3He
= p - 11B
= p - 7Li
10-22
<σv> (m3/s)
10-23
10-24
10-25
10-26
100 101 102 103
T (keV)
The reactivity is
ξ12 =< σv >12 [m3 s−1 ] (1.18)
which we may regard as the efficiency of a reaction. The reactivities of some reactions as a
function of plasma temperature are shown in figure 1.15. The assumption has been made
that the plasma is thermalized, meaning that the particles are Maxwell distributed. A high
temperature plasma, just as an ordinary gas, relaxes rapidly into a Maxwell distribution
(the time is of the same order of magnitude as the Coulomb collision time). This is for
example relevant when the plasma has been heated by a current flowing through it. But,
as we shall see, other methods of heating, for example ion beams and micro waves, must
also be used to reach fusion temperatures. This disturbs the velocity distribution and a
35
Figure 1.15: Fusion reactivity < σv > for some reactions as function of the plasma temperature
in unit keV. (Figure credit: By permission of Oxford University Press)
high energetic non-Maxwellian component must be taken into account, making it more
complicated to calculate the reactivity.
36
which can be rewritten as
1 2 3neT
n ξDT eEDT > (1.24)
4 τE
The reactivity ξDT is, within 10 % proportional to T 2 for fusion relevant temperatures
T = 10 − 20 keV, which means ξDT /T 2 is nearly constant, and the value is approximately
ξDT
≈ 1.1 × 10−30 m3 s−1 (eV)−2 (1.25)
T2
We can then write the break-even condition as
12
nT τE > (1.26)
EDT (ξDT /T 2 )
in a stationary state
Pin = Ploss (1.30)
leading to the Lawson criterion,
or
η
Ploss < Pf us (1.32)
1−η
37
1021
nτΕ
(m-3 s)
1020
10 100
T (keV)
Figure 1.16: Lawson diagram for the DT reaction, calculated with a reactor efficiency of
η = 1/3. (Figure credit: By permission of Oxford University Press)
We note incidentally, that this criterion gives the condition Q > 1 for a theoretical round
trip efficiency of η = 1/2. However, for a realistic smaller value for η < 1/2, Lawson’s
analysis gives a more strict criterion, as we anticipated.
The fusion power density Pf us for an optimum fuel mixture n1 = n2 = n/2. This gives
us
1
Pf us = n2 ξ12 eE12 (1.33)
4
where E12 is the fusion energy per reaction in unit eV. The electromagnetic radiation loss
power density in form of bremsstrahlung emission is
where C = 1.7 × 10−38 W m3 (eV)−1/2 . Lawson writes the loss power density as
W
Ploss = Prad + . (1.35)
τE
This requires some explanation. In the fusion research literature, the energy confinement
time τE is sometimes used for characterizing the total loss power from the plasma, and
sometimes it is used for characterizing only the non-radiative part. In Lawson’s criterion,
τE refers to the non-radiative losses. The second term is the loss power density due to
all other processes, and it is expressed as the plasma kinetic energy density W = 3neT
divided by the energy confinement time τE .
The Lawson criterion is finally written as
3eT
nτE > √ . (1.36)
η 1
ξ12 eE12 − C T
1−η4
where T and E12 are in unit eV. Observe that the right hand side is function of the plasma
temperature only (for a given value of η). The marginal condition, the break-even case,
is illustrated in figure 1.16 for η = 1/3. The curve has a minimum at T ≈ 25 keV, where
nτE = 6 × 1019 m−3 s. We have seen that the Lawson criterion sets the required values for
38
Figure 1.17: Ignition by heating from α-particles. (Figure credit: By permission of Oxford
University Press)
plasma density, temperature and confinement time if the fusion reactions are to heat the
plasma. The plasma temperature should be about 10-40 keV for the DT and DD reactions
since minima for the product nτE as a function of T is in this range. The Lawson criteria
for the two types of fusion reaction are
DT-reaction nτE > 6 × 1019 m−3 s for T ∼ 25 keV
DD-reaction nτE > 4 × 1021 m−3 s for T ∼ 30 keV.
In a reactor, where we wish to produce for example 1000 MWel in excess of the
electricity used for plasma heating, we need to position ourselves a bit above the curve.
Furthermore, Lawson’s analysis tells us what temperatures, densities and confinement
times we must achieve in the different fusion experiments. In magnetic confinement,
densities should be n = 1020 –1021 m−3 and confinement times τE = 0.1–10 s.
1.6.8 Ignition
Since one of the fusion products always is an ion, it is possible to improve on the Lawson
criterion. If these fast ions can be made to remain in the plasma they will contribute to
the heating by colliding with the plasma particles, transfering their energy directly back
to the plasma particles, thus heating the plasma. This is a much more efficient way of
heating than, as the Lawson criterion assumed, letting the ions hit the wall causing heat
that is transformed to electrical energy. It is even possible that the fast ions will produce
so much heat that all external heating of the plasma can be turned off. This is called
ignition of the plasma.
In the DT case, it is possible to confine the fast α-particle during its deceleration.
This will cause α-particle heating with a power per unit volume of
1
Pα = n2 ξDT eEα (1.37)
4
where Eα = 3.5 × 106 eV is the kinetic energy of the α-particle. If we express the total
power loss density as
3neT
Ploss = , (1.38)
τE
the ignition condition becomes Pα > Ploss , or
12T
nτe > . (1.39)
ξDT Eα
39
where T and Eα both are in unit eV.
This is, of course, a stricter criterion than that of Lawson. But if we take into account
that α-particle heating is becoming significant at 5 keV, ignition may be achieved already
at 10 keV, as illustrated in figure 1.17.
However, a plasma does not have to ignite to be interesting for fusion. A driven
plasma can also be used in a fusion reactor, but, as mentioned above, the economy will
be much better if α-particles heat the plasma directly, instead of giving off heat which,
at a limited efficiency, is transformed to electricity.
An ignited plasma cannot become uncontrolled in the same way as in a fission reactor.
The slightest disturbance of the plasma equilibrium will cause energy losses and cooling.
To control a plasma which has reached ignition, a burning plasma, is however not without
problems.
The European programme also involves the development of a fusion materials test facility,
the International Fusion Materials Irradiation Facility - Demo Oriented Neutron Source
(IFMIF-DONES), to be constructed in Spain.35
Currently, the central fusion research facility in Europe is the JET device, at the
Culham Centre for Fusion Energy (CCFE) in UK. It is the largest and most successful
magnetic confinement fusion experiment in the world. JET is collectively used under
EUROfusion management, and over 350 scientists and engineers from all over Europe
31
https://energy.ec.europa.eu/topics/research-and-technology/strategic-energy-technology-plan en
32
https://www.iter.org/
33
https://euro-fusion.org/eurofusion/roadmap/
34
https://euro-fusion.org/
35
https://ifmif-dones.es/
40
Figure 1.18: Schematic drawing of the JET magnetic confinement device in UK. (Figure
credit: EUROfusion-JET).
contributes to the JET programme. In the beginning of 2022, JET scientists announced
a breakthrough result, in which roughly 10 MW of fusion power was produced for a time
period of around 5 seconds.36 The record is considered a clear demonstration of the
potential for fusion energy, and is fully in line with the predictions for ITER.
Later in the same year 2022, the Lawrence Livermore National Laboratory’s National
Ignition Facility (LLNL-NIF) in the US, announced another breakthrough result, reaching
scientific energy break even, using a different approach to fusion energy - inertial confine-
ment. Break even is considered as a very significant milestone in fusion research, meaning
that more energy was produced from the fusion process than the energy put in to drive
it. 37 These recent achievements in magnetic confinement fusion and inertial confinement
fusion research contributes to the current optimism about fusion and the belief that fusion
power plants some day will be realized and deliver electricity to the grid.
41
the 1.7 cm long muon traces corresponded to the 5.5 MeV released by fusion of a proton
and a deuterium nucleus,
µ + p + D → 3 He + µ + 5.5 MeV
What happened was that pDµ-molecules were formed, and the muon catalyzed the fusion
reaction by a drastic reduction of the distance between the proton and the deuterium
nucleus, thus allowing quantum mechanical tunnelling. Early hopes about fusion energy
were dashed, however, by calculations (which actually already existed) showing that each
muon on average catalyzed less than one fusion reaction. The reason for this is the short
life time of muons, about 2.2 µs. Another reaction, the DT reaction, was more promising
and the focus of research until the late 1980s. In this case, the muon will catalyze more
than a hundred fusion reactions during its life time.
A major problem with muon catalyzed fusion is the probability for the muon ”sticking”
to one of the produced α-particles. If the muon sticks to the α-particle, it will not catalyze
more fusion reactions, and the so called catalytic chain will be broken. The energy cost for
producing muons by accelerators is very large, in the GeV range, requiring that each muon
catalyze of the order of one thousand fusion reaction in order for muon-catalyzed fusion
to reach energy break-even. The required number of reactions per muon is larger than
what currently is observed in experiments. In order for muon catalyzed fusion to become
a viable fusion energy source, a novel, less energy-expensive method for muon-production
is needed.
42
Chapter 2
Thermonuclear fusion has proved a successful route to fusion. The principle is that a gas
mixture of for example D and T is heated to the plasma state and further to very high
temperatures and is confined in a suitable magnetic bottle. We shall now consider the
requirements on the plasma and confinement that must be fulfilled if reactor conditions
are to be reached. We shall also discuss the plasma models that we need in our further
analysis of the fusion plasma.
Before moving into the large and complex area of plasma modelling. We will need
some common physics constants, such as elementary particle masses and charges. Special
attention has to be paid to the energy and temperature units used in the formulas in this
book. Analysis of a plasma fluid model involves vector algebra, and differentiating vector
fields, using the nabla symbol. We often use the cylindrical coordinate system, which
requires some care when differentiating scalar or vector fields. Further some formulas
from electromagnetics is often employed, for example in order to compute magnetic fields
from currents, and vice versa, in some simple symmetric geometries. We have collected
all this information in an appendix at the end of this book, so it is easy to find.
3
W = (ne eTe + ni eTi ) = 3neT
2
It is here assumed a plasma with equal electron and ion density n and temperature T , and
with the plasma temperature T in unit eV. The volume average < . . . > is taken over the
plasma volume Vp and Ploss is the total power loss from the plasma volume. The energy
confinement time can thus be interpreted as a typical ”energy leakage time” or ”cool off”
time. Let us now assume the plasma to be confined in a tokamak, see figure 2.1. The
importance of τE follows from the Lawson criterion for break-even in a fusion reactor. A
large tokamak such as the JET device has τE ∼ 1 s.
43
Figure 2.1: Principal drawing of plasma and coils in a tokamak . The direction of the ring is
called toroidal and the direction along the edge of the plasma cross section is called poloidal .
(Figure credit: EUROfusion-JET)
The second important quality parameter for the magnetically confined plasma is the β
value, defined as the ratio of plasma kinetic pressure p to magnetic field pressure B 2 /2µ0 :
<p> < 2neT >
β= = . (2.2)
< B 2 /2µ0 > < B 2 /2µ0 >
44
product (eq. 1.28), we find the required parameter combination nτE T ∼ 1021 m−3 s keV.
Typical parameters for a large tokamak such as the JET device, with capacity to reach
near break-even condition, are
which are line with the triple product required. The stored plasma energy density for
these parameters is
The plasma volume of JET is of the order of Vp ∼ 90 m3 , which together with the energy
confinement time τE determines the necessary power input to keep the plasma at this
energy density
W 4.8 × 105
Pin = Vp ∼ × 90 ∼ 4.3 × 107 W ≈ 40 MW (2.6)
τE 1
The plasma pressure is
A typical beta value is β ∼ 5 %, and together with the plasma pressure, this determines
the required confinement magnetic field strength. Inserting parameter values we obtain
2µ0 p
B2 = ∼ 2 × 4π × 10−7 × 3.2 × 105 /0.05 = 16 T2 (2.8)
β
giving a required tokamak magnetic field of order B ∼ 4 T.
45
y
B
rL
vi
x
ve
during this long journey. Considering this, it seem like a good idea to investigate how
individual particles move in stationary magnetic and electric fields. In order to have a
convenient way of studying interactions in a system of particles, we subsequently introduce
particle distributions and kinetic theory. To be able to describe complicated geometries
in a reasonable fashion, we will need the so called fluid model .
46
components are
vx = −vθ sin ωt
vy = vθ cos ωt
vz = const (2.12)
which easily can be verified by insertion in 2.10. It follows that the direction of motion for
positive charges (ions) is in the negative θ-direction, or clock-wise in figure 2.2, whereas
the direction of motion for negative charges (electrons) is in the positive θ-direction, anti-
clockwise. We introduce v⊥ for the thermal velocity in the plane perpendicular to B, and
write (
−v⊥ for positive charges (ions)
vθ = (2.13)
v⊥ for negative charges (electrons)
It is convenient to introduce the following relation between vθ and v⊥
q
vθ = − v⊥ (2.14)
|q|
We integrate eq. 2.12 to get the position of the particle (xp (t), yp (t)) as
xp = rL cos ωt
yp = rL sin ωt (2.15)
The radius of the circular orbit in the xy-plane is
vθ vθ
rL = =− (2.16)
ω ωc
which can be written as
mv⊥
rL = . (2.17)
|q|B
In a thermal plasma, we can set v⊥ = vT , where the thermal velocity vT is obtained as
1
2
mvT2 = eT . The Larmor radius is then
√
2meT
rL = (2.18)
|q|B
Now, how large is the Larmor radius? This is obviously an important question for mag-
netic confinement. Definitely, rL /a 1 must be fulfilled, let’s say rL /a = 0.01 - 0.1.
From (2.17) it is given that the ion Larmor radius is dominating the electron Larmor
radius by a factor of (mi /me )1/2 . For deuterium, we have
√ √
2meT −4 T
rL /a = = 2.0 × 10 (2.19)
eaB aB
which gives rL /a = 2.0 × 10−3 for a = 1 m, T = 10 keV and B0 = 10 T, being comforting.
We also realize that rL /a ≈ 0.07 for the fast α particles created in the DT reaction.
The cyclotron frequency for the thermal deuterons is ω ≈ 5 × 108 s−1 for the above
parameters. Since the ions travel 1 m in about 10−6 s, they gyrate about 500 times around
the field line during this distance! (The magnetic field is not homogeneous in magnetic
bottles, but the purpose of the investigation above is to give understanding of principles
and orders of magnitude).
In conclusion, a particle rotates rapidly at a distance rL from the centre of gyration
of the field line in a homogeneous magnetic field. It is also drifting at constant velocity
along the field line.
47
z
B
E y
+ -
Vxd
x
Figure 2.3: Guiding center drift due to crossed electric and magnetic fields.
E × B drift
When electric fields affect the motion of the particle, it will cause a slow drift of the
orbiting particles guiding center (gyro-center) in a direction perpendicular to the mag-
netic field. Does this mean that confinement is lost? Let us investigate this. Assume
perpendicular homogeneous magnetic and electric fields B = Bz ẑz and E = Ey ŷ
y , as shown
in figure 2.3. The equation of motion then becomes
dv
m = q(E + v × B). (2.20)
dt
Writing out the components
dvx
m = qvy Bz
dt
dvy
m = −qvx Bz + qEy (2.21)
dt
dvz
m =0
dt
and introducing ωc , we get
dvx
= ωc vy
dt
dvy Ey
= −ωc vx + ωc (2.22)
dt Bz
which gives an equation for vy
d2 v y
= −ωc2 vy
dt2
The expression for the velocity components is finally obtained as
Ey
vx = −vθ sin ωt +
Bz
vy = vθ cos ωt (2.23)
vz = const
The guiding center drift velocity vE = Ey /Bz in the x direction is superimposed on the
gyration motion. We can write in more general form the expression for the E × B drift
in vector notation as
E×B
vE = . (2.24)
B2
48
z
B
Δ
i B
Vxd y
-
x + e
Vxd
This guiding center drift is equal for ions and electrons in both absolute value and direc-
tion!
∇B drift
For a homogeneous magnetic field, the Lorentz force averaged over one gyration period
< F > =< qv × B > is zero. Consider now an inhomogeneous magnetic field in the z
direction with variation of B in the y direction, as is shown in figure 2.4. In this case the
average Lorentz force will not be exactly zero, which will lead to a drift motion, called
∇B drift. In order to calculate this drift, we will consider the drift as a small perturbation
to the gyration motion, and regard the average Lorentz force experienced by the particle
over one gyration period as an external force. This external force is then inserted in the
general expression for the guiding center drift. Assume the magnetic field in the positive
z-direction.
< Fx > = q < vy Bz >
< Fy > = −q < vx Bz >
The particle velocity is approximately same as for homogeneous field:
vx = −vθ sin ωt, vy = vθ cos ωt
Taylor expand the Bz field around the centre of the gyration:
dBz
Bz ≈ B0z + y , y = rL sin ωt
dy
49
Fcf
VII
i e
Vdz R c Vdz
B
z
Figure 2.5: Illustration of the curvature drift for ions and electrons.
dBz
< Fx >= qvθ rL < cos(ωt) sin(ωt) >= 0
dy
dBz 1 dBz
< Fy >= qvθ rL < sin2 (ωt) >= qvθ rL
dy 2 dy
We note that
2
mv⊥
qvθ rL = −
B
and write the expression for the average Lorentz force
2
1 mv⊥
< F >= − ∇B (2.27)
2 B
Inserting this force in the general expression eq. 2.25 for the guiding-centre drift gives the
∇B drift
2
< F > ×B 1 mv⊥ ∇B × B
v∇B = 2
= − (2.28)
qB 2 qB B2
Curvature drift
When the magnetic field is inhomogeneous, the field lines are generally also curved.
Assume cylinder symmetry, and introduce cylindrical coordinates (R, φ, Z). The magnetic
field is in the azimuthal direction B = Bφ φ̂
φ. Particles are moving at constant speed along
field lines with a velocity v = vk φ̂
φ as shown in figure 2.5. The curvature of the field line
is given by the radial vector Rc = R from the center of curvature to the field line. In
a frame moving with the particle, the particle experiences a centrifugal force opposite to
the inertial force
dv dφ̂
φ
Fc = −m = −mvk (2.29)
dt dt
The particle angular position is
vk
φp (t) = ωt = t (2.30)
R
and we write
dφ̂
φ ∂ φ̂
φ dφp vk
= = − R̂ (2.31)
dt ∂φp dt R
The expression for the centrifugal force is
mvk2
Fc = R̂ (2.32)
R
50
z Δ
φ B
Bφ
-I
- Bφ + Fcf
-
R
Rc
The drift is obtained by inserting the centrifugal force in the general expression for the
guiding center drift, eq.2.25, which gives
2
Fc × B 1 mvk mvk2
vc = = B R̂ × φ̂
φ= Ẑ (2.33)
qB 2 qB 2 R qRB
The general expression for the curvature drift is
mvk2 Rc × B
vc = (2.34)
qB 2 Rc2
Particle motion in the toroidal direction along the curved field lines gives a vertical curva-
ture drift. The radius-of-curvature vector of the toroidal field is Rc = R, which is inserted
in the equation for the curvature drift.
In a practical situation, the thermal motion of the particles generally has components in
directions both parallel and perpendicular to the field lines, so both types of drift are
present simultaneously. Unfortunately, the gradient drift and the curvature drift are in
51
z
φ
Bφ
Bθ +
Bθ B
θ
Bφ Bφ
Bθ -I
-φ
+
the same direction, and they add up, giving the total vertical guiding center drift in a
toroidal magnetic field
m 1 2 mvT2 v2
vd = v∇B + vc = (vk2 + v⊥ )Ẑ ≈ Ẑ = T Ẑ (2.35)
qRBφ 2 qRBφ ωc R
Ions and electrons have guiding-centre drift in opposite vertical direction (positive ions
drift upward and negative electrons drift downward for the magnetic field direction shown
in the figure. Since electrons and ions are drifting in the opposite direction a vertical elec-
trical field is built up. From eq. (2.24) it can be seen that an E × B drift in the outward
radial direction results for both ions and electrons! The vertical guiding-centre drift in
a purely toroidal magnetic field prevents confinement. This is the reason why toroidal
plasmas must always have a significant field component in the poloidal direction.
52
- - - PLASMA + +
- - -- - - - + + +
+ +
-
- - + -- -
- - - --
- -
+ +
++ - ++
+ +
λ D -- - -- λD
+
++ +
+
φ
φo
x
φ = φo e
- λ
D
λD
given by me ve2 /2 = eT /2 (one degree of freedom) their thermal energy density is ne eTe /2.
Obviously, the separation can only be maintained on a distance less than the distance
that results if the kinetic energy of the electrons is less than the electrostatic energy. If
the two energy densities are equal, we get the so called Debye length
r
ε0 T
λD = . (2.36)
ne
For a fusion plasma, λD ≈ 10−5 – 10−4 m. This can be compared to the average distance
between the particles, which is about n−1/3 ≈ 10−7 m. Thus, every particle interacts with
(103 )3 = 109 other particles within the Debye sphere. Particles outside this volume are in
principle shielded. Not until we reach distances much larger than the Debye length is it
relevant to use fluid models, which do not include effects of charge separation. In these
models, the relation ni = ne is used and the plasma is considered charge neutral.
df ∂f ∂q ∂p
≡ + w · ∇f = −f ∇ · w = −f (∇q · + ∇p · ). (2.38)
dt ∂t ∂t ∂t
On the left hand side, we have
d ∂ ∂ ∂q ∂p
≡ +w·∇= + · ∇q + · ∇p . (2.39)
dt ∂t ∂t ∂t ∂t
The right hand side is equal to zero. This is obvious from Hamilton’s equations
∂q
= ∇p H
∂t
∂p
= −∇q H (2.40)
∂t
which transforms the right hand side to −f (∇q · ∇p H − ∇p · ∇q H) = 0. Thus, what
remains is the very general and very important kinetic equation
∂f ∂q ∂p
+ · ∇q f + · ∇p f = 0. (2.41)
∂t ∂t ∂t
Now, let us rewrite (2.41) for a plasma with low density, which means that collisions
between particles can be neglected on the time scale we are interested in. The equation
of motion becomes
∂p ∂v
=m = q(E + v × B) (2.42)
∂t ∂t
and (2.41) becomes the so called Vlasov equation
∂f q
+ v · ∇f + (E + v × B) · ∇v f = 0. (2.43)
∂t m
This is a very useful relation when collisions between particles can be neglected. The
applications are fast phenomena, instabilities and so on. Note that each kind of particle
species in a plasma is described by its own Vlasov equation!
is relatively easy to calculate from statistical Markov chain analysis and will contain con-
tributions from both dynamic friction and diffusion in velocity space such as angular
scattering. However, we shall not give the expression here. As a rule, the Fokker-Planck
equation must be solved numerically. Semi-analytical solutions can be achieved for ex-
ample from calculation of the deceleration of a particle beam in the plasma.
54
-e, m e
- v χ
r
+ Ze
Figure 2.8: Deflection of an electron due to the Coulomb force from an ion.
55
where
Ze2
r0 = (2.48)
4πε0 me ve2
The incidence parameter r = r0 corresponds to a deflection angle π/2. The cross-section
for these ”large-angle” collisions is σ = πr02 , and the corresponding collision frequency in
a plasma with ion density ni is then ν = ni σve = ni πr02 ve .
However, the actual collision frequency is much higher, since the contribution from
”small-angle” collisions at larger distances where r r0 dominates the summation over
field particles. The deflection angle vary with the incidence parameter as χ ∝ r−1 and at
the same time, the number of field ions N , within a circle of radius r, increase
P as N ∝ r2 .
It follows that adding contributions from small angle collisions gives χ ∝ r and, the
sum diverges.
It has be shown that the electrical field from a charged particle is shielded at a distance
greater than the Debye length λD (eq. 2.21). Thus, only field particles within a sphere
of radius λD around the test particle interact in the collision process. The collision time
is calculated by summing up contributions from field particles within this Debye sphere
under the assumption that the motion of field particles are uncorrelated.
Consider the change of the electron forward momentum p = me ve , as it moves through
a thin layer dx of field ions:
Z λD
dp = − ni ∆p(r) 2πrdr dx (2.49)
0
where ∆p(r) is the forward momentum lost by interaction with a field ion at the incidence
parameter r. The rate of change with time is obtained if we assume a test electron with
velocity ve passing the layer in a time interval dt = dx/ve , giving
Z λD
dp
= −2πve ni ∆p(r) rdr (2.50)
dt 0
56
For electrons and ions in a plasma, the velocity distribution of the electrons has to be
considered, and for a thermal distribution we get scaling of the collision frequency with
plasma parameters as
νei =∝ ni Z 2 ln ΛTe−3/2 (2.56)
ln Λ ≈ 15. (2.62)
From this, it can be seen that in a D plasma, the ion-ion collision time is about 60
times longer than the electron-electron collision time. For a D fusion plasma with density
n = 1020 m−3 and temperature T = 10 keV we get τee = 0.23 ms τii = 14 ms
τex = 420 ms. All these times are well below the lifetime of the plasma. In a fusion
experiment with T = 200 eV we get τee = 0.65 µs τii = 54 µs τeq = 850 µs instead.
57
2.6.3 Plasma resistivity
Now a brief study of the resistivity of the plasma, which is determined by collisions
between electrons and ions. Let us first consider the simplest model possible. Assume
that an external electric field E generates a constant current through the plasma. The
electrons with charge −e and mass me have reached a drift velocity ve , and a force balance
is achieved between the driving electric field and the friction caused by ion collisions with
characteristic collision frequency νei
If we further assume Ohm’s law in the simple form E = ηj ≈ −ηne eve , we have
me νei
η= . (2.64)
ne e 2
Inserting the expression for the collision frequency νei , eq. 2.55 we get
ni Z 2 e2 ln Λ
η= (2.65)
4πε20 ne me ve3
For electrons and ions in a plasma, the velocity distribution of the electrons has to be
considered, and for a thermal distribution we get
ni Z 2
η∝ ln Λ Te−3/2 (2.66)
ne
The plasma charge neutrality condition requires that the number of positive and negative
charges in the plasma balances, and it is expressed in term of the densities as ne = Zni .
Using this condition, we obtain
η ∝ Z ln Λ Te−3/2 (2.67)
We now extend this simple analysis to a plasma with multiple ion species. The force
balance equation is in this case written as
X
−eE = me ve νei . (2.68)
i
where the summation is over different ion species enumerated i = 1, 2, 3, . . .. The resis-
tivity for multiple ions is
me X
η= νei . (2.69)
ne e2 i
Inserting eq. 2.55 for the collision frequency we obtain
e2 ln Λ 1 X
η= ni Zi2 (2.70)
4πε20 me ve3 ne i
58
which is inserted in the expression for the resistivity
e2 Zeff ln Λ
η= (2.73)
4πε20 me ve3
For electrons and ions in multi-species plasma, the scaling of the resistivity with plasma
parameters is then obtained as
Spitzer performed a careful Fokker-Planck analysis and found, in the case of singly charged
ions, a more accurate expression for the resistivity
me
ηS = 0.51 . (2.75)
ne2 τe
Inserting the expression for the electron collision time, eq. 2.57, the so called classical
Spitzer resistivity is obtained as
The resistivity is obtained in unit ohm-meter, and the value for the electron temperature
is inserted in unit electron volt eV ( or simply volt).
The resistivity of the plasma is comparable to that of copper at T = 1.4 keV and is one
order of magnitude smaller for thermonuclear temperatures. We will find that this fact
makes it difficult to heat the plasma to thermonuclear temperatures by drawing a current
through it (“ohmic heating”). In a magnetized plasma, the Spitzer resistivity is the
resistivity in the direction along the magnetic field. It is noteworthy that perpendicular
to the magnetic field, the resistivity increases with almost a factor of two; η⊥ = 1.96ηk .
59
Parallel to the field lines, fluid models are formally not valid! However, kinetic effects can
often damp the deviations from the MHD model, in particular the fast dynamics along
the field lines. But there are still a number of phenomena that the fluid model is unable
to handle.
Here P is the pressure tensor, which can describe also an anisotropic plasma pressure
and non-diagonal contributions from viscosity. For an isotropic distribution function the
tensor reduces to a scalar pressure and is given by
Z
1
p = m f (r, v0 , t)(v0 − v)2 dv0 (2.79)
3
We begin the derivation by integrating eq.2.77 over velocity space. The first term becomes
Z Z
∂f 3 0 ∂ ∂n
dv = f d3 v 0 = (2.80)
∂t ∂t ∂t
where n = n(r, t) is the plasma density. The second term becomes
Z Z
v · ∇f d v = ∇ · v0 f d3 v 0 = ∇ · (nv)
0 3 0
(2.81)
We will now show that the third term of the left side of (2.77) is equal to zero. With
F = q(E + v × B) we have
Z
F
· ∇v0 f d3 v 0 = {∇ · (φA) = φ∇ · A + A · ∇φ} =
m
Z Z
F F
= ∇v0 · ( f )d v − f ∇v0 · d3 v 0 =
3 0
(2.82)
m m
Z I I
F
= { ∇ · AdV = A · dS} = f · dS − 0 = 0
V S m
60
Here, we have used
F q q
∇v 0 · = ∇v0 · (E + v0 × B) = (0 + B · (∇v0 × v0 )) = 0 (2.83)
m m m
and the fact that f → 0 on the surface of the volume in velocity space as it expands to
contain all particles (the surface integral in (2.82)). The result is the so called equation
of continuity
∂n
+ ∇ · (nv) = 0 (2.84)
∂t
The collision term on the right side of (2.77) does not contribute, since collisions cannot
change the number of particles in phase space, and (2.84) describes a conservation, or
continuity, of particles in the phase space. This can be seen if we integrate (2.84) over
space; Z I
∂
ndV = − nv · dS. (2.85)
∂t V S
We may interpret (2.85) as that the time evolution of the total number of particles in a
volume is given only by the particle flow in and out of the surface of the volume. In a
similar way, we can multiply by v0 and integrate. After some calculations, which are not
given here, the result is the equation of motion
∂v
nm( + (v · ∇)v) = −∇ · P + nq(E + v × B) + R (2.86)
∂t
The term R is the change of the momentum caused by collisions with all kinds of particles.
Note that for each kind of particles (different hydrogen ions, electrons, impurities) there
is a pair of equations (2.84) and (2.86). Next we may calculate an energy equation which
describes energy changes in the plasma (heat generation, heat convection and so on) by
multiplying (2.44) by (v0 − v)(v0 − v) and integrating over velocity space.
∂ρ
+ ∇ · (ρv) = 0 (2.88)
∂t
61
By adding the equations of motion for electrons and ions we obtain the equation of motion
for the fluid
∂v
ρ + (v · ∇)v = j × B − ∇p (2.89)
∂t
Where we have used j = e(ne ve + ni vi ) ≈ enve and assumed isotropic plasma pressure.
Now, we need an equation for p. In principle, this can be found from the equation of
energy discussed earlier. Even early researchers, however, realized the importance of the
special case which is also present in classical thermodynamics, i.e. adiabatic behaviour.
Here, the dissipation is neglected, and the time is assumed too short for the small plasma
elements to be able to exchange heat with their surrounding. This gives the relation (the
adiabatic index Γ = Cp /CV , which is often assumed to be 5/3)
d p
= 0. (2.90)
dt ρΓ
Ohm’s law (here in simplified form) for a plasma is;
E + v × B = ηj. (2.91)
In ideal MHD the resistivity is equal to zero. Equation (2.91) then expresses that the
electric field in a (Lorentz) system which follows the motion always is zero because of the
infinite conductivity of the plasma. To close the system of equations, we add Maxwell’s
equations (the displacement current is neglected):
∇·B=0 (2.92)
∂B
∇×E=− (2.93)
∂t
∇ × B = µ0 j. (2.94)
Let us summarize the ideal MHD equations:
∂ρ
+ ∇ · (ρv) = 0 E+v×B=0
∂t
∂v ∂B
ρ + (v · ∇)v = j × B − ∇p ∇×E=− (2.95)
∂t ∂t
d p
=0 ∇ × B = µ0 j
dt ρΓ
This is a simplified model, but it still cover many important phenomena in a fusion plasma.
Even in this simple model, we have to handle 14 scalar coupled nonlinear time dependent
differential equations in three space dimensions!
We can expand the continuity equation to obtain an equation for the time rate of
change of the density as follows:
∂ρ
+ (v · ∇)ρ = −ρ∇ · v (2.96)
∂t
The left hand side of the equation is the convective time derivative of the density, which
represents the time rate of change of the density of a fluid element that follows the flow of
the fluid. The term on the right hand side represents the time rate change of the density
in that fluid element due to compression or expansion of the fluid. We may more simply
write the continuity equation as follows:
dρ
= −ρ∇ · v (2.97)
dt
62
The equation for the pressure may be expanded as follows
dp Γ dρ
=p (2.98)
dt ρ dt
In the equation for the pressure, the time derivative must be understood as the convective
time derivative. The equation can be combined with the density equation to obtain
dp
= −Γp∇ · v (2.99)
dt
Again the right hand side of the equation represents the time rate of change of the pressure
in a fluid element due to compression or expansion. In an incompressible fluid, the terms
on the right and side of the density and pressure equations would be zero. One important
simplification can, in some cases, be the assumption of incompressibility.
The electric field can be eliminated from the system using Faraday’s law. Combining
Faraday’s law
∂B
∇×E=− (2.100)
∂t
and the ideal-MHD Ohm’s law with zero resistivity,
E+v×B=0 (2.101)
we have
∂B
= ∇ × (v × B). (2.102)
∂t
The ideal MHD equations can then be written in the alternative form
∂ρ
+ (v · ∇)ρ = −ρ∇ · v
∂t
∂v
ρ + (v · ∇)v = j × B − ∇p
∂t
∂p
+ (v · ∇)p = −Γp∇ · v
∂t
µ0 j = ∇ × B
∂B
= ∇ × (v × B) (2.103)
∂t
63
Chapter 3
EQUILIBRIUM
While the current heats the plasma through resistive losses, it also counteracts the ex-
pansive force of the plasma column. Here, as elsewhere, we will assume that there is at
least one state available for the plasma, which does not change over time. We shall call
this the equilibrium state. Bennett’s equilibrium equations can be obtained from ideal
MHD if we assume a stationary state, i.e. if all time derivatives are set to zero. If we, like
Bennett, also assume that the equilibrium is static (v = 0) we determine the equilibrium
equation
∇p = j × B. (3.1)
The current j can be eliminated using Ampére’s law
∇ × B = µ0 j. (3.2)
∇·B=0 (3.3)
does the job. Equations (3.1)-(3.3) describe the equilibrium of a fusion plasma, and their
validity can be shown to be very general. Even in more advanced plasma models, the
64
terms in the system (3.1)-(3.3) are dominating. Thus, there is ample reason to study the
different solutions to the equilibrium equations.
Bennett’s pinch (named after the confining j × B force) has a circular cylindrical
geometry. In cylindrical coordinates, (r, θ, z) we have ∂/∂θ = 0 = ∂/∂z. In equation
(3.1) only the r component remains non-trivial:
dp
= jθ Bz − jz Bθ . (3.4)
dr
Further (3.2) gives
1 dBz
jθ = −
µ0 dr
1 d
jz = (rBθ ) (3.5)
µ0 r dr
Combining (3.4) and (3.5) we obtain the equilibrium equation for a circular cylinder ;
dp Bθ d Bz dBz
=− (rBθ ) − . (3.6)
dr µ0 r dr µ0 dr
Bennett studied a z pinch, in which Bz = 0. By multiplying (3.6) by r2 and integrating
from r = 0 to r = a, assuming p(a) = 0, he obtained the Bennett relation
µ0 Ip2 = 16πN eT̄ (3.7)
where the line density is
Za
N= 2πrndr (3.8)
0
and the Bennett temperature is
Za
1
T̄ = 2πrnT dr. (3.9)
N
0
Here the ideal gas law, combining electrons and ions, p = 2neT has been used. Equation
(3.7) is an integral form of the equilibrium equation, which in a clear way shows the
relation between a plasma current, density and temperature in a z-pinch. We see, for
example, that for N = 1020 m−1 and T = 10 keV the current I = 2.5 MA is required.
Is it difficult to solve the equilibrium equation (3.6)? It depends. If the magnetic
fields are specified as functions of r, only one integration is needed to calculate p(r) and
that would be it. If one of the magnetic fields must be calculated, the analytic integration
usually becomes less straightforward. It is clear, however, that two functions must be
given if the equilibrium equation is to be solved. Now, if (3.6) has several solutions,
which one will the plasma “choose”? The answer will come from a discussion of the
stability of the fusion plasma.
Circular cylindrical fusion configurations are problematic. Since the plasma is “open”
at the ends, the plasma suffers large losses of particles and energy. Further, the electrodes
at the edge usually cause vast amounts of contamination entering the plasma. Finally,
the z pinch and most other 1-D “screw pinches” (when Bz in included, the field lines form
a screw in the axial direction) are unstable. Thus, there is reason to close the ends of the
plasma to form a circular configuration, a torus. Consequently, we must also generalize
the equilibrium equation to two dimensions.
65
3.2 Beta value
The beta parameter can be defined in different ways. The two basic definitions are poloidal
β and toroidal β. The poloidal β is defined as
<p>
βp ≡ (3.10)
Bθ2 (a)/(2µ0 )
R
where the plasma pressure p is the average over the plasma volume; < p >= pdV /V .
For a straight circular cylindrical equilibrium, we can calculate
R Bθ (a) = µ02I/(2πa) and p
can be calculated by integrating eq. 3.6. If we write < p >= 2πrpdr/πa and integrate
by parts, we have
Za
1 dBφ2 2
βp = 1 + 2 r dr. (3.11)
a Bθ (a)2 dr
0
is instead often used, where B0 is the toroidal field on the magnetic axis r = 0.
66
8R
F = µo /(ln + � - 1) . (7.221)
l
2 2
FB2 =
h a
The outward force FP exerted by the plasma pressure is The f
F1 =
•
B Balan
111 l)9 µo l/
�
+ =
and t
Fig. 7.10
Poloidal magnetic field due to the combined plasma current and vertical magnetic field. B .1=
Figure 3.2: Plasma expansion
106
force balanced by an external vertical field. (From Plasma
Plasma Equilibria
Physics for Nuclear Fusion, K. Miyamoto, The MIT Press, 1980, p. 186)
�6
I
R
�------
67
I constant, the change in magnetic energy is
dW = Fh dR0 (3.15)
Here L is the sum of the internal and external self inductance of the current ring:
L = Le + Li (3.17)
It is easily verified that a constant current density corresponds to li = 1/2. The external
self-inductance of a current ring with minor radius a and major radius R0 is obtained by
integrating the poloidal magnetic flux that passes through the center of the ring, crossing
the circular area with radius R = R0 − a
8R0
Le = µ0 R0 ln −2 (3.21)
a
and the total inductance is thus
8R0 li
L = µ0 R0 ln + −2 (3.22)
a 2
Accordingly, the hoop force is
µ0 I 2
8R0 li
Fh = ln + −1 (3.23)
2 a 2
Plasma pressure
There is an additional outward force Fp exerted by the plasma pressure. Consider a thin
sector of the toroidal plasma ring, as shown in figure 3.4: The plasma pressure p on both
sides of the sector will result in a net radial outward force. For a plasma ring with minor
radius a the total force on one side of the sector (in the toroidal direction) is
68
Figure 3.4: Plasma expansion force balanced by an external vertical field. (From Plasma
Physics for Nuclear Fusion, K. Miyamoto, The MIT Press, 1980, p. 187)
where < p > is the average plasma pressure. If the sector has a toroidal angular extension
dφ, the resultant outward radial component of the forces from both sides of the sector is
dφ
dFp = 2Fpφ sin ≈ Fpφ dφ (3.25)
2
where < Bφ2 > is the average toroidal field in the plasma. This results in an inward total
radial force
< Bφ2 >
FB1 = −2π 2 a2 (3.28)
2µ0
69
The toroidal external field produced by the toroidal field windings is stronger on the
inboard side than the outboard side, varying with R as
R0
Bφv = Bφv0 (3.32)
R
which results in a net outward force on the plasma sector. By writing R = R0 + a cos θ
(figure 3.3), and using the Taylor expansion
1 1 a
≈ 1− cos θ (3.33)
R R0 R0
the outward force on the plasma ring becomes
2
B2
Z 2π
Bφv0
a 2 2 φv0
FB2 = − 2π R0 1 − cos θ a cos θ dθ = 2π a (3.34)
2µ0 0 R0 2µ0
Stability condition
The additional requirement of stability of the plasma equilibrium to vertical and horizontal
70
Z
B
dFZ =Ip×BR dl
BR Ip
dFZ
dFZ R
Ip
BR
Figure 3.5: Illustration of vertical restoring force from concave vertical equilibrium field.
displacements gives restrictions on the radial variation of the external vertical field. If the
vertical field lines are slightly concave towards the major axis, the interaction of the R
component of the field and the plasma produce a vertical restoring force, see figure 3.5.
The curvature of the external vertical field is described by the field index n, defined as
R ∂Bzv
n=− v (3.42)
Bz ∂R
The requirement that the field be concave for vertical stability gives the condition that
the field index must be positive
n>0 (3.43)
Considering radial displacements around the equilibrium major radius R0 , stability
requires that the net radial force is inward for an outward displacement, and vice versa.
This condition can be written as
∂FR
< 0. (3.44)
∂R0
It can be shown that this condition leads to the requirement
3
n< . (3.45)
2
In summary, the field index of the applied external vertical field has to fulfil
3
0<n< , (3.46)
2
to ensure simultaneous stability to both vertical and radial displacements.
71
Figure 3.6: Magnetic flux surfaces.
B · ∇p = 0
j · ∇p = 0 (3.47)
Thus, the current is flowing along the magnetic flux surfaces. These are defined by the
magnetic flux due to the poloidal magnetic field Bp , confined from the magnetic axis r = 0
and outwards in the radial direction, calculated as
ZZ ZZ I
Φ= Bp · dS = (∇ × A)p · dS = {Stoke’s law } = A · dl = 2πRAϕ (3.48)
S S Γ
We have introduced coordinates (R, ϕ, z) according to figure 3.3. The magnetic flux
inside each flux surface must be the same everywhere along the surface, and may thus be
characterised using Φ.
72
z
R
ϕ
Figure 3.7: Coordinate system for 2-D geometry in R, z. Axial symmetry (∂/∂ϕ = 0) is
assumed.
since only the magnetic field component Bφ is associated with a poloidal current. Hence,
RBϕ
χ= . (3.53)
µ0
It can be shown that χ, just like p, is a function only of ψ.
jp × eϕ Bϕ + jϕ eϕ × Bp = ∇p (3.54)
where eϕ is the unity vector in the toroidal direction. Then, we use the flux functions to
eliminate Bp and jp ;
1 ∂ψ ∂ψ 1
Bp = − , 0, = (∇ψ × eϕ ) (3.55)
R ∂z ∂R R
1 ∂χ ∂χ 1
jp = − , 0, = (∇χ × eϕ ) (3.56)
R ∂z ∂R R
Introducing this into (3.54) gives
1 1 Bϕ jϕ
(∇χ × eϕ ) × eϕ Bϕ + jϕ eϕ × (∇ψ × eϕ ) = − ∇χ + ∇ψ = ∇p (3.57)
R R R R
73
since eϕ · ∇ψ = 0 = eϕ · ∇χ. The chain rule then gives
dp
∇p(ψ) = ∇ψ
dψ
dχ
∇χ(ψ) = ∇ψ (3.58)
dψ
Inserting this into (3.57), we can separate jϕ ;
dp dχ dp µ0 dχ
jϕ = R + Bϕ =R + χ . (3.59)
dψ dψ dψ R dψ
We now rewrite jϕ using Ampére’s law (3.2):
∂BpR ∂Bpz 1 ∂ 2ψ ∂ 1 ∂ψ
µ0 jϕ = − =− − . (3.60)
∂z ∂R R ∂z 2 ∂R R ∂R
Finally, combining Eqs.(3.59) and (3.60) we obtain the Grad-Shafranov equation
∂ 2ψ 1 ∂ψ ∂ 2 ψ dp dχ
2
− + 2 = −µ0 R2 − µ20 χ . (3.61)
∂R R ∂R ∂z dψ dψ
This is a so called generalized partial differential equation; ψ is present both as a dependent
variable on the left hand side and as a independent variable on the right hand side! In
order to solve the Grad-Shafranov equation, the functions p = p(ψ) and χ = χ(ψ) are
first specified. Generally, the Grad-Shafranov equation is difficult to solve analytically.
Because of this, a large number of numerical codes have been developed.
An analytical equilibrium solution ψ(R, z) has been found by Solov’ev. By assuming
the following dependencies for the functions p = p(ψ) and χ = χ(ψ);
dp
µ0 = −a
dψ
dχ
µ20 χ = −bR02 (3.62)
dψ
where a, b and R0 are constants, the Grad-Shafranov equation can be solved analytically,
resulting in the following expression for ψ(R, z):
ψ(R, z) = 1
2
(bR02 + cR2 )z 2 + 81 (a − c)(R2 − R02 )2 . (3.63)
Here, R0 is the major radius of the magnetic axis, and c is a constant. Solov’ev’s equilib-
rium is mainly useful for testing equilibrium codes, due to its restricted modelling of the
current density jφ .
Then, why is the Grad-Shafranov equation important? In order to answer that, it is
useful to think about the meaning of magnetic confinement. The toroidal plasma must be
completely balanced by a magnetic field, so it does not drift towards a wall and in a way
which conserves its stability. As we will find in the next chapter, stability problems have
been dominating fusion research over the years. The Grad-Shafranov equation can be
used to calculate the design of the poloidal and toroidal coils balancing the plasma, and
the currents required in these. Just like the 1-D case, two functions (besides geometry
and boundary conditions) must be specified if the equation shall be solved. By studying
different dependencies for p(ψ) and χ(ψ) the stability can be investigated and optimized.
Then, when the plasma current, temperature, density and so on can be measured, more
precise solutions to the Grad-Shafranov equation can be calculated and used for a stability
analysis of the equilibrium.
74
Chapter 4
75
magnetic pressure reacting back on the plasma. There is also a ”tension” in the magnetic
field, which can create a force on the plasma. If the field lines are stretched out by
the plasma, the field will produce a force that reacts back on the plasma, opposing the
stretching. It is sometimes helpful to think of the magnetic field lines as elastic bands.
The forces acting on the plasma in the cylindrical equilibrium can be described as
follows: The radial variation (gradient) of the combined kinetic and magnetic pressures
produces ”pressure gradient forces” which, in case the plasma pressure is peaked on axis,
and the magnetic pressure is higher on the outside, can provide the required force equi-
librium. One such example is a plasma cylinder in an external axial magnetic field Bz . In
a cylinder with an axial plasma current, there is an azimuthal field Bθ produced by the
current. In this case the plasma and magnetic field pressures are not perfectly balanced,
and there is a net radially outward force. The force equilibrium is in this case main-
tained by a radially inward force produced by the ”tension” in the azimuthal magnetic
field lines. One can imagine the azimuthal magnetic field lines as rubber bands that are
wound around the plasma cylinder and acts to compress it radially.
If the plasma would have moved towards a surrounding wall, repelling currents may
have been generated, which would also be stabilizing. Intuitively, this can be seen as if the
magnetic fields lines between the plasma and the wall are compressed when the plasma is
dislocated towards the wall, thus yielding an opposing magnetic pressure. The magnetic
field may also be stabilizing if it has a convex curvature in the direction of the plasma;
we will discuss this later, in connection to pressure driven instabilities. But let us first
return to the important stabilisation due to field line bending.
For sake of simplicity, we consider a case where the equilibrium magnetic field lines
are straight. If the plasma, for some reason, would be disturbed from the equilibrium
state, and the perturbation would cause a bending of the magnetic field lines, a force
arises due to the stretching of the field, which is stabilizing, acting to straighten the field
lines, returning them back to the original state. Think of the magnetic field lines as
elastic bands, striving to be as straight as possible, and to return to their equilibrium
state after a disturbance. Imagine a plasma, originally at equilibrium, being displaced
somewhat. Consider a small volume of plasma with the equilibrium position r0 . Let
us call the displacement ξ ≡ r − r0 , where r is the new position. It is interesting to
study the change of the magnetic field due to the plasma displacement. We now write
B = B0 + B1 , where index 0 represents the equilibrium field and index 1 represents the
field at the perturbed position. Thus, we have linearized the magnetic field B about the
equilibrium position. In the same manner we may linearize the velocity field about the
static equilibrium position; v = v0 + v1 = v1 . If we introduce these expressions in eq.
2.102, 0:th order terms yield a trivial relation, and first order terms yield
∂B1
= ∇ × (v1 × B0 ). (4.5)
∂t
Since we are interested in B1 , and since v1 = ∂ξξ /∂t, we can integrate eq. 4.5 over time
to obtain
B1 = ∇ × (ξξ × B0 ). (4.6)
In order to illustrate the stabilization effect of field line bending, we can assume a homo-
geneous equilibrium field B0 . The, nabla operator in eg. 4.6 then only operates on the
perturbation ξ , and the equation can be written as
76
Then we make a plane wave ansatz, which introduces the wave vector k of the perturba-
tion, pointing in the direction of the wave propagation:
We may then replace ∇ → ik, transforming the partial differential equation to a simpler
algebraic equation:
B1 = i(k · B0 )ξξ − i(k · ξ )B0 . (4.9)
This equation is interesting. It shows that if the direction of the wave propagation is
perpendicular to the equilibrium magnetic field (k · B0 = 0), the perturbed field must
be parallel to the equilibrium field (B1 kB0 ). Hence, the perturbation does not bend the
magnetic field, and there is no reacting magnetic force due to stretching of the field. As
a result, the magnetic field’s defence against instability is weakened! It turns out that
perturbations with k · B0 = 0 are the most dangerous from the stability point of view.
With this notation, we have k ≡ (0, m/r, −n/R) and r ≡ (0, rθ, Rϕ) ; the periodicity
requires 2πr = m(2π/kθ ), or kθ = m/r and similarly for kϕ . Then, we can write
m n
k · B0 = 0 ⇔ Bθ − Bϕ = 0 (4.11)
r R
We now introduce a very important concept, the safety factor q:
rBϕ
q≡ . (4.12)
RBθ
The relation 4.11 tells us where in the plasma localized instabilities might arise; on flux
surfaces where q is a rational number.
m
q= . (4.13)
n
This rational number, m/n, is defined by the rotation of the perturbation along the
magnetic field. Thus, the condition (4.11) can be thought of as a resonance between the
perturbation and the equilibrium magnetic field.
But where in the plasma are the flux surfaces described in (4.13)? It turns out that
in most of the tokamak plasma q > 1, increasing monotonically radially outwards. This
means that the most global perturbation (m, n) = (1, 1) does not satisfy (4.13), which is
very advantageous from a stability point of view. However, close to r = 0, there is a small
region where we may have q < 1, and consequently there is a flux surface where q = 1 and
a resonance for the (m, n) = (1, 1) instability is possible. This instability indeed appears,
and it has a periodic growth and decay manner, which gives it it’s name, the sawtooth
instability.
77
It can further be concluded that no (1, n) mode, where n ≥ 2, can be resonant in the
tokamak plasma. This leads to good stability, since instabilities with low m, n usually
have the highest energy and thus are the most dangerous (the stabilizing effect of the
magnetic field increases with m and n). We will shortly study the instabilities which are
relevant for the tokamak. First, however, we will study q in a geometric context.
Differential geometry tells us that the equation for a field line can be written as
ds Bθ
= , (4.14)
Rdϕ Bϕ
where ds and Rdϕ are the poloidal and toroidal distances of an infinitesimal segment along
the magnetic field. The total change in ϕ as we make a complete turn in the poloidal
direction is I
Bϕ ds
∆ϕ = , (4.15)
RBθ
which, for a circular cylindrical plasma with R/a 1 is
rBϕ
∆ϕ = 2π = 2πq. (4.16)
RBθ
It is obvious that q = ∆ϕ/(2π) can be interpreted as the number of turns along a magnetic
field line in the toroidal direction (the ”long way”) for each revolution in the poloidal
direction (the ”short way”).
Now it is easier to realize why (4.11) is considered a resonance condition; when the
rotation of the perturbation exactly matches the rotation of the magnetic field, q becomes
a rational number. A physical interpretation of the instability risk associated with the
resonance is that the field lines on the resonance surface close after a few turns (fewer turns
if m and n are small). This means that any electrical field that arises because of some
perturbation is difficult to extinguish; the electrons are tied to field lines which do not
cover the flux surface. The electrical field can then act as a driving force for the unwanted
plasma movement we call instability. If q had been a non-rational real number, the field
lines, and thus the electrons, could have covered the whole flux surface ergodically.
In a general geometry, the safety factor should be calculated from
I
1 Bϕ ds
q= . (4.17)
2π RBθ
Since the integration follows a poloidal turn along a magnetic field line, q is a flux surface
function. A high q is desired from a stability point of view, and can be achieved for
example if the toroidal field is high or the major radius is small.
In some configurations, such as the stellarator, the term rotation transform ι, which
is related to the safety factor by ι ≡ 2π/q is often used.
78
• Stabilizing effects come from the plasma itself (magnetic field bending, curvature
and compression) and from a conducting wall.
• Ideal instabilities arise even if the magnetic field is perfectly conducting and not
destroyed anywhere.
• Resistive instabilities (or modes) can arise when the plasma is not perfectly con-
ducting.
• Rational surfaces for q are the flux surfaces which cause most problems for the
stability. This is true both for ideal and resistive modes.
We want to use a plasma model which includes the effects discussed above and at
the same time is simple. The simplest, physically reasonable, theory is the one-fluid
model ideal MHD (and resistive MHD if resistivity is included in Ohm’s law). But the 14
coupled nonlinear differential equations in three space dimensions and one time dimension
are usually far too complicated to handle. In this case, we can make it easy for ourselves
by studying the linear stability.
In linear stability analysis, the plasma is assumed to have attained some equilibrium
state. As shown earlier, a small displacement ξ (r) of the plasma, away from equilibrium,
is performed. In a stable configuration, the plasma returns to equilibrium. In other cases,
the equilibrium is unstable. The different cases in figure 4.1 provide intuitive illustrations
of this.
Figure 4.1: A ball in different stable/unstable situations:(a) stable, (b) unstable, (c)
marginally stable, (d) stable with dissipation (e) unstable with dissipation, (f) linearly
stable, non-linearly unstable, (g) linearly unstable, non-linearly stable
Four methods are used to study the stability. The first three can often be performed
analytically.
1. The energy principle is a way to investigate how the potential energy of the plasma
is affected by a perturbation.
79
Figure 4.2: Unstable motion in a tokamak, calculated from the eigenvalue equations for ideal
MHD. (Figure credit: By permission of Oxford University Press)
2. The eigenvalue problem for a certain perturbation can be expressed and solved using
the MHD equations. This results in eigenfrequencies, which tell us the growth rate
of the instability and also the eigenfunctions. One example of such a solution for
tokamaks can be seen in figure 4.2.
3. Marginal stability can be investigated by assuming that the growth rate is zero for
a given mode (m, n).
4. The initial value problem can be solved by numerical integration of the time devel-
opment after a given perturbation of the plasma.
Finally, it is important to note that ideal MHD stability almost always must be fulfilled,
while fluctuations caused by other instabilities often are present in the plasma.
4.4.1 Linearization
Let us exemplify linearization by continuing the discussion on charge separation. We
showed earlier that the Debye length λD is the largest distance in a plasma at which
charge separation can be maintained; at larger distances, the thermal energy is too large.
It turns out that the separation corresponds to a motion with a characteristic frequency.
We can derive this so called plasma frequency through linearization.
Let us assume that, for short times, the ions are stationary and the motion of the
electrons are parallel to the magnetic fields described by the fluid equations. We will
80
begin by considering the continuity equation;
∂ne
= −∇ · (ne ve ). (4.18)
∂t
If we now include the first order terms from a Taylor expansion around the equilibrium
state which existed just before the perturbation, we have
ne = ne0 + ne1
ve = 0 + ve1 (4.19)
where we have assumed that the fluid motion of electrons in equilibrium can be neglected.
All perturbed (first order) quantities are now expanded as in (4.10). Since the equations
we will study are linear, we may consider one mode at the time (algebra tells us they are
independent). Further, we see that the gradient and the time derivative of the perturbed
quantities can be replaced according to
∇ → ik
∂
→ −iω (4.20)
∂t
Equation (4.18) then gives one 0:th and one 1:st order equation:
∂ne0
=0
∂t
−iωne1 = −i ne0 k · ve1 (4.21)
The first equation provides no new information. Let us then linearize the electron equation
of motion and Poisson’s equation to obtain a closed system of equations:
dve
me = −eE ⇒ me (−i ω ve1 ) = −eE1
dt
ne e ne1 e
∇·E=− ⇒ ik · E = − (4.22)
ε0 ε0
Combining the algebraic equations (4.21) and (4.22) we find
s
ne0 e2
ω = ωp = . (4.23)
ε0 me
This characteristic frequency for plasma oscillations is named, not surprisingly, the plasma
frequency and has been measured experimentally. The example shows how powerful the
linearization tool can be in combination with Fourier methods. Complicated differential
equations are reduced to easily solvable algebraic expressions.
81
B and electric field E. There are four differential equations that couple these four vector
fields as follows:
dv
ρ =j×B
dt
j = ∇ × B/µ0
∂B
= −∇ × E
∂t
E = −v × B. (4.24)
The non-perturbed state is an infinite, static (v0 = 0), current-less (j0 = 0) plasma in a
homogeneous magnetic field B0 . We can then make a plane wave ansatz, assuming that
the perturbed quantities vary with space and time as
which converts the differential equations to algebraic equations. We linearize, and write
down the 1:st order equations:
−ρiωv1 = j1 × B0
j1 = ik × B1 /µ0
−iωB1 = −ik × E1
E1 = −v1 × B0 . (4.26)
We use a cartesian (x, y, z) coordinate system with the magnetic field B0 in the z-direction,
and assume an Alfvén wave that propagates along the magnetic field direction (k k B0 ).
From the equations we find the directions of the vector fields as follows: v1 ⊥ B0 , j1 ⊥ k,
and v1 ⊥ j1 . We can then orient the system so that v1 is in the y-direction, and j1 is in
the x-direction. Further E1 ⊥ B0 , E1 ⊥ v1 , and B1 ⊥ k, B1 ⊥ E1 , from which follows
that B1 is in the y-direction and E1 is in the x-direction. It turns out that this wave is a
transverse wave, with all perturbation vectors in a plane perpendicular to the wave vector.
The problem is then simplified to four scalar equations for four unknown perturbation
fields vy1 , jx1 , By1 , Ex1 ,
ρiωvy1 = jx1 B0
jx1 = −ikBy1 /µ0
−iωBy1 = −ikEx1
Ex1 = −vy1 B0 . (4.27)
and elimination gives the dispersion relation for the Alfvén wave
B0
ω = kvA , vA = √ . (4.28)
µ0 ρ
Already in 1942 Alfvén described the wave in a short paper in Nature. The discovery
contributed to Alfvén being rewarded the Nobel Price in physics in 1970. This wave is
very important in fusion plasma physics.
If the plasma pressure is included into the equation of motion, two more MHD waves
can be derived using the dispersion relation:
ω2 p
2
= 21 [c2s + va2 ± (c2s + va2 )2 + 4c2s vA2 cos2 θ]. (4.29)
k
82
Figure 4.3: Alfvén wave with kkB0 . (Figure credit: By permission of Oxford University Press)
Figure 4.4: General Alfvén wave. k and B0 are in the plane of the figure, while the velocity
fluctuations are perpendicular to this plane. The wave is propagating in the x direction with
velocity vA . (Figure credit: By permission of Oxford University Press)
Figure 4.5: Magnetoacoustic wave. The velocity fluctuations are in the same plane as k and
B0 . (Figure credit: By permission of Oxford University Press)
p
Here, cs = γp/ρ is the sound velocity and θ is the angle between k and B0 . The two
waves are called the fast and the slow magnetoacoustic wave, respectively. In the case of
the plasma being incompressible, only the Alfvén wave results. Notice the direction of the
three MHD waves in figures 4.3-4.5. The Alfvén wave is propagating by bending the field
lines (B1 ⊥B0 ) while the magnetoacoustic waves are mainly compressing the field lines.
Of course, a large number of other waves can also exist in a plasma, such as electro-
magnetic waves. These cannot be described using the MHD models, however, which only
cover low frequency waves with ω < ωci .
One may wonder why the above perturbations caused waves and not instabilities. The
answer is that there are no driving forces for instabilities in the homogeneous plasma we
discussed.
But how come ideal MHD can describe waves and instabilities, but not damped or
growing oscillations? If we analyze the linearized MHD equations, we realize that the
eigenvalue ω 2 always is real. This is a result of the hermitian character of the system,
which we will not discuss further here.
83
4.5 Ideal MHD stability
The starting point is the ideal MHD equations in the following form, see egn. 2.103;
∂ρ
+ (v · ∇)ρ = −ρ∇ · v
∂t
∂v
ρ + (v · ∇)v = j × B − ∇p
∂t
∂p
+ (v · ∇)p = −Γp∇ · v
∂t
µ0 j = ∇ × B
∂B
= ∇ × (v × B) (4.30)
∂t
Assume a static equilibrium v = 0, and linearize the equations. The density equation
is not needed, since the density does not appear in any of the other equations. The 1:st
order equations can then be written as
∂v1
ρ0 = j1 × B0 + j0 × B1 − ∇p1
∂t
∂p1
+ (v1 · ∇)p0 = −Γp0 ∇ · v1
∂t
µ0 j1 = ∇ × B1
∂B1
= ∇ × (v1 × B0 ) (4.31)
∂t
Here subscript ‘0’ is used for equilibrium and ‘1’ for perturbed quantities. We now intro-
duce the plasma displacement ξ (r) as dξξ /dt = v, and time integrate the equations for the
pressure and the magnetic field:
∂ 2ξ
ρ0 = j1 × B0 + j0 × B1 − ∇p1
∂t2
p1 = −ξξ · ∇p0 − Γp0 ∇ · ξ
µ0 j 1 = ∇ × B1
B1 = ∇ × (ξξ × B0 ) (4.32)
The linearized equations can then be combined into a single equation for ξ :
∂ 2ξ
ρ0 = F(ξξ ) (4.33)
∂t2
where
1
F(ξξ ) = (∇ × (∇ × (ξξ × B0 ))) × B0 + (∇ × B0 ) × (∇ × (ξξ × B0 ))
µ0
(4.34)
+∇ ξ · ∇p0 + Γp0 ∇ · ξ
84
4.5.1 The energy principle
We shall now turn to the energy principle, which is the most important method to study
the stability of the ideal MHD plasma. The method is frequently used, since it is compact
and relatively simple to handle analytically. Sometimes, the result is a so called stability
criterion, which tells us what is required from the equilibrium in order to stabilize a given
perturbation. A disadvantage is that the growth rate of the unstable modes discovered
cannot be calculated; we only know if we have stability or not. For growth rates, we have
to solve the eigenvalue equation, which is usually done numerically. Let us now see how
the energy principle works.
To simplify the notation we introduce
∂ξξ ∂ 2ξ
ξ̇ξ = , ξ̈ξ = (4.35)
∂t ∂t2
and write the equation of motion 4.33 as
We then perform a scalar multiplication of the equation with ξ̇ξ and integrate over the
plasma volume to obtain Z Z
ρ0 ξ̇ξ · ξ̈ξ dV = ξ̇ξ · F(ξξ )dV (4.37)
V V
The system may be isolated by using using the boundary condition ξ ⊥ = 0. The left hand
side of the equation can be rewritten as
Z Z Z
∂ 1 ∂ 1 ˙2
ρ0 ξ̇ξ · ξ̈ξ dV = ρ0 ξ̇ξ · ξ̇ξ dV = ρ0 ξ dV (4.38)
∂t 2 ∂t 2
V V V
It can be shown that this particular form F(ξξ ) is self-adjoint, which means that
Z Z
η · F(ξξ )dV = ξ · F(ηη )dV (4.39)
V V
for any choices of the vector fields η and ξ . We can use this property to show that
Z Z Z
∂
ξ̇ξ · F(ξξ )dV = ξ · F(ξ̇ξ )dV = ξ · F(ξξ ) dV (4.40)
∂t
V V V
The second equality follows from the fact that F is a linear function of ξ . The right hand
side of 4.37 can be rewritten as
Z Z Z Z
1 1∂
ξ̇ξ · F(ξξ )dV = ξ̇ξ · F(ξξ )dV + ξ · F(ξ̇ξ )dV = ξ · F(ξξ )dV (4.41)
2 2 ∂t
V V V V
We can now rewrite eq. 4.37, inserting eq. 4.38 and 4.41, as
Z Z
∂ 1 ˙2 1
ρ0 ξ dV − ξ · F(ξξ )dV = 0 (4.42)
∂t 2 2
V V
85
We interpret this equation as follows: The first term is the perturbed kinetic energy of the
system Z
1 ˙2
δK = ρ0 ξ dV (4.43)
2
V
and the second term is the perturbed potential energy of the system
Z
1
δW = − ξ · F(ξξ )dV (4.44)
2
V
The equation,
∂
δK + δW = 0 (4.45)
∂t
states that the total (kinetic + potential) perturbed energy is constant in time,
δK + δW = const. (4.46)
It then follows that any perturbation that decreases the potential energy (δW < 0)
produces a corresponding increase in the kinetic energy (δK > 0), which means that the
initial velocity increases, and the system is unstable. Alternatively, if all perturbations
lead to an increase in the potential energy (δW > 0), the system is stable.
Now, the idea is the following. First we define the geometry (e.g. a straight circular
cylinder). Then, we try all possible small displacements ξ of the plasma. If the displace-
ment that minimizes δW gives δW > 0 the plasma is stable, otherwise it is unstable. The
energy principle is a compact and neat way of testing for stability.
Test functions for ξ are very useful when it is possible to show that δW < 0 for a
certain perturbation. But how can we know if a sufficient number of perturbations ξ
have been investigated? There are two good approaches: either it is possible to see how
ξ affects δW in each case, or ξ must be expanded in a complete series of eigenfunctions
(usually truncated after a few terms).
In order to see what affects the stability, we write the perturbed potential energy as
follows:
Z
1
δW = − ξ · j1 × B0 + j0 × B1 − ∇p1 dV (4.47)
2
V
86
The first term in the integral can be rewritten as
Z
ξ · j1 × B0 dV
V
Z
1
= ξ· (∇ × B1 ) × B0 dV
µ0
V
Z
1
= ∇ · (B1 × B0 ) × ξ − B1 · ∇ × (ξξ × B0 ) dV
µ0
V
I Z (4.48)
1 1
= (B1 × B0 ) × ξ · n̂ dS − B1 · B1 dV
µ0 µ0
S V
I Z
1 1 2
= (B1 · ξ )(B0 · n̂) − (B1 · B0 )(ξξ · n̂) dS − B dV
µ0 µ0 1
S V
I Z
1 1 2
=− (B1 · B0 )(ξξ · n̂) dS − B dV
µ0 µ0 1
S V
We have here used the Gauss’ theorem to transform part of the the volume integral to a
surface integral over the plasma boundary. Further, we have used B0 · n̂ = 0 since the
equilibrium magnetic field is assumed to be aligned with the plasma boundary. Similarly,
we can rewrite the third term in the integral as
Z
ξ · ∇p1 dV
V
Z
= ∇ · (p1 ξ ) − p1 (∇ · ξ ) dV
V
I Z (4.49)
= p1 (ξξ · n̂) dS − p1 (∇ · ξ ) dV
S V
I Z
2
= p1 (ξξ · n̂) dS + (ξξ · ∇p0 )(∇ · ξ ) + Γp0 (∇ · ξ ) dV
S V
We now insert eq. 4.48 and eq. 4.49 in eq. 4.47, giving the following expression for the
perturbed potential energy:
1 B12
Z
δW = − j0 · (B1 × ξ ) + (ξξ · ∇p0 )(∇ · ξ ) + Γp0 (∇ · ξ )2 dV
2 µ0
V
I (4.50)
1 1
+ B1 · B0 + p1 (ξξ · n̂) dS
2 µ0
S
87
This term represents the stabilizing effect that arises from compression and bending of
magnetic field lines. The second term is a potentially destabilizing term.
Z
− j0 · (B1 × ξ )dV
V
It involves the equilibrium current density j0 . An instability that is caused by this term
being negative is referred to as a current-driven instability. The third term is also a
potentially destabilizing term.
Z
(ξξ · ∇p0 )(∇ · ξ )dV
V
It involves the equilibrium pressure gradient ∇p0 . An instability that is caused by this
term being negative is referred to as a pressure-driven instability. The fourth term is a
positive, always stabilizing term.
Z
Γp0 (∇ · ξ )2 dV
V
B12
Z
1
δWV = dV (4.54)
2 µ0
vacuum
The brackets in the expression for WS describes the discontinuous ”jump” in the quantity
inside the bracket across the boundary, calculated as the quantity just outside the bound-
ary minus the quantity just inside the boundary. In the surface integral, n̂ is the normal
unit vector in the outward direction from the plasma’s confining surface. The integra-
tion must be performed over the entire unperturbed plasma volume (F ), over the plasma
88
edge surface (S) and over the outer vacuum volume (V ) between the plasma and the sur-
rounding wall. The transformation between surface and outer volume integrals have been
performed using Gauss’ theorem, and the boundary conditions at the plasma-vacuum
boundary for magnetic fields and pressure.
4.5.3 Eigenmodes∗
Since the linearized equation of motion has constant coefficients in time, it can be writ-
ten as a complex valued eigenvalue equation. The displacement vector is obtained from
the real part of the complex valued eigenfunction ξ (r) and the corresponding complex
eigenvalue ω as
ξ (r, t) = <{ξξ (r) exp(iω)} (4.55)
Inserting this expression in eq. 4.33 gives the following form of the equation of motion
By scalar multiplication by ξ ∗ and integration over the entire plasma volume, we obtain
the potential energy change associated with ξ (r). It is better to multiply with the complex
conjugate ξ ∗ rather than ξ since ξ ∗ · ξ is real.
Z Z
− ρ0 ω ξ · ξ dV = ξ ∗ · F(ξξ )dV.
2 ∗
(4.57)
V V
The left hand side is the perturbed kinetic energy δK produced by the eigenfunction, and
the right hand side is the corresponding perturbed potential energy δW . If follows from
the self-adjoint property of F(ξξ ) that ω 2 is real-valued. If we solve for ω 2 from (4.57) we
find;
− 12 V ξ ∗ · F(ξξ )dV
R
2 δW
ω = R 1 2
≡ . (4.58)
ρ |ξξ | dV
V 2 0
K
Since K > 0, the sign of ω 2 is given by δW . If δW < 0 we have ω = ±i|γ|, where γ is
the growth rate of the mode. There will be a growing eigenmode, growing as e|γ|t , as well
a damped eigenmode e−|γ|t , which quickly becomes dominated by the growing mode. For
δW > 0 we get oscillating, undamped eigenmodes(waves) with time dependence e±i|ω|t .
89
Figure 4.6: Toroidal instability in an early fusion experiment. (Figure: EUROfusion-JET)
90
Figure 4.8: Diagram of unstable modes in a tokamak. (Figure credit: By permission of Oxford
University Press)
It was early realized that a powerful axial magnetic field and a conducting wall close
to the plasma had stabilizing effects, but it would turn out that several other instabilities
can be significant. We will now familiarize ourselves with so called macroscopic insta-
bilities which may be studied in fluid models. There are also a great number of more
or less important microscopic instabilities which are caused by, for example, the charge
distribution or non-Maxwellian velocity distributions. These instabilities must be studied
within kinetic theory.
In principle, we can say that the price paid for the high stability of the tokamak is the
high toroidal magnetic field. This causes the safety factor to assume the value q(0) ≈ 1
increasing steadily to q(a) = 3–4. Because of this, the number of resonant modes is
greatly reduced in the plasma, see figure 4.8. Since q = m/n on the resonance surfaces,
the mode (m, n) = (1, 1) is destabilized close to r = 0. This mode is called internal since
it is resonant in the plasma. Further, it is obvious that no modes (0, n); for all n or (1, n);
n > 1 can be resonant. However, modes like (2,1), (3,1), (3,2) and so on may be. These
are called kinks since the instability is kinking.
91
We will now demonstrate use of the energy principle. If b = a, i.e. the plasma extends
to the conducting wall, we have ξa = 0. This leaves us with only positive terms in the
integrand in (4.59). Thus, δW > 0 in this case, and the plasma is stable to perturbations.
If the conducting wall is further away, the plasma is still stable to all modes (m, n) where
qa > m/n. With m = 2, n = 1, 2, 3 . . ., we see that the plasma is stable if qa = q(a) > 2,
which is often the case in a tokamak.
Then, which modes are unstable? Let’s assume q(a) = 2.5. Outside the plasma edge
(but inside the conducting wall), q increases as rBϕ /(RBθ ) ∝ (Bϕ /R)[r/(1/r)] which is
∝ r2 since Bϕ ∝ 1/R in a tokamak. Therefore q may equal to 3 between the plasma and
the wall, making the perturbation resonant there and qa < m/n is a possibility. There
is nothing mysterious about this; the energy principle treats the stability of the entire
system plasma + vacuum field. Modes that are resonant between the plasma and the
wall are known as external modes. They are stabilized by, as discussed, having q(a) > 3.
Then, only modes with m > 3 are resonant, and these are generally stable at reasonable
beta values since they cause a considerable bending of the magnetic field.
But, you may ask, does the fact that one term in (4.59) is negative mean that δW < 0?
It does. The advantage with the energy principle is precisely this. We may choose any
perturbation ξ(r) when we minimize δW in search of an instability. We may, for example,
choose ξ(r) very small in the plasma volume, so that the positive contribution of the
integral may be neglected. Then, the negative term after the integration will dominate
close to the resonant surface and we will have an instability.
Figure 4.9: Stability diagram for external kinks in a tokamak. (Figure credit: By permission
of Oxford University Press)
The method for keeping q(a) > 3 is to “peak” the current density profile close to
r = 0. Figure 4.9 shows how the stability against external kinks is improved by increasing
v in the current density profile jϕ = j0 [1 − (r/a)2 ]v , which also increases q(a)/q(0). The
diagram has been calculated using a numerical code, applying the energy principle to the
assumed current profile. We see that the plasma is completely stable to external kinks
for q(a)/q(0) > 3.5, i.e. v ≥ 3.
In conclusion, the tokamak plasma is usually stable to external kinks if q(a) > 3.
This is important since unstable MHD modes grow very fast, on the time scale of one
Alfvén time (less than one µs). Recently, experiments have been carried out in tokamaks
with “reversed shear” where the q profile is permitted to have a minimum in the plasma
and small values on the edge. The significance of having the conducting wall close to
92
q
1
r
ξ δ
ξ0
rs r
Figure 4.10: Test function ξ(r) for internal (m, n)=(1,1) internal kink instability
the plasma is then increased. The reason for this is that these “advanced” tokamak
configurations allow a big part of the plasma current, called bootstrap current, to be driven
by the plasma itself (more on this in the chapter on heating). Further, this configuration
allows the reduction of the tokamak’s dimensions.
where ε = a/R. Let us investigate if δW < 0 for any perturbation when (m, n) = (1, 1).
It seems like δWmin = 0 if we use the test function ξ(r) = ξ0 , where ξ0 is constant. This
test function is not, however, a possible solution since it does not fulfil ξ(a) = 0 except
for the trivial case ξ0 = 0.
Let us instead use the following test function (see figure 4.10)
ξ0
0 ≤ r < rs − δ
ξ(r) = ξ0 (rs − r)/δ rs − δ ≤ r ≤ rs (4.61)
0 r > rs
If we use the Taylor expansion around q(r) at distance x = r − rs from the resonance
surface; q ≈ qs + q 0 x = m/n + q 0 x, we have n/m − 1/q = (nq − m)/(mq) ≈ (n2 /m2 )q 0 x.
For m = 1, the second term in the square bracket in the integral is zero and we have
Z0
ξ0 2 0 2 ξ2
δW ∝ − (q x) dx ∝ 02 q 02 δ 3 = O(δ) (4.62)
δ δ
−δ
which approaches zero when δ → 0. So, there are allowed test functions resulting in
δWmin = 0. Thus, the mode is marginally stability to this order, and stability depends
on the terms O(ε2 ). But if terms O(ε2 ) are included we find, after some calculations,
that δWmin < 0 for realistic beta values. This mode, which is independent on an ideally
conducting wall, is called an internal kink .
93
Figure 4.11: Combined stability diagram for external and internal kinks. (Figure credit: By
permission of Oxford University Press)
Experimentally, internal kinks are manifested as sawteeth in measured signals of, for
example the temperature close to r = 0. What happens is this: since the resistivity
is proportional to T 3/2 and the temperature profile decreases outwards from r = 0, the
current density tends to increase in time close to r = 0. This causes a local decrease
of q ∝ rBϕ /(RBθ ) ∝ 1/jϕ . As q approaches 1, the (1,1) mode is destabilized and the
resistivity close to r = 0 increases. This causes a decreasing jϕ and increasing q. The
cycle is then repeated and sawtooth shaped oscillations arise. In this way, q(0) is fixed to
a value close to one in a tokamak.
In general, the (1,1) kinks negatively affect the flux surfaces close to r = 0, causing
flat profiles and poor confinement locally. Resistive effects do have a limiting effect on
the instability, but the best thing would be to eliminate the mode completely.
In the complemented stability diagram in figure 4.11 the area of internal kink modes
has been added.
94
high aspect ratio tokamak, which allow as to approximate the plasma geometry with
a circular cylinder of length L = 2πR, and cylindrical co-ordinates (r, θ, ϕ), ϕ = z/R.
The eigenmodes in this geometry are Fourier modes, and the eigenvalues are the inverse
characteristic growth times of the modes, the so called growth rates γ. The ansats for the
eigenmodes are thus
The eigenvalue problem is reduced to solving the one dimensional differential equations for
the complex valued Fourier mode coefficients Q(r) with appropriate boundary conditions.
If a helical eigenmode has the same pitch as the equilibrium magnetic field on some
flux surface, the perturbation is said to be resonant. In a toroidal geometry, this resonance
condition can be expressed using the safety factor as q = m/n. Magnetic field surfaces
where this condition can be fulfilled are called resonant surfaces. It is the region around
the resonant surface that plasma resistivity matters. Outside this reqion, the ideal MHD
model works well. The procedure to calculate the growth rate of the tearing mode is
a type of boundary layer analysis. An internal boundary layer is assumed around the
resonant radius. Outside this boundary layer the resistivity can be neglected, and the
ideal MHD model is used. Further, resistive instabilities grow slowly, so the inertial term
can be neglected in the other, using a perturbed force equilibrium equation. The starting
point of the analysis is thus the force equilibrium equation:
j × B − ∇p = 0 (4.64)
The aim is to obtain a differential equation for the radial magnetic field perturbation
B1r . We start by taking the rotation of the force balance equation, thus eliminating ∇p:
∇ × (j × B) = 0 (4.65)
Using the low-aspect ratio, tokamak approximation, if follows that B · ∇jϕ j · ∇Bϕ ,
and, we can write the axial component of the equation as
B · ∇jϕ = 0 (4.68)
95
term approaches zero and (4.70) is no longer valid. We introduce an inner boundary layer
at the resonance surface. In the boundary layer, rs − δ < r < rs + δ, resistivity and inertia
are important, and we have to use the resistive MHD equations, with the inertial term in
the equation of motion included. First we use Ohm’s, Faraday’s, and Ampere’s laws
ηj =E + v × B (4.71)
∂B
− =∇ × E (4.72)
∂t
µ0 j =∇ × B
to write an equation for the time rate of change of the magnetic field, as follows
∂B1 η
− + ∇ × (v1 × B0 ) = ∇ × (∇ × B1 ). (4.73)
∂t µ0
We can rewrite the second term as
∇ × (v1 × B0 ) = (B0 · ∇)v1 − B0 )(∇ · v1 ) (4.74)
and assume incompressible perturbations (∇ · v1 ). We then obtain an equation for the
radial component of the perturbed magnetic field
η
γBr1 − B0 · ∇v1r = ∆Br1 . (4.75)
µ0
We need another equation for the velocity perturbation v1r . The equation of motion,
including the inertial term, but omitting the term 2:nd order in the velocity is used
∂v
ρ = j × B − ∇p. (4.76)
∂t
Taking the rotation gives the following equation
ργ∇ × v1 = ∇ × (j × B) (4.77)
which can be expanded to yield an equation for v1r . We combine the two differential
equations, eqs. 4.75 and 4.77 to obtain a differential equation for B1r in the boundary
layer. The growth rate γ is obtained by matching the solutions for B1r in the outer
ideal MHD region to the solution in the inner resistive boundary layer at the boundaries
between the two regions. The matching procedure introduces the matching parameter ∆0 ,
which has given name to this analysis method which is sometimes called the ”delta-prime
method”. "
0
out #r=rs +δ " 0 in #r=rs +δ
B B1r
∆0 = 1r
= (4.78)
B1r B1r
r=rs −δ r=rs −δ
The bracket indicates the difference of the values taken at the boundaries located at
r = rs + δ and r = rs − δ. We now introduce the characteristic time for magnetic field
diffusion in the plasma, the resistive time constant
µ 0 a2
τR = (4.79)
η
and the inertial time scale, the Alfvén transit time
√
a µ0 ρ
τA = (4.80)
Bϕ
96
Figure 4.12: Cross section of a) normal flux surfaces b) flux surfaces with m=1,2,3 islands
a aq 0 −3/5 −2/5
γ = 0.55 n τR τA (a∆0 )4/5 (4.81)
R q
where n is the toroidal mode number of the eigenmode. The quantities in the expression
for γ are evaluated at the position of the resonant layer r = rs . The characteristic growth
time for the tearing mode is intermediate between the resistive and inertial time scales.
Let us illustrate with a numerical example: Consider a tokamak with magnetic field
Bϕ = 5 T, deuterium plasma ion density n = 1020 m−3 , and plasma temperature T = 1
keV, in a torus with minor radius a=0.5 m. The resistivity η = 2.5 × 10−8 Ωm is estimated
from Spitzer’s formula with these parameters, giving a resistive time constant of τR = 13
s, in the order of seconds. The Alfvén transit time τA = 0.065 µs, which characterize the
growth of ideal MHD modes, is less than microseconds. The tearing mode growth time for
3/5 2/5
these parameter is of the order of milliseconds γ −1 ∼ τR τA ∼ 6 ms. Resistive tearing
modes thus grow considerably slower than ideal MHD modes at high plasma temperature.
For these “pure” tearing modes, γ increases as η 3/5 with plasma resistivity. The growth
rate is approximately inversely proportional to plasma temperature since η ∝ T −3/2 and
γ ∝ η 3/5 ∝ T −9/10 . When realistic beta values are used in the analysis, the dependence of
the growth rate on resistivity can instead become the so called resistive g-mode scaling
η 1/3 .
97
Figure 4.13: Transformation of the toroidal geometry to cartesian coordinates.
to the resonant surface will then vary as in figure 4.13d. At the resonant surface r = rs
we have B0 = 0. (Alternatively B0 can be seen as the poloidal equilibrium field after
substraction of a field with constant pitch equal to the equilibrium field on the resonant
surface.) Assume that a resonant perturbation with a radial magnetic field component
arises. For simplicity, we assume the perturbation to be sinusoidal in poloidal and toroidal
directions. Some distance away from the resonant surface the effect of the perturbation
will be small, but close to the surface it will be comparable to the equilibrium magnetic
field in this projection. The result is that local flux surfaces, i.e. magnetic islands are
created as in figure 4.13e. The boundary between local flux surfaces in a magnetic is-
land and the surrounding global flux surfaces is called the separatrix . The width of the
separatrix at the widest point on the island is called the island width, denoted W .
An expression for the island width W for a given amplitude of a magnetic perturbation
can be obtained in the following way. Introduce cartesian coordinates (x, y, z) and let
the x direction correspond to radial direction, with x = 0 at r = rs , and B0 is in y
direction (figure 4.13 e). Assume a resonant magnetic field perturbation B1 arises. The
perturbation will be periodic in y direction, and its x component can be written
where k = 2π/λ = m/rs is the perturbations wave vector in poloidal direction, λ is the
wave length and m is the poloidal mode number. Close to the resonant surface, the y
98
component of B0 can be approximated by Taylor expansion around x = 0:
B0y ≈ xB00 (4.83)
where
dB0y
B00 = . (4.84)
dx x=0
In this 2-D geometry it is helpful to introduce a flux function to describe B. The flux
function Ψ is chosen such as
1. Ψ is constant on magnetic flux surfaces: ∇Ψ · B = 0
2. B fulfils the condition: ∇ · B = 0.
These two condition are fulfilled by the following definition for Ψ:
∂Ψ
Bx = −
∂y
∂Ψ
By = (4.85)
∂x
By introducing eq. 5.99 in eq. 5.102 we obtain the perturbation part of the flux function;
B1
Ψ1 = cos ky (4.86)
k
and by inserting eq. 5.100 we get the equilibrium part;
x2 0
Ψ0 = B. (4.87)
2 0
The total flux function is
Ψ = Ψ0 + Ψ1 . (4.88)
The flux function is also constant on the separatrix, i.e. the line separating the area with
local flux surfaces (the magnetic island) from the area with global flux surfaces (the rest
of the plasma). Let the flux function assume the value Ψ = Ψs on the separatrix. An
expression for the island width can be reached as follows: The separatrix passes x = 0 at
ky = p 2π, p = ±0, 1, 2, . . .. Inserting this in eq. 5.105 gives
B1
Ψs = . (4.89)
k
The separatrix also defines the island width. One special case is W = 2x for ky =
π + p 2π, p = ±0, 1, 2, . . ., which are the y values when the separatrix passes the widest
part of the island. Is we insert this in eq. 5.105 we get a second expression for Ψs :
2
1 W B1
Ψs = B00 − . (4.90)
2 2 k
Combining eq. 5.106 and 5.107 gives us the following expression for the island width W :
s
B1
W =4 . (4.91)
kB00
In cylinder geometry, it can be shown that the expression for the island width W has a
similar form (using the safety factor q):
s
Br1 q
W = 4r . (4.92)
mBθ0 rq 0
99
Figure 4.14: Positive curvature (the inside of the torus) is, on average, stabilizing. (Figure
credit: By permission of Oxford University Press)
100
An apt name for these instabilities is ballooning modes. Since these strongly limit the
possible beta values in tokamak, different ways of handling them have been explored. One
interesting possibility has turned out to be the positive effect from the inside of the torus
when the plasma is given a bean like shape. This allows the so called second stability
region to be reached.
4.6.7 Disruptions
An Achilles’ heel of the tokamak is disruptions, and limit operational parameters of the
tokamak such as the maximum current and density. Disruptions are fast events when
plasma current and confinement drastically declines. The sudden termination of the
plasma cause high heat loads and strong mechanical strains. Future tokamak experiments
are designed to withstand a number of smaller disruptions, but in a reactor, it would be
preferred to get rid of them completely. Not all tokamaks suffer from disruptions. At low
aspect ratio, spherical tokamaks, disruptions do not occur.
The cause of disruptions is only partly known. We will now discuss two important
factors; too low qa value and too high density. Tearing modes with high amplitude have
proved to be connected to these cases. When the qa value is too low (≤ 2), the current
density at the edge tends to be steep, which causes strong tearing modes. Thus, it is
important not to have too high plasma current for a given toroidal field.
Disruptions at high density is usually caused by an increase in radiation due to edge
impurities. When the plasma edge thus is cooled, the current density profile it contracted
resulting in a destabilizing gradient. An empirical law for the upper density limit has
been found by Greenwald; n < 1 × 1020 Ip (MA)/(πa2 ). The relation has turned out to
be valid for many experiments. It gives the maximum allowed density 6 × 1020 m−3 in a
reactor with a = 1 m and Ip = 20 MA.
101
Figure 4.15: Comparison between experimental results for beta and equation (4.96). (Figure
credit: By permission of Oxford University Press)
Figure 4.15 shows comparisons with experiments where the Troyon factor 2.8 above
has been replaced by li (the internal inductance of the plasma), which is a measure of
the current profile’s influence on the βϕ limit. The inductance increases when current is
concentrated near the magnetic axis r = 0. The parameter βN denotes the “normalized
beta value”; βN = β(%)/[Ip (MA)/(aBϕ )], where β is the experimental value.
102
Figure 4.16: Comparison of L-mode and H-mode temperature profiles showing the temperature
”pedestal” at the edge of the plasma. (Figure credit: EUROfusion-JET)
lead to ELM:s.
ELMs appear in different forms. In their weakest form they occur as frequent small
pulses, producing time traces that look like grass, giving them the name ”grassy” ELMs.
With higher power heating they become more separated in time and larger in size. When
they have become ”giant” ELMs they present a serious problem, the associated sudden
energy release produces a high heat loading on material surfaces, which leads to reactor
technological constraints.
103
with the field lines of neighbouring magnetic surfaces. This gives rise to magnetic short-
circuiting. Since now fast energy exchange can also occur transversally to the surfaces,
the plasma temperature and plasma pressure drop right across the island. This restricts
the plasma pressure attainable.
The formation of such magnetic islands can be hampered, however, by accurately
beaming microwaves into the centre of an emerging island. This generates a local electric
current that makes the island vanish. The magnetic field perturbation is suppressed
and the plasma pressure can rise again. It is presently possible to eradicate an island
completely. One may stabilise a particularly disturbing tearing mode that can lead to
disruption of a discharge; after improvement of the beaming method this could be achieved
with quite low microwave power. For stabilisation, less than ten per cent of the total
heating power applied is sufficient when precisely beamed to the right spot.
104
Chapter 5
TRANSPORT
105
n
r no
r
ro
t
τ
x
Δx
plasma and in neutral gas is that the particles in the plasma do not move in straight
paths between the collisions, but in curved paths given by the gyration in the magnetic
field. Another important difference is that the collisions are not close collisions, but the
result from interactions by many particles over long distances. Despite these differences,
estimations of the order of magnitude of diffusion coefficients can be made with models of
the type used for neutral gases, so called random walk models, which is the starting point
for this description of classical transport theory. The purpose of these estimations is to, in
the easiest possible way, illustrate the principle of the particle transport mechanism. The
analysis can in principle be generalized to include the more important heat diffusivity.
At the end of this section, we will summarize the relations between particle and heat
diffusivity which are useful for an analysis of the plasma energy balance.
3. particles move in positive and negative direction with equal probability after a col-
lision
106
+ ρ
-
Δx
Z0 !
1 dn
Γ+ = n0 + x dx
2τ dr r0
−∆x
(∆x)2 dn
Γ = Γ+ − Γ− = −
2τ dr r0
(∆x)2
D= (5.2)
2τ
This important result, which can be used to estimate the diffusion coefficient, tells us that
the coefficient is proportional to the step length squared divided by the collision time.
This is the diffusion we call random walk process.
107
+
z
j
V Δ
p
B
collisions between particles of different charge (see figure 5.4). The collision time used for
the estimate is thus the one for Coulomb collisions between electrons and ions.
Collisions between particles of different charge cause a movement of the guiding centre
of the particles of the same order of magnitude as the Larmor radius
∆x ≈ rL (5.3)
Inserting this in the expression for the diffusion coefficient in the random walk model
gives
r2
D≈ L (5.4)
2τ
Including the Larmor radius and collision time τe of the electrons, the expression for the
particle diffusion coefficient D becomes
me T
D≈ (5.5)
eB 2 τe
If we insert the Larmor radius and collision time τi of the ions, we find that the
diffusion coefficients for ions and electrons are equal, i.e. De = Di = D. This means that
the classical diffusion is naturally ambipolar , i.e. it conserves the charge neutrality of the
plasma.
108
n
a
a
r
expression for the radial fluid velocity perpendicular to the magnetic field:
1 η
v⊥ = 2
E × B − 2 ∇p (5.6)
B B
The pressure can be written as p = 2neT . Now, assume that E = 0 and that T = const.
The radial particle flow is Γ = nv⊥ , or
2neT
Γ = −η ∇n (5.7)
B2
Comparing the above expression with Fick’s law gives an expression for the diffusion
coefficient D:
2neT
D=η 2 (5.8)
B
Inserting the expression for plasma resistivity η gives the diffusion coefficient as
2me T
D= (5.9)
eB 2 τ
It is obvious that both resistive MHD and random walk gives similar expressions for the
classical diffusion coefficient in a straight cylinder.
109
The constant τp in known as the particle confinement time. In this simple model, τp is
the time constant for the decrease in density at a given radius due to the radial particle
flow. The diffusion equation in cylindrical coordinates is of the form of a Bessel equation:
1 ∂ ∂n n
r + =0 (5.14)
r ∂r ∂r τp D
The general solution is the zero’th order Bessel function J0 (r). The condition n(a) = 0
at the plasma edge r = a finally gives
r t
n(r, t) = n0 J0 2.405 exp − (5.15)
a τp
where n0 in the plasma density at the cylinder axis r = 0 at time t = 0 and x = 2.405 is
the first zero of J0 (x). Inserting the solution gives the relation between τp And D:
a2 a2
τp = ≈ 0.17 (5.16)
2.4052 D D
This is a very useful expression for a good approximation of the confinement time when
the diffusion coefficient is given.
χe ≈ D (5.22)
r
mi
χi ≈ D (5.23)
me
The energy confinement times τE,e and τE,i are given from the heat diffusivities χe andχi .
110
5.1.6 Energy confinement time for ohmically heated plasma∗
Consider a tokamak or pinch in which the toroidal current provides both heating and
confinement. The energy confinement time τE in stationary state (Pin = Ploss ) is
3
2
< 2neT >
τE = Vp . (5.24)
Pin
where Vp is the toroidal plasma volume, approximated as Vp = 2π 2 Ra2 . The ohmic heating
power is given by the product of the plasma current Ip and the loop voltage Uloop .
Poloidal beta βθ is exactly equal to unity for a z-pinch. The z-pinch, however, is unstable.
In the more attractive reversed-field pinch, βθ ≈ 0.1. In tokamaks, βθ can reach values of 1-
3 depending on the shape of the plasma cross section. Use Ampere’s law 2πaBθ (a) = µ0 Ip
and rewrite the expression for the energy confinement time as
3 Ip
τE = µ0 Rβθ , (5.27)
8 Uloop
3βθ µ0 a2 3βθ
τE = = τR (5.29)
16fG η 16fG
where τR = µ0 a2 /η is the plasma resistive time, characterizing the magnetic field diffusion
in the plasma. If we introduce the Spitzer resistivity, eq. 2.76 and combine equations we
get
2.8 × 1014 2 2.5
τE = a βθ (Ip /N )1.5 Ip1.5 (5.30)
fG Zeff ln Λ
The plasma density has been included in the so called line density N = nπa2 . This is the
highest possible energy confinement time for an ohmically heated plasma. The Coulomb
logarithm ln Λ can be set to ≈ 15 for a fusion plasma. With parameters fG = 5, a = 1 m,
Zeff = 1, βθ = 1, n = 1.0×1021 m−3 and Ip = 10 MA we find τE = 21 s. If experiments had
given this value, fusion energy would already be reality. Unfortunately, real life tokamaks
have energy confinement times one or two orders of magnitude smaller. One of the most
important problems in fusion research is to try to explain and diminish this difference.
111
z
φ
Bφ
Bθ +
Bθ B
θ
Bφ Bφ
Bθ -I
-φ
+
112
Driftytor
z Projektion av driftbanor
på vertikalplan
φ
VII Ω
r
B
θ
Rc R
d d
VII
Fluxyta
the torus, since the vertical drift changes direction relative to the minor axis making the
average drift motion away from the minor axis small on average. This means that the
particle confinement depends on the particles ability to move freely along the field lines
around the tokamak. Particle collisions, which can be regarded as a friction force slowing
down the motion along the field line, cause increased transport, as we will show in the
next section.
In a tokamak with circular magnetic flux surfaces the projection of the guiding-centre
drift orbits on a vertical plane will be circles with the centre displaced somewhat in
horizontal direction relative the flux surface (figure 5.8). The displacement will be inwards
or outwards depending on the direction of the particle’s velocity along the field line.
The size of the displacement of the guiding-centre orbit relative the flux surface can
be estimated in the following way: Consider a magnetic flux surface with a circular cross
section, its centre at R = Rc and a radius r = r0 . Introduce polar coordinates (r, θ) with
origin at R = Rc to describe the projection of a guiding-centre orbit on a vertical plane
through the torus’ symmetry axis. The guiding-centre drift velocity is split into two parts;
motion along the field line and guiding-centre drift across the field: v = vk + vd .
The parallel velocity has one component in toroidal direction and one in the vertical
plane. When the field line is turning around the minor axis due to the rotation transform,
the velocity component in the vertical plane will be in the poloidal direction; the parallel
velocity vk corresponds to a poloidal rotation with angular frequency Ω;
Bθ
vk,θ = vk = rΩ (5.31)
B
where Bθ in the poloidal component of the magnetic field and B is the absolute value of
the magnetic field. The perpendicular velocity component is caused by the guiding-centre
drift in vertical direction vd = vd Ẑ, given by eq. 2.20:
vT2
vd ≈ (5.32)
ωc R
where R is the major radius and vT is the thermal velocity. The cyclotron angular
frequency have different sign depending on the charge of the particle (positive ion or
negative electron). The perpendicular angular velocity is divided into polar components:
vd,r = vd sin θ
vd,θ = vd cos θ (5.33)
In the vertical plane, the particle position is described by radius r as a function of the
poloidal angle; rp = rp (θp ), where θp = θp (t). The expressions for the radial and poloidal
113
velocities are:
drp
= vd,r
dt
dθp
rp = vk,θ + vd,θ (5.34)
dt
drp
= vd sin θp
dt
dθp vd
= Ω + cos θp (5.35)
dt rp
r0 B vT2
d≈ (5.41)
vk Bθ ωc R
The sign of d depend on the sign of ωc and vk , i.e. the charge of the particle and the di-
rection of its motion parallel to the magnetic field. The expression for d may be simplified
by assuming vk ≈ v⊥ ≈ vT and introducing the Larmor radius rL and the safety factor q
114
z
Magn. fluxytor
<
θ
<
z
the particles, at each collision, can move one step length, corresponding to the distance
between the flux surfaces, which is of the same order of magnitude as the displacement d
of the guiding-centre orbit relative the flux surface (figure 5.9).
Inserting the distance d in the random walk model expression for the particle diffusion
coefficient gives the diffusion coefficient for this transport process.
d2
DP S = (5.43)
2τ
The expression for the distance d is
DP S ≈ q 2 DC (5.44)
Fz = −qvθ Br (5.45)
Use ∇ · B = 0 to rewrite as
1 ∂Bz
Fz = qvθ rL (5.46)
2 ∂z
115
Fm Fm
θ
B min z
B
Figure 5.11: Schematic drawing of trapping of a particle between two magnetic mirror
fields
1 mv 2
qvθ rL = − ⊥ = −µ (5.47)
2 2B
where the magnetic moment µ has been introduced. Finally, we have the following ex-
pression for the magnetic mirror force Fz
∂Bz
Fz = −µ . (5.48)
∂z
A particle gyrating with a given velocity v⊥ in a homogeneous magnetic field B will have
a constant magnetic moment µ. Since the magnetic field varies very little over space, the
motion will be almost periodic (gyration around a magnetic field line and a slow drift).
If the variation of the magnetic field is negligible over a Larmor radius, the magnetic
moment will remain constant. (The magnetic moment is a so called adiabatic invariant.)
This can be used to construct a condition for reflection of a charged particle against a
magnetic mirror field. Consider a particle with velocity vk0 and v⊥0 in an area with weak
magnetic field Bmin surrounded by two magnetic mirror fields (figure 5.11). The magnetic
momentum µ is assumed to be constant when the particle is moving along the magnetic
field line:
2
v⊥ v2
= ⊥0 . (5.49)
B Bmin
Further, the kinetic energy of the particle is constant, i.e.:
1 1
m(vk2 + v⊥
2 2
) = m(vk0 2
+ v⊥0 ). (5.50)
2 2
At the reflection point, vk = 0 which gives a condition on the magnetic field strength:
2
B 2
v⊥ vk0
= 2 =1+ 2 . (5.51)
Bmin v⊥0 v⊥0
To describe trapping of particles in the outer part of a tokamak, we assume the magnetic
field strength to be equal to the magnitude of the toroidal magnetic field, which is the
dominating magnetic field component in a tokamak (figure 5.12). The toroidal field is
approximately given by the expression
R0 1
B ≈ B0 = B0 (5.52)
R 1 + ε cos θ
where the inverse aspect ratio ε has been introduced:
r
ε= . (5.53)
R
116
Δ
z F= μ B
φ II
B max
B min
B r θ VII o , VI o
I
Ro R
Δ
F= μ II
B
The parameter will be of importance for the neoclassical transport of trapped particles
in the tokamak. The variation of the magnetic field strength along field lines decides
whether particles can be trapped:
Bmax
≈ 1 + 2ε. (5.54)
Bmin
The condition for trapping of particles in the outer part of the tokamak is given by
vk0 √
< 2ε. (5.55)
v⊥0
It can be shown that the ratio of trapped particles is
117
z
φ
"Banana orbit"
r θ
θb R
Δr=Wb
Figure 5.13: Projection of a guiding-centre orbit in the vertical plane - a banana orbit
118
The particle orbit is given by the differential equation
drp vd θp
≈ s 2 . (5.68)
dθp θb ωb θp
1−
θb
For a particle orbit coinciding with the magnetic flux surface r = r0 at the turning point
θ = θb the expression for the guiding-centre orbit is
2 !
θp
(rp − r0 )2 = wb2 1 − (5.69)
θb
where
vd θb
wb = . (5.70)
ωb
The expression describes a “banana shaped” drift orbit with half width wb in the midplane
(figure 5.13). Introducing the approximations v⊥0 ≈ vT and vk0 /v⊥0 ≈ ε1/2 , then gives a
simple expression for the banana-orbit half-width:
wb2
DB = ftrap . (5.72)
τeff
Here, ftrap is the ratio of trapped particles and τeff is an effective collision time corre-
sponding to the particle confinement time in the magnetic mirror field. It can be shown
that
τeff ≈ ε τ (5.73)
119
z Magn. flux surfaces
φ
Drift orbit
r
Ro R
Wb
Figure 5.14: Somewhat simplified illustration of the collisions causing banana orbit diffu-
sion.
D 2
q Dc
D PS
D
c
2
/
q
2
3
ε-
rL2
Dc τ
B
D
Classical
3/
2 Vth Vth n
ε qR qR
where τ is the usual confinement time (defined as the time required for 90◦ deflection
a test particle’s velocity due to the total effect of all collisions). The banana diffusion
coefficient is finally
DB = ε−3/2 q 2 DC (5.74)
where DC is the classical diffusion coefficient.
120
z
φ Δ
B
B + + +
+ Vion +
B + + V--I =
1
2 Ex B
B
θ
r
R
j ll j ll
Eθ B E θ j ll
- Velec. -
- - -
- -
Pfirsch-Schlüter requires that the particles collide during a poloidal rotation along the
field line. This condition can be written as
vT
ν> (5.77)
qR
for Pfirsch-Schlüter diffusion. The neoclassical diffusion coefficient as a function of the
collision frequency is shown in figure 5.15. There is an area between banana orbit diffusion
and Pfirsch-Schlüter diffusion which is usually called the plateau regime .
121
Radiation
3 from edge
Radial
2 transport
Ohmic
1 heating
ELECTRONS Radiation
Additional Energy
1 heating exchange
Charge
IONS exchange
2 Radial
transport
Charge exchange
3 from edge
60
τE (Neoclassical)
40
τE
(ms)
20
τE (Experimental)
0
0 1 2 3 4 5 6
20
ne (10 m-3 )
Figure 5.18: Comparison between experimental and neoclassical confinement time in the
Alcator tokamak at MIT (1978).
122
Table 5.1: Scaling laws for τE in an ohmically heated tokamak
123
Eθ
Eθ
z Vr = B
θ B
Vr
Figure 5.19: Cylinder model for calculation of fluctuation-driven radial particle flux.
some kind of plasma instability and is usually observed as, for example, fluctuations in
local plasma density. These fluctuations are called plasma turbulence. The fluctuation
driven transport is usually much larger than the neoclassical. In fusion experiments, it is
desirable to suppress the fluctuations in order to minimize anomalous transport.
The anomalous transport processes can roughly be divided into two groups. Transport
driven by electrostatic fluctuations (1) and transport driven by electromagnetic fluctua-
tions (2). The electrostatic fluctuations can give a moderate increase of transport. This
kind is present in tokamaks and is thought to cause the anomalous transport of electron
energy which can be observed experimentally. If magnetic fluctuations are present they
can cause the magnetic field to become stochastic (this is discussed later), and the result-
ing energy flow along magnetic field lines is one kind of electromagnetic turbulence driven
transport.
124
Δ
y n0 n
+
2π
E1 E1 λ= ky
V1= B -
B
x
z
or, alternatively
1
Γ= < n1 E1 > . (5.83)
B0
Transport will only occur if the density fluctuation and the velocity fluctuation are cor-
related and “in phase”, i.e. the time average of the fluctuations must not be zero.
125
The continuity equation provides a correlation between v1 and n1 (after linearisation
where terms of 2nd order are neglected).
dn0
v1 − iωn1 = 0. (5.88)
dx
The phase velocity of the wave can be written as
ω Te dn0
vDe = =− . (5.89)
k B0 n0 dx
The phase velocity for drift waves vDe is equal to the so called diamagnetic drift velocity
for electrons. In this approximative model, the drift wave is oscillating with a constant
amplitude. The density and velocity perturbations are 90◦ out of phase, which means no
radial particle flux. In a realistic plasma model with finite ion mass, the inertia of the ions
causes a phase difference between density and velocity perturbations. E.g., the velocity
perturbation v1 is in the outwards direction when the density perturbation is already
displaced outwards. Thus, a so called drift instability arises. Regard the perturbation
as a plane wave with increasing amplitude. The phase difference between n1 and ϕ1 can
be expressed as a function of the growth rate γ in the following way: All perturbation
quantities are assumed to vary in time and space as
126
n0
n
A diffusion coefficient describing drift wave transport can finally be obtained using Fick’s
law 2
1 2 n1
DD ≈ γδn . (5.100)
2 n0
A very approximate, but simple, estimate of DD is reached if we assume that the amplitude
of the density perturbation is of the same order of magnitude as the gradient of the
equilibrium density
n0
kn1 ≈ . (5.101)
δn
This should be a upper limit for the perturbation amplitude, which in this case cancel
the gradient of the equilibrium density, which is the driving force for the instability, see
figure 5.21. The drift wave diffusion coefficient is finally
γ
DD ≈ . (5.102)
k2
In the literature, this approximation is called the mixing length approximation of the
turbulence driven diffusion coefficient.
127
Figure 5.22: Numerical simulation of gradually increasing magnetic perturbation ampli-
tude.
perturbations at the two surfaces (figure 5.22b), but this does not influence the rest of the
plasma. When the amplitude of the perturbations increases an interesting phenomenon
occurs. Magnetic field lines leave the flux surfaces and are visible as individual points at a
certain distance away from these (figure 5.22c). When the amplitude of the perturbation
is sufficiently large, the magnetic field will look like in figure 5.22d. The magnetic field
lines wander randomly between the resonant surfaces. This is seen in the graph as an
even pattern of points between the surfaces. The magnetic field in this area is stochastic.
Since transport along field lines is much larger than across the magnetic field, plasma
confinement is heavily affected when the magnetic field becomes stochastic. An impor-
tant criterion telling us whether a magnetic field is stochastic or not is if overlapping of
magnetic islands occurs, i. e. if the radial extent (island width W ) of islands on adjacent
magnetic flux surfaces is larger than the radial distance between flux surfaces.
128
distance δr radially while followed a distance lc in the toroidal direction:
δr δBr
= . (5.103)
lc B
The correlation length in this case has the same meaning as the collision time in the
random walk model of the particle diffusion coefficient:
(δr)2
DM ≈ . (5.104)
2lc
It is obvious that the dimension of DM is meter and it may thus only be regarded as
a diffusion coefficient in a mathematical sense. For an order of magnitude estimate of
DM , δBr can be put equal to the radial amplitude of the magnetic field perturbation
δBr ≈ B1r . The expression for the magnetic field line diffusion coefficient is
2
B1r
DM ≈ lc . (5.105)
B0
Consider the energy transport for electrons in a stochastic magnetic field. Depending on
the collision frequency of the electrons ν there are two cases. These are characterized by
the relation between the electron mean free path length along the field lines λ and the
correlation length of the magnetic field lc :
In the collision free area (λ > lc ) the electrons are assumed to move radially the same
distance δr as magnetic field lines when they move along the field lines. The radial heat
diffusivity perpendicular to the magnetic field χe⊥ is estimated in the same way as the
particle diffusion using the random walk model:
(δr)2
χe⊥ ≈ . (5.106)
2τ
In this collision free model, the collision time is replaced by the time interval τ during
which electrons move the distance lc ;
lc
τ= (5.107)
vT
where vT is the thermal velocity of the electrons, and
2
B1r
χe⊥ ≈ DM vT = lc vT . (5.108)
B0
129
The heat diffusivity along field lines can also be estimated using the random walk model:
λ2
χek ≈ . (5.110)
2τ
Insertion in the random walk estimate of the radial heat diffusivity gives, using the ap-
proximation δBr ≈ B1r :
2 2
λ2 δBr
B1r
χe⊥ ≈ = χek . (5.111)
2τ B0 B0
In both cases (collision free and collision dominated) the heat diffusivity is proportional
to the magnetic fluctuation amplitude squared. (Compare to the diffusion coefficient for
electrostatic fluctuations) When the parallel heat transport χek is very high, the radial heat
transport across the magnetic field can be significant already at low levels of fluctuation.
130
Chapter 6
RADIATION, BOUNDARY
LAYER, HEATING
One of the most important processes causing energy losses in a fusion plasma is emission
of electromagnetic radiation. In the central plasma, which is fully ionized, bremsstrahlung
is dominating. In a magnetized plasma there is also cyclotron radiation caused by the
gyrating motion of the electrons in the magnetic field. This kind of radiation is, however,
mainly reabsorbed in the plasma.
In the boundary layer of the plasma there are ions which are not fully ionized. These
cause line radiation in this region. This is also the case for the central plasma if the
temperature is low. Line radiation in the boundary layer is often desirable, since energetic
ions and atoms from the central plasma may transform their kinetic energy to radiation.
The advantage is that they hit the wall with lower energy, which limits the problems
associated with plasma wall interaction. By plasma wall interaction, we mean a number of
processes taking place at the plasma boundary. This often involves energetic ions or atoms
from the plasma hitting the surrounding walls. It can be sputtering or melting of the wall
material, which then cause atoms from the wall to enter the plasma as contaminations.
These contaminations increase the radiation losses from the plasma, and for this reason, it
is desirable to avoid intense plasma wall interactions. Devices designed to control plasma
wall interaction in fusion reactors are so called limiters or magnetic diverters, discussed
later in this chapter.
To reach the temperatures required for a fusion reactor energy must be supplied to the
plasma somehow. There are three main methods used: 1) Resistive, or ohmic heating. 2)
Injection of high energy beams of neutral atoms. 3) Heating with electromagnetic waves
in the radio frequency range.
6.1.1 Bremsstrahlung
Coulomb collisions between electrons and ions cause acceleration of charged particles
which leads to so called bremsstrahlung. Since the electrons are much lighter than the
131
d
-e, v d
d
+Ze
Figure 6.1: Schematic illustrating the approximation used in the Bremsstrahlung power
estimate.
ions, the electrons will be more accelerated and thus make the dominating contribution
to the bremsstrahlung emission from the plasma.
A simple semi-classical estimate of the emitted bremsstrahlung power can be made
in the following way: The bremsstrahlung power from an electron in a Coulomb collision
is given by a general expression from electromagnetics for the electromagnetic radiated
power emitted from an electron with charge −e and acceleration a.
e2 a2
P = . (6.1)
6πε0 c3
Here c is the velocity of light. The acceleration a is given by the expression for the
Coulomb force
F Ze2
a= = (6.2)
me 4πε0 me r2
where r is the distance between the charged particles. The radiated power increases
rapidly with decreasing distance, P ∼ r−4 , which means that so called close collisions
are dominating. (See the section on Coulomb collisions). Quantum mechanics give the
following estimate of the smallest possible distance d = rmin (see Heisenberg’s uncertainty
relation);
h
d= (6.3)
2πme v
where h is Plank’s constant and v is the electrons velocity. An estimate of the emitted
energy δE from a Coulomb collision can be made using the following simplified model of
the collision (see figure 6.1). Assume that the electron passes on a distance d away from
the ion and that the interaction occurs during a time t corresponding to the time needed
for the electron to travel a distance 2d along it’s path, or
2d
t= . (6.4)
v
Emitted energy per collision is given by eq. 6.1, 6.2 and 6.4 as
Z 2 e6
δE = P · t = . (6.5)
6(2πε0 cd)3 m2e v
σ = πd2 (6.6)
132
2
Pbr
Pα
1
0
0 5 10 15
T(keV)
Figure 6.2: The relation between bremsstrahlung power and alpha particle heating.
ν = ni σv. (6.7)
The total emitted bremsstrahlung power density for electrons of density ne and velocity
v is finally given by eq. 6.5–6.7
ne ni Z 2 e 6 v −38 2
p
PB = ne νδE = = 1.2 × 10 n n
e i Z Te [W/m3 ]. (6.8)
24πε30 c3 hme
In a multi-species plasma, we may use the definition of the effective charge, eq. 2.72, and
write p
PB = 1.7 × 10−38 n2e Zeff Te [W/m3 ]. (6.10)
The temperature is inserted in unit electron volt (or simply volt).
The bremsstrahlung loss provides a criterion for the lowest possible electron temper-
ature in a fusion reactor (approx. 5 keV). Since the heating power from alpha particles
increase faster with higher temperature than the bremsstrahlung loss, a positive net effect
is possible when the electron temperature is higher than this value (figure 6.2).
133
Inserting eq. 6.10 in eq. 6.1 gives the following expression for emitted power per unit
volume for a given electron density ne :
e5 B 2 Te
PC ≈ n e = 6.1 × 10−20 ne B 2 Te [W/m3 ]. (6.12)
3πε0 c3 m3e
The cyclotron radiation frequency spectrum shows that the radiation is emitted at
certain frequencies corresponding to the electron cyclotron frequency and its higher har-
monics.
ωce
f =n (6.13)
2π
where n = 1, 2, 3 . . .. In a fusion reactor plasma it is likely that re-absorption of radiation
at these frequencies is very large and the power loss from the plasma surface will be
small. For this reason we often neglect this power loss when analysing the energy balance
of the fusion plasma. Electron cyclotron radiation, however, has great importance as a
diagnostic used to measure the electron temperature in the plasma.
134
Charge neutrality, which is a fundamental property of a plasma, then gives the electron
density as X
ne = Z i ni . (6.16)
i
The distribution of density in the different charge states is mainly dependent on the
plasma electron temperature.
The distribution of charge state densities can be estimated with a model for the ion-
ization and recombination in the plasma in the following manner: The density in a given
charge state z, for a given ion type i, is described by the following system of ionization
balance equations
dni,z
= ni,z+1 ne αi,z+1 − ni,z ne αi,z + ni,z−1 ne Si,z−1 − ni,z ne Si,z . (6.17)
dt
In equilibrium, when ni,z = const, the density in charge state z is given by the balance
between recombination and ionization to and from adjacent states z − 1 and z + 1. These
two atomic processes are described by the recombination coefficient αi,z = αi,z (Te ) and
the ionization coefficient Si,z = Si,z (Te ). These are mainly dependent on the electron
temperature. Solving the system of equations (eq. 6.15) then gives the charge distribution
ni,z = ni,z (Te ) and average charge Z i = Z i (Te ) expressed in the electron temperature Te .
The density in a given excited state k, for a given charge number z is described by similar
excitation balance equations :
dni,z,k X
= {ni,z,l Ai,z,lk − ni,z,k Ai,z,kl + ni,z,l ne Ci,z,lk − ni,z,k ne Ci,z,kl } . (6.18)
dt l
The density of the excited state in equilibrium is given by the balance between excitation
and de-excitation processes; absorption and emission of radiation and electron collisions.
Here, the constants Ai,z,kl and Ai,z,lk are probabilities for absorption and emission; Ci,z,kl =
Ci,z,kl (Te ) and Ci,z,lk = Ci,z,lk (Te ) are collision coefficients which are functions of the
electron temperature. Under certain conditions (especially low electron density), the
excitation balance is mainly controlled by balance between excitation from the ground
state l = 1, caused by electron collisions and de-excitation by spontaneous emission
of radiation. This approximation is called corona equilibrium. The excitation balance
equation (eq. 6.16) for a state k can then be rewritten as
X
{−ni,z,k Ai,z,kl + ni,z,1 ne Ci,z,1k } = 0 (6.19)
l<k
Ci,z,1k
ni,z,k = ni,z,1 ne P . (6.20)
l<k Ai,z,kl
The density and thus the radiated power is proportional to the electron density ne and
the density in the ground state ni,z,1 , which is given by the ionization balance (eq. 6.15).
The total emitted power PR for a given ion type can be calculated as a function of the
electron temperature. The temperature dependence of the radiated power is described by
the radiation parameter Ri = Ri (Te );
PR = ne ni Ri (6.21)
135
10-30
Tungsten 8 Oxygen
10-31 z
4 Carbon
-32
R 3 10 Iron 0
(Wm ) 60
10-33 z
Oxygen 40 Tungsten
-34
10 20 Iron
Carbon
(a) (b)
10-35 0
1 10 100 1 10 100 1 10 100 1 10 100
eV eV eV keV keV keV eV eV eV keV keV keV
Te Te
Figure 6.3: Radiation parameters R and average charge number Z as a function of electron
temperature
The above expression give the radiated power per unit volume. The radiation parameter
R also include the bremsstrahlung so that PR is the total emitted power for a given
impurity. The left graph in figure 6.3 shows the radiation parameter R for some common
impurity ions as a function of Te . The right graph in figure 6.3 shows the average charge
number Z, also as a function of Te . The figure shows that low Z and high Z impurities
have different effect on the plasma. Low Z impurities cause maximum radiation power at
relatively low temperatures (Te ≈ 10 − 100 eV). At high temperatures (Te > 1 keV), they
are fully ionized and will not cause any line radiation. High Z contaminations have their
maximum radiated power at higher temperatures (Te ≈ 100 eV − 1 keV) and the radiated
power at reactor temperatures (Te ≈ 10 keV) is high, due to the higher bremsstrahlung
(PB ∼ Z 2 ).
For a given impurity, the densities in the different charge states, and thus the emitted
power from the corresponding lines, will vary with the electron temperature. At higher
electron temperatures, the higher charge states will become dominating. This can be
observed by a spectroscopic measurement of the emitted plasma radiation. Figure 6.4
shows an experimental measurement of the radiative power (in relative units) for different
charge states in oxygen (O1+ , O2+ , . . .) as a function of time. In this pulsed experiment
the electron temperature is increasing over time from Te ≈ 0 at time t = 0 to Te ≈ 250 eV
at t = 0.5 ms. Thus, the graph can also be seen as radiation power as a function of
electron temperature and it clearly shows that higher charge states have their maximum
radiative power at higher electron temperatures.
136
0.3
0.2
o 1+
0.1
0.0
0.4
o 2+
0.2
0.0
0.6
o 3+
0.4
0.2
0.0
0.4
o 4+
0.2
0.0
0.6
o 5+
0.4
0.2
0.0
6+
0.05
o
0.00
0.0 0.1 0.2 0.3 0.4 0.5 0.6
t(msec) time
Figure 6.4: Experimentally measured radiated power for different charge states in oxygen.
The electron temperature increases from Te ≈ 0 at time t = 0 to Te ≈ 250 eV at t = 0.5 ms.
excitation and charge exchange collisions. The flux of impurities is caused by different
processes releasing atoms and molecules from the wall, e.g. desorption, sputtering, erosion
and vaporization.
The flow of impurities from the wall increases radiation losses from the plasma. To
decrease the negative effects of plasma wall interaction the adsorbed impurity atoms are
removed from the wall. This is done using by conditioning the vacuum vessel through
“baking” and glow discharge cleaning. Other methods used to decrease the effects of
plasma wall interaction in tokamaks are the use of limiters and divertors.
Here, m1 , m2 and Z1 , Z2 are mass number and charge number for the incident ion (1) and
the wall material (2), and E is the ion energy in keV. A corresponding universal curve
also exists for the energy backscattering coefficient RE (figure 6.5).
137
0
10
Energy
Backscattering coefficient
reflection
-1
10
Particle
-2 reflection
10
-3
10
-3 -2 -1 0 1 2
10 10 10 10 10 10
Reduced energy ε ∝ E
The rate coefficient (< σv >) for these processes as a function of the particle energy is
shown in figure 6.6.
The physically adsorbed particles are loosely bound to the wall surface, i.e. the so called
binding energy is around 0.3 eV. The chemically adsorbed particles are more tightly
138
-13 H + Charge exchange on H
10
-14
-15
10
Charge exchange of
+
-16 H 2 on H
10
Proton ionization of H
-17
10
10 0 10 1 10
2
10 3 104 10
5
10 6
Particle energy (eV)
Figure 6.6: The reactivity for different atomic and molecular processes in hydrogen gas.
bound with a binding energy of around 3 eV. The flux of plasma ions to the wall surface
during a plasma discharge causes desorption of absorbed particles and thus a flux of low
Z impurities into the plasma.
To decrease the influx of impurities from the absorbed gas, the gas is removed as much
as possible before the fusion plasma discharge is started in the vacuum vessel. There are
two main methods for this:
1. Heating of the vacuum vessel, so called “baking” of the vessel. This method is used
to remove physically bound adsorbed gas.
2. “Cleaning” of the wall surface using a continuous, weakly ionized plasma discharge
(glow discharge) in the vacuum vessel (known as Glow Discharge Cleaning). Glow
discharges without magnetic field in hydrogen or helium are used. Plasma ions hit
the wall surface and chemically adsorbed gas is released.
139
Carbon
Stainless steel He
-1 Molybdenum
10
He
H
-3
10
-4
10
-2 -1
10 10 1 10
Ion energy (keV)
Figure 6.7: Sputtering yield for different combinations of particles and materials.
λD
V
-V s
Figure 6.8: Schematic picture of the space charge sheath at the wall.
140
3
10
Cr
-6
Ti
10
-9
Ni
10
Figure 6.9: Vapour pressure as a function of temperature for different wall materials.
6.2.7 Limiters
A plasma limiter is a material surface extending some distance from the surrounding wall
(figure 6.10). Limiters can have different geometries, either local or extended in poloidal
or toroidal direction.
The purpose of the limiters are mainly to:
141
Wall
Lc
Boundary layer Γ =nV
Γ =-D dn
Limiter dr
Confined plasma Separatrix
a
B
142
which have the solution x
n(x) = n0 exp − (6.28)
λ
where n0 is the density at r = a and the characteristic length for the density variation in
the boundary layer is p
λ = D⊥ τk . (6.29)
The confinement time τk is calculated by dividing the total number of particles in the
boundary layer N = ALc n by the total flow along field lines to the side surfaces of the
limiter 2A
N Lc n
τk = = . (6.30)
2AΓk 2Γk
The flux to the wall surface is controlled by a negative space charge layer at the wall
(figure 6.8) limiting the electron flux so that ion and electron fluxes are equal
Γk = 0.5cs n (6.31)
where r
2eTe
cs = . (6.32)
mi
Inserting eq. 6.29 in eq. 6.28 finally gives
Lc
τk = . (6.33)
cs
1. Impurities are pumped from the target chamber and is not brought into the plasma.
2. Lower power flux to the wall with a so called gas target decreases the production of
impurities.
143
z
φ
Main chamber
Separatrix
ID Divertor
By increasing the neutral gas density in the target chamber the power released as elec-
tromagnetic radiation from the plasma in the divertor is increased. The power flux to the
wall in the divertor is thus decreased correspondingly and the production of impurities
at the wall surface is reduced. The density gradient along field lines in the outer mag-
netic field area required for this scenario (which is called gas target) can in principle be
maintained by the friction force from the plasma flux into the divertor.
PΩ = ηj 2 (6.34)
where η is the plasma resistivity and j is the plasma electrical current density. The
plasma resistivity depends on Coulomb collisions and decreases with increased plasma
temperature (η ∼ T −3/2 ). This means that the ohmic heating is less efficient at high
plasma temperatures and other heating methods are needed as well.
In a tokamak, the toroidal plasma current is limited by the condition for MHD stability.
MHD stability in the tokamak requires that the safety factor q > 1 in the centre (see
chapter 4). The current density in the centre is in the toroidal direction, i.e. j = jϕ .
Ampére’s law provides a relation between the poloidal magnetic field Bθ and the toroidal
current density in the centre (where jϕ is assumed constant):
1
Bθ = µ0 rjϕ . (6.35)
2
The stability condition q > 1 then gives the following condition for the current density:
2Bϕ
jϕ < (6.36)
µ0 R
144
where R is the major radius and Bϕ is the toroidal magnetic field. Inserting eq. 2.76
gives the following condition for PΩ :
2
−5 −3/2 2Bϕ
PΩ < 5.22times10 Zeff ln ΛT (6.37)
µ0 R
The energy confinement time for ohmically heated tokamaks has been experimentally
measured and can, in a first approximation, be described by the following simple empirical
scaling law (see Table 5.1):
τE ≈ 5 × 10−21 na2 [s] (6.39)
where n is the plasma density and a is the minor radius. Since the plasma energy density
is W = 3neT the power loss per unit volume will be
3neT T
Ploss = = 9.6 × 101 2 [W/m3 ]. (6.40)
τE a
The energy balance equation (PΩ = Ploss ), with eq 6.36 and eq. 6.38 inserted, gives the
following condition for the plasma temperature
a 4/5
T < 8.5 × 102 Bϕ [eV]. (6.41)
R
For a tokamak assumed to have aspect ratio R/a = 3 and the magnetic field Bϕ = 10 T
we have that T < 2.2 keV, which indicates that the plasma temperature possible to reach
by ohmic heating in this case is limited to 2–3 keV.
145
Energy of hydrogen beam atoms
-19 10 keV 100 keV
10
Charge
Exchange (σch )
Cross-
section
(m2 ) Ionization
byions (σi )
-20
10
T
e =1
keV
Electron
ionization T =
e 10
keV
-21
10
10 keV 100 keV 1 MeV
Energy of deuterium beam atoms
Figure 6.12: Cross section for some atomic processes in hydrogen gas.
3. Charge exchange collision between hydrogen atom and hydrogen ion (σch )
Hbeam + H+ → H+beam + Hplas
The cross sections for these processes as a function of the energy of the injected neutral
atoms are shown in figure 6.12.
An expression for the absorption of the neutral beam due to the mentioned processes
can be obtained in the following way: The intensity of the beam I(x) as a function of the
distance x the beam has passed through the plasma is written as
where vb is the velocity of the beam particles and Nb (x) is the number of particles per
unit length in the beam. The intensity is the number of beam particles incident on a
surface in the plasma perpendicular to the beam per unit time. Consider a disc of plasma
with area A equal to the cross section of the beam and with thickness dx. The number
of plasma particles in the disc is
Np = nAdx (6.43)
where n is the plasma density. The probability p for collision between an incident beam
particle and a plasma particle in the disc is described by the collision cross section σ:
σ
p= Np = σndx. (6.44)
A
The number of particles in the beam dI absorbed (per unit time) in the disc due to
collisions is then
dI = −pI = −σnIdx (6.45)
146
1.0 0
0 1.0
1 10 100
E b/ Te
Figure 6.13: The relation between heating of electrons and heating of ions.
thus
dI
= −σnI. (6.46)
dx
This differential equation has the solution
x
I(x) = I0 exp − (6.47)
λ
where
1
λ= (6.48)
σn
is the characteristic penetration depth for the neutral beam in the plasma and I0 is the
incident intensity at the plasma surface. The collision cross section is the sum of the cross
sections for the different processes described above (as shown in figure 6.12), i.e.
1
σ = σch + σi + < ve σe > . (6.49)
vb
In the case of ionization by electron collisions, where the electron velocities are higher than
the neutral atom velocity, an integration (< . . . >) is made over the velocity distribution
of the electrons.
The energetic ions created by absorption of the neutral beam transfer their energy
to the plasma particles when they are decelerated by Coulomb collisions. The fast ions
in the beam, with mass mb and velocity vb , experience a friction force Fbe from electron
collisions
mb vb
Fbe = (6.50)
τie
and a corresponding friction force from ion collisions
mb vb
Fbi = (6.51)
τii
where τie and τii are collision times for electron and ion collisions respectively. The heating
power transferred to each particle type can be written as
Pe = vb Fbe (6.52)
and
mb
Pi = vb Fbi (6.53)
mb + mi
147
Accelerator
Ion Deflector
source Neutralizer magnet
Neutrals
D+
D+ D0
D0
Plasma
Gas D+
Ions
1 2 3 4
where the energy from ion collisions is divided between the beam ion and the plasma
ion. (In electron-ion collisions, the electron is assumed to absorb all energy, because of its
smaller mass.) The collision time increases with increased relative velocity v between the
colliding particles (τ ∼ v 3 ) (eq. 2.31). This means that, according to eq. 6.48–6.51, at low
neutral beam energies the ions will absorb most of the heating power, and at high neutral
beam energies the electrons will absorb the most. The ratio of the electron heating power
Pe and the ion heating power Pi to the total heating power Pe + Pi as a function of neutral
beam energy Eb is shown in figure 6.13.
Charged particles are produced in an ion source. After the acceleration, the charged
particles pass through the neutralization chamber. Since the neutralization is not perfect,
a deflection magnet is used to remove the remaining energetic ions from the beam.
The neutralization chamber contains neutral gas and the charged particles are neu-
tralized by charge exchange collisions with neutral atoms. The energetic ions created
can be ionized again by collisions with neutral atoms. The relation between neutral and
ion densities, the equilibrium neutral fraction, will be determined by the ratio of the
cross-sections for the two inverse processes:
N0 σch
= (6.54)
N+ σion
The neutral fraction as a function of the beam energy is shown in figure 6.15. The figure
also shows how the penetration depth (eq. 6.46) depends on the beam energy.
148
1.0
0.8 Neutral
fraction
0.6
0.4
Penetration
0.2 distance(m)
0
10 20 40 60 80 100 200 300
Beam energy (keV)
Figure 6.15: Neutral fraction and penetration depth as a function of beam energy.
z
Bo
k
θ
y
x
Figure 6.16: Diagram showing different directions assumed in the wave model
which describes plane waves with angular frequency ω propagating in a direction given
by the wave vector k (figure 6.16). A wave is described by the relation k(ω), which is
called the dispersion relation. In this section we use E, B, . . . to describe the complex
149
amplitude of a perturbation quantity. Introducing the plane wave anzats in Maxwell’s
equations gives
ik × E = iωB
ik × B = µ0 j − µ0 ε0 iω.E (6.56)
The relation between the vectors j and E in the plasma is described by the conductivity
↔
tensor σ:
↔
j = σ · E. (6.57)
↔
The plasma effect on the wave propagation is contained in the conductivity tensor σ. The
conductivity tensor is calculated from the linearised fluid equations (see chapter 4). In
this model the plasma pressure is assumed small and collision terms are neglected. The
↔
force balance equation and Ohm’s law gives an expression for σ. The refractive index n
is defined as
c
n = k. (6.58)
ω
Eq. 6.54 and 6.55 gives the wave equation
↔
n2 E − n(n · E) − ε · E = 0 (6.59)
↔
where ε is the dielectric tensor ,
i
εkl = σkl + δkl . (6.60)
ωε0
The dielectric tensor has the following form:
ε⊥ εxy 0
εkl = −εxy ε⊥ 0 (6.61)
0 0 εk
where
2
X ωpj
ε⊥ = 1 − 2
(6.62)
j
ω 2 − ωcj
2
X ωpj
εk = 1 − (6.63)
j
ω2
2
X ωpj ω
cj
εxy = i 2
. (6.64)
j
ω 2 − ωcj ω
Summation over j is over particle types (ions, electrons). Here, ωcj is ion and electron
cyclotron frequencies. respectively (eq. 2.9) and ωpj is ion and electron plasma frequencies,
respectively (eq. 4.18). Assume that the wave vector is in the yz plane (see figure 6.16)
where θ is the angle between k and B. Introduce
nk = n cos θ
n⊥ = n sin θ. (6.65)
The wave equation (eq. 6.57) can then be written in matrix form
ε⊥ − n2k − n2⊥
εxy 0 Ex
2
−εxy ε⊥ − n k n k n ⊥ Ey = 0. (6.66)
0 nk n⊥ εk − n2⊥ Ez
150
E
z ωt y
Nontrivial solutions to the homogeneous system of equations (eq. 6.64) are obtained when
the determinant of the matrix is equal to zero. This condition finally gives the plasma
dispersion relation, expressed using n2k and n2⊥ ;
ε⊥ n4⊥ − (ε⊥ − n2k )(ε⊥ + εk ) + ε2xy n2⊥ + εk (ε⊥ − n2k )2 + ε2xy = 0.
(6.67)
The waves propagating in the plasma are described by n2k and n2⊥ that are solutions to
eq. 6.65.
Consider the special case k k B0 , i.e. n⊥ = 0. The following equation is obtained for
2
nk ;
(ε⊥ − n2k )2 + ε2xy = 0. (6.68)
Introduce
R = ε⊥ + iεxy
L = ε⊥ − iεxy (6.69)
which gives two solutions
“R − wave00
R
n2k = (6.70)
L “L − wave00 .
Consider the first component of the wave equation (eq. 6.57), which gives a relation
between the perpendicular components of the wave field E:
(ε⊥ − n2k )Ex + εxy Ey = 0 (6.71)
from which is obtained
−i “R − wave00
Ex
= (6.72)
Ey
i “L − wave00 .
The two perpendicular components of the wave are π/2 out of phase, which means that
the wave is circular polarized, and the E field vector will rotate with angular frequency ω.
If the wave is viewed along z direction, in direction of decreasing z (figure 6.17), the “R”
wave rotation is counter clock-wise, i.e. the “R” wave is right circular polarized (RCP ).
For the “L” wave the rotation is clock-wise, i.e. the “L” wave is left circular polarized
(LCP ).
In the other special case k⊥B0 , i.e. nk = 0, the following equation is obtained for n2⊥ ;
ε⊥ n4⊥ − (ε⊥ εk + RL)n2⊥ + εk RL = 0 (6.73)
which has two solutions
ε “O − wave00
k
n2⊥ = (6.74)
RL
“X − wave . 00
ε⊥
151
Δ
n n 0
152
V
E
ω
k
which, after inserting eq. 6.63 gives the condition for resonance
εk
= − tan2 θ. (6.77)
ε⊥
In the case k k B0 , i.e. θ = 0 there will be resonances according to eq. 6.75 when ε⊥ → ∞.
Eq. 6.60 shows that this occurs when
2
ωci “ion cyclotron resonance”
ω2 = (6.78)
2
ωce “electron cyclotron resonance”
For wave propagation perpendicular to the magnetic field, when k⊥B0 , i.e. θ = π/2,
resonances is given by eq. 6.75 with the condition ε⊥ = 0. Eq. 6.60 gives, with reasonable
2 2 2 2 2 2
approximations ωpi ωpe , ωci ωce , ωci ωpi the following equation for ω 2 ;
ω 4 − ω 2 (ωce
2 2
+ ωpe 2 2
) + ωpi ωce ≈ 0 (6.79)
153
ω ci
+
z
k BO y
ωt
E
ω
propagating along the magnetic field (k k B0 ) (figure 6.21). Seen from the positive z axis,
the wave E field vector is rotating clock-wise with angular frequency ω. Ions in the plasma
also gyrate clock-wise around magnetic field lines with the ion cyclotron frequency ωci .
When the wave frequency fulfils the resonance condition ω ≈ ωci , the ion will experience
an electric field mainly perpendicular to the gyration movement during a longer time.
This will accelerate the ion. The acceleration cause energy to be transferred from the
wave to the ion, thus damping the wave. This process is called cyclotron damping. In
this case as well, wave-particle resonance is necessary, since if ω 6≈ ωci the E field will no
longer be in phase with the particle. Thus, the acceleration is successively in different
directions, preventing an effective energy transfer.
Through wave-particle interaction described above (Landau damping and cyclotron
damping) energy will be transferred from the wave to the particles. The wave is damped,
or absorbed. At the same time, the energy transfer increase the particle’s kinetic energy,
thus heating the plasma. This is the effect used for heating of the plasma with radio
frequency waves (RF heating).
2. Lower Hybrid Resonance Heating (LHRH). The typical frequency is f ≈ 1–5 GHz.
The resonance condition depends, through the plasma frequency, on the plasma
density. This means that the absorption region in the plasma varies with the plasma
154
density. LHRH heats the ions through damping at the lower hybrid resonance, but
can also heat the electrons by electron-Landau damping.
155
Chapter 7
PLASMA DIAGNOSTICS
Plasma diagnostics is the common name for different types of measurement equipment
used to study, i.e. diagnose, fusion experiments. The purpose is to obtain information
about plasma quantities such as temperature, density and so on. A space resolved mea-
surement is usually necessary, since it allows reconstruction of radial profiles of these
quantities. Sometimes dynamic events are studied, and this requires measurements with
high time resolution. These instruments can both be simple, e.g. a magnetic probe mea-
suring the magnetic field outside the plasma, and complicated, e.g. equipment used to
measure the electron temperature by looking at the scattering of electromagnetic waves
in the plasma. A large tokamak experiment usually has several diagnostic systems si-
multaneously measuring important plasma parameters. Quite often, several methods are
used to measure a parameter, e.g. one diagnostic used to measure a parameter locally
in different positions in the plasma at a certain time combined with another measuring
the same parameter in a larger area over time. One such example is the measurement of
the electron temperature, which is done locally at a specific time by scattering of elec-
tromagnetic waves and also continuously by a measurement of bremsstrahlung radiation
emitted from a larger part of the plasma. An example of the use of information from
several plasma diagnostics is the estimation of the confinement time τE in an ohmically
heated tokamak. This requires a simultaneous measurement of electron temperature Te ,
ion temperature Ti , electron density ne , average ion charge (which gives the ion density),
toroidal electric current in the plasma I, toroidal loop voltage V (which gives the ohmic
power input PΩ = V · I).
Methods for plasma diagnostics can be sorted by the physical property of the plasma
on which the measurement is based:
1. Magnetic field. Measurements of the magnetic field inside or outside the plasma
using different types of probes or measurement coils.
2. Plasmas particle flux. Different types of probes in direct contact with the plasma
measuring the flux of plasma particles.
156
5. Electromagnetic radiation from bound electrons. Spectroscopic measurement of line
radiation emitted from, for example, not fully ionized impurities in the plasma.
Some examples of diagnostic methods and plasma quantities measured are shown in Table
7.1. These methods are only a small selection of those used at large fusion experiments.
Some of the methods are described in more detail in this chapter.
7. Ion processes
neutral particle analysis neutral density, ion temperature na , Ti
neutron spectroscopy ion temperature Ti
157
Figure 7.1: Schematic picture of Rogowski coil used to measure plasma current.
of a long solenoid coil encircling the plasma in the poloidal direction. The coil provides a
measurement of the line integral of the poloidal magnetic field component along poloidal
direction, and Ampére’s law gives the toroidal current passing inside the coil as
I
1
I= B · dl. (7.1)
µ0
Consider a coil with N loops and length L. The voltage from loop i is equal to the time
derivative of the magnetic flux through the loop (figure 7.2).
d ˆ i = − δA d (Bi · dli ) .
Vi = −δA Bi · ds (7.2)
dt dl dt
The total voltage from the coil is
N N
! I
X N δA d X N δA d N δA dI
V = Vi = − Bi · dli ≈− B · dl =− µ0 . (7.3)
i=1
L dt i=1
L dt L dt
dA
ds i
B
L
dl = N
158
7.2.1 Principle of interferometry measurement
An expression for the phase change of an electromagnetic wave passing through a plasma
can be described in the following way. Consider an electromagnetic wave in vacuum
ω 2π
= k0 = . (7.4)
c λ0
Here, k0 is the wave-number and λ0 is the wavelength in vacuum. When the electromag-
netic wave propagates in the plasma, it is in the form of an “O” wave described by the
plasma dispersion relation (see the section on wave propagation in the plasma) given by
eq. 6.56, 6.61 and 6.72: 2 2
k ωpe
= k = 1 − 2 (7.5)
k0 ω
where k = 2π/λ and λ is the wavelength in the plasma. We have assumed that ω ωpi .
Since we also have the condition ω ωpe we have
r
2 2
k ωpe 1 ωpe
= 1− 2 ≈1− . (7.6)
k0 ω 2 ω2
The phase for a measurement beam passing a distance L in the plasma can be written as
Z
φ= k dl. (7.7)
L
Plasma
dl
O
ref. O
Ο
159
7.2.2 A laser interferometer system
An interferometer for continuous measurement of the electron density often uses modu-
lation of the incident electromagnetic wave. This is desirable because it helps to avoid
problems when interpreting the measured phase shift. Is some cases it is impossible to
decide if the plasma density is increasing or decreasing when the measured phase shift
changes. This problem can be avoided by modulation. An example of an interferometer
1 2 ω
ω
I
ω BC l + ω Light detector
BC
OSC.
40 MHz to quad. unit ω 4
3 BC
Quadrature unit
with modulation technique for measuring a time dependent plasma density is shown in
figure 7.4. This interferometer uses a carbon dioxide laser as wave source. The interfer-
ometer consists of the following components:
2. Acoustic-optical modulator, a so called Bragg cell . The Bragg cell divides the laser
beam into two parts by deflecting part of the incident beam. The deflected beam
part is frequency shifted ωBC /2π = 40 MHz, which is the carrier frequency of
the Bragg cell. The deflected beam is the interferometer reference beam IR ∼
cos(ωl t + ωBC t). The beam which is not deflected is the interferometer measurement
beam IM . It passes through the plasma, is reflected by mirror M 1 , passes through
the plasma again and is reflected against the exit window of the Bragg cell in the
same direction as the reference beam. The measurement beam that passes through
the plasma will have a time dependent phase shift ∆φ(t) relative the reference beam;
IM ∼ cos(ωl t + ∆φ(t)).
3. Detector. The interference between the measurement beam and the reference beam
takes place in the detector. When the beams are added the result is a signal which,
after low-pass filtering, has a frequency equal to that of the difference of the two
beam frequencies, ID ∼ cos(ωBC t − ∆φ(t)).
4. Quadrature unit. Quadrature detection means that the detector signal is multiplied
(mixed) with the signal from the RF-oscillator driving the Bragg cell. After low-pass
filtering the result is both a cosine signal IC ∼ cos(∆φ(t)) and, by phase shifting
the carrier frequency π/2, a sine signal IS ∼ sin(∆φ(t)). (figure 7.5)
160
Mixer
Phase shifter LP-filter
cos ωBC t
osc. ωBC o
0 C cos φ
o
90 sin ω t
BC
Digitized
U= ωBC t - φ 0
o
cos U S sin φ
o
0
Finally, we get the phase shift as a function of time ∆φ(t), and thus the time variation
of the line density by digitizing the cosine and sine signals from the quadrature detector
and calculate ∆φ(t) = arctan(IS /IC ).
a R
n.
θ
ikt
s.r
ob
An electron being accelerated will emit radiation which can be described using the
expression for radiation from a dipole. The radiated power per unit solid angle in a given
direction of observation is
dP e2
= a2 sin2 θ (7.12)
dΩ 16π 2 0 c3
where a is the electron acceleration and θ is the angle between the observation direction
and the acceleration vector a (figure 7.6). Assume that a plane electromagnetic wave is
incident on the plasma
Ei = Ei0 exp i(ki · r − ωi t) . (7.13)
161
The electron is accelerated in the electric field of the incident wave
eEi
a= . (7.14)
me
The average power of the scattered wave is then
dP e2 Ei02
= sin2 θ. (7.15)
dΩ 16π 2 0 c3 m2e
ki
Ei
v
ks
vt
R2
Es R1
R0
The Doppler effect can be easily understood from figure 7.7. The scattered wave is
assumed to be generated by an electron at a point given by the vector R1 . The electron is
assumed to have velocity v, and the electron position at time t is given by the expression
R1 = R0 + vt. (7.18)
162
In the observation point, R2 , the scattered wave is
Es ∼ Ei (R1 , t) exp iks · (R2 − R1 )
∼ exp i(ki · R1 − ωi t + ks · (R2 − R1 ))
∼ exp i(ks · R2 − ωi t − (ks − ki ) · R1 ) . (7.20)
Inserting eq. 7.18 gives
Es ∼ exp i(ks · R2 − ωi t − (ks − ki ) · R0 − (ks − ki ) · vt) . (7.21)
A comparison between eq. 7.19 and eq. 7.21 gives the following relation between ωs and
ωi
ωs = ωi + (ks − ki ) · v. (7.22)
The frequency shift ω and the wave vector k is introduced
ω = ωs − ωi
k = ks − ki (7.23)
which means that eq. 7.22 can be rewritten in the following simple way
ω = k · v. (7.24)
This relation shows that the frequency shift is proportional to the electron velocity com-
ponent in the direction of the wave vector k. Assume that the electrons have a velocity
distribution f (v) along k. A thermal distribution (Maxwell distribution) is described as
me v 2
r
me
f (v)dv = exp − dv. (7.25)
2πeTe 2eTe
It is usually assumed that individual electrons move and radiate independent of each other
(i.e. the electrons are assumed uncorrelated). With this assumption, we get the spectral
function S(ω) for the scattered wave by inserting v = ω/k in eq. 7.25;
me ω 2
r
me 1
S(ω)dω = exp − dω. (7.26)
2πeTe k 2eTe k 2
The radiative power per unit solid angle, angular frequency and volume is finally (for
θ = π/2) :
dP
= ne r02 Pi S(ω). (7.27)
dΩdωdV
A measurement of the electron temperature using Thomson scattering is done by shining
a laser beam through the plasma and measure the spectral function S(λs ) of the scattered
light using spectrometry. Because of this, the spectral function (eq. 7.26) is written as a
function of λs , the wavelength of the scattered light. Assume that the angle between ki
and ks is π/2;
ω 2 (ω − ω )2 2 2
s i 2 (λi − λs ) 2 λ
= =c ≈c (7.28)
k ks2 + ki2 λ2i + λ2s 2λ2i
where λ = λi − λs has been used. The spectral function S(λ) is given by eq. 7.26 as
2 !
m e c2 λ
S(λ) ∼ exp − (7.29)
4eTe λi
from which the electron temperature is estimated
√ by measuring the spectral function
half-maximum width This is proportional to Te (figure 7.8).
163
S(ω)
ω
ω α Te
ωo ωs
1. Laser. Usually a pulsed laser with high power, ∼ 1 GW, is used. A common laser
type uses a ruby rod and can give one or several pulses. The beam is focused by
a lens to a selected point in the plasma. To decrease the scattering of light from
the vessel wall, the beam is lead through a number of apertures on its way to the
plasma. A beam dump stops light from being reflected back after passing through
the plasma.
2. Spectrometer. The light scattered an angle π/2 from the incident laser beam is
observed using a lens from a selected region, the scattering volume. The spectrum
of the scattered light is analysed with a spectrometer using filters for the wavelength
intervals measured in the spectrum. The light passing a filter is detected with a
photo multiplier. The measured scattered radiation power (∼ 1 mW) is very small
compared to the incident radiation power.
164
7.4 Measurement of electron cyclotron radiation
The electromagnetic radiation emitted by electrons when they gyrate in the magnetic
field, the electron cyclotron radiation, can be used to measure the plasma electron tem-
perature. In general, if the electrons emit radiation of a given frequency ω, they will also
absorb radiation of the same frequency. The electrons can be assumed to be in thermal
equilibrium, i.e. they have a thermal velocity distribution given by an electron tempera-
ture Te (Maxwell distribution). In this case, the emittance j(ω) has a simple relation to
the absorption coefficient α(ω) (Kirchoff’s law);
where IB (ω) is the radiation intensity for a black body . The absorption coefficient is in
the range 0 < α(ω) < 1. A black body is defined as a medium absorbing all radiation
which is incident on it. Its radiation intensity is given by Planck’s radiation law
−1
h̄ω 3
h̄ω
IB (ω) = 3 2 exp −1 . (7.31)
8π c eTe
Here, h̄ = h/2π, and h is Planck’s constant. The intensity is defined as radiation power
per unit area, solid angle and angular frequency. The electron cyclotron radiation can in
this case be considered low frequency, h̄ω eTe , and an approximative expression for IB
can be used
ω 2 eTe
IB (ω) ≈ . (7.32)
8π 3 c2
The radiation intensity observed in a given direction outside the plasma is caused by
emission and absorption of the beam along its path through the plasma.
The beam path can be considered a straight line given by the coordinate s where
s1 < s < s2 . The plasma far boundary is given by s = s1 (where the beam “starts”)
and the plasma near boundary is given as s = s2 (where it “leaves the plasma”). The
intensity of the beam along the path is given by the beam transport equation which gives
the change dI of the intensity when the beam travels a distance ds;
dI = j ds − αI ds. (7.33)
dτ = −α ds. (7.34)
The quantity τ (s) is the optical depth for the plasma layer s < s0 < s2 . Inserting eq. 7.30
and 7.34 in eq. 7.33 gives
dI
= −IB + I (7.37)
dτ
165
a) 3 ω ce 5.0
2 b)
I α
B
ω
4.0
2ω ce
7s
3.0
Te (keV)
3s
I B(I - e-τ )
I
2.0
ω ce
ω 2 Te (R) 1.0
0.5s
0
0 50 100 150 200 250 300 350 3.0 3.5 4.0
Frequency (GHz) Major radius (m)
The intensity incident on the far plasma boundary can be neglected, i.e. I(s1 ) = 0. The
following expression for the intensity observed from a plasma with optical depth τ can
then be obtained from eq. 7.39;
Two limiting cases for the emitted radiation intensity, depending on the plasma re-
absorption, exists:
1. Optically thin plasma, τ 1: In this case, the intensity of a beam with length L in
the plasma is equal to I ≈ IB τ ≈ IB αL = jL. The plasma absorption of radiation
is negligible.
2. Optically thick plasma, τ 1: In this case, I ≈ IB , i.e. the plasma radiates with
the same intensity as a black body. This is the highest possible radiation intensity.
In a plasma where the magnetic field varies over space, electron cyclotron radiation can
give information on the local electron temperature, i.e. the plasma temperature profile
Te (R) can be calculated. Emission of electron cyclotron radiation of a given frequency
ω occurs at a resonance layer in the plasma where the resonance condition ω ≈ nωce is
fulfilled. According to Kirchoff’s law absorption of radiation of this frequency will mainly
take place in the resonance layer. Outside the resonance layer we can assume that the
radiation passes with negligible absorption, i.e. the plasma is considered optically thin.
(τ (ω) 1). If the resonance layer is optically thick (τ (ω) 1), the emitted electron
cyclotron radiation is given directly by IB (ω) and thus the local electron temperature Te
166
Figure 7.11: Example of measured X-ray radiation spectrum.
in the resonance layer. The radiation intensity is in general given by eq. 7.40 using τ (ω),
the optical depth for the resonance layer.
The magnetic field in a tokamak is, in a first approximation, assumed equal to the
vacuum toroidal field, B ∼ 1/R, where R is the major radius. The resonance condition
for electron cyclotron radiation then gives
R0
ω = nωce0 (7.41)
R
where R0 is the plasma major radius in the centre and ωce0 is the electron cyclotron
frequency at R0 . The resonance frequency thus varies with position in the plasma as
ω ∼ 1/R. The observed broadening of the emission line for a given harmonic n, i.e. the
continuous radiation spectrum around ω ≈ nωce0 (figure 7.10a), thus gives the electron
temperature as a function of the frequency
8π 3 c2
Te (ω) = I(ω) (7.42)
eω 2
from which the profile of Te (R) is given if we use the relation between ω and R (eq. 7.41)
(figure 7.10b).
n2e Zeff
dN E
E ∼ √ exp − (7.43)
dE eTe eTe
167
Torus wall Liquid
N 2tank
Aperture and Be
foil box window
Plasma 2
1
Flight Si(Li)
Aperture box tube crystal
3
To multichannel Connection
analyser to electronics
X-ray pulse Amplifier
Figure 7.13: Schematic picture of the signal from the silicon detector.
where dN/dE is the total number of photons with energy E emitted per unit energy, time
and volume. The bremsstrahlung gives a continuous spectrum. Measuring the absolute
value of the intensity in a given energy interval can in principle give an estimate of the
effective ion charge Zeff . The quantity E(dN/dE) has a simple exponential dependence
on the photon energy E. This can be used to estimate the electron temperature Te (figure
7.11). Besides the continuous spectrum there is also a discrete spectrum, caused by line
radiation which complicates the measurement. Usually an energy interval ∼ 1–10 keV of
the spectrum is used, i.e. soft X-ray radiation. An example of an experimental setup for
measurement of the soft X-ray spectrum is shown in figure 7.12.
The diagnostic consists of the following components:
1. Absorption foil. A thin disc of a suitable material (e.g. beryllium) absorbs the
low-energy part of the spectrum since it should not be analysed.
3. Multichannel pulse height analyser. Pulse height analysis, i.e. counting the num-
ber of voltage pulses in a certain voltage interval, finally gives the photon energy
spectrum.
168
-18
σ 0l -13
10 10
σ 0l υ i (m 3 s -1)
-19 -14
10 10
σ 0l (m 2 )
σ 0l υ i
-20 -15
10 10
-21 -16
10 10
1 10 10
2
10
3
10 4 10 5
1
2m i νi2(eV)
Figure 7.14: Cross section for charge exchange as a function of ion energy.
neutral hydrogen atom is transferred to the ion. Thus, a new neutral atom is formed.
This new neutral particle will have high energy since the ion keeps its kinetic energy when
capturing the electron. In this way, a new population of “fast” neutral particles is created
in the plasma.
+ +
Hslow + Hfast → Hslow + Hfast (7.44)
These high energy neutral particles leave the plasma since they are not confined by the
magnetic field. The “fast” neutral particles will have the same velocity distribution as the
ions. The ions usually have a thermal velocity distribution given by the ion temperature
Ti . The cross section for the charge exchange collision is known (figure 7.14). Thus, an
analysis of the energy spectrum of the flux of fast neutral particles outside the plasma is a
measurement of the ion temperature Ti . Measuring the absolute value of the fast neutral
flux also tells us the density of the “slow” neutral component in the plasma.
An expression for the spectral density of the fast neutral component is obtained in
the following way: Consider the normal case that the electrons have a thermal velocity
distribution. The distribution function f (v)dv gives the number of ions with velocity v
in the interval dv per unit volume in the plasma
3/2
mi v 2
mi
f (v)dv = n exp − 4πv 2 dv (7.45)
2πeTi 2eTi
where n is the ion density and mi is the ion mass. The reactivity for charge exchange
collisions for a given velocity v is
ξ(v) =< σch (|v − va |)|v − va | >va (7.46)
where σch is the cross section for charge exchange collisions (which depends on the rel-
ative velocity of the colliding particles), v is the ion velocity, va is the velocity of the
“slow” neutral atoms, and < . . . >va denotes an averaging over the slow neutral velocity
distribution. Since |v| = v va , the following expression for the reactivity can be used
as an approximation
ξ(v) ≈ σch (v)v (7.47)
The number of fast neutral particles created by charge exchange collisions with velocity
v in an interval dv per unit time and unit volume is
S(v)dv = na ξ(v)f (v)dv = na σch (v)vf (v)dv (7.48)
169
Neutral
flux
35
30
25
0 10 20 30 40
E 0(keV)
Energy
where na is the density of slow neutral particles. Assume that the plasma volume observed
can be considered a cylinder with area A and length L, and that the detector sees a flux
within a solid angle Ωs . The number of particles detected per unit time can then be
written as
Ωs
F (v)dv = A L S(v)dv. (7.49)
4π
Inserting eq. 7.45 and 7.48 in 7.49 gives the following expression for the measured neutral
particle flux;
3/2
mi v 2
mi
F (v)dv = Ωs ALnna σch (v) exp − v 3 dv. (7.50)
2πeTi 2eTi
Alternatively, the expression can be written as a function of the ion energy E (figure
7.15);
3/2
1 E E
F (E)dE = Ωs ALnna σch (E) √ exp − dE. (7.51)
πeTi 2mi eTi
The spectral density F (E) divided by σch (E) E has a simple exponential dependence
on the ion energy E
F E
∼ exp − (7.52)
σch E eTi
which gives the ion temperature Ti .
170
2 Magnetic
pole piece
+
H trajectory B
Channeltron
detectors
171
Chapter 8
The tokamak is the most developed magnetic confinement configuration. The reason for
this is the good confinement properties and relatively simple geometry of the tokamak.
The resources for tokamak research, in terms of capital investment and man-years, have
been incomparably stronger than those for other confinement schemes.
The diagram in figure 8.1 summarizes results from some of the largest tokamak exper-
iments JET (Europe), TFTR (USA), JT-60U (Japan) and DIII-D (USA). Data from the
compact high-field tokamak Alcator (USA) is also included. By the Q factor we denote
the relation between injected energy and generated fusion energy. Thus, Q = 1 stands for
“break even”. It is clear that enormous progress has been achieved. The triple product
nτE T has increased by a factor 1000, and only a factor of five remains to the region, where
ignition is obtained.
In this chapter we will familiarize ourselves with some of the most successful alter-
natives to the conventional tokamak. This includes both variations of the tokamak, e.g.
the advanced, the spherical and the compact tokamaks, and quite different configurations
like the reversed-field pinch (RFP), the stellarator and the high density pinch. We will
also briefly discuss the principles for inertial confinement. Finally, we will touch upon the
hybrid reactor, which combines characteristics of fusion and fission reactors.
8.1 Tokamak
But let us start with the future development of the conventional tokamak. A schematic
drawing of a conventional tokamak is shown in figure 8.2.
The next step device in the international tokamak research line is the ITER device,
currently under construction in France.
8.1.1 ITER
ITER is latin for ”the way”. It is the name of the next-step device in the tokamak research
line, aimed to ”lead the way to fusion energy”. It is a scaled-up version of present-day
tokamaks. Based on extrapolations of the fusion plasma performance in these devices,
ITER is projected to generate 500 MW of fusion power. It will be the next major step in
fusion research. ITER is a collaboration between EU, Japan, Russia, China, South Korea,
India and USA. Hence nations, representing more than half of the world’s population are
engaged in the project.
172
Figure 8.1: nτE T diagram (Figure credit: EUROfusion-JET).
173
ITER is currently (2022) under construction at the Cadarache site in France, and is
expected to come into operation around 2025. ITER is planned to be the final step before a
demonstration reactor aimed at electricity production is constructed. The most important
task for ITER is to generate significant self-heating from nuclear fusion reactions - at the
same time providing a good enough plasma energy confinement - so that the tokamak can
be maintained in a nearly self-sustained ”burning” state. This has never been achieved in
any tokamak, and reaching this milestone would take fusion research a big step forward to
the realisation of fusion reactors. ITER will also serve as a test bench for the technology
of an integrated reactor system. This includes effects of high heat and neutron flux in
addition to studies of the functionality of reactor components.
ITER is a large tokamak with superconducting magnets generating about 500 MW
of fusion power for at least 400 seconds. The plasma volume will be ten times that of
JET and will be close to the size of future commercial reactors. The ITER project will,
for the first time, enable the study of plasmas heated by the hot alpha particles (helium
nuclei) produced in the fusion reaction. ITER will generate ten times more power than the
externally input power required to heat the plasma, equivalent to an energy multiplication
factor of Q ∼ 10.
For more information on ITER, see www.iter.org.
174
Figure 8.3: The riggatron.
Several of these are correlated. A high beta value and magnetic field, for example,
results in higher density and smaller dimensions. Let us classify today’s most important
research paths in alternative magnetic confinement by beta value and magnetic field.
High beta value, low magnetic field : Reversed-field pinch (RFP), spherical tokamak.
Low beta value, high magnetic field : Compact tokamak (Ignitor, Riggatron).
Low beta value, low magnetic field: Stellarator.
175
Figure 8.4: The stellarator Wendelstein VII-AS in Garching, Germany. (Figure credit: EURO-
fusion)
A toroidal beta of 20 % has been measured in the first experimental device START
(Small Tight Aspect Ratio Tokamak) in Culham, England. Some other data for this
experiment are R/a = 0.3 m/0.25 m, Ip ≤ 250 kA, Bϕ ≤ 0.5 T, τ ≈ 50 ms and τE ≈ 2
ms.
The high energy confinement and the high beta value seem very promising and larger
machines like MAST (Culham) and (Princeton, USA) are now evaluating possible plasma
parameters and scalings. The problems associated with spherical tokamaks are mainly
design of blanket and cooling of the central conductor.
8.2.3 Stellarator
One of the earliest principles for magnetic confinement is the stellarator, which was men-
tioned already in the introduction. Since the plasma shape does not have toroidal sym-
metry, it has recently become possible to computationally design the magnetic field coils
with sufficiently low field errors and high confinement. The stellarator is now considered
the most important alternative to the tokamak.
A major advantage is that no toroidal plasma current is needed. This enables the con-
figuration to be driven continuously using negative ion beam heating, possibly combined
with ECH. The modular coils are superconducting. The configuration also features low
orbit losses and low neoclassical diffusion.
Like the tokamak, the geometry is toroidal. The rotational transform is maintained
by a deformation fo the torus (early experiments showed poor confinement) or by having
the external coils form a spiral. Almost the entire confining magnetic field is generated
by external coils. The resulting plasma cross section will not be axisymmetric. The
stellarator is thus a typical three-dimensional configuration (see figure 8.4).
Since the magnetic field in the stellarator is mainly maintained by external coils, the
beta value may not be too high, in which case it will cause internal currents break-
ing the topology. Consequently, a stellarator must be physically large. The dimensions
R/a = 20 m/(1–2) m are realistic, with plasma parameters βϕ = 0.02 och Bϕ = 5 T.
176
Particle motion in a stellarator is divided into three cathegories:
* circulating orbits, where the particles travel toroidally without being mirrored
* “helically trapped” particles, mirror reflected in the spiraling field
* “toroidally traped” particles, moving in banana orbits as in tokamaks
Depending on how the rotational transform is created the stellarator configurations
may be divided into different sub-types. Examples are the classical stellarator, the tor-
satron, the helias (like the Wendelstein VII-X experiment in Greifswald, Germany) and
the heliotron (like the large active Japanese experiment Large Helical Device; LHD).
Much stellarator research is focused on the possibility of reaching high beta values,
hopefully 5–10 %. Transport studies, however, are difficult to perform because of the 3-D
geometry. Particle and energy transport are dominated by the way particles are helically
or toroidally captured (and by unwanted variations, so called ripples) in the magnetic field
and by the collision frequency. The latter influences how the particles are scattered in
the different orbits, which are characterized by loss frequencies. In axially asymmetrical
systems electrons and ions may diffuse independently of each other.
In general it can be said that the energy confinement time has a similar scaling de-
pendence as for the tokamak, but with a certain deviation in the low collision range,
where the scattering of helically trapped particles enhances the diffusion. The so-called
advanced stellarators aim to suppress this mechanism.
Summarizing, we conclude that the stellarator has certain advantages compared to
the tokamak. The stationary magnetic field simplifies the design of the magnetic coils
(the superconducting coils do not require pulsing, and thus their energy does not need to
be stored). The high aspect ratio in combination with the absence of a transformer core
and the modular design simplifies construction and makes the stellarator interior easily
accessible. The continuous operation allows simpler blanket design (and less wear). The
large wall surface results in less wall load and longer operation time before exchange of
wall material.
177
Figure 8.5: Experimentally measured magnetic fields in an RFP.
Figure 8.6: The toroidal magnetic field in an RFP changes direction in the outer region.
and thus more economical fusion devices. In addition, already early experiments showed
that high beta values could be sustained in the configuration.
The main problems for the conventional RFP is, at present, the quite low energy
confinement and the need for wall stabilization. We will now discuss these issues in some
detail.
MHD instabilities in an RFP can be divided into five main categories:
• m = 0 modes; are resonant at the reversal surface, where the axial field changes
sign.
• Internally resonant m = 1 modes; have n < 0 and are resonant inside the reversal
surface.
• Externally resonant m = 1 modes; have n > 0 and are resonant outside the reversal
surface.
• Internally non-resonant m = 1 modes; have n < 0 and a pitch similar to that of the
equilibrium magnetic field at the magnetic axis.
• Externally non-resonant m = 1 modes; have n > 0 and a pitch similar to that of
the magnetic field just outside the plasma boundary.
In a conventional RFP, the plasma is surrounded by an adjacent conducting wall with an
L/R time (“shell time”) being longer than the lifetime of the plasma pulse. This stabilizes
all modes except those of category 2) above, the internally resonant modes. The latter
modes are also called dynamo modes, since they are responsible for the plasma’s constant
relaxation towards the minimum energy state.
The pulse length in an RFP reactor must be greater than the shell time. Because of
this, it is important to experimentally demonstrate that the RFP configuration can be
maintained with a resistive shell. This is a field of intense experimental and theoretical
research.
178
Figure 8.7: The RFP Extrap T2R at the Alfvén laboratory, KTH.
One possibility is to rotate the plasma, so that the magnetic field constantly “sees” a
new wall to penetrate. Unfortunately, this requires a rotation velocity close to the Alfvén
velocity, which is unrealistic. So, even if the plasma is rotating, it is unstable to dynamo
modes, as well as to internal and external non-resonant modes (growing on a characteristic
L/R time scale). These modes are called resistive shell modes.
In the Extrap T2R RFP at the Alfvén laboratory at KTH in Stockholm (see figure
8.7) the effect of resistive shell modes on the plasma dynamics is investigated. Previous
experiments at the RFP HBTX-1C in Culham, England, showed that non-resonant, non-
rotating resistive shell modes have a strong, undesired effect on the plasma behaviour.
The Extrap T2R experiment at KTH is characterized by the following parameters:
R/a = 1.24 m/0.183 m, Ip = 100–150 kA, βp = 0.05–0.10 and < n > = 0.5–1.0×1019
m−3 .
From the discussion above it follows that the resistive shell modes must be reduced
through some sort of active feedback. This is similar to the case for the advanced tokamak,
in which the plasma current is continuously generated from an external source in combi-
nation with internal, pressure gradient driven currents. The enhanced current results in
higher beta value and more compact dimensions that those of the conventional tokamak.
In the advanced tokamak, there is only one resistive shell mode, which makes the design
of the feedback coils relatively simple. In the RFP, on the other hand, a large number of
modes with different pitch must be stabilized, which gives a more difficult problem, being
intensively studied at present.
In the conventional RFP, the confinement time is limited due to the magnetic field
fluctuations caused by the dynamo. A stochastization of the magnetic field in the inner
region of the plasma results, with accompanying radial heat transport along the field lines.
At the boundary, heat is transported by perpendicular ion heat conduction along a rela-
tively sharp temperature gradient. This also occurs in the presence of a nearby conducting
179
Figure 8.8: The principle for inertial fusion. (Figure credit: http://physci.llnl.gov/divisions
/fep/ifep/IFEP.html)
wall. Computer simulations have given the relations βp ≈ 3.4 × 10−5 a−0.2 (I 2 /N )−0.40 and
τE ≈ 0.77a1.4 (I 2 /N )0.34 , where N is the line density. These approximatively correspond
to optimized transport in the conventional RFP. The earlier reactor study TITAN uses
the parameters R/a = 3.8 m / 0.6 m, Ip = 18 MA, βp = 0.2 and < n >= 9 × 1020 m−3 ,
resulting in a too optimistic beta value according to the scaling above. In the advanced
RFP, in which the dynamo activity has been reduced by the method of pulsed poloidal
current drive (PPCD), confinement is strongly enhanced. Experiments and computer sim-
ulations show that the energy confinement can be improved several times by this method.
Experimental and theoretical RFP research is mainly carried out in Italy, USA, Japan,
China and Sweden. It is still too early to evaluate the reactor potential of the RFP. Topics
that must be addressed are, as mentioned, resistive wall modes, current profile control
techniques to master high magnetic fluctuations and weak energy confinement, as well as
heating and current drive.
180
ns. Thus τE ≈ τ ≈ R/vT , where R is the grain radius and vT is the thermal velocity. The
radius of the grain is then 0.01 m.
In reality, the DT grains are not made of massive fuel. In the above example a massive
grain would require around 100 MJ to heat up to 10 keV, an energy input which is not
possible with the lasers available today.
The grain is instead made as a hollow sphere. The innermost layer is frozen DT.
Outside is added a so called pusher and the outermost layer is some easily evaporated
material. When the grain is irradiated, the surface is evaporated. The reaction force
compresses the DT fuel to around 1000 times the normal density. The laser energy 1–10
MJ is incident during 1–10 ns, corresponding to an enormous power. The associated
Lawson criterion for inertial fusion is usually written ρR ≥ 10 kg/m2 .
It is not straightforward to carry out the process described above. The major obstacle
has been Rayleigh-Taylor instabilities which easily arise during the compression. By
fabricating very smooth grains and developing homogeneous laser and ion beams, most
stability problems now seem to be solved. The figure below shows a computer simulation
of the time development of Rayleigh-Taylor instabilities.
Figure 8.9: Contours of density and temperature in a computer simulation of inertial fusion.
(Figure credit: FUSION EXPO-EUROfusion)
To achieve homogeneous radiation of the grain, researchers do not only use “direct
drive” as above. The grain may also be surrounded by a gold chamber. When the walls
of the chamber are radiated, they emit high energy bremsstrahlung which compresses and
heats the grain in a process called “indirect drive”. The process is imaged in figure 8.10.
A reactor based on inertial fusion will necessarily be large, because of the large beam
sources. The reactor chamber will have a diameter of about 10 m. The inner walls will
be covered with a fluid Li curtain which captures neutrons and transforms Li to T. The
reactions must occur with a frequency of about 10 Hz. So far, it has not been decided
if direct or indirect radiation will be used, or if laser or ion beams are the final solution.
181
Figure 8.10: The principle for “indirect drive”. (Figure credit: Lawrence Livermore National
Security, LLC, Lawrence Livermore National Laboratory, and the Department of Energy, USA)
(He++ ion beams give higher efficiency in the beam source and better coupling to the
plasma). Reactor studies indicate a rather low COE; 0.05–0.07 Euro/kWh(electricity).
182
Chapter 9
In the European “fast-track schedule” for realization of fusion power, the goal is to have a
demonstration reactor (DEMO), demonstrating electricity to the grid, by the mid-century.
In a reactor, a larger plasma volume than that in today’s tokamaks is required, to provide
the increase in confinement time necessary for power break-even. Further, a reactor must
generate it’s magnetic field for much longer time periods than current tokamaks, which
necessitates the use of superconducting magnetic field coils.
In this chapter we will just scratch the surface of fusion technology. As fusion physics
problems successively become solved or effectively handled, fusion technology issues be-
comes increasingly important in the research and development of fusion. ITER will cer-
tainly provide an important test bench for developing many technological solutions before
we may embark on DEMO and commercial fusion. A new fusion materials test facility
(IFMIF-DONES), will help to optimise the materials that should be used for plasma-
facing components, blanket and the reactor structure.
The next step tokamak device, ITER, will integrate some of the technologies essential
for a fusion reactor, such as superconducting magnets, and test tritium breeding blanket
designs. Many new technologies, from superconducting magnets to novel materials are
being developed for ITER. The ITER device comprises the tokamak itself and associated
systems for heating, fuelling, exhaust of waste and gas, control, and diagnostic measure-
ments. ITER, like JET has a vertical D-shaped plasma and a lower divertor system. The
divertor is a critical component and the main area where the plasma will contact the
material wall. ITER will use low temperature (4 K) superconducting magnets for both
its 18 toroidal and 6 poloidal field coils as well as the central solenoid coil. The toroidal
field coils generate a magnetic field of 5.3 T in the plasma region. The whole ITER device
is enclosed in a 70 K cryostat which helps insulate the superconducting magnetic coils.
External heating systems provide input heating power of about 50 MW using neutral
beam injection and radio-frequency electromagnetic wave heating. The inner surfaces of
the vacuum vessel are covered with blanket modules. These modules will provide shield-
ing from the high-energy neutrons produced by the fusion reactions and will be used for
testing tritium breeding concepts.
The ITER plasma containment vessel will be more than twice as large and 16 times
as heavy as any previously manufactured fusion vessel. The internal structure of a fusion
reactor will become radioactive during operation due to neutron radiation and the presence
of tritium. Remote-handling systems are therefore vital to be able to replace components,
such as the divertor and eventually breeder blanket modules, inside the machine. In a
fusion power plant, cryogenic systems are used to remove impurities from the plasma, cool
the superconducting coils to allow them to operate, separate the waste gases into their
183
1 5
6
2
7
3
Figure 9.1: ITER. 1) Central solenoid coil, 2) Toroidal field coil, 3) Poloidal field coil, 4)
Diagnostics, 5) Blanket module, 6) Vacuum vessel, 7) Cryostat, 8) External heating sys-
tems 9) Divertor (From ITER - Uniting science today global energy tomorrow, European
Commission, 2007, p. 14)
different individual components for disposal or recycling for fuel, provide the cooling for
radio-frequency heating sources, and control the gas pressure of neutral beam systems.
Large-scale vacuum systems are required to ensure an ultra-high vacuum in the large
reactor vessels and to maintain the vacuum surrounding the superconducting magnets in
the cryostat.
Plasma heating systems are essential for obtaining a high-temperature plasma. For
ITER, the fusion reaction would not continue if the plasma was not heated by an external
source, and for both ITER and future power plant operation it is likely that the heating
systems will be an essential tool to ensure stability and control of the plasma. Initially,
three main types of heating systems will be deployed for ITER:
• Ion Cyclotron Resonance Heating - in this system, ions in the plasma are heated
by electromagnetic waves with a frequency of 30-50 MHz. The main issues concern
how to couple the intense radiation to the plasma and what effect this has on the
performance of the plasma.
• Electron Cyclotron Resonance Heating - here the electrons in the plasma are heated
by electromagnetic waves with a frequency of 100-200 MHz. This radiation has the
advantage that it can be transmitted through air which simplifies the design and
means that the source can be far from the plasma.
• Neutral beam injection - in this system, charged fusion fuel particles are accelerated
to very high speed corresponding to a kinetic energy of 1 MeV and neutralized so
184
Figure 9.2: The reactor. (Figure credit: EUROfusion-JET)
that high-energy neutral particles can pass through the magnetic field and enter the
fusion plasma. As a result, the plasma is heated by the transfer of kinetic energy.
185
Figure 9.3: Limiter and divertor. (Figure credit: EUROfusion-JET)
186
9.1.2 First wall
The part of the wall being closest to the plasma is called the first wall. It has a protecting
layer of low Z material and inside this some material with very high heat conductivity.
Copper and steel are less suitable materials in a reactor because of their relatively long
lived radioactive isotopes. Preferable materials are low-activity materials such as silicon
carbide, titanium or aluminium. The first wall must be able to withstand high neutron
fluxes. Most materials are destroyed with time, primarily through embrittlement, swelling
and endurance failure (especially due to disruptions). A rule of thumb is that if the
number of dpa (displacements per atom) 1, the crystal structure may become severely
damaged. In a fusion reactor, certain components might be exposed to 200 dpa. Thus,
the neutron flux is a measure of the lifetime of the first wall and the blanket. Stainless
steel is considered to withstand 15 MWy/m2 , but the requirement is somewhat higher.
The conventional tokamak has a relatively large wall surface / plasma volume ratio, which
prolongs the time before a change of the first wall/blanket is needed. In the ARIES-I study,
the calculated neutron load is 2.5 MW/m2 , which corresponds to 6 years of operation using
stainless steel walls. A more compact reactor would require more frequent changes of wall
material. A solution is here to reduce the cost of the FPC to about 10 % of the reactor’s
total cost, making it economically viable to exchange the entire FPC.
187
The first reaction is exothermal and gives a positive contribution to the energy balance
in the reactor. The neutron balance is also important, of course; each fusion reaction must
result in at least one T atom through these reactions, or the supply of T will eventually
run out.
No experimental device is running yet, but the most popular design solution is to
either let the cooling agent in part or in full be a Li salt, or to use solid Li2 O. The total
inventory of T in a reactor is 1–40 kg, depending on the type of reaction.
9.1.4 Shield
The magnetic field coils, generally superconducting, must be protected from neutrons
and heat. The shield, which is just outside the blanket, takes care of this while also
contributing to T containment.
The combined thickness of the blanket and the shield must be 1.5 m for the DT reaction
and 1.0 m for the DD reaction. The thickness only weakly depends on the material; it is
the same for 316 SS as for V and Ti.
A result of this is that the coils must be placed some distance away from the plasma,
which results in a lower field strength in the plasma. Apparently, the reactor becomes
larger and more costly when we take neutrons into account. Compact reactors, where
water cooled copper coils are inside the blanket have been proposed, but this causes
problems for the neutron balance (absorption in the coils). It is obvious that the use of
advanced fuels can mean a lot for reactor economy.
188
region where the strength of the magnetic field is low enough for this material to be used.
However, the position also means that replacing poloidal field coils will be very difficult,
so each coil will be designed with redundant turns so that any faults can be isolated to
ensure that operation of ITER continues unhindered. Once energized, the magnets can
work continuously with very high efficiency. As these magnets run at liquid helium tem-
perature it is necessary to operate them in vacuum to prevent heat in the atmosphere
from boiling off the helium.
9.2.1 Tritium
Tritium beta decays into a neutrino and a neutron (activity: 9600 Ci/g) with a half-life of
12.3 years.. The radiation can cause damage in biological tissue if T enters into the tissue.
The most common form of T when entering tissue is T2 O and THO, which replaces H2 O.
At present, three different reactor options are studied in the EU, based on different
construction materials. The fusion power is 3000 MW. In all these cases, passive safety is
important. No serious accidents are allowed to happen if the safety systems fail for some
reason. The requirement is that the general public should not be exposed to a dose higher
189
than 50mSv/y. Thus, the most serious accident must lead to no more than 50mSv doses,
and the present designs fulfils this requirement.
Since the T inventory (except in DHe3 reactors) is more than 50 g, a special confine-
ment of the entire reactor is needed. The vacuum vessel is separated from the surroundings
by both a vacuum wall and an outer wall. There is also a requirement on a minimum dis-
tance from densely populated areas. It is estimated that tritium confinement is disturbed
10−7 times/y, and this must be handled by the design. Evacuation is not considered to
become necessary.
We should note that also DD and DHe3 reactions generate T (around 0.35 and 0.005
kg/day, respectively), but this is mainly consumed in the plasma.
The handling of T dominates the external costs for fusion, which are low in comparison
with other energy sources. The external costs may be further reduced if the production
of C14 is avoided. This radioactive isotope (half-life 5730 y) stems from reactions with
neutrons and nitrogen, the latter being present in certain types of steel. To avoid the
production of C14 , other materials with less content of natural nitrogen could be used.
Alternatively, the materials could be enriched with respect to the isotope N15 (reduced
activation ferritic steels), resulting in lower levels of C14 .
If air gets into the vacuum vessel, an explosive mix of hydrogen and oxygen can be
formed, and cause release of tritium. The risk is, however, considered very small. Some
designs use an inert gas between the two outer protective walls.
190
Figure 9.4: Remote handling robot inside JET’s vacuum vessel. (Figure: EUROfusion-JET)
caused by helium formation during radiation. Disadvantages are that the material can
become brittle at high working temperatures, 550–6500 C is the upper limit, and that it
affects the coil’s magnetic fields. The best alloy at the moment is 9CrWTaV, which also
is resistant to corrosion. From a radiological point of view, F/M steel is not an optimum
material, as mentioned earlier.
Vanadium alloys. These (e.g. V4Ti4Cr) fulfil most of the requirements above, but are
sensitive to reactions with oxygen. Thus, only liquid lithium can be used as coolant, and
this strongly limits the interest for these materials.
Titanium and chromium alloys are mostly studied for their low neutron activation.
Their brittlement limits the possibilities, but models using powdered metal are promising.
Silicon carbide.This composite material (SiC/SiC) is definitely the most promising
when it comes to compromises between heat resistivity (withstands up to 11000 C), dura-
bility and waste. The research is focused on making silicon carbide a better thermal and
electrical conductor, less sensitive to radiation and cheaper.
The case for all these materials is that the testing is difficult since powerful enough
radiation sources are lacking.
191
Appendix A
192
The conversion from the energy unit eV to Joule is
1 eV = 1.602177 × 10−19 J
Temperature in the SI-system is measured in Kelvin (K). In the field of plasma physics
we refer to temperatures by their equivalent energies in eV. In one species of a plasma at
a temperature T , for example ions, the ion energy density is (3/2)nkB T , where n is the
ion particle density, and kB is Boltzmann’s constant. We associate an average energy of
(1/2)kB T to each of the three degrees of freedom for a plasma particle. This equivalence
of temperature and energy is used in plasma physics to define the temperature in electron
volts. A temperature in electron volts of TeV = 1 eV is equal to the temperature in Kelvin
(TK ), for which the energy kB TK is 1 eV. The conversion factor is thus obtained from the
relation
eTeV = kB TK
The conversion from the temperature unit eV to Kelvin is
In this book we use the symbol T referring to TeV . We have also used the convention to
write out the electron charge e explicitly in the expression eT , which means that T has
the SI-unit Volt (V).
a · (b
b × c) = b · (cc × a) = c · (a
a × b)
a × (bb × c) = (a
a · c)bb − (a
a · b)cc
∇ × ∇Ψ = 0
∇ × (∇ × a) = ∇(∇ · a) − ∇2 a
∇ · (aa × b) = b · (∇ × a) − a · (∇ × b)
∇ × (a a × b) = a(∇ · b) − b(∇ · a) + (b b · ∇)a
a − (a
a · ∇)b
b
193
Figure A.1: Cylindrical coordinate system
Line element
dl = drr̂r + rdθθ̂θ + dzẑz
Volume element
dV = rdrdθdz
Gradient, divergence, curl
∂Ψ 1 ∂Ψ ∂Ψ
∇Ψ = r̂r + θ̂θ + ẑz
∂r r ∂θ ∂z
1 ∂ 1 ∂aθ ∂az
∇·a = (rar ) + +
r ∂r r ∂θ ∂z
1 ∂az ∂aθ ∂ar ∂az 1 ∂ 1 ∂ar
∇×a = − r̂r + − θ̂θ + (raθ ) − ẑz
r ∂θ ∂z ∂z ∂r r ∂r r ∂θ
A.6 Electromagnetics
Low frequency approximation of Maxwell’s equations in differential form
∇·B=0
1
∇·E= ρ
ε0
∂B
∇×E=−
∂t
∇ × B = µ0 j
Ampere’s law in integral form is useful for calculating the magnetic field from a current
where there is for example cylinder symmetry. The following form is obtained by applying
Stoke’s theorem I Z
B · dl = µ0 j · n̂
ndS = µ0 Iencl
C S
194
where Iencl is the current enclosed by the curve C. One example is the calculation of the
azimuthal magnetic field from the axial plasma current in the case of circular cylinder
geometry, which then is direcly obtained from Ampere’s law as Bθ = µ0 Iz /(2πr). When
modelling high frequency plasma waves, and electromagnetic waves in vacuum, the exact
form of Ampere’s law is used
1 ∂E
∇ × B = µ0 j +
c2 ∂t
195
Index
196
coal, 4 drift wave instability, 91
cold fusion, 21, 29 dynamo, 145
collision cross section, 115 dynamo mode, 146
collision frequency, 23, 86
collision time, 39, 74 earth, 21
compact tokamak, 148 effective charge, 33
conductivity tensor, 119 efficiency, 25
confinement, 68, 70 eigenvalue problem, 54
continuity equation, 41 Einstein, 19
continuous spectrum, 137 electricity, 4
conventional RFP, 145 electromagnetic fluctuation, 89
conventional tokamak, 144 electromagnetic radiation, 100, 125
cooling water, 25 electron cyclotron radiation, 103, 134
corona equilibrium, 104 electron cyclotron resonance, 122
Coulomb collision, 73 electron cyclotron resonance heating, 124
Coulomb force, 19 electron temperature, 130
Coulomb logarithm, 32, 39 electrostatic deflection, 140
critical disturbance, 159 electrostatic fluctuation, 89
cross section, 21 electrostatic turbulence, 90
curvature drift, 37 elongation, 66
cut-off, 121 embrittlement, 156
cyclotron damping, 123 emission, 100
cyclotron frequency, 35 emittance, 134
cyclotron radiation, 102 endurance failure, 156
energy backscattering coefficient, 106
Debye length, 37 energy confinement time, 25, 31, 33
deceleration time, 39 energy consumption, 4
deexcitation, 104 energy density, 21
deflection angle, 74 energy equation, 42
deflection magnet, 117 energy exchange time, 40
DEMO, 143 energy principle, 54, 58
density fluctuation, 90 energy situation, 3, 4
desorption, 107 energy sources, 13, 17
deuterium, 3, 19, 21 energy spectrum, 136
diagnostic method, 126 environment, 5
diamagnetic drift velocity, 92 equation of motion, 42
dielectric tensor, 119 equilibrium, 45
diffusion coefficient, 73 equilibrium equation, 45
direct transfer, 23 evaporation, 110
discharge, 45 excitation, 104
displacement, 52 excitation balance, 104
disruption, 66, 157 exotic fusion reactions, 27
distribution function, 38 exotic reaction, 22, 27
disturbance of cooling, 157, 160 external heating, 33
disturbance of superconducting coils, 157 external kink, 61
divertor, 154, 155 external mode, 62
Doppler effect, 130
fast magnetoacoustic wave, 57
drift, 36
fast neutral particle, 139
drift wave, 91
Fick’s law, 73
197
field emission, 109 hydrogen, 7, 14
first wall, 154, 155 hydropower, 11
fission, 16
fission energy, 19 ideal instability, 54
fission reactor, 4, 15 ideal MHD, 43
flood, 5 ignition, 26
fluid model, 34, 41 impurity ions, 103
fluid velocity, 42 incidence parameter, 74
flux function, 47 incompressibility, 44
flux surface, 47 inertial confinement, 20, 142
Fokker-Planck equation, 39 inertial fusion, 150
Forsmark, 16 initial value problem, 55
fossil fuel, 3, 4 instability, 54
fossil fuel resources, 4 interferometer, 127
fossil gas, 5, 6 intermediate level waste, 15
FPC, 148 internal kink, 63
free electron, 100, 125 internal mode, 61
fuel cell, 14 ion cyclotron resonance, 122
fusion, 19 ion cyclotron resonance heating, 124
fusion cross section, 23 ion process, 126
fusion energy, 18 ion source, 117
fusion power, 18, 50 ion temperature, 139
fusion product, 18 ionization, 104, 107
fusion reaction, 19, 21 ionization coefficient, 104
fusion reactivity, 24 island width, 94
fusion reactor, 4, 19, 154 ITER, 27, 142
fusion-fission, 151 kinetic model, 38
gamma radiation, 21 kinetic theory, 34
Gamow, 24 kink, 61
gas, 4 kink instability, 60
gas target, 112 Kirchoff’s law, 134
geothermal energy, 9
L wave, 120
geothermal power, 12
Landau damping, 122
glow discharge cleaning, 108
Larmor radius, 35
Grad-Shafranov equation, 48, 49
Lawson, 23, 24
gravitational force, 36
Lawson criterion, 24, 31
greenhouse effect, 5
Lawson diagram, 26
Greenwald, 67
left circular polarized wave, 120
growth rate, 59
LHD, 149
guiding-centre drift, 78
limiter, 110, 154, 155
guiding-centre orbit, 78
line density, 33, 129
heat conductivity, 77 line radiation, 103
heat diffusivity, 77 linear stability, 54
heating, 26 linearize, 52
high energy neutral particle, 138 LMJ, 27, 150
high level waste, 15 low level waste, 15
hybrid reactor, 4, 142, 151 lower hybrid resonance, 122
hydrodynamic equation, 42 lower hybrid resonance heating, 124
198
magnetic coils, 154, 157 Ohmös law, 43
magnetic confinement, 20 ohmic heating, 113
magnetic deflection, 140 ohmicly heated plasma, 33
magnetic divertor, 112 Ohms law, 40
magnetic field surface, 47 oil, 4
magnetic island, 93 one-fluid velocity, 43
magnetic mirror field, 82 optical depth, 134
magnetic mirror force, 82 optically thick plasma, 135
magnetic moment, 82 optically thin plasma, 135
major radius, 32 oscillatory motion, 84
marginal stability, 55 overlapping, 97
martensitic steel, 158
mass density, 43 particle, 34
MAST, 148 particle confinement time, 76
Maxwell distribution, 23, 24 particle drift, 36
Maxwellös equations, 43 particle dynamics, 35
mean free path, 23 particle flux, 125
Mercier, 65 Pfirsch-Schlüter transport, 81
MHD wave, 57 phase shift, 127
MHD waves, 56 pinch, 2
mixing length approximation, 93 Planck’s radiation law, 134
mode conversion, 118 plasma, 23
modulation, 129 plasma current, 32
multichannel pulse height analyzer, 137 plasma diagnostics, 125
muon-catalyzed fusion, 21, 28 plasma dispersion relation, 118
plasma electron density, 127
NbSn super conductor, 157 plasma frequency, 56
NbTi super conductor, 157 plasma model, 34
neoclassical diffusion coefficient, 87 plasma oscillation, 56
neoclassical transport, 71, 78 plasma parameter, 34
neptunium, 15 plasma wall interaction, 105
neutral beam heating, 114 plateau regime, 87
neutral beam injection, 114, 117 plutonium, 15
neutral fraction, 117 poloidal β, 50
neutral particle, 138 poloidal beta, 32
neutral particle analysis, 137 poloidal beta value, 31
neutral particle flux, 139 poloidal direction, 31
neutralization chamber, 117 pressure driven instability, 65
neutralize, 117 pressure gradient, 65
neutron, 21 pressure tensor, 42
neutron flux, 156
neutron free reaction, 22 quadrature detection, 129
NIF, 27, 150 quality parameters, 31
NSTX, 148 quantum mechanics, 24
nuclear accidents, 4
R wave, 120
nuclear energy, 15
radiation parameter, 104
nuclear power, 18
radiative power, 132
O-wave, 120 radio-frequency heating, 117
ocean wave, 11 radio-frequency waves, 117
199
radioactivity, 19 scattering volume, 133
random walk model, 72 separatrix, 94, 155
rate coefficient, 107 shear, 65
rational surfaces, 54 shield, 154, 157
Rayleigh-Taylor instability, 60, 151 silicon carbide, 159
reaction kinematics, 23 silicon detector, 137
reaction rate, 23 slow magnetoacoustic wave, 57
reaction velocity, 21 soft X-ray radiation, 137
reactivity, 23 solar cells, 10
recombination, 104 solar energy, 9
recombination coefficient, 104 solar power, 9
recycling, 106 solenoid coil, 126
recycling coefficient, 106 Solov’ev, 49
reduced ion energy, 106 spallation, 16
reference beam, 127 spectral function, 132
reflection, 121 spectral function half-maximum width, 133
refractive index, 119, 125 spectroscopic measurement, 126
relaxation time, 39 spherical tokamak, 147
remote collision, 74 Spitzer, 1
renewable energy, 9, 12 sputtering, 108
renewable energy source, 4, 12 sputtering yield, 108
resistive g-mode, 64 stability, 58
resistive instability, 54, 64 stability diagram, 64
resistive MHD, 75 stabilizing effect, 54
resistive mode, 65 stabilizing magnetic fields, 51
resistive shell mode, 147 stabilizing term, 60
resistivity, 40 stainless steel 316 SS, 156
resonance, 52, 121 star, 19
resonance layer, 135 START, 147
resonance surface, 121 static equilibrium, 45
resonant surface, 93 stationary state, 45
rest gas, 103 stellarator, 1, 2, 37, 142, 148
reversed toroidal field, 145 step length, 73
reversed-field pinch, 50, 142, 144 stochastic magnetic field, 97
RFP, 50, 144 storage of carbon dioxide, 5
Riggatron, 148 strong interaction, 25
right circular polarized, 120 sun, 1, 3, 21
rippling mode, 64 Suydam, 65
Rogowski coil, 126 swelling, 156
rotation transform, 53 Szilard, 1
Rutherford’s scattering law, 74
target chamber, 112
safety factor, 52 tearing mode, 64, 67
safety issue, 157 temperature, I
sausage instability, 60 test function, 59
sawtooth, 53, 63 thermal emission, 109
scattered wave, 126 thermalized plasma, 23
scattering, 126 thermodynamic equilibrium, 51
scattering time, 39 thermonuclear fusion, 1, 20, 31
200
Thomson scattering, 130
titanium alloy, 159
tokamak, 2, 32
tokamak stability, 61
toroidal beta value, 50
toroidal direction, 31
transmission, 126
transmutation, 15
trapped particle, 82
trapping, 114
triple product, 141
tritium, 19, 21, 157
tritium handling, 157
tritium., 3
Troyon factor, 68
turning point, 85
wave energy, 11
wave power plant, 11
wave propagation, 118
wave source, 129
wave-particle interaction, 118
wave-particle resonance, 122
wavelength shift, 130
Wendelstein, 149
wind power, 10
wind power plant, 10
z-pinch, 46, 50
201