Three-wave interactions in magnetized warm-fluid plasmas: general theory with evaluable coupling coefficient

TL;DR

This paper derives a general, evaluable formula for three-wave coupling coefficients in magnetized warm-fluid plasmas, enabling precise analysis of wave interactions in contexts like fusion and solar physics.

Contribution

It introduces a novel, general formula for coupling coefficients in magnetized warm-fluid plasmas, applicable to arbitrary wave angles and interactions.

Findings

Coupling coefficients for laser scattering in magnetized plasmas are calculated.

Growth rates of laser scattering are affected by magnetic fields and oblique angles.

Decay of Alfvén waves shows slight preference for backward geometry in solar conditions.

Abstract

Resonant three-wave coupling is an important mechanism via which waves interact in a nonlinear medium. When the medium is a magnetized warm-fluid plasma, a previously-unknown formula for the coupling coefficients is derived by solving the fluid-Maxwell's equations to second order using multiscale perturbative expansions. The formula is not only general but also evaluable, whereby numerical values of the coupling coefficient can be determined for any three resonantly interacting waves propagating at arbitrary angles. As one example, coupling coefficient governing laser scattering is evaluated. In conditions relevant to magnetized inertial confinement fusion experiments, lasers scatter from magnetized plasma waves and the growth rates are modified at oblique angles. As another example, coupling coefficient between two Alfv\'en waves via a sound wave is evaluated. In conditions relevant to solar corona, the decay of a parallel Alfv\'en wave only slightly prefers exact backward geometry.

Figures

Equations

$$\frac{\partial n_s}{\partial t}+\nabla \cdot (n_s{\mathbf{v}}_s)=0,$$

$$m_{s}n_{s}\Big{(}\frac{\partial\mathbf{v}_{s}}{\partial t}+\mathbf{v}_{s}\cdot\nabla\mathbf{v}_{s}\Big{)}=-\nabla p_{s}+e_{s}n_{s}(\mathbf{E}+\mathbf{v}_{s}\times\mathbf{B}),$$

$$n_s{d}_t{p}_s={\xi }_s{p}_s{d}_t{n}_s,$$

$$\frac{\partial \mathbf{B}}{\partial t}=-\nabla \times \mathbf{E},$$

$$\frac{\partial \mathbf{E}}{\partial t}=c^2\nabla \times \mathbf{B}-\frac{1}{{ϵ}_0}\sum _s{e}_s{n}_s{\mathbf{v}}_s,$$

$$U=\frac{1}{2}\epsilon_{0}\mathbf{E}^{2}+\frac{1}{2\mu_{0}}\mathbf{B}^{2}+\sum_{s}\frac{1}{2}\Big{(}m_{s}n_{s}\mathbf{v}_{s}^{2}+\frac{p_{s}}{\xi_{s}-1}\Big{)},$$

$$\mathbf{S}=\frac{1}{\mu_{0}}\mathbf{E}\times\mathbf{B}+\sum_{s}\frac{1}{2}\Big{(}m_{s}n_{s}\mathbf{v}_{s}^{2}\mathbf{v}_{s}+\frac{\xi_{s}}{\xi_{s}-1}p_{s}\mathbf{v}_{s}\Big{)},$$

$${\partial }_{t}U+\nabla \cdot \mathbf{S}=0,$$

$$\mathbf{Z}$$

$${\partial }_t$$

$${\partial }_i$$

$${\partial }_{t_0}{n}_{s1}$$

$$m_s{n}_{s0}{\partial }_{t_0}{\mathbf{v}}_{s1}$$

$$n_{s0}{\partial }_{t_0}{p}_{s1}$$

$${\partial }_{t_0}{\mathbf{B}}_1$$

$${\partial }_{t_0}{\mathbf{E}}_1$$

$${\partial }_{t_0}{n}_{s2}+n_{s0}{\nabla }_0\cdot {\mathbf{v}}_{s2}$$

$$m_s{n}_{s0}{\partial }_{t_0}{\mathbf{v}}_{s2}+{\nabla }_0{p}_{s2}-e_s{n}_{s0}({\mathbf{E}}_2+{\mathbf{v}}_{s2}\times {\mathbf{B}}_0)$$

$$n_{s0}{\partial }_{t_0}{p}_{s2}-{\xi }_s{p}_{s0}{\partial }_{t_0}{n}_{s2}$$

$${\partial }_{t_0}{\mathbf{B}}_2+{\nabla }_0\times {\mathbf{E}}_2=-{\partial }_{t_1}{\mathbf{B}}_1-{\nabla }_1\times {\mathbf{E}}_1,$$

$${\partial }_{t_0}{\mathbf{E}}_2-c^2{\nabla }_0\times {\mathbf{B}}_2+\frac{1}{{ϵ}_0}\sum _s{e}_s{n}_{s0}{\mathbf{v}}_{s2}$$

$${\mathbf{E}}_1$$

$${\mathbf{B}}_1$$

$$p_{s1}={\epsilon }_s{n}_{s1},$$

$${\mathbf{v}}_{s1}$$

$$n_{s1}$$

$${\hat{\mathbb{F}}}_{s,\mathbf{k}}\mathbf{Z}=\mathbf{Z}+i{\beta }_{s,\mathbf{k}}({\hat{\mathbb{F}}}_{s,\mathbf{k}}\mathbf{Z})\times \mathbf{b}+\frac{u_s^2}{{\omega }_{\mathbf{k}}^2}\mathbf{k}(\mathbf{k}\cdot {\hat{\mathbb{F}}}_{s,\mathbf{k}}\mathbf{Z}),$$

$$\hat{\mathbb{F}}=\mathbb{FP}={\mathbb{P}}^{†}\mathbb{F},$$

$$\mathbb{F}\mathbf{Z}={\gamma }^2[\mathbf{Z}+i\beta \mathbf{Z}\times \mathbf{b}-{\beta }^2(\mathbf{Z}\cdot \mathbf{b})\mathbf{b}],$$

$$\mathbb{P}\mathbf{Z}=\mathbf{Z}+{\hat{\gamma }}^2\frac{u^2}{{\omega }^2}\mathbf{k}(\mathbf{k}\cdot \mathbb{F}\mathbf{Z}),$$

$${\mathbb{D}}_{\mathbf{k}}^{ij}=({\omega }_{\mathbf{k}}^2-c^2{\mathbf{k}}^2){\delta }^{ij}+c^2{k}^i{k}^j-\sum _s{\omega }_{ps}^2{\hat{\mathbb{F}}}_{s,\mathbf{k}}^{ij}.$$

$${\chi }_{s,\mathbf{k}}=-\frac{{\omega }_{ps}^2}{{\omega }_{\mathbf{k}}^2}{\hat{\mathbb{F}}}_{s,\mathbf{k}}.$$

$$\displaystyle\hat{\mathcal{F}}\!=\!\left(\begin{array}[]{ccc}\gamma^{2}(1\!+\!\gamma^{2}\rho^{2}\mathrm{s}^{2}_{\theta})\!&\!i\beta\gamma^{2}(1\!+\!\gamma^{2}\rho^{2}\mathrm{s}^{2}_{\theta})\!&\!\gamma^{2}\rho^{2}\mathrm{s}_{\theta}\mathrm{c}_{\theta}\\ -i\beta\gamma^{2}(1\!+\!\gamma^{2}\rho^{2}\mathrm{s}^{2}_{\theta})\!&\!\gamma^{2}(1\!+\!\beta^{2}\gamma^{2}\rho^{2}\mathrm{s}^{2}_{\theta})\!&\!-i\beta\gamma^{2}\rho^{2}\mathrm{s}_{\theta}\mathrm{c}_{\theta}\\ \gamma^{2}\rho^{2}\mathrm{s}_{\theta}\mathrm{c}_{\theta}\!&\!i\beta\gamma^{2}\rho^{2}\mathrm{s}_{\theta}\mathrm{c}_{\theta}\!&\!1\!+\!\rho^{2}\mathrm{c}^{2}_{\theta}\end{array}\right)\!,$$

Full text

I Introduction

Plasmas are dielectric media wherein waves can interact nonlinearly. Unlike crystals whose optical properties may have limited range of tunability, plasma parameters can vary by many orders of magnitude. In particular, an adjustable optical axis can be introduced by applying an external magnetic field. Thereby, all nonlinear optical phenomena seen in crystals Bloembergen (1996) also occur in magnetized plasmas with ample flexibility. In addition to hosting optical phenomena, magnetized plasmas also support a zoo of other waves. These additional waves, such as the Alfvén wave, Bernstein waves, and hybrid waves, not only mediate new interactions between light waves, but also couples nonlinearly among themselves. For example, interactions between Alfvén waves is thought to be a major mechanism for anisotropic turbulence Sridhar and Goldreich (1994); Ng and Bhattacharjee (1996); Goldreich and Sridhar (1997); Chandran (2008); Chaston et al. (2008); Gogoberidze (2007); Schekochihin et al. (2012) and particle heating Voitenko and Goossens (2005); Araneda et al. (2008) in astrophysical plasmas.

While nonlinear wave coupling occurs in any dielectric medium, what makes the biggest difference is perhaps the coupling coefficient. When the coupling is weak, very large amplitude waves are needed in order to cause sizable effects. On the contrary, when the coupling is strong, even small amplitude waves can lead to observable consequences. Since plasma parameters span many orders of magnitude, it is impractical to exhaust the multidimensional parameter space by experiments and first-principle simulations. An analytical formula, which can be used to determine numerical values of the coupling coefficient, is thereof invaluable for mapping out wave-wave coupling behaviors in magnetized plasmas.

For over half a century, numerous attempts are made to calculate wave coupling in magnetized plasmas due to three-wave interactions, which are the leading-order terms of the nonlinear response tensor Sagdeev and Galeev (1969); Davidson (1972). However, the presence of a background magnetic field significantly complicates the calculation, and most attempts start by restricting to a particular set of wave triad in some special geometry. For example, theories have been developed when waves propagate perpendicular to the magnetic field Platzman et al. (1968); Stenflo (1972), and explicit results have been obtained when the pump is the extraordinary wave Grebogi and Liu (1980a, b); Barr et al. (1984); Boyd and Rankin (1985); Simon (1995); Mourenas (1996); Vyas et al. (2016); Dodin and Arefiev (2017), the ordinary wave Purohit et al. (2010); Shi et al. (2017a), the upper-hybrid wave Ram (1982), and the lower-hybrid wave Sanuki and Schmidt (1977). For wave propagation nearly parallel to the magnetic field, transverse and longitudinal modes decouple Sjölund and Stenflo (1967a); Stenflo (1970), and results have been obtained when all waves are electrostatic Shivamoggi (1982), when the pump wave is a circularly polarized laser Laham et al. (1998), the whistler wave Kumar and Tripathi (2011), the fast wave Voitenko and Goossens (2002); Modi and Sharma (2013), and the Alfvén waves Hasegawa and Chen (1976a, b); Erokhin et al. (1978); Derby Jr (1978); Goldstein (1978); Wong and Goldstein (1986); Jayanti and Hollweg (1993a); Matsukiyo and Hada (2003); Sweeney and Stewart (1978); Brodin and Stenflo (1988, 1990); Viñas and Goldstein (1991a); Jayanti and Hollweg (1993b); Hollweg (1994); Voitenko (1998); Voitenko and Goossens (2000); Shukla et al. (2004); Ruderman and Simpson (2004); Brodin et al. (2006); Nariyuki and Hada (2007). Although more general theories exist Galloway and Kim (1971); Boyd and Turner (1978); Liu and Tripathi (1986); Stenflo (1994); Brodin and Stenflo (2012), the formal expressions of the coupling coefficient are too cumbersome to be useful and are rarely evaluated in practice Viñas and Goldstein (1991b). Moreover, in order to simplify results, increasing number of assumptions are usually made as the discussion progresses, and conflicting assumptions have led to numerous disputes in the literature.

In order to obtain a formula for the coupling coefficient that is not only general but also evaluable, a mathematically robust approach is necessary. In a previous paper Shi et al. (2017b), an approached based on multiscale perturbative solution has been demonstrated for magnetized cold-fluid plasmas. The key to simplifying the general result is not to make additional assumptions, but to package seemingly complex terms into well-motivated operators. By studying properties of these mathematical operators, profound simplifications can then be unveiled, which would otherwise be buried under tedious arithmetics. This approach is not only useful for obtaining analytical expressions, but also necessary to avoid brute-force manipulation of large matrices during numerical evaluations.

In this paper, I will further demonstrate the operator approach by considering three-wave interactions in magnetized warm-fluid plasmas. The ideal warm-fluid model is applicable when the wavelengths of interest are much longer than the Debye length, while much shorter than the collisional mean free path. In this regime, plasma particles respond to perturbations collectively and dissipative effects are small. For fusion and astrophysical plasmas, the fluid-Maxwell model has a reasonable range of applicability. For example, in inertial confinement fusion conditions, the plasma density $n\sim 10^{20}\,\text{cm}^{-3}$ and temperature $T\sim 1$ keV. Correspondingly, the Debye length $\lambda_{D}\sim 10^{-2}\,\mu\text{m}\,(T/n)^{1/2}$ is usually much smaller than the laser wavelength, which is in turn much smaller than the collisional mean free path $\lambda_{\text{mfp}}\sim 10\,\mu\text{m}\,(nZ^{2})^{-1}$ in low-$Z$ plasmas.

