## Abstract

The nonlinear coupling term in the Gross-Pitiaevski equation which describes a
Bose-Einstein condensate (BEC) can cause four-wave mixing (4WM) if three BEC
wavepackets with momenta **k**
_{1}, **k**
_{2},
and **k**
_{3} interact. The interaction will produce a fourth
wavepacket with momentum **k**
_{4} = **k**
_{1}
+ **k**
_{2} - **k**
_{3}. We study
this process using numerical models and suggest that experiments are feasible.
Conservation of energy and momentum have different consequences for 4WM with
massive particles than in the nonlinear optics case because of the different
energy-momentum dispersion relations.

© 1998 Optical Society of America

## 1. Introduction

Interference of matter waves formed from Bose-Einstein condensates (BECs) [1, 2, 3, 4, 5] has been demonstrated experimentally [1, 2]. By virtue of the nonlinear nature of the self-interaction
term in the Gross-Pitaevskii(GP) equation, which describes the dynamics of such
systems at zero temperature, one may expect nonlinear phenomena to occur in BEC
dynamics. The equivalent of the self-focusing nonlinearity in optical Kerr media [6, 7, 8] is actually self-defocusing for the case of positive
scattering length. This example of nonlinear behavior has been observed in the
expansion of the condensate due to the mean field energy when the trap is released [9]. Goldstein *et al*. [10, 11] have proposed that phase conjugation of matter waves should
be possible in analogy to this phenomenon in nonlinear optics, including the case of
multiple spin-component condensates [12]. They consider the case where a
“probe” BEC wavepacket interacts with two counterpropagating
“pump” wavepackets to generate a fourth that is phase
conjugate to the “probe”; the “probe” is
weak and causes negligible depletion of the “pump”. Law et al. [13] also suggest analogies between interactions in multiple
spin-component condensates and four-wave mixing. Here we examine the four-wave
mixing (4WM) in a single spin-component condensate that occurs as a result of the
nonlinear self-interaction term in the GP equation when three BEC wavepackets with
momenta **k**
_{1}, **k**
_{2}, and
**k**
_{3} collide and interact. Nonlinear 4WM can generate a
new BEC wavepacket with a new momentum
**k**
_{1}+**k**
_{2}-**k**
_{3}.
Our assumptions on geometry and number of atoms in the wavepackets are less
restrictive than those of Goldstein et al. We suggest that experiments with such
wavepackets should be feasible, for example, using Raman output coupling techniques
which have been demonstrated experimentally by the NIST group [14, 15]

## 2. Theory of four-wave mixing

The nature of 4WM in BEC collisions is unlike 4WM for optical pulse collisions in
dispersive media [8, 16, 17], since the momentum and energy constraints imposed are
different in the two cases. This is because the energy-momentum dispersion relation
for massive particles is quadratic in *k*, whereas it is linear in
*k* for the case of light. Moreover, in dispersive optical media,
the momentum of light waves is proportional to the product of the frequency of the
light and the refractive index, and the refractive index depends upon frequency (and
the propagation direction in non-isotropic media -hence conservation of energy does
not in general guarantee conservation of momentum in optical 4WM). This complication
involving the properties of the medium does not arise in the BEC case. For 4WM in
BEC collisions, when the magnitude of the momenta of each of the wavepackets are
identical, i.e., |**k**
_{i}| = |**k**
_{j}| for *i* ≠ *j*,
conservation of momentum and energy. However, in general, when |**k**
_{i}| ≠ |**k**
_{j}| for *i* ≠ *j*,
conservation of momentum does not imply conservation of energy for 4WM in BEC
collisions. Clearly, creation of new BEC wavepackets in 4WM is limited to cases when
momentum and energy conservation are simultaneously satisfied.

Here we will study 4WM of BECs by numerical calculations on a model with three BEC
wavepackets. There are two possible choices of initial conditions: (1) a
“whole collision” in which initial spatially separated
wavepackets come together at the same time, or (2) a “half
collision” in which the wavepackets are initially formed in the same
condensate at (nearly) the same time. Although we will assume initial condition (1),
similar to optical 4WM experiments with light pulses, the “half
collision” version (2) will be feasible exeprimentally using the methods
of [14, 15]; similar conclusions regarding 4WM will apply to such a
case. We assume the initial condensates are comprised of magnetically confined atoms
in the same *F, M*_{F}
state, and any trapping potentials are
turned off before propagation begins. The three condensate wavepackets are given
initial momenta and positions as in Fig. 1 so that they collide at a given point. We have carried
out calculations in 1D, 2D, and 3D.

