Quantum critical region of ultracold Bose gases exhibiting universal density-probability distribution after free expansion

We report the observation of a universal behavior of ultracold quantum critical (QC) Bose gases trapped in a one-dimensional optical lattice. We extract the density probability distributions from the measured atomic density of the rubidium Bose gases after a free expansion time. The density probability is found to follow a simple exponential law once the lattice system has entered the QC region above the Berezinskii–Kosterlitz–Thouless transition. We show that, in addition to relative phase fluctuations between the subcondensates in different lattice wells, there also exist spatial phase fluctuations within single lattice wells in this QC region. The universal density-probability distribution can thus be well understood by a simple theoretical model taking into account these two kinds of phase fluctuations.


Introduction
The nature of quantum criticality [1][2][3] driven by quantum fluctuations is still a great puzzle, despite of the remarkable advances in heavy-fermion metals and rare earth-based intermetallic compounds, etc [4]. New understanding of quantum criticality is widely believed to be a key to resolving open questions in metal-insulator transitions [5], high temperature superconductivity [6] and novel material design, etc. Cold atoms in optical lattices provide a unique chance to not only simulate other strongly correlated systems [7], but also study some models unaccessible in solid state systems, particularly for the Bose-Hubbard model [8,9].
Ultracold Bose atoms in one-dimensional (1D), two-dimensional (2D), three-dimensional (3D) optical lattices have been widely studied with the Bose-Hubbard model [10]. At zero temperature, there is a continuous quantum phase transition from superfluid (SF) to Mott insulator (MI). Due to the strongly correlated quantum behavior, the quantum critical point (QCP) at absolute zero temperature distorts strongly the structure of the phase diagram at finite temperatures, leading to the emergence of an unconventional 'V-shaped' quantum critical (QC) region [2,[11][12][13][14][15][16][17]. In analogy with a black hole, the crossover to QC gas involves crossing a 'material event horizon', which implies strongly that the QC gas has a simple and universal behavior [3]. In a recent work [18], a scaling behavior of the equation of state for a 2D optical lattice was observed based on in situ density measurements.
In the commonly used 3D optical lattice systems, the atoms in a lattice well are strongly confined in all directions, with a mean occupation number per lattice well of about ∼ 1 3 [10]. Thus the spatial phase fluctuations for the atoms in a single lattice well can be omitted, and the relevant theoretical calculations based on the Bose-Hubbard model [8,9] and Wannier function [19] can give quantitative descriptions for almost all experimental phenomena [7]. In contrast, for ultracold Bose atoms in a 1D optical lattice, unique characteristics may arise due to the quasi-2D nature of the Bose gases as well as the Berezinskii-Kosterlitz-Thouless (BKT) transition [20][21][22][23] associated with such a quasi-2D gas.
In this paper we report the observation of universal behavior for ultracold QC Bose gases in a 1D optical lattice. Different from the usual method, the density probability distributions are extracted from the measured density distributions of the released gases after sufficiently long free expansion time, for different depths of the lattice potential. It was found that the density probability follows a simple exponential law when the Bose gases reach the QC region above the BKT transition. Our theoretical model reveals that this universal densityprobability distribution is attributed to the existence of two kinds of phase fluctuations, i.e., relative phase fluctuations between different lattice wells and spatial phase fluctuations within single lattice wells. These intersite and intra-site phase fluctuations in the QC region will lead to significant density fluctuations for the released gases after the free expansion, and then result in the observed exponentially-decaying density probability.
The exponential density probability distribution is universal in a sense because for a system with completely random phase fluctuations in two or three dimensions, there is always this unique behavior. It has been shown that the light can have exponential density probability distribution if there exist phase fluctuations with two freedoms [24]. It has also been shown that for different quantum chaotic systems of matter wave in two dimensions, there is always exponential density probability distribution after sufficiently long dynamic evolution [25]. Our system is not a quantum chaotic system. However, due to the existence of two kinds of phase fluctuations, we will give the first experimental demonstration of this exponential density probability distribution with matter waves. Our work paves the way to further demonstrate the similar exponential behavior for fermionic SF, molecular condensate and other novel systems with large phase fluctuations. It also provides a new tool in analyzing the property of the QC region for other cold atomic systems.
The paper is organized as follows. In section 2, we describe the experiment and the measurement of density probability. In section 3, we plot the phase diagram for our 1D optical lattice system. It is followed by section 4, where a simple theoretical model is presented to explain the observed exponentially decaying densityprobability. Numerical simulation is also presented in this section. Finally, we conclude the paper in section 5.

Experiment
The experiment starts with a nearly pure 87 Rb condensate containing ∼ × 1.3 10 5 atoms. The condensate is created in a magnetic trap with axial and transverse trapping frequencies of ω ω π = ⊥ { , } 2 {18.7, 205} z Hz. The magnetic trap is then adiabatically relaxed to lower trapping frequencies of π 2 {8.2, 90.4} Hz over a time of 500 ms. Accordingly, the condensate attains a larger spatial extension and is thus able to occupy more lattice wells. The 1D optical lattice is formed by a retroreflected laser beam of λ = 1064 nm along the axis (z direction) of the condensate. The laser beam is focused to a waist of μ 140 m at the position of the condensate. It is ramped up to a given intensity over a time of 10 ms, yielding a lattice potential π λ = V sE z sin (2 ) r opt 2 , with E r being the recoil energy of an atom absorbing one lattice photon and s the potential depth in units of E r . After a holding time of 5 ms, we suddenly switch off the combined potential and allow the cold atomic cloud to expand freely for a time of 30 ms. The expanded atomic gases are probed using the conventional absorption imaging technique. Since the probe beam is applied along the x direction, what an absorption image records is the 2D (y z , ) column density profile of the atomic cloud. We denote by N y z ( , ) The potential depth of the optical lattice is calibrated using the method of Kapitza-Dirac scattering as in [26]. The temperature of an ultracold gas is determined by the method of bimodal-profile fit as in our previous work [27]. After the time of flight of 30 ms, the normal gas in the optical lattice has a Gaussian distribution in the transverse (y) direction while the condensed gas has a Thomas-Fermi profile. The temperature can be deduced from the rms size of the normal gas. In our experiment, the measured temperatures are in the order of tens of nano-Kelvins.
The left column of figures 1(a1)-(a6) displays the spatial distributions of the released atomic gas for different lattice depth s. From the top three images ((a1)-(a3)), we see increasing side peaks in the interference patterns as s is increased. Further increasing s, the interference fringes become blurred and disappear eventually (see (a4)-(a6)), similar to the experimental observation in [28]. The gradual disappearance of the interference fringes can be explained by the increase of the random relative phase during the SF-QC transition. As the lattice system moves away from the SF state with increased lattice depth, density fluctuation along the z direction becomes large which is attributed to the random relative phase between the subcondensates in different lattice wells [28]. In addition, we find that significant density fluctuations also exist along the y direction. Figures 2(c) and (d) demonstrate a typical case of simultaneous existence of density fluctuations in these two directions. In contrast, for a shallow depth, the released atomic gas shows much smaller density fluctuations (see figures 2(a) and (b)).
To qualitatively analyze the statistical property of the atomic density, we define a density probability as follows: total where δn is the width of the density interval. S total is the total area of the region where the atomic gas could reach, which is shown by a red dashed line. The red dashed line is chosen so that all region effectively occupied by atoms is included. S n ( )is the area of the region where the density is larger than n. The averaged density is = n N S s total with N being the total atom number within the red dashed line, measured for each run of the experiment. For convenience, we prefer a dimensionless density of = n n n s 0 . Then we can use a dimensionless density The blue circles in the right column of figure 1 are density probabilities calculated from the corresponding density distributions in the left column. For each image, the calculation is performed only for the region enclosed by a red dashed line. The boundary of the atomic cloud has been set in such a way that the atomic density outside the enclosed area has dropped to zero within the measurement error. Note that the horizontal spacing of the discrete data points is actually the small interval of n 0 used in the calculation of P n ( ) 0 0 . Due to the optical noise of the probe beam, P n ( ) 0 0 can be nonzero even for negative n 0 . This can be understood from the expression of n y z ( , )mentioned above. If n b 0 denotes the observed maximum negative density, optical noises in each P n ( ) 0 0 plot are concentrated within a small range of − < < n n n b b 0 0 0 . Thus, the range of > n n b 0 0 is regarded as the effective region where the optical noise is negligible and P n ( ) 0 0 truly reflects the density probability distribution. In our experiment, n b 0 is typically less than 0.5. The red solid lines in figure 1 correspond to an ideal exponentially-decaying density probability in the following form: e n 0 0 0 They are plotted for the convenience of comparison with the observed density probability values. It is obvious that, with increasing s, P 0 has a tendency to P e 0 . In figures 1(b4)-(b6), P 0 agrees with P e 0 in the whole effective region of n 0 .
We should point out that a large quantity of pixels are favorable for analyzing the statistical property of a density distribution. However, limited by the pixel size of the CCD camera in the imaging system, the total number of the occupied pixels in an absorption image is typically less than 3000. It will cause noticeable error in determining the density probability. Therefore, the relative error in figure 1 is noticeably higher for very large n 0 due to the small quantity of pixels involved. In the numerical simulation in section 4, we deliberately use a smaller pixel size in order to increase the pixel number, and hence to improve the accuracy of the density probability.
To quantify the discrepancy between the measured density probability and the exponentially decaying curve of P e 0 , we first define a standard deviation as follows:   Figure 3 displays the standard deviation for different lattice depths. Apparently, the deviation is very large for low lattice depth. However, it drops quickly with increased lattice depth, until a depth of ∼ E 40 r is reached. Beyond this point, the deviation is small and rather flat, indicating that the density probability distribution is almost identical to P e 0 . As will be discussed in section 4 the exponential distribution is an indicator of the QC region above the BKT transition.

Phase diagram
Before the free expansion, the Rb gas is trapped in a combined potential consisting both a 1D optical lattice and a harmonic trap. The total trap potential is given by The two trap frequencies ω ⊥ and ω z depend mainly on the magnetic trap since the optical trap is much weaker in both the axial and radial directions. For atoms in a lattice well, they experience an effective harmonic potential along the z direction with an angular frequency ω ≃ ℏ ∼ s E 2 z r . There can be hundreds of atoms in a single lattice well. The Bose gas in each lattice well becomes quasi-2D [29][30][31][32], as long as ω ℏ ≫ ∼ k T z B is satisfied as in our case.
For such a quasi-2D Bose gas, the critical temperature in the occupied lattice well is given by , with N l being the atomic number in the lattice well and ζ the Reimann zeta function. Well below T 2D , the whole system becomes a chain of subcondensates. In QC and MI region, there is no correlation between the subcondensates in different lattice wells, and a quasi-2D gas should undergo the BKT transition [20,21] at a critical temperature T BKT . By minimizing the free energy = − F E TS (E and S are the energy and entropy of the Bose gas in a lattice well when vortices are considered), one can deduce the BKT transition temperature: = T T 4 BKT 2D [32]. Above the critical value of T BKT , it is favorable to create vortices in a quasi-2D subcondensate, and the unbinding of bound vortices will lead to strong spatial phase fluctuations within the single subcondensate.
The situation is more complicated in the SF region where the Bose gases are highly correlated. For a series of M highly correlated quasi-2D Bose gases, the free energy becomes = − F ME TS. Since the possibly existing vortices in different lattice wells are highly correlated in spatial locations, S can be approximated as the entropy of a single subcondensate. Then the BKT transition temperature becomes 2D . This means that the temperature of a quasi-2D condensate is always much below the BKT transition of the entire lattice system. Hence, it is not necessary to consider BKT transition in the SF region. In figure 4 we plot the phase diagram for our 1D optical lattice. The QCP is calculated based on the usual Bose-Hubbard model for 1D optical lattice [7,33]. Note that, phase fluctuations are not considered in the calculation of QCP. Phase fluctuations in each lattice site will not cause significant change in U and J, and thus lead only to high order correction to the QCP. This high order correction will not cause serious problem to our work, because what we study is the QC region, and we cannot access the QCP due to the finite-temperature effect. The crossover from the SF to the QC region is obtained from the phase transition by considering the effective mass due to the presence of the optical lattice [11]. Similar to [11], the dashed line in figure 4 between the QC region and the MI region is only a schematic because there is no clear boundary between the QC region and the MI region.
Typical parameters of the central lattice well have been adopted in the calculation of T 2D , T BKT , tunneling energy J, as well as the on-site repulsion U. Phase transition to SF state occurs at a critical temperature T C , which has been calculated by using the analytic estimation derived in our previous work [27]. One see from figure 4 that, at zero temperature, a continuous phase transition from SF to MI will occur as the lattice depth is increased. Whereas, for finite temperature, a V-shaped QC region appears between the SF and MI regions. The line of T BKT divides the QC region into two parts (type I and type II). Our experiments focused on the SF region and the QC region above the BKT transition. We were unable to access the QC region below the BKT transition due to the relatively high temperature.

Density probability distribution
We here present the analytical derivation of the exponentially decaying density probability found in the type-I QC region. The 1D lattice system consists of a chain of subcondensates, with + k 2 1 M occupied wells (−k M to k M ). The atomic number in each lattice well can be calculated for given trap parameters and a total atomic number [27]. The initial wave function in the kth lattice well can be written as φ = x y z t t ( , , , ) 0 . After switching off the combined potential, the atomic gas is allowed to expand over a time of t so that the expanded subcondensates overlap completely with each other. The 3D density distribution is then written as: can be obtained directly from the free expansion along z direction, while ϕ k r is a completely random relative phase which accounts for the loss of coherence between the subcondensates. Without the random phase ϕ k r , the phase ϕ k 0 will lead to clear interference fringes along z direction with a period of π λ ℏt m 4 [29]. In the QC region, ϕ k s is completely random, and there would be no interference fringes along z direction, just as shown in figures 1(a5)-(a6). . Finite-temperature phase diagram of an ultracold rubidium gas in a 1D optical lattice. U and J refer to the mean field energy per particle and the tunneling energy, respectively. T C is the critical temperatures of Bose-Einstein condensation for the whole combined trap, while T har is the critical temperature for a purely harmonic trap in absence of the optical lattice. By increasing the depth of an optical lattice, U/J can be increased by several orders of magnitude. All the phase boundary lines, except for the dashed line between QC and MI regions, are calculated with a total atom number of × 1.3 10 5 and the same trap parameters as in our experiment. The red solid circles (a2)-(a6) correspond to the absorption images in figures 1(a2)-(a6).
vanishes, and ϕ k s is also very small due to the nearly perfect coherence between the subcondensates. On the contrary, in the QC region above the BKT transition, both of them are completely random. For the QC region below the BKT transition, ϕ k s is completely random, while ϕ ⊥ k is negligible. With numerically calculated density distribution, we further analyze the corresponding density probability so as to make a comparison with the experimental results. In our calculation all parameters have been chosen to match our experiment except for the pixel size. A smaller pixel size has been used to reduce the error of density probability. The boundary of an atomic cloud is chosen in a similar way as for the experiment, so that the expanded wave packets of the subcondensates are almost completely enclosed. Figure 5 displays the simulation results for three typical cases in which an atomic gas initially stays in SF, type-I QC and type-II QC regions, respectively. For the case of SF state (see (a1) and (b1)), the density probability distribution does not match the exponential curve at all. For the state in type-I QC region ((a2) and (b2) in figure 5), the density probability distribution shows a good match to the exponential curve, in agreement with our experimental results. Significant density fluctuations appear in both the y and z directions, which has also been verified in our experiment (figure 2).
In our simulation shown in figure 5(a2), there is a slight difference compared with the experimental result given by figure 1(a5). In our experiment, there is a broad density peak near the origin. A possible reason for this difference lies in that in our simulation both kinds of phase fluctuations have been assumed to be completely random, but the lattice depth is not high enough to satisfy these assumptions. For a 1D optical lattice with large filling factors in the order of several hundreds as in our case, the relative phases between neighboring lattice sites are indeed completely random in the MI region. However, it is very difficult for such a 1D optical lattice to reach the MI region which requires extremely high depth levels. Figures (a3) and (b3) were obtained using the same parameters as in (a2) and (b2) except that ϕ ⊥ k is set to zero. This is of course a situation of QC region below the BKT transition. Obviously, the discrepancy between the computed density probability and the exponential curve is noticeably larger than that in the type-I QC region. Besides, the interference pattern is quite different from that observed in our experiment. High density fluctuations appear in the z direction while a roughly smooth density profile is seen in the y direction.
There are two regions in the phase diagram which have not been studied experimentally. Here we give a brief discussion for their density probability distributions. (i) In the MI region, the random phase between neighboring lattice site is completely random. Hence, if the temperature is higher than T BKT but much lower than T 2D , we expect the exponential form similar to the case of type-I QC region. In contrast, if the temperature is lower than T BKT , we expect the same density probability distribution as that of the type-II QC region. (ii) As the temperature increases toward T 2D , the density probability distribution should approach the case of normal gas, which can be approximated as a Gaussian distribution [25].
We stress that from the derivation of this section and [24,25], the exponential density probability distribution is universal because it originates from the phase fluctuations and does not depends on the details of the system. Once there are appropriate completely random phases for either light or matter waves, there would be exponential density probability distribution. This exponential law is unique since the normal gas will show a completely different Gaussian distribution [25]. Hence, the density probability distribution provides a powerful tool to reveal the phase fluctuations in different systems.

Discussion and conclusion
The coexistence of two random phases makes the lattice system similar to the situation of the interference of randomly scattered waves [24], where an exponential distribution is also found for the field intensity. The exponential law is due to the fact that the intensity distribution is always exponential for completely random superposition of quantum states [25]. Our results have revealed a simple universal behavior despite of the extremely complex many-body state in the QC region. It seems that density probability distribution is a useful tool to identify some quantum phases in optical lattice systems.
The exponential density-probability distribution implies that a system in the QC region above the BKT transition cannot be regarded as a classical thermal cloud, although the temperature is higher than the critical temperature T c . A classical thermal cloud features density fluctuations as δ ≈ n n 2 . In contrast, the exponential density probability distribution corresponds to much larger density fluctuations of δ ≈ n n 2 2 [25]. In fact, we had also measured the atomic densities for a purely thermal gas which was initially confined in the magnetic trap with a temperature well above T c . Smaller density fluctuations were observed, and the density probability showed larger deviation from the exponential curve than a QC gas.
The present work supports the long-standing belief that the quantum states in the QC region are strongly correlated, and there should be a simple universal behavior, once an appropriate physical quantity is found [2][3][4]. In the future, we wish to reach the type-II QC region by further cooling the atomic gas. The BKT transition is expected to manifest itself in the statistic property of the atomic density as the temperature is decreased.