Within the range of its applicability, the fluid model may then be solve perturbatively when amplitudes of fluctuations are small. To obtain solutions beyond the linear order, special procedures are necessary in order to avoid secular behaviors. A well-suited procedure is multiscale expansion, which expand spatial and temporal scales in addition to expanding the amplitudes. By renormalizing the spacetime, well-behaved high-order perturbative solutions can then be obtained. In this weak-coupling regime, no further assumption is needed, and the model equations can be solved using rigorous procedures to study interactions between arbitrary waves in the most general geometry under a wide variety of plasma conditions.

This paper is organized as followes. In Sec. II, the fluid model and the multiscale method will be reviewed. In sec. III, I will introduce important operators and review linear waves from the operator perspective. In Sec. IV, I will derive the coupling coefficient by solving the second-order equations. In Sec. V, known results in the literature will be recovered, and evaluation of the general formula will be demonstrated using two examples. Discussion is made in Sec. VI followed by a summary. Supplemental details are provided in the Appendix.

II Warm-fluid model

The fluid model describes plasma species as charged gases, which couple with self-consistent electromagnetic fields through the Lorentz force law and the Maxwell’s equations. The multi-fluids model can be regarded as moments of the kinetic model, and can be used to obtain magnetohydrodynhamics (MHD) models after further simplifications.

II.1 Fluid-Maxwell’s equations

For each plasma species, its density evolves according to the continuity equation. In the absence of ionization and recombination, the number of particles is conserved, and the continuity equation is

[TABLE]

where $n_{s}$ is the density of species $s$, whose fluid velocity is $\mathbf{v}_{s}$. The continuity equation contains a nonlinear term $n_{s}\mathbf{v}_{s}$, which will contribute to wave-wave couplings.

The fluid velocity evolves according to the momentum equation. Using the continuity equation, the nonrelativistic momentum equation can be written as

[TABLE]

where $m_{s}$ and $e_{s}$ are the mass and charge of each particle of species $s$, whose thermal motion leads to a pressure $p_{s}$. The above is the simplest momentum equation for warm plasmas, assuming collisions play negligible role, and the internal stress tensor $\tau_{ij}=-p\delta_{ij}$ remains isotropic despite of external forces.

To close the infinite hierarchy of fluid equations, we can express the pressure in terms of density and velocity. For simplicity, consider polytropic process for which $pV^{\xi}$ is a constant, where $\xi\geq 0$ is the polytropic index, $p$ is the pressure, and $V$ is the volume of the fluid element. Suppose the number of particles in the fluid element is constant, then the polytropic condition relates changes of pressure and density by

[TABLE]

where $d_{t}=\partial_{t}+\mathbf{v}_{s}\cdot\nabla$ is the convective derivative. The polytropic process assumes that the heat to work ratio is a constant. In particular, the polytropic process recovers the isobaric process when $\xi=0$; the isothermal process when $\xi=1$; the isochoric process when $\xi=\infty$; and the adiabatic process when $\xi=C_{p}/C_{v}$, where $C_{p}$ and $C_{v}$ are heat capacities at constant pressure and volume.

To model plasmas with self-consistent electric and magnetic fields, we can couple the fluid equations with the Maxwell’s equations. The time evolution of the magnetic field is given by the Faraday’s law:

[TABLE]

which is independent of plasma dynamics. In comparison, the time evolution of the electric field is given by the Ampère’s law:

[TABLE]

where $\epsilon_{0}$ is the vacuum permittivity. The other two Maxwell’s equations $\nabla\cdot\mathbf{E}=\sum_{s}e_{s}n_{s}/\epsilon_{0}$ and $\nabla\cdot\mathbf{B}=0$ are guaranteed once they are satisfied at the initial time.

The fluid-Maxwell system satisfies local energy-momentum conservation laws. The total energy density of the system is

[TABLE]

which is the sum of the field, the kinetic, and the thermal energy densities. Similarly, the energy flux is

[TABLE]

which is constituted of the Poynting flux, the kinetic flux, and the thermal flux. The local conservation law is

[TABLE]

which can be verified by straightforward calculations using Eqs. (1)-(5). Analogously, one can show that the local momentum is also conserved: $\partial_{t}\Pi_{j}+\partial_{i}\sigma_{ij}=0$, where $\Pi_{i}$ is the momentum density, and $\sigma_{ij}$ is the stress tensor.

II.2 Multiscale perturbative expansions

The fluid-Maxwell’s equations are a set of nonlinear partial differential equations. The equations self-consistently determine the fluid variables $n_{s}$, $\mathbf{v}_{s}$, and $p_{s}$, as well as the field variables $\mathbf{E}$ and $\mathbf{B}$. Although calculating the general solution is difficult, perturbative solutions may be obtained when fluctuations have small amplitudes. Since the equations are nonlinear, special procedures are needed in order to remove secular behaviors beyond the leading order. Once secular behaviors are removed, the perturbative solutions are well-behaved without violating the small-amplitude assumption.

Consider perturbations from an equilibrium state of the plasma. Then, a generic fluid or field variable $\mathbf{Z}$ can be expanded as

[TABLE]

Here, the equilibrium state is labeled by the subscript “0”, and $\lambda$ is an auxiliary smallness parameter. Notice that at this step, it is not necessary to assumed any property of higher order terms $\mathbf{Z}_{j}$. In particular, the average $\langle\mathbf{Z}_{j}\rangle$ needs not be zero, and the equilibrium state $\mathbf{Z}_{0}$ is not necessarily the averaged quantity.

One way of removing secular behavior from the perturbative solution is to also expand temporal and spatial scales. Using multiscale expansions, the time and space derivatives are

[TABLE]

where $\lambda$ is the same auxiliary expansion parameter. The multiscale expansion assumes that weaker interactions occur on slower time scales and larger spatial scales. This is intuitive because weaker couplings require further accumulations before their effects become appreciable. In the multiscale expansion, $t_{1}$ is a time scale slower than $t_{0}$ by a factor of $\lambda$, and processes that occur on $t_{1}$ scale is assumed to be well separated from processes on $t_{0}$ scale. Similarly, other temporal and spatial scales are ordered by $\lambda$ and are assumed to be independent.

Consider the simplest equilibrium where the plasma is uniform and stationary under a constant background magnetic field. In this case, $\mathbf{E}_{0}$ and $\mathbf{v}_{s0}$ are zero, whereas $\mathbf{B}_{0}$, $n_{s0}$, and $p_{s0}$ are nonzero but constant. Although this simple situation is rarely encountered in realistic plasmas, it provides a reasonable simplification when the scales of inhomogeneities are well-separated from characteristic scales of three-wave interactions.

When perturbed from the simple equilibrium, the $\lambda$-order equations are homogeneous linear partial differential equations with constant coefficients. The linearized fluid equations are

[TABLE]

which couple $n_{s1}$, $\mathbf{v}_{s1}$, and $p_{s1}$ in pairs. The linearized Maxwell’s equations are

[TABLE]

These linear partial differential equations can be easily solves in the Fourier space to obtain the full spectrum of linear waves in magnetized warm-fluid plasmas.

To compute three-wave coupling between linear waves, we need to solve the equation to the next order. The $\lambda^{2}$-order continuity equation is

[TABLE]

the $\lambda^{2}$-order momentum equation is

[TABLE]

and the $\lambda^{2}$-order pressure equation is

[TABLE]

The above equations may be simplified using $\lambda$-order equations, as will be done later in Sec. IV.1. Similarly, we can write down the $\lambda^{2}$-order Maxwell’s equation. The second-order Faraday’s law is

[TABLE]

and the second-order Ampère’s law is

[TABLE]

Although these equations may look complicated, they are in fact linear equations for fluid variables $n_{s2}$, $\mathbf{v}_{s2}$, and $p_{s2}$, as well as field variables $\mathbf{E}_{2}$ and $\mathbf{B}_{2}$. Moreover, these second order variables couple in exactly the same way as in the $\lambda$-order equations. The only difference is the presence of source terms, which I have arranged to appear on the right-hand-sides (RHS) of the above equations. Once the first-order variables are solved from the $\lambda$-order equations, these source terms can be regarded as known. The above $\lambda^{2}$-order equations are then a system of inhomogeneous linear partial differential equations, which can be solved again in the Fourier space.

III Magnetized linear waves

Before discussing three-wave interactions, it is useful to familiarize with linear waves in magnetized warm-fluid plasmas. In this section, I will review the eigenvalues, the eigenvectors, and the eigenenergies of linear waves. During this review, I will also introduce important concepts that will become indispensable in the next section.

III.1 First-order equations

Now let us solve the first-order equations. Since the equations are linear, the general solution is a superposition of plane waves. In particular, the first-order electric and magnetic fields are given by

[TABLE]

where Faraday’s law has been used to related $\mathbf{B}_{1}$ to $\mathbf{E}_{1}$. In the above spectral expansion, $\theta_{\mathbf{k}}=\mathbf{k}\cdot\mathbf{x}_{0}-\omega_{\mathbf{k}}t_{0}$ is the fast varying phase, $\mathbf{\mathbfcal{E}}_{1,\mathbf{k}}(\mathbf{x}_{1},t_{1},\dots)$ is the slowly varying amplitude, and the summation is over a discrete spectrum $\mathbb{K}_{1}$. Since the electric field is real-valued, whenever $\mathbf{k}\in\mathbb{K}_{1}$, we must also have $-\mathbf{k}\in\mathbb{K}_{1}$. Moreover, we need the reality conditions $\omega_{-\mathbf{k}}=-\omega_{\mathbf{k}}$ and $\mathbf{\mathbfcal{E}}_{-\mathbf{k}}=\mathbf{\mathbfcal{E}}_{\mathbf{k}}^{*}$, where the star denotes complex conjugation. It is easy to check that once these conditions are satisfied, $\mathbf{B}_{1}$ is also real-valued.

The three fluid variables can also be expressed in terms of the electric field. The pressure equation is easy to solve, which gives a simple linear relation

[TABLE]

where the constant $\varepsilon_{s}:=\xi_{s}p_{s0}/n_{s0}$ has the unit of energy. Assuming ideal gas law $p_{0}=n_{0}k_{B}T_{0}$, then $\varepsilon=\xi k_{B}T_{0}$ is proportional to the temperature. Substituting the above relation into the momentum equation, $\mathbf{v}_{s1}$ and $n_{s1}$ can be solved in conjunction with the continuity equation:

[TABLE]

Here, the solution is expressed in terms the warm forcing operator $\hat{\mathbb{F}}_{s,\mathbf{k}}:\mathbb{C}^{3}\rightarrow\mathbb{C}^{3}$, which is a linear operator and is specific to each species and wave vector.

To convert the above symbolic expressions to actual solutions, we need to find an explicit expression for the warm forcing operator. Using the momentum equation, the warm forcing operator satisfies

[TABLE]

for any $\mathbf{Z}\in\mathbb{C}^{3}$. In the above equation, $\beta_{s,\mathbf{k}}=\Omega_{s}/\omega_{\mathbf{k}}$ is the magnetization ratio, where $\Omega_{s}=e_{s}B_{0}/m_{s}$ is the gyrofrequency; $\mathbf{b}$ is the unit vector along $\mathbf{B}_{0}$; and $u_{s}^{2}:=\varepsilon_{s}/m_{s}=\xi_{s}k_{B}T_{s0}/m_{s}$ is the thermal speed. It is easy to see that the inverse operator satisfies $\hat{\mathbb{F}}^{-1}_{ij}=\delta_{ij}-i\beta\epsilon_{ijl}b_{l}-u^{2}k_{i}k_{j}/\omega^{2}$, where $\epsilon_{ijl}$ is the Levi-Civita symbol and I have abbreviated all subscripts for simplicity. Inverting $\hat{\mathbb{F}}^{-1}$, the forcing operator can be expressed as the composite:

[TABLE]

where $\mathbb{F}$ is the cold forcing operator and $\mathbb{P}$ is the pressure operator. The cold forcing operator acts on any complex vector by Shi et al. (2017b)

[TABLE]

where $\gamma^{2}=1/(1-\beta^{2})$ is the magnetization factor. It is easy to check that $\mathbb{F}$ recovers the identity operator in the unmagnetized limit. The pressure operator acts on any complex vector by

[TABLE]

