Four-wave mixing in spin–orbit coupled Bose–Einstein condensates

Traditionally, four-wave mixing (FWM) processes are associated with photon interactions via a non-linear polarization. This is a third-order parametric process in which two particles (from two so-called writing or pump beams—one particle from each beam) are annihilated when passing through a non-linear medium, and at the same time two new particles (constituting probe and signal beams) are generated. In optics, FWM is commonly associated with the third order Kerr nonlinearity. The phenomenon is ubiquitous (see e.g., [1, 2]) and its applications are very widespread, including: fiber optic communication (very often not welcome), wavelength conversion, parametric amplification, optical regeneration, optic phase conjugation and correction of the aberration of images.

FWM can be observed also for massive particles as it was predicted and observed in cold atomic gasses two decades ago [3–5] (see also [6]). In this case resonantly interacting particles are neutral atoms rather than photons. Flexibility of control of trapping potentials, as well as of nonlinear interactions, in atomic systems open interesting perspectives of managing both momentum conservation and energy conservation laws through the interplay of additional linear and nonlinear potentials. This issue was already explored in references [7–11]. In this paper we explore a similar idea of controlling FWM processes through artificially created gauge potentials.

In order to consider FWM in a specific medium, one has to identify the characteristic eigenmodes: the elementary solutions to the linearized equations of motion in the form of plane waves, identified by their wavevectors and frequencies which satisfy the so-called phase matching conditions, that are equivalent to momentum and energy conservation laws. These are often quite demanding constraints depending on the particular form of dispersion relation characteristic for the system under investigation. For instance, in one dimension they cannot be satisfied for a system of cold atoms obeying parabolic dispersion relation and confined to (quasi-)one dimension. The situation can be improved by artificial change of the dispersion law using linear optical lattices [7–9] or by the manipulation of the wavenumbers of the matter waves involved in the process by means of nonlinear lattices [10, 11]. These modifications introduce the internal texture to the propagation medium, making it inherently inhomogeneous.

If a system has a spinor nature, i.e., consists of two subsystems, an alternative way to manipulate the linear properties of the medium, even preserving homogeneity, is to employ coherent coupling of the constituents. In optics, for example, one can satisfy the matching conditions for the FWM of light propagating in homogeneous parity-time ($\mathcal{PT}$)-symmetric coupled waveguides with gain and losses [12]. One can also consider the matching conditions, and thus observation of the FWM in $\mathcal{PT}$-symmetric optical lattices which are available experimentally [13]. And yet another interesting application was demonstrated in multi-component vector solitons consisting of two perpendicular FWM dipole components created by electromagnetically induced gratings [14].

A similar situation naturally occurs for spinor Bose–Einstein condensate, where coupling between two atomic states by means of the spin–orbit coupling (SOC) allows one to manipulate the dispersion relation in the presence of external potential. This idea becomes attractive, since using various experimental techniques, spin–orbit-coupled Bose–Einstein condensates (SOC-BECs) of hyperfine states of ⁸⁷Rb atoms has been created [15–20]. Notice that atoms of s and p bands of the static lattice were considered as pseudospins [21].

The main goal of this paper is to show that with properly adjusted SOC, one can satisfy the phase matching conditions for homogeneous one-dimensional SOC-BECs. Importantly, the SOC properties in atomic systems are highly adjustable [22–25] and the matching conditions reported below are experimentally feasible. Inter-atomic interactions are also tunable, most commonly by Feshbach resonance; see reference [26] for the observation of Feshbach resonances for SOC fermions, or reference [27] for observation of partial waves, with nonlinear interactions controlled by the dressing technique. Additional degrees of freedom in manipulating the effective dispersion relation, and thus, of the matching conditions (see (8) and (9) below), may be reached by using moving lattices [17, 21].

Understanding the dynamics of the FWM process is important in the context of creation of pairs of correlated particles, which are crucial in the ultra-precise metrology. Particles undergoing elastic collisions may create entanglement and it can be of practical use in the experiment. Here we mention the experiments [30–32], study of Cauchy–Schwarz inequality [33, 34], and twin beam generation [35]. Also, recently there was a prediction [36] and its experimental confirmation [37], of the Bell inequality is violated in a many-body system of massive particles with spinor condensates.

