bluemira.plasma_physics.collisions
Tokamak plasma collision formulae.
Functions
|
Debye length |
|
Calculate the reduced mass of a two-particle system |
|
|
|
Calculate the de Broglie wavelength |
|
Calculate the perpendicular impact parameter, a.k.a. b90 |
|
Calculate the value of the Coulomb logarithm for an electron hitting a proton. |
|
Formula for electrical conductivity in a plasma as per L. Spitzer. |
Module Contents
- bluemira.plasma_physics.collisions.debye_length(temperature: float, density: float) float
Debye length
- Parameters:
temperature (float) – Temperature [K]
density (float) – Density [m^-3]
- Returns:
Debye length [m]
- Return type:
float
Notes
Debye length is given by the formula:
\[\lambda_D = \sqrt{\frac{\varepsilon_0 k_B T}{n e^2}}\]where:
\(\varepsilon_0\) is the vacuum permittivity,
\(k_B\) is the Boltzmann constant,
\(T\) is the temperature in Kelvin,
\(n\) is the number density of particles,
\(e\) is the elementary charge.
The conversion of the temperature from Kelvin to Joules implicitly includes the Boltzmann constant (\(k_B\)).
- bluemira.plasma_physics.collisions.reduced_mass(mass_1: float, mass_2: float) float
Calculate the reduced mass of a two-particle system
- Parameters:
mass_1 (float) – Mass of the first particle
mass_2 (float) – Mass of the second particle
- Returns:
Reduced mass
- Return type:
float
Notes
Reduced mass of a two-particle system :
\[\mu_{AB} = \frac{m_A m_B}{m_A + m_B}\]where \(m_A\) and \(m_B\) are the masses of the particles.
- bluemira.plasma_physics.collisions.thermal_velocity(temperature: float, mass: float) float
- Parameters:
temperature (float) – Temperature [K]
mass (float) – Mass of the particle [kg]
- Returns:
Thermal velocity [m/s]
- Return type:
float
Notes
The thremal velocity is calculated as
\[\sqrt{\frac{2 k_B T}{m}}\]The conversion of the temperature from Kelvin to Joules implicitly includes the Boltzmann constant ( \(k_B\) ).
The \(\sqrt 2\) term is for a 3-dimensional system and the most probable velocity in the particle velocity distribution.
- bluemira.plasma_physics.collisions.de_broglie_length(velocity: float, mu_12: float) float
Calculate the de Broglie wavelength
- Parameters:
velocity (float) – Velocity [m/s]
mu_12 (float) – Reduced mass [kg]
- Returns:
De Broglie wavelength [m]
- Return type:
float
Notes
The de Broglie length is given by
\[\lambda_{th} = \frac{h}{2 \cdot \mu_{12} \cdot velocity}\]where \(h\) is the Planck Constant.
- bluemira.plasma_physics.collisions.impact_parameter_perp(velocity: float, mu_12: float) float
Calculate the perpendicular impact parameter, a.k.a. b90
- Parameters:
velocity (float) – Velocity [m/s]
mu_12 (float) – Reduced mass [kg]
- Returns:
Perpendicular impact parameter [m]
- Return type:
float
Notes
\[b_{90} = \frac{e^2}{4 \pi \epsilon_0 \mu_{12} v^2}\]where:
\(e\) is the elementary charge (absolute charge of an electron)
\(\epsilon_0\) is the vacuum permittivity,
\(\mu_{12}\) is the reduced mass,
\(v\) is the relative velocity.
- bluemira.plasma_physics.collisions.coulomb_logarithm(temperature: float, density: float) float
Calculate the value of the Coulomb logarithm for an electron hitting a proton.
- Parameters:
temperature (float) – Temperature [K]
density (float) – Density [1/m^3]
- Returns:
Coulomb logarithm value
- Return type:
float
Notes
The Coulomb logarithm is calculated using the formula:
\[\ln{\Lambda} = \ln{\left(1 + \left(\frac{\lambda_{Debye}} {b_{min}}\right)^2\right)^{1/2}}\]where:
\(\lambda_{Debye}\) is the Debye length,
\(b_{min}\) is the minimum impact parameter, which is defined as the maximum value between \(b_{90}\) (the perpendicular impact parameter) and the de Broglie wavelegth \(\lambda_{th}\)
- bluemira.plasma_physics.collisions.spitzer_conductivity(Z_eff: float, T_e: float, ln_lambda: float) float
Formula for electrical conductivity in a plasma as per L. Spitzer.
- Parameters:
Z_eff (float) – Effective charge [dimensionless]
T_e (float) – Electron temperature on axis [keV] The equation takes in temperature as [eV], so an in-line conversion is used here.
ln_lambda (float) – Coulomb logarithm value
- Returns:
Plasma resistivity [1/Ohm/m]
- Return type:
float
Notes
Spitzer and Haerm, 1953
\(\sigma = 1.92e4 (2-Z_{eff}^{-1/3}) \dfrac{T_{e}^{3/2}}{Z_{eff}ln\Lambda}\)