Physics background

This page documents the physical model behind StarkZee. The equations follow the notation of Ferri, Peyrusse & Calisti (2022) [Ferri2022].

Liouville space representation

The spectral line-shape intensity \(I(\omega)\) is the Fourier transform of the dipole autocorrelation function:

\[I(\omega) = \frac{1}{\pi}\,\mathrm{Re} \int_0^\infty C(t)\,e^{i\omega t}\,dt\]

In Liouville space notation:

\[C(t) = \langle\langle \vec{d}^{\,*} | U(t) | \vec{d}\,\rho_0 \rangle\rangle\]

where \(\vec{d}\) is the electric transition dipole operator, \(\rho_0\) the equilibrium density matrix of the initial manifold, and \(U(t)\) the time-evolution propagator. In the anti-symmetric subspace the Liouvillian is

\[L = \frac{1}{\hbar}(H_u \otimes I^d - I \otimes H_l^d)\]

whose diagonal elements are the transition frequencies and whose off-diagonal elements encode the Stark couplings.

Radiator Hamiltonian

For a hydrogen-like radiator with nuclear charge \(Z\) and principal quantum number \(n\) in an external field \(\vec{B} = B\hat{z}\), the atomic Hamiltonian in the uncoupled \(|n,l,m_l,s,m_s\rangle\) basis is

\[H_A = H_0 + V_\text{SO} + H_Z^{(1)} + H_Z^{(2)}\]

Unperturbed energy (degenerate across the shell):

\[H_0 = -\frac{Z^2\,\text{Ry}}{n^2}\]

Spin-orbit coupling \(V_\text{SO} = \xi\,\vec{L}\cdot\vec{S}\):

\[\xi = \frac{Z^4\,\alpha^2\,\text{Ry}}{n^3\,l\,(l+\tfrac{1}{2})\,(l+1)}\]

where \(\alpha\) is the fine-structure constant. Together with the mass-velocity and Darwin corrections (also included when fine_structure=True), this reproduces the Dirac fine-structure splitting and restores the \(2s_{1/2} = 2p_{1/2}\) degeneracy.

Linear Zeeman (diagonal):

\[H_Z^{(1)} = \mu_B B\,(m_l + g_s\,m_s), \qquad g_s \approx 2.0023192\]

Quadratic (diamagnetic) Zeeman:

\[H_Z^{(2)} = \frac{e^2 B^2}{8 m_e}\,r^2\sin^2\theta\]

\(H_Z^{(2)}\) is non-diagonal in \(l\) (connects \(\Delta l = 0\) and \(\Delta l = \pm 2\) states with the same \(m_l\), \(m_s\)). The off-diagonal radial integrals

\[\langle n,l_1 | r^2 | n,l_2 \rangle, \qquad |l_1 - l_2| = 2\]

do not reduce to a simple product of diagonal elements. Using the geometric-mean approximation \(\sqrt{\langle r^2\rangle_{l_1}\langle r^2\rangle_{l_2}}\) overestimates the true value by 12 % (\(n=3\)), 21 % (\(n=4\)), and 41 % (\(n=5\)). StarkZee computes these integrals by direct numerical quadrature (scipy.integrate.quad) over the analytical hydrogenic radial wavefunctions

\[R_{nl}(r) = \sqrt{\left(\frac{2Z}{n}\right)^{\!3} \frac{(n-l-1)!}{2n\,(n+l)!}} \;e^{-Zr/n}\!\left(\frac{2Zr}{n}\right)^{\!l} L_{n-l-1}^{2l+1}\!\!\left(\frac{2Zr}{n}\right)\]

with results cached after the first call.

Stark perturbation

Under the quasi-static ion approximation the radiator sits in a constant microfield \(\vec{F}\) at angle \(\theta\) to \(\vec{B}\):

\[F_z = F\cos\theta, \qquad F_x = F\sin\theta\]

The Stark interaction is

\[V_E = -e\,(z\,F_z + x\,F_x)\]

Within the \(n\)-shell the radial element is analytic:

