Nonequilibrium thermodynamics and optimal cooling of a dilute atomic gas

TL;DR

This paper develops a theoretical and experimental framework to analyze and optimize the nonequilibrium thermodynamics of a few-atom gas during cooling, focusing on entropy production and thermalization processes.

Contribution

It introduces a method to optimize cooling sequences by minimizing entropy production and verifies a refined second law in a few-particle system.

Findings

Single Raman pulse creates a nonequilibrium state that does not thermalize alone.

Combining free evolution and cooling pulses achieves thermalization.

Optimized pulse spacing minimizes entropy production and enhances cooling efficiency.

Abstract

Characterizing and optimizing thermodynamic processes far from equilibrium is a challenge. This is especially true for nanoscopic systems made of few particles. We here theoretically and experimentally investigate the nonequilibrium dynamics of a gas of few noninteracting Cesium atoms confined in a nonharmonic optical dipole trap and exposed to degenerate Raman sideband cooling pulses. We determine the axial phase-space distribution of the atoms after each Raman cooling pulse by tracing the evolution of the gas with position-resolved fluorescence imaging. We evaluate from it the entropy production and the statistical length between each cooling steps. A single Raman pulse leads to a nonequilibrium state that does not thermalize on its own, due to the absence of interparticle collisions. Thermalization may be achieved by combining free phase-space evolution and trains of cooling pulses. We minimize the entropy production to a target thermal state to specify the optimal spacing between a sequence of equally spaced pulses and achieve in this way optimal thermalization. We finally use the statistical length to verify a refined version of the second law of thermodynamics. Altogether, these findings provide a general, theoretical and experimental, framework to analyze and optimize far-from-equilibrium processes of few-particle systems.

Equations

$${\rho }_f(z,p_z)\propto \exp \left(-\frac{V(z)}{k{T}_{\text{R}}}-\frac{p_z^2}{2mk{T}_{\text{R}}}\right)$$

$$D({\rho }_i||{\rho }_f)=\int dzd{p}_z{\rho }_i(z,p_z)\ln \left(\frac{{\rho }_i(z,p_z)}{{\rho }_f(z,p_z)}\right).$$

$$\Delta S=S_f-S_0={\beta }_{f}Q+\Sigma .$$

$$L_i=\sqrt{2{\Sigma }_i}=\sqrt{2[D\left({\rho }_{i-1}||{\rho }_f\right)-D\left({\rho }_i||{\rho }_f\right)]},$$

$$\Sigma \ge \frac{L^2}{2n}$$

$$\left(\frac{\partial }{\partial t}+\frac{p_z}{m}\frac{\partial }{\partial z}+F(z)\frac{\partial }{\partial p_z}\right)\rho (z,p_z,t)=I_{\mathrm{coll}}(\rho ,p_z),$$

$$\rho (z,p_z,t)=A\exp \left[-a(t)z^2-b(t)p_z^2-c(t)z{p}_z\right].$$

$$\dot{a}$$

$$\dot{b}$$

$$\dot{c}$$

$$D\left({\rho }_i||{\rho }_f\right)$$

$$=D(f_i(z)||f_f(z)).$$

$$K({\rho }_a||{\rho }_b)=D({\rho }_a||({\rho }_a+{\rho }_b)/2).$$

$$\Sigma =D({\rho }_0||{\rho }_{\text{eq}})=\int dzd{p}_z({\rho }_0\ln {\rho }_0-{\rho }_0\ln {\rho }_{\text{eq}}).$$

$${\Sigma }_0=D({\rho }_0||{\rho }_{\text{eq}})-D({\rho }_1||{\rho }_{\text{eq}}).$$

$${\Sigma }_i=D({\rho }_i||{\rho }_{\text{eq}})-D({\rho }_{i+1}||{\rho }_{\text{eq}}).$$

$$f(z)=\int \rho (z,p_z)\mathrm{d}{p}_z$$

$$\tilde{f}(p_z)=\int \rho (z,p_z)\mathrm{d}z$$

$$D_{z{p}_z}({\rho }_1||{\rho }_2)$$

$$+\frac{2a_1{b}_2+2a_2{b}_1-c_1{c}_2}{(4a_1{b}_1-c_1^2{)}^2}$$

$$D_z(f_1||f_2)$$

$$a_2$$

$$b_2$$

$$c_2$$

$${\chi }^2=\sum _{z_i}{\left(\frac{f_{\mathrm{sim}}(z_i)-f_{\exp }(z_i)}{\Delta f_{\mathrm{sim}}(z_i)+\Delta f_{\exp }(z_i)}\right)}^2$$

$$D(f_1||f_2)$$

$$=$$

$$K(f_1||f_2)$$

$$=$$

Full text

I Introduction

Nonequilibrium processes are omnipresent in nature. Owing to their complexity and diversity, their description far away from thermal equilibrium is nontrivial leb08 . A defining property of out-of-equilibrium systems is that they dissipate energy in the form of heat, leading to an irreversible increase of their entropy. The irreversible entropy production is thus a central quantity of nonequilibrium thermodynamics, the same way that entropy is a central quantity of equilibrium physics leb08 . In the past decades, the laws of thermodynamics have been successfully extended to small nonequilibrium systems bus05 ; sek10 ; sei12 ; jar11 ; cil13 . In these systems, thermal fluctuations can no longer be neglected and thermodynamic variables are therefore random. In particular, the second law has been generalized in the form of fluctuation theorems that quantify the occurrence of negative entropy production events bus05 ; sek10 ; sei12 ; jar11 ; cil13 . The stochastic properties of the nonequilibrium entropy production have been extensively investigated, both theoretically and experimentally, for microscopic systems such as colloidal particles bus05 ; sek10 ; sei12 ; jar11 ; cil13 ; wan02 ; car04 ; sch05 ; bli06 ; dou06 ; tie06 . On the other hand, only few experiments have probed nonequilibrium thermodynamics in nanoscopic systems so far. These include one-particle systems, such as a single spin-1/2 bat14 ; bat15 or a single harmonic oscillator an15 ; bru18 , two-spin systems mic19 ; pal19 , and many-particle systems, such as cold-atomic gases bru18 ; kin06 ; gri12 ; cer17 . However, to our knowledge, no such nonequilibrium thermodynamic experiment has been realized in the intermediate regime of few-particle systems.