where $\hat{\gamma}^{2}=1/(1-\hat{\beta}^{2})$ is the thermal factor and $\hat{\beta}^{2}=u^{2}(\mathbf{k}\cdot\mathbb{F}\mathbf{k})/\omega^{2}$ is the thermal ratio. It is easy to check that $\mathbb{P}$ recovers the identity operator in the cold limit. Using $\mathbf{k}\cdot\mathbb{F}\mathbf{k}=\gamma^{2}[\mathbf{k}^{2}-\beta^{2}(\mathbf{k}\cdot\mathbf{b})^{2}]$, it is a straightforward calculation to verify that $\hat{\mathbb{F}}$ given by the above formulas satisfies Eq. (27). The warm forcing operator inherits a number of properties from $\mathbb{F}$ and $\mathbb{P}$. First, since $\mathbb{F}^{\dagger}=\mathbb{F}$ is self-adjoint with respect to vector inner products, the warm forcing operator $\hat{\mathbb{F}}^{\dagger}=\hat{\mathbb{F}}$ is also self-adjoint, although $\mathbb{P}^{\dagger}\neq\mathbb{P}$ is not. Second, since $\mathbb{F}_{-\mathbf{k}}=\mathbb{F}^{*}_{\mathbf{k}}$, the pressure operator $\mathbb{P}_{-\mathbf{k}}=\mathbb{P}^{*}_{\mathbf{k}}$ and the warm forcing operator $\hat{\mathbb{F}}_{-\mathbf{k}}=\hat{\mathbb{F}}^{*}_{\mathbf{k}}$. It is then easy to see that the fluid variables given by Eqs.(24)-(26) are real-valued.

Having expressed all fluctuations in terms of $\mathbfcal{E}_{1,\mathbf{k}}$, the only remaining equation is the Ampère’s law, which can be use to constrain the electric field. Substituting Eqs. (23) and (25) into Eq. (16), each Fourier component satisfies the matrix equation $\mathbb{D}_{\mathbf{k}}\mathbfcal{E}_{1,\mathbf{k}}=\mathbf{0}$, where the dispersion tensor

[TABLE]

Here, $\omega_{ps}^{2}=e_{s}^{2}n_{s0}/\epsilon_{0}m_{s}$ is the plasma frequency of species $s$. From the above first-order electric-field equation, it is easy to see that the forcing operator $\hat{\mathbb{F}}_{s,\mathbf{k}}$ is related to the linear susceptibility by

[TABLE]

Although the susceptibility is commonly used in linear wave theories, the forcing operator is more convenient when discussing nonlinear wave-wave couplings.

III.2 Dispersion relations

The first-order electric-field equation has nonzero solutions if and only if the dispersion tensor is degenerate. The degeneracy condition gives the dispersion relation $\det\mathbb{D}_{\mathbf{k}}=\mathbf{0}$, which constrains the wave frequency $\omega_{\mathbf{k}}$ as a function of the wavevector. For each wavevector $\mathbf{k}$, there can be multiple solutions of $\omega_{\mathbf{k}}$, each living on a separate dispersion branch.

When evaluating determinant of the dispersion tensor, it is convenient to use its matrix representations. A particularly convenient coordinate is the field coordinate $(\hat{\mathbf{x}},\hat{\mathbf{y}},\hat{\mathbf{z}})$, in which $\mathbf{B}_{0}=(0,0,B_{0})$ is aligned with $\hat{\mathbf{z}}$ and $\mathbf{k}=k(\sin\theta,0,\cos\theta)$ is in the $x$-$z$ plane. In this coordinate, the warm forcing operator $\hat{\mathbb{F}}$ is represented by the Hermitian matrix

[TABLE]

where the reduced thermal factor $\rho^{2}=\hat{\gamma}^{2}u^{2}k^{2}/\omega^{2}$, and I have abbreviated $\mathrm{s}_{\theta}:=\sin\theta$ and $\mathrm{c}_{\theta}:=\cos\theta$. Summing over responses of all plasma species, the dispersion tensor is represented by the matrix

[TABLE]

where $n=ck/\omega$ is the refractive index, and components of the dielectric tensor are related to

[TABLE]

The above expressions recover the standard Stix symbols in the cold limit, where $\rho^{2}$ becomes zero and the $k$-dependence of the dielectric tensor vanishes.

Taking determinant of the dispersion matrix, the wave dispersion relation can be written in the form

[TABLE]

Coefficients in the above equation depend on $n^{2}$ as well as $\omega$ due to thermal effects:

[TABLE]

Here, $S^{\prime}=S-T\mathrm{c}^{2}_{\theta}$, $D^{\prime}=D+E\mathrm{c}^{2}_{\theta}$, and $P^{\prime}=P-T\mathrm{s}^{2}_{\theta}$. Analogous to the cold case, $R^{\prime}=S^{\prime}+D^{\prime}$ and $L^{\prime}=S^{\prime}-D^{\prime}$. If we formally treat Eq. (46) as a quadratic equation for $n^{2}$, the determinant $F^{2}=B^{2}-4AC=(R^{\prime}L^{\prime}-S^{\prime}P^{\prime})^{2}\mathrm{s}^{4}_{\theta}+4(D^{\prime}P^{\prime}-EA)^{2}\mathrm{c}^{2}_{\theta}\geq 0$. The formal solutions $n^{2}=(B\pm F)/2A$ then give two implicit equations for $n^{2}$ as a function of the wave frequency. In the limit $ck\rightarrow 0$, thermal effects vanish, and the asymptotic dispersion relation is discussed in Appendix A.

A numerically robust procedure for evaluating the dispersion relation is to converted it to a polynomial equation for $\omega^{2}$, using which wave frequencies can be solved as functions of the wavevector (Fig. 1a). To see what multiplicative prefactor is needed, notice that the rational functions $A$, $B$, and $C$ have a pole at $\omega^{2}=0$. In addition, each warm species contribute two poles at $\omega^{4}-(\Omega_{s}^{2}+u_{s}^{2}k^{2})\omega^{2}+\Omega_{s}^{2}u_{s}^{2}k^{2}\mathrm{c}^{2}_{\theta}=0$. One of these poles becomes degenerate with the $\omega^{2}=0$ pole either when the species is cold, in which case the other pole becomes the magnetic pole $\omega^{2}-\Omega_{s}^{2}=0$, or when the species is unmagnetized, in which case the other pole becomes the thermal pole $\omega^{2}-u_{s}^{2}k^{2}=0$. For parallel wave propagation, $\mathrm{c}^{2}_{\theta}=1$, so the magnetic and the thermal poles decouple; for perpendicular wave propagation, $\mathrm{c}^{2}_{\theta}=0$, so one pole becomes $\omega^{2}=0$ while the other pole becomes the hybrid pole $\omega^{2}-\Omega_{s}^{2}-u_{s}^{2}k^{2}=0$; the two poles are otherwise mixed at general angles of propagation. After multiplying the minimal pole-removing prefactor on both sides of Eq. (46), it becomes a polynomial equation for $\omega^{2}$ of degree $N$. For an unmagnetized plasma $N=3+N_{c}$, where $N_{c}=N_{t}+\text{sgn}(N_{s}-N_{t})-1$ is the number of sound waves. Here, $N_{s}$ is the total number of plasma species, $N_{t}$ is the number of warm species, and sgn is the sign function. When the plasma becomes magnetized, $N=3+N_{c}+N_{s}$, because each species contributes an additional cyclotron resonance. At general propagation angles, the dispersion relation is constituted of $N$ separate branches with hybrid characteristics.

III.3 Polarization of eigenmodes

Once the dispersion relation is satisfied, the first-order electric-field equation has nontrivial solutions. The solution space is a one-dimensional vector space when $\omega_{\mathbf{k}}$ is nondegenerate. In this case, the vector space is $\mathbfcal{E}_{1,\mathbf{k}}=\mathcal{E}_{1,\mathbf{k}}\mathbf{e}_{\mathbf{k}}$, where $\mathcal{E}_{1,\mathbf{k}}\in\mathbb{C}$ is an arbitrary complex scalar and the unit polarization vector $\mathbf{e}^{\dagger}_{\mathbf{k}}\mathbf{e}_{\mathbf{k}}=1$ is completely specified, up to the $U(1)$ symmetry, by two polarization angles on the unit sphere.

It is physically meaningful to specify the two polarization angles in relation to the wavevector and the magnetic field. When $\mathbf{k}$ and $\mathbf{B}_{0}$ are not aligned, the unit vector can be decomposed as $\mathbf{e}=\hat{\mathbf{k}}e_{k}-i\hat{\mathbf{y}}e_{y}+(\hat{\mathbf{k}}\times\hat{\mathbf{y}})e_{\times}$, where the unit vector $\hat{\mathbf{y}}\mathbin{\!/\mkern-5.0mu/\!}\mathbf{B}_{0}\times\mathbf{k}$. In spherical coordinate, components of $\mathbf{e}$ can be written as $e_{k}=\cos\phi$, $e_{y}=\sin\phi\sin\psi$, and $e_{\times}=\sin\phi\cos\psi$. The wave is longitudinal when $\phi=0^{\circ}$ and transverse when $\phi=90^{\circ}$; the wave electric field is in the $\mathbf{k}$-$\mathbf{B}_{0}$ plane when $\psi=0^{\circ}$ and at maximum angle with the plane when $\psi=90^{\circ}$. Since arbitrary scaling is allowed, the unit vector $\mathbf{e}$ is defined on the projective space. In terms of $\phi$ and $\psi$, the wave polarization is invariant under transformations $\Psi_{\pm}:(\phi,\psi)\rightarrow(-\phi,\psi\pm 180^{\circ})$ and $\Phi_{\pm}:(\phi,\psi)\rightarrow(\phi\pm 180^{\circ},\psi)$. since the polarization angles are periodic in $360^{\circ}$, the wave polarization is invariant under actions of $\Psi_{\pm}$ and $\Phi_{\pm}$ in arbitrary compositions.

To compute polarization angles for each eigenmode, it is more convenient to use the wave coordinate, which is related to the field coordinate by $(\hat{\mathbf{k}},-i\hat{\mathbf{y}},\hat{\mathbf{k}}\times\hat{\mathbf{y}})=(\hat{\mathbf{x}},\hat{\mathbf{y}},\hat{\mathbf{z}})\mathcal{L}_{y}(\theta)$, where $\mathcal{L}_{y}(\theta)$ is a left-handed rotation around $\hat{\mathbf{y}}\rightarrow-i\hat{\mathbf{y}}$ by angle $\theta$. In the wave coordinate, the dispersion tensor is represented by a different matrix $\mathcal{D}^{\prime}=\mathcal{L}^{-1}\mathcal{D}\mathcal{L}$, which can be written explicitly as

[TABLE]

The degenerate matrix equation $\mathcal{D}^{\prime}\mathbf{e}=0$ is solved when the polarization angles satisfy

[TABLE]

These polarization angles can be computed after solving $\omega$ as a function of $ck$, using which the refractive index and the dispersion symbols can be evaluated.

A numerically robust procedure for computing the unit polarization vector is by directly solving the degenerate matrix equation. Denoting $\mathbf{d}^{\prime}_{i}$ the $i$-the row vector of the matrix $\mathcal{D}^{\prime}$, then $\mathbf{e}\propto\alpha_{1}\mathbf{d}^{\prime}_{2}\times\mathbf{d}^{\prime}_{3}+\alpha_{2}\mathbf{d}^{\prime}_{3}\times\mathbf{d}^{\prime}_{1}+\alpha_{3}\mathbf{d}^{\prime}_{1}\times\mathbf{d}^{\prime}_{2}$, where $\alpha_{i}$ is an arbitrary parameter. Since $\mathcal{D}^{\prime}$ is a rank-2 matrix when the dispersion relation is satisfied, the three vectors on the RHS are parallel, and at most two of them can be simultaneously zero when the plasma is magnetized. By summing up the three terms, it is guaranteed that the RHS is never a zero vector as long as special values of $\alpha_{i}$ are avoided. Then, the unit polarization vector can be determined after normalization. Having obtained $\mathbf{e}$ in Cartesian coordinate, the polarization angles in spherical coordinate can be easily determined. An example is shown in Fig. 1b, where $\phi$ (solid lines) and $\psi$ (dashed lines) are plotted as functions of wave frequency.

While the above procedure is generally applicable, it is instructive to note two special propagation angles. When $\theta=0^{\circ}$ or $90^{\circ}$, Eq. (54) and (55) become indeterminate, even though the polarization vector remains well defined. When $\theta=0^{\circ}$, the longitudinal electrostatic modes, which satisfy $P=0$, decouple with the transverse modes. One set of transverse modes satisfy $n^{2}=R$ and are right-handed (R) circularly polarized with $\tan\psi=1$; the other set of transverse modes satisfy $n^{2}=L$ and are left-handed (L) circularly polarized with $\tan\psi=-1$. When $\theta=90^{\circ}$, the ordinary (O) wave decouples. The O wave is unmagnetized EM wave, which satisfies $n^{2}=P$ with the wave electric field along $\mathbf{B}_{0}$. The remaining modes are the extraordinary (X) wave hybridized with plasma waves, which satisfies $n^{2}=RL/S+T$ with $\tan\phi=S/D$ neither transverse nor longitudinal. For these modes $\cos\psi=0$, and the wave electric field is always perpendicular to the background magnetic field.

III.4 Energy of linear waves

Inherited from the nonlinear fluid-Maxwell equations, the linear system also conserves energy locally. The energy density of the linear system is of $\lambda^{2}$ order:

[TABLE]

where the last term comes from $p_{s2}$ as will become clear later when I discuss $\lambda^{2}$-order equations. The energy flux of the linear system is