In the context of present considerations we would like to mention two recent experimental achievements. In the first a stripe phase with supersolid properties in SOC-BECs has been observed [28]. In this case spin flip process with a momentum transfer has been realized and observed using Bragg scattering. In another experiment radio-frequency (rf) photons were dressed with tunable recoil momentum by combining rf pulses with an oscillating magnetic force. This leads to a new application of Floquet engineering: periodically driven systems can have time-averaged properties which cannot be achieved with constant fields and in our opinion hold a promise of mixing different waves [29].

Due to the spinor nature of the one-dimensional (1D) SOC-BEC, it is characterized by two branches of the dispersion relation. As a result, the matching conditions can be readily satisfied, as we shall see below. Moreover, unlike in the case BECs without SOC, where additional external potentials, like linear [7–9] or nonlinear [10, 11] lattices are employed for ensuring matching conditions, now one can find a diversity of distinct FWM processes, where the interacting waves, as well as waves generated may represent different spinor states at the same values of parameters of the system (similarly to the FWM with laser pulses reported in [12]). It is a goal of the present study to identify the FWM processes available in the FWM in SOC-BECs and show they can be efficient enough to be observed.

The organization of the paper is as follows: in section 2 we discuss phase matching conditions in a quasi-1D SOC-BEC, and we focus on the degenerate case, where the central pump wave serves as a source for stimulated enhancement of two probe waves. Here we borrow the terminology from optics but we make no distinction between signal and probe beams. Next we identify four possible configurations of FWM and in section 3 we perform feasibility study using real time simulations for all predicted configurations. The outcomes are summarized in conclusion.

2.1. General relations

Let us consider a quasi-1D SOC-BEC which is described by a two-component order parameter $\mathbf{\Psi }\left(x,t\right)={\left({{\Psi}}_{1}\left(x,t\right),{{\Psi}}_{2}\left(x,t\right)\right)}^{\mathrm{T}}$ (hereafter T stands for transposition). The dynamics of the spinor Ψ(x, t) is governed by the coupled Gross–Pitaevskii equations (GPEs):

$$\begin{equation}\mathrm{i}{\partial }_{t}\mathbf{\Psi }=H\mathbf{\Psi }+\frac{1}{2}G\left(\mathbf{\Psi }\right)\mathbf{\Psi },\end{equation} \tag{ 1 }$$

where

$$\begin{equation}H=\frac{1}{2}\left(-{\partial }_{x}^{2}+\mathbf{\Omega }\cdot \boldsymbol{\sigma }-\mathrm{i}\alpha {\sigma }_{x}{\partial }_{x}\right)\end{equation} \tag{ 2 }$$

is the linear mean field Hamiltonian of the two component BEC, α is the SOC strength, Ω is the vector of the Zeeman coupling (we admit external magnetic field), σ = (σ_x, σ_y, σ_z) is the vector of Pauli matrices σ_x,y,z. Physical interpretation of these parameters depends on the particular realization. For instance in the experiment performed in an ⁸⁷Rb Bose–Einstein condensate, a pair of Raman lasers created a momentum sensitive coupling between two internal atomic states [15]. This SOC was equivalent to that of an electronic system with equal contributions of Rashba and Dresselhaus couplings, and with a uniform magnetic field B in the (x, z)-plane. In materials the SOC is due to intrinsic properties, which are largely determined by the specific material and the details of its growth. In these and other proposed schemes Ω may have different physical interpretation, including Rabi frequencies of the dressing laser fields. It has to be taken into account when one consider tuning and range of vector Ω. Here we will call it Zeeman coupling, and Ω_z Zeeman splitting.

The nonlinearity 2 × 2 matrix is given by

$$\begin{equation}G\left(\mathbf{\Psi }\right)=\left(\begin{matrix}{cc}{g}_{1}\vert {{\Psi}}_{1}{\vert }^{2}+g\vert {{\Psi}}_{2}{\vert }^{2}& 0\\ 0& g\vert {{\Psi}}_{1}{\vert }^{2}+{g}_{2}\vert {{\Psi}}_{2}{\vert }^{2}\end{matrix}\right)\end{equation} \tag{ 3 }$$

