ASCOT5
Loading...
Searching...
No Matches
physlib.h File Reference

Methods to evaluate elementary physical quantities. More...

#include <math.h>
#include "consts.h"

Go to the source code of this file.

Macros

#define physlib_gamma_vnorm(v)
 Evaluate Lorentz factor from velocity norm.
 
#define physlib_vnorm_gamma(gamma)
 Evaluate velocity norm from Lorentz factor.
 
#define physlib_gamma_pnorm(m, p)
 Evaluate Lorentz factor from momentum norm.
 
#define physlib_gamma_vpar(m, mu, vpar, B)
 Evaluate Lorentz factor from parallel velocity.
 
#define physlib_gamma_ppar(m, mu, ppar, B)
 Evaluate Lorentz factor from parallel momentum.
 
#define physlib_Ekin_gamma(m, gamma)
 Evaluate kinetic energy [J] from Lorentz factor.
 
#define physlib_gamma_Ekin(m, ekin)
 Evaluate Lorentz factor from kinetic energy [J].
 
#define physlib_Ekin_pnorm(m, p)
 Evaluate kinetic energy [J] from momentum norm.
 
#define physlib_Ekin_ppar(m, mu, ppar, B)
 Evaluate kinetic energy [J] from parallel momentum.
 
#define physlib_vnorm_pnorm(m, p)
 Evaluate velocity norm [m/s] from momentum norm.
 
#define physlib_pnorm_vnorm(m, v)
 Evaluate momentum norm [kg m/s] from velocity norm.
 
#define physlib_gc_ppar(p, xi)
 Evaluate guiding center parallel momentum [kg m/s] from momentum norm and pitch.
 
#define physlib_gc_mu(m, p, xi, B)
 Evaluate guiding center magnetic moment [J/T] from momentum norm and pitch.
 
#define physlib_gc_p(m, mu, ppar, B)
 Evaluate guiding center momentum norm [kg m/s] from parallel momentum and magnetic moment.
 
#define physlib_gc_xi(m, mu, ppar, B)
 Evaluate guiding center pitch from parallel momentum and magnetic moment.
 
#define physlib_gyrolength_p(q, p, B)
 Evaluate gyroradius [m] from momentum vector.
 
#define phys_gyrolength_ppar(m, q, mu, ppar, B)
 Evaluate gyroradius [m] from parallel momentum and magnetic moment.
 
#define phys_gyrofreq_pnorm(m, q, p, B)
 Evaluate gyrofrequency [rad/s] from momentum norm.
 
#define phys_gyrofreq_ppar(m, q, mu, ppar, B)
 Evaluate gyrofrequency [rad/s] from parallel momentum and magnetic moment.
 
#define phys_ptoroid_gc(q, R, ppar, psi, B, Bphi)
 Evaluate toroidal canonical momentum for particle.
 
#define phys_vperp_gc(v, vpar)
 Evaluate perpendicular speed from velocity norm and parallel component of velocity for guiding centre.
 
#define phys_pperp_gc(p, ppar)
 Evaluate perpendicular momentum from momentum norm and parallel component of momentum for guiding centre.
 
#define phys_ppar_Ekin(m, ekin, mu, B)
 Evaluate magnitude of parallel momentum from kinetic energy.
 
#define phys_ppar_mu(m, mu, B, q, gyrof)
 Evaluates magnitude of perpendicular velocity based on the magnetic moment.
 

Detailed Description

Methods to evaluate elementary physical quantities.

Definition in file physlib.h.

Macro Definition Documentation

◆ physlib_gamma_vnorm

#define physlib_gamma_vnorm ( v)
Value:
( \
sqrt( 1.0 / ( (1.0 - v / CONST_C) * (1.0 + v / CONST_C) ) ) )
#define CONST_C
Speed of light [m/s].
Definition consts.h:23

Evaluate Lorentz factor from velocity norm.

$ \gamma = \sqrt{\frac{1}{1-v^2/c^2}}$

where

  • $v$ is velocity norm [m/s]