[TABLE]

where the last term is due to thermal flux and kinetic flux does not contribute at $\lambda^{2}$ order. Using the first-order equations (12)-(16), it is a straightforward calculation to show the energy conservation law on fast scales:

[TABLE]

Notice that the conservation law mixes all linear waves that are present in the system. In other words, not only do linear waves contribute individually, but their interference also contributes to the total energy.

Now let us focus on the energy of a single linear wave. Using the Fourier expansion $\mathbf{E}_{1}=(\mathbfcal{E}e^{i\theta}+\mathbfcal{E}^{*}e^{-i\theta})/2$, $\mathbf{E}_{1}^{2}=(\mathbfcal{E}^{2}e^{2i\theta}+2\mathbfcal{E}\mathbfcal{E}^{*}+\mathbfcal{E}^{*2}e^{-2i\theta})/4$. We see the electric energy has a rapidly-oscillating part and a slowly-varying part. Only the later remains after averaging on $t_{0}$ and $\mathbf{x}_{0}$ scales, namely, $\langle\mathbf{E}_{1}^{2}\rangle_{0}=\mathbfcal{E}\mathbfcal{E}^{*}/2$. Following similar arguments, the averaged magnetic energy can be computed using Fourier expansion Eq. (23). Then, the averaged field energy density is $\epsilon_{0}\mathcal{E}_{i}^{\dagger}(2\delta^{ij}-\sum_{s}\omega_{ps}^{2}\hat{\mathbb{F}}_{s}^{ij}/\omega^{2})\mathcal{E}_{j}/4$, where I have used $\mathbb{D}_{\mathbf{k}}\mathbfcal{E}_{1,\mathbf{k}}=\mathbf{0}$ and Eq. (31). Similarly, the averaged kinetic energy density is $\epsilon_{0}\sum_{s}\omega_{ps}^{2}\mathcal{E}_{i}^{\dagger}(\hat{\mathbb{F}}^{2}_{s})^{ij}\mathcal{E}_{j}/4\omega^{2}$, where I have used the self-adjoint property of $\hat{\mathbb{F}}_{s}$. Finally, the averaged thermal energy density is $\epsilon_{0}\sum_{s}\omega_{ps}^{2}u_{s}^{2}\mathcal{E}_{i}^{\dagger}\hat{\mathbb{F}}_{s}^{ia}k_{a}k_{b}\hat{\mathbb{F}}_{s}^{bj}\mathcal{E}_{j}/4\omega^{4}$. Summing up the three terms, the total energy density of the linear wave can be written as

[TABLE]

The Hamiltonian $\mathbb{H}$ of the linear wave is related to the dielectric tensor $\epsilon=\mathbb{I}+\sum_{s}\chi_{s}$ by the usual relation $\omega\mathbb{H}=\partial(\omega^{2}\epsilon)/\partial\omega$, where the partial derivative is at fixed $\mathbf{k}$. Since the susceptibility $\chi_{s}$ is related to $\hat{\mathbb{F}}_{s}$ by Eq. (32), the wave energy operator

[TABLE]

where $\mathbb{I}$ is the idensity operator and the partial derivative is again at fixed $\mathbf{k}$. To see the connection between Eqs. (59) and (60), note the following identity

[TABLE]

An easy way to show the above identity is to take patial derivative on both sides of $\hat{\mathbb{F}}\hat{\mathbb{F}}^{-1}=\mathbb{I}$, then $\partial\hat{\mathbb{F}}/\partial\omega=-\hat{\mathbb{F}}(\partial\hat{\mathbb{F}}^{-1}/\partial\omega)\hat{\mathbb{F}}$. Using the expression of the inverse operator, it is easy to see $\omega(\partial\hat{\mathbb{F}}^{-1}/\partial\omega)\mathbf{Z}=i\beta\mathbf{Z}\times\mathbf{b}+2u^{2}\mathbf{k}(\mathbf{k}\cdot\mathbf{Z})/\omega^{2}$. Replacing $\mathbf{Z}$ by $\hat{\mathbb{F}}\mathbf{Z}$, and using the Eq. (27), the above identity is then obvious.

The averaged wave energy depends on the wave envelope. To leading order, the wave envelope is a function of $\mathbf{x}_{1}-\mathbf{v}_{g}t_{1}$, where $v_{g}^{i}=\partial\omega/\partial k_{i}$ is the group velocity of the linear wave. Consequently, the averaged energy satisfies the advection equation on slow scales:

[TABLE]

This equation is consistent with the $\lambda^{3}$-order conservation law if and only if

[TABLE]

Now let me show this is indeed the case by direct calculations. The averaged Poynting flux is $\epsilon_{0}c^{2}\mathcal{E}_{a}^{\dagger}(2k^{i}\delta^{ab}-k^{a}\delta^{ib}-k^{b}\delta^{ia})\mathcal{E}_{b}/4\omega$, and the averaged thermal flux is $\epsilon_{0}\sum_{s}\omega_{ps}^{2}\mathcal{E}_{a}^{\dagger}(\partial\hat{\mathbb{F}}_{s}^{ab}/\partial k_{i})\mathcal{E}_{b}/4\omega$. For the thermal flux, I have used the identity that the partial derivative of $\hat{\mathbb{F}}$ at fixed $\omega$ is given by

[TABLE]

which can be shown similarly to Eq. (61). Summing the Poynting and the thermal fluxes, Eq. (63) is satisfied whenever $\omega v_{g}^{i}\mathcal{E}_{a}^{*}\mathbb{H}^{ab}\mathcal{E}_{b}=\mathcal{E}_{a}^{*}[c^{2}(2k^{i}\delta^{ab}-k^{a}\delta^{ib}-k^{b}\delta^{ia})+\sum_{s}\omega_{ps}^{2}\partial\hat{\mathbb{F}}_{s}^{ab}/\partial k_{i}]\mathcal{E}_{b}.$ This equation is nothing other than $\mathcal{E}_{a}^{*}d_{k_{i}}\mathbb{D}^{ab}\mathcal{E}_{b}=0$, which is trivially satisfied as a consequence of $\mathbb{D}_{\mathbf{k}}\mathbfcal{E}_{1,\mathbf{k}}=\mathbf{0}$ for all linear eigenmodes, whose frequency satisfies the dispersion relation and polarization solves the first-order electric-field equation. I have thus verified that the envelope of a single linear wave advects at the wave group velocity as expected.

IV Magnetized three-wave interactions

Building upon a thorough understanding of linear waves, we are now ready to study their interactions. Due to these interactions, waves become coupled. Consequently, instead of passing through each other uneventfully with only linear superpositions, waves now actually “collide” and exchange energy and momentum. In this section, I will investigate couplings mediated by three-wave interactions. These lowest-order interactions are usually the strongest whenever resonance conditions can be satisfied.

IV.1 Second-order equations

The $\lambda^{2}$-order fluid-Maxwell’s equations (17)-(21) are linear partial differential equations for $\mathbf{E}_{2}$, $\mathbf{B}_{2}$, $p_{s2}$, $n_{s2}$, and $\mathbf{v}_{s2}$ with source terms. The general solution is again a superposition of plane waves, whose spectrum is completely determined by existing linear waves in the system.

Let us express all second-order fluctuations in terms of electric-field fluctuations:

[TABLE]

where the $\lambda^{2}$-order spectrum $\mathbb{K}_{2}$ and amplitudes $\mathbfcal{E}_{2,\mathbf{k}}$ will be determined later. Using the second-order Faraday’s law [Eq. (20)], the $\lambda^{2}$-order magnetic field is

[TABLE]

The second line involves slow derivatives of linear wave amplitudes, which are unknown at this point.

To solve the fluid equations, let us first express $p_{s2}$ in terms of $n_{s2}$. Using $p_{s1}=\varepsilon_{s}n_{s1}$ [Eq. (24)], many terms in the second-order pressure equation cancels. Integrating Eq. (19) on $t_{0}$ time scale, the second-order pressure is

[TABLE]

Notice that due to quadratic nonlinearities, the average $\langle p_{s2}\rangle_{0}$ is in general nonzero. We see that $p_{s2}$ is related to $n_{s1}^{2}$, as anticipated from the energy of linear waves.

Next, let us express $n_{s2}$ in terms of $\mathbf{v}_{s2}$. Suppose the Fourier expansion of the second-order velocity is $\mathbf{v}_{s2}=\sum_{\mathbf{k}}\exp(i\theta_{\mathbf{k}})\mathbfcal{V}_{s2,\mathbf{k}}/2$. Then, substituting $\mathbf{v}_{s1}$ and $n_{s1}$ [Eqs. (25) and (26)] into the second-order continuity equation [Eq. (17)], the $\lambda^{2}$-order density is

[TABLE]

where summation on the first line is over the spectrum of $\mathbf{v}_{s2}$. Since thermal effect does not enter through the continuity equation directly, the above expression is identical to the cold-fluid case, where the current beating is

[TABLE]

The current beating comes from the nonlinearty $n_{s1}\mathbf{v}_{s1}$, and thermal effects enter only indirectly through the warm forcing operator.

Eliminating $p_{s2}$ and $n_{s2}$, we can now solve for $\mathbf{v}_{s2}$. Using Eq. (13) to simplify Eq. (18), the equation is of the form $\sum_{\mathbf{k}}(\mathbfcal{V}_{2,\mathbf{k}}-i\beta\mathbfcal{V}_{2,\mathbf{k}}\times\mathbf{b}-u^{2}\mathbf{k}\mathbf{k}\cdot\mathbfcal{V}_{2,\mathbf{k}}/\omega_{\mathbf{k}}^{2})\exp(i\theta_{\mathbf{k}})=\sum_{\mathbf{k}}\mathbf{Z}_{\mathbf{k}}\exp(i\theta_{\mathbf{k}})$, where I have suppressed the species index for simplicity. Then, using Eq. (27) of the warm forcing operator, the solution is

[TABLE]

The spectrum $\mathbfcal{V}_{2,\mathbf{k}}$ can now be read out from the above equation, and explicit expressions of $n_{s2}$ and $p_{s2}$ can then be obtained. On the third line of Eq. (IV.1), the first term $\mathbf{L}_{\mathbf{p},\mathbf{q}}$ is the longitudinal beating introduced by the $\mathbf{v}_{1}\times\mathbf{B}_{1}$ nonlinearity:

[TABLE]

The second term $\mathbf{T}_{\mathbf{p},\mathbf{q}}$ is the turbulent beating introduced by the $\mathbf{v}_{1}\cdot\nabla_{0}\mathbf{v}_{1}$ nonlinearity:

[TABLE]

These two terms are the same as in the cold case, except that the cold forcing operator is now replaced by the warm forcing operator. Additionally, the thermal effect introduces thermal beating as a third term:

[TABLE]

which is caused by nonlinearities in $n_{s2}$ and $p_{s2}$. The turbulent and thermal beatings can be rewritten in terms of the velocity perturbation $\mathbfcal{V}_{1,\mathbf{k}}$. These beatings are purely fluid effects, which exist even when the fluid is neutral (Appendix B). On the other hand, the longitudinal and current beatings are genuine electromagnetic nonlinearities, whereby transverse EM waves in the vacuum become mixed with the otherwise longitudinal motion of the plasma.

Having expressed all fluctuations in terms of electric-field fluctuations, we can now solve for the electric field. Substituting in $\mathbf{B}_{2}$ and $\mathbf{v}_{s2}$ into Eq. (21), the $\lambda^{2}$-order electric-field equation can be grouped into four sets of terms, involving $\mathbfcal{E}_{2,\mathbf{k}}$, $\partial_{t_{1}}\mathbfcal{E}_{1,\mathbf{k}}$, $\nabla_{1}\mathbfcal{E}_{1,\mathbf{k}}$, and $\mathbfcal{E}_{1,\mathbf{p}}\mathbfcal{E}_{1,\mathbf{q}}$. Differentiating on $t_{0}$ scale, the first set simplifies to $\mathbb{D}_{\mathbf{k}}\mathbfcal{E}_{2,\mathbf{k}}$ using Eq. (31). The second set simplifies to $(\partial\mathbb{D}_{\mathbf{k}}/\partial\omega_{\mathbf{k}})\partial_{t_{1}}\mathbfcal{E}_{1,\mathbf{k}}$, using $\mathbb{D}_{\mathbf{k}}\mathbfcal{E}_{1,\mathbf{k}}=\mathbf{0}$ and Eq. (61). The third set simplifies to $-(\partial\mathbb{D}_{\mathbf{k}}/\partial\mathbf{k})\cdot\nabla_{1}\mathbfcal{E}_{1,\mathbf{k}}$, using Eq. (64). Finally, the fourth set encapsulates all beatings. With all these simplifications, the second-order electric-field equation is then

[TABLE]

where the scattering strength summed over all plasma species is

[TABLE]

The above equation is formally identical to the cold case, except that the dispersion tensor $\mathbb{D}$ now contains thermal modifications. In addition, thermal effects directly enter the quadratic response of each species:

[TABLE]

