Absorption Calculation

The absorption coefficient \(\alpha\) is the ratio of the power absorbed per unit volume to the power flux carried by the wave. Raytrax evaluates it at every point of the ray trajectory from three ingredients:

1. the polarization vector \(\hat{\boldsymbol{e}}\) of the wave, obtained from the weakly relativistic dielectric tensor,
2. the anti-Hermitian part of the dielectric tensor, obtained from a fully relativistic resonance integral over the electron distribution function, and
3. the power flux vector \(\boldsymbol{F}\), obtained from the cold dielectric tensor.

The sections below derive the absorption coefficient and describe each ingredient.

Loss Density

We start from Poynting's theorem in the absence of free currents,

where \(\boldsymbol{S}\) is the power flux density, \(u\) the energy density, \(\boldsymbol{D}=\varepsilon_0 \boldsymbol{\varepsilon_r}\boldsymbol{E}\), and \(\boldsymbol{B}=\mu_0 \boldsymbol{H}\).

For harmonic fields \(\boldsymbol{E}(\boldsymbol{r},t)=\mathrm{Re}\left[\hat{\boldsymbol{E}}(\boldsymbol{r}) e^{-i \omega t}\right]\), the time average of this equation can be written as

Splitting the relative dielectric tensor into its Hermitian and anti-Hermitian parts, \(\boldsymbol{\varepsilon_r}=\boldsymbol{\varepsilon}_r^H + i \boldsymbol{\varepsilon}_r^A\), where

the time-averaged power loss density (units W/m³) becomes

Both \(\boldsymbol{\varepsilon}_r^H\) and \(\boldsymbol{\varepsilon}_r^A\) are Hermitian matrices, so the quadratic form on the right is real. It is positive in an absorbing plasma.

Power Flux Density

The time average of the power flux density (Poynting vector) is given by the energy density times the group velocity,

Using the Hermitian part of the dispersion tensor \(\boldsymbol{\mathsf{D}}^H=\frac{1}{2}(\boldsymbol{\mathsf{D}} + \boldsymbol{\mathsf{D}}^\dagger)\) (see overview), it can be shown that this quantity can also be written as

Writing \(\hat{\boldsymbol{E}} = \hat E \hat{\boldsymbol{e}}\), where \(\hat{\boldsymbol{e}}\) is the unit polarization vector, and using \(\partial/\partial\boldsymbol{k} = (c/\omega)\,\partial/\partial\boldsymbol{N}\) with the refractive index vector \(\boldsymbol{N} = c\boldsymbol{k}/\omega\), this becomes

The dimensionless vector \(\boldsymbol{F}\) is the power flux vector. In vacuum, \(\hat{\boldsymbol{e}}^* \cdot \boldsymbol{\mathsf{D}} \cdot \hat{\boldsymbol{e}} = 1 - N^2\) for a transverse wave, so \(\boldsymbol{F} = \boldsymbol{N}\) with \(|\boldsymbol{F}| = 1\). This recovers the vacuum Poynting flux \(\varepsilon_0 c |\hat E|^2/2\).

Absorption Coefficient

The power carried by the ray decays as \(dP/ds = -\alpha P\), so the absorption coefficient (units 1/m) is the ratio of the power loss density (units W/m³) to the power flux density (units W/m²),

The magnitude of the electric field drops out, and we obtain

It is implemented in absorption_coefficient.

Optical Depth

The absorption coefficient \(\alpha\) is related to the (dimensionless) optical depth \(\tau\) as

where \(s\) is the arc length, and the linear absorption power density (units W/m) along the ray trajectory is given as

where \(P_0\) is the initial power of the ray.

Polarization Vector

All tensors in this section are expressed in Stix coordinates: the \(z\) axis is aligned with the magnetic field \(\boldsymbol{B}_0\) and the \(x\) axis is chosen such that the refractive index vector lies in the \(x\)–\(z\) plane, \(\boldsymbol{N} = (N_\perp, 0, N_\parallel)\). In these coordinates the dispersion tensor reads

The polarization vector is a solution of \(\boldsymbol{\mathsf{D}} \cdot \hat{\boldsymbol{e}} = 0\). Raytrax fixes the component that dominates for quasi-perpendicular propagation and solves a \(2\times 2\) subsystem for the other two:

• O mode: \(e_z = 1\); \(e_x\) and \(e_y\) follow from the first and second rows.
• X mode: \(e_y = 1\); \(e_x\) and \(e_z\) follow from the first and third rows.

The result is then normalized to \(|\hat{\boldsymbol{e}}| = 1\).

The refractive index \(\boldsymbol{N}\) at which \(\boldsymbol{\mathsf{D}}\) is evaluated comes from the ray trajectory, which is traced with the cold plasma Hamiltonian (see ray tracing). It therefore does not satisfy the dispersion relation of the warm tensor exactly, and the equation that is left out is not satisfied exactly either. Selecting a \(2\times 2\) subsystem gives a well-defined polarization vector in this situation.

It is implemented in polarization.

Weakly Relativistic Dielectric Tensor

Absorption depends on how \(\hat{\boldsymbol{e}}\) projects onto the gyration of resonant electrons (see below). This projection can be small, e.g. for the X mode at the fundamental harmonic. It is then sensitive to thermal corrections to the Hermitian part of the dielectric tensor close to the resonance. For this reason, raytrax follows Travis¹ and computes \(\hat{\boldsymbol{e}}\) from the weakly relativistic dielectric tensor of Krivenski and Orefice² rather than from the cold tensor.