Definition at line 21 of file physlib.h.

◆ physlib_vnorm_gamma

#define physlib_vnorm_gamma ( gamma)
Value:
( \
sqrt( 1.0 - 1.0 / ( gamma * gamma ) ) * CONST_C )

Evaluate velocity norm from Lorentz factor.

$ v = \sqrt{1 - \frac{1}{\gamma^2}}c$

where

  • $v$ is velocity norm [m/s]

Definition at line 33 of file physlib.h.

◆ physlib_gamma_pnorm

#define physlib_gamma_pnorm ( m,
p )
Value:
( \
sqrt(1.0 + ( p * p ) / ( m * m * CONST_C2 ) ) )
#define CONST_C2
Speed of light squared [m^2/s^2].
Definition consts.h:26

Evaluate Lorentz factor from momentum norm.

$\gamma = \sqrt{1 + \left(\frac{p}{mc}\right)^2}$

where

  • $m$ is mass [kg]
  • $p$ is momentum norm [kg m/s]

Definition at line 46 of file physlib.h.

◆ physlib_gamma_vpar

#define physlib_gamma_vpar ( m,
mu,
vpar,
B )
Value:
( \
sqrt( ( 1.0 + (2.0 * mu * B) / ( m * CONST_C2 ) ) / \
( (1.0 - vpar / CONST_C) * (1.0 + vpar / CONST_C) ) ) )

Evaluate Lorentz factor from parallel velocity.

$\gamma = \sqrt{\frac{1 + (2\mu B/mc^2)}{1 - v_\parallel^2/c^2}}$

where

  • $m$ is mass [kg]
  • $\mu$ is magnetic moment [J/T]
  • $v_\parallel$ is parallel velocity [m/s]
  • $B$ is magnetic field norm [T]

Definition at line 61 of file physlib.h.

◆ physlib_gamma_ppar

#define physlib_gamma_ppar ( m,
mu,
ppar,
B )
Value:
( \
sqrt( 1.0 + 2 * mu * B / ( m * CONST_C2 ) + \
ppar * ppar / ( m * m * CONST_C2 ) ) )

Evaluate Lorentz factor from parallel momentum.

$\gamma = \sqrt{1 + 2\mu B/mc^2 + (p_\parallel/mc)^2}$

where

  • $m$ is mass [kg]
  • $\mu$ is magnetic moment [J/T]
  • $p_\parallel$ is parallel momentum [kg m/s]
  • $B$ is magnetic field norm [T]

Definition at line 77 of file physlib.h.

◆ physlib_Ekin_gamma

#define physlib_Ekin_gamma ( m,
gamma )
Value:
( ( gamma - 1.0 ) * m * CONST_C2 )

Evaluate kinetic energy [J] from Lorentz factor.

$E_\mathrm{kin}=(\gamma - 1) * m c^2$

where

  • $m$ is mass [kg]
  • $\gamma$ is the Lorentz factor

Definition at line 91 of file physlib.h.

◆ physlib_gamma_Ekin

#define physlib_gamma_Ekin ( m,
ekin )
Value:
( ekin / ( m * CONST_C2 ) + 1.0 )

Evaluate Lorentz factor from kinetic energy [J].

$\gamma = \frac{E_\mathrm{kin}}{m c^2} + 1$

where

  • $m$ is mass [kg]
  • $\gamma$ is the Lorentz factor

Definition at line 103 of file physlib.h.

◆ physlib_Ekin_pnorm

#define physlib_Ekin_pnorm ( m,
p )
Value:
( \
( physlib_gamma_pnorm(m, p) - 1.0 ) * m * CONST_C2 )
#define physlib_gamma_pnorm(m, p)
Evaluate Lorentz factor from momentum norm.
Definition physlib.h:46

Evaluate kinetic energy [J] from momentum norm.

$E_\mathrm{kin}=(\gamma(p) - 1) * m c^2$

where

  • $m$ is mass [kg]
  • $p$ is momentum norm [kg m/s]