with intra- and inter-component interactions g_1,2 and g, respectively, and we use the dimensionless units with ℏ = m = 1.

To address the matching conditions for the FWM we start with the eigenmodes of the linear spectral problem, representing them in the form of the plane waves

$$\begin{equation}{\mathbf{\Psi }}_{{\pm}}\left(x,t\right)={\mathrm{e}}^{\mathrm{i}kx-\mathrm{i}{\mu }_{{\pm}}\left(k\right)t}{\boldsymbol{\psi }}_{{\pm}}\left(k\right),\end{equation} \tag{ 4 }$$

where ${\boldsymbol{\psi }}_{{\pm}}\left(k\right)={\left({\psi }_{{\pm}}^{\left(1\right)}\left(k\right),{\psi }_{{\pm}}^{\left(2\right)}\left(k\right)\right)}^{\mathrm{T}}$ is a constant (i.e. x- and t-independent) spinor, k is a mode wavenumber, μ_±(k) is its frequency, and ± indicate the upper ('+') and lower ('−') branches of the spectrum (i.e., μ₋(k) ⩽ μ₊(k)). We will concentrate in the case of the magnetic field in the (x, z) plane, i.e. Ω = (Ω_x, 0, Ω_z) (without loss of generality we fix Ω_x, Ω_z ⩾ 0), for which we compute the two branches of the dispersion relation

$$\begin{equation}{\mu }_{{\pm}}\left(k\right)=\frac{{k}^{2}}{2}{\pm}\frac{\varepsilon \left(k\right)}{2}.\end{equation} \tag{ 5 }$$

Here

$$\begin{equation}\varepsilon \left(k\right)=\sqrt{{{\Omega}}_{z}^{2}+{{\sim}{{\Omega}}}^{2}\left(k\right)}\end{equation} \tag{ 6 }$$

with ${\sim}{{\Omega}}\left(k\right)=\alpha k+{{\Omega}}_{x}$, is the gap between the spectral branches at a given k. The lower (−) and the upper (+) branches of the respective eigenvectors are defined as

$$\begin{equation}{\boldsymbol{\psi }}_{{\pm}}\left(k\right)=\frac{1}{\sqrt{{{\sim}{{\Omega}}}^{2}+{\left(\varepsilon \left(k\right)\mp {{\Omega}}_{z}\right)}^{2}}}\left(\begin{pmatrix}{c}{\sim}{{\Omega}}\hfill \\ \hfill {\pm}\varepsilon \left(k\right)-{{\Omega}}_{z}\hfill \end{pmatrix}\right).\end{equation} \tag{ 7 }$$

It is to be emphasized that (4), as well as the matching conditions (8), (9), strictly speaking are valid only in linear regime, or approximately valid during the initial stages of the evolution. If however, the pulses with newly generated central frequencies separate in space from the initial pulse due to different group velocities, further generation of new harmonics is suppressed and one can observe propagation of quasi-monochromatic pulses. Confirmation of the described scenario requires numerical simulation of the full nonlinear problem. This approach is adopted below.

To reduce the number of parameters, here we investigate the degenerate FWM process, where two input spinor states are identical. We label their wavevectors by k₁ and using the optical terminology, we call them pump waves. Two spinors that are created in the FWM process with central wavevectors k₂ and k₃, will be referred to as probe waves. Respectively, the conservation of wavenumbers and frequencies of the pump and probe waves, are expressed in the from of phase-matching conditions

$$\begin{equation}2{k}_{1}={k}_{2}+{k}_{3},\end{equation} \tag{ 8 }$$

$$\begin{equation}2{\mu }_{{\nu }_{1}}\left({k}_{1}\right)={\mu }_{{\nu }_{2}}\left({k}_{2}\right)+{\mu }_{{\nu }_{3}}\left({k}_{3}\right).\end{equation} \tag{ 9 }$$

The indexes ν_j (j = 1, 2, 3) refer to either '+' or '−' branch of the spectrum.

Below we consider only the cases k_2,3 ≠ k₁, excluding the trivial case of self-phase modulation where all wavenumbers are equal. We note that the system (1) obeys gauge, rather than Gallilean, invariance, at which the generalized momentum