The laws of thermodynamics are not only of fundamental but also of practical importance. A primary objective of thermodynamics is thus to optimize processes. Optimization goals vary depending on the application, ranging from the minimization of dissipation to the maximization of work output or of cooling power bej06 . For macroscopic systems, the properties of optimal transformations have been studied within finite-time thermodynamics sal83 ; and84 ; nul85 ; and11 . The two central quantities of this approach are the entropy production, that characterizes energy dissipation, and the thermodynamic length, that measures the distance from equilibrium at which a system operates. Both are commonly calculated in the linear response regime by expanding thermodynamic potentials, such as entropy or internal energy, to second order around equilibrium sal83 ; and84 ; nul85 ; and11 . Optimization schemes are usually developed by minimizing one of the two. These techniques have been employed to optimize fractional distillation and other processes sal83 ; and84 ; nul85 ; and11 ; sal98 ; sal01 ; nul02 ; che99 . On the other hand, for microscopic systems, where thermal fluctuations are sizable, this optimization framework has been extended to the level of single trajectories within stochastic thermodynamics for linear sch07 and nonlinear aur11 systems. Methods to theoretically compute and experimentally evaluate the thermodynamic length have been proposed cro07 ; fen09 ; siv12 ; gin16 . However, despite these theoretical studies, such nonequilibrium optimization schemes have still to be demonstrated experimentally. In particular, thermodynamic distances have not been measured yet.

A further complication arises in atomic systems. A central assumption of finite-time thermodynamics and stochastic thermodynamics is indeed that systems are coupled to ideal heat baths that induce full phase-space thermalization, that is, of both position and momentum degrees of freedom. However, this hypothesis is often not fulfilled at the atomic level. A prominent instance is provided by laser cooling of atoms which plays an essential role in the study of new states of matter and high-resolution spectroscopy coh11 . Most laser cooling schemes only induce thermalization of the momentum degrees of freedom met99 . In dense atomic samples, frequent atomic collisions redistribute the energy and establish thermal equilibrium. By contrast, in dilute gases with rare interparticle collisions, these nonideal reservoirs lead to far from equilibrium states that do not thermalize on their own. Their description thus lies outside the currently existing framework. New experimental and theoretical tools are hence required to achieve their thermalization.

We here report the theoretical and experimental investigation of the nonequilibrium dynamics and the thermalization of a dilute gas of Cesium atoms confined in an optical dipole trap met99 , and illuminated by laser pulses for degenerate Raman sideband cooling (DRSC) vul98 ; ker00 ; han00 . This technique is a standard subrecoil cooling scheme for a variety of atomic systems tre01 ; web03 ; mon95 ; des04 ; lee96 ; che15 ; par15 ; gro17 ; kau12 ; tho13 . The present study of a few-particle system coupled to an engineered bath allows us to experimentally access key nonequilibrium quantities in a well-controlled atomic setup. It further gives us the opportunity to illustrate and validate our general nonequilibrium optimization approach with a common laser cooling example. We determine in particular, for each thermalization step, the nonequilibrium entropy production and a generalized thermodynamic length appropriate for these nonideal reservoirs. We use the former quantity to optimize the cooling of the few-atom system and the latter one to gain physical insight into the optimal cooling process and verify a refined version of the second law of thermodynamics known as the horse-carrot theorem and11 ; sal98 .

In our experiment, short pulses of Raman cooling lasers are applied to an initially thermal cloud along the axial direction of our nonharmonic trap. Axial and radial directions are only weakly coupled, rendering the problem essentially one-dimensional. The Raman pulses thermalize the atomic momentum distribution to the Raman temperature, thus cooling the system, but leave the position distribution unchanged. They hence create for most initial conditions a nonequilibrium state that does not thermalize on its own, due to the absence of interparticle collisions. In order to realize complete phase-space thermalization at the Raman temperature, we devise protocols consisting of a train of Raman pulses separated by intervals of free evolution (Fig. 1). For concreteness, we consider a sequence of three equally spaced pulses. The first Raman pulse (RP1) decreases the energy of the gas and moves it out of equilibrium. The second and third pulses (RP2 and RP3) drive the gas back towards a thermal state while cooling it further. For quasi harmonic trapping potentials, thermalization is routinely established by using a pulse spacing of a quarter of the trap period dep99 . This method has, for instance, recently led to the all-optical Bose-condensation of Rb atoms without an evaporative cooling stage hu17 . For the strongly anharmonic potential of our experiment, the trap leads to nontrivial dynamics of the nonequilibrium states and raises the question of the choice of the pulse spacing in this case. We seek the optimal pulse spacing $\tau$ by minimizing the entropic distance to the equilibrium target state $\rho_{f}$ at the Raman temperature, employing both a static and a dynamical criterion, which lead to the same result. An analysis of the nonequilibrium statistical length furthermore reveals that optimal thermalization is mainly reached during the first two cooling stages, with nearly equal statistical distances.

The outline of the paper is as follows. We begin in Sect. II by deriving the nonequilibrium entropy production and the statistical length for the nonideal reservoirs occurring in the experiment. We further present the horse-carrot theorem and the two criteria used to optimize the thermalization. In Sect. III we illustrate the physical meaning of the static and dynamical optimization criteria for the analytically solvable case of a harmonic trapping potential. We additionally present the experimental setup in Sect. IV and the numerical phase-space reconstruction of the phase-space distributions in Sect. V. Finally, in Sect. VI, we demonstrate optimal thermalization in a strongly nonharmonic trap and an experimental verification of the horse-carrot theorem.

II Nonequilibrium quantities and optimization criteria

As illustrated in Fig. 1, the goal of the DRSC protocol is to reach the final thermal state

[TABLE]

at the DRSC temperature $T_{\text{R}}$, where $\rho(z,p_{z})$ is the projected phase-space density onto the $(z,p_{z})$ plane, which is the relevant subspace for our experiment, and $V(z)$ denotes the axial potential. In order to quantify the approach of a nonthermal state $\rho_{i}$ produced by the DRSC scheme to the final target state $\rho_{f}$, we employ the relative entropy between these two states cov06

[TABLE]

Similarly, the corresponding quantities $D(f_{i}||f_{f})=\int dz\,f_{i}\ln(f_{i}/f_{f})$ and $D(\tilde{f}_{i}||\tilde{f}_{f})=\int dp_{z}\,\tilde{f}_{i}\ln(\tilde{f}_{i}/\tilde{f}_{f})$ can be defined for the respective position and momentum projections, $f(z)$ and $\tilde{f}(p_{z})$, of the phase-space distribution. The relative entropy is an information-theoretic quantity that satisfies the important property that $D(\rho_{i}||\rho_{f})\geq 0$, equality being only achieved for $\rho_{i}=\rho_{f}$  cov06 . This renders the relative entropy a useful indicator for the approach to the final target state.