Definition at line 115 of file physlib.h.

◆ physlib_Ekin_ppar

#define physlib_Ekin_ppar ( m,
mu,
ppar,
B )
Value:
( \
( physlib_gamma_ppar(m, mu, ppar, B) - 1.0 ) * m * CONST_C2 )
#define physlib_gamma_ppar(m, mu, ppar, B)
Evaluate Lorentz factor from parallel momentum.
Definition physlib.h:77

Evaluate kinetic energy [J] from parallel momentum.

$E_\mathrm{kin}=(\gamma(m, \mu, p_\parallel, B) - 1) * m c^2$

where

  • $m$ is mass [kg]
  • $p_\parallel$ is parallel momentum [kg m/s]

Definition at line 128 of file physlib.h.

◆ physlib_vnorm_pnorm

#define physlib_vnorm_pnorm ( m,
p )
Value:
( \
p / sqrt(m * m + ( p * p ) / CONST_C2 ) )

Evaluate velocity norm [m/s] from momentum norm.

$v = p/\gamma(p)m$

where

  • $m$ is mass [kg]
  • $p$ is momentum norm [kg m/s]

Definition at line 141 of file physlib.h.

◆ physlib_pnorm_vnorm

#define physlib_pnorm_vnorm ( m,
v )
Value:
( m * v * physlib_gamma_vnorm(v) )
#define physlib_gamma_vnorm(v)
Evaluate Lorentz factor from velocity norm.
Definition physlib.h:21

Evaluate momentum norm [kg m/s] from velocity norm.

$p = \gamma(v)mv$

where

  • $m$ is mass [kg]
  • $v$ is velocity norm [m/s]

Definition at line 154 of file physlib.h.

◆ physlib_gc_ppar

#define physlib_gc_ppar ( p,
xi )
Value:
( p * xi )

Evaluate guiding center parallel momentum [kg m/s] from momentum norm and pitch.

$p_\parallel = \xi p$

where

  • $p$ is momentum norm [kg m/s]
  • $\xi$ is pitch

Definition at line 167 of file physlib.h.

◆ physlib_gc_mu

#define physlib_gc_mu ( m,
p,
xi,
B )
Value:
( \
p * p * ( 1.0 - xi * xi ) / ( 2 * B * m ) )

Evaluate guiding center magnetic moment [J/T] from momentum norm and pitch.

$\mu = (1-\xi^2)p/(2mB)$

where

  • $m$ is mass [kg]
  • $p$ is momentum norm [kg m/s]
  • $\xi$ is pitch
  • $B$ is magnetic field norm [T]

Definition at line 182 of file physlib.h.

◆ physlib_gc_p

#define physlib_gc_p ( m,
mu,
ppar,
B )
Value:
( \
m * CONST_C * sqrt( pow( physlib_gamma_ppar(m, mu, ppar, B), 2) - 1 ) )

Evaluate guiding center momentum norm [kg m/s] from parallel momentum and magnetic moment.

$p = \sqrt{\gamma(m,\mu,p_\parallel,B)^2 - 1} m c$

where

  • $m$ is mass [kg]
  • $\mu$ is magnetic moment [J/T]
  • $p_\parallel$ is parallel momentum [kg m/s]
  • $B$ is magnetic field norm [T]

Definition at line 198 of file physlib.h.

◆ physlib_gc_xi

#define physlib_gc_xi ( m,
mu,
ppar,
B )
Value:
( \
ppar / physlib_gc_p(m, mu, ppar, B) )
#define physlib_gc_p(m, mu, ppar, B)
Evaluate guiding center momentum norm [kg m/s] from parallel momentum and magnetic moment.
Definition physlib.h:198

Evaluate guiding center pitch from parallel momentum and magnetic moment.

$\xi = p_\parallel / p(m,\mu,p_\parallel,B) $

where

  • $m$ is mass [kg]
  • $\mu$ is magnetic moment [J/T]
  • $p_\parallel$ is parallel momentum [kg m/s]
  • $B$ is magnetic field [T]