$$\begin{equation}{\Pi}={\int }_{-\infty }^{\infty }{\mathbf{\Psi }}^{{\dagger}}\left(-\mathrm{i}{\partial }_{x}+\frac{\alpha }{2}{\sigma }_{x}\right)\mathbf{\Psi }\mathrm{d}x\end{equation} \tag{ 10 }$$

is conserved: dΠ/dt = 0. For our consideration this means that the input wavenumber k₁ cannot be set arbitrarily to zero without changing the spinor eigenstates. It also means that at zero Zeeman field Ω = 0, the linear Hamiltonian H is gauge equivalent to the usual one-dimensional Schrödinger Hamiltonian ${H}_{0}=-{\partial }_{x}^{2}$ which does not support the matching conditions (8) and (9). In other words, while SOC controls the waves involved in resonant processes, the FWM itself requires nonzero Zeeman field.

For the following consideration it is convenient to rewrite (8) as

$$\begin{equation}{k}_{2}={k}_{1}+q,\quad {k}_{3}={k}_{1}-q.\end{equation} \tag{ 11 }$$

Now matching condition for frequencies (9) can be rewritten in the form

$$\begin{equation}2{q}^{2}=2{s}_{1}\varepsilon \left({k}_{1}\right)-{s}_{2}\varepsilon \left({k}_{1}+q\right)-{s}_{3}\varepsilon \left({k}_{1}-q\right),\end{equation} \tag{ 12 }$$

where s_j = ±1. Since each wave belongs to either upper or lower branch, these are eight different equations for given k₁ and q. However only four of them have nontrivial solutions. To justify this we first notice that if s₂ = s₃, the condition (12) is symmetric under q ↔ −q exchange. If however s₂ ≠ s₃, then (12) is symmetric with respect to simultaneous change (s₂, q) ↔ (s₃, −q). This allows one to restrict the analysis to the case q > 0. Next we use ɛ(k) > 0 property and conclude that the case (s₁, s₂, s₃) = (−1, 1, 1) does not have solutions since the right-hand side of equation (12) becomes negative. Let us now consider k₁ ⩾ 0 (the case k₁ < 0 is fully analogous). For non-negative values of k₁ we find the following inequalities

$$\begin{equation}0{\leqslant}{q}^{2}{\leqslant}{\varepsilon }_{+}\left({k}_{1}-q\right){\leqslant}{\varepsilon }_{+}\left({k}_{1}+q\right),\end{equation} \tag{ 13 }$$

excluding the cases (−1, 1, −1) and (−1, −1, 1). Hence, the initial pulse from the lower brunch of the spectrum may originate degenerate FWM in processes involving modes from the lowest branch only (this is the configuration 4 in the table 1 below). Finally, using inequalities (13) one can exclude also the case (1, 1, 1), leaving only four possible configurations summarized in the table 1.

Table 1. Possible configurations of degenerate FWM processes (positive and negative q are included). The first column are numbers identifying configurations, which corresponds to the specific choice of the spectrum branches for the pump and probe waves, indicated in second, third, and fourth columns, respectively. In the last column we show the maximal number, N_max, of q values solving equation (12).

Configuration	s1	s2	s3	Nmax
1	1	−1	−1	2
2	1	1	−1	2
3	1	−1	1	2
4	−1	−1	−1	4

In the last column of the table 1 we list the maximal number of solutions for particular configuration. In what follows we present analytical and graphic considerations that led us to these counts.

2.2. Analysis of possible configurations

By straightforward algebraic manipulations we can eliminate square-root terms in the equation (12) (simple sequence of transfers and squaring). As a result all four phase matched processes listed in the table 1 are determined by the following cubic equations

$$\begin{align}\hfill {Q}^{3}& -\left(1+4{s}_{1}\sqrt{{{\sim}{\omega }}^{2}+{\omega }_{z}^{2}}\right){Q}^{2}+\left(2{s}_{1}\sqrt{{{\sim}{\omega }}^{2}+{\omega }_{z}^{2}}+5{{\sim}{\omega }}^{2}+5{\omega }_{z}^{2}\right)Q\hfill \\ \hfill & -{\omega }_{z}^{2}-2{s}_{1}\sqrt{{{\sim}{\omega }}^{2}+{\omega }_{z}^{2}}\left({{\sim}{\omega }}^{2}+{\omega }_{z}^{2}\right)=0\hfill \end{align} \tag{ 14 }$$