The relative entropy also possesses a simple thermodynamic interpretation pro76 ; sch80 ; esp10 ; def11 . For a nonequilibrium process from an initial thermal state $\rho_{0}$, at inverse temperature $\beta_{0}=(kT_{0})^{-1}$, to a final thermal state $\rho_{f}$, at inverse temperature $\beta_{f}=(kT_{f})^{-1}$, the (axial) Gibbs-Shannon entropy, $S=-\int dzdp_{z}\,\rho\ln\rho$ satisfies pro76 ; sch80 ; esp10 ; def11

[TABLE]

Here, $Q=\int dzdp_{z}(\rho_{f}-\rho_{0})H$ is the heat absorbed by the system, $H$ its Hamiltonian and $\Sigma=D(\rho_{0}||\rho_{f})$ the nonequilibrium entropy production given as the relative entropy between initial and final states.

For a discrete sequence of nonthermal intermediate states $\rho_{i}$, $(i=1,2,3)$, as created after each Raman cooling pulse in our experiment, the entropy production associated with each step reads $\Sigma_{i}=D\left(\rho_{i-1}||\rho_{f}\right)-D\left(\rho_{i}||\rho_{f}\right)$ (Appendix A). The statistical length defined as,

[TABLE]

then quantifies the distance from equilibrium at which the system operates. It vanishes when $\rho_{i}=\rho_{f}$ for all $i$. Equation (4) reduces to the usual thermodynamic length in the limit of quasistatic processes where all the intermediate states are close to thermal nul85 ; and11 ; sal98 ; sal01 ; nul02 . The above nonequilibrium quantities allow the investigation of not only the final state reached after the application of the DRSC protocol but also of the cooling process itself, by providing direct information on the intermediate states.

The total entropy production $\Sigma=\sum_{i=1}^{n}\Sigma_{i}$ (multiplied by $T_{f}$) is a measure of the amount of energy that is irreversibly extracted from the system during thermalization pro76 ; sch80 ; esp10 ; def11 . It is bounded from below by the square of the total statistical length $L=\sum_{i=1}^{n}L_{i}$ divided by twice the number of steps (Appendix A)

[TABLE]

in analogy to the horse-carrot theorem and11 ; sal98 ; sal01 . The name horse-carrot process finds its origin in the analogy with a system (the horse) which is coaxed along a sequence of states by controlling the state of its environment (the carrot). The importance of the horse-carrot theorem stems from the fact that it provides a sharper lower bound to the nonequilibrium entropy production than the second law of thermodynamics which only requires $\Sigma\geq 0$. It additionally implies that optimal quasistatic horse-carrot processes (for which inequality is replaced by an equality) correspond to steps of equal thermodynamic length sal83 ; and84 ; nul85 ; and11 ; sal98 ; sal01 ; nul02 . We shall find that this also holds exactly for a harmonic confining potential and approximately for a nonharmonic trap for the generalized nonequilibrium statistical length (4) (Sect. VI).

Commonly considered optimization schemes minimize the nonequilibrium entropy production with fixed initial and final states nul85 ; and11 ; sal98 ; sal01 ; nul02 . By contrast, the state $\rho_{i}$ produced by the DRSC protocol depends on the entire cooling sequence and is hence not fixed. Our strategy is therefore to minimize the entropic distance to the target thermal state $\rho_{f}$ and identify the final temperature with the Raman temperature, $\beta_{f}=\beta_{\mathrm{R}}$. We concretely consider two optimization criteria:

(1) Static criterion: the first condition minimizes the relative entropy between $\rho_{i}$ and the target state $\rho_{f}$, $D\left(\rho_{i}||\rho_{f}\right)$. This corresponds to minimizing the entropy production $\Sigma_{i}$ [Eq. (4)] in the case of successful thermalization.

(2) Dynamical criterion: the second condition minimizes the amplitude of oscillations of the positional relative entropy, $\Delta D=\max_{t}D(f_{i}(t)||f_{f})-\min_{t}D(f_{i}(t)||f_{f})$, during the free time evolution of the atomic cloud after the Raman pulse. This criterion is based on the stationarity of a thermal state: for an equilibrium state, the distribution $f_{i}$ is constant in time and hence $\Delta D=0$. The closer the state is to equilibrium, the smaller the oscillation amplitude $\Delta D$.

The application of both optimization criteria requires to extract the relative entropy $D$ from measured data.

III Harmonic Case

In order to better understand the physical meaning of the above optimization criteria, we first consider the problem of a harmonic potential which is analytically solvable. In this case, the optimal pulse spacing is given by a quarter of the oscillation period dep99 ; hu17 . This result is intuitively clear as it corresponds to the time needed to switch position and momentum axes in phase-space during free evolution. Phase-space compression, and hence cooling and thermalization, is therefore optimal.

We analyze the phase-space dynamics by solving the Boltzmann equation for the density $\rho(z,p_{z},t)$ hua87

[TABLE]

where $m$ is the atomic mass, $F(z)=-m\omega_{z}^{2}x$ the force acting on the atom and $I_{\mathrm{coll}}$ the collision integral which takes in to account atomic interactions. For the few-atom samples that we consider, atomic interactions are negligible and hence $I_{\mathrm{coll}}=0$. Equation (6) can then be solved exactly with the Gaussian ansatz gue14

[TABLE]

This leads to a system of three coupled linear differential equations of first-order for the time-dependent coefficients $a(t)$, $b(t)$ and $c(t)$

[TABLE]

Equations (7)-(10) can be used to compute analytical expressions for the relative entropies (2) (Appendix B).

