"""Data for different fusion reactions."""
from __future__ import annotations
from typing import ClassVar
import numpy as np
from cfspopcon.unit_handling import Quantity, Unitfull, convert_units, ureg, wraps_ufunc
[docs]
class FusionReaction:
"""A base class for different fusion reactions."""
instances: ClassVar[dict[str, FusionReaction]] = dict()
def __init__(self) -> None:
"""Records children classes in instances."""
self.instances[self.__class__.__name__] = self
[docs]
class DTFusionBoschHale(FusionReaction):
"""Deuterium-Tritium reaction using Bosch-Hale cross-section."""
def __init__(self) -> None:
"""Sets the reaction energies for the Deuterium-Tritium reaction."""
super().__init__()
self.energy_per_reaction = Quantity(17.6, ureg.MeV)
self.energy_to_neutrals_per_reaction = self.energy_per_reaction * 4.0 / 5.0
self.energy_to_charged_per_reaction = self.energy_per_reaction * 1.0 / 5.0
[docs]
@staticmethod
@wraps_ufunc(input_units=dict(ion_temp=ureg.keV), return_units=dict(sigmav=ureg.cm**3 / ureg.s))
def calc_rate_coefficient(ion_temp: float) -> float:
r"""Calculate :math:`\\langle \\sigma v \rangle` for a given ion temperature.
Cross-section from :cite:`bosch_improved_1992`, equation 12, using coefficients from table VII, page 625 (first column is DT reaction)
Args:
ion_temp: [keV]
Returns:
:math:`\\langle \\sigma v \rangle` [cm^3/s]
"""
C = [0.0, 1.17302e-9, 1.51361e-2, 7.51886e-2, 4.60643e-3, 1.35000e-2, -1.06750e-4, 1.36600e-5]
B_G = 34.3827
mr_c2 = 1124656
theta = ion_temp / (
1 - (ion_temp * (C[2] + ion_temp * (C[4] + ion_temp * C[6]))) / (1 + ion_temp * (C[3] + ion_temp * (C[5] + ion_temp * C[7])))
)
eta = (B_G**2 / (4 * theta)) ** (1 / 3)
sigmav: float = C[1] * theta * np.sqrt(eta / (mr_c2 * ion_temp**3)) * np.exp(-3 * eta)
return sigmav
[docs]
def calc_average_ion_mass(self, heavier_fuel_species_fraction: Unitfull) -> Unitfull:
"""Calculate the average mass of the fuel ions.
Args:
heavier_fuel_species_fraction: n_heavier / (n_heavier + n_lighter) number fraction.
Returns:
:term:`average_ion_mass` [amu]
"""
average_ion_mass = 2.0 * (1 - heavier_fuel_species_fraction) + 3.0 * heavier_fuel_species_fraction
return average_ion_mass * ureg.amu
[docs]
def calc_energy_per_reaction(self) -> Unitfull:
"""Returns the total energy per reaction."""
return convert_units(self.energy_per_reaction, ureg.MJ)
[docs]
def calc_energy_to_neutrals_per_reaction(self) -> Unitfull:
"""Returns the energy going to uncharged species per reaction."""
return convert_units(self.energy_to_neutrals_per_reaction, ureg.MJ)
[docs]
def calc_energy_to_charged_per_reaction(self) -> Unitfull:
"""Returns the energy going to charged species per reaction."""
return convert_units(self.energy_to_charged_per_reaction, ureg.MJ)
[docs]
def calc_power_density(self, ion_temp: Unitfull, heavier_fuel_species_fraction: Unitfull) -> Unitfull:
"""Returns the total power density."""
sigmav = self.calc_rate_coefficient(ion_temp)
fuel_ratio = heavier_fuel_species_fraction * (1.0 - heavier_fuel_species_fraction)
power_density = sigmav * self.energy_per_reaction * fuel_ratio
return convert_units(power_density, ureg.MW * ureg.m**3)
[docs]
def calc_power_density_to_neutrals(self, ion_temp: Unitfull, heavier_fuel_species_fraction: Unitfull) -> Unitfull:
"""Returns the power density going to uncharged species."""
return (4.0 / 5.0) * self.calc_power_density(ion_temp, heavier_fuel_species_fraction)
[docs]
def calc_power_density_to_charged(self, ion_temp: Unitfull, heavier_fuel_species_fraction: Unitfull) -> Unitfull:
"""Returns the power density going to charged species."""
return (1.0 / 5.0) * self.calc_power_density(ion_temp, heavier_fuel_species_fraction)
[docs]
class DTFusionHively(DTFusionBoschHale):
"""Deuterium-Tritium reaction using Hively cross-section."""
[docs]
@staticmethod
@wraps_ufunc(input_units=dict(ion_temp=ureg.keV), return_units=dict(sigmav=ureg.cm**3 / ureg.s))
def calc_rate_coefficient(ion_temp: float) -> float:
r"""Calculate :math:`\\langle \\sigma v \rangle` for a given ion temperature.
Cross-section from :cite:`hively_convenient_1977`
Args:
ion_temp: [keV]
Returns:
:math:`\\langle \\sigma v \rangle` [cm^3/s]
"""
A = [-21.377692, -25.204054, -7.1013427 * 1e-2, 1.9375451 * 1e-4, 4.9246592 * 1e-6, -3.9836572 * 1e-8]
r = 0.2935
sigmav: float = np.exp(
A[0] / ion_temp**r + A[1] + A[2] * ion_temp + A[3] * ion_temp**2.0 + A[4] * ion_temp**3.0 + A[5] * ion_temp**4.0
)
return sigmav
[docs]
class DDFusionBoschHale(FusionReaction):
"""Deuterium-Deuterium reaction using Bosch-Hale cross-section."""
def __init__(self) -> None:
"""Sets the reaction energies for the Deuterium-Deuterium reaction."""
super().__init__()
self.energy_per_DD_to_pT_reaction = Quantity(1.01 + 3.02, ureg.MeV)
self.energy_per_DD_to_nHe3_reaction = Quantity(0.82 + 2.45, ureg.MeV)
self.energy_to_neutrals_per_DD_to_nHe3_reaction = Quantity(2.45, ureg.MeV)
self.energy_to_charged_per_DD_to_nHe3_reaction = Quantity(0.82, ureg.MeV)
self.energy_to_neutrals_per_DD_to_pT_reaction = Quantity(0.0, ureg.MeV)
self.energy_to_charged_per_DD_to_pT_reaction = Quantity(1.01 + 3.02, ureg.MeV)
[docs]
@staticmethod
@wraps_ufunc(
input_units=dict(ion_temp=ureg.keV),
return_units=dict(sigmav_combined=ureg.cm**3 / ureg.s, sigmav_DD_to_pT=ureg.cm**3 / ureg.s, sigmav_DD_to_nHe3=ureg.cm**3 / ureg.s),
output_core_dims=[(), (), ()],
)
def calc_rate_coefficient(ion_temp: float) -> tuple[float, float, float]:
r"""Calculate :math:`\\langle \\sigma v \rangle` for a given ion temperature.
Cross-section from :cite:`bosch_improved_1992`
Args:
ion_temp: [keV]
Returns:
:math:`\\langle \\sigma v \rangle` [cm^3/s]
"""
# For D(d,n)3He
cBH_1 = [((31.3970**2) / 4.0) ** (1.0 / 3.0), 5.65718e-12, 3.41e-03, 1.99e-03, 0, 1.05e-05, 0, 0] # 3.72e-16,
mc2_1 = 937814.0
# For D(d,p)T
cBH_2 = [((31.3970**2) / 4.0) ** (1.0 / 3.0), 5.43360e-12, 5.86e-03, 7.68e-03, 0, -2.96e-06, 0, 0] # 3.57e-16,
mc2_2 = 937814.0
thetaBH_1 = ion_temp / (
1
- (
(cBH_1[2] * ion_temp + cBH_1[4] * ion_temp**2 + cBH_1[6] * ion_temp**3)
/ (1 + cBH_1[3] * ion_temp + cBH_1[5] * ion_temp**2 + cBH_1[7] * ion_temp**3)
)
)
thetaBH_2 = ion_temp / (
1
- (
(cBH_2[2] * ion_temp + cBH_2[4] * ion_temp**2 + cBH_2[6] * ion_temp**3)
/ (1 + cBH_2[3] * ion_temp + cBH_2[5] * ion_temp**2 + cBH_2[7] * ion_temp**3)
)
)
etaBH_1 = cBH_1[0] / (thetaBH_1 ** (1.0 / 3.0))
etaBH_2 = cBH_2[0] / (thetaBH_2 ** (1.0 / 3.0))
sigmav_DD_to_nHe3 = cBH_1[1] * thetaBH_1 * np.sqrt(etaBH_1 / (mc2_1 * (ion_temp**3.0))) * np.exp(-3.0 * etaBH_1)
sigmav_DD_to_pT = cBH_2[1] * thetaBH_2 * np.sqrt(etaBH_2 / (mc2_2 * (ion_temp**3.0))) * np.exp(-3.0 * etaBH_2)
sigmav_combined = sigmav_DD_to_pT + sigmav_DD_to_nHe3
if not np.isreal(sigmav_combined):
return np.nan, np.nan, np.nan
else:
return sigmav_combined, sigmav_DD_to_pT, sigmav_DD_to_nHe3
[docs]
def calc_average_ion_mass(self) -> Unitfull:
"""Returns the average mass of the fuel ions.
Returns:
:term:`average_ion_mass` [amu]
"""
return 2.0 * ureg.amu
[docs]
def calc_energy_per_reaction(self, ion_temp: Unitfull) -> Unitfull:
"""Returns the total energy per reaction."""
sigmav_combined, sigmav_DD_to_pT, sigmav_DD_to_nHe3 = self.calc_rate_coefficient(ion_temp)
energy_per_reaction = (
self.energy_per_DD_to_pT_reaction * sigmav_DD_to_pT + self.energy_per_DD_to_nHe3_reaction * sigmav_DD_to_nHe3
) / sigmav_combined
return convert_units(energy_per_reaction, ureg.MJ)
[docs]
def calc_energy_to_neutrals_per_reaction(self, ion_temp: Unitfull) -> Unitfull:
"""Returns the energy going to uncharged species per reaction."""
sigmav_combined, sigmav_DD_to_pT, sigmav_DD_to_nHe3 = self.calc_rate_coefficient(ion_temp)
energy_to_neutrals_per_reaction = (
self.energy_to_neutrals_per_DD_to_pT_reaction * sigmav_DD_to_pT
+ self.energy_to_neutrals_per_DD_to_nHe3_reaction * sigmav_DD_to_nHe3
) / sigmav_combined
return convert_units(energy_to_neutrals_per_reaction, ureg.MJ)
[docs]
def calc_energy_to_charged_per_reaction(self, ion_temp: Unitfull) -> Unitfull:
"""Returns the energy going to charged species per reaction."""
sigmav_combined, sigmav_DD_to_pT, sigmav_DD_to_nHe3 = self.calc_rate_coefficient(ion_temp)
energy_to_charged_per_reaction = (
self.energy_to_charged_per_DD_to_pT_reaction * sigmav_DD_to_pT
+ self.energy_to_charged_per_DD_to_nHe3_reaction * sigmav_DD_to_nHe3
) / sigmav_combined
return convert_units(energy_to_charged_per_reaction, ureg.MJ)
[docs]
def calc_power_density(self, ion_temp: Unitfull) -> Unitfull:
"""Returns the total power density."""
_, sigmav_DD_to_pT, sigmav_DD_to_nHe3 = self.calc_rate_coefficient(ion_temp)
power_density = sigmav_DD_to_pT * self.energy_per_DD_to_pT_reaction + sigmav_DD_to_nHe3 * self.energy_per_DD_to_nHe3_reaction
return convert_units(power_density, ureg.MW * ureg.m**3)
[docs]
def calc_power_density_to_neutrals(self, ion_temp: Unitfull) -> Unitfull:
"""Returns the power density going to uncharged species."""
_, sigmav_DD_to_pT, sigmav_DD_to_nHe3 = self.calc_rate_coefficient(ion_temp)
power_density = (
sigmav_DD_to_pT * self.energy_to_neutrals_per_DD_to_pT_reaction
+ sigmav_DD_to_nHe3 * self.energy_to_neutrals_per_DD_to_nHe3_reaction
)
return convert_units(power_density, ureg.MW * ureg.m**3)
[docs]
def calc_power_density_to_charged(self, ion_temp: Unitfull) -> Unitfull:
"""Returns the power density going to charged species."""
_, sigmav_DD_to_pT, sigmav_DD_to_nHe3 = self.calc_rate_coefficient(ion_temp)
power_density = (
sigmav_DD_to_pT * self.energy_to_charged_per_DD_to_pT_reaction
+ sigmav_DD_to_nHe3 * self.energy_to_charged_per_DD_to_nHe3_reaction
)
return convert_units(power_density, ureg.MW * ureg.m**3)
[docs]
class DDFusionHively(DDFusionBoschHale):
"""Deuterium-Deuterium reaction using Hively cross-section."""
[docs]
@staticmethod
@wraps_ufunc(
input_units=dict(ion_temp=ureg.keV),
return_units=dict(sigmav_combined=ureg.cm**3 / ureg.s, sigmav_DD_to_pT=ureg.cm**3 / ureg.s, sigmav_DD_to_nHe3=ureg.cm**3 / ureg.s),
output_core_dims=[(), (), ()],
)
def calc_rate_coefficient(ion_temp: float) -> tuple[float, float, float]:
r"""Calculate :math:`\\langle \\sigma v \rangle` for a given ion temperature.
Cross-section from :cite:`hively_convenient_1977`
Args:
ion_temp: [keV]
Returns:
:math:`\\langle \\sigma v \rangle` [cm^3/s]
"""
a_1 = [
-15.511891,
-35.318711,
-1.2904737 * 1e-2,
2.6797766 * 1e-4,
-2.9198685 * 1e-6,
1.2748415 * 1e-8,
] # For D(d,p)T
r_1 = 0.3735
a_2 = [
-15.993842,
-35.017640,
-1.3689787 * 1e-2,
2.7089621 * 1e-4,
-2.9441547 * 1e-6,
1.2841202 * 1e-8,
] # For D(d,n)3He
r_2 = 0.3725
# Ti in units of keV, sigmav in units of cm^3/s
sigmav_DD_to_pT = np.exp(
a_1[0] / ion_temp**r_1 + a_1[1] + a_1[2] * ion_temp + a_1[3] * ion_temp**2.0 + a_1[4] * ion_temp**3.0 + a_1[5] * ion_temp**4.0
)
sigmav_DD_to_nHe3 = np.exp(
a_2[0] / ion_temp**r_2 + a_2[1] + a_2[2] * ion_temp + a_2[3] * ion_temp**2.0 + a_2[4] * ion_temp**3.0 + a_2[5] * ion_temp**4.0
)
sigmav_combined = sigmav_DD_to_pT + sigmav_DD_to_nHe3
return sigmav_combined, sigmav_DD_to_pT, sigmav_DD_to_nHe3
[docs]
class DHe3Fusion(FusionReaction):
"""Deuterium-Helium-3 reaction (Bosch-Hale cross-section)."""
def __init__(self) -> None:
"""Sets the reaction energies for the Deuterium-Tritium reaction."""
super().__init__()
self.energy_per_reaction = Quantity(18.3, ureg.MeV)
self.energy_to_neutrals_per_reaction = Quantity(0.0, ureg.MeV)
self.energy_to_charged_per_reaction = self.energy_per_reaction
[docs]
@staticmethod
@wraps_ufunc(input_units=dict(ion_temp=ureg.keV), return_units=dict(sigmav=ureg.cm**3 / ureg.s))
def calc_rate_coefficient(ion_temp: float) -> float:
r"""Deuterium-Helium-3 reaction.
Calculate :math:`\langle \sigma v \rangle` for a given characteristic ion energy.
Function tested on available data at [1, 2, 5, 10, 20, 50, 100] keV.
Maximum error = 8.4% within range 2-100 keV and should not be used outside range [2, 100] keV.
Uses DD cross section formulation :cite:`bosch_improved_1992`.
Args:
ion_temp: [keV]
Returns:
:math:`\langle \sigma v \rangle` in cm^3/s.
"""
# For He3(d,p)4He
cBH_1 = [
((68.7508**2) / 4.0) ** (1.0 / 3.0),
5.51036e-10, # 3.72e-16,
6.41918e-03,
-2.02896e-03,
-1.91080e-05,
1.35776e-04,
0,
0,
]
mc2_1 = 1124572.0
thetaBH_1 = ion_temp / (
1
- (
(cBH_1[2] * ion_temp + cBH_1[4] * ion_temp**2 + cBH_1[6] * ion_temp**3)
/ (1 + cBH_1[3] * ion_temp + cBH_1[5] * ion_temp**2 + cBH_1[7] * ion_temp**3.0)
)
)
etaBH_1 = cBH_1[0] / (thetaBH_1 ** (1.0 / 3.0))
sigmav = cBH_1[1] * thetaBH_1 * np.sqrt(etaBH_1 / (mc2_1 * (ion_temp**3.0))) * np.exp(-3.0 * etaBH_1)
return float(sigmav)
[docs]
def calc_average_ion_mass(self, heavier_fuel_species_fraction: Unitfull) -> Unitfull:
"""Calculate the average mass of the fuel ions.
Args:
heavier_fuel_species_fraction: n_heavier / (n_heavier + n_lighter) number fraction.
Returns:
:term:`average_ion_mass` [amu]
"""
average_fuel_ion_mass = 2.0 * (1 - heavier_fuel_species_fraction) + 3.0 * heavier_fuel_species_fraction
return average_fuel_ion_mass * ureg.amu
[docs]
def calc_energy_per_reaction(self) -> Unitfull:
"""Returns the total energy per reaction."""
return convert_units(self.energy_per_reaction, ureg.MJ)
[docs]
def calc_energy_to_neutrals_per_reaction(self) -> Unitfull:
"""Returns the energy going to uncharged species per reaction."""
return convert_units(self.energy_to_neutrals_per_reaction, ureg.MJ)
[docs]
def calc_energy_to_charged_per_reaction(self) -> Unitfull:
"""Returns the energy going to charged species per reaction."""
return convert_units(self.energy_to_charged_per_reaction, ureg.MJ)
[docs]
def calc_power_density(self, ion_temp: Unitfull, heavier_fuel_species_fraction: Unitfull) -> Unitfull:
"""Returns the total power density."""
sigmav = self.calc_rate_coefficient(ion_temp)
fuel_ratio = heavier_fuel_species_fraction * (1.0 - heavier_fuel_species_fraction)
power_density = sigmav * self.energy_per_reaction * fuel_ratio
return convert_units(power_density, ureg.MW * ureg.m**3)
[docs]
def calc_power_density_to_neutrals(self) -> Unitfull:
"""Returns the power density going to uncharged species."""
return Quantity(0.0, ureg.MW * ureg.m**3)
[docs]
def calc_power_density_to_charged(self, ion_temp: Unitfull, heavier_fuel_species_fraction: Unitfull) -> Unitfull:
"""Returns the power density going to charged species."""
return self.calc_power_density(ion_temp, heavier_fuel_species_fraction)
[docs]
class pB11Fusion(FusionReaction):
"""Proton-Boron-11 reaction (Nevins-Swain cross-section)."""
def __init__(self) -> None:
"""Sets the reaction energies for the Deuterium-Tritium reaction."""
super().__init__()
self.energy_per_reaction = Quantity(8.7, ureg.MeV)
self.energy_to_neutrals_per_reaction = Quantity(0.0, ureg.MeV)
self.energy_to_charged_per_reaction = self.energy_per_reaction
[docs]
@staticmethod
@wraps_ufunc(input_units=dict(ion_temp=ureg.keV), return_units=dict(sigmav=ureg.cm**3 / ureg.s))
def calc_rate_coefficient(ion_temp: float) -> float:
r"""Proton (hydrogen)-Boron11 reaction.
Calculate :math:`\langle \sigma v \rangle` for a given characteristic ion energy.
Uses cross section from Nevins and Swain :cite:`nevins_thermonuclear_2000`.
Updated cross sections in :cite:`sikora_new_2016`, and :cite:`putvinski_fusion_2019` are not in analytic form.
Args:
ion_temp: [keV]
Returns:
:math:`\langle \sigma v \rangle` in cm^3/s.
"""
# High temperature (T>60 keV)
# For B11(p,alpha)alpha,alpha
cBH_1 = [
((22589.0) / 4.0) ** (1.0 / 3.0),
4.4467e-14,
-5.9357e-02,
2.0165e-01,
1.0404e-03,
2.7621e-03,
-9.1653e-06,
9.8305e-07,
]
mc2_1 = 859526.0
thetaBH_1 = ion_temp / (
1
- (
(cBH_1[2] * ion_temp + cBH_1[4] * ion_temp**2 + cBH_1[6] * ion_temp**3)
/ (1 + cBH_1[3] * ion_temp + cBH_1[5] * ion_temp**2 + cBH_1[7] * ion_temp**3)
)
)
etaBH_1 = cBH_1[0] / (thetaBH_1 ** (1.0 / 3.0))
sigmavNRhigh = (cBH_1[1] * thetaBH_1 * np.sqrt(etaBH_1 / (mc2_1 * (ion_temp**3.0))) * np.exp(-3.0 * etaBH_1)) * 1e6 # m3 to cm3
# Low temperature (T<60 keV)
E0 = ((17.81) ** (1.0 / 3.0)) * (ion_temp ** (2.0 / 3.0))
deltaE0 = 4.0 * np.sqrt(ion_temp * E0 / 3.0)
tau = (3.0 * E0) / ion_temp
Mp = 1.0 # *1.67e-27
MB = 11.0 # *1.67e-27
Mr = (Mp * MB) / (Mp + MB)
C0 = 197.000 * 1e-25 # MeVb to kev/m^2
C1 = 0.240 * 1e-25 # MeVb to kev/m^2
C2 = 0.000231 * 1e-25 # MeVb to kev/m^2
Seff = C0 * (1 + (5.0 / (12.0 * tau))) + C1 * (E0 + (35.0 / 36.0) * ion_temp) + C2 * (E0**2.0 + (89.0 / 36.0) * E0 * ion_temp)
sigmavNRlow = (np.sqrt(2 * ion_temp / Mr) * ((deltaE0 * Seff) / (ion_temp ** (2.0))) * np.exp(-tau)) * 1e6 # m3 to cm3
# 148 keV resonance
sigmavR = ((5.41e-21) * ((1.0 / ion_temp) ** (3.0 / 2.0)) * np.exp(-148.0 / ion_temp)) * 1e6 # m3 to cm3
sigmav = sigmavNRhigh
if ion_temp < 60.0: # keV
sigmav = sigmavNRlow + sigmavR
elif (ion_temp > 60.0) and (ion_temp < 130): # keV
sigmav = sigmavNRhigh + sigmavR
return float(sigmav)
[docs]
def calc_average_ion_mass(self, heavier_fuel_species_fraction: Unitfull) -> Unitfull:
"""Calculate the average mass of the fuel ions.
Args:
heavier_fuel_species_fraction: n_heavier / (n_heavier + n_lighter) number fraction.
Returns:
:term:`average_ion_mass` [amu]
"""
average_fuel_ion_mass = 1.0 * (1 - heavier_fuel_species_fraction) + 11.0 * heavier_fuel_species_fraction
return average_fuel_ion_mass * ureg.amu
[docs]
def calc_energy_per_reaction(self) -> Unitfull:
"""Returns the total energy per reaction."""
return convert_units(self.energy_per_reaction, ureg.MJ)
[docs]
def calc_energy_to_neutrals_per_reaction(self) -> Unitfull:
"""Returns the energy going to uncharged species per reaction."""
return convert_units(self.energy_to_neutrals_per_reaction, ureg.MJ)
[docs]
def calc_energy_to_charged_per_reaction(self) -> Unitfull:
"""Returns the energy going to charged species per reaction."""
return convert_units(self.energy_to_charged_per_reaction, ureg.MJ)
[docs]
def calc_power_density(self, ion_temp: Unitfull, heavier_fuel_species_fraction: Unitfull) -> Unitfull:
"""Returns the total power density."""
sigmav = self.calc_rate_coefficient(ion_temp)
fuel_ratio = heavier_fuel_species_fraction * (1.0 - heavier_fuel_species_fraction)
power_density = sigmav * self.energy_per_reaction * fuel_ratio
return convert_units(power_density, ureg.MW * ureg.m**3)
[docs]
def calc_power_density_to_neutrals(self) -> Unitfull:
"""Returns the power density going to uncharged species."""
return Quantity(0.0, ureg.MW * ureg.m**3)
[docs]
def calc_power_density_to_charged(self, ion_temp: Unitfull, heavier_fuel_species_fraction: Unitfull) -> Unitfull:
"""Returns the power density going to charged species."""
return self.calc_power_density(ion_temp, heavier_fuel_species_fraction)
REACTIONS = dict(
DTFusionBoschHale=DTFusionBoschHale(),
DTFusionHively=DTFusionHively(),
DDFusionBoschHale=DDFusionBoschHale(),
DDFusionHively=DDFusionHively(),
DHe3Fusion=DHe3Fusion(),
pB11Fusion=pB11Fusion(),
)