where ${\sim}{\omega }={\sim}{{\Omega}}/{\alpha }^{2}$, ω_z = Ω_z/α² and Q = q²/α². For obvious reasons we are interested only in positive roots of (14) and exclude the root Q = 0 (i.e. q = 0) which does not correspond to FWM but to the self-phase modulation.

A number of real roots of equation (14) is determined by the sign of the discriminant

$$\begin{equation}{{\Delta}}_{{s}_{1}}={{\sim}{\omega }}^{2}\left[15{{\sim}{\omega }}^{2}+4{s}_{1}\sqrt{{{\sim}{\omega }}^{2}+{\omega }_{z}^{2}}\left({{\sim}{\omega }}^{2}+{\omega }_{z}^{2}+3\right)-12{\omega }_{z}^{2}-4\right].\end{equation} \tag{ 15 }$$

If ${{\Delta}}_{{s}_{1}}{ >}0$, there exists one real root. Three distinct real solutions exist if ${{\Delta}}_{{s}_{1}}{< }0$. At ${{\Delta}}_{{s}_{1}}=0$ all roots are real and at least one is multiple [38].

Now we inspect systematically configurations listed in table 1 which manifest qualitatively different types of dynamics. Starting with the last one, we set s₁ = −1. Now the discriminant (15) depends on the two parameters $\left\{{\sim}{\omega },{\omega }_{z}\right\}$ for different values of which ${{\Delta}}_{{s}_{1}}$ can acquire any sign or be zero. Analyzing Vieta formulae one can exclude the possibility of all three real roots being positive [38]. It means that in the configuration 4 (see the table 1) there may exist either one or two real positive roots of equation (14). Taking into account that roots appear in pairs, ±q, this corresponds to at most four possible arrangements allowed by the phase matching condition (12).

Next we turn to the configurations 1, 2, and 3 in the table 1, by setting s₁ = +1 in (14) and in (15). In this case we find that ${{\Delta}}_{{s}_{1}}{< }0$ for all values of the parameters $\left\{{\sim}{\omega },{\omega }_{z}\right\}$, i.e. all roots are real. Moreover, they are all positive as follows from the Vieta formulae. The three real positive roots of Q, i.e. six roots q, correspond to the three different configurations (notice that the difference among these configurations was removed upon squaring in obtaining equation (14)).

It is easy to establish one-to-one correspondence of the roots and the configurations. For if one root is common for two configurations we use equation (12) to show that one of the terms ɛ_s(k₁ ± q) is equal to zero. Obviously, for each of the configurations with s₁ = +1 always there is at least one positive root (all roots cannot be negative, because their product is positive). Let 0 < q₁, q₂, q₃ denote the positive roots, so that the configuration 1 has two symmetric roots, which we denote by ±q₁ while the roots of the configurations 2 and 3 are given by (q₂, −q₃) and (−q₂, q₃), respectively.

These arguments are exemplified in figure 1 where we show graphical solutions of phase matching equation. In the panels (a), obtained for α = 2, and (b), obtained for α = 10, the dashed blue curve is the plot of the left-hand side (LHS) of equation (12). Right-hand side (RHS) of this equation is represented by red, pink and black lines correspond to configurations 1, 2 and 3 in table 1. Crossings of solid and dashed lines yield the roots of the phase matching equations. Comparing the values of the LHS and RHS of equation (12) at q = 0 and at q → ∞, we conclude that each of the configurations 1, 2 and 3 must have at least one root for q > 0.The fourth configuration (s₁ = s₂ = s₃ = −1) is illustrated in panel (c). Here we observe three different possibilities: no positive solutions of the phase matching condition (12) (the light green curve does not cross the dash blue curve at α = 5); one positive root (the deep green and dash blue curves tangent to each other when α ≈ 7.67), and two positive solutions (the brown curve crosses the dash blue curve in two points, when α = 9).

Zoom In Zoom Out Reset image size