In this approximation the dielectric tensor takes the form

where \(X = \omega_p^2/\omega^2\), \(\mu = m_e c^2 / T_e\), and the elements of \(\boldsymbol{\mathcal{D}}\) are finite sums over harmonics \(|s| \le s_\text{max}\) and finite Larmor radius (FLR) orders \(k \le k_\text{max}\). Each term is a power \(\lambda^{|s|+k-1}\) of the FLR parameter

multiplied by combinations of the Shkarofsky functions \(\mathcal{F}_{q+1/2}\)³, which contain the weakly relativistic resonance denominators. As in Krivenski and Orefice, the electron cyclotron frequency enters with its sign, \(\omega_c < 0\).

Raytrax sets \(s_\text{max} = k_\text{max}\) to Beam.max_harmonic. The FLR expansion order thus is at least as large as the highest resonant harmonic, so the polarization vector contains the leading-order FLR correction at each resonance that is included.

It is implemented in weakly_relativistic_dielectric_tensor.

Relativistic Resonance Integral

The weakly relativistic tensor expands the relativistic factor of the electrons to first order. The anti-Hermitian part, which controls the absorbed power, is sensitive to the shape of the distribution function in the resonant region. Raytrax therefore evaluates it from the fully relativistic expression^{4 5}.

Resonance Condition

We use the normalized momentum \(\boldsymbol{u} = \boldsymbol{p}/(m_e c)\), with components \(u_\parallel\) and \(u_\perp\) relative to \(\boldsymbol{B}_0\), the Lorentz factor \(\gamma = \sqrt{1 + u_\parallel^2 + u_\perp^2}\), and \(Y = |\omega_c|/\omega\). An electron interacts resonantly with the wave at harmonic \(s\) if

Inserting this into the definition of \(\gamma\) gives

For \(|N_\parallel| < 1\) this is an ellipse in the \((u_\parallel, u_\perp)\) plane. It intersects the \(u_\parallel\) axis at

A resonance exists only if \(s^2 Y^2 + N_\parallel^2 > 1\). Since \(\gamma \ge 1\), absorption therefore takes place only where \(s|\omega_c| > \omega\sqrt{1 - N_\parallel^2}\), i.e. on the high-field side of the cold resonance for perpendicular propagation (relativistic downshift). For \(s \le 0\) and \(|N_\parallel| < 1\) the resonance condition has no solution, so only positive harmonics \(s = 1, \dots, s_\text{max}\) contribute.

Anti-Hermitian Part

For an isotropic distribution function \(f(\gamma)\), normalized as \(\int f\, d^3u = 1\), the anti-Hermitian part contracted with the polarization vector is

where \(u_\perp\) is given by the resonance ellipse and

Here \(b = k_\perp v_\perp \gamma / |\omega_c|\) is the Larmor radius argument and \(J_s\) is the Bessel function of the first kind. The product \(\hat{\boldsymbol{e}} \cdot \boldsymbol{\Pi}_s = \sum_i e_i \Pi_{s,i}\) is taken without complex conjugation. Up to normalization, \(|\hat{\boldsymbol{e}} \cdot \boldsymbol{\Pi}_s|^2\) is the quasilinear diffusion coefficient of the resonant electrons⁴. For a general, anisotropic distribution, \(\partial f/\partial\gamma\) is replaced by \(\partial f/\partial\gamma + N_\parallel\, \partial f/\partial u_\parallel|_\gamma\), but the second term vanishes for isotropic \(f\).

Raytrax uses the relativistic Maxwell–Jüttner distribution

where \(K_2\) is the modified Bessel function of the second kind. Since \(\partial f/\partial \gamma < 0\), the result is positive, as required for an absorbing plasma.

In the weakly relativistic limit, the expression above agrees with the anti-Hermitian part of the Krivenski–Orefice tensor. The unit tests check this.

Numerical Evaluation

The integral over \(u_\parallel\) is evaluated with the trapezoidal rule on 1000 points. The interval \([u_-, u_+]\) is restricted to \(|u_\parallel| \le 7\,v_\text{th}/c\), where the Maxwell–Jüttner distribution has decayed by about \(10^{-11}\) relative to its peak. A harmonic whose resonance ellipse does not reach into this region contributes zero.

It is implemented in anti_hermitian_dielectric_form.

Power Flux Vector

The power flux vector \(\boldsymbol{F}\) in the denominator of the absorption coefficient is evaluated with the cold dielectric tensor, and the polarization vector inside the gradient is computed from the cold tensor as well. Differentiating the weakly relativistic tensor with respect to \(\boldsymbol{N}\) near a cyclotron harmonic leads to resonant singularities, which the cold tensor avoids. The thermal corrections to \(\boldsymbol{F}\) are of order \(v_\text{th}/c\). This approximation is least accurate at the fundamental harmonic, in particular for X mode, where the thermal corrections to the Hermitian part of the dielectric tensor are largest.

It is implemented in cold_power_flux_vector_stix.

Summary of Plasma Models

Different parts of the calculation use different plasma models:

Because the trajectory is traced with the cold Hamiltonian, \(\boldsymbol{N}\) does not include the thermal modification of the refractive index near the resonance. See limitations.