Abstract

This paper studies electromagnetic-wave propagation in anisotropic dipole-active crystals under nonstationary stimulated Raman scattering (SRS) by polaritons when interacting waves have different polarizations. We derive a coupled system of nonlinear differential equations for the laser pump, Stokes-wave components, polariton-wave components, and dipole-active phonon vibrations. The model combines shortened Maxwell’s equations for waves of different polarizations and frequencies with a forced-oscillator equation for the dipole-active phonons. We then verify the resulting system by deriving the corresponding Manley-Rowe relation, which provides an energy-consistency criterion for the proposed nonstationary SRS model. The main contribution of this paper is the derivation of an energy-consistent nonstationary SRS model for anisotropic dipole-active crystals in which the interacting electromagnetic waves may have different polarizations. The Manley-Rowe relation is used not as an auxiliary result but as a consistency test showing that the coupled pump, Stokes, polariton, and phonon equations preserve the expected balance of energy exchange.

Keywords

Nonstationary SRS, Polaritons, Phonons, Polarization, Manley-Rowe

Share and Cite:

1. Introduction

In recent years, nonlinear optics has advanced from studies of bulk crystals to micro- and nanoscale photonic systems [1]-[16]. Work on soliton frequency-comb generation in X-cut ferroelectric materials, thin-film lithium niobate, and lithium tantalate microresonators shows that Raman effects, photorefractive effects, anisotropy, and nonlinear coupling are central to modern integrated photonic platforms [5]-[9]. These developments motivate further research on nonstationary stimulated Raman scattering by polaritons in nonlinear anisotropic crystals [17]-[22]. Building on our earlier studies of nonstationary SRS by polaritons for a stationary pump field [23] and soliton solutions for the pump, Stokes, and polariton waves [24], this paper derives a more general nonlinear system for anisotropic dipole-active crystals. We also show that the system satisfies the Manley-Rowe relation, which confirms its energy-consistent structure [25]-[30]. The paper is organized as follows. Section 2 states the wave, phonon, and polarization assumptions and derives the coupled nonlinear equations for the pump, Stokes, polariton, and phonon amplitudes. Section 3 introduces dimensionless amplitudes and uses them to derive the Manley-Rowe relation for waves with different polarizations. Section 4 summarizes the result.

2. Basic Principles and Equations

We assume that the pump field is a linearly polarized plane electromagnetic wave and that the nonlinear medium is a crystal layer bounded by the planes z = 0 and z = L. The crystal is taken to have cubic symmetry. As shown in [22], a pumped cubic crystal becomes anisotropic because the linearly polarized pump wave deforms the dielectric constant. Under these assumptions, the pump wave is written as

${\overrightarrow{E}}_l\left(\overrightarrow{r},t\right)={\widehat{e}}_l{A}_l\left(z,t\right)\exp \left[i\left({\overrightarrow{k}}_l\cdot \overrightarrow{r}-{\omega }_{l}t\right)\right]+c.c.$ (1)

The subscripts l, s, p, and f denote the pump (laser), Stokes, polariton wave fields, and phonon vibration; $\omega$ are the frequencies, n and $\overrightarrow{k}$ are the refractive indices and the wave vectors in the unpumped medium, and $\widehat{e}$ are the real unit vectors of electromagnetic fields. The medium is nonmagnetic and transparent at frequencies ${\omega }_{l,s}$ . We use the Stokes and polariton fields as

$\begin{array}{l}{\overrightarrow{E}}_s\left(\overrightarrow{r},t\right)={\displaystyle \sum _{\mu =1,2}{\widehat{e}}_s^{\left(\mu \right)}{A}_s^{\left(\mu \right)}\left(z,t\right)\exp \left[i\left({\overrightarrow{k}}_s\cdot \overrightarrow{r}-{\omega }_{s}t\right)\right]}+c.c.,\\ {\overrightarrow{E}}_p\left(\overrightarrow{r},t\right)={\displaystyle \sum _{\sigma =1,2,3}{\widehat{e}}_p^{\left(\sigma \right)}{A}_p^{\left(\sigma \right)}\left(z,t\right)\exp \left[i\left(\overrightarrow{W}\cdot \overrightarrow{r}-{\omega }_{p}t\right)\right]}+c.c.,\end{array}$ (2)

where: ${\widehat{e}}_s^{\left(\mu \right)}\perp {\overrightarrow{k}}_s$ , ${\widehat{e}}_s^{\left(1\right)}\perp {\widehat{e}}_s^{\left(2\right)}$ , $k_s=q_s{n}_s$ , $q_s={\omega }_s/c$ , $\overrightarrow{W}={\overrightarrow{k}}_l-{\overrightarrow{k}}_s$ , ${\widehat{e}}_p^{\left(1,2\right)}\perp \overrightarrow{W}$ , ${\widehat{e}}_p^{\left(1\right)}\perp {\widehat{e}}_p^{\left(2\right)}$ , ${\widehat{e}}_p^{\left(3\right)}=\frac{\overrightarrow{W}}{W}$ , ${\omega }_p={\omega }_l-{\omega }_s$ .

The longitudinal component of the Stokes wave can be neglected, but this cannot be done for the polariton wave in the phonon region. As will be seen from formulas that follow, with further advance into this region all three polariton amplitudes become comparable. The phase shift of the polariton wave is determined by the vector $\overrightarrow{W}$ and not by ${\overrightarrow{k}}_p$ ($k_p=q_p\sqrt{{\epsilon }_p}$ , $q_p={\omega }_p/c$ , ${\epsilon }_p={{\epsilon }^{\prime }}_p+i{{\epsilon }^{″}}_p$ is the dielectric constant at the frequency ${\omega }_p$ ).

The expression for the phonon vibration is:

${\overrightarrow{Q}}_f\left(\overrightarrow{r},t\right)={\displaystyle \sum _{\sigma =1,2}{\widehat{e}}_f^{\left(\sigma \right)}{u}_f^{\left(\sigma \right)}\left(z,t\right)\exp \left[i\left(\overrightarrow{W}\cdot \overrightarrow{r}-{\omega }_{f}t\right)\right]}+c.c.,$ (3)