Finally, in the panel (d) we varied k₁ (while holding α = 9) and shown the regions where the phase matching equation of the fourth configuration supports one positive (the black curves) and two positive solutions (the pink curves). In all panels (a)–(d) we fixed Ω_x = 2.5 and Ω_z = 8.

2.3. Matching of group velocities

While no matching conditions on the group velocities is imposed, for practical observation of different scenarios of FWM in numerical simulations, the issue of the group velocities (GVs)

$$\begin{equation}{v}_{{\pm}}\left(k\right)\equiv \frac{\partial {\mu }_{{\pm}}\left(k\right)}{\partial k}=k{\pm}\frac{{\sim}{{\Omega}}}{2\sqrt{{{\sim}{{\Omega}}}^{2}+{{\Omega}}_{z}^{2}}},\end{equation} \tag{ 16 }$$

becomes relevant. On the one hand where all wavepackets move with respect to each other, it is important that the spinor involved in the process have similar values of GVs: otherwise fast separation of wavepackets in space may drastically reduce the conversion efficiency. On the other hand, GVs should have sufficient difference in order to observe spatial separation of the probe wavepackets. Thus in addition to solving the matching conditions we set a task of finding optimal conditions in the context of FWM numerical simulations (they are presented in the next section).

For each configuration listed in table 1, at a given k₁ one can determine q, i.e., the wavenumbers k₂ and k₃, and consequently their GVs.

In figure 2(a) we observe that the GV of the pump wavepacket, v₊(k₁), is close to either v₋(k₂) or v₋(k₃), almost for all k₁ except the vicinity of k₁ = 0. This means that to obtain clear separation of the wavepackets, generated in the FWM at relatively short time intervals, k₁ should be chosen close to zero. Then, the separation between velocities grows rapidly enough allowing direct observation of separated pulses. On the other hand, the time that pulses overlap is still long enough to generate a substantial four wave mixing signal.

Zoom In Zoom Out Reset image size

Figure 2(b) shows the dependence of GVs of phase matched wavepackets versus pump momentum for configurations 2 and 3 from the table 1. We again observe that in some regions GVs are close to each other or even coincide what does not allow observation of separation of the generated wavepackets from the initial one. However, an interesting situation occurs in the vicinity of k₁ = 0. Here GVs of the second and third waves have bigger absolute values than v₊(k₁) and the same sign. In this case both created waves move faster that the initial wavepacket. We should emphasize that this is not common in the usual realizations of FWM and this is solely due to the SOC coupling. Note that in these configurations it is also possible to initiate FWM process where one of the velocities of the probe waves is smaller and one bigger than that of the pump.The situation is more complicated for the case of the fourth configuration in table 1, as shown in figure 2(c). The (thick) red curve represents GV of the pump wave of the negative branch v₋(k₁) (see equation (16)). The other (black) curves, that have forms of three loops, represent GVs of generated waves [v₋(k₂) and v₋(k₃)] that correspond to other (non-trivial) solutions. Like in the previous cases, to reach significant separation of the pulses in the real space we choose k₁ in a region far from the crossing of the curves.

Equipped with the solutions of matching conditions and with the ideas of optimization the conversion efficiency in terms of the GVs we now turn to direct numerical simulations of the configurations of the FWM processes summarized in table 1. In order to find favorable conditions to observe clear evidence of specific FWM process, first one has to select proper momentum k₁. Note that phase matching will automatically determine all participating wavepackets GVs as explained in the previous section. Then, appropriate initial widths and amplitudes of the pump and probe waves need to be adjusted to ensure long enough and strong enough nonlinear interaction.

In all simulations we use the wavepackets having equal widths and completely overlapping at t = 0, i.e.,

$$\begin{equation}\mathbf{\Psi }\left(x,t=0\right)={\mathrm{e}}^{-{x}^{2}/{w}^{2}}\sum _{j=1}^{3}{A}_{j}{\boldsymbol{\psi }}_{{s}_{j}}\left({k}_{j}\right){\mathrm{e}}^{\mathrm{i}{k}_{j}x}.\end{equation} \tag{ 17 }$$

Here A_j are initial amplitudes of the wave-packets with the central wavevectors k_j of the spinors defined in accordance with equation (7).