The thermal beating $\mathbf{U}^{s}_{\mathbf{p},\mathbf{q}}$ vanishes when the thermal speed $u_{s}^{2}\rightarrow 0$, while the three other beatings remain finite when the species becomes cold.

Since the second-order electric-field equation must be satisfied for each Fourier component, the equation can be split into two sets. The first set involves only $\lambda$-order spectrum $\mathbb{K}_{1}$. Suppose within the spectral bandwidths, $\mathbf{k},\mathbf{p},\mathbf{q}\in\mathbb{K}_{1}$ satisfy the three-wave resonance conditions $\mathbf{k}=\mathbf{p}+\mathbf{q}$ and $\omega_{\mathbf{k}}=\omega_{\mathbf{p}}+\omega_{\mathbf{q}}$, then the on-shell equation is of the form

[TABLE]

The RHS is simply zero if $\mathbf{k}\in\mathbb{K}_{1}$ is not in resonance with other waves. The second set of equations generate $\mathbb{K}_{2}$ from $\mathbb{K}_{1}$. To aviod secular behavior, we can demand that $\mathbb{K}_{2}\cap\mathbb{K}_{1}$ is the empty set, so that all $\mathbf{k}\in\mathbb{K}_{2}$ is generated by off-resonance beatings. Then, each off-shell equation is of the form

[TABLE]

where $\mathbf{k}=\mathbf{p}+\mathbf{q}$ and $\omega_{\mathbf{k}}=\omega_{\mathbf{p}}+\omega_{\mathbf{q}}$ are still satisfied, except that $\mathbf{k}\notin\mathbb{K}_{1}$ is no longer a linear eigenmode. In other words, the wave dispersion relation is not satisfied for off-shell waves, so the operator $\mathbb{D}_{\mathbf{k}}$ can be inverted to give $\mathbfcal{E}_{2,\mathbf{k}}=i\mathbb{D}_{\mathbf{k}}^{-1}\mathbf{S}_{\mathbf{p},\mathbf{q}}$. Thus, the $\lambda^{2}$-order spectrum is completely determined by the $\lambda$-order spectrum. After solving both the on-shell and the off-shell equations, $\mathbf{B}_{2},\mathbf{v}_{s2},n_{s2},p_{s2}$ can be determined and the second-order equations are then solved.

Now that the $\lambda^{2}$-order equations have been formally solved, it is important to note $\mathbf{S}_{\mathbf{p},\mathbf{q}}$ satisfies a number of identities, which are required in order for the solutions to be valid. First, from its definition, it is obvious that $\mathbf{S}_{\mathbf{p},\mathbf{q}}=\mathbf{S}_{\mathbf{q},\mathbf{p}}$ is symmetric, which is expected from the symmetry of three-wave interactions. Second, it is a straightforward calculation to check that $\mathbf{S}_{\mathbf{p},\mathbf{q}}^{*}=-\mathbf{S}_{-\mathbf{q},-\mathbf{p}}$ satisfies the reality condition. Consequently, all second-order fluctuations are real-valued. Finally, $\mathbf{S}_{\mathbf{p},-\mathbf{p}}=\mathbf{0}$ is secular-free. In other words, a wave does not beat with itself to generate a mode with both $\omega=0$ and $\mathbf{k}=\mathbf{0}$. To show the last identity, notice that for each species $\mathbf{U}_{\mathbf{p},-\mathbf{p}}+\mathbf{U}_{-\mathbf{p},\mathbf{p}}=\mathbf{0}$. Moreover, since $\omega_{-\mathbf{p}}=-\omega_{\mathbf{p}}$, we have $\mathbf{R}_{\mathbf{p},-\mathbf{p}}+\mathbf{R}_{-\mathbf{p},\mathbf{p}}=\hat{\mathbb{F}}_{\mathbf{0}}(\mathbf{L}_{\mathbf{p},-\mathbf{p}}+\mathbf{L}_{-\mathbf{p},\mathbf{p}}+\mathbf{T}_{\mathbf{p},-\mathbf{p}}+\mathbf{T}_{-\mathbf{p},\mathbf{p}})=-i\beta_{\mathbf{p}}\hat{\mathbb{F}}_{\mathbf{0}}[(\mathbf{p}\cdot\hat{\mathbb{F}}_{\mathbf{p}}^{*}\mathbfcal{E}_{\mathbf{p}}^{*})(\hat{\mathbb{F}}_{\mathbf{p}}\mathbfcal{E}_{\mathbf{p}}\times\mathbf{b})+\text{c.c.}]/\omega_{\mathbf{p}}^{2}$. To see what $\hat{\mathbb{F}}_{\mathbf{0}}$ is, we can take limits $\omega\rightarrow 0$ and $\mathbf{k}\rightarrow\mathbf{0}$. Although $\hat{\mathbb{F}}_{\mathbf{0}}$ depends on the how these two limits are taken, $\mathbf{R}_{\mathbf{p},-\mathbf{p}}+\mathbf{R}_{-\mathbf{p},\mathbf{p}}=\mathbf{0}$ is independent of the limiting procedure. Therefore, $\mathbf{S}_{\mathbf{p},-\mathbf{p}}=\mathbf{0}$ is always satisfied.

IV.2 On-shell equations and action conservation

While the off-shell equations are easy to solve, the on-shell equations are where the nontrivial dynamics is contained. These equations are nonlinearly-coupled advection equations. Due to the vector nature of these equations, not only does wave amplitude change, but the wave polarization can also evolve. Moreover, the wave phase, trajectory, and angular momentum are usually altered as well due to three-wave interactions.

Since the wave polarization can change in general, it is important to verify whether the on-shell equation is compatible with the first-order equation. Suppose $\mathbb{D}\mathbfcal{E}=\mathbf{0}$ is satisfied over the entire space before the waves encounter, then $\mathbb{D}\mathbfcal{E}=\mathbf{0}$ will be satisfied for all time if $\mathbb{D}\partial_{t_{1}}\mathbfcal{E}=\mathbf{0}$ for all $\mathbf{x}_{1}$. Using the on-shell equation (77), this compatibility condition is satisfied if

[TABLE]

where I have used the fact that the Hamiltonian $\omega\mathbb{H}=\partial\mathbb{D}/\partial\omega$ is invertible. Since the dispersion operator $\mathbb{D}$ is degenerate for linear eigenmodes, the above condition only requires that $(\partial\mathbb{D}/\partial\mathbf{k})\cdot\nabla_{1}\mathbfcal{E}+\mathbf{S}$ is in the null space of $\mathbb{D}(\partial\mathbb{D}/\partial\omega)^{-1}$. Notice that the rank of $\mathbb{D}$ is at most two, so the above condition imposes at most two constraints. Therefore, there always remains degree of freedom allowing $\mathbfcal{E}$ to evolve in time.

The compatibility condition can be used to remove the redundant degree of freedom of the on-shell equation. Taking total $\mathbf{k}$ derivative on both sides of $\mathbb{D}\mathbfcal{E}=\mathbf{0}$, where $\omega$ and $\mathbfcal{E}$ are now regarded as functions of $\mathbf{k}$, we have $(d\mathbb{D}/d\mathbf{k})\mathbfcal{E}+\mathbb{D}d\mathbfcal{E}/d\mathbf{k}=\mathbf{0}$. Here, the total derivative of the dispersion tensor is $d\mathbb{D}/d\mathbf{k}=\mathbf{v}_{g}\partial\mathbb{D}/\partial\omega+\partial\mathbb{D}/\partial\mathbf{k}$, where $\mathbf{v}_{g}=\partial\omega/\partial\mathbf{k}$ is the wave group velocity. Using the wave energy operator [Eq. (60)], the on-shell equation (77) becomes $\omega\mathbb{H}(\partial_{t_{1}}+\mathbf{v}_{g}\cdot\nabla_{1})\mathbfcal{E}+\mathbb{D}\partial_{l}d\mathbfcal{E}/dk_{l}=\mathbf{S}$. The compatibility condition Eq. (79) then guarantees that the advection keeps $\mathbfcal{E}$ inside the eigenspace. To be more specific, denoting $\Pi$ the projection operator into the null space of $\mathbb{D}$, then after applying the compatibility condition, the on-shell equation is reduced to

[TABLE]

where $\mathbf{S}^{\pi}=\mathbb{H}\,\Pi\,\mathbb{H}^{-1}\mathbf{S}$ is the eigen projection of $\mathbf{S}$. This somewhat abstract notation is illustrated using cold unmagnetized plasma as an example in Appendix C. Since $\mathbfcal{E}^{\dagger}\mathbb{D}=\mathbf{0}$ for eigenmodes, $\mathbfcal{E}^{\dagger}\mathbf{S}^{\pi}=\mathbfcal{E}^{\dagger}\mathbf{S}$. We see only the eigen projection of the scattering strength affects the evolution of the linear wave.

Now let us focus on the simplest nontrivial case, namely, the resonant interaction between three on-shell waves. Without loss of generality, the resonance conditions can be written as

[TABLE]

where $\omega_{i}$ are positive. For simplicity, I will abbreviate $\mathbfcal{E}_{j}:=\mathbfcal{E}_{1,\mathbf{k}_{j}}$. Moreover, since the slow dynamics is on $t_{1}$ and $\mathbf{x}_{1}$ scales only, I will suppress the index of the temporal and spatial scales. Then, the three on-shell equations can be written as

[TABLE]

where $\bar{j}:=-j$, and $d:=\partial_{t}+\mathbf{v}_{g}\cdot\nabla$ is the convective derivative at respective wave group velocities. Notice that the group velocity $\mathbf{v}_{g}$ is in general not aligned with $\mathbf{k}$ when the plasma is magnetized.

What is nontrivial about the on-shell equations is that they guarantee action conservation for resonant three-wave interactions. The conservation laws are consequences of the action identity

[TABLE]

Using this identity, $\mathbfcal{E}^{\dagger}\mathbf{S}^{\pi}=\mathbfcal{E}^{\dagger}\mathbf{S}$, and $\mathbf{S}_{-\mathbf{p},-\mathbf{q}}=-\mathbf{S}^{*}_{\mathbf{p},\mathbf{q}}$, it is easy to show the action conservation laws

[TABLE]

where $\langle U_{j}\rangle$ is the energy of wave $j$ averaged on fast scales [Eq. (59)]. The action conservation laws are manifestations of the Feynman rules of three-wave interactions Shi et al. (2016): each quanta of wave “1” is converted to a quanta of wave “2” and a quanta of wave “3”, or vice versa. Using the action conservation laws and $\omega_{1}=\omega_{2}+\omega_{3}$, the total wave energy is also conserved:

[TABLE]

Notice that the above conservation laws hold only when the three waves are in resonance.

The action identity can be shown by direct calculations. During the calculation, one will encounter terms like $\mathbfcal{F}_{1}^{*}\cdot\mathbfcal{F}_{2}$, where $\mathbfcal{F}_{j}=\hat{\mathbb{F}}_{j}\mathbfcal{E}_{j}$. Such terms can be simplified using the following quadratic identity of the Forcing operator:

[TABLE]

which can be obtained from Eq. (28) using property of the cold forcing operator Shi et al. (2017b): $(\beta_{1}-\beta_{2})\mathbb{F}_{1}\mathbb{F}_{2}=\beta_{1}\mathbb{F}_{1}-\beta_{2}\mathbb{F}_{2}$. A suite of similar identities can be obtained using $\hat{\mathbb{F}}_{\bar{j}}=\hat{\mathbb{F}}_{j}^{*}$ and $\hat{\mathbb{F}}^{\dagger}=\hat{\mathbb{F}}$. The product is then

[TABLE]

where I have used property of the pressure operator: $\mathbb{P}\mathbf{Z}=\mathbf{Z}+u^{2}(\mathbf{k}\cdot\hat{\mathbb{F}}\mathbf{Z})\mathbf{k}/\omega^{2}$. The action identity [Eq. (86)] can then be shown by straightforward calculation of $\mathbfcal{E}_{1}\cdot(\mathbf{R}_{\bar{2},\bar{3}}+\mathbf{R}_{\bar{3},\bar{2}})/\omega_{1}^{2}$, and comparing it with the other two terms of the same structure.

In fact, terms in the action identity can be organized into a very simple and intuitive form:

[TABLE]

Here, $\mathcal{E}_{i}$ is the scalar amplitude such that $\mathbfcal{E}_{i}=\mathcal{E}_{i}\mathbf{e}_{i}$, where $\mathbf{e}_{i}$ is the unit polarization vector. In the above expression, $\Theta^{s}$ and $\Phi^{s}$ are the nondimensionalized electromagnetic and the thermal scattering strengths. The electromagnetic scattering $\Theta^{s}$ is due to the $P_{s}^{i}(\partial_{i}A_{l})J^{l}_{s}$ coupling in the Lagrangian Shi et al. (2017b), where $\mathbf{P}_{s}$ is the displacement and $\mathbf{J}_{s}$ is the current of species $s$ in response to perturbations of the gauge field $\mathbf{A}$. Same as the cold case, the electromagnetic scattering contains six permutations:

[TABLE]