The evolution of the phase-space density $\rho(z,p_{z},t)$, together with the relative entropy $D_{zp_{z}}(t)=D(\rho(t)||\rho_{f})$ and the positional relative entropy $D_{z}(t)=D(f_{i}(t)||f_{f})$ are shown in Fig. 2a as a function of time. The phase-space distribution $\rho(z,p_{z},t)$ is circular (equilibrium) for the initial thermal state. It is elliptic (nonequilibrium) after the first Raman pulse applied at $t=0$ and rotates with period $t_{\text{osc}}=2\pi/\omega_{z}$. It is again circular (equilibrium) for the thermal state attained after the second Raman pulse applied at $t=(3/4)t_{\text{osc}}$. The relative entropy $D_{zp_{z}}(t)$ is constant during free evolution. Its value is halved after each Raman pulse (dots) until it vanishes once the target state $\rho_{f}$ is reached. This is the effect captured by the first (static) thermalization criterion. On the other hand, the positional relative entropy $D_{z}(t)$ is only constant for equilibrium states and oscillates for nonequilibrium distributions, reflecting the rotation of the phase-space density $\rho(z,p_{z},t)$. The amplitude of these oscillations vanishes once the target state $\rho_{f}$ is reached. This is the physical content of the second (dynamical) thermalization criterion. Note that we have applied the second Raman pulse at $t=(3/4)t_{\text{osc}}$ in this example only to display the oscillations of the intermediate nonequilibrium state. Optimal thermalization can already be achieved at $t=(1/4)t_{\text{osc}}$. We may further characterize the cooling process in one single diagram by combining the relevant quantities for the static and dynamical criteria, $D_{zp_{z}}$ and $\Delta D$, (Fig. 2b), in a schematic representation that qualitatively resembles Fig. 1.

For the ideal harmonic trap, full phase-space thermalization at the Raman cooling temperature $T_{\text{R}}$ is already obtained after the second DRSC pulse. The situation is more involved for nonharmonic potentials. Owing to the nonlinearity of the trapping force, each atom has a different period which depends on the oscillation amplitude. Determining the optimal pulse spacing for arbitrary nonharmonic potentials is therefore a highly nontrivial task. We will next show that our optimization strategy successfully works for arbitrary potentials.

IV Experimental Setup

