"""Functions related to gyro-kinetic diamagnetic fluxes."""
import xarray as xr
from torx.arrays import make_xarray
from torx.grid import Grid2D
from torx.equilibrium import EquilibriumBaseClass
from torx.units import Normalization
from torx.decorators import autodoc_function
from torx.operators import grad, rot
from torx.vector import vector_dot, vector_cross
[docs]
@autodoc_function
def diamagnetic_particle_flux(
grid: Grid2D,
equi: EquilibriumBaseClass,
norm: Normalization,
W_par: xr.DataArray,
W_perp: xr.DataArray,
charge: float=1.00
):
"""Calculate the diamagnetic particle flux using 2nd order moments."""
assert W_par.dims == W_perp.dims
r_norm, z_norm = grid.coords_like(W_par)
absB = equi.magfield_absolute(r_norm, z_norm)
normb = equi.magfield_vector(r_norm, z_norm) / absB
er = equi.magfield_vector_radial(r_norm, z_norm, normalize=True)
gradB = grad(absB, grid=grid)
curlb = rot(normb, grid=grid)
er_dot_curlb = vector_dot(er, curlb) / absB
er_dot_gradB = vector_dot(er, vector_cross(normb, gradB)) / absB**2
prefac = 1 / charge
return make_xarray(
prefac * 2 * W_par * er_dot_curlb \
+ prefac * W_perp * er_dot_gradB,
norm=norm.n0 * norm.c_s0 * norm.rho_s0.m / norm.R0.m,
)
[docs]
@autodoc_function
def diamagnetic_particle_flux_maxw(
grid: Grid2D,
equi: EquilibriumBaseClass,
norm: Normalization,
W_tot: xr.DataArray,
charge: float=1.00
):
"""Calculate the diamagnetic particle flux for an equilibrium distribution."""
r_norm, z_norm = grid.coords_like(W_tot)
absB = equi.magfield_absolute(r_norm, z_norm)
normb = equi.magfield_vector(r_norm, z_norm) / absB
er = equi.magfield_vector_radial(r_norm, z_norm, normalize=True)
gradB = grad(absB, grid=grid)
curlb = rot(normb, grid=grid)
er_dot_curlb = vector_dot(er, curlb) / absB
er_dot_gradB = vector_dot(er, vector_cross(normb, gradB)) / absB / absB
prefac = 1 / charge
return make_xarray(
prefac * W_tot * (er_dot_curlb + er_dot_gradB),
norm=norm.n0 * norm.c_s0 * norm.rho_s0 / norm.R0,
)
[docs]
@autodoc_function
def diamagnetic_heat_flux(
grid: Grid2D,
equi: EquilibriumBaseClass,
norm: Normalization,
K_par: xr.DataArray,
K_perp: xr.DataArray
):
"""Calculate the diamagnetic heat flux using 4th order moments."""
assert K_par.dims == K_perp.dims
r_norm, z_norm = grid.coords_like(K_par)
absB = equi.magfield_absolute(r_norm, z_norm)
normb = equi.magfield_vector(r_norm, z_norm) / absB
er = equi.magfield_vector_radial(r_norm, z_norm,
normalize=True)
gradB = grad(absB, grid=grid)
curlb = rot(normb, grid=grid)
er_dot_curlb = vector_dot(er, curlb)
er_dot_gradB = vector_dot(er, vector_cross(normb, gradB))
return make_xarray(
er_dot_curlb * K_par \
+ er_dot_gradB * K_perp,
norm=norm.n0 * norm.c_s0 * norm.Te0,
)
[docs]
@autodoc_function
def diamagnetic_heat_flux_maxw(
grid: Grid2D,
equi: EquilibriumBaseClass,
norm: Normalization,
dens: xr.DataArray,
temp: xr.DataArray,
charge: float=1.00
):
"""Calculate the diamagnetic heat flux for an equilibrium distribution."""
assert dens.dims == temp.dims
r_norm, z_norm = grid.coords_like(dens)
absB = equi.magfield_absolute(r_norm, z_norm)
normb = equi.magfield_vector(r_norm, z_norm) / absB
er = equi.magfield_vector_radial(r_norm, z_norm, normalize=True)
gradB = grad(absB, grid=grid)
curlb = rot(normb, grid=grid)
er_dot_curlb = vector_dot(er, curlb) / absB
er_dot_gradB = vector_dot(er, vector_cross(normb, gradB)) \
/ absB / absB
prefac = 2.5 * dens * temp * temp / charge
return make_xarray(
prefac * er_dot_curlb + prefac * er_dot_gradB,
norm=norm.n0 * norm.c_s0 * norm.Te0 * norm.rho_s0 / norm.R0,
)