The Gross-Pitaevskii(GP) equation for a single component BEC can be written as [18],

where ${T}_{x}=\frac{-{\mathit{\u0127}}^{2}}{2m}\left(\frac{{\partial}^{2}}{\partial {x}^{2}}+\frac{{\partial}^{2}}{\partial {y}^{2}}+\frac{{\partial}^{2}}{\partial {z}^{2}}\right)$ is the kinetic energy operator, *V*(**x**,
*t*) is the potential imposed on the atoms and ${U}_{0}=\frac{4\pi {a}_{0}{\mathit{\u0127}}^{2}}{m}{N}_{T}$ is the atom-atom interaction strength, proportional to the ** s**-wave scattering length,

*a*

_{0}, atomic mass,

*m*, and the total number of atoms in all the wavepackets,

*N*

_{T}. The initial wavefunction is comprised of three BEC wavepackets,

where *Ψ*(**x** - **x**
_{i}) is the solution to the GP equation with a locally harmonic potential
centered around **x** = **x**
_{i}, *i* = 1,2,3; the normalization constant N chosen so that the
norm of Ψ is unity. We assume the three inital positions **x**
_{i} are spatially separated so the initial wavepackets are non-overlapping (one
could also consider the “half-collision” case where the three
**x**
_{i} are the same and the wavepackets are generated *in situ* from
the same initial condensate). Although the initial wavepackets can have arbitrary
phases multiplying the amplitudes *Ψ*(**x** -
**x**
_{i}) [3, 5], for simplicity we take the initial relative phase between
wavepackets to be zero, since 4WM will occur for an arbitrary set of initial
relative phases.

The nonlinear term in the GP equation will have terms with the factor
exp[*i*(**k**
_{i} + **k**
_{j} - **k**
_{l}) · **x**] where *i, j* and
*l* can be 1, 2 or 3 respectively. These terms can generate a
wavepacket with a central momentum that is not in the initial wavefunction
Ψ(**x**, *t* = 0). For example, if
**k**
_{2} = -**k**
_{1} (see Fig. (2)), then it is possible to produce a wavepacket with
central momentum **k**
_{4} = **k**
_{1} -
**k**
_{2} + **k**
_{3} =
2**k**
_{1} + **k**
_{3}.