We initialize our system by trapping an average of 7 Cs atoms from background vapor in a magneto-optical trap and transfer them into a crossed optical-dipole trap which creates a conservative potential (Fig. 3a). The trap is formed by a horizontal laser beam propagating along the $z$-axis with a beam waist of $21\text{\,}\mathrm{\SIUnitSymbolMicro m}$ and power of $0.25\text{\,}\mathrm{W}$, and a second crossed, vertical beam pointing in $-x$-direction with a waist of $165\text{\,}\mathrm{\SIUnitSymbolMicro m}$ and power of $3.5\text{\,}\mathrm{W}$. The atomic collision rate of $36\text{\,}\mathrm{Hz}$ at peak density is smaller than the inverse evolution time used in the experiment. The cloud is thus effectively noninteracting. We extract the atomic positions along the axial $z$-direction by employing fluorescence imaging in a 1D optical lattice sch16 and obtain the experimental position distribution $f(z)$ after binning (Fig. 3b). Every measurement is repeated several hundred times with identical parameters to get sufficient statistics. The dipole trap potential is approximately harmonic in radial direction $(x,y)$ with trap frequency $\omega_{r}=2\pi\times$1.1\text{,}\mathrm{kHz}. The initial thermal state at temperature $T_{0}$ is prepared by applying a sufficiently long optical molasses pulse [met99 ](#bib.bib38). The potential is markedly anharmonic in the axial $z$-direction and the position distribution features pronounced wings. At the center of the trap, the harmonic approximation yields an axial frequency of $\omega_{a}=2\pi\times$60\text{\,}\mathrm{Hz}. We extract the initial temperature of the gas, $T_{0}=$12.1(11)\text{,}\mathrm{\SIUnitSymbolMicro}\mathrm{K}$$, by comparing the measured position distributions $f(z)$ to numerical simulations of the three-dimensional trapping potential for atoms at various temperatures in a $\chi^{2}-$analysis (Appendix C).

We cool the initial state of the atomic cloud by applying a train of DRSC pulses following the scheme of Ref. ker00 . Details on the experiment may be found in Refs. hoh16 ; may18 . The setup comprises four DRSC lattice beams and a pump beam, as illustrated in Fig. 3a. At a detuning of $-6\text{\,}\mathrm{MHz}$ from the Cs $\mathrm{D}_{2}$-transition $\ket{F=4}\rightarrow\ket{F^{\prime}=4}$, the DRSC lattice lasers create an interference pattern with lattice sites at a trap depth of $U_{\mathrm{Raman}}=k\times$44\text{,}\mathrm{\SIUnitSymbolMicro K} and trap frequencies of $\omega_{\mathrm{trap}}=2\pi\times(71,29,28)\,$\mathrm{kHz} along the principal axis of the trap minimum. During a DRSC pulse the Cs atoms are tightly confined in a lattice site. The magnetic background field of $100\text{\,}\mathrm{m}\mathrm{G}$ applied along the $x$-$y$-diagonal is chosen such that neighboring Zeeman and vibrational states $\ket{m_{F}+1,\nu}$ and $\ket{m_{F},\nu-1}$ are energetically degenerate. A Raman coupling induced by the DRSC lasers leads to the exchange of population between these degenerate states and thereby facilitates the transfer of vibrational energy to Zeeman energy. An additional DRSC pumping beam which drives mainly $\sigma^{+}$-transitions at a detuning of $12\text{\,}\mathrm{MHz}$ to the $\ket{F=3}\rightarrow\ket{F^{\prime}=2}$-transition dissipates the Zeeman-energy, while preserving the vibrational state during the absorption and subsequent emission of the pump-photons. The Lamb-Dicke factors along the three principal axis of the Raman lattice sites are $\eta=(0.17,0.27,0.27)$. This leads on average to a reduction of the vibrational quantum number $\nu$, and thereby to a decrease of the kinetic energy of the atoms.

We apply a train of three such DRSC pulses with duration of $10\text{\,}\mathrm{ms}$ each and equal spacing of $\tau$ to the atomic sample as illustrated in Fig. 4a. The state $\rho_{3}$ resulting from this protocol is imaged after a variable evolution time $t$, which allows to record the time evolution of the state. Interrupting the DRSC protocol at any intermediate step $i$ by disabling subsequent pulses provides experimental access to the intermediate states $\rho_{i}$ (Fig. 4b).

V Numerical Phase-Space Reconstruction

The specific properties of the DRSC interaction facilitate a simple effective description of the cooling effect. First, the tight confinement of the Cs atoms in the 3D DRSC lattice potential pins the atomic position to a specific DRSC lattice site. The lattice spacing of the DRSC lattice is of the order of $1\text{\,}\mathrm{\SIUnitSymbolMicro m}$, which is much less than the typical dimension of the atomic sample in the optical dipole trap. The positions of the atoms in the dipole trap are therefore effectively frozen and the position distribution $f(z)$ does not change during the DRSC. Second, the cooling effect of the DRSC imposes a new distribution of atomic momenta to the sample. This distribution can be described in good approximation by a Maxwell-Boltzmann distribution. The temperature $T_{\text{R}}$ that characterizes the momentum distribution $\tilde{f}(p_{z})$ will be referred to as Raman temperature. Since the potential energy in the crossed dipole trap is unchanged, the DRSC pulse creates in general a nonthermal state.

The validity of this effective description of the DRSC is confirmed by experimental data in Fig. 5a for the first Raman cooling pulse. The measurements of the position distribution before and after the pulse verify that it remains unchanged during the DRSC pulse. The effect of the DRSC in the momentum distribution can be studied by observing the free evolution of the system. The state $\rho_{1}$ after the first pulse is not thermal since position and momentum distributions correspond to different temperatures, $T_{0}$ and $T_{\text{R}}$, respectively. This imbalance gives rise to the phase-space dynamics shown in Fig. 5b and can be employed to extract the Raman temperature $T_{\text{R}}$. We compare the measured evolution of the position distribution $f_{1}(z)$ to numerical simulations of the three-dimensional trapping potential with the temperature $T_{\text{R}}$ being the only free parameter. We obtain a Raman temperature of $T_{\text{R}}=$2.9(2)\text{,}\mathrm{\SIUnitSymbolMicro}\mathrm{K}$$ in a $\chi^{2}$-analysis (Appendix C). The simulation data can additionally be used as an efficient way to extract the full phase-space information as shown in the insets of Fig. 5b. This information is commonly only available at the price of additional technical effort or much larger atom number than used here afe17 ; ber18 . While the axial phase-space distribution $\rho_{1}$ would simply freely rotate in the two-dimensional space $(z,p_{z})$ for a harmonic potential, we here observe the creation of whorls induced by the nonlinearity of the trap mil86 . The projection onto the position axis shows excellent agreement between numerics and experimental data at all times. We may thus conclude that the effective model for the DRSC captures all the relevant features of the phase-space evolution. We can further simulate the full cooling protocol without any free parameters, once we have determined the initial and Raman temperatures.

The properties of the DRSC also enable the evaluation of the relative entropy (2) right after a DRSC pulse (at evolution time $t=0$ in Fig. 5b). Since the momentum distribution is randomized to the same Maxwell distribution, $\tilde{f}_{i}(p_{z})=\tilde{f}_{f}(p_{z})$, characterized by only the Raman temperature $T_{\text{R}}$ during each Raman pulse, it is independent from the position distribution $f_{i}(z)$. As a result, the phase-space distribution factorizes $\rho_{i}(z,p_{z})=f_{i}(z)\tilde{f}_{f}(p_{z})$ directly after a Raman pulse. We can thus determine the full axial phase-space distribution $\rho_{i}(z,p_{z})$ immediately after each cooling pulse. Since the factorization property also holds true for a thermal state, we have for the final state $\rho_{f}(z,p_{z})=f_{f}(z)\tilde{f}_{f}(p_{z})$. The additivity of the relative entropy for independent distributions cov06 then implies that the entropic distance between $\rho_{i}$ and the target state $\rho_{f}$ simplifies to

[TABLE]

The full relative entropy can hence be determined from the measured position distribution $f_{i}(z)$. We next discuss how this central quantity of nonequilibrium thermodynamics can be evaluated from experimental data in order to optimize the cooling of the few-particle gas.

VI Application to Experimental Data

VI.1 Optimal thermalization

The practical implementation of the two optimization criteria based on the total and positional relative entropies faces the problem that the relative entropy is only well-defined for probability distributions that are absolutely continuous with respect to one another, that is, there exists no point in phase space where one distribution vanishes, while the other one does not lin91 . Any occurrence of zero bins, due to finite statistics, in the experimentally measured or in the numerically simulated distribution in the denominator will thus result, for a nonvanishing numerator, in a division by zero (Appendix D). This issue does not seem to have been noticed in the nonequilibrium thermodynamics literature so far bus05 ; sek10 ; sei12 ; jar11 ; cil13 . We solve it by replacing the relative entropy by the closely related $K$ directed divergence, well-known in engineering, and defined as lin91

[TABLE]

It satisfies $K(\rho_{a}||\rho_{b})\geq 0$ and $K(\rho_{a}||\rho_{b})=0$ if and only if $\rho_{a}=\rho_{b}$, like the relative entropy (2). It is always well defined irrespective of $\rho_{a}$ and $\rho_{b}$. It is further bounded by the relative entropy, $K(\rho_{a}||\rho_{b})\leq D(\rho_{a}||\rho_{b})/2$ lin91 . It thus provides a lower bound to the energy irreversibly dissipated from the system during thermalization. We shall see below that the use of the $K$ directed divergence allows the optimal thermalization of the atomic gas.

Figure 6a presents the implementation of the first (static) optimization criterion for the state $\rho_{3}$. The $K$ directed divergence $K\left(\rho_{3}||\rho_{f}\right)$ is shown for various pulse spacings: the triangles correspond to numerical simulations for a harmonic trap, while the large dots are the experimental results for the nonharmonic trap. The small dots are the related simulations (Appendix D). We observe a vanishing minimum in the harmonic case at $4.2\text{\,}\mathrm{ms}$ which corresponds to a quarter of a trap period. The state $\rho_{3}$ after the last Raman pulse is here equal to the target thermal state $\rho_{f}$, revealing perfect thermalization. We experimentally find a minimum for the nonharmonic case at $6.3\text{\,}\mathrm{ms}$, in good agreement with the numerical simulations. The entropic distance to the target state $\rho_{f}$ is reduced by almost a factor two at this point compared to the nonoptimal protocols.

Figure 6b displays the results of the second (dynamical) optimization criterion for the state $\rho_{3}$. The oscillation amplitude $\Delta K=\max_{t}K(f_{3}(t)||f_{f})-\min_{t}K(f_{3}(t)||f_{f})$ after the last cooling pulse for a free evolution up to $9\text{\,}\mathrm{ms}$ is shown for different pulse spacings, both for the harmonic (triangles) and anharmonic (dots) potentials. We again observe a minimum at $4.2\text{\,}\mathrm{ms}$ for the simulated harmonic case and at $6.3\text{\,}\mathrm{ms}$ for the experimental nonharmonic potential, thus confirming the findings obtained with the first, static condition. Figures 6c-f show the time evolution of the $K$ directed divergence $K(f_{3}(t)||f_{f})$ after each cooling pulse for the optimal spacing. No oscillations are seen for the initial thermal state $\rho_{0}$ (black). These oscillations strongly increase after the first cooling pulse (red), revealing the nonthermal nature of state $\rho_{1}$, before decreasing again for the states $\rho_{2}$ and $\rho_{3}$ after the application of each additional Raman pulse (orange and blue). Finally, the oscillation amplitude reaches a minimum for $\rho_{3}$.

Both criteria may be combined to draw a map (Figs. 7a-b) of the cooling process in the plane $(K\left(\rho_{i}||\rho_{f}\right),\Delta K(f_{i}))$, similar to Fig. 1a. Figures 7c-d further show the overlap between the measured (blue bars) and simulated (blue lines) axial distributions after the last pulse, as well as the simulated target distribution (green lines) for $\tau=$ $2.1\text{\,}\mathrm{ms}6.3\text{\,}\mathrm{ms}10.5\text{\,}\mathrm{ms}$. We observe an overlap of $75\%$ for the optimal spacing of $6.3\text{\,}\mathrm{ms}$, twice the value for the other two times (Appendix F). This offers an additional confirmation of the validity of the two thermodynamic optimization criteria. We note, however, that the overlap integral does not possess any simple thermodynamic interpretation in contrast to the relative entropy or the $K$ directed divergence.

VI.2 Statistical length and horse-carrot theorem

The experimental reconstruction of the axial phase-space distribution after each Raman pulse allows us to analyze the whole cooling process by evaluating the statistical length $L_{i}$ of each cooling step. Figure 8a presents the simulated lengths $L_{i}^{K},(i=1,2,3)$ based on the $K$ directed divergence for the harmonic trap. We note that the first two steps (red and orange) have equal length for the optimal spacing of $4.2\text{\,}\mathrm{ms}$, while the length of the last step vanishes. Optimal thermalization thus occurs during the first two Raman pulses with identical entropy production. This picture is still approximately true for the nonharmonic potential (Fig. 8b): the first two statistical lengths are nearly equal for the optimal pulse spacing of $6.3\text{\,}\mathrm{ms}$, while the third one is much smaller. This is a nontrivial result: it was originally theoretically derived for close-to-equilibrium quasistatic processes sal83 ; and84 ; nul85 ; and11 ; sal98 ; sal01 ; nul02 (see also Refs. ton87 ; spi95 ; dio96 ) and has never been confirmed experimentally to our knowledge. The fact that it also holds true (exactly for the harmonic case and approximately for the anharmonic trap) for the generalized statistical length (4) (even when the relative entropy is replaced by the $K$ directed divergence) is remarkable. It suggests a quite general range of validity of the principle of equal statistical distances (or of equipartition of entropy production, as it is sometimes called ton87 ; spi95 ; dio96 ) for optimal nonequilibrium processes. Figure 9 additionally shows an experimental verification of the generalized horse-carrot theorem (5), $\Sigma^{K}\geq(L^{K})^{2}/(2n)$ for the $K$ directed divergence, as a function of the pulse spacing. It shows that the entropy production is maximal for the optimal pulse spacing, corresponding to maximal heat extraction from the system. This situation is somewhat different from the usual one, where the system of interest is continuously coupled to an ideal heat bath. In the present experiment, the phase-space evolution is mostly nondissipative as the system is only punctually coupled to a nonideal reservoir that thermalize the momentum degree of freedom. Figure 9 confirms the validity of a sharpened second law in this nonequilibrium situation.

VII Conclusion

We have experimentally studied the nonequilibrium thermodynamics of a few-particle system consisting of a gas of noninteracting Cesium atoms driven by Raman laser cooling pulses. Tracing the evolution of the gas with position-resolved fluorescence imaging enabled us to access the full phase-space density of the effectively one-dimensional system. We have used this distribution to evaluate the nonequilibrium entropy production and the statistical length based on the $K$ directed divergence. The latter quantity is always defined, in contrast to the usual relative entropy, and provides a lower bound to it. It further belongs to the family of $f$-divergences and shares their properties csi67 . As a first application, we have optimized the thermalization of the atomic gas and determined the optimal Raman pulse spacing for a nonharmonic trap potential by minimizing the entropy production to a final target state. We have additionally verified a horse-carrot theorem and analyzed the entire cooling process with the help of the statistical length. We have found that optimal thermalization is mainly achieved during the first two cooling stages, corresponding to nearly equal statistical lengths. Our findings demonstrate an effective, theoretical and experimental, method to characterize and optimize general nonequilibrium processes of few-particle systems. They further highlight the practical usefulness of nonequilibrium concepts such as entropy production and statistical lengths down to the atomic level. While we have validated our generic approach with the example of laser cooling of noninteracting atoms, the same theoretical and experimental techniques can be straightforwardly employed to include external time-dependent drivings, tunable interactions or dissipation effects. Our results thus provide a versatile platform to engineer nonequilibrium states and investigate complex far-from-equilibrium-optimization protocols for driven-dissipative interacting particles lab16 , as well as for power output mechanisms and thermal machines ros16 , both in the classical and quantum regimes.

Acknowledgements

We acknowledge financial support from the German Science Foundation (DFG) under Project No. 277625399 - TRR 185 and Grant No. FOR 2724.

Appendix A: Entropy production and statistical length

We begin by reminding the derivation of the entropy production for a single equilibration step def11 . We consider a system with Hamiltonian $H$ in an initial state $\rho_{0}$ that thermalizes to the equilibrium state $\rho_{\text{eq}}$ with inverse temperature $\beta$. The entropy production is defined as $\Sigma=\Delta S-\beta Q$, where $\Delta S=-\int dzdp_{z}\,(\rho_{\text{eq}}\ln\rho_{\text{eq}}-\rho_{0}\ln\rho_{0})$ is the entropy difference between final and initial states and $Q=\int dzdp_{z}\,H(\rho_{\text{eq}}-\rho_{0})$ the corresponding heat. Using $\rho_{\text{eq}}=\exp(-\beta H)/Z$, one readily finds sch80 ; pro76 ; def11 ,

[TABLE]

Expression (14) is the maximal amount of work that can be extracted during thermalization pro76 ; def11 . Let us now consider a multistep equilibration process with one intermediate (nonthermal) state $\rho_{1}$. The entropy production between this state and the equilibrium state $\rho_{\text{eq}}$ is $\Sigma_{1}=D(\rho_{1}||\rho_{\text{eq}})$. Using the additivity of the entropy production, $\Sigma=\Sigma_{0}+\Sigma_{1}$, we obtain the entropy production between state $\rho_{0}$ and $\rho_{1}$ as esp10 ; def11 (see also Refs. cus18 ; man18 ),

[TABLE]

Equation (15) can be generalized to an arbitrary number of nonthermal intermediate steps by recursion, yielding,

[TABLE]

Common optimization schemes consider equilibrium intermediate states generated by coupling the system to different baths at (slightly) different temperatures $T_{i}$ nul85 ; and11 ; sal98 ; sal01 ; nul02 . In this quasistatic case, $\Sigma_{i}^{\text{qs}}=D(\rho_{i}||\rho_{\text{eq},i})$, where $\rho_{\text{eq},i}$ is a thermal state at inverse temperature $\beta_{i}$. The square root, $L_{i}^{\text{qs}}=\sqrt{2\Sigma_{i}^{\text{qs}}}$, defines a statistical length in thermodynamic space nul85 ; and11 ; sal98 ; sal01 ; nul02 . It is a proper (Riemannian) distance in contrast to the relative entropy that does not satisfy the triangle inequality. Interestingly, the total entropy production, $\Sigma^{\text{qs}}=\sum_{i=0}^{n}\Sigma_{i}^{\text{qs}}$, is bounded from below by the square of the total length $L^{\text{qs}}=\sum_{i=0}^{n}L_{i}^{\text{qs}}$, that is, $\Sigma^{\text{qs}}\geq(L^{\text{qs}})^{2}/(2n)$. This result, which follows from the Cauchy-Schwarz inequality, is often referred to as the horse-carrot theorem and11 ; sal98 . It is significant because it provides a sharper lower bound to the entropy production than the second law of thermodynamics, which only states that the entropy production is non-negative. The lower bound can actually be reached, showing that dissipation can be reduced by coaxing the system along the desired path, much like guiding a horse along by waving a carrot in front of it and11 ; sal98 .

Similarly, the square root of Eq. (16), $L_{i}=\sqrt{2\Sigma_{i}}$, defines a statistical length, which reduces to the usual thermodynamic length for quasistatic processes nul85 ; and11 ; sal98 ; sal01 ; nul02 . The total entropy production is still bounded from below by the square of the total statistical length divided by twice the number of steps, $\Sigma\geq L^{2}/(2n)$, generalizing the horse-carrot theorem to nonthermal intermediate states.

Appendix B: Analytical Solution in the Harmonic Case

Using the analytical expression of the Gaussian phase-space density given in the main text, the position and momentum projections are easily integrated to

[TABLE]

where the time-dependent parameters $a^{*}=a-{c^{2}}/{4b}$ and $b^{*}=b-{c^{2}}/{4a}$ are the corresponding projected variables. The relative entropies follow as,

[TABLE]

The above expressions are employed for the calculations presented in Fig. 2 of the main text. The effect of the DRSC pulse is taken into account by setting the parameters $a_{1}$, $b_{1}$, and $c_{1}$ before the DRSC pulse to new values $a_{2}$, $b_{2}$, and $c_{2}$, determined by incorporating the constrains arising from the characteristics of the cooling: First, the DRSC erases all correlations of the state, implying $c_{2}=0$. Second, the velocity distribution is given by a Maxwellian at the Raman cooling temperature $T_{\text{R}}$. And third, the position distribution is not influenced by the DRSC. Summing up these conditions yields the parameters after the Raman cooling pulse to be

[TABLE]

Appendix C: Numerical Simulations

For the numerical simulation of the phase-space dynamics in the DRSC protocols, the atomic motion in the trap is modeled with a Monte-Carlo approach where full three-dimensional trajectories of $N=10^{5}$ atoms are calculated. This simulation only features two free parameters: First, the initial temperature of the atomic cloud $T_{0}$ determines the initial, thermal phase-space distribution, which sets the starting point for the simulation. Second, the Raman cooling temperature $T_{\text{R}}$ is employed to model the effect of the DRSC by resetting the atomic velocities to a Maxwell-Boltzmann distribution corresponding to $T_{\text{R}}$, whenever a DRSC pulse is applied. Using these two temperatures together with precise information on the trap, which was specified by independent trap frequency and beam shape measurements, the effect of arbitrary pulse sequences on the phase-space distribution and the ensuing dynamics can be computed. In this section, we show how the experimental value for $T_{0}$ is extracted from the measured initial distribution $f_{0}(z)$ and the Raman cooling temperature $T_{\text{R}}$ is determined from the measured evolution after the first Raman cooling pulse $f_{1}(t,z)$.

In order to model the position distribution $f_{0}(z)$, we employ a simulation scenario, where atoms are initially located at the trap center. A heat bath at temperature $T_{0}$ is emulated by resetting the atomic velocities repeatedly to random velocities corresponding to the desired initial temperature $T_{0}$. Due to the resulting damped motion of the atoms in the trap, the atomic position distribution approaches a thermal distribution at $T_{0}$ hoh16 . We compare the simulated position distributions $f_{\mathrm{sim}}(z)$ for various temperatures $T_{0}$ to the experimentally measured initial position distribution $f_{\mathrm{exp}}(z)$ shown in Fig. 10a by calculating the $\chi^{2}-$value,

[TABLE]

for the binned data as a measure for the goodness of the fit bev03 ($\Delta f_{\mathrm{sim}}$ and $\Delta f_{\mathrm{exp}}$ are the statistical uncertainties of $f_{\mathrm{sim}}$ and $f_{\mathrm{exp}}$). The $\chi^{2}-$value for simulations at various temperatures is shown in Fig. 10b, where we use a polynomial fit to the data in order to extract the initial temperature $T_{0}=$12.1(11)\text{,}\mathrm{\SIUnitSymbolMicro}\mathrm{K}$$.

The final temperature $T_{\text{R}}$ which corresponds to the DRSC temperature is not visible in the position distribution directly after a DRSC-pulse. However, as illustrated in Fig. 5b of the main text, the evolution in the trapping potential after the first DRSC pulse shows clear evidence of the cooling effect by featuring a breathing behavior. In order to extract the value of $T_{\text{R}}$, we simulate the time evolution of atomic samples which are prepared at the initial temperature $T_{0}$ and then reset the atomic velocities to values corresponding to different Raman cooling temperatures $T_{\text{R}}$. For every evolution time $t$, we extract a Raman cooling temperature $T_{\text{R}}$ with a $\chi^{2}-$analysis, analogous to the strategy employed for the initial distribution, by comparing the simulations for different Raman cooling temperatures to the experimental distribution (Fig. 10c). We combine the results of all measured evolution times shown in Fig. 10d by a weighted constant fit to the data, thereby extracting the DRSC temperature $T_{\text{R}}=$2.9(2)\text{,}\mathrm{\SIUnitSymbolMicro}\mathrm{K}. The red shaded area in the plot indicates small evolution times $t$ where the $\chi^{2}-$analysis fails, because the information about the velocity distribution is not yet transformed into the position distribution. This behavior is also visible in the size of the error bars, which first decreases until $t=$6\text{\,}\mathrm{ms} and then increases again. Combining the extracted values for $T_{0}$ and $T_{\text{R}}$, the simulation dataset corresponding to $T_{0}=$12\text{,}\mathrm{\SIUnitSymbolMicro K} and $T_{\text{R}}=$3\text{\,}\mathrm{\SIUnitSymbolMicro K} is the best fitting simulation. Therefore, this dataset is employed for the calculation of the simulation data points presented in the main text. Accordingly, the final thermal state $\rho$ is also represented by the simulation data for a thermal state at temperature $T_{\text{R}}=$3\text{,}\mathrm{\SIUnitSymbolMicro K}$$.

DRSC is in general a subrecoil cooling scheme, because it fundamentally allows to reach subrecoil temperatures. The temperature of $2.9\text{\,}\mathrm{\SIUnitSymbolMicro K}$ observed in the experiment is clearly above the recoil temperature of Cs which is $0.1\text{\,}\mathrm{\SIUnitSymbolMicro K}$ for the DRSC laser light. This optimum is not reached due to technical limitations like laser power noise, laser linewidths, off-resonant photon scattering.

Appendix D: Numerical calculation of the relative entropy

For the analysis of our data, we typically bin the atomic positions from an experiment or a Monte-Carlo simulation in order to create a numerical representation of the density distribution. The integral for the relative entropy then corresponds to a sum over all bins $z_{i}$, where

[TABLE]

As discussed in Ref. lin91 , typical data with finite statistics may exhibit bins where $f_{2}(z_{i})=0$, meaning that no atom has been observed in this bin. However, this corresponds to a division by zero in Eq. (23), rendering the calculation of the integrand value $\mathcal{D}_{i}$ impossible for this specific bin. In contrast, the directed divergence lin91 ,

[TABLE]

can be evaluated even at bins where $f_{2}(z_{i})=0$. It is therefore much more robust especially when analyzing experimental data, where statistical errors are usually even more pronounced.

In order to illustrate the problem, we employ the data set used for the harmonic approximations shown in Figs. 6 and 8 of the main text. The corresponding initial temperature for the simulation is $T_{0}=$1\text{,}\mathrm{\SIUnitSymbolMicro K} and the final temperature is $T_{\text{R}}=$0.25\text{\,}\mathrm{\SIUnitSymbolMicro K}. While these values are more than one order of magnitude colder than the experimental parameters, the ratio of the two temperatures is the same as in the experiment, thereby providing a comparable cooling process. Nevertheless, at these low temperatures, the harmonic approximation of the trapping potential holds also in axial direction. In fact, the density distributions of the Monte-Carlo simulation (bars) shown in Fig. 11a fit very well to Gaussian distributions (solid lines) which are expected for the harmonic case. Figures 11b and c show the integrands $\mathcal{K}_{i}$ and $\mathcal{D}_{i}$, where again the bars correspond to the numerical data and the solid lines show the Gaussian fit. The missing bars seen in Fig. 11b clearly indicate the numerical problem connected to the relative entropy. By contrast, the integral for the directed divergence in Fig. 11c can be evaluated in the whole range.

Appendix E: Contributions of the momentum distribution

For factorized distributions $\rho_{i}=f_{i}(z)\tilde{f}_{i}(p_{z})$ the relative entropy $D$ can be split into two contributions $D\left(\rho_{i}||\rho_{f}\right)=D(f_{i}(z)||f_{f}(z))+D(\tilde{f}_{i}(p_{z})||\tilde{f}_{f}(p_{z}))$, where the first term accounts for the position distribution and the second takes into account the momentum distribution. After a Raman cooling pulse (for $i=1,2,3$), the contribution of the momentum distributions is zero, because $f_{i}(p_{z})$ and $f_{f}(p_{z})$ are identical. For the initial distribution ($i=0$), however, this contribution is not zero, as here the momentum distributions are not equal. In the measured position distributions $f_{0}(p_{z})$ at $t=0$, this contribution is not visible. However, as the initial ($T_{0}$) and final ($T_{\text{R}}$) temperatures are known, the contribution $D(\tilde{f}_{0}(p_{z})||\tilde{f}_{f}(p_{z}))$ can be calculated from the thermal momentum distributions $\tilde{f}(p_{z},T)=1/\sqrt{2\pi mkT}\cdot\exp(-p_{z}^{2}/(2mkT))$. We find for the directed divergence employed in Figs. 6 and 7 this contribution of the momentum distribution to be $K(\tilde{f}(p_{z},T_{0})||\tilde{f}(p_{z},T_{\text{R}}))=0.071$ by solving the integral numerically. The hollow experimental points in Figs. 6 and 7 are thus a combination of the measured contribution to the directed divergence from the position distribution and the numerically deduced contribution from the momentum distribution.

Appendix F: Overlap calculation

The overlap of the distribution after the last cooling pulse $f_{3}(0)$ with the final distribution $f_{f}$ is evaluated in the following way. We first renormalize the final distribution (green line) to the maximum of the experimental data (blue bars). Integration of the renormalized final distribution then yields the amount of atoms in the experimental distribution that match the final distribution. We identify this value with the overlap. We find the largest overlap at a pulse spacing of $6.3\text{\,}\mathrm{ms}$ which corresponds to our optimization result (Fig. 12).