For the FWM process corresponding to configuration 1, the initial state is formed with A₁ = 1, A₂ = 0.2, A₃ = 0 and s_j are chosen according to table 1: s₁ = 1, s₂ = −1 and s₃ = −1. In figure 3 we show an example of the FWM for this configuration. Due to the FWM process, by the expense of the highly populated initial state A₁ we observe strong amplification of the seed state and growth of the third matter wave with phase matched momentum k₃.

Zoom In Zoom Out Reset image size

This process is depicted with snapshots at the beginning and end (i.e., at t = 0 and t = 300) of the simulations in figure 3, where main panels (a) and (c) [(b) and (d)] refer to the first [second] spinor component. In particular blue contours in panels (a) and (b) represent initially overlapping pump (k₁) and probe (k₂) waves in the configuration space. They are fully separated after evolution time (t = 300) due to the difference in GVs and new, clearly visible, wave of central momentum k₃ is generated. Panels (c) and (d) show the corresponding features in the Fourier space. Here we distinct two waves as narrow blue peaks, at initial time and again three waves at the end of evolution. The most explicit feature is substantial broadening of all participating matter waves during the evolution. The inset in each panel shows full time evolution of modulus of the spinor components—(a) and (b) in the real space and (c) and (d) in the momentum space.

In the next two figures we illustrate the FWM process corresponding the second and third configuration from table 1, with GV configuration shown in panel (b) of the figure 2. As mentioned above in these two configurations the exist two different roots (q₂ ≠ q₃) of phase matching condition (12). As one can see directly in panels (a) and (b) of figure 1 these configurations are related by the transformation q₁ → −q₂ and q₂ → −q₁, i.e. the analysis of second and third configuration are analogous.

Interestingly, in figure 4 both probe and created waves are generated in the same side of the pump wave. The evolution of the probe wavepackets in figure 5 looks qualitatively similar to that one shown in figure 3. However, since each newborn wavepacket bears a quasi-spin, the emergent spinors (more precisely the left propagating waves) are different in these cases.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Turning to the fourth in the table (1), we recall that in this case one can obtain up to five solutions from the phase matching condition (including the trivial case of q = 0). We start with figure 6, where initial group velocities can be identified in the panel (c) of the figure 2 and are marked as dots on red and black curves. In principle, the dynamics presented in this case is very similar to that shown for the configurations 1 and 3, except that now different spinor states are involved (respectively the wavepackets bear different quasi-spins). Also closer look at figures 6(c) and (d), showing the spectra of the components reveals an interesting feature. Namely one can spot oscillations of the amplitude of pump wavepacket and we attribute them to the self phase modulation which was mentioned above as the trivial solution Q = 0 (or q = 0) of the equation (14).

Zoom In Zoom Out Reset image size

For completeness we present the fourth configuration in the case when the phase matching allows simultaneously for two FWM processes. In this last simulation we used parameters corresponding to the region, where the phase matching equation supports four different root: q₁, q₂, −q₁ and −q₂. Having these roots in hand, we define momenta of the two sets of probe waves in the following way:

$$\begin{align}\hfill {k}_{2}={k}_{1}+{q}_{1},\quad {k}_{3}={k}_{1}-{q}_{1}\\ \hfill {k}_{4}={k}_{1}+{q}_{2},\quad {k}_{5}={k}_{1}-{q}_{2}.\end{align}$$

In present case dynamics of the FWM process involves waves with the amplitudes: A₁ as a pump, and two pairs of probes A₂ and A₃, as well as A₄ and A₅ (by convention a mode with A_j has momentum k_j, where j = 1, ..., 5). In numerical simulations we initiate the dynamics by putting A₁ = 1, A₃ = A₅ = 0.2, and A₂ = A₄ = 0 (see equation (17); the values of the rest of parameters are listed in the figure 7 caption). In the nonlinear evolution the pump is interacting with both sets of the probe wavepackets creating two new waves A₂ and A₄ with momenta k₂ and k₄ respectively, in two simultaneous FWM processes. Figure 7 shows two components of the spinor wavefunctions in the configuration [(a) and (b)] and momentum [(c) and (d)] spaces, where we can clearly see two sets of probe waves. This time we can not refer to figure 2, since the values of parameters were slightly different from those used for panel (c), but for the sake of clarifications we propose to look at panel (d) of the figure 1. Due to the strong spreading there is substantial overlap of wavepackets, but when we look at them in the momentum space all peaks can be clearly identified. Close inspection reveals also oscillations on the pump wave due to the self-phase modulation mentioned above.

