Reactions

Reactions#

The following content gives some background to the reactions implemented in Hermes-3. Please see the Post-processing section for related diagnostics.

See also

For a more detailed description of how reaction source terms are handled in the code, please see the ‘Reactions’ section of the developer manual.

Reaction basics#

The formula for the reaction is used as the name of the component. This makes writing the input file harder, since the formula must be in the exact same format (e.g. h + e and e + h won’t be recognised as being the same thing), but makes reading and understanding the file easier.

To include a set of reactions, it is probably easiest to group them, and then include the group name in the components list

[hermes]
components = ..., reactions

[reactions]
type = (
        h + e -> h+ + 2e,  # ionisation
        h+ + e -> h,    # Radiative + 3-body recombination
       )

Note that brackets can be used to split the list of reactions over multiple lines, and trailing commas are ignored. Comments can be used if needed to add explanation. The name of the section does not need to be “reactions”, and multiple components could be created with different reaction sets. Be careful not to include the same reaction twice.

When reactions are added, all the species involved must be included, or an exception should be thrown.

Transfer channels#

Reactions typically convert species from one to another, leading to a transfer of mass momentum and energy. For a reaction converting species \(a\) to species \(b\) at rate \(S\) (units of events per second per volume) we have transfers:

\[\begin{split}\begin{aligned} \frac{\partial}{\partial t} n_a =& \ldots - S \\ \frac{\partial}{\partial t} n_b =& \ldots + S \\ \frac{\partial}{\partial t}\left( m n_a u_a\right) =& \ldots + F_{ab} \\ \frac{\partial}{\partial t}\left( m n_a u_a\right) =& \ldots + F_{ba} \\ \frac{\partial}{\partial t}\left( \frac{3}{2} p_a \right) =& \ldots - F_{ab}u_a + E_{ab} - \frac{1}{2}mSu_a^2 \\ \frac{\partial}{\partial t}\left( \frac{3}{2} p_b \right) =& \ldots - F_{ba}u_b + E_{ba} + \frac{1}{2}mSu_b^2 \end{aligned}\end{split}\]

where both species have the same mass: \(m_a = m_b = m\). In the pressure equations the \(-F_{ab}u_a\) comes from splitting the kinetic and thermal energies; \(E_{ab}=-E_{ba}\) is the energy transfer term that we need to find; The final term balances the loss of kinetic energy at fixed momentum due to a particle source or sink.

The momentum transfer \(F_{ab}=-F_{ba}\) is the momentum carried by the converted ions: \(F_{ab}=-m S u_a\). To find \(E_{ab}\) we note that for \(p_a = 0\) the change in pressure must go to zero: \(-F_{ab}u_a + E_{ab} -\frac{1}{2}mSu_a^2 = 0\).

\[\begin{split}\begin{aligned} E_{ab} =& F_{ab}u_a + \frac{1}{2}mSu_a^2 \\ =& - mS u_a^2 + \frac{1}{2}mSu_a^2\\ =& -\frac{1}{2}mSu_a^2 \end{aligned}\end{split}\]

Substituting into the above gives:

\[\begin{split}\begin{aligned} \frac{\partial}{\partial t}\left( \frac{3}{2} p_b \right) =& \ldots - F_{ba}u_b + E_{ba} + \frac{1}{2}mSu_b^2 \\ =& \ldots - mSu_au_b + \frac{1}{2}mSu_a^2 + \frac{1}{2}mSu_a^2 \\ =& \ldots + \frac{1}{2}mS\left(u_a - u_b\right)^2 \end{aligned}\end{split}\]

This has the property that the change in pressure of both species is Galilean invariant. This transfer term is included in the Amjuel reactions and hydrogen charge exchange.

Adjusting reactions#

The reaction rates can be adjusted by a user-specified arbitrary multiplier. This can be useful for the analysis of the impact of individual reactions. The multiplier setting must be placed under the neutral species corresponding to the reaction, e.g. under [d] when adjusting deuterium ionisation, recombination or charge exchange. The multiplier for the fixed fraction impurity radiation must be placed under the impurity species header, e.g. under [ar] for argon. This functionality is not yet currently implemented for helium or neon reactions.

Setting

Specified under

Reaction