Each scattering channel, which satisfies $\Theta^{s}_{i,\bar{j}\bar{l}}=\Theta^{s*}_{\bar{i},jl}$ is given by the simple formula

[TABLE]

where $\mathbf{f}_{s,j}:=\hat{\mathbb{F}}_{s,j}\mathbf{e}_{j}$. The thermal scattering $\Phi^{s}$ is due to warm-fluid effects, which is present even in neutral fluids (Appendix B). The thermal scattering contains four contributions

[TABLE]

The symmetric thermal scattering is formed by contracting $\mathbf{f}$ with its own wave vector:

[TABLE]

On the other hand, the skewed symmetric thermal scattering, which satisfies $\Phi^{s}_{\bar{j}}=-\Phi^{s}_{j}$, is formed by contracting $\mathbf{f}$ with a common wave vector:

[TABLE]

Since $\Phi^{s}$ is proportional to $u_{s}^{2}/c^{2}$, it is usually very small, otherwise a relativistic plasma model is required in the first place. It is obvious that when the species is cold, the thermal scattering vanishes and the three-wave scattering reduces to purely electromagnetic scattering in the cold-fluid case.

IV.3 Coupling coefficient and growth rate

When polarization is not of concern, the on-shell equations can be reduced to scalar-amplitude equations called the three-wave equations, which contain a single essential parameter: the coupling coefficient. Denoting $\mathbfcal{E}=\mathcal{E}\mathbf{e}$, we can define the wave energy coefficient

[TABLE]

such that the averaged wave energy $\langle U\rangle=\epsilon_{0}u|\mathcal{E}|^{2}/2$. It is then convenient to normalize the scalar amplitude by

[TABLE]

where $e$ and $m_{e}$ are the charge and mass of electrons. With this normalization, the wave energy $\langle U\rangle\propto\omega^{2}|a|^{2}$. Notice that the decomposition $\mathbfcal{E}=\mathcal{E}\mathbf{e}$ is not unique, and we can always perform $U(1)$ rotations $\mathcal{E}\rightarrow\mathcal{E}e^{i\alpha}$ and $\mathbf{e}\rightarrow\mathbf{e}e^{-i\alpha}$ such that the vector amplitude $\mathbfcal{E}$ is invariant. To remove this arbitrariness, we can impose the condition that $\mathcal{E}\in\mathbb{R}$ is real-valued. Then, the decomposition is unique up to the $\mathbb{Z}_{2}$ symmetry $\mathcal{E}\rightarrow-\mathcal{E}$ and $\mathbf{e}\rightarrow-\mathbf{e}$. With this reduced symmetry, the normalized scalar amplitude is also real-valued.

The equation for the normalized scalar amplitude can be derived from the on-shell equations. Allowing the polarization to evolve, the derivative $d(\mathcal{E}u^{1/2})=[\mathbf{e}^{\dagger}\mathbb{H}d(\mathcal{E}\mathbf{e})+\text{c.c.}]/4u^{1/2}$. Then, the real-valued amplitude evolves according to $da_{1}=e(\mathbf{e}_{1}^{\dagger}\mathbf{S}_{2,3}+\text{c.c.})/(4m_{e}c\omega_{1}^{2}u_{1}^{1/2})$. Using Eq. (75) for $\mathbf{S}_{2,3}$ and Eq. (90) for the inner product, the normalized scalar amplitudes satisfy the following three-wave equations:

[TABLE]

where the convective derivatives are at respective group velocities of the three waves. Due to the $\mathbb{Z}_{2}$ symmetry, only the relative signs of the above equations are of importance. The essential parameter of the three-wave equations is the coupling coefficient $\Gamma_{r}$, which is the real part of the complex-valued coupling coefficient

[TABLE]

Here, $Z_{s}:=e_{s}/e$ and $M_{s}:=m_{s}/m_{e}$ are the normalized charge and mass of species $s$. As a consequence of wave interference, contributions of different plasma species add up in the complex plane. Moreover, the three waves also interfere. The three-wave interference depends on the relative wave phase, whose change corresponds to a rotation of $\Gamma$ in the complex plane. When the three waves are phase locked, $\Gamma_{r}=|\Gamma|$ is maximized. In this case, the beat wave of $a_{1}$ and $a_{2}$ are in phase with $a_{3}$, so the plasma responses constructively interfere. The above formula is formally identical to the cold case, except for the extra $\Phi$ term due to thermal scattering.

The three-wave coupling may be small for three distinct reasons Shi et al. (2018). First, the coupling coefficient $\Gamma$ may be interference-suppressed because terms in the summation cancel one another. In this case, although scattering due to each species is appreciable, the nonlinear responses are of opposite phases and destructively interfere. Second, $\Gamma$ may be polarization-suppressed because the numerator of each terms is small. In this case, the wave vectors and polarization vectors are at orthogonal angles, so that vector inner products in $\Theta^{s}$ and $\Phi^{s}$ are small. Finally, $\Gamma$ may be energy-suppressed because its denominator is large. In this case, a large fraction of the wave energy is kinetic or thermal, so the electromagnetic field amplitude is small for given wave energy.

A consequence of the three-wave interactions is the parametric decay instability. During the instability, a large amplitude pump wave $a_{1}$ decays to the frequency-downshifted daughter waves $a_{2}$ and $a_{3}$, whose relative phases are automatically locked. In the linear regime of parametric interaction, $a_{1}$ barely changes and $a_{2}$ and $a_{3}$ grow almost exponentially with growth rate

[TABLE]

The actual growth rate observed in experiments are likely influenced by wave damping, which includes both collisional and collisionless damping. Damping effects can be important and may be introduced phenomenologically in the three-wave equations. However, the ideal-fluid model does not captured damping self-consistently.

The growth rate may be compared to that of Raman backscattering $\gamma_{R}=\sqrt{\omega_{1}\omega_{p}}|a_{1}|/2$ in cold unmagnetized plasmas of the same density, where $\omega_{p}^{2}=\sum_{s}\omega_{ps}^{2}$ is the total plasma frequency. We can write the growth rate $\gamma_{0}=\gamma_{R}\mathcal{M}$, then

[TABLE]

The normalized growth rate is now symmetric with respect to the three waves, and is proportional to the coupling coefficient up to some kinematic factors. Since $\Gamma\sim\omega_{p}^{2}$, the normalized growth rate is zero in the limit $\omega_{p}\rightarrow 0$. This is expected because there is no three-wave coupling in the vacuum.

V Examples

The above general theory is applicable to a discrete spectrum of weakly-damped and weakly-coupled waves in magnetized warm-fluid plasmas. The waves can propagate in any directions and have arbitrary frequencies. An example of resonant interaction is shown in Fig. 2, where $a_{1}$ is on the L branch, $a_{2}$ is on the U branch, and $a_{3}$ is on the F branch. The matching of resonance conditions in the Fourier space is shown in the inset for collimated scattering. Also shown in the inset is the interaction geometry in the configuration space, where $\langle\hat{\mathbf{k}}_{1},\mathbf{B}_{0}\rangle$ is fixed at $30^{\circ}$, while $\hat{\mathbf{k}}_{2}$ has polar angle $\theta_{2}$ and azimuthal angle $\phi_{2}$. Due to mirror symmetry, the frequency downshift $\Delta\omega=\omega_{2}-\omega_{1}$ is plotted only for the western hemisphere (Fig. 2a), while the normalized growth rate $\mathcal{M}$ is shown only for the eastern hemisphere (Fig. 2b). Plasma parameters used in this example are the same as in Fig. 1. The pump wave (marked by green dot) has frequency $\omega_{1}=75$ Trad/s, which corresponds to $ck_{1}\approx 51.54$ Trad/s. Due to the presence of the magnetic field, the coupling has intricate angular dependence.

A numerically robust procedure for evaluating the coupling coefficient and the growth rate is as follows. First, imagine we launch a pump wave with frequency $\omega_{1}$ in direction $\hat{\mathbf{k}}_{1}$ on a given branch. Then, the wave number $k_{1}$ can be solved from the dispersion relation using procedures described in Sec. III.2, and the unit polarization vector $\mathbf{e}_{1}$ can be computed using procedures described in Sec. III.3. The matrix representation of the warm forcing operator $\hat{\mathbb{F}}_{1}$ can be computed, based on which the wave energy operator $\mathbb{H}_{1}$ and the wave energy coefficient $u_{1}$ can be evaluated. Second, suppose we place a detector along $\hat{\mathbf{k}}_{2}$, we can in principle detect all waves whose wave number is such that $\omega_{1}-\omega_{2}(\mathbf{k}_{2})=\omega_{3}(\mathbf{k}_{1}-\mathbf{k}_{2})$ is in resonance with a third wave. For each pair of wave branches, solving the above resonance condition gives $k_{2}$, from which $\omega_{2}$ and $\mathbf{e}_{2}$ can be determined and $\hat{\mathbb{F}}_{2}$ and $u_{2}$ can be evaluated. At the same time, $\mathbf{k}_{3}$ also becomes known, from which $\omega_{3}$, $\mathbf{e}_{3}$, $\hat{\mathbb{F}}_{3}$ and $u_{3}$ can be evaluated. Finally, once these quantities are determined, the scattering strengths can be evaluated, so are the coupling coefficient [Eq. (101)] and the growth rate [Eq. (103)]. This procedure is used to obtain Fig. 2, and will be further demonstrated below using two examples.

V.1 Scattering of high-frequency lasers

When the pump $a_{1}$ and the seed $a_{2}$ are high-frequency EM waves, whose frequencies $\omega_{1},\omega_{2}\gg\omega_{ps},|\Omega_{s}|$, they asymptote to vacuum EM waves. In the high-frequency limit, $\beta,\hat{\beta}\rightarrow 0$, $\gamma,\hat{\gamma}\rightarrow 1$, and the forcing operators $\hat{\mathbb{F}}_{1},\hat{\mathbb{F}}_{2}\sim\mathbb{I}$. Since $a_{3}$ is a plasma wave with much lower frequency, the electromagnetic scattering is dominantly $\Theta^{s}\simeq-(c\mathbf{k}_{3}\cdot\mathbf{f}^{*}_{s,3})(\mathbf{e}_{1}\cdot\mathbf{e}_{2}^{*})/\omega_{3}$ whenever $\mathbf{e}_{1}\cdot\mathbf{e}_{2}^{*}$ is of order unity. Moreover, since vacuum EM waves are transverse, the thermal scattering becomes $\Phi^{s}\simeq u_{s}^{2}k_{3}^{2}(c\mathbf{k}_{3}\cdot\mathbf{f}^{*}_{s,3})(\hat{\mathbf{k}}_{3}\cdot\mathbf{e}_{1})(\hat{\mathbf{k}}_{3}\cdot\mathbf{e}_{2}^{*})/\omega_{1}\omega_{2}\omega_{3}$. When the angle between $\mathbf{k}_{1}$ and $\mathbf{k}_{2}$ is not too small, $ck_{3}\sim\omega_{1},\omega_{2}$ is always large, so $\mathbf{e}_{3}\sim\hat{\mathbf{k}}_{3}$ is approximately longitudinal. Then, $\mathbb{P}_{3}\mathbf{e}_{3}\simeq\hat{\gamma}^{2}_{3}\mathbf{e}_{3}$ and $\mathbf{e}_{3}\cdot\mathbf{f}_{3}^{*}\simeq\gamma^{2}_{3}\hat{\gamma}^{2}_{3}(1-\beta_{3}^{2}\cos^{2}\theta_{3})$, where $\theta_{3}=\langle\mathbf{k}_{3},\mathbf{B}_{0}\rangle$ and the thermal ratio $\hat{\beta}_{3}^{2}\simeq u^{2}k_{3}^{2}\gamma_{3}^{2}(1-\beta_{3}^{2}\cos^{2}\theta_{3})/\omega_{3}^{2}$. Finally, the wave energy coefficients $u_{1},u_{2}\simeq 1$, and $u_{3}$ can be evaluated with $\mathbf{f}_{3}\!\cdot\!\mathbf{f}_{3}^{*}=\hat{\gamma}_{3}^{4}[\cos^{2}\theta_{3}+\gamma_{3}^{4}(1+\beta_{3}^{2})\sin^{2}\theta_{3}]$. Using these asymptotics, the coupling coefficient and the parametric growth rate can be approximated.