\[\langle n,l | r | n,l-1 \rangle = \frac{3n}{2Z}\sqrt{n^2 - l^2} \quad [a_0]\]

The combined Hamiltonian \(H = H_A + V_E\) is diagonalized by numpy.linalg.eigh at every quadrature point to obtain the Stark-dressed transition frequencies \(\omega_k\) and dipole weights \(|d_k|^2\).

Plasma microfield distributions

The ion microfield is averaged over a probability distribution \(W(F)\).

Holtsmark (unscreened) [Holtsmark1919]:

\[W_H(\beta) = \frac{2\beta}{\pi} \int_0^\infty y\,\sin(\beta y)\,e^{-y^{3/2}}\,dy, \qquad \beta = F/F_0\]

where \(F_0 = e/(4\pi\varepsilon_0 r_e^2)\) is the Holtsmark normal field and \(r_e = (3/4\pi N_e)^{1/3}\) is the mean inter-particle distance.

Hooper (Debye-screened, default) [Hooper1968]:

\[W(\beta, a) = \frac{2\beta}{\pi} \int_0^\infty y\,\sin(\beta y)\, e^{-y^{3/2}\,S(y,a)}\,dy\]

where the screening function is

\[S(y, a) = \left(1 + \frac{1.5\,a^2}{y^2}\right)^{-3/4}, \qquad a = r_e / \lambda_D\]

and \(\lambda_D = \sqrt{\varepsilon_0 T_e / N_e e^2}\) is the electron Debye length. Setting \(a = 0\) recovers the Holtsmark distribution.

The profile integrator uses a uniform grid of num_f points in \([0, 10\,F_0]\) and num_mu Gauss-Legendre points for the angle \(\mu = \cos\theta \in [0, 1]\).

Electron impact broadening (GBK)

Fast electrons are treated in the impact (completed-collision) approximation [Griem1997]. Their contribution is a homogeneous Lorentzian broadening of half-width (HWHM):

\[W_e(\Delta\omega) = W_\text{pref}\,\langle r^2\rangle_n \bigl[C_n + G(\Delta\omega)\bigr]\]

where the prefactor is

\[W_\text{pref} = \frac{4\pi}{3}\,N_e \sqrt{\frac{2 m_e}{\pi k_B T_e}}\, \left(\frac{\hbar}{m_e}\right)^{\!2}\]

and the shell-averaged mean-square radius is

\[\langle r^2\rangle_n = \frac{1}{n^2} \sum_{l=0}^{n-1}(2l+1)\,\frac{n^2}{2Z^2} \bigl[5n^2+1-3l(l+1)\bigr]\,a_0^2 \;\propto\; \frac{n^4}{Z^2}\]

The GBK dynamical factor

\[G(\Delta\omega) = \tfrac{1}{2}\,E_1(y), \qquad y = \left(\frac{n^2}{2Z}\right)^{\!2} \frac{\Delta\omega^2 + \omega_c^2}{2\,\text{Ry}\cdot T_e}\]

uses the exponential integral \(E_1\) and the cutoff frequency

\[\omega_c = \max(\omega_p,\;\omega_L,\;\omega_e)\]

where \(\omega_p = \sqrt{N_e e^2/\varepsilon_0 m_e}\) is the plasma frequency, \(\omega_L = eB/m_e\) the electron Larmor frequency (dominant at high \(B\)), and \(\omega_e = 2\pi v_\text{th}/r_e\) the configuration-change rate.

Each Stark-dressed transition component \((i \to j, q)\) is broadened by a Lorentzian of half-width \(\gamma_e\):

\[L_{ij}(\omega) = \frac{\gamma_e/\pi}{(\omega - \omega_{ij})^2 + \gamma_e^2}\]

and its contribution is weighted by \(|d_q(i,j)|^2\) and the microfield quadrature weight.

By default (frequency_dependent_width=True) \(\gamma_e\) is evaluated at the actual detuning of each component. Setting frequency_dependent_width=False fixes it at the line-center value \(\gamma_e(0)\), which is faster but less accurate in the far wings.

The strong-collision constants \(C_n\) from Ferri et al. (2022):

n

C_n

≤ 2