We set up a numerical calculation where the initial state, Eq. (2), evolves according to the GP equation (1). The time evolution is carried out using a split-operator
Fourier transform method [19, 20]. We have verfied numerically that energy and momentum are
conserved during the course of our calculation. Thus,
*dE*(*t*)/*dt* = 0, where $E\left(t\right)=\u3008\Psi \left(t\right)\mid \left({T}_{x}+\frac{1}{2}{U}_{0}{\mid \Psi \mid}^{2}\right)\mid \Psi \left(t\right)\u3009$ is the energy per particle in the BEC, and
*d*
**P**(*t*)/*dt* = 0,
where **P**(*t*) =
-*iħ*⟨Ψ(*t*)|∇|Φ(*t*)⟩
is the momentum per particle.

In order to estimate the importance of the various terms in the GP equation, it can
be written in terms of characteristic time scales
*t*_{DF}*,t*_{NL}
in the following manner [19, 20]:

Here the diffraction time *t*_{DF}
=
2${\mathit{\text{mw}}}_{0}^{2}$/*ħ*, and
the nonlinear interaction time *t*_{NL}
=
(*U*
_{0}|Ψ_{m}|^{2}/*ħ*)^{-l}, where ${\mid {\Psi}_{m}\mid}^{2}={\left(\frac{4}{3}\pi {w}_{0}^{3}\right)}^{-1}$ approximates the mean value of
|Φ(**x**)|^{2}, and
*w*
_{0} stands for initial halfwidth of the colliding
wavepackets. The smaller the characteristic time, the more important the
corresponding term in the GP equation. We also define the collision duration time
*t*_{col}
=
(2*w*
_{0})/*v*, where *v* =
*k*
_{1}/*m* is the initial group velocity
of a wavepacket. The ratio
*t*_{col}
/*t*_{NL}
gives an indication
of the strength of the nonlinearity during the collision. The larger the ratio of
*t*_{col}
/*t*_{NL}
, the
stronger the effects of the nonlinearity during the overlap of the wavepackets.
These characteristic times stand in the ratios ${t}_{\mathit{DF}}:{t}_{\mathit{col}}:{t}_{\mathit{NL}}=1:\frac{\lambda}{2\pi {w}_{0}}:\frac{{w}_{0}}{6{a}_{0}{N}_{T}}$, where *λ* is the De Broglie wavelength
associated with the wavepacket velocity *v*. Experimental condensates
with *t*_{col}
/*t*_{NL}
≫ 1
can be readily achieved. Thus, the nonlinear term will have time to act while the
BEC wavepackets remain physically overlapped during a collision.

## 3. Numerical simulations

We solve Eq. (3) numerically in reduced form by chosing the units of length
*x*
_{0} and time *t*
_{0} so that
(*t*
_{0}/*t*_{DF}
)(*w*
_{0}/*x*
_{0})^{2}
= 1/2; once *x*
_{0} is chosen, *t*
_{0}
is given by *t*
_{0} =
${\mathit{\text{mx}}}_{0}^{2}$/*ħ*. Here
we choose the unit of length *x*
_{0} to be
*x*
_{0} = 10*μ*m, so that
*t*
_{0} = 36.2 ms for ^{23}Na atoms.
Consequently, the unit of energy *E*
_{0} =
*m*(*x*
_{0}/*t*
_{0})^{2}
= *ħ*/*t*
_{0} is
*E*
_{0} = 2.91 × 10^{-33} J =
*h*(4.39 Hz), and the unit of momentum
*p*
_{0} =
*m*(*x*
_{0}/*t*
_{0})
= *ħ*/*x*
_{0} is
*p*
_{0} = 1-05 × 10^{-29} Kg m/s. For
comparison purposes, the recoil energy and momentum for a 589 nm photon (the Na
resonance transition) are 5690*E*
_{0} and
107*p*
_{0} respectively.

In 1D, only head-on collisions of the condensates are possible. With the
energy-momentum dispersion relation, *E* =
*ħ*
^{2}
*k*
^{2}
/(2*m*), the constraints imposed by conservation of energy and
momentum during the collision do not permit additional wavepackets to be created in
1D. Let us consider the following initial conditions: *x*
_{1}
= -*x*
_{2}, and *x*
_{3} = 0, and
*k*
_{2} = -*k*
_{1}, and
*k*
_{3} = 0. Two wavepackets move symmetrically towards
the central wavepacket whose center is at *x*
_{3} = 0. Since
the nonlinear term in the propagation equation is of third order in Ψ,
and Ψ is a superposition of condensates with momenta
*k*
_{1}, -*k*
_{1} and 0, the
nonlinear term could become a source of wavepackets propagating with momentum 0,
±*k*
_{1},
±2*k*
_{1} and
±3*k*
_{1}; in addition the collision could
transfer population between the condensate wavepackets. We carried out numerical
experiments with different values of nonlinearity *U*
_{0}, up
to the value of the
*τ*_{col}
/*τ*_{NL}
= 10. As expected, the only effect observed in 1D simulations was a slight delay of
the maximum of the moving wavepacket peaks. No transfer of population, or additional
peak creation was present, i.e., no wavepackets of momentum
±2*k*
_{1} and
±3*k*
_{1} were created. Even if we had taken
*k*
_{3} ≠ 0, or
|*k*
_{1}| ≠
|*k*
_{2}|, no new wavepacket would
appear in a 1D calculation.

Next we consider the two dimensional case with initial configurations such that
**x**
_{1} = -**x**
_{2} = (20,0),
**k**
_{1} = -**k**
_{2} = (10,0), and
**x**
_{3} =
(20*α*,-40*α*),
**k**
_{3} = (-10*α*,
20*α*), with several different values of the parameter
*α*; *α* = 0.5, 0.7, 1, 1.25
and 1.5. We take *N*_{T}
= 1.4 × 10^{6}
^{23}Na atoms equally partitioned between wavepackets with initial
*w*
_{0} = 10 *μ*m, for a
typical tight trap with a mean trap frequency of 200 Hz. For this case
*t*_{DF}
= 72 ms, *t*_{col}
=
7.2 ms, and *t*_{NL}
= 0.031 ms. The wavepacket momenta
**k**
_{i} are significantly larger than the internal momentum. Fig. 2a shows the initial configuration of the colliding
wavepackets. The initial conditions were selected such that the three wavepackets
collide at *t* = 0 at the origin of the reference frame. Fig. 2 shows results of our 2D calculations for the case of
*α* = 1.0. Fig. 2a–c shows the probability distribution
|ΨΦ(*x,y,t*)|^{2} before,
during, and after the collision, and Fig. 2 d–f shows its Fourier transform, the
momentum distribution |ΨΦ(*k*_{x}*,
k*_{y}*,t*)|^{2}. Fig. 2b illustrates the wavepacket interference during the
collision.

The striking feature in Figs. 2c and 2f is an additional 4WM wavepacket created in the collision
with momentum **k**
_{4} = **k**
_{1} -
**k**
_{2} + **k**
_{3}. The three
other wavepackets seen in Fig. 2c are the ones that pass through the collision region
without changing their central momenta. Redistribution of population between
different wavepackets must satisfy conservation of energy and momentum. Before or
after the collision, when the wavepackets are separated in space, the population of
the i-th one is given by the integral *N*_{i}
=
*N*_{T}
∫_{Vi}
*d*
**x**⟨Ψ_{i}|Ψ_{i}⟩ where the integration region *V*_{i}
is
selected to include the region around the i-th wavepacket. If we use
*N*_{i}
and *N′*_{i}
to
denote inital and final populations, conservation of energy and momentum show that
*N′*
_{4} =
Δ*N*
_{3} =
Δ*N*
_{1} =
-Δ*N*
_{2}, where
Δ*N*_{i}
=
*N′*_{i}
- *N*_{i}
. These
relationships are satisfied in our numerical simulations.

## 4. Interpretation

A physical interpretation of our results can be made in analogy with 4WM in nonlinear
optics. We can regard the collision of the condensates as producing a grating formed
due to the nonlinear term in Eq. (1). If Bragg conditions are satisfied (i.e., if energy and
momentum are conserved in the formation of new peaks), new 4WM induced wavepackets
may be formed in addition to the wavepackets present initially. Consider a nonlinear
term in Eq. (1) and assume that Φ(**x**,
*t*) consists of three wavepackets: Ψ(**x**,
*t*) = ${\mathrm{\sum}}_{i=1}^{3}$ Ψ_{i}(**x**, *t*), as in Eq. (2). Hence, the nonlinear term
|Ψ(**x**,*t*)|^{2}Ψ(**x**,*t*)
in Eq. (1) becomes a sum of nine contributions. Terms homogeneous in
the index i describe self phase modulation (self-focusing or actually,
self-defocusing for the case of positive scattering length), in analogy with
nonlinear optics, and these terms can not be a source of new wavepackets. Moreover,
crossed phase modulation terms of the form |Ψ_{i}(**x**,*t*)|^{2}Ψ_{j}(**x**,*t*) for *i,j* = 1,2,3 and
*j* ≠ *i* also can not contribute to
the formation of new peaks. Only mixed terms, containing different indices may be a
source of new wavepackets, but only when the Bragg conditions (i.e., conservation of
momentum and energy) are satisfied. Specifically, these conditions are: (a)
**k**
_{4} = **k**
_{i} - **k**
_{j} + **k**
_{l}, and (b) conservation of energy, which for our initial configuration (see Fig. 2) is equivilant to
|**k**
_{4}| =
|**k**
_{3}| where *i* = 3
and *j* = 1, *l* = 2. Referring again to the analogy
with nonlinear optics, one can say that wavepackets with indices *j*
and *l* create a grating during the collision whose grating vector is
**K** = -**k**
_{j} + **k**
_{l};, and the wavepacket
*i* = 3 scatters off of this grating. This is confirmed by our
numerical simulations. Only when *α* = 1 does the central
**K** vector satisfy the Bragg conditions stated above. However, since
the wavepackets contain momentum components spread around the central **K**
vector, formation of the additional wavepacket may still occur for
*α* ≠ 1 due to internal momentum
compensating the momentum mismatch of the central **K** vectors. As
|*α* - 1| increases, the number
of atoms in the additional wavepacket decreases rapidly.

It is important to note that the new 4WM wavepacket with momentum
**k**
_{4} is *not* a phase conjugate of the
wavepacket with momentum **k**
_{3}, since the wavepackets with
momenta **k**
_{1} and **k**
_{2} do not form a
static grating, but rather, form a dynamic grating which changes in time due to the
evolution and interaction of the wavepackets throught the course of the collision.
The new **k**
_{4} wavepacket is cigar-shaped, whereas the initial
**k**
_{3} wavepacket is spherical. The final
**k**
_{1} and **k**
_{3} wavepackets are mirror
images of one another (about a plane containing the centers of wavepackets 2 and 4).
Their final shape is distorted relative to their initial shape; a bite has been
removed from wavepackets 1 and 3 by the 4WM process which (1) created the new
wavepacket with momentum **k**
_{4} and (2) added Bose atoms to the
wavepacket with momentum **k**
_{2}. Fig. 3 is an enlarged view of the final wavepackets in
coordinate space with the details of their shapes shown more clearly. If we were to
view the final wavepackets in the reference frame in which the new wavepacket is
stationary, the trailing edges of wavepackets 1 and 3 would be reduced and the
*trailing* edge of wavepacket 2 would be enhanced. The structure
of the wavepackets clearly shows that the nature of the matter waves obtained after
the collision is sensitive to the details of the collision dynamics and the
properties of the initial wavepackets. A static grating picture is not sufficient to
explain the results.

Our 3D calculations show similar results to the 2D ones. The quantity

which indicates the z-averaged distribution, is similar to the 2D distribution and to
the probability distribution cut at the collision plane *z* = 0.
These different distributions show only a few percent difference in the ratio of
number of atoms in the four final wavepackets.

The formulation developed here assumes that the BECs are at zero temperature and can
therefore be treated in the mean-field approximation as given by the GP equation. At
finite temperatures, the collision must include the above-the-mean-field part of the
wavefunction (order parameter), and a density matrix treatment is required.
Moreover, even at zero temperature, there may be effects due to non-vanishing
expectation values of the above-the-mean-field part of the wavefunction involving
the creation of excitations due to the collision [21]. Moreover, in our treatment we have assumed that only one
value of *M*_{F}
is present for atoms in the condensates
(e.g., *F* = 1, *M*_{F}
= -1) with the z-axis
perpendicular to the scattering plane. Another means of carrying out 4WM experiments
is to use far off-resonance light traps to confine separate BECs [22], and impart momentum boosts to the BECs. In this case, the
BECs would be multi-component in nature with all values of
*M*_{F}
being present in the condensates. A multi-component
GP equation could be used to describe such experiments, as for example, the
phase-conjugation experiment proposed by [12].

## 5. Acknowledgments

This work was supported in part by grants from the US-Israel Binational Science Foundation, the James Franck Binational German-Israel Program in Laser-Matter Interaction (YBB) and the U.S. Office of Naval Research (PSJ). We are grateful to Eduard Merzlyakov for assisting with the 3D computations at the Israel Supercomputer Center.

## References

**1. **M. R. Andrews, C. G. Townsend, H.-J. Miesner, D. S. Durfee, D. M. Kurn, and W. Ketterle, “Observation of interference between two Bose-Einstein condensates,” Science **275**, 637 (1997). [CrossRef] [PubMed]

**2. **D. S. Hall, M. R. Matthews, C. E. Wieman, and E. A. Cornell, “Measurements of relative phase in two-component Bose-Einstein condensates,” Phys. Rev. Lett. **81**, 1543 (1998). [CrossRef]

**3. **H. Wallis, A. Rohrl, M. Naraschewski, and A. Schenzle, “Phase-space dynamics of Bose condensates: Interference versus interaction,” Phys. Rev. **A55**, 2109 (1997).

**4. **A. Rohrl, M. Naraschewski, A. Schenzle, and H. Wallis, “Transition from phase locking to the interference of independent Bose condensates: Theory versus experiment,” Phys. Rev. Lett. **78**, 4143 (1997). [CrossRef]

**5. **J. Javanainen and M. Wilkens, “Phase and phase diffusion of a split Bose-Einstein condensate,” Phys. Rev. Lett. **78**, 4675 (1997). [CrossRef]

**6. **R. Y. Chiao, E. Garmire, and C. H. Townes, “Self-trapping of optical beams,” Phys. Rev. Lett. **13**, 479 (1964). [CrossRef]

**7. **O. Svelto, “Self-focusing, self-trapping, and self-phase modulation of laser beams,” Prog. Opt. **12**, 3 (1973).

**8. **R. W. Hellwarth, “Third-order optical susceptibilities of liquids and solids,” Prog. Quant. Electr. **5**, 1 (1977). [CrossRef]

**9. **M.-O. Mewes, M R. Andrews, N. J. van Druten, D. M. Kurn, D. S. Durfee, and W. Ketterle, “Bose-Einstein Condensation in a Tightly Confining dc Magnetic Trap,” Phys. Rev. Lett. **77**, 416 (1996). [CrossRef] [PubMed]

**10. **E. Goldstein, K. Plättner, and P. Meystre, “Atomic phase conjugation,” Quantum Semiclass. Opt. **7**, 743 (1995). [CrossRef]

**11. **E. Goldstein, K. Plättner, and P. Meystre, “Atomic phase conjugation from a Bose condensate,” J. Res. Nat. Inst. Stand. Technol. **101**, 583 (1996). [CrossRef]

**12. **E. Goldstein and P. Meystre, “Phase conjugation of trapped Bose-Einstein condensates,” LANL Preprint Archive Cond-mat/9806165 (1998).

**13. **C. K. Law, H. Pu, and N. P. Bigelow, “Quantum spins mixing of spinor Bose-Einstein condensates,” LANL Preprint Archive Cond-mat/9807258 (1998).

**14. **M. Kozuma, L. Deng, E. Hagley, J. Wen, R. Lutwak, K. Helmerson, S. L. Rolston, and W. D. Phillips, “Coherent splitting of Bose-Einstein condensed atoms with optically induced Bragg diffraction,” Phys. Rev. Lett., in press (1998).

**15. **E. Hagley, L. Deng, M. Kozuma, J. Wen, K. Helmerson, S. L. Rolston, and W. D. Phillips, “A well-collimated quasi-continuous atom laser,” preprint (1998).

**16. **P. D. Maker and R. W. Terhune, “Study of optical effects due to induced polarization third order in the electric field strength,” Phys. Rev. **A137**, 801 (1965). [CrossRef]

**17. **A. Yariv and D. M. Pepper, “Amplified reflection , phase-conjugate, and oscillation in degenerate four-wave mixing,” Opt. Lett. **1**, 16 (1977). [CrossRef] [PubMed]

**18. **P. Noziéres and D. Pines, *The Theory of Quantum Liquids, Vol. II*, (Addison-Wesley, Redwood City, 1990).

**19. **M. Trippenbach and Y. B. Band, “Dynamics of short pulse splitting in dispersive nonlinear media,” Phys. Rev. **A56**, 4242 (1997)

**20. **M. Trippenbach and Y. B. Band, “Effects of Self-steepening and self-frequency shifting on short-pulse splitting in dispersive nonlinear media,” Phys. Rev. **A57**, 4791 (1998).

**21. **Y. Japha, S. Choi, K. Burnett, and Y. B. Band, “Coherent output, stimulated quantum evaporation and pair breaking in a trapped atomic Bose gas,” Phys. Rev. Lett. (in press).

**22. **J. Stenger, S. Inouye, D. M. Stamper-Kurn, H. J. Miesner, A. P. Chikkatur, and W. Ketterle, “Spin domains in ground-state Bose-Einstein condensates,” Nature **396**, 345 (1998). [CrossRef]