The above approximations clearly recover the cold magnetized case Shi et al. (2017b), and they also recover the well-known Raman and Brillouin scatterings in warm unmagnetized plasmas. In the unmagnetized limit, $\beta=0$, $\gamma=1$, and $\hat{\beta}_{3}^{2}=u^{2}k_{3}^{2}/\omega_{3}^{2}$. Moreover, the plasma waves are purely longitudinal with the dispersion relation $\omega_{3}^{2}=\sum_{s}\omega_{ps}^{2}\hat{\gamma}_{s,3}^{2}$. Then, $\mathbf{e}_{3}\cdot\mathbf{f}_{3}^{*}=\hat{\gamma}^{2}_{3}$, $\mathbf{f}_{3}\!\cdot\!\mathbf{f}_{3}^{*}=\hat{\gamma}_{3}^{4}$, and $u_{3}=\sum_{s}\omega_{ps}^{2}\hat{\gamma}_{s,3}^{4}/\omega_{3}^{2}$. In most cases $\Phi^{s}\ll\Theta^{s}$, because $\omega_{3}\ll\omega_{1},\omega_{2}\sim ck_{3}$. For the same reason, $k_{3}\simeq 2k_{1}\sin(\alpha/2)$, where $\alpha$ is the angle between $\mathbf{k}_{1}$ and $\mathbf{k}_{2}$. Now let us focus on quasi-neutral electron-ion plasma with $Z_{i}=1$. In this two-species plasma, there are two longitudinal waves. The high-frequency wave is the Langmuir wave, whose mediation gives rise to the Raman scattering. Since $\lambda_{D}k_{3}\ll 1$ is required for weak collisionless damping, the dispersion relation is $\omega_{3}^{2}\simeq\omega_{p}^{2}$, so $\hat{\beta}_{i,3}^{2},\hat{\beta}_{e,3}^{2}\ll 1$, and $\hat{\gamma}_{i,3}^{2}\simeq\hat{\gamma}_{e,3}^{2}\simeq 1$. The normalized unmagnetized Raman growth rate is then

[TABLE]

We see responses by the two species destructively interfere, and exact cancellation occurs in electron-positron plasma where $M_{i}=1$. On the contrary, the responses constructively interfere for the low-frequency sound wave, whose mediation gives rise to the Brillouin scattering. For sound wave, the dispersion relation is $\omega_{3}^{2}\simeq c_{s}^{2}k_{3}^{2}$, where the sound speed $c_{s}^{2}=2M_{i}u_{i}^{2}/(M_{i}+1)$ assuming $u_{e}^{2}=M_{i}u_{i}^{2}$. Then, $\hat{\gamma}_{e}^{2}\simeq 2/(1-M_{i})$ and $\hat{\gamma}_{i}^{2}\simeq 2M_{i}/(M_{i}-1)$ are of opposite signs. The wave energy coefficient $u\simeq 4\omega_{p}^{2}M_{i}/\omega^{2}(M_{i}-1)^{2}$, and the unmagnetized Brillouin growth rate is then

[TABLE]

for both electron-ion and electron-positron plasmas when temperature is not too high. The above recovers the weak-coupling results in the literature Sjölund and Stenflo (1967b); Gorbunov (1969); Litvak and Trakhtengerts (1971); Forslund et al. (1973); Drake et al. (1974); Forslund et al. (1975); Edwards et al. (2016), which were derived for unmagnetized plasmas in the parametric interaction picture.

Without any approximation, the growth rate can be evaluated numerically to determine collective laser scattering in magnetized plasmas. Consider an example relevant to laser-driven magnetized liner fusion Davies et al. (2017); Barnak et al. (2017). In the experimental design, a preheat laser with 351-nm wavelength propagates along $\mathbf{B}_{0}$ of about 30 T. The $\text{D}_{2}$ plasma, for which $Z_{i}=1$ and $M_{i}\approx 3671$, has density $\sim$1.5 mg/$\text{cm}^{3}$. After fully ionized, the number density is about $n_{e}=n_{i}=4.5\times 10^{20}\,\text{cm}^{-3}$, and the plasma temperature is about $T_{e}=400$ eV and $T_{i}=150$ eV. In this example, the Debye length $\sim 10^{-2}\,\mu$m is much smaller than the laser wavelength $\lambda_{0}$, and $\lambda_{0}$ is much smaller than the collisional mean free path $\sim 10\,\mu$m, so the ideal fluid model is applicable. Moreover, since $\lambda_{0}$ is much smaller than the ion gyro radius, ions are essentially unmagnetized and the A branch has minuscule contribution. Therefore, effects of magnetization are mainly due to electrons.

Due to cylindrical symmetry, the scattering only depends on the polar angle $\theta_{2}$, which is $0^{\circ}$ for forward scattering and $180^{\circ}$ for backward scattering. The growth rates, in units of Raman backscattering, are shown in Fig. 3, where the curves are color-coded by frequency downshifts. The growth rates are polarization-dependent, and the eigenmodes are elliptically polarized, except when $\theta=90^{\circ}$ where they become the linearly polarized X and O wave. For the R-wave pump, scattering to the R branch (a, c) is polarization-suppressed for near backward scattering, while scattering to the L branch (b, d) is polarization-suppressed for near forward scattering. When $B_{0}=30$ T (a, b), $\Omega_{e}\approx 5.3$ Trad/s is much smaller than $\omega_{p}\approx 1.2\times 10^{3}$ Trad/s, so scattering from the U branch is close to Raman. Similarly, since $v_{A}/c\approx 7\times 10^{-5}$ is smaller than $c_{s}/c\approx 9\times 10^{-4}$, the sound wave is little modified, and scattering from the S branch is close to Brillouin. Other than modifying Raman and Brillouin, the magnetic field introduces additional modes from which the laser can scatter. However, in weak magnetic fields, scattering from the F branch is energy-suppressed, because the F branch is dominated by electron cyclotron motion. In larger magnetic fields, for example $B_{0}=300$ T (c, d), effects of magnetization then become larger.

V.2 Scattering of MHD waves

To illustrate that the general formula is applicable to any wave triplets, let us consider scattering of MHD waves as another example. In this case, the wave frequency $\omega\ll\Omega_{i}$, and the fluid-Maxwell’s equations asymptote to MHD equations. Consequently, the wave dispersion relation also asymptotes to that of the MHD waves (Appendix A). The asymptotics are particularly simple for wave propagation parallel to $\mathbf{B}_{0}$, where the forcing operator becomes $\hat{\mathbb{F}}_{\parallel}\mathbf{Z}=\gamma^{2}(\mathbf{Z}+i\beta\mathbf{Z}\times\mathbf{b})+(\hat{\gamma}^{2}-\gamma^{2})(\mathbf{Z}\cdot\mathbf{b})\mathbf{b}$. It may be tempting to already take the $\omega/\Omega\rightarrow 0$ limit for $\hat{\mathbb{F}}_{\parallel}$. However, the limit should be taken only after summations over species is carried out, because leading terms may turn out to cancel.

First, let us determine the approximate wave dispersion relations. Consider two species plasmas with $Z_{i}=1$, then the sum in $\mathbb{D}_{11}$ is $\omega_{pe}^{2}\gamma_{e}^{2}+\omega_{pi}^{2}\gamma_{i}^{2}\simeq-c^{2}\omega^{2}/v_{A}^{2}$, where $v_{A}^{2}=c^{2}M_{i}\Omega_{i}^{2}/\omega_{p}^{2}$ is the Alfvén speed. The sum in $\mathbb{D}_{12}$ is $\beta_{e}\omega_{pe}^{2}\gamma_{e}^{2}+\beta_{i}\omega_{pi}^{2}\gamma_{i}^{2}\simeq-(1-1/M_{i})c^{2}\omega^{3}/v_{A}^{2}\Omega_{i}$. Finally, the sum in $\mathbb{D}_{33}$ is $\omega_{pe}^{2}\hat{\gamma}_{e}^{2}+\omega_{pi}^{2}\hat{\gamma}_{i}^{2}\simeq\omega^{2}\omega_{p}^{2}(\omega^{2}-c_{s}^{2}k^{2})/(\omega^{2}-u_{e}^{2}k^{2})(\omega^{2}-u_{i}^{2}k^{2})$. The dispersion tensor $\mathbb{D}$ for parallel wave propagation can then be easily determined in the field coordinate. The longitudinal wave satisfies $\omega^{2}\simeq c_{s}^{2}k^{2}$, which is essentially the unmagnetized sound wave. The transverse waves satisfy $(1+v_{A}^{2}/c^{2})\omega^{2}=v_{A}^{2}k^{2}\pm(1-1/M_{i})\omega^{3}/\Omega_{i}$. The “+” branch has higher phase velocity and is right-handed circularly polarized with $\mathbf{e}\propto(1,i,0)$; the “–” branch has lower phase velocity and is left-handed circularly polarized with $\mathbf{e}\propto(1,-i,0)$. To lowest order in $\omega/\Omega_{i}$, the two branches merge into dispersionless Alfvén waves $\omega^{2}=c_{A}^{2}k^{2}$, where $c_{A}^{2}=v_{A}^{2}/(1+v_{A}^{2}/c^{2})$. The energy coefficient of the Alfvén waves is $u\simeq c^{2}/c_{A}^{2}$, which is usually very large because most wave energy is contained in magnetic and fluid fluctuations instead of the wave electric field. In the MHD limit $v_{A}^{2}\ll c^{2}$, the above results recover the dispersion relations of parallel-propagating MHD waves.

Now we can compute three-wave coupling between MHD waves. For parallel wave propagation, other than the coupling between three sound waves, which is discussed in Appendix B, the only nonzero coupling is between two Alfvén waves of the same polarization ($a_{1},a_{2}$) and a sound wave ($a_{3}$). Since the waves are dispersionless, the resonance conditions can be satisfied only when $a_{1}$ and $a_{2}$ are counter propagating. In this geometry, the resonant wave vectors are $k_{2}/k_{1}=|c_{s}-c_{A}|/(c_{s}+c_{A})$ and $k_{3}/k_{1}=2c_{A}/(c_{s}+c_{A})$. To compute the the scattering between these waves, notice that for the sound wave $\mathbf{f}=\hat{\gamma}^{2}\mathbf{b}$, while for the Alfvén wave $\mathbf{f}=\mathbf{e}/(1\pm\beta)$ where $\pm$ corresponds to the L and R polarizations. The electromagnetic scattering due to each species is thereof $\Theta_{s}\simeq-[1/(\omega_{1}+\Omega_{s})+1/(\omega_{2}+\Omega_{s})]c\omega_{1}\omega_{2}\omega_{3}/c_{A}(\omega_{3}^{2}-u_{s}^{2}k_{3}^{2})$, and the thermal scattering is $\Phi_{s}\simeq 0$ because two waves are transverse. Summing over species and then take the limit $\omega/\Omega_{i}\rightarrow 0$, the coupling coefficient is

[TABLE]

The coupling can also be expressed in terms of the magnetic field $a_{1}=M_{i}\Omega_{i}\mathcal{B}_{1}/\omega_{1}B_{0}$ whereby the growth rate [Eq. (102)] can be readily evaluated, which agrees with the weak-coupling result in the literature Jayanti and Hollweg (1993a); Matsukiyo and Hada (2003); Wong and Goldstein (1986); Derby Jr (1978); Goldstein (1978) for both electron-ion and electron-positron plasmas.

Without any approximation, the exact formula of the resonant coupling coefficient can be evaluated numerically. Let us consider an example relevant for solar corona at a height comparable to the solar radius Phillips et al. (1995); Guhathakurta et al. (1996); Aschwanden (2006); Tomczyk et al. (2007). There, the plasma is mostly hydrogen with $M_{i}\approx 1837$. The plasma density $n_{e}\sim n_{i}\sim 10^{7}\,\text{cm}^{-3}$, the plasma temperature $T_{e}\sim T_{i}\sim 100$ eV, and the magnetic field $B_{0}\sim$ 1 G. Correspondingly, $\omega_{p}\approx 1.8\times 10^{8}$ rad/s, $c_{s}/c\approx 8.0\times 10^{-4}$, and $v_{A}/c\approx 2.3\times 10^{-3}$. The ion cyclotron frequency $\Omega_{i}\approx 10^{4}$ rad/s is much higher than the observed Alfvén wave frequency, which is in the mHz band. In this frequency range, the low-frequency waves are well-approximated by ideal MHD waves. Consider the coupling between two Alfvén waves via the sound wave. The daughter wave frequency $f_{2}/f_{1}$ is shown in Fig. 4a, and the coupling coefficient $\Gamma/\Gamma_{\parallel}$ is shown in Fig. 4b, where $\Gamma_{\parallel}$ is given by Eq. (106). While the coupling has strong dependence on $c_{s}/v_{A}$, it has very weak dependence on the frequency of the parallel pump Alfvén wave. Moreover, the dependence on $\theta_{2}$, the angle of the daughter wave with respect to the local magnetic field, is also weak. Consequently, the decay of the parallel pump wave only slight prefers exact backward geometry.

VI Discussion and summary

Beyond linear waves, this paper treats coherent three-wave interactions in magnetized plasma by solving the warm-fluid model to second order. Unlike previous attempts, which were specialized for each wave triad in restricted geometry, here, the systematic treatment using perturbation theory offers a unified description of all possible interactions at arbitrary angles. This methodology, first introduced for magnetized cold-fluid plasma Shi et al. (2017b), is extended to incorporate thermal effects. The ideal-fluid model is applicable when all wavelengths are much larger than the Debye length, while much shorter than the collisional mean free path. In this regime, thermal effects enter indirectly through the forcing operator [Eq. (27)], as well as directly in the quadratic response [Eq. (76)]. Nevertheless, the second-order electric-field equation [Eq. (IV.1)] remains formally unchanged.

The formalism developed in this paper is not only general, but also practical, whereby numerical values of the coupling coefficient can be obtained. The coupling coefficient is an essential parameter in the commonly-used three-wave equations [Eqs. (98)-(100)]. Previously, little is known about the numerical value of the coupling coefficient when the plasma becomes magnetized. Now, a general formula has been provided [Eq. (101)], which can be evaluated for any three resonantly interacting waves. To demonstrate the powerfulness of the general formula, the coupling between high-frequency lasers via Raman [Eq. (104)] and Brillouin [Eq. (105)] scatterings in unmagnetized plasmas are recovered as special cases. Moreover, the same formula also recovers coupling between two Alfvén waves and a sound wave [Eq. (106)], which are at the other extreme of the wave spectrum. While asymptotic expressions of the general formula may be found for special cases, the exact formula can always be evaluated using numerical procedures demonstrated in this paper. Based on nontrivial analytic simplifications, numerical evaluations the of coupling coefficient can now be made efficient and robust.

In summary, this paper derives a general formula governing resonant three-wave interactions in magnetized warm-fluid plasmas in the weak-coupling regime. Applying the formula to magnetized inertial confinement fusion conditions, the magnetic field is found to modify Raman and Brillouin scatterings of lasers, as well as introduce additional scattering modes at oblique angles. For parameters relevant to solar corona, the formula for parallel coupling between two Alfvén waves via the sound wave is found to give good approximations also at oblique angles. Due to weak angular dependence, exact backscattering is only slightly preferred over oblique decays.

Acknowledgements.

The author thanks John D. Moody, Bradley B. Pollock, David J. Strozzi, Pierre A. Michel, Nathaniel J. Fisch, Matthew R. Edwards, and Matthew W. Kunz for helpful discussions. This work was performed under the auspices of the U.S. Department of Energy by Lawrence Livermore National Laboratory under Contract DE-AC52-07NA27344 and was supported by the Lawrence Fellowship through LLNL-LDRD Program under Project No. 19-ERD-038.

Appendix A Asymptotic dispersion $ck\rightarrow 0$

The dispersion relation contains gapped and gapless modes. For gapped modes, the wave frequency $\omega\rightarrow\omega_{c}$ when $ck\rightarrow 0$, where $\omega_{c}$ is some finite cutoff frequency. For gapless modes, $\omega\rightarrow 0$ when $ck\rightarrow 0$, but the refractive index $n=ck/\omega$ approaches some finite constant. The asymptotic dispersion relation is useful for analytic approximations, and may be used as initial guesses for numerical root finding.

For gapped modes, the cutoff frequencies are solutions of $C(\omega_{c})=0$, where $C$ is given by Eq. (49). Since thermal effects vanish, one cutoff frequency is always $\omega_{p}$. The other cutoff frequencies are solutions of $R(\omega_{c})=L(-\omega_{c})=0$. In a magnetized plasma of $N_{s}$ species, there are $N_{s}+1$ non-negative solutions, which becomes strictly positive when the plasma is not quasi-neutral. For finite but small $ck$, we can expand near $\omega_{c}$. The asymptotic dispersion relation is quadratic: $\omega\simeq\omega_{c}^{2}+\delta\omega^{2}$, where $\delta\omega^{2}=2Bc^{2}k^{2}/\omega_{c}\partial_{\omega}C$. Here, $B$ and $\partial_{\omega}C$ are evaluated at $\omega=\omega_{c}$ and $ck=0$. The analytic expression is simple, since thermal effects vanish.

On the other hand, thermal effects are important for gapless modes. To obtain asymptotic dispersion relation when $\omega\rightarrow 0$, we can expand using Laurent series. After tedius but otherwise straightforward expansions, the leading terms in a quasi-neutral plasma are

[TABLE]

assuming $c_{\theta}^{2}\gg\omega^{2}/\Omega^{2}$. The dispersion coefficients $I_{0}=1+\sum_{s}\omega_{ps}^{2}\eta_{s}^{2}/\Omega_{s}^{2}$, $I_{1}=\sum_{s}\omega_{ps}^{2}\eta_{s}^{2}/\Omega_{s}$, and $I_{2}=\sum_{s}\omega_{ps}^{2}\eta_{s}^{2}$, where $\eta_{s}^{2}=1/(1-n^{2}c_{\theta}^{2}u_{s}^{2}/c^{2})$. In the cold limit, $I_{0}\rightarrow c^{2}/c_{A}^{2}$, $I_{1}\rightarrow 0$, and $I_{2}\rightarrow\omega_{p}^{2}$. Substituting Eqs. (107)-(109) into the dispersion relation [Eq. (46)], we obtain an equation for $n^{2}$. The Alfvén wave decouples with the dispersion relation

[TABLE]

What remains are the fast wave mixed with the sound waves, which is given by

[TABLE]

A special case is when all species are cold. Then, the sound wave vanishes and the above recovers the cold fast wave $\omega^{2}=c_{A}^{2}k^{2}$. In more general cases, a numerically robust procedure for solving the dispersion relation is to remove poles of Eq. (111), and convert it to a polynomial of $n^{2}c_{\theta}^{2}$ of degree $N_{c}$+1. When $N_{t}\geq N_{s}-1$, the leading coefficient is $\sum_{s}\omega_{ps}^{2}\prod_{s^{\prime}\neq s}(-\mu_{s^{\prime}}^{2})$, otherwise the leading coefficient is $\prod_{\mu_{s^{\prime}}\neq 0}(-\mu_{s^{\prime}}^{2})\cdot\sum_{\mu_{s^{\prime}}=0}\omega_{ps}^{2}$. The polynomial equation has exactly $N_{c}$+1 real and positive roots, which can be found by standard numerical methods. The above is a multi-fluid extension of MHD, which retains only one sound wave.

Appendix B Three-wave in neutral fluid

To illustrate that that turbulent and thermal beatings are originated from fluid nonlinearities, let us consider three-wave interactions in neutral fluid, which is described by

[TABLE]

where $\rho$ is the mass density and $d_{t}=\partial_{t}+\mathbf{v}\cdot\nabla$ is the convective derivative at the fluid velocity.

The linearized fluid equations describe sound waves. Suppose we weakly perturbe the equilibrium with constant $\rho_{0}$, $p_{0}$, and $\mathbf{v}_{0}=\mathbf{0}$, the first-order fluid velocity is $\mathbf{v}_{1}=\frac{1}{2}\sum_{\mathbf{k}\in\mathbb{K}_{1}}\exp(i\theta_{\mathbf{k}})\mathbfcal{V}_{1,\mathbf{k}}$. The continuity equation then gives $\rho_{1}/\rho_{0}=\frac{1}{2}\sum_{\mathbf{k}\in\mathbb{K}_{1}}\exp(i\theta_{\mathbf{k}})\mathbf{k}\cdot\mathbfcal{V}_{1,\mathbf{k}}/\omega_{\mathbf{k}}$, and the pressure equation gives $p_{1}=u^{2}\rho_{1}$, where $u^{2}=\xi p_{0}/\rho_{0}$ is the thermal speed. Substituting these into the momentum equation, each Fourier amplitude satisfies

[TABLE]

The dispersion operator is now $\bar{\mathbb{D}}_{\mathbf{k}}=\omega_{\mathbf{k}}^{2}-u^{2}\mathbf{k}\mathbf{k}$. The eigenmode satisfies the dispersion relation $\omega^{2}=u^{2}\mathbf{k}^{2}$, and is the longitudinally polarized sound wave.

To second order in multiscale perturbative analysis, the equations can be obtained from Eqs. (17)-(19) by setting the electromagnetic contributions to zero. The pressure equation gives $p_{2}$ by Eq. (67) after replacements $mn\rightarrow\rho$ and $\varepsilon/m\rightarrow u^{2}$. Expanding the second-order velocity as $\mathbf{v}_{2}=\sum_{\mathbf{k}}\exp(i\theta_{\mathbf{k}})\mathbfcal{V}_{2,\mathbf{k}}/2$, the continuity equation gives $\rho_{2}$ by Eq. (68) after replacing $ie\hat{\mathbb{F}}\mathbfcal{E}/m\omega\rightarrow\mathbfcal{V}$. Using $\bar{\mathbb{D}}_{\mathbf{k}}\mathbfcal{V}_{1,\mathbf{k}}=\mathbf{0}$, the second-order momentum equation can then be written as

[TABLE]

which is formally identical to Eq. (IV.1) if the later is written in terms of velocity perturbations. Analogously, $\bar{\mathbf{S}}_{\mathbf{p},\mathbf{q}}=(\bar{\mathbf{R}}_{\mathbf{p},\mathbf{q}}+\bar{\mathbf{R}}_{\mathbf{q},\mathbf{p}})/2$, and the only difference is that now $\bar{\mathbf{R}}_{\mathbf{p},\mathbf{q}}=\bar{\mathbf{T}}_{\mathbf{p},\mathbf{q}}+\bar{\mathbf{U}}_{\mathbf{p},\mathbf{q}}$. Here, $\bar{\mathbf{T}}_{\mathbf{p},\mathbf{q}}$ and $\bar{\mathbf{U}}_{\mathbf{p},\mathbf{q}}$ can be obtained from Eqs. (72) and (73) by replacing $\hat{\mathbb{F}}\mathbfcal{E}/\omega\rightarrow\mathbfcal{V}$. We see that turbulent and thermal beatings are intrinsically fluid nonlinearities.

The second-order velocity equation can be split into off-shell and on-shell equations. The off-shell equations can be solved by inverting the nondegenerate $\bar{\mathbb{D}}$, and the on-shell equations can be simplified using eigen projections. Suppose the resonance conditions are of the form “$p=q+l$”, the three-wave amplitude equations are

[TABLE]

where $v_{l}=v_{-l}^{*}$ is the complex amplitude such that $\mathbfcal{V}_{1,\mathbf{k}_{l}}=v_{l}\hat{\mathbf{k}}_{l}$. In the parametric decay picture, the growth rate $\gamma_{0}=(1+\xi)(k_{2}k_{3})^{1/2}|v_{1}|/4$. To obtain the above three-wave equations, I have used the fact that resonance conditions can be satisfied only when $\hat{\mathbf{k}}_{1}=\hat{\mathbf{k}}_{2}=\hat{\mathbf{k}}_{3}$, because the sound waves are dispersionless. Due to the special dispersion relation, three-wave interactions are one dimensional, along which any two copropagating waves can resonantly interact.

Appendix C Eigen projection

To illustrate how the compatibility condition [Eq. (79)] can be used in conjunction with the on-shell equation [Eq. (77)], let us consider unmagnetized cold plasma as an example. In this case, $\hat{\mathbb{F}}=\mathbb{I}$ is the identity operator and the dispersion tensor is

[TABLE]

The partial derivatives are $\partial\mathbb{D}/\partial\omega=2\omega\mathbb{I}$ and $\partial\mathbb{D}^{ij}/\partial k_{l}=c^{2}(k_{i}\delta_{jl}+k_{j}\delta_{il}-2k_{l}\delta_{ij})$. The on-shell equation is then of the form

[TABLE]

where I have omitted the subscripts of $t_{1}$ and $\mathbf{x}_{1}$. Notice that the Eq. (119) has redundant degrees of freedom, because the spatial derivatives originate from the projection operator $\mathbb{I}-\hat{\mathbf{k}}\hat{\mathbf{k}}$, which has a nontrivial kernel.

For electromagnetic waves, the dispersion relation is $\omega^{2}=\omega_{p}^{2}+c^{2}k^{2}$. The dispersion tensor then becomes $\mathbb{D}=c^{2}\mathbf{k}\mathbf{k}$, which is a rank-1 operator. The null space is two dimensional, and the eigenmodes are transverse, which satisfy $\mathbf{k}\cdot\mathbfcal{E}=0$. The compatibility condition is satisfied if and only if

[TABLE]

Substituting the solution of $\nabla\cdot\mathbfcal{E}$ into Eq. (119), the on-shell equation becomes

[TABLE]

where $c^{2}\mathbf{k}/\omega$ is nothing other than the group velocity, and $\mathbf{S}^{\perp}=(\mathbb{I}-\hat{\mathbf{k}}\hat{\mathbf{k}})\mathbf{S}$ is the transverse projection. While the transverse projection is typically put in “by-hand” when studying unmagnetized three-wave interactions Michel et al. (2014), here I have shown why the projection necessarily arises.

For the cold Langmuir waves, $\omega^{2}=\omega_{p}^{2}$, and $\mathbb{D}=c^{2}(\mathbf{k}\mathbf{k}-k^{2}\mathbb{I})$ is a rank-2 operator. The null space is thereof one dimensional, and the eigenmode is longitudinal, which satisfies $\mathbfcal{E}\mathbin{\!/\mkern-5.0mu/\!}\mathbf{k}$. The compatibility condition is satisfied if and only if

[TABLE]

Substituting this into Eq. (119), which can be separated into parallel and perpendicular components, the on-shell equation becomes

[TABLE]

As expected, the group velocity of the cold Langmuir wave is zero, and only the longitudinal component $\mathbf{S}^{\parallel}=\hat{\mathbf{k}}(\hat{\mathbf{k}}\cdot\mathbf{S})$ affects the wave evolution.