K_iz_multiplier

Neutral species

Ionisation rate

R_ex_multiplier

Neutral species

Ionisation (excitation) radiation rate

K_rec_multiplier

Neutral species

Recombination rate

R_rec_multiplier

Neutral species

Recombination radiation rate

K_cx_multiplier

Neutral species

Charge exchange rate

R_multiplier

Impurity species

Fixed frac. impurity radiation rate

The charge exchange reaction can also be modified so that the momentum transfer channel is disabled. This can be useful when testing the impact of the full neutral momentum equation equation compared to purely diffusive neutrals. A diffusive only model leads to all of the ion momentum being lost during charge exchange due to the lack of a neutral momentum equation. Enabling neutral momentum introduces a more accurate transport model but also prevents CX momentum from being lost, which can have a significant impact on the solution and may be difficult to analyse. Disabling the momentum transfer channel allows you to study the impact of the improved transport only and is set as:

[hermes]
components = ..., c, ...

[reactions]
no_neutral_cx_mom_gain = true

Fixed fraction radiation model#

In the fixed fraction radiation model, the impurity density is assumed to be a constant fraction of the electron density. The impurity only affects the plasma solution through radiation, which is calculated using a “cooling curve”, which is a function of the radiation in terms of temperature in units of \(Wm^{-3}\). More information on cooling curves can be found in literature, e.g. A. Kallenbach PPCF 55(12) (2013).

To use this component you can just add it to the list of components and then configure the impurity fraction:

[hermes]
components = ..., c, ...

[c]
type = fixed_fraction_carbon
fraction = 0.05   # 5% of electron density
diagnose = true   # Saves Rc (R + section name)

This will create a diagnostic variable named Rc containing the radiation source in \(Wm^{-3}\) where c is the species name.

Several ADAS rates are provided: nitrogen, neon, argon, krypton, xenon and tungsten. The component names are fixed_fraction_carbon, fixed_fraction_nitrogen, fixed_fraction_neon, fixed_fraction_argon, fixed_fraction_krypton, fixed_fraction_xenon and fixed_fraction_tungsten.

Each rate is in the form of a 10 coefficient log-log polynomial fit of data obtained using the open source tool radas, except xenon and tungsten that use 15 and 20 coefficients respectively. The \(n {\tau}\) parameter representing the density and residence time assumed in the radas collisional-radiative model has been set to \(1\times 10^{20} \times 0.5ms\) based on David Moulton et al 2017 PPCF 59(6).

Each rate has an upper and lower bound beyond which the rate remains constant. Please refer to the source code in fixed_fraction_radiation.hxx for the coefficients and bounds used for each rate.

In addition to the above rates, there are three simplified cooling curves for Argon: fixed_fraction_argon_simplified1, fixed_fraction_argon_simplified2 and fixed_fraction_argon_simplified3. They progressively reduce the nonlinearity in the rate by taking out the curvature from the slopes, taking out the RHS shoulder and taking out the LHS-RHS asymmetry, respectively. These rates may be useful in investigating the impact of the different kinds of curve nonlinearities on the solution.

There is also a very simple carbon radiation function fixed_fraction_hutchinson_carbon which is in coronal equilibrium, using a simple formula from I.H.Hutchinson Nucl. Fusion 34 (10) 1337 - 1348 (1994):

\[L\left(T_e\right) = 2\times 10^{-31} \frac{\left(T_e/10\right)^3}{1 + \left(T_e / 10\right)^{4.5}}\]

which has units of \(Wm^3\) with \(T_e\) in eV.

NOTE:

By default, fixed fraction radiation is disabled in the core region. This represents the fact that the impurity will likely be coronal on closed field lines and feature reduced radiation. This can prevent unphysical MARFE-like behaviour in deep detachment. This behaviour can be disabled by setting no_core_radiation=false in the impurity options block.

Implemented reactions#

Hydrogen and Helium#

Reionisation and recombination reactions for Hydrogen isotopes and Helium are handled by the IznReaction and RecReaction classes. Reaction rates themselves are evaluated by the RateHelper class, using Amjuel datasets by default.

Multiple isotopes of Hydrogen can be evolved and are kept track of using species labels h, d and t. The following could therefore be used to enable ionisation of Deuterium and Tritium in the configuration file:

[hermes]
components = d, t, reactions

[reactions]
type = (
        d + e -> d+ + 2e,  # Deuterium ionisation
        t + e -> t+ + 2e,  # Tritium ionisation
       )

The full set of ionisation and recombination reactions currently available for Hydrogen isotopes and Helium are:

Table 1 Ionisation and recombination reactions for Hydrogen isotopes and Helium#

Reaction

Description
(default data source)

Reference / Notes

h + e -> h+ + 2e

Hydrogen ionisation
(Amjuel H.4 2.1.5)

D.Reiter, K.Sawada and T.Fujimoto (2016)

d + e -> d+ + 2e

Deuterium ionisation
(Amjuel H.4 2.1.5)

t + e -> t+ + 2e

Tritium ionisation
(Amjuel H.4 2.1.5)

h+ + e -> h

Hydrogen recombination
(Amjuel H.4 2.1.8)

Effective rates, including both radiative and 3-body contributions.

d+ + e -> d

Deuterium recombination
(Amjuel H.4 2.1.8)

t+ + e -> t

Tritium recombination
(Amjuel H.4 2.1.8)

he + e -> he+ + 2e

He ionisation
(Amjuel H.4 2.3.9a)

Unresolved metastables

he+ + e -> he

He+ recombination
(Amjuel H.4 2.3.13a)

Unresolved metastables

Energy losses associated with both the ionisation potential energy cost and photon emission during excitation and de-excitation during multi-step ionisation are also included. The default datasets for these calculations are Amjuel reaction H.10 2.1.5 for ionisation and Amjuel reaction H.10 2.1.8 for recombination.

Charge exchange reactions between Hydrogen isotopes are handled by the CXReaction class. By default, rates are calculated using Amjuel H.3 3.1.8, parameterised by an effective temperature scaled to different isotope masses and finite neutral particle temperatures according to (Amjuel p43):

\[T_{eff} = \frac{M}{M_1}T_1 + \frac{M}{M_2}T_2\]

The available charge exchange reactions for Hydrogen isotopes are:

Table 2 Charge exchange reactions for Hydrogen isotopes#

Reaction

Description (default data source)

h + h+ -> h+ + h

Hydrogen charge exchange (Amjuel H.3 3.1.8)

d + d+ -> d+ + d

Deuterium charge exchange (Amjuel H.3 3.1.8)

t + t+ -> t+ + t

Tritium charge exchange (Amjuel H.3 3.1.8)

h + d+ -> h+ + d

Mixed Hydrogen isotope CX (Amjuel H.3 3.1.8)

d + h+ -> d+ + h

“

h + t+ -> h+ + t

“

t + h+ -> t+ + h

“

d + t+ -> d+ + t

“

t + d+ -> t+ + d

“

Lithium#

These rates are taken from ADAS (‘96 and ‘89)

Reaction

Description

li + e -> li+ + 2e

Lithium ionisation

li+ + e -> li+2 + 2e

li+2 + e -> li+3 + 2e

li+ + e -> li

Lithium recombination

li+2 + e -> li+

li+3 + e -> li+2

li+ + h -> li + h+

Charge exchange with hydrogen

li+2 + h -> li+ + h+

li+3 + h -> li+2 + h+

li+ + d -> li + d+

Charge exchange with deuterium

li+2 + d -> li+ + d+

li+3 + d -> li+2 + d+

li+ + t -> li + t+

Charge exchange with tritium

li+2 + t -> li+ + t+

li+3 + t -> li+2 + t+

The implementation of these rates is in ADASLithiumIonisation, ADASLithiumRecombination and ADASLithiumCX template classes:

template<int level>
struct ADASLithiumIonisation : public OpenADAS

ADAS effective ionisation (ADF11)

Template Parameters:

level – The ionisation level of the ion on the left of the reaction

template<int level>
struct ADASLithiumRecombination : public OpenADAS

ADAS effective recombination coefficients (ADF11)

Template Parameters:

level – The ionisation level of the ion on the right of the reaction

Public Functions

inline ADASLithiumRecombination(std::string, Options &alloptions, Solver*)
Parameters:

alloptions – The top-level options. Only uses the [“units”] subsection.

template<int level, char Hisotope>
struct ADASLithiumCX : public OpenADASChargeExchange
Template Parameters:

  • level – The ionisation level of the ion on the right of the reaction

  • Hisotope – The hydrogen isotope (‘h’, ‘d’ or ‘t’)

Public Functions

inline ADASLithiumCX(std::string, Options &alloptions, Solver*)
Parameters:

alloptions – The top-level options. Only uses the [“units”] subsection.

Neon#

These rates are taken from ADAS (96): SCD and PLT are used for the ionisation rate and radiation energy loss; ACD and PRB for the recombination rate and radiation energy loss; and CCD (89) for the charge exchange coupling to hydrogen. The ionisation potential is also included as a source or sink of energy for the electrons.

Reaction

Description

ne + e -> ne+ + 2e

Neon ionisation

ne+ + e -> ne+2 + 2e

ne+2 + e -> ne+3 + 2e

ne+3 + e -> ne+4 + 2e

ne+4 + e -> ne+5 + 2e

ne+5 + e -> ne+6 + 2e

ne+6 + e -> ne+7 + 2e

ne+7 + e -> ne+8 + 2e

ne+8 + e -> ne+9 + 2e

ne+9 + e -> ne+10 + 2e

ne+ + e -> ne

Neon recombination

ne+2 + e -> ne+

ne+3 + e -> ne+2

ne+4 + e -> ne+3

ne+5 + e -> ne+4

ne+6 + e -> ne+5

ne+7 + e -> ne+6

ne+8 + e -> ne+7

ne+9 + e -> ne+8

ne+10 + e -> ne+9

ne+ + h -> ne + h+

Charge exchange with hydrogen

ne+2 + h -> ne+ + h+

ne+3 + h -> ne+2 + h+

ne+4 + h -> ne+3 + h+

ne+5 + h -> ne+4 + h+

ne+6 + h -> ne+5 + h+

ne+7 + h -> ne+6 + h+

ne+8 + h -> ne+7 + h+

ne+9 + h -> ne+8 + h+

ne+10 + h -> ne+9 + h+

ne+ + d -> ne + d+

Charge exchange with deuterium

ne+2 + d -> ne+ + d+

ne+3 + d -> ne+2 + d+

ne+4 + d -> ne+3 + d+

ne+5 + d -> ne+4 + d+

ne+6 + d -> ne+5 + d+

ne+7 + d -> ne+6 + d+

ne+8 + d -> ne+7 + d+

ne+9 + d -> ne+8 + d+

ne+10 + d -> ne+9 + d+

ne+ + t -> ne + t+

Charge exchange with tritium

ne+2 + t -> ne+ + t+

ne+3 + t -> ne+2 + t+

ne+4 + t -> ne+3 + t+

ne+5 + t -> ne+4 + t+

ne+6 + t -> ne+5 + t+

ne+7 + t -> ne+6 + t+

ne+8 + t -> ne+7 + t+

ne+9 + t -> ne+8 + t+

ne+10 + t -> ne+9 + t+

The implementation of these rates is in ADASNeonIonisation, ADASNeonRecombination and ADASNeonCX template classes:

template<std::size_t level>
struct ADASNeonIonisation : public OpenADAS

ADAS effective ionisation (ADF11)

Template Parameters:

level – The ionisation level of the ion on the left of the reaction

template<std::size_t level>
struct ADASNeonRecombination : public OpenADAS

ADAS effective recombination coefficients (ADF11)

Template Parameters:

level – The ionisation level of the ion on the right of the reaction

Public Functions

inline ADASNeonRecombination(std::string, Options &alloptions, Solver*)
Parameters:

alloptions – The top-level options. Only uses the [“units”] subsection.

template<std::size_t level, char Hisotope>
struct ADASNeonCX : public OpenADASChargeExchange
Template Parameters:

  • level – The ionisation level of the ion on the right of the reaction

  • Hisotope – The hydrogen isotope (‘h’, ‘d’ or ‘t’)

Public Functions

inline ADASNeonCX(std::string, Options &alloptions, Solver*)
Parameters:

alloptions – The top-level options. Only uses the [“units”] subsection.