Thermalization and Bose–Einstein condensation of quantum light in bulk nonlinear media
We study the thermalization and the Bose–Einstein condensation of a paraxial, spectrally narrow beam of quantum light propagating in a lossless bulk Kerr medium. The spatiotemporal evolution of the quantum optical field is ruled by a Heisenberg equation analogous to the quantum nonlinear Schrödinger equation of dilute atomic Bose gases. Correspondingly, in the weak-nonlinearity regime, the phase-space density evolves according to the Boltzmann equation. Expressions for the thermalization time and for the temperature and the chemical potential of the eventual Bose–Einstein distribution are found. After discussing experimental issues, we introduce an optical setup allowing the evaporative cooling of a guided beam of light towards Bose–Einstein condensation. This might serve as a novel source of coherent light.
In the last few years, many-body physics has embraced a novel class of systems, the so-called quantum fluids of light . In these optical systems, light and matter combine to generate new photonlike particles that, differently from vacuum photons, are characterized by sizeable effective masses and mutual interactions and, therefore, may give rise to novel states of matter.
One of the most used platforms to study the physics of quantum fluids of light is the semiconductor planar microcavity, in which the cavity photons and the quantum-well excitons strongly couple to form mixed light-matter interacting bosonic quasiparticles called exciton polaritons . Numerous quantum-hydrodynamics collective phenomena have been investigated theoretically and successfully observed experimentally in such exciton-polariton fluids . Nevertheless, fluids of light in cavity-based systems are inevitably subject to losses, which is typically detrimental for the experimental observation of coherent quantum dynamical features. A more promising configuration for the study of quantum phenomena in fluids of light consists in the paraxial propagation of a quasimonochromatic beam of light in a nonabsorbing bulk nonlinear medium of Kerr type.
It is well known [3, 4, 5] that in such a cavityless, propagating, geometry the complex amplitude of the classical optical field is a slowly varying function of space and time which satisfies a nonlinear wave equation formally identical to the Gross–Pitaevskii (GP) equation of dilute Bose–Einstein (BE) condensates  after exchanging the roles of the propagation coordinate and of the time parameter. This classical paraxial bulk dynamics may be regarded as the emerging mean-field description of an underlying quantum nonlinear Schrödinger dynamics, as formalized in full generality in a recent work by two of us .
In a recent experimental study , C. Sun et al. have provided the first observation of classical-wave condensation using a beam of classical monochromatic light propagating in a nonlinear photorefractive crystal. The mechanism underlying this condensation of classical light finds its origin in the thermalization of the classical optical field [9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19] towards an equilibrium state whose statistics obeys the Rayleigh–Jeans (RJ) thermal law, which corresponds to the classical (high-temperature and/or long-wavelength) limit of the BE distribution.
In this letter, we push this research line forward by investigating the very quantum aspects of the thermalization dynamics of the propagating fluid of light. Making use of the fully quantum theory developed in ref. , we discuss the possibility of measuring the Boltzmann tails of the eventual BE distribution, which constitutes the hallmark of the particlelike, quantum, nature of the paraxial beam of light at thermal equilibrium. Inspired by recent advances towards atom-laser devices based on in-waveguide evaporative-cooling schemes [20, 21, 22, 23], we finally propose a mechanism leading to a complete BE condensation in the quantum fluid of light. If realized, such a process would offer a novel route to generate spontaneous optical coherence in a novel concept of coherent-light source.
We consider the propagation in the positive- direction of a paraxial, spectrally narrow beam of light of central angular frequency in a bulk, electrically neutral, nonmagnetic, nonabsorbing, nonlinear medium of real-valued intensity-dependent refractive index . Here, is the background refractive index, describes the spatial profile of the refractive index, quantifies the strength of the —spatially local and instantaneous— Kerr nonlinearity of the medium and is the slowly varying [3, 4, 5] envelope of the light wave’s electric field of propagation constant in the increasing- direction, where denotes the vacuum speed of light. For simplicity’s sake, we neglect light polarization and we assume that Raman and Brillouin light-scattering processes on phonons in the optical medium occur at a negligible rate.
Following ref. , it is possible to map the quantum propagation of the beam of light in the positive- direction onto a quantum nonlinear Schrödinger evolution of a closed system of many interacting photons in a three-dimensional space spanned by the two-dimensional transverse position vector and by the physical time parameter . Introducing the time parameter and the three-dimensional position vector , where denotes the inverse of the group velocity of the photons in the medium at , the quantum mechanical propagation equation of the light beam may be reformulated in the Heisenberg form , where the quantum field operator is the second-quantized slowly varying envelope of the electric field, normalized ( is the vacuum permittivity) in a way to satisfy the usual equal- Bose commutation relations and , and where
is the many-body Hamiltonian operator of the system.
In eq. ( ‣ Quantum formalism), is the external potential experienced by the photons, due to the spatial variation of the refractive index, and is the strength of the effective photon-photon interactions induced by the Kerr nonlinearity.
Even more importantly, and are the effective masses of the paraxial photons in, respectively, the transverse plane and the direction. In generic media, the values of are typically very different, as they have completely different physical origins: the former originates from paraxial diffraction in the transverse plane while the latter, inversely proportional to the group-velocity-dispersion parameter of the medium at , starts playing a crucial role for nonmonochromatic optical fields having a nontrivial time dependence. Unless the carrier frequency lies in the neighborhood of some optical resonance where dispersion is strong, is generally much smaller than ; as an example, using tabulated data for fused silica  around (), one obtains a ratio ().
As the Hamiltonian ( ‣ Quantum formalism) is only valid within a limited angular-frequency and wavevector range around , one has to ensure that photon-photon scattering induces no sizeable photon population outside this paraxial region. Thanks to the conservation of the energy ( ‣ Quantum formalism), a necessary and —unless the chromatic dispersion has an unusually complex shape— sufficient condition is that the two masses have the same sign. The robustness of a coherent photon wave against modulational instabilities imposes further conditions that the longitudinal mass be positive, , and the photon-photon interactions be repulsive, ; by definition, this amounts to assume that the dielectric is characterized by an anomalous group-velocity dispersion, , and a self-defocusing Kerr nonlinearity, .
In this section, we provide an analytical estimate of the time —that is, of the propagation distance along the Kerr medium— that is necessary for the isolated quantum fluid of light described by the Hamiltonian ( ‣ Quantum formalism) to thermalize. It is worth stressing that the thermalization process is here assumed to occur via photon-photon collisions within the fluid only, and not to involve any thermal equilibration with the underlying optical medium, e.g., by photon-phonon scattering or repeated absorption-emission cycles as it was instead the case in the experiment of refs. [16, 17].
Assuming for the sake of simplicity that the dielectric is spatially homogeneous, , i.e., in eq. ( ‣ Quantum formalism), and that the total interaction energy is small with respect to the total kinetic one in the eventual thermal-equilibrium state, the latter has to be characterized by an occupation number in the plane-wave state of wavevector and energy of the BE form
( is the Boltzmann constant), where and are respectively the temperature and the chemical potential of the thermalized quantum fluid of light. As we have assumed there is no thermal contact with the underlying optical medium, is not related to the temperature of the latter as in refs. [16, 17], and both and are fully determined as functions of the energy and number densities of the photon fluid entering the medium, as detailed in the next section.
A simple model —based on the quantum nonlinear Schrödinger formalism ( ‣ Quantum formalism)— to investigate the relaxation dynamics of the initial state of the photon fluid, at (i.e., ), towards thermal equilibrium, at (i.e., ), is provided by the homogeneous [as ] Boltzmann kinetic equation 
for the uniform phase-space density of the paraxial photons occupying the plane-wave state of wavevector and energy at the propagation time . At long times, i.e., when , the solution of eq. ( ‣ Thermalization time) approaches the stationary BE distribution ( ‣ Thermalization time). Equation ( ‣ Thermalization time) is valid (i) in the absence of condensate and (ii) in the weak-interaction regime. The constraint (i) is satisfied as long as one considers energies and densities yielding noncondensed equilibrium states; otherwise, one has to include the coherent dynamics of the condensate’s order parameter in eq. ( ‣ Thermalization time) . The condition (ii) may be checked a posteriori by requiring that, in the eventual thermal state, the total interaction energy is small compared to the total kinetic energy, as already supposed in the second paragraph of the present section.
To estimate the effective relaxation time towards thermal equilibrium, we are going to mutuate well-known results from the theory of weakly interacting atomic Bose gases. A numerical study  demonstrated that the thermalization time of weakly interacting bosonic atoms not too far from thermal equilibrium is typically of the order of , where denotes the average collision rate. This means that about three collisions per particle are sufficient to make the system thermalize. This collision rate may be expressed as  , where denotes the mean number density of the gas, is the average norm of the velocity of ideal classical bosons at temperature and is the low-energy boson-boson-scattering cross section. This expression for holds in the case of an isotropic three-dimensional system. In the present optical case, as highlighted in the previous section, the system is characterized by an anisotropic mass tensor. As a result, the above-given formula for cannot be applied directly to estimate the time for the quantum fluid of light to relax towards thermal equilibrium.
In order to be able to safely use it, one has to make the kinetic contribution to the Hamiltonian ( ‣ Quantum formalism) isotropic with a common mass in all the , , directions. To do so, we introduce the mass parameter —that corresponds to the geometric mean of the paraxial-photon effective masses in the transverse , and longitudinal directions— and the rescaled position vector . As an inversed rescaling holds in momentum space, one readily verifies that such a transformation preserves the spatial as well as the phase-space densities. In the isotropic space, we are then allowed to use the estimate for the thermalization time in terms of the mean number density , the Boltzmann-averaged velocity and the scattering cross section , where is the -wave scattering length written as a function of the two-body interaction parameter in the isotropic space . As our coordinate change preserves both the spatial and the phase-space densities, it is immediate to check that is equal to the original photon-photon coupling constant in the anisotropic space, .
Combining the results of the previous paragraph, one eventually gets an explicit formula for the thermalization time ,
that corresponds to the usual expression of the thermalization time of a three-dimensional weakly interacting atomic Bose gas, with a mass given by the geometric average of the masses in the transverse , and longitudinal directions.
Temperature and chemical potential at thermal equilibrium
As the kinetic-energy density and the photon number density are quantities conserved during the evolution of the quantum fluid of light described by eq. ( ‣ Thermalization time), the temperature and the chemical potential characterizing the thermal-equilibrium, at , BE distribution ( ‣ Thermalization time) may be fixed by the initial, at , values of and : and , where the left-hand sides depend on and and the right-hand sides are functions of the parameters of the incoming electromagnetic field.
In the final equilibrium state (, i.e., ), using eq. ( ‣ Thermalization time), one readily gets and , where is the fugacity and refers to the Bose integral, with the Euler gamma function. These equations are similar to the well-known results of the ideal Bose gas, with the difference that, in the present optical case, there are two different thermal de Broglie wavelengths due to the anisotropy of the kinetic energy in eq. ( ‣ Quantum formalism). By means of the equation for , one finds that , where is the mean-field interaction-energy density  of the fluid of light at equilibrium. Thus, as the Bose integrals and are of the same order [ for ], the weak-interaction condition required for eq. ( ‣ Thermalization time) to be valid reads , which may be reexpressed in terms of the -wave scattering length in the isotropic space as . Note that this constraint directly implies the usual diluteness condition when one enters the quantum-degeneracy regime, .
We assume that the initial fluid of light (), i.e., the incident beam of light (), is characterized by the following Gaussian distribution in real space:
with finite correlation lengths and in the transverse plane and the longitudinal direction, respectively. In an actual experiment, the input density of the quantum fluid of light is tuned by varying the intensity of the incoming light beam, the transverse correlation length may be tuned by processing the input beam through spatial light modulators  and the longitudinal correlation length may in principle be varied by modifying the coherence time of the incident beam. Note that must be larger than the wavelength of the carrier wave to ensure the paraxiality of the beam of light in the medium and must be smaller than the frequency range within which the quadratic approximation of the dispersion relation of the medium is valid. Fourier transforming eq. ( ‣ Temperature and chemical potential at thermal equilibrium) yields the expression of the initial occupation number at . From this, one obtains, at , and .
Making use of the conservation laws and , one eventually gets the following 2-by-2 system:
the resolution of which makes it possible to obtain and in the final thermal-equilibrium state in terms of , , , and . Introducing the effective temperatures , the first of eqs. ( ‣ Temperature and chemical potential at thermal equilibrium) may be rewritten as , which shows that the transverse and longitudinal modes, initially distributed at different temperatures , eventually equilibrate at the same temperature .
Reminding the definition of the spatial coordinate , the third component of the paraxial-photon wavevector may be expressed as [7, 27] , where is the detuning from the angular frequency of the pump. As a result, the measurement of the BE distribution ( ‣ Thermalization time) as a function of requires a good angular resolution to isolate the light deflected with a transverse wavevector as well as a good spectral resolution to isolate the angular-frequency component of the transmitted light at .
On the other hand, to have access to the large-momentum, Boltzmann, tails of the BE distribution —and so, in turn, to the whole BE distribution as a function of — at the exit face of the nonlinear dielectric where the fluid of light is imaged, some conditions have to be satisfied.
The inverse of the de Broglie wavelengths and being the natural scales of variation of as a function of and , a first condition for detecting the whole BE distribution in the transmitted beam of light is that and must verify the constraints satisfied respectively by and (see the third paragraph of the previous section).
A second, perhaps more challenging, condition concerns the length of the bulk nonlinear medium, which has to be at least of the order of the distance necessary for the quantum fluid of light to fully relax towards thermal equilibrium. Making use of the analytical result ( ‣ Thermalization time) and of the first of eqs. ( ‣ Temperature and chemical potential at thermal equilibrium) with the reasonable estimate for , one finds that must behave at a given carrier wave at as
where depends on , , , and on , , . As a most important contribution, it is immediate to see that the stronger the Kerr nonlinearity is, the shorter the thermalization distance is. Plugging explicit values into ( ‣ Experimental considerations), we estimate for a light beam of wavelength, intensity and initial coherence lengths propagating in bulk silica  an unreasonably long , i.e., of the order of the estimated radius of the solar system…
While an experiment using such standard bulk nonlinear media looks clearly unfeasable, very promising alternatives are offered by resonant media where photons are strongly mixed with matter excitations. In this way, very strong effective photon-photon interactions may be obtained, e.g., for polaritons in bulk semiconducting materials showing narrow exciton lines such as GaAs or ZnSe . This effect can be further reinforced by many orders of magnitude if the chosen material excitation involves spatially wide (even almost micron-sized) Rydberg states, either in optically dressed atomic gases in the so-called Rydberg-EIT regime  or in highest-quality solid-state samples . A further advantage of resonant media is the wide tunability of the optical parameters simply by changing the carrier frequency , which is of a great utility to ensure the dynamical stability of the photon fluid.
Discussion of a recent experiment
In ref. , C. Sun et al. reported having experimentally observed the relaxation of a classical, i.e., not quantum, fluid of interacting photons towards a thermal-equilibrium state. A beam of classical monochromatic light, initially prepared in a nonthermal state via a suitable tayloring of the incident phase profile, was made to propagate in a photorefractive crystal whose optical nonlinearity was strong enough to make the transverse angular distribution of the beam of light fastly evolve towards a RJ-type, i.e., classical, thermal law. For small enough initial kinetic energies, a marked peak around was observed in the transverse-momentum- distribution, which was interpreted as a signature of the occurrence of a kinetic condensation of classical waves.
In order to fully understand the analogies and the differences with our quantum study, we can start by noting that a key conceptual assumption of the experiment  is that the light beam remains perfectly monochromatic all along its propagation across the nonlinear crystal. Under a mean-field approximation and provided no spontaneous temporal modulations such as self-pulsing  occur, monochromaticity at all distances is a trivial consequence of the classical GP form of the nonlinear Schrödinger field equation corresponding to the quantum Hamiltonian ( ‣ Quantum formalism).
On the other hand, monochromaticity corresponds within the framework of our quantum theory to having at all propagation times a factorized momentum distribution , where the transverse-momentum distribution evolves with while the longitudinal one remains constant and proportional to the Dirac function at all ’s. Monochromaticity at all ’s then requires that no scattering process can change the ’s of the colliding paraxial photons.
Most remarkably, the specific form of the optical nonlinearity of the photorefractive crystal used in the experiment  automatically serves this purpose, as its slow response involves the time- average of the optical intensity and —in many-body terms— corresponds to infinite-range interactions along the axis. As a result, all processes that would generate frequencies different from the incident one are suppressed. Keeping in mind that the population is sharply peaked on the only occupied states with , it is then straightforward to see that the kinetics will eventually relax to the classical RJ distribution rather than to the BE one ( ‣ Thermalization time): because of the -shaped factor in the ’s, all the quantum “” terms in the Boltzmann equation ( ‣ Thermalization time) are in fact irrelevant, so that the quantum kinetics reduces to a classical one.
The situation is of course completely different if a local and instantaneous nonlinearity is used in an experiment. Within our theory , this corresponds to a local interaction in the three-dimensional , , space. As a result, wave-mixing processes can mix all the three components of the momentum, therefore allowing for a full three-dimensional thermalization of the photon gas in both its transverse-momentum- distribution and its physical-frequency- distribution, where is measured from the carrier wave at . Given the quantum nature of our model, the final result of this thermalization process will be a BE distribution of the form ( ‣ Thermalization time), which automatically solves all the ultraviolet black-body catastrophes that infest classical theories such as the one used in ref. . As a final point, it is worth highlighting that thermalization to a quantum distribution is based on the quantum “” terms in the Boltzmann equation and thus does not benefit from the large Bose stimulation factor involved in the thermalization of classical waves. Together with the typically weaker Kerr optical nonlinearity of fast media, this explains why our prediction for is dramatically longer than the experimental one of ref. .
Evaporative cooling and BE condensation of a beam of light
An interesting consequence of the above-investigated thermalization process appears when the quantum fluid of light enters the BE-condensed phase. From the theory of the ideal Bose gas , the critical line for BE condensation in the plane may be obtained by imposing in the second of the thermal-state equations ( ‣ Temperature and chemical potential at thermal equilibrium), which yields the usual formula for the BE-condensation critical temperature,
in terms of the Riemann zeta function at , , and the before-introduced geometric mean of the paraxial-photon effective masses. To realize a BE condensate of light in a bulk geometry, the experimentalist has to choose the rescaled intensity and the correlation lengths and of the incident beam in such a way that the temperature in the thermal state, solution of eqs. ( ‣ Temperature and chemical potential at thermal equilibrium), is smaller than given by eq. ( ‣ Evaporative cooling and BE condensation of a beam of light).
Following the theoretical and experimental investigations [20, 21, 22, 23] of the evaporative cooling of an atomic beam propagating in a magnetic trap, a promising way to facilitate BE condensation in the quantum fluid of light consists in progressively making the photon beam evaporate in the transverse directions.
This can be obtained by introducing a time-dependent trapping potential into the Hamiltonian ( ‣ Quantum formalism), for instance of truncated harmonic form for and for , where the radius is a decreasing function of the propagation time . Based on the relation between the spatial profile of the refractive index and the effective potential in eq. ( ‣ Quantum formalism), this truncated harmonic trap may be realized by means of a conically tapered multimode optical waveguide, as pictorially sketched in fig. 1. The core  is taken to have an inverse-parabolic refractive-index profile, while the cladding  is homogeneous with a refractive index smoothly connecting the one of the core’s edge.
As a result of this tapering, the maximum value of the trapping potential decreases as increases, so that the large-momentum (or large-energy) tails of the photon distribution are progressively removed. At the same time, the remaining photons keep reequilibrating to lower and lower temperatures under the effect of collisions, until the fluid of light eventually crosses the critical temperature for BE condensation.
Upon the mapping, BE condensation from an initially thermal photon gas corresponds to the appearance of spontaneous optical coherence when an initially incoherent beam of light is injected into the nonlinear medium: the long-range order of the BE condensate of light reflects into optical coherence extending for macroscopically long times and distances , . In contrast to trivial angular- and frequency-filtering processes, a key element of our proposal are the collisions between the photons, that allow the fluid of light to reestablish thermal equilibrium at lower and lower temperatures while the most energetic photons keep being removed.
In this letter, we have investigated the relaxation dynamics of a paraxial, quasimonochromatic beam of quantum light towards thermal equilibrium in a lossless bulk Kerr medium. Following ref. , the propagation of the quantum light field has been mapped onto a quantum nonlinear Schrödinger evolution of a conservative quantum fluid of many interacting bosons. Correspondingly, in the weak-interaction regime, the evolution of the momentum distribution from an arbitrary nonthermal state towards a thermal state with a BE form can be modeled by the Boltzmann kinetic equation, which offers analytical formulas for the thermalization time and for the final temperature and chemical potential in terms of the parameters of the input beam and of the medium.
In addition to extending the concept of classical-light-wave condensation  to a fully quantum level and solving well-known ultraviolet pathologies of existing classical theories, our results suggest an intriguing long-term application as a novel source of coherent light: taking inspiration from related advances in atom-laser devices [20, 21, 22, 23], we have pointed out a novel in-waveguide evaporative-cooling scheme to obtain spontaneous macroscopic optical coherence from an initially incoherent beam of light. As our proposal does not rely on population-inverted atomic transitions, it holds the promise of being implemented in an arbitrary domain of the electromagnetic spectrum.
We are grateful to G. Ferrari, A. Gambassi, A. Minguzzi, C. Miniatura and M. Richard for helpful and stimulating discussions. AC acknowledges the kind hospitality of the BEC Center, where this work was partially done. P-ÉL and IC were funded by the EU-FET Proactive Grant AQuS, Project No. 640800, and by the Provincia Autonoma di Trento, partly through the Call “Grandi Progetti 2012”, Project SiQuro.
-  Carusotto I. Ciuti C., Rev. Mod. Phys., 85 (2013) 299
-  Deng H. et al., Rev. Mod. Phys., 82 (2010) 1489
-  Boyd R. P., Nonlinear Optics (Academic Press, San Diego) 1992
-  Agrawal G. P., Nonlinear Fiber Optics (Academic Press, San Diego) 1995
-  Rosanov N. N., Spatial Hysteresis and Optical Patterns (Springer, Berlin) 2002
-  Pitaevskii L. P. Stringari S., Bose–Einstein Condensation and Superfluidity (Oxford Science Publications, Oxford) 2016
-  Larré P.-É. Carusotto I., Phys. Rev. A, 92 (2015) 043802
-  Sun C. et al., Nat. Phys., 8 (2012) 470
-  Picozzi A., Opt. Express, 15 (2007) 9063
-  Lagrange S. et al., EPL, 79 (2007) 64001
-  Picozzi A., Opt. Express, 16 (2008) 17171
-  Picozzi A. Rica S., EPL, 84 (2008) 34004
-  Barviau B. et al., Opt. Lett., 33 (2008) 2833
-  Barviau B. et al., Opt. Express, 17 (2009) 7392
-  Suret P. et al., Phys. Rev. Lett., 104 (2010) 054101
-  Klaers J. et al., Nat. Phys., 6 (2010) 512
-  Klaers J. et al., Nature, 468 (2010) 545
-  Aschieri P. et al., Phys. Rev. A, 83 (2011) 033838
-  Michel C. et al., Phys. Rev. A, 84 (2011) 033848
-  Mandonnet E. et al., Eur. Phys. J. D, 10 (2000) 9
-  Castin Y. et al., J. Mod. Opt., 47 (2000) 2671
-  Lahaye T. et al., Phys. Rev. Lett., 93 (2004) 093003
-  Lahaye T. et al., Phys. Rev. A, 72 (2005) 033411
-  Malitson I. H., J. Opt. Soc. Am., 55 (1965) 1205
-  Griffin A. et al., Bose-Condensed Gases at Finite Temperature (Cambridge University Press, Cambridge) 2009
-  Wu H. Foot C. J., J. Phys. B: At. Mol. Opt. Phys., 29 (1996) L321
-  Larré P.-É. Carusotto I., Eur. Phys. J. D, 70 (2016)
-  Peyronel T. et al., Nature, 488 (2012) 57
-  Jang J. I., Optoelectronics: Materials and Techniques (InTech) 2011, Chap. 5