Definition at line 214 of file physlib.h.

◆ physlib_gyrolength_p

#define physlib_gyrolength_p ( q,
p,
B )
Value:
( \
math_dot(p, B) / ( fabs(q) * math_dot(B, B) ) )
#define math_dot(a, b)
Calculate dot product a[3] dot b[3].
Definition math.h:31

Evaluate gyroradius [m] from momentum vector.

$\rho_g = \frac{\mathbf{p}\cdot\mathbf{B}}{|q|B^2}$

where

  • $q$ is charge [C]
  • $\mathbf{v}$ is momentum vector [kg m/s]
  • $\mathbf{B}$ is magnetic field vector [T]

Definition at line 228 of file physlib.h.

◆ phys_gyrolength_ppar

#define phys_gyrolength_ppar ( m,
q,
mu,
ppar,
B )
Value:
( \
sqrt( 2 * m * mu * \
physlib_gamma_ppar(m, mu, ppar, B) / B ) / fabs(q) )

Evaluate gyroradius [m] from parallel momentum and magnetic moment.

$\rho_g = \frac{1}{|q|}\sqrt{\frac{2\gamma(m, \mu, p_\parallel, B) m \mu}
            {B}}$

where

  • $m$ is mass [kg]
  • $q$ is charge [C]
  • $\mu$ is magnetic moment [J/T]
  • $p_\parallel$ is parallel momentum [kg m/s]
  • $B$ is magnetic field norm [T]

Definition at line 245 of file physlib.h.

◆ phys_gyrofreq_pnorm

#define phys_gyrofreq_pnorm ( m,
q,
p,
B )
Value:
( \
fabs(q) * B / ( m * physlib_gamma_pnorm(m, p) ) )

Evaluate gyrofrequency [rad/s] from momentum norm.

$\omega_g = \frac{q B}{\gamma(p) m}$

where

  • $m$ is mass [kg]
  • $q$ is charge [C]
  • $p$ is momentum norm [kg m/s]
  • $B$ is magnetic field norm [T]

Definition at line 261 of file physlib.h.

◆ phys_gyrofreq_ppar

#define phys_gyrofreq_ppar ( m,
q,
mu,
ppar,
B )
Value:
( \
fabs(q) * B / ( m * physlib_gamma_ppar(m, mu, ppar, B) ) )

Evaluate gyrofrequency [rad/s] from parallel momentum and magnetic moment.

$\omega_g = \frac{q B}{\gamma(m, \mu, p_\parallel, B) m}$

where

  • $m$ is mass [kg]
  • $q$ is charge [C]
  • $\mu$ is magnetic moment [J/T]
  • $p_\parallel$ is parallel momentum [kg m/s]
  • $B$ is magnetic field norm [T]

Definition at line 278 of file physlib.h.

◆ phys_ptoroid_gc

#define phys_ptoroid_gc ( q,
R,
ppar,
psi,
B,
Bphi )
Value:
( \
ppar * R * (Bphi / B) + q * psi )
@ R
Definition hist.h:18

Evaluate toroidal canonical momentum for particle.

$\Pi_{\phi,fo} = Rp_{\phi} + q\psi$

where

$\Pi_{\phi,fo}$ is the canonical momentum conjugate to (the toroidal angle) $\phi$ [J s/rad] $R$ is the major radius [m] $p_{\phi}$ is the toroidal momentum [kg m/s] \fq $ is the charge of the particle [C]
$\psi\f$ is the poloidal magnetic flux [Tm^2] */ #define phys_ptoroid_fo(q, R, pphi, psi) ( \ R * pphi + q * psi )