1.50

3, 4

0.75

≥ 5

0.40

Thermal Doppler broadening

Thermal motion of the radiating ions Doppler-shifts each photon frequency by \(\delta\omega = \omega_0\,v_z/c\). Averaging over a Maxwell–Boltzmann velocity distribution gives a Gaussian with \(1/e\) half-width

\[\Delta E_D = E_0\sqrt{\frac{2\,T_i}{m_\text{ion}\,c^2/e}}\]

For H Balmer-\(\alpha\) at \(T_i = 5\) eV this gives \(\Delta E_D \approx 0.062\) meV — comparable to the Zeeman splitting at 1 T and negligible relative to the Stark width at \(N_e \sim 10^{23}\) m-3.

Doppler broadening is not applied automatically. The static profile solver returns the pure Stark-Zeeman result; convolutions are applied afterwards via starkzee.convolutions:

\[I_\text{obs}(\lambda) = \bigl[I_\text{SZ}(\lambda) * G_D(\lambda)\bigr] * G_\text{inst}(\lambda)\]

where \(G_D\) is the Doppler Gaussian and \(G_\text{inst}\) the instrumental slit function. Both operate on a uniform wavelength grid; convert energy grids to wavelength before calling.

Frequency Fluctuation Model

When ion dynamics are important [Talin1995] (high \(N_e\), low \(B\), or high-\(n\) lines) the static-ion approximation overestimates the central peak height. The FFM treats the microfield as a Markov jump process switching between field configurations at rate \(\nu_i\).

The line profile is (Sherman–Morrison form):

\[I(\omega) = \frac{r^2}{\pi}\,\mathrm{Re} \frac{S(\omega)}{1 - \nu_i\,S(\omega)}, \qquad S(\omega) = \sum_k \frac{p_k}{\nu_i + \gamma_k + i(\omega - \omega_k)}\]

where \(p_k = |d_k|^2 / r^2\) are the normalized SDT weights, \(\gamma_k\) the electron-impact half-width, and \(r^2 = \sum_k |d_k|^2\). The ion fluctuation rate is

\[\nu_i = \frac{v_\text{th}}{r_i}\,\frac{\hbar}{e}, \qquad v_\text{th} = \sqrt{\frac{2 T_i}{m_i}}\]

with \(r_i = (3/4\pi N_i)^{1/3}\) the mean ion-sphere radius.

The limit \(\nu_i \to 0\) recovers the static profile; the limit \(\nu_i \to \infty\) gives a single Lorentzian (motional narrowing).

Observation geometry

The intensity observed at angle \(\theta\) to \(\vec{B}\) is

\[I(\theta) = I_\pi\sin^2\theta + \tfrac{1}{2}(I_{\sigma+} + I_{\sigma-})(1+\cos^2\theta)\]

Useful special cases:

  • Transverse (\(\theta = 90°\)): \(I_\pi + \tfrac{1}{2}(I_{\sigma+} + I_{\sigma-})\)

  • Along B (\(\theta = 0°\)): \(I_{\sigma+} + I_{\sigma-}\)

  • Angle-averaged: \(\tfrac{2}{3}I_\pi + \tfrac{1}{3}(I_{\sigma+} + I_{\sigma-})\)

At \(B = 0\) the quantization axis is undefined; all three polarization components are equal by spherical symmetry.

Oscillator strengths and line strengths

The line strength summed over all polarizations and substates is

\[S_{ul} = \sum_{q,i,j} |\langle l_j | r_q | u_i \rangle|^2 \quad [a_0^2]\]

The weighted absorption oscillator strength is

\[gf = \frac{2}{3}\,\frac{\Delta E}{E_\text{H}}\,S_{ul}\]

and the Einstein A coefficient

\[A_{ul} = \frac{4\,\alpha^3}{3}\, \left(\frac{\Delta E}{E_\text{H}}\right)^{\!3} \frac{S_{ul}}{g_u\,\tau_\text{au}}\]