Zoom In Zoom Out Reset image size

In general finding optimal conditions to observe effective FWM process is not an easy task especially taking into account that the problem is multi-parametric and the values of parameters have to be carefully selected. In order to support this statement we have chosen the Zeeman splitting Ω_z as a control parameter and in figure 8 we present the efficiency (unit of per cent) of FWM defined as

$$\begin{equation}\eta \left({k}_{3}\right)=\frac{{\sim}{N}\left({k}_{3}\right)}{N}\cdot 100\%,\end{equation} \tag{ 18 }$$

where N is total number of atoms in the system and the number of atoms corresponding to the generated wave ${\sim}{N}\left({k}_{3}\right)$ is evaluated after long enough evolution time, when the wave packets are well separated. The efficiency is calculated for the FWM conversion in the process identified above as configuration 1 for different for different types of the nonlinearity parameterized by the quantities $\overline{g}=\left({g}_{1}+{g}_{2}\right)/2$ and Δg = (g₁ − g₂)/2. In all of the calculations presented in the figure 8 we took $g=\overline{g}$.

Zoom In Zoom Out Reset image size

In both panels of figure 8 we vary Ω_z while keeping the rest of the parameters fixed (see figure 3 for the exact values of all parameters). In the left panel different curves correspond to different values of the ratio Δg/g with g = 0.8, while in the right panel Δg is equal to zero and different curves were obtained choosing different values of the interaction strength g. In both panels we observe non-monotonic behavior of the efficiency on Ω_z. At g₁ = g₂ = g (the right panel and the red curve in the left panel) we find that the efficiency goes to zero in the limits Ω_z → 0 and Ω_z → ∞ having maximum at some finite value of the Zeeman splitting (notice that both these limits correspond to the integrable cases [39]). Such behavior may be explained by the fact that the case of all equal nonlinearities corresponds to the SU(2) symmetric (also known as Manakov [40]) nonlinearity. Thus at Ω_z = 0 the linear Hamiltonian (2) can be diagonalized by the global rotation. An approximate diagonalization can be also performed in the formal limit Ω_z → ∞ at Ω_x fixed, which after rescaling can be viewed as the limit of Ω_x → 0 at Ω_z fixed. Next we recall, that for a 1D NLS equation the matching conditions cannot be satisfied, this meaning that no FWM processes involving eigenstates from the same branch can be observed. On the other hand, the SU(2) nonlinearity does not support FWM processes where states from the different branches are involved. This can be easily seen if we rewrite the Manakov nonlinearity in the form G(Ψ) ≡ Ψ^†Ψ. Indeed, since the states belonging to different branches of the spectrum are mutually orthogonal, this nonlinearity does not support mixing of mentioned states, i.e. the respective frequency conversion processes is not phase matched. On the other hand, when the nonlinear coefficients are not exactly equal (blue and green lines in the left panel of figure 8) we observe that even at Ω_z → 0 the efficiency does not vanish, although becomes very small, as shown in the inset in the left panel of figure 8.

In our study we analyzed the four-wave mixing process in Bose–Einstein condensates with spin–orbit coupling. We found all phase matched configurations for degenerate case where two identical initial states interact with two probe ones. We performed numerical simulations to illustrate the dynamics in which we seeded one of the probe and observed stimulated growth of the latter combined with resonant generation of extra waves. We found unique conditions where both probe waves have smaller group velocity than pump wave, and also reported the case when two FWM process can occur simultaneously.

The work was supported by the Polish National Science Center 2016/22/M/ST2/00261 (M.T.), Vietnam National Foundation for Science and Technology Development (NAFOSTED) under Grant number 103.01-2017.55 (N.V.H.), and Portuguese Foundation for Science and Technology (FCT) under Contract no. UIDB/00618/2020 (V.V.K).