/**

Evaluate toroidal canonical momentum for guiding center

$\Pi_{\phi,gc} = p_{\parallel}R\frac{B_{\phi}}{B} + q\psi$

where

$\Pi_{\phi,gc}$ is the canonical momentum conjugate to (the toroidal angle) $\phi$ [J s/rad] $p_{\parallel}$ is the parallel momentum [kg m/s] $R$ is the major radius [m] $B_{\phi}$ is the toroidal component of the magnetic flux density [T] $B$ is the magnitude of the magnetic flux density [T] \fq $ is the charge of the particle [C]
$\psi\f$ is the poloidal magnetic flux [Tm^2]

Definition at line 314 of file physlib.h.

◆ phys_vperp_gc

#define phys_vperp_gc ( v,
vpar )
Value:
sqrt(pow(v,2) - pow(vpar,2))

Evaluate perpendicular speed from velocity norm and parallel component of velocity for guiding centre.

$v_{\perp} = \sqrt{v^2 - v_{\parallel}^2}$

where

  • $v_{\perp}$ is perpendicular speed [m/s]
  • $v$ is speed [m/s]
  • $v_{\parallel}$ is parallel velocity component [m/s]

Definition at line 329 of file physlib.h.

◆ phys_pperp_gc

#define phys_pperp_gc ( p,
ppar )
Value:
sqrt(pow(p,2) - pow(ppar,2))

Evaluate perpendicular momentum from momentum norm and parallel component of momentum for guiding centre.

$p_{\perp} = \sqrt{p^2 - p_{\parallel}^2}$

where

  • $p_{\perp}$ is perpendicular momentum [kg m/s]
  • $p$ is momentum [kg m/s]
  • $p_{\parallel}$ is parallel momentum component [kg m/s]

Definition at line 343 of file physlib.h.

◆ phys_ppar_Ekin

#define phys_ppar_Ekin ( m,
ekin,
mu,
B )
Value:
(sqrt((physlib_gamma_Ekin(m, ekin)*\
physlib_gamma_Ekin(m, ekin) - 1.0)*m*m*CONST_C2 - 2.0*m*mu*B))
#define physlib_gamma_Ekin(m, ekin)
Evaluate Lorentz factor from kinetic energy [J].
Definition physlib.h:103

Evaluate magnitude of parallel momentum from kinetic energy.

$abs(p_{\parallel}) = \sqrt{(\gamma^2 - 1)m^2c^2 - 2m\mu B}$

where

$abs(p_{\parallel})$ is the magnitude of the parallel component of momentum $\gammaf$ is the Lorentz factor [-]
$@_fakenlmf$ is the mass of the particle [kg] $c$ is the speed of light in vacuum [m/s] $\mu$ is the magnetic moment of the particle [J/T] $B$ is the magnitude of magnetic flux density [T]

Definition at line 360 of file physlib.h.

◆ phys_ppar_mu

#define phys_ppar_mu ( m,
mu,
B,
q,
gyrof )
Value:
(m*CONST_C*sqrt((q*B/(gyrof*m)) * \
(q*B/(gyrof*m)) - 2*mu*B/(m*CONST_C2) - 1))

Evaluates magnitude of perpendicular velocity based on the magnetic moment.

$v_{\perp} = \sqrt{\frac{2\mu B}{m}}$

where

$v_{\perp}$ is the perpendicular speed [m/S] $\mu$ is the magnetic moment [J/T] \fB $ is the magnitude of the magnetic flux density [T]
$\m\f$ is the mass of the particle [kg] */ #define phys_vperp_mu(m, mu, B) sqrt(2*mu*B/m)

/**

Evaluates magnitude of parallel momentum from magnetic moment and gyrofrequency.

$abs(p_{\parallel}) = mc\sqrt{\left(\frac{qB}{\Omega m}\right)^2 - \frac{2
                         \mu B}{mc^2} - 1}$

where

$abs(p_{\parallel})$ is the magnitude of the parallel component of momentum \fm $ is the mass of the particle [kg]
$@_fakenlc $ is the speed of light in vacuum [m/s]
$@_fakenlq $ is the charge of the particle [C]
$@_fakenlB $ is the magnitude of the magnetic flux density [T]
$\Omega\f$ is the gyrofrequency [rad/s] $\mu$ is the magnetic moment [J/T]

Definition at line 396 of file physlib.h.