where \(E_\text{H} = 2\,\text{Ry}\) is the Hartree energy, \(g_u = 2n_u^2\) the upper-shell statistical weight, and \(\tau_\text{au} = \hbar / E_\text{H} \approx 2.419 \times 10^{-17}\) s the atomic unit of time. Both \(gf\) and \(A_{ul}\) are validated against NIST ASD tabulated values to 0.5 % and 1 % respectively.

Physical features and validation

Quadratic Zeeman polarization wings

At \(B \geq 500\) T the diagonal QZ shifts within \(n=3\) differ by up to 5 meV, separating the \(3p\,(m_l=\pm 1)\to 2s\) transitions (which carry ≈ 21 % of the H\(\alpha\) oscillator strength and the largest \(n=3\) QZ shift, \(+8.87\) meV at \(B=1000\) T) from the dominant \(3d\to 2p\) cluster by ≈ 7 meV. The resulting wings appear only in quadratic_zeeman=True profiles and are validated in test_12_halpha_qz_wings.py.

±2μBB Stark-Zeeman satellite features

A transverse microfield \(F_x\) couples adjacent-\(m_l\) states within the same principal shell via \(\Delta l = \pm 1,\;\Delta m_l = \pm 1\) matrix elements. The upper eigenstate near \(E_{0,n} + 2\mu_B B\) acquires a small admixture of the neighbouring \(np\) state, producing a \(\sigma^+\) transition to the zero-Zeeman lower state at photon energy \(E_0 + 2\mu_B B\).

H\(\beta\) (n=4→2). The satellite is a distinct peak at +117 meV (≈ 2 % of the main \(\sigma^+\) intensity), because \(n=4\) has both \(|4d,\,m_l{=}2\rangle\) and \(|4f,\,m_l{=}2\rangle\) degenerate at \(+2\mu_B B\), greatly amplifying Stark mixing.

H\(\alpha\) (n=3→2). Only \(|3d,\,m_l{=}2\rangle\) sits at \(+2\mu_B B\) (no \(l=3\) substates exist for \(n=3\)). The satellite amplitude (≈ 0.07 % of main peak) is buried in the Lorentzian tail of the \(\sigma^+\) main peak, so no distinct local maximum appears.

The mixing coefficient \(\beta = F_x\langle np|r|nd\rangle/(\mu_B B)\) scales as \(B^{-1}\), so \(\beta\) is larger at lower \(B\). However, the satellite position \(2\mu_B B\) also shrinks with \(B\). At \(B = 100\) T, \(2\mu_B B = 11.6\) meV falls inside the Stark-broadened \(\sigma^+\) cluster; the satellite merges into the tail and is not resolved. High \(B\) is required for the satellite to appear as a separated peak. Validated in test_13_stark_zeeman_satellites.py.

Numerical simplifications

  • Within-shell Stark matrix — the quadratic Stark effect (coupling to \(n \pm 1\) shells) is neglected. This is valid when the Stark shift \(\ll Z^2\,\text{Ry}(1/n^2 - 1/(n+1)^2)\).

  • Resonance-center impact width (optional) — setting frequency_dependent_width=False fixes \(\gamma_e\) at its line-center value \(\gamma_e(0)\), avoiding repeated \(E_1\) evaluations. Faster but less accurate in the wings.

  • Gauss-Legendre angle quadrature — the integration over microfield direction uses num_mu Gauss-Legendre points on \(\mu = \cos\theta \in [0, 1]\) (default 6).

  • Cached radial integrals\(\langle r\rangle\) and \(\langle r^2\rangle\) matrix elements are computed once and cached with functools.lru_cache.

References

[Ferri2022]

S. Ferri, O. Peyrusse, A. Calisti, Matter and Radiation at Extremes 7, 015901 (2022).

[Talin1995]

B. Talin, A. Calisti, L. Godbert, R. Stamm, R. W. Lee, Phys. Rev. A 51, 1918 (1995).

[Holtsmark1919]
  1. Holtsmark, Ann. Phys. 58, 577 (1919).

[Hooper1968]
    1. Hooper, Phys. Rev. 165, 215 (1968).

[Griem1997]

H. R. Griem, Principles of Plasma Spectroscopy, Cambridge University Press (1997).