where ${\widehat{e}}_f$ is collinear with ${\widehat{e}}_{\text{'}p}$ . This vibration is supposed to be combinational active and solitary.

The fields ${\overrightarrow{E}}_{l,s,p}$ are related through the nonlinear part of the polarization ${\overrightarrow{P}}_{l,s,p}\left(\overrightarrow{r},t\right)$. The sources of transient anisotropy are the processes of SRS where the polarization is proportional to the cubic power of the amplitude [22]. In this case, the strength of the nonlinear polarization along a specific crystal axis scales with the projection of the electric field strength onto that exact axis ($P_l\propto \left({\gamma }_{ll}{\left|A_l\right|}^2{A}_l+{\gamma }_{ls}^{\left(\mu \right)}{\left|A_s^{\left(\mu \right)}\right|}^2{A}_l\right)$ and $P_s\propto \left({\gamma }_{sl}^{\left(\mu \right)}{\left|A_l\right|}^2{A}_s^{\left(\mu \right)}+{\gamma }_{ss}^{\left(\mu \right)}{\left|A_s^{\left(\mu \right)}\right|}^2{A}_s^{\left(\mu \right)}\right)$ ).

Here we consider the system of right-handed unit vectors with the z-axis along which the electromagnetic waves propagate (z-axis is to the right), the x-axis is upwards (${\widehat{e}}_l={\widehat{e}}_s^{\left(1\right)}={\widehat{e}}_s^{\left(x\right)}=\widehat{i}$ ), and the y-axis is out of the page (${\widehat{e}}_s^{\left(2\right)}={\widehat{e}}_s^{\left(y\right)}=\widehat{j}$ ). The corresponding polarizations can be presented as follows:

$\begin{array}{l}{P}_{l,x}\propto ({\gamma }_{xxxx}{e}_l^{\left(x\right)}{e}_l^{\left(x\right)}{e}_l^{\left(x\right)}{e}_l^{\left(x\right)}{\left|A_l\right|}^2{A}_l+{\gamma }_{xxxx}{e}_l^{\left(x\right)}{e}_s^{\left(x\right)}{e}_s^{\left(x\right)}{e}_l^{\left(x\right)}{\left|A_s^{\left(x\right)}\right|}^2{A}_l\\ +{\gamma }_{xyyx}{e}_l^{\left(x\right)}{e}_s^{\left(y\right)}{e}_s^{\left(y\right)}{e}_l^{\left(x\right)}{\left|A_s^{\left(y\right)}\right|}^2{A}_l),\end{array}$ (4)

where:

$\begin{array}{l}{\gamma }_{ll}\equiv {\gamma }_{xxxx}{e}_l^{\left(x\right)}{e}_l^{\left(x\right)}{e}_l^{\left(x\right)}{e}_l^{\left(x\right)},\\ {\gamma }_{ls}^{\left(x\right)}\equiv {\gamma }_{xxxx}{e}_l^{\left(x\right)}{e}_s^{\left(x\right)}{e}_s^{\left(x\right)}{e}_l^{\left(x\right)}\equiv {\gamma }_{ls}^{\left(1\right)},\\ {\gamma }_{ls}^{\left(y\right)}\equiv {\gamma }_{xyyx}{e}_l^{\left(x\right)}{e}_s^{\left(y\right)}{e}_s^{\left(y\right)}{e}_l^{\left(x\right)}\equiv {\gamma }_{ls}^{\left(2\right)}.\end{array}$ (5)

$\begin{array}{l}{P}_{s,x}\propto \left({\gamma }_{xxxx}{e}_s^{\left(x\right)}{e}_l^{\left(x\right)}{e}_l^{\left(x\right)}{e}_s^{\left(x\right)}{\left|A_l\right|}^2{A}_s^{\left(x\right)}+{\gamma }_{xxxx}{e}_s^{\left(x\right)}{e}_s^{\left(x\right)}{e}_s^{\left(x\right)}{e}_s^{\left(x\right)}{\left|A_s^{\left(x\right)}\right|}^2{A}_s^{\left(x\right)}\right),\\ P_{s,y}\propto \left({\gamma }_{yxxy}{e}_s^{\left(y\right)}{e}_l^{\left(x\right)}{e}_l^{\left(x\right)}{e}_s^{\left(y\right)}{\left|A_l\right|}^2{A}_s^{\left(y\right)}+{\gamma }_{yyyy}{e}_s^{\left(y\right)}{e}_s^{\left(y\right)}{e}_s^{\left(y\right)}{e}_s^{\left(y\right)}{\left|A_s^{\left(y\right)}\right|}^2{A}_s^{\left(y\right)}\right),\end{array}$ (6)

where:

$\begin{array}{l}{\gamma }_{sl}^{\left(x\right)}\equiv {\gamma }_{sl}^{\left(1\right)}\equiv {\gamma }_{xxxx}{e}_s^{\left(x\right)}{e}_l^{\left(x\right)}{e}_l^{\left(x\right)}{e}_s^{\left(x\right)},{\gamma }_{ss}^{\left(x\right)}\equiv {\gamma }_{ss}^{\left(1\right)}\equiv {\gamma }_{xxxx}{e}_s^{\left(x\right)}{e}_s^{\left(x\right)}{e}_s^{\left(x\right)}{e}_s^{\left(x\right)},\\ {\gamma }_{sl}^{\left(y\right)}\equiv {\gamma }_{sl}^{\left(2\right)}\equiv {\gamma }_{yxxy}{e}_s^{\left(y\right)}{e}_l^{\left(x\right)}{e}_l^{\left(x\right)}{e}_s^{\left(y\right)},{\gamma }_{ss}^{\left(y\right)}\equiv {\gamma }_{ss}^{\left(2\right)}\equiv {\gamma }_{yyyy}{e}_s^{\left(y\right)}{e}_s^{\left(y\right)}{e}_s^{\left(y\right)}{e}_s^{\left(y\right)}.\end{array}$ (7)

Thus, because of the presence of the components of the fourth-rank susceptibility tensor ${\gamma }_{xxxx}$ , ${\gamma }_{xyyx}$ , and ${\gamma }_{yxxy}={\gamma }_{xyxy}$ , the third-order nonlinear anisotropy parameter $\sigma$ is defined as $\sigma ={\gamma }_{xxxx}-\left({\gamma }_{xxyy}+{\gamma }_{xyxy}+{\gamma }_{xyyx}\right)$ , will not be zero, which means that the cubic crystal becomes anisotropic.

The shortened equations for the amplitudes $A_{l,s,p}^{\left(\mu ,\sigma \right)}$ are obtained from Maxwell’s equations by the standard procedure [31]. The phonon equation is obtained by treating the vibration as an oscillator under an external force ($\sigma ,\mu =1,2$ ):

$\begin{array}{c}\frac{\partial A_l^{*}}{\partial z}+\frac{1}{v_l^z}\frac{\partial A_l^{*}}{\partial t}=-i\frac{2\pi {\omega }_l}{c{n}_l\cos {\theta }_l^z}({\chi }_0{A}_s^{\left(\mu \right)*}{A}_p^{\left(\sigma \right)*}+\frac{{\eta }_f{\omega }_f^2{\alpha }^{\left(f\right)}}{g^{\left(f\right)}}{A}_s^{\left(\mu \right)*}{\tilde{u}}_f^{\left(\sigma \right)*}\\ \begin{array}{c}\\ \end{array}+{\gamma }_{ll}{\left|A_l\right|}^2{A}_l^{*}+{\gamma }_{ls}^{\left(\mu \right)}{\left|A_s^{\left(\mu \right)}\right|}^2{A}_l^{*}),\end{array}$ (8)

$\begin{array}{c}\frac{\partial A_s^{\left(\mu \right)*}}{\partial z}+\frac{1}{v_s^{z\left(\mu \right)}}\frac{\partial A_s^{\left(\mu \right)*}}{\partial t}=-i\frac{2\pi {\omega }_s}{c{n}_s^{\left(\mu \right)}\cos {\theta }_s^{z\left(\mu \right)}}({\chi }_0{A}_l^{*}{A}_p^{\left(\sigma \right)}+\frac{{\eta }_f{\omega }_f^2{\alpha }^{\left(f\right)}}{g^{\left(f\right)}}{A}_l^{*}{\tilde{u}}_f^{\left(\sigma \right)}\\ \begin{array}{c}\\ \end{array}+{\gamma }_{sl}^{\left(\mu \right)}{\left|A_l\right|}^2{A}_s^{\left(\mu \right)*}+{\gamma }_{ss}^{\left(\mu \right)}{\left|A_s^{\left(\mu \right)}\right|}^2{A}_s^{\left(\mu \right)*}),\end{array}$ (9)

$\begin{array}{l}\frac{\partial A_p^{\left(\sigma \right)\ast }}{\partial z}+\frac{1}{v_p^{z\left(\sigma \right)}}\frac{\partial A_p^{\left(\sigma \right)\ast }}{\partial t}-\frac{W{e}_p^{\left(\sigma \right)z}}{2W^z}\frac{\partial A_p^{\left(3\right)*}}{\partial z}\\ =i\frac{q_p^2{\epsilon }_p^{\infty }}{2W^z}\left(\frac{W^2}{q_p^2{\epsilon }_p^{\infty }}-1\right)A_p^{\left(\sigma \right)\ast }-i\frac{s_f{\omega }_f^2{\omega }_p^2}{2g^{\left(f\right)}{W}^z{c}^2}{\tilde{u}}_f^{\left(\sigma \right)*}\\ -i\frac{2\pi q_p^2}{W^z}\left({\gamma }_{pl}^{\left(\sigma \right)}{\left|A_l\right|}^2{A}_p^{\left(\sigma \right)\ast }+{\gamma }_{ps}^{\left(\mu \right)}{\left|A_s^{\left(\mu \right)}\right|}^2{A}_p^{\left(\sigma \right)\ast }\right),\end{array}$ (10)

$\begin{array}{l}-iW{e}_p^{\left(\sigma \right)z}\frac{\partial A_p^{\left(\sigma \right)\ast }}{\partial z}+i\frac{\partial A_p^{\left(3\right)*}}{\partial z}\left(W^z-W{e}_p^{\left(3\right)z}\right)+i\frac{2{\omega }_p{\epsilon }_p^{\infty }}{c^2}\frac{\partial A_p^{\left(3\right)\ast }}{\partial t}-k_p^2{A}_p^{\left(3\right)*}\\ -\frac{s_f{\omega }_f^2{\omega }_p^2}{g^{\left(f\right)}{c}^2}{\tilde{u}}_f^{\left(3\right)*}=4\pi q_p^2\left({\gamma }_{pl}^{\left(\sigma \right)}{\left|A_l\right|}^2{A}_p^{\left(\sigma \right)*}+{\gamma }_{ps}^{\left(\mu \right)}{\left|A_s^{\left(\mu \right)}\right|}^2{A}_p^{\left(\sigma \right)*}\right),\end{array}$ (11)

$\begin{array}{c}\frac{\partial {\tilde{u}}_f^{\left(\sigma \right)*}}{\partial t}+\frac{{\gamma }_f}{2}{\stackrel{⌢}{u}}_f^{\left(\sigma \right)*}=i\Delta \omega {\tilde{u}}_f^{\left(\sigma \right)*}-i\frac{g^{\left(f\right)}}{2{\omega }_f}{A}_p^{\left(\sigma \right)*}\\ -i{\alpha }^{\left(f\right)}{\alpha }^{\left(f\right)}{}^{\prime }2\pi /\left(h{\upsilon }_0\right)\left({\left|A_l\right|}^2+{\left|A_s^{\left(\mu \right)}\right|}^2\right){\stackrel{⌢}{u}}_f^{\left(\sigma \right)*},\end{array}$ (12)

where: ${\chi }_0$ and ${\gamma }_{ll},{\gamma }_{ls}^{\left(\mu \right)},{\gamma }_{sl}^{\left(\mu \right)},{\gamma }_{ss}^{\left(\mu \right)},{\gamma }_{pl}^{\left(\sigma \right)},{\gamma }_{ps}^{\left(\mu \text{,}\sigma \right)}$ are the corresponding tensor contractions of non-resonant quadratic and cubic nonlinear polarizabilities with unit vectors of polarization of interacting waves; ${\eta }_f={\left(s_f/\left(2\pi h{\omega }_f{v}_0\right)\right)}^{1/2}$ ; $v_0$ is the volume of the elementary cell; $s_f$ is the oscillator strength of the oscillation ${\omega }_f$ ; $g^{\left(f\right)}$ is the ratio of effective charge to effective mass for a given vibration; ${\epsilon }_p^{\infty }$ is the non-resonant part of dielectric permeability at frequency ${\omega }_p$ ; $q_{l,s,p}={\omega }_{l,s,p}/c$ ; ${\alpha }^{\left(f\right)}=e_s^i{e}_l^j{\alpha }_{ij}^{\left(f\right)}$ ; ${\alpha }_{ij}^{\left(f\right)}$ is the tensor of combinational scattering per one elementary cell; ${\tilde{u}}_f^{\left(\sigma \right)}=u_f^{\left(\sigma \right)}\cdot \exp \left(-i\Delta \omega t\right)$ ; $\Delta \omega ={\omega }_f-{\omega }_p$ ; ${\gamma }_f$ is the constant of phonon attenuation; ${\upsilon }_{l,s,p}^z$ -z- components of velocities of waves on ${\omega }_{l,s,p}$ . Here we consider the processes in the vicinity of the phonon resonance ($\Delta \omega ={\omega }_f-{\omega }_p$ ). Consequently, we omit the equation simulating the evolution of the difference between the populations of the upper and lower levels characterized by T₁ (the time of longitudinal relaxation) while keeping the equation for phonon vibration. Nonstationary processes caused by the ultrashort pulses result in minimization of dissipation since the duration of the ultrashort pulses $\tau$ obeys $\tau \ll {\gamma }_f^{-1}$ (${\gamma }_f^{-1}$ after averaging over all vibrations would become T₂ (the time of transverse relaxation)). That is why we omit the term ${\gamma }_f{\stackrel{⌢}{u}}_f^{\left(\sigma \right)*}/2$ and keep $i\Delta \omega {\tilde{u}}_f^{\left(\sigma \right)*}$ , which makes sense in the non-resonant range of polariton frequencies $\left|\Delta \omega \right|\gg {\gamma }_f/2$ . As was shown in [22], in the case of non-resonant polariton frequencies, the contribution of $A_p^{\left(3\right)}$ can be neglected. To illustrate it, let’s consider the case of strong absorption:

$\left|\frac{\partial A_p^{\left(\sigma \right)*}}{\partial z}+\frac{1}{v_p^{z\left(\sigma \right)}}\frac{\partial A_p^{\left(\sigma \right)*}}{\partial t}\right|\ll \frac{q_p^2{\epsilon }_p^{\infty }}{2W^z}\left(\frac{W^2}{q_p^2{\epsilon }_p^{\infty }}-1\right)A_p^{\left(\sigma \right)*},$ (13)

when we can neglect in (10) and (11) the terms with derivatives, after which these equations yield

$\begin{array}{l}\left(W^2-q_p^2{\epsilon }_p^{\infty }\right)A_p^{\left(\sigma \right)*}=\frac{s_f{\omega }_f^2{\omega }_p^2}{g^{\left(f\right)}{c}^2}{\tilde{u}}_f^{\left(\sigma \right)*}+4\pi q_p^2\left({\gamma }_{pl}^{\left(\sigma \right)}{\left|A_l\right|}^2{A}_p^{\left(\sigma \right)*}+{\gamma }_{ps}^{\left(\mu \right)}{\left|A_s^{\left(\mu \right)}\right|}^2{A}_p^{\left(\sigma \right)*}\right),\\ -k_p^2{A}_p^{\left(3\right)*}=\frac{s_f{\omega }_f^2{\omega }_p^2}{g^{\left(f\right)}{c}^2}{\tilde{u}}_f^{\left(3\right)*}+4\pi q_p^2\left({\gamma }_{pl}^{\left(\sigma \right)}{\left|A_l\right|}^2{A}_p^{\left(\sigma \right)*}+{\gamma }_{ps}^{\left(\mu \right)}{\left|A_s^{\left(\mu \right)}\right|}^2{A}_p^{\left(\sigma \right)*}\right).\end{array}$

Since the right parts are approximately the same, we equate the left parts as well.

$\left|\left(W^2-q_p^2{\epsilon }_p^{\infty }\right)A_p^{\left(\sigma \right)*}\right|\approx \left|-k_p^2{A}_p^{\left(3\right)*}\right|$ (14)

We assume that W corresponds to the polariton dispersion as $W^2=q_p^2{{\epsilon }^{\prime }}_p$ $\left({{\epsilon }^{\prime }}_p=\text{Re}\left\{{\epsilon }_p\right\}\right)$ , $k_p^2=q_p^2{{\epsilon }^{\prime }}_p$ . This allows us to express the ratio of amplitudes as

$\left|\frac{A_p^{\left(3\right)*}}{A_p^{\left(\sigma \right)*}}\right|\approx \left|\frac{W^2-q_p^2{\epsilon }_p^{\infty }}{k_p^2}\right|=\left|\frac{{{\epsilon }^{\prime }}_p-{\epsilon }_p^{\infty }}{{{\epsilon }^{\prime }}_p}\right|=\left|1-\frac{{\epsilon }_p^{\infty }}{{{\epsilon }^{\prime }}_p}\right|$ (15)

This ratio can be generalized for both the lower polariton branch (${\omega }_p<{\omega }_f$ ) and the upper one (${\omega }_p>{\omega }_f$ ) as

$\left|\frac{A_p^{\left(3\right)*}}{A_p^{\left(\sigma \right)*}}\right|\approx \left|1\pm \frac{{\epsilon }_p^{\infty }}{{{\epsilon }^{\prime }}_p}\right|,$ (16)

Where the sign + corresponds to the case of the upper branch in which the polariton frequency approaches the resonance from the right. In this case, ${{\epsilon }^{\prime }}_p$ is negative and becoming even more negative (due to the normal dispersion), resulting in increasing (16) toward unity. And conversely, if the polariton frequency were to go to ${\omega }_f$ from the left, ${{\epsilon }^{\prime }}_p$ is positive and becoming more positive (it corresponds to the – sign in (16)).

After applying these approximations, the final system of equations is

$\begin{array}{c}\frac{\partial A_l^{*}}{\partial z}+\frac{1}{v_l^z}\frac{\partial A_l^{*}}{\partial t}=-i\frac{2\pi {\omega }_l}{c{n}_l\cos {\theta }_l^z}({\chi }_0{A}_s^{\left(\mu \right)*}{A}_p^{\left(\sigma \right)*}+\frac{{\eta }_f{\omega }_f^2{\alpha }^{\left(f\right)}}{g^{\left(f\right)}}{A}_s^{\left(\mu \right)*}{\tilde{u}}_f^{\left(\sigma \right)*}\\ \begin{array}{c}\\ \end{array}+{\gamma }_{ll}{\left|A_l\right|}^2{A}_l^{*}+{\gamma }_{ls}^{\left(\mu \right)}{\left|A_s^{\left(\mu \right)}\right|}^2{A}_l^{*}),\end{array}$ (17)

$\begin{array}{c}\frac{\partial A_s^{\left(\mu \right)*}}{\partial z}+\frac{1}{v_s^{z\left(\mu \right)}}\frac{\partial A_s^{\left(\mu \right)*}}{\partial t}=-i\frac{2\pi {\omega }_s}{c{n}_s^{\left(\mu \right)}\cos {\theta }_s^{z\left(\mu \right)}}({\chi }_0{A}_l^{*}{A}_p^{\left(\sigma \right)}+\frac{{\eta }_f{\omega }_f^2{\alpha }^{\left(f\right)}}{g^{\left(f\right)}}{A}_l^{*}{\tilde{u}}_f^{\left(\sigma \right)}\\ \begin{array}{c}\\ \end{array}+{\gamma }_{sl}^{\left(\mu \right)}{\left|A_l\right|}^2{A}_s^{\left(\mu \right)*}+{\gamma }_{ss}^{\left(\mu \right)}{\left|A_s^{\left(\mu \right)}\right|}^2{A}_s^{\left(\mu \right)*}),\end{array}$ (18)

$\begin{array}{c}\frac{\partial A_p^{\left(\sigma \right)\ast }}{\partial z}+\frac{1}{v_p^{z\left(\sigma \right)}}\frac{\partial A_p^{\left(\sigma \right)\ast }}{\partial t}=i\frac{q_p^2{\epsilon }_p^{\infty }}{2W^z}\left(\frac{W^2}{q_p^2{\epsilon }_p^{\infty }}-1\right)A_p^{\left(\sigma \right)\ast }-i\frac{s_f{\omega }_f^2{\omega }_p^2}{2g^{\left(f\right)}{W}^z{c}^2}{\tilde{u}}_f^{\left(\sigma \right)*}\\ -i\frac{2\pi q_p^2}{W^z}\left({\gamma }_{pl}^{\left(\sigma \right)}{\left|A_l\right|}^2{A}_p^{\left(\sigma \right)*}+{\gamma }_{ps}^{\left(\mu ,\sigma \right)}{\left|A_s^{\left(\mu \right)}\right|}^2{A}_p^{\left(\sigma \right)*}\right),\end{array}$ (19)

$\frac{\partial {\tilde{u}}_f^{\left(\sigma \right)*}}{\partial t}=i\Delta \omega {\tilde{u}}_f^{\left(\sigma \right)*}-i\frac{g^{\left(f\right)}}{2{\omega }_f}{A}_p^{\left(\sigma \right)*}-i{\alpha }^{\left(f\right)}{\alpha }^{\left(f\right)}{}^{\prime }2\pi /\left(h{\upsilon }_0\right)\left({\left|A_l\right|}^2+{\left|A_s^{\left(\mu \right)}\right|}^2\right){\stackrel{⌢}{u}}_f^{\left(\sigma \right)*},$ (20)

$\sigma ,\mu =1,2.$

The system’s structure reflects the different temporal responses of non-resonant electronic and resonant phonon components of the quadratic and cubic polarizabilities to ultrashort laser pulses. The non-resonant components respond almost instantaneously, whereas the resonant components respond with a delay arising from inertia and the mass difference between electrons and ions. In the equations, the terms associated with non-resonant polarizabilities are multiplied by the corresponding electromagnetic field strengths, while the resonant terms include ${\tilde{u}}_f^{\left(\sigma \right)*}$ . In the stationary case ($\frac{\text{d}{\tilde{u}}_f^{\left(\sigma \right)*}}{\text{d}t}=0$) in (20), substituting the solution for ${\tilde{u}}_f^{\left(\sigma \right)*}$ into (17) - (19) and collecting like terms yields the standard expressions for the polarizabilities, with their non-resonant and resonant components combined. This system generalizes the nonlinear systems considered earlier in [32]-[34] by allowing waves with different polarizations to be coupled through the polariton and phonon amplitudes.

3. The Manley-Rowe Relation for Waves Propagating with Different Polarizations at Frequencies ${\omega }_{l,s,p,f}$

We consider nonstationary SRS by polaritons when the polariton frequencies lie close to the phonon resonance, $\Delta \omega \gg {\gamma }_f/2$ . In this regime, the phonon contribution to SRS can be described through the phonon displacement, which forms the basis for the transverse polarization with the characteristic transverse relaxation time T₂, rather than through the population difference between the upper and lower levels, as in the resonant case. Because the polarization-induction and electromagnetic-wave re-emission processes are very fast ($\tau \ll T_2$ ), energy exchange between the propagating electromagnetic waves and the medium is expected to be negligible. Consequently, the Manley-Rowe relation is expected to hold only for the average characteristics of the electromagnetic waves.

To derive the Manley-Rowe relation, we first introduce dimensionless amplitudes that normalize the pump, Stokes, polariton, and phonon variables with respect to the characteristic laser-pulse amplitude and duration:

$\begin{array}{l}{\tilde{M}}_l\equiv {\left(\frac{c{n}_l\cos {\theta }_l^z}{2\pi {\omega }_l{z}_0}\right)}^{1/2}{\tilde{A}}_l,{\tilde{M}}_s^{\left(\mu \right)}\equiv {\left(\frac{c{n}_s^{\left(\mu \right)}\cos {\theta }_s^{z\left(\mu \right)}}{2\pi {\omega }_s{z}_0}\right)}^{1/2}{\tilde{A}}_s^{\left(\mu \right)},\\ {\tilde{M}}_p^{\left(\sigma \right)}\equiv {\left(\frac{c{n}_p^{\left(\sigma \right)}\cos {\theta }_p^{z\left(\sigma \right)}}{2\pi {\omega }_p{z}_0}\right)}^{1/2}{\tilde{A}}_p^{\left(\sigma \right)},{\tilde{M}}_f^{\left(\sigma \right)}\equiv {\left(\frac{c{n}_p^{\left(\sigma \right)}\cos {\theta }_p^{z\left(\sigma \right)}}{2\pi {\omega }_p{z}_0}\right)}^{1/2}{\tilde{\tilde{u}}}_f^{\left(\sigma \right)},\end{array}$

where ${\tilde{A}}_l,{\tilde{A}}_s^{\left(\mu \right)},{\tilde{A}}_p^{\left(\sigma \right)}\equiv A_l/A_0,A_s^{\left(\mu \right)}/A_0,A_p^{\left(\sigma \right)}/A_0$ ; ${\tilde{\tilde{u}}}_f^{\left(\sigma \right)}\equiv {\tilde{u}}_f^{\left(\sigma \right)}/z_0$ ; $z_0=c{\tau }_0$ , $A_0$ and ${\tau }_0$ are the peak amplitude and characteristic duration of the laser pump pulse. The amplitudes ${\tilde{A}}_{l,s,p}^{\left(\mu ,\sigma \right)}$ and ${\tilde{u}}_f^{\left(\sigma \right)}$ are complex because they have phase differences relative to the real-valued quantity $A_0$ . By contrast, intensities and photon numbers are real. We therefore use the absolute values of ${\tilde{M}}_{l,s,p}^{\left(\mu ,\sigma \right)}$ . The quantities ${\left|{\tilde{M}}_{l,s,p}^{\left(\mu ,\sigma \right)}\right|}^2$ are proportional to the instantaneous number of photons generated at a certain frequency per unit time at a certain position relative to the total energy of pump photons passed the boundary. For phonons, we do it formally, keeping in mind that the group velocity of phonons is small ($\left|\frac{\text{d}{\tilde{u}}_f^{\left(\sigma \right)}}{\text{d}z}\right|\approx 0$).With these definitions, system (17)-(20) can be written in dimensionless form as follows:

$\begin{array}{c}{\tilde{M}}_l\left(\frac{\partial {\tilde{M}}_l^{*}}{\partial \tilde{z}}+\frac{1}{{\tilde{v}}_l^z}\frac{\partial {\tilde{M}}_l^{*}}{\partial \tilde{t}}\right)=-i(C1{\tilde{M}}_l{\tilde{M}}_s^{\left(\mu \right)*}{\tilde{M}}_p^{\left(\sigma \right)*}+C2{\tilde{M}}_l{\tilde{M}}_s^{\left(\mu \right)*}{\tilde{M}}_f^{\left(\sigma \right)*}\\ +C3{\left|{\tilde{M}}_l\right|}^4+C{4}^{\left(\mu \right)}{\left|{\tilde{M}}_l\right|}^2{\left|{\tilde{M}}_s^{\left(\mu \right)}\right|}^2),\end{array}$ (21)

$\begin{array}{c}{\tilde{M}}_l^{*}\left(\frac{\partial {\tilde{M}}_l}{\partial \tilde{z}}+\frac{1}{{\tilde{v}}_l^z}\frac{\partial {\tilde{M}}_l}{\partial \tilde{t}}\right)=i(C1{\tilde{M}}_l^{*}{\tilde{M}}_s^{\left(\mu \right)}{\tilde{M}}_p^{\left(\sigma \right)}+C2{\tilde{M}}_l^{*}{\tilde{M}}_s^{\left(\mu \right)}{\tilde{M}}_f^{\left(\sigma \right)}\\ +C3{\left|{\tilde{M}}_l\right|}^4+C{4}^{\left(\mu \right)}{\left|{\tilde{M}}_l\right|}^2{\left|{\tilde{M}}_s^{\left(\mu \right)}\right|}^2),\end{array}$ (22)

$\begin{array}{c}{\tilde{M}}_s^{\left(\mu \right)}\left(\frac{\partial {\tilde{M}}_s^{\left(\mu \right)*}}{\partial \tilde{z}}+\frac{1}{{\tilde{v}}_s^{\left(\mu \right)z}}\frac{\partial {\tilde{M}}_s^{\left(\mu \right)*}}{\partial \tilde{t}}\right)=-i(C1{\tilde{M}}_l^{*}{\tilde{M}}_s^{\left(\mu \right)}{\tilde{M}}_p^{\left(\sigma \right)}+C2{\tilde{M}}_l^{*}{\tilde{M}}_s^{\left(\mu \right)}{\tilde{M}}_f^{\left(\sigma \right)}\\ +C{4}^{\left(\mu \right)}{\left|{\tilde{M}}_l\right|}^2{\left|{\tilde{M}}_s^{\left(\mu \right)}\right|}^2+C{5}^{\left(\mu \right)}{\left|{\tilde{M}}_s^{\left(\mu \right)}\right|}^4),\end{array}$ (23)

$\begin{array}{c}{\tilde{M}}_s^{\left(\mu \right)*}\left(\frac{\partial {\tilde{M}}_s^{\left(\mu \right)}}{\partial \tilde{z}}+\frac{1}{{\tilde{v}}_s^{\left(\mu \right)z}}\frac{\partial {\tilde{M}}_s^{\left(\mu \right)}}{\partial \tilde{t}}\right)=i(C1{\tilde{M}}_l{\tilde{M}}_s^{\left(\mu \right)*}{\tilde{M}}_p^{\left(\sigma \right)*}+C2{\tilde{M}}_l{\tilde{M}}_s^{\left(\mu \right)*}{\tilde{M}}_f^{\left(\sigma \right)*}\\ +C{4}^{\left(\mu \right)}{\left|{\tilde{M}}_l\right|}^2{\left|{\tilde{M}}_s^{\left(\mu \right)}\right|}^2+C{5}^{\left(\mu \right)}{\left|{\tilde{M}}_s^{\left(\mu \right)}\right|}^4),\end{array}$ (24)

$\begin{array}{c}{\tilde{M}}_p^{\left(\sigma \right)}\left(\frac{\partial {\tilde{M}}_p^{\left(\sigma \right)*}}{\partial \tilde{z}}+\frac{1}{{\tilde{v}}_p^{\left(\sigma \right)z}}\frac{\partial {\tilde{M}}_p^{\left(\sigma \right)*}}{\partial \tilde{t}}\right)=iC6{\left|{\tilde{M}}_p^{\left(\sigma \right)}\right|}^2-iC7{\tilde{M}}_p^{\left(\sigma \right)}{\tilde{M}}_f^{\left(\sigma \right)*}\\ -iC{8}^{\left(\sigma \right)}{\left|{\tilde{M}}_l\right|}^2{\left|{\tilde{M}}_p^{\left(\sigma \right)}\right|}^2\\ -iC{9}^{\left(\mu ,\sigma \right)}{\left|{\tilde{M}}_S^{\left(\mu \right)}\right|}^2{\left|{\tilde{M}}_P^{\left(\sigma \right)}\right|}^2,\end{array}$ (25)

$\begin{array}{c}{\tilde{M}}_p^{\left(\sigma \right)*}\left(\frac{\partial {\tilde{M}}_p^{\left(\sigma \right)}}{\partial \tilde{z}}+\frac{1}{{\tilde{v}}_p^{\left(\sigma \right)z}}\frac{\partial {\tilde{M}}_p^{\left(\sigma \right)}}{\partial \tilde{t}}\right)=-iC6{\left|{\tilde{M}}_p^{\left(\sigma \right)}\right|}^2+iC7{\tilde{M}}_p^{\left(\sigma \right)*}{\tilde{M}}_f^{\left(\sigma \right)}\\ +iC{8}^{\left(\sigma \right)}{\left|{\tilde{M}}_l\right|}^2{\left|{\tilde{M}}_p^{\left(\sigma \right)}\right|}^2\\ +iC{9}^{\left(\mu ,\sigma \right)}{\left|{\tilde{M}}_S^{\left(\mu \right)}\right|}^2{\left|{\tilde{M}}_P^{\left(\sigma \right)}\right|}^2,\end{array}$ (26)

$\begin{array}{c}\left(\frac{s_f{\omega }_f^4{z}_0^2}{g^{{\left(f\right)}^2}{A}_0^2}\right)\times {\tilde{M}}_f^{\left(\sigma \right)}\frac{\partial {\tilde{M}}_f^{\left(\sigma \right)*}}{\partial \tilde{t}}=\left(\frac{s_f{\omega }_f^4{z}_0^2}{g^{{\left(f\right)}^2}{A}_0^2}\right)(iC10{\left|{\tilde{M}}_f^{\left(\sigma \right)}\right|}^2-iC11{\tilde{M}}_f^{\left(\sigma \right)}{\tilde{M}}_p^{\left(\sigma \right)*}\\ -iC12{\left|{\tilde{M}}_l\right|}^2{\left|{\tilde{M}}_f^{\left(\sigma \right)}\right|}^2-iC13{\left|{\tilde{M}}_s^{\left(\mu \right)}\right|}^2{\left|{\tilde{M}}_f^{\left(\sigma \right)}\right|}^2),\end{array}$ (27)

$\begin{array}{c}\left(\frac{s_f{\omega }_f^4{z}_0^2}{g^{{\left(f\right)}^2}{A}_0^2}\right)\times {\tilde{M}}_f^{\left(\sigma \right)*}\frac{\partial {\tilde{M}}_f^{\left(\sigma \right)}}{\partial \tilde{t}}=\left(\frac{s_f{\omega }_f^4{z}_0^2}{g^{{\left(f\right)}^2}{A}_0^2}\right)\times (-iC10{\left|{\tilde{M}}_f^{\left(\sigma \right)}\right|}^2+iC11{\tilde{M}}_f^{\left(\sigma \right)*}{\tilde{M}}_p^{\left(\sigma \right)}\\ +iC12{\left|{\tilde{M}}_l\right|}^2{\left|{\tilde{M}}_f^{\left(\sigma \right)}\right|}^2+iC13{\left|{\tilde{M}}_s^{\left(\mu \right)}\right|}^2{\left|{\tilde{M}}_f^{\left(\sigma \right)}\right|}^2),\end{array}$ (28)

where

$\tilde{z}=\frac{z}{z_0},\tilde{t}=\frac{t}{{\tau }_0},{\tilde{v}}_l^z=\frac{v_l^{\text{z}}}{c},{\tilde{v}}_s^{(\mu )z}=\frac{v_s^{(\mu )z}}{c},{\tilde{v}}_p^{(\sigma )z}=\frac{v_p^{(\sigma )z}}{c}$

$C1={\chi }_0{A}_0{\left(\frac{8{\pi }^3{z}_0^3{\omega }_l{\omega }_s{\omega }_p}{c^3{n}_l{n}_s^{\left(\mu \right)}{n}_p^{\left(\sigma \right)}\cos {\theta }_l^z\cos {\theta }_s^{z\left(\mu \right)}\cos {\theta }_p^{z\left(\sigma \right)}}\right)}^{1/2},$

$C2=\frac{{\eta }_f{\omega }_f^2{\alpha }^{\left(f\right)}{z}_0}{g^{\left(f\right)}}{\left(\frac{8{\pi }^3{z}_0^3{\omega }_l{\omega }_s{\omega }_p}{c^3{n}_l{n}_s^{\left(\mu \right)}{n}_p^{\left(\sigma \right)}\cos {\theta }_l^z\cos {\theta }_s^{z\left(\mu \right)}\cos {\theta }_p^{z\left(\sigma \right)}}\right)}^{1/2},$

$C3={\gamma }_{ll}{A}_0^2{\left(\frac{2\pi {\omega }_l{z}_0}{c{n}_l\cos {\theta }_l^z}\right)}^2,$

$C{4}^{\left(\mu \right)}={\gamma }_{ls}^{\left(\mu \right)}{A}_0^2\left(\frac{4{\pi }^2{\omega }_l{\omega }_s{z}_0^2}{c^2{n}_l{n}_s^{\left(\mu \right)}\cos {\theta }_l^z\cos {\theta }_s^{z\left(\mu \right)}}\right),$ $\left({\gamma }_{ls}^{\left(\mu \right)}\approx {\gamma }_{sl}^{\left(\mu \right)}\right)$

$C{5}^{\left(\mu \right)}={\gamma }_{ss}^{\left(\mu \right)}{A}_0^2{\left(\frac{2\pi {\omega }_s{z}_0}{c{n}_s^{\left(\mu \right)}\cos {\theta }_s^{z\left(\mu \right)}}\right)}^2,$

$C6=\frac{q_p^2{\epsilon }_p^{\infty }{z}_0}{2W^z}\left(\frac{W^2}{q_p^2{\epsilon }_p^{\infty }}-1\right),$

$C7=\frac{s_f{\omega }_f^2{\omega }_p^2{z}_0^2}{2g^{\left(f\right)}{W}^z{c}^2{A}_0},$

$C{8}^{\left(\sigma \right)}=\frac{2\pi q_p^2{z}_0}{W^z}{\gamma }_{pl}^{\left(\sigma \right)}{A}_0^2\left(\frac{2\pi {\omega }_l{z}_0}{c{n}_l\cos {\theta }_l^z}\right),$

$C{9}^{\left(\mu ,\sigma \right)}=\frac{2\pi q_p^2{z}_0}{W^z}{\gamma }_{ps}^{\left(\mu ,\sigma \right)}{A}_0^2\left(\frac{2\pi {\omega }_s{z}_0}{c{n}_s^{\left(\mu \right)}\cos {\theta }_s^{z\left(\mu \right)}}\right),$

$C10=\Delta \omega {\tau }_0,$

$C11=\frac{g^{\left(f\right)}{A}_0{\tau }_0}{2{\omega }_f{z}_0},$

$C12=\frac{2\pi A_0^2{\tau }_0{\alpha }^{\left(f\right)}{\alpha }^{\left(f\right)}{}^{\prime }}{h{\upsilon }_0}\left(\frac{2\pi {\omega }_l{z}_0}{c{n}_l\cos {\theta }_l^z}\right),$

$C13=\frac{2\pi A_0^2{\tau }_0{\alpha }^{\left(f\right)}{\alpha }^{\left(f\right)}{}^{\prime }}{h{\upsilon }_0}\left(\frac{2\pi {\omega }_s{z}_0}{c{n}_s^{\left(\mu \right)}\cos {\theta }_s^{z\left(\mu \right)}}\right).$

Adding the left-hand sides of Equations (21)-(28) and then combining the corresponding right-hand sides gives the balance relation:

$\begin{array}{l}\left(\frac{\partial }{\partial \tilde{z}}+\frac{1}{{\tilde{\nu }}_l^z}\frac{\partial }{\partial \tilde{t}}\right){\left|{\tilde{M}}_l\left(\tilde{z},\tilde{t}\right)\right|}^2+\left(\frac{\partial }{\partial \tilde{z}}+\frac{1}{{\tilde{v}}_s^{z\left(\mu \right)}}\frac{\partial }{\partial \tilde{t}}\right){\left|{\tilde{M}}_s^{\left(\mu \right)}\left(\tilde{z},\tilde{t}\right)\right|}^2\\ +\left(\frac{\partial }{\partial \tilde{z}}+\frac{1}{{\tilde{v}}_p^{z\left(\sigma \right)}}\frac{\partial }{\partial \tilde{t}}\right){\left|{\tilde{M}}_p^{\left(\sigma \right)}\left(\tilde{z},\tilde{t}\right)\right|}^2+\left(\frac{s_f{\omega }_f^4{z}_0^2}{g^{{\left(f\right)}^2}{A}_0^2}\right)\frac{\partial }{\partial \tilde{t}}{\left|{\tilde{M}}_f^{\left(\sigma \right)}\left(\tilde{t}\right)\right|}^2=0,\end{array}$ (29)