text:
# Bose-Einstein condensation and superfluidity of dilute Bose gas in a random potential

## I Introduction

Bose-Einstein condensation (BEC) and superfluidity of liquid 4He in random environments including Aerogel and Vycor glass are active problems in quantum fluid research. In particular, finding out how spatial confinement affects the Bose fluid has stimulated both experimental and theoretical studies.

Below the $`\lambda `$ temperature of 2.17 K, liquid 4He enters the superfluid state and behaves as though it has no viscosity. Superfluidity is a macroscopic quantum phenomenon as well as superconductivity, and understanding both have been one of the major goals of quantum statistical physics. The various observations of superfluidity was successfully explained by the phenomenological two fluid model Tilley , which is based on the idea that the system consists of an inviscid superfluid and a viscous normal fluid. On the other hand, the lambda transition had been thought to be caused by BEC, which was confirmed by neutron scattering experiments Sokol . With a BEC, a macroscopic number of particles occupies a single particle ground state and is described by a macroscopic wave function. The inviscid superflow can be described by this wave function Huang1 . However the relation between BEC and superfluidity is not completely understood. Although superfluidity and BEC are closely related to each other, one is not necessary or sufficient for the other. For example, in a two-dimensional Bose system, Kosterlitz and Thouless proved that superfluidity can exist even without BEC Kosterlitz , and the superfluidity was actually observed in 4He films Bishop . The Bose system in a random environment might be another good example for studying the relation between BEC and superfluidity Trivedi ; Huang2 . This system has received considerable attention because localization effects allow condensed particles to belong to the normal fluid rather than the superfluid and BEC can separate from the superfluid. The phase diagram of this system has been discussed, showing a specific non-superfluid phase. Thus, studies of this system can reveal the relation between BEC and superfluidity.

Porous glass such as Vycor is often used as a random media in experimental studies. Vycor glass is a from 30% to 70% porous glass containing wormholelike pores, the characteristic diameters of which vary from 30Å to 100Å. By adjusting the pore size and the adsorbed 4He coverage, we can change the density of 4He and the superfluid transition. By using torsional oscillators, Reppy et al. Reppy ; Crooker  interesting features at various pore sizes of Vycor glass or coverages were observed, particularly the behavior of the superfluid critical temperature and the temperature dependence of the superfluid density. The superfluid component has two-dimensional behavior when the pore size of Vycor glass is large, and becomes three-dimensional as the pore size is reduced. The superfluid density in such porous glass is smaller than that of bulk 4He and its critical temperature decreases with the coverage. Below a certain coverage, the superfluid density can no longer exist, even near 0 K. These results show that superfluidity is broken by the random environment. It is also important to find out how disorder affects the BEC. BEC and its elementary excitation in liquid 4He can be observed by neutron scattering. Bulk 4He has the usual excitations such as phonons, maxons and rotons Tilley . They are collective excitations on BEC, the excitations except for phonons being absent above the critical temperature. Dimeo et al. Dimeo  and Plantevin et al. Plantevin  used neutron scattering and a torsional oscillator to measure the elementary excitations and the superfluid transition, respectively, of 4He in porous glass. Surprisingly, the dispersion curve in porous glass was the same as that in the bulk, which means that the disorder does not affect the elementary excitations. Furthermore, these elementary excitations were observed even above the superfluid critical temperature. Hence, the BEC might persist in the disordered system even above the superfluid critical temperature. These results have not been understood completely yet; thus, the exact relationship between BEC and superfluidity remains puzzling.

This problem is also interesting from the following theoretical standpoint. In a Bose system confined in a random media, the long-range-order correlation due to BEC can compete with the disorder, so that the BEC critical temperature can be reduced. Huang and Meng proposed a model of the three-dimensional dilute Bose gas in a random potential Huang3  that assumed a small coverage of 4He in Vycor glass . Because it is difficult to formulate the random potential for the porous glass, they used a delta-functional impurity potential and analyzed their model using the Bogoliubov transformation and taking an ensemble average. They found that both BEC and superfluidity are depressed by the random potential and the superfluidity disappears below a critical density, even at 0 K, which is qualitatively consistent with the observations by Reppy et al. They also predicted a reentrant transition at low densities; that is, the superfluid phase enters a normal phase again with decreasing temperatures. However, the random potential of their model does not include the pore size, and thus it is difficult to quantitatively compare to experimental results for a range of pore sizes. Another model is the Bose Hubbard model with the random potential. By considering the transfer energy, the on-site repulsion, and the random potential, Fisher et al. Fisher  found that the Bose glass phase can exist with the superfluid phase and the Mott insulating phase. The Bose glass phase is similar to the Anderson insulating phase Anderson  in metal. In the Bose glass phase, the condensed particles are localized and thus do not contribute to superfluidity. Thus, the Bose glass phase could influence the collective excitations even above the superfluid critical temperature. However, the theoretical excitation energy Krauth  for the Bose glass phase disagrees with measurements by neutron scattering experiment Plantevin , so it is not yet clear whether the Bose glass phase has actually been detected. Finally, it should be noted that Huang and Meng’s model cannot describe the Bose glass phase because the ensemble average makes the system uniform.

Few theoretical studies of this random system are quantitatively compare to the experiment. Thus, in this work, we improve Huang and Meng’s model Huang3  by adding the size dependence of the random potential instead of using their the delta-functional potentials. The strength of the random potential can be estimated by comparing calculated and experimental critical coverages below which the superfluid density disappears, even at 0 K. As a result, our model has no free parameters and can be used for quantitative comparisons to experimental data. This enables us to determine whether or not our picture of the three-dimensional dilute Bose gas in random potential is applicable to a real system. Our formulation cannot address this question at high temperatures due to the high number of thermally excited quasi particles. As far as the condensate density is almost independent of temperature at low temperatures, however, our formulation works well, leading to the following results. (1) The specific heat agrees quantitatively with experimental data at low temperatures. (2) Because of the random potential, the specific heat is not proportional to $`T^3`$, as occurs for phonons, but to $`T`$. Furthermore, by obtaining the condensate density and the superfluid density, we found the following. (3) When the total density is sufficiently low, BEC can persist even when the superfluid density disappears below that critical coverage. (4) The random potential causes a reentrant transition of the superfluid phase. Finally, we show why decreasing the open pore density of the Vycor glass should allow the reentrant phase to be detected experimentally

A brief summary of our paper is as follows. In Sec. II, we describe our model of the dilute Bose gas in a random potential and derive the partition function. Section III tests our model by quantitatively comparing calculated to experimental specific heats. In Sec. IV, the BEC density and the superfluid density are obtained and their characteristics are discussed. Section V is discussion and conclusions.

## II Model

Superfluid 4He adsorbed in Vycor glass can be modeled by a three-dimensional dilute Bose gas in a random external potential Huang3 . The grand canonical Hamiltonian is

$`\widehat{H}\mu \widehat{N}`$ $``$ $`\widehat{K}`$ (1)
$`=`$ $`{\displaystyle 𝑑x^3\widehat{\mathrm{\Psi }}^{}(𝒙)\left[\frac{\mathrm{}^2}{2m}\mathbf{}^2+U(𝒙)\mu \right]\widehat{\mathrm{\Psi }}(𝒙)}`$
$`+{\displaystyle \frac{v_0}{2}}{\displaystyle 𝑑x^3\widehat{\mathrm{\Psi }}^{}(𝒙)\widehat{\mathrm{\Psi }}^{}(𝒙)\widehat{\mathrm{\Psi }}(𝒙)\widehat{\mathrm{\Psi }}(𝒙)},`$
where $`\widehat{\mathrm{\Psi }}(𝒙)`$ is the field operator for Bose particles of mass $`m`$, $`\widehat{N}=𝑑x^3\widehat{\mathrm{\Psi }}^{}(𝒙)\widehat{\mathrm{\Psi }}(𝒙)`$ is the number operator, $`\mu `$ is the chemical potential, and $`U(𝒙)`$ is the external random potential that represents the effect of Vycor glass. The first term of the Hamiltonian is the kinetic energy and the external potential, whereas the second term refers to the hard-sphere interaction between particles with $`v_0=4\pi a\mathrm{}^2/m`$ being the coupling constant with the s-wave scattering length $`a`$. This repulsive interaction prevents all particles from being localized at the minimum of $`U(𝒙)`$. This has similarities to the Fermi system with disorder Anderson ; Langer ; Gor'kov ; for example, Fermions cannot localize in a single orbital in space due to the Pauli exclusion principle. Therefore, the Fermion system is stable even if it is free from the repulsive interaction. On the other hand, to prevent the system from collapsing into the minimum of $`U(𝒙)`$, the Bose system should include a repulsive interaction. This makes the problem more complicated than that of the Fermi system.

Proceeding in a standard fashion, we introduce the free particle annihilation and creation operators $`\widehat{a}_𝒌`$ and $`\widehat{a}_𝒌^{}`$. We assume that the level with $`𝒌=\mathrm{𝟎}`$ is macroscopically occupied with occupation number $`N_0`$, so $`\widehat{a}_0`$ and $`\widehat{a}_0^{}`$ are replaced by a c-number $`\sqrt{N_0}`$. By making a Fourier transformation and neglecting all off-diagonal terms $`U_𝒌\widehat{a}_𝒌^{}^{}\widehat{a}_{𝒌^{\prime \prime }}`$ and $`v_0\widehat{a}_𝒌^{}\widehat{a}_𝒌^{}^{}\widehat{a}_{𝒌^{\prime \prime }}\widehat{a}_{𝒌^{\prime \prime \prime }}`$, we obtain

$`\widehat{K}_{\mathrm{eff}}`$ $`=`$ $`V\left(\mu n_0+{\displaystyle \frac{1}{2}}v_0n_0^2+{\displaystyle \frac{n}{V}}U_0\right)+{\displaystyle \underset{𝒌0}{}}\left[{\displaystyle \frac{\mathrm{}^2k^2}{2m}}\mu +v_0(n+n_0)\right]\widehat{a}_𝒌^{}\widehat{a}_𝒌+\sqrt{{\displaystyle \frac{n_0}{V}}}{\displaystyle \underset{𝒌0}{}}\left(U_𝒌\widehat{a}_𝒌^{}+U_𝒌\widehat{a}_𝒌\right)`$ (2)
$`+{\displaystyle \frac{1}{2}}v_0n_0{\displaystyle \underset{𝒌0}{}}\left(\widehat{a}_𝒌^{}\widehat{a}_𝒌^{}+\widehat{a}_𝒌\widehat{a}_𝒌\right),`$
where $`V`$ is the volume of the system, $`n_0=N_0/V`$ is the number density of condensate, and $`U_𝒌`$ is the Fourier transformation of $`U(𝒙)`$. By neglecting the off diagonal terms, we are neglecting the interactions between the excited particle and the random potential and that between pairs of excited particles; these become important as the temperature rises and the condensate density decreases. Hence, this approximation is poor when many particles are thermally excited. Nevertheless, this approximation is useful at low temperatures where the condensate density is almost independent of temperature. All results here are obtained for these low temperatures.

This Hamiltonian can be diagonalized by the Bogoliubov transformation

$$\widehat{a}_𝒌=\frac{\widehat{c}_𝒌+\gamma _𝒌\widehat{c}_𝒌^{}}{\sqrt{1\gamma _𝒌^2}}+g_𝒌.$$
(3)
Then the coefficients $`\gamma _𝒌`$, $`g_𝒌`$ and the quasiparticle spectrum $`\omega _𝒌`$ are given by

$`\gamma _𝒌`$ $`=`$ $`\xi 1+\sqrt{\xi (\xi +2)},`$ (4a)
$`g_𝒌`$ $`=`$ $`\sqrt{{\displaystyle \frac{n_0}{V}}}{\displaystyle \frac{U_𝒌}{(\xi +2)v_0n_0}},`$ (4b)
$`\omega _𝒌`$ $`=`$ $`v_0n_0\sqrt{\xi (\xi +2)},`$ (4c)
$`\xi `$ $`=`$ $`{\displaystyle \frac{\mathrm{}^2k^2}{2mv_0n_0}}+\mathrm{\Delta },\mathrm{\Delta }={\displaystyle \frac{v_0n_0\mu }{v_0n_0}}.`$ (4d)
Next, we take an ensemble-average to quench the random potential. The random potential simulates Vycor glass with a characteristic pore size $`r_p`$ as follows. The quenched potential $`U_𝒌`$ may decay above the characteristic wave number $`k_p=2\pi /r_p`$. Thus, we assume the averaged potential

$$\frac{1}{V}U_𝒌U_𝒌_{\mathrm{av}}=R_0\mathrm{exp}\left[\frac{k^2}{2k_p^2}\right],$$
(5)
where av denotes the ensemble-average. $`R_0`$, with dimension $`(energy)^2(length)^3`$, is the characteristic strength of the random potential. Equation (5) makes our model completely different from Huang and Meng’s, and we will show that the results are also different. The coherence length of the BEC is thought to be from hundreds to thousands of Å, whereas the spatial scale of disorder is the pore size in the glass, which is dozens of Å. Hence, the macroscopic wave function of BEC is not sensitive to disorder in and between pores but instead depend on the disorder averaged over the coherence length. Hence, the ensemble-averaged system can become nearly uniform. For a uniform Bose system, it has been proven that the elementary excitation spectrum becomes the gapless Goldstone mode Stringari . Thus, we set $`\mathrm{\Delta }=0`$ in Eq. (4d).

The resultant diagonalized and ensemble-averaged Hamiltonian is

$`\widehat{K}_{\mathrm{eff}}`$ $`=`$ $`V(\mu n_0+{\displaystyle \frac{n}{V}}U_0+ϵ_1+ϵ_R)`$ (6a)
$`+{\displaystyle \underset{𝒌0}{}}\mathrm{}\omega _k\widehat{c}_𝒌^{}\widehat{c}_𝒌,`$
$`\omega _k`$ $`=`$ $`{\displaystyle \frac{\mathrm{}}{2m}}k\sqrt{k^2+16\pi an_0},`$ (6b)
$`ϵ_1`$ $`=`$ $`{\displaystyle \frac{2\pi an_0^2\mathrm{}^2}{m}}\left[1+{\displaystyle \frac{128}{15\sqrt{\pi }}}\sqrt{n_0a^3}\right],`$ (6c)
$`ϵ_R`$ $`=`$ $`{\displaystyle \frac{m\sqrt{a}n_0^{3/2}R_0}{4\sqrt{\pi }\mathrm{}^2}}[\mathrm{e}^{2\alpha }(5+4\alpha )\{1\mathrm{erf}(\sqrt{2\alpha }\}`$ (6d)
$`+\sqrt{{\displaystyle \frac{2}{\pi \alpha }}}(1+\alpha )],`$
$`\alpha `$ $`=`$ $`{\displaystyle \frac{4\pi an_0}{k_p^2}},`$ (6e)
where $`ϵ_1`$ is the hard sphere interaction energy at 0 K, similarly, $`ϵ_R`$ is that for the random potential. The quasiparticle spectrum $`\omega _k`$ is the same as that in the hard sphere Bose gas model Huang2  and is independent of the random potential. This independence is confirmed by neutron scattering experiments Plantevin , which justifies the above assumption of $`\mathrm{\Delta }=0`$; conversely, if $`\mathrm{\Delta }0`$, the spectrum would depend on the random potential.

This Hamiltonian enables us to obtain the grand partition function $`Q=\mathrm{Tr}\{\mathrm{exp}(\beta \widehat{K})\}`$ and various physical quantities. The condensate density is defined by the following relation:

$$n_0=n\frac{1}{V}\underset{𝒌0}{}\widehat{a}_𝒌^{}\widehat{a}_𝒌,$$
(7)
where $`n`$ is the particle number density. The second term represents the noncondensate particle number as

$`{\displaystyle \frac{1}{V}}{\displaystyle \underset{𝒌0}{}}\widehat{a}_𝒌^{}\widehat{a}_𝒌=n_1+n_R,`$ (8a)
$`n_1`$ $`=`$ $`{\displaystyle \frac{8}{3\sqrt{\pi }}}(n_0a)^{3/2}`$ (8b)
$`+{\displaystyle \frac{4}{\sqrt{\pi }\lambda ^3}}{\displaystyle _0^{\mathrm{}}}𝑑t{\displaystyle \frac{t(t^2+\theta /2)}{\sqrt{t^2+\theta }\{\mathrm{e}^{t\sqrt{t^2+\theta }}1\}}},`$
$`n_R`$ $`=`$ $`{\displaystyle \frac{m^2R_0}{8\pi ^{3/2}\mathrm{}^4}}\sqrt{{\displaystyle \frac{n_0}{a}}}[\mathrm{e}^{2\alpha }(1+4\alpha )\{1\mathrm{erf}(\sqrt{2\alpha })\}`$ (8c)
$`2\sqrt{{\displaystyle \frac{2\alpha }{\pi }}}],`$
$`\lambda `$ $`=`$ $`\sqrt{{\displaystyle \frac{2\pi \beta \mathrm{}^2}{m}}},\theta ={\displaystyle \frac{8\pi a\mathrm{}^2\beta n_0}{m}},t^2={\displaystyle \frac{\mathrm{}^2\beta }{2m}}k^2.`$ (8d)
Here $`n_1`$ is the noncondensate density excited by the hard sphere interaction, $`n_R`$ is the density due to the scattering of condensate particles with the random potential, and $`\lambda `$ is the thermal de Broglie wave length. When $`a`$ vanishes, $`n_R`$ becomes infinite. This means that the system would collapse if there were no repulsive interactions between particles.

Because superfluidity is described by the two fluid model, the particle density $`n`$ consists of the normal fluid density $`n_n`$ and the superfluid density $`n_s`$. The superfluid density $`n_s`$ can be calculated by linear response theory Hohenberg . Because of its viscosity, only the normal fluid responds to a small, applied velocity field. Thus the normal fluid density can be defined by the response of the momentum density $`j_i(𝒙,t)`$ to the external velocity field $`v_i(𝒙,t)`$. Linear response theory gives the following relation:

$`j_i(𝒙,t)`$ $`=`$ $`\chi _{ij}(𝒙,t)v_j(𝒙,t),`$ (9a)
$`\chi _{ij}(𝒙,t)`$ $`=`$ $`[j_i(𝒙,t),j_j(0,0)],`$ (9b)
$`j_i(𝒙,t)`$ $`=`$ $`{\displaystyle \frac{\mathrm{}}{2i}}\{\widehat{\mathrm{\Psi }}^{}(𝒙,t){\displaystyle \frac{\widehat{\mathrm{\Psi }}(𝒙,t)}{x_i}}`$ (9c)
$`{\displaystyle \frac{\widehat{\mathrm{\Psi }}^{}(𝒙,t)}{x_i}}\widehat{\mathrm{\Psi }}(𝒙,t)\},`$
$`\widehat{\mathrm{\Psi }}(𝒙,t)`$ $`=`$ $`\mathrm{e}^{\mathrm{i}(\widehat{H}\mu \widehat{N})t/\mathrm{}}\widehat{\mathrm{\Psi }}(𝒙)\mathrm{e}^{\mathrm{i}(\widehat{H}\mu \widehat{N})t/\mathrm{}},`$ (9d)
where $`\widehat{\mathrm{\Psi }}(𝒙,t)`$ is the Heisenberg field operator. The static susceptibility $`\chi _{ij}(𝒌)`$ is defined as

$`\chi _{ij}(𝒙,t)`$ $`=`$ $`{\displaystyle \frac{d\omega }{2\pi }\frac{d^3k}{(2\pi )^3}\mathrm{e}^{\mathrm{i}\omega t}\mathrm{e}^{\mathrm{i}𝒌𝒙}\chi _{ij}(𝒌,\omega )},`$ (10a)
$`\chi _{ij}(𝒌)`$ $`=`$ $`\underset{\omega 0}{lim}\chi _{ij}(𝒌,\omega ).`$ (10b)
Because of the rotational invariance, the static susceptibility $`\chi _{ij}(𝒌)`$ can be written

$$\chi _{ij}(𝒌)=\frac{k_ik_j}{k^2}A(𝒌)+\left(\delta _{ij}\frac{k_ik_j}{k^2}\right)B(𝒌),$$
(11)
where $`A(𝒌)`$ and $`B(𝒌)`$ are the longitudinal and transverse parts, respectively. The transverse susceptibility $`B(0)`$ is the normal fluid mass density. The superfluid number density $`n_s`$ is $`nB(0)/m`$. The susceptibility $`B(0)`$ can be calculated by the Bogoliubov transformation in Eq. (3). After some tedious calculations, the resultant superfluid density is given by

$`n_s=nn_{n1}n_{nR},`$ (12a)
$`n_{n1}={\displaystyle \frac{8}{3\sqrt{\pi }\lambda ^3}}{\displaystyle _0^{\mathrm{}}}𝑑t{\displaystyle \frac{t^4\mathrm{e}^{t\sqrt{t^2+\theta }}}{(1\mathrm{e}^{t\sqrt{t^2+\theta }})^2}},`$ (12b)
$`n_{nR}={\displaystyle \frac{4}{3}}n_R,`$ (12c)
where $`n_{n1}`$ is the normal fluid density due to the elementary excitations, and $`n_{nR}`$ is that due to scattering with the random potential. The density $`n_{n1}`$ can be also obtained using Khalatnikov’s method that is based on Galilean invariance Khalatnikov . The relation $`n_{nR}=4/3n_R=n_R+1/3n_R`$ shows that the random potential causes the larger normal fluid density than the noncondensate density; some condensate particles are captured by the random potential to participate in the normal fluid. This makes it possible to destroy superfluidity even at 0 K when $`n_{nR}`$ becomes comparable to $`n`$.

This formulation can be used to obtain various physical quantities including the condensate density, the superfluid density, and the specific heat.

## III Comparison with experiments

In this section, we compare the calculated specific heat and the superfluid density with experimental results. Quantitative agreement is shown to be good at low temperatures, which supports our assumption of a dilute Bose gas in the random potential. Furthermore, we show that the random potential leads to as-yet unobserved behavior of the specific heat.

To make a quantitative comparison, we give the following numerical values to the parameters: $`m6.6\times 10^{27}\mathrm{kg}`$ and $`a5\times 10^{10}\mathrm{m}`$ are the mass and the s-wave scattering length of a 4He atom. Other parameters are from the experiments of Reppy et al. Reppy . The volume $`V`$ of open pores in the Vycor glass (about 40% of the total volume of the Vycor glass) is about 1 $`\mathrm{cm}^3`$. The particle density $`n`$ of 4He inside the Vycor glass is estimated as follows. In Vycor glass, the atoms are adsorbed and fully cover the surfaces of the open pores due to the van der waals attraction. The pore area is about $`10^8`$ $`\mathrm{m}^2/\mathrm{m}^3`$. The rest of atoms, which do not participate in the first-layer solid, can behave as a dilute gas inside the pores. The particle density $`n`$ of the dilute gas is obtained by subtracting the adsorbed amount from the total amount. This density is estimated to be from 0.001% to 70% of the density of bulk liquid 4He $`n_{bulk}2.1\times 10^{28}/\mathrm{m}^3`$. Because the first layer of 4He adsorbed on the surfaces cannot move and behaves as a solid, we assume the pore size 30Å of Vycor glass is effectively reduced by $`2a`$. Thus, $`r_p`$ is estimated to be 20Å. The last parameter $`R_0`$, which is the strength of the random potential, can be fixed by comparing to experiment. Shown in Fig. 1(a) are the data of zero temperature superfluid signals taken in an experiment that used a torsional oscillator (Fig. 12 of reference Reppy ). Because the superfluid component does not contribute to the moment of inertia, the resonant frequency and the period of oscillation differ from those without superfluid. The period difference $`\mathrm{\Delta }P`$ is approximately proportional to the superfluid component.

Here, the superfluid density is nearly proportional to $`\mathrm{\Delta }P`$ and disappears at a coverage of 17.5 mg. Figure 1(b) shows the superfluid density at 0 K from Eq. (12). As in the experiment, the superfluid density becomes zero at a certain coverage that depends on $`R_0`$. Thus, the value of $`R_0`$ can be fixed using the comparison with Fig. 1(a); i.e., $`R_0=5\times 10^{75}\mathrm{J}^2\mathrm{m}^3`$. Here, we define $`R_w\sqrt{R_0n}`$, which is the single particle energy converted from $`R_0`$. In the Vycor glass, $`R_w/k_B`$ is about $`0.0011`$ K. Just above the critical coverage, the superfluid density increases linearly for both the experiment and the calculation; however, their slopes cannot be compared because the amplitude of $`n_s`$ is unknown in the experiment.

Because all parameters are now fixed, we will quantitatively compare calculations to experiments. The specific heat can be obtained from temperature differentiation of the free energy:

$`\mathrm{\Omega }`$ $`=`$ $`\mathrm{\Omega }_1+\mathrm{\Omega }_R,`$ (13a)
$`\mathrm{\Omega }_1`$ $`=`$ $`V(\mu n_0+ϵ_1)`$ (13b)
$`+{\displaystyle \frac{4V}{\sqrt{\pi }\beta \lambda ^3}}{\displaystyle _0^{\mathrm{}}}𝑑t\{t^2\mathrm{log}[1\mathrm{e}^{t\sqrt{t^2+\theta }}]\},`$
$`\mathrm{\Omega }_R`$ $`=`$ $`V(ϵ_R+{\displaystyle \frac{n}{V}}U_0).`$ (13c)
where $`\mathrm{\Omega }_1`$ is the free energy of the elementary excitation and the hard sphere interaction, and $`\mathrm{\Omega }_R`$ is the free energy from the random potential. Figure 2 compares our results to data of low temperature specific heat taken from Fig. 1 of reference Reppy . In Fig. 2(a), which shows the data for high density, the density $`n`$ is fixed from the experimental coverage, whereas we fix the density from the superfluid critical temperature in Fig. 2(b) at low density. This is because we have no information about the data on the coverage. Figure 2(b) also shows the superfluid density.

The theoretical results agree quantitatively with experiment without using free parameters. Above 1.0 K in Fig. 2(a), the calculated condensate density begins to decrease rapidly; here, our criterion of constant condensate density fails, which likely causes the discrepancy with experiment. However, Fig. 2(b) shows that the calculated specific heat agrees with experiment up to temperatures near the superfluid critical temperature; in this temperature region, the calculated condensate density hardly decreases. This means that the system is more dilute than that of Fig. 2(a) and thus is affected by the random potential rather than the elementary excitations. We discuss this in the next section. Nevertheless, these comparisons show that our model is accurate at low temperatures.

Our model predicts an effect at low temperature that is due to the random potential but has not yet been observed. This is shown in Fig. 3, which is the log-log plot of Fig. 2(a). When the system is free from the random potential, the specific heat should increase a $`T^3`$ because of the contribution from the phonons. However, with a random potential, the dependence is linear in $`T`$ at low temperatures. This means that at low temperatures, the contribution from the random potential is larger than that from the elementary excitations (phonons). The free energy $`\mathrm{\Omega }_R`$ in Eq. (13) depends on the temperature only through the condensate density $`n_0`$. $`\mathrm{\Omega }_R`$ is the energy from the scattering of the condensate particles with the random potential, and the resultant specific heat is given by the energy that the condensate particles need to slip out of the random potential. Experimental observation of this $`T`$-linear dependence might clearly identify the influence from the random potential.

## IV Condensate density and superfluid density

This section describes some characteristic behavior of the condensate density $`n_0`$ and the superfluid density $`n_s`$ derived from our model.

Figure. 4 shows the dependence of $`n_0`$ and $`n_s`$ on temperature and density. Both $`n_0`$ and $`n_s`$ decrease with decreasing density, even at 0 K. This means that the effect of the random potential on $`n_0`$ and $`n_s`$ becomes larger as the density is reduced. Figure. 4 shows clearly the difference between $`n_s`$ and $`n_0`$. Below the critical density, the superfluid density disappears, although the condensate persists.

This situation indicates that the condensate particles cannot move as a superfluid because they are trapped by the random potential. We expect that this theoretical result will be confirmed by measurements of condensate density.

Figure. 5 shows the temperature dependence of the superfluid density just before superfluidity disappears. This figure shows the reentrant transition at which the superfluid density $`n_s`$ goes to zero with a decrease of temperature. In this temperature region, the condensate density is almost constant; hence, our formulation should work well in accordance with the criterion described in Sec. II. This reentrant transition is understood as follows.

The condensate depletion $`n_R`$ of Eq. (8c) comes from the scattering of the condensate particles with the random potential. This decreases with $`n_0`$ as the temperature rises, so that $`n_{nR}=4/3n_R`$ of Eq. (12c) also decreases. As shown in Fig. 6, the magnitude of this decrease in $`n_{nR}`$ exceeds the increase in $`n_{n1}`$ of Eq. (12b), which is the normal fluid density due to the elementary excitations, in the very low temperature region $`T0.5`$ mK. In other words, condensate particles that are trapped by the random potential at lower temperatures can escape at higher temperatures and thus participate in superfluidity.

This reentrant transition has not been observed experimentally, probably because it should only occur at very low densities and low temperatures. However, large values of $`R_0`$ can make the reentrant transition observable as follows. We define the temperature $`T_{max}`$ as that which maximizes the superfluid density $`n_s`$, and define $`\mathrm{\Delta }n_sn_s(T=T_{max})n_s(T=0)`$. Fig. 7 shows the density dependence of $`T_{max}`$ and $`\mathrm{\Delta }n_s`$. Both variables increase with $`R_0`$.

Therefore, the reentrant transition is more likely to be measured at larger $`R_0`$. The parameter $`R_0`$ is the strength of the random potential over the entire space, and one way to increase $`R_0`$ is to decrease the open pore density of the Vycor glass.

## V Conclusions

The present paper describes the dilute Bose gas system in a random potential. The outcomes of our studies are as follows.

By including the pore size dependence of Vycor glass in the random potential, our model could closely match the experimental conditions of liquid 4He in Vycor glass. We fixed the strength of the random potential by equating the theoretical and experimental critical coverages below which the superfluid density at 0 K vanishes. No other parameters could be adjusted, and thus we could quantitatively compare theory to experiment for other physical quantities.

First, we showed that the calculated specific heat for Vycor glass quantitatively agrees with measurements. This indicates that liquid 4He in Vycor glass behaves as a dilute Bose gas in a random potential. For low temperatures, the calculated specific heat was linear in $`T`$ because of the random potential. Second, the BEC was shown to persist even when superfluidity disappears below the critical density. Finally, we showed that the reentrant transition of the superfluid phase is more likely to be observed experimentally by increasing the strength of the random potential.

Because we neglected interactions between pairs of excited particles and between excited particles and the random potential, this model does not apply to systems at high temperatures. To overcome this limitation, we are improving the model to include these interactions and will report on this more general model in the near future.

###### Acknowledgements.

MT acknowledges support by a Grant-in-Aid for Scientific Research (Grant No.12640357) by Japan Society for the Promotion of Science.


format: text

text:
# On H. Friedrich’s formulation of the Einstein equations with fluid sources.

## 1 Introduction.

Cattaneo \[C\] and Ferrarese \[F1\], have stressed long ago that the fundamental physical data are those defined on time lines instead of those defined on spacelike hypersurfaces. In the case of a fluid, timelike world lines have a material reality given by the flow, the unknowns associated to them correspond to a lagrangian picture, in contrast with the usual eulerian one. In many situations, particularly in cosmology, the lagrangian description seems more fundamental. However, one cannot be too confident because General Relativity holds many surprises. For example, consider a fluid in internal equilibrium and in equilibrium with a stationary black hole. Under some conditions (cf. \[BY\]) the thermodynamic temperature is found by certain eulerian observers whose world lines have vanishing vorticity. This corresponds to our ”well posedness” result given in section 11.4. We believe that eulerian and lagrangian descriptions are both usefull.

Ferrarese \[F2\], \[F3\] has given a complete 3+1 decomposition of the connection and Riemann tensor in a frame whose time axis is tangent to the time lines and the space axes orthogonal to it. Using this decomposition he has written the Einstein- dust equations as a (infinite dimensional) dynamical system. Recently H. Friedrich \[HF\] has used a lagangian description of the flow together with a tetrad formalism to write the Einstein equations coupled with a perfect fluid as a symmetric first order system, formally hyperbolic. His treatment includes dust as a particular case.

The Einstein- dust system, as well as the Einstein- Euler system for a perfect fluid, have been proved long ago \[C.B1\] to form a well posed hyperbolic system in the sense of Leray. The relativistic Euler system has been put in first order symmetric hyperbolic form (FOSH) first by K.O. Friedrichs himself \[KOF\] , then by general methods by Anile, Boillat and Ruggeri (see references in \[A\]). Rendall \[R\] has obtained, for well chosen equations of state, FOSH systems for the Euler equations which extend to vacuum, but do not yet permit, as he points out himself, the solution of realistic problems of motion of compact fluid bodies in otherwise empty space. The system obtained by H. Friedrich can be important to solve such realistic problems.

In this article we establish a variant of Friedrich’s system for a lagrangian and tetrad description of the flow and the spacetime geometry. Contrary to Friedrich we will not work with the Weyl tensor but with the Riemann tensor, as in \[AC.BY\]. It avoids the introduction of unphysical characteristics and seems to us more natural in this non vacuum case as well as somewhat simpler. The systems we obtain for vacuum, dust or perfect fluid are symmetric and hyperbolic in the same sense as Friedrich’s system. They are hyperbolic in K.O. Friedrichs’s sense, which leads to existence theorems for solutions in Sobolev spaces, only if the time lines admit global spacelike sections, and the initial data are given on such a section. These existence theorems are needed even if one wants to consider the evolution equations as a well posed dynamical system.

## 2 Metric and coframe.

We write the metric in an orthonormal frame, i.e.,

$$g=(\theta ^0)^2+\underset{i=1}{\overset{3}{}}(\theta ^i)^2.$$
(2.1)
We choose the time axis to be tangent to the time lines, i.e. the cobasis $`\theta `$ is such that $`\theta ^i`$ does not contain $`dx^0.`$ We set

$$\theta ^i=a_j^idx^j,\text{ }\theta ^0=Udx^0+b_idx^i.$$
(2.2)
We will call such a frame a C.F. (Cattaneo- Ferrarese) frame. We denote by $`A_i^j`$ the matrix inverse of $`a_i^j`$. It holds that

$$\text{ }dx^i=A_j^i\theta ^j,\text{ }dx^0=U^1(\theta ^0A_j^ib_i\theta ^j)\text{ }.$$

The Pfaff derivatives $`_\alpha `$ in the C.F. frame are linked to the partial derivatives $`/x^\alpha `$ by the relations

$$\frac{}{x^0}=U_0,\text{ }\frac{}{x^i}=a_i^j_j+b_i_0.$$

$$_0=U^1\frac{}{x^0},\text{ }_i=A_i^j[\frac{}{x^j}U^1b_j\frac{}{x^0}].$$

## 3  Bianchi equations

Whatever the coframe, the components of the Riemann tensor satisfy the identities

$$_\alpha R_{\beta \gamma ,\lambda \mu }+_\beta R_{\gamma \alpha ,\lambda \mu }+_\gamma R_{\alpha \beta ,\lambda \mu }0,$$
(3.1)
hence, if the Ricci tensor $`R_{\alpha \beta }`$ satisfies the Einstein equations

$$R_{\alpha \beta }=\rho _{\alpha \beta },$$
(3.2)
it holds that

$$_\alpha R_{\beta \gamma ::\mu ,}^{\mathrm{}..\alpha }=_\beta \rho _{\gamma \mu }_\gamma \rho _{\beta \mu }$$
(3.3)
Equations 3.1 and 3.3 imply the following ones (cf. an analogous system in \[AC.BY\])

$$_0R_{hk,0j}+_kR_{0h,0j}_hR_{0k,0j}=0,$$
(3.4)
and, in using a symmetry of the Riemann tensor,

$$_0R_{:::i,0j}^0+_hR_{:::i,0j}^h=J_{i,0j}_0\rho _{ji}_j\rho _{0i}.$$
(3.5)
The equations 3.4 and 3.5 are for each given pair ($`0j)`$ a first order system for the components $`R_{hk,0j\text{ }}`$ and $`R_{0h,0j}.`$ The principal operator is a symmetric 6 by 6 matrix:

$`\left(\begin{array}{cccccc}_0& \text{0}& \text{0}& _2& \text{-}_1& \text{0}\\ \text{0}& _0& \text{0}& \text{0}& _3& \text{-}_2\\ \text{0}& \text{0}& _0& \text{-}_3& \text{0}& _1\\ _2& \text{0}& \text{-}_3& _0& \text{0}& \text{0}\\ \text{-}_1& _3& \text{0}& \text{0}& _0& \text{0}\\ \text{0}& \text{-}_2& _1& \text{0}& \text{0}& _0\end{array}\right)`$

Analogous equations and result hold with the pair $`(0j)`$ replaced by $`(lm)`$, $`l<m`$. The system finally obtained has a principal matrix consisting of 6 identical 6 by 6 blocks around the diagonal, it is symmetric.

However, the system obtained for the components of the Riemann tensor in the coframe with time axis $`\theta ^0`$ and space axis $`\theta ^i,`$ with principal operator the matrix $`M^\alpha _\alpha ,`$ cannot be said to be hyperbolic in the usual sense: the principal matrix $`M^0`$ is the unit matrix hence positive definite, but the operators $`_\alpha `$ are not partial derivatives. We say that the system is a quasi FOSH (First Order Symmetric Hyperbolic) system (Cf. section 11).

###### Remark 1

One can associate with the Riemann tensor, as in the usual case of the 3+1 splitting, two pairs of “electric” and “magnetic” 2- tensors whose non zero components in the frame are:

$`E_{ij}`$ $`R_{0i,0j},\text{ }D_{ij}{\displaystyle \frac{1}{4}}\epsilon _{ihk}\epsilon _{jlm}R^{hk,lm}`$
$`H_{ij}`$ $`{\displaystyle \frac{1}{2}}\epsilon _{ihk}R_{:::::,0j}^{hk},\text{ }B_{ji}{\displaystyle \frac{1}{2}}\epsilon _{ihk}R_{0j}^{:::::hk},`$
where $`\epsilon _{ijk}`$ is the totally antisymmetric Kronecker tensor. These definitions imply here

$$R_{hk,0j}\epsilon _{::hk}^iH_{ij},\text{ }R_{hk,lm}\epsilon _{::hk}^i\epsilon _{::lm}^jD_{ij},\text{ }R_{0j,hk}=\epsilon _{::hk}^iB_{ji}.$$

The previous system is a quasi FOSH system for these pairs of tensors.

###### Remark 2

The system 3.1, 3.2 contains also the equations, which we do not use for evolution,

$$_hR_{ij,\lambda \mu }+_jR_{hi,\lambda \mu }+_iR_{jh,\lambda \mu }=0$$
(3.6)
and

$$_hR_{\mathrm{..0},\lambda \mu }^h=_\lambda \rho _{\mu 0}_\mu \rho _{\lambda 0}.$$
(3.7)
These equations are not usual constraints on a submanifold $`t=cons\mathrm{tan}t`$ because $`_i`$ contains the transversal derivative $`/t.`$ We call them Bianchi quasi-constraints.

## 4 Coframe structure coefficients.

The structure coefficients $`c`$ of the coframe are defined by

$$d\theta ^\alpha \frac{1}{2}c_{\beta \gamma }^\alpha \theta ^\beta \theta ^\gamma .$$

We have

$$d\theta ^0=_iU\theta ^idx^0+_\alpha b\theta ^\alpha dx^i$$

that is,

$$d\theta ^0=U^1_iU\theta ^i(\theta ^0A_j^hb_h\theta ^j)+A_j^h_\alpha b_h\theta ^\alpha \theta ^j,$$

therefore

$$c_{0i}^0=U^1_iUA_i^j_0b_j=c_{i0}^0$$

and, with $`f_{[ij]}f_{ij}f_{ji},`$

$$c_{ij}^0=U^1b_hA_{[j}^h_{i]}UA_{[j}^h_{i]}b_h$$

that is

$$c_{ij}^0=UA_{[i}^h_{j]}(U^1b_h)$$

also

$$d\theta ^i=A_k^j(_ha_j^i\theta ^h\theta ^k+_Oa_j^i\theta ^0\theta ^k)$$

hence

$$c_{k0}^i=A_k^j_0a_j^i$$

$$c_{hk}^i=A_{[h}^j_{k]}a_j^i$$

###### Remark 3

If the time lines are not hypersurface orthogonal (i.e. if $`b_i0)`$ the coefficients $`c_{ij}^h`$ are different from the structure coefficients of the space frame $`\theta ^i.`$

## 5 Connection.

The connection is such that

$$\omega _{\beta \gamma }^\alpha \omega _{\gamma \beta }^\alpha =c_{\beta \gamma }^\alpha $$

and

$$\omega _{\gamma \beta ,\lambda }+\omega _{\gamma \lambda ,\beta }=0$$

with ($`\eta `$ is the Minkowski metric)

$$\omega _{\gamma \beta ,\lambda }\eta _{\alpha \lambda }\omega _{\gamma \beta }^\alpha $$

We set

$$c_{\gamma \beta ,\lambda }\eta _{\alpha \lambda }c_{\gamma \beta }^\alpha ,$$

and we have (Cf. C.B-DeWitt I p.308$`)`$

$$\omega _{\alpha \lambda ,\mu }=\frac{1}{2}(c_{\alpha \lambda ,\mu }c_{\alpha \mu ,\lambda }c_{\lambda \mu ,\alpha }).$$

We find that, as foreseen from antisymmetry

$$\omega _{00}^0=\omega _{i0}^0$$

We set

$$Y_i\omega _{00,i},$$

and we have

$$Y_i\omega _{00,i}=\omega _{0i}^0=\omega _{00}^i=c_{0i,0}=c_{0i}^0=U^1_iUA_i^j_0b_j.$$

We know that $`\omega _{0i}^j=\omega _{0i,j}`$ is antisymmetric in $`i`$ and $`j.`$ We set

$$\omega _{0i,j}f_{ij}=\frac{1}{2}\{A_j^h_0a_h^iA_i^h_0a_i^j+A_{[i}^h_{j]}b_h\}.$$

###### Remark 4

Let $`e_{(\alpha )}_\alpha `$ be the frame dual to $`\theta ^\alpha `$, i.e. with components $`\delta _{(\alpha )}^\lambda .`$ Then

$$_\beta e_{(\alpha )}^\lambda =\omega _{\beta \alpha }^\lambda ,$$

in particular $`\omega _{0i,j}`$ is the projection on $`e_{(j)}`$ of the derivative of $`e_{(i)}`$ in the direction of $`e_{(0)}.`$ We have $`f_{ij}=0`$ if the frame is Fermi transported along the time line.

The connection coefficient $`\omega _{i0,j}`$ is the sum of a term symmetric in $`i`$ and $`j`$ and an antisymmetric one, we have

$$X_{ij}\omega _{i0,j}=\omega _{i0}^j=\frac{1}{2}\{A_j^h_0a_h^i+A_i^h_0a_i^j+A_{[i}^h_{j]}b_h\}.$$

The antisymmetric term vanishes if the time lines are hypersurface orthogonal ($`b_i=0).`$

The coefficients $`\omega _{ij}^h`$ are also linear expressions in terms of the first derivatives of the frame coefficients, they are identical to the connection constructed with the $`a_i^j`$ at fixed $`x^0`$ if $`b_i=0.`$

## 6 Frame evolution.

The relations between the connection coefficients $`\omega _{ij}^h`$ and the frame coefficients $`a_i^j`$ give equations for the dragging of these coefficients along the time lines when the connection is known. Indeed we have found

$$_0a_j^i=a_j^kc_{0k}^i=2a_j^k(\omega _{0k}^i\omega _{k0}^i)2a_j^k(f_k^iX_k^i)$$
(6.1)
$$_0b_i=a_i^h(Y_h+U^1_hU$$
(6.2)
No equation gives the evolution of $`U`$. It can be considered as a gauge variable fixing the time parameter. H. Friedrich chooses the coordinate $`x^0`$ to be the proper time along the given time lines, i.e. $`U=1`$

The quantity $`f_j^i`$ is also a gauge variable fixing the evolution of the space frame. If we choose it to be Fermi transported like Friedrich, then $`f_j^i=0.`$

## 7 Curvature.

### 7.1 Evolution of connection.

Using the general formula (cf. C.B- DeWitt I p. 306, with opposite sign convention)

$$R_{\lambda \mu \mathrm{}.\beta }^\mathrm{}\mathrm{}\alpha _\lambda \omega _{\mu \beta }^\alpha _\mu \omega _{\lambda \beta }^\alpha \omega _{\rho \beta }^\alpha (\omega _{\lambda \mu }^\rho \omega _{\mu \lambda }^\rho )+\omega _{\lambda \rho }^\alpha \omega _{\mu \beta }^\rho \omega _{\mu \rho }^\alpha \omega _{\lambda \beta }^\rho $$
(7.1)
$$R_{0h..j}^{\mathrm{}..i}_0\omega _{hj}^i\stackrel{~}{}_hf_j^iY_hf_j^i(f_h^kX_h^{..k})\omega _{kj}^i+Y^iX_{hj}^{..}Y_jX_h^i$$
(7.2)
where $`\stackrel{~}{}`$ is the pseudo covariant derivative constructed with $`_i`$ and $`\omega _{ij}^h`$ (Cattaneo- Ferrarese transversal derivative). In the Fermi gauge $`R_{0h..j}^{\mathrm{}..i}`$ reduces to

$$R_{0h..j}^{\mathrm{}..i}_0\omega _{hj}^i+X_h^{..k}\omega _{kj}^i+Y_iX_h^\mathrm{}jY_jX_h{}_{}{}^{i}.$$
(7.3)
The principal operator in these equations is the dragging along the time lines of $`\omega _{ij}^h.`$

We have:

$$R_{h\mathrm{0..0}}^{\mathrm{}..i}\stackrel{~}{}_hY_i_0X_h^{..i}+Y^iY_hX_h^{..j}X_j^{..i}+f_h^{..j}X_j^{..i}f_j^{..i}X_h^{..j}.$$
(7.4)
We deduce from this identity, in the Fermi gauge,

$$R_{00}\stackrel{~}{}_iY^i_0X_i^i+Y^iY_iX_h^{..j}X_j^{..h}.$$
(7.5)
In the next sections we will write other equations to achieve the determination of both $`Y`$ and $`X.`$

### 7.2 Quasi-constraints.

The other components of the Riemann tensor are given in the C-F frame by:

$$R_{hk..j}^{\mathrm{}..i}\stackrel{~}{R}_{hk\mathrm{}j}^\mathrm{}\mathrm{}^if_j^i(X_{hk}X_{kh})+X_{.k}^{.i}X_{jh}X_{jk}^.X_{..h}^i,$$
(7.6)
where $`\stackrel{~}{R}_{hk..j}^{\mathrm{}..^i}`$ denotes the expression formally constructed as a Riemann tensor with the coefficients $`\omega _{ij}^h.`$

$$R_{kh\mathrm{}.j}^\mathrm{}\mathrm{}0_k\omega _{hj}^0_h\omega _{kj}^0\omega _{\rho j}^0(\omega _{kh}^\rho \omega _{hk}^\rho )+\omega _{k\rho }^0\omega _{hj}^\rho \omega _{h\rho }^0\omega _{kj}^\rho ,$$
(7.7)
that is ,with previous notations,

$$R_{kh\mathrm{}.j}^\mathrm{}\mathrm{}0\stackrel{~}{}_kX_{hj}\stackrel{~}{}_hX_{kj}Y_j(X_{kh}X_{hk}).$$
(7.8)
These equations do not enter in the evolution system for the connection that we are considering, they do not contain the derivative $`_0,`$ we call them connection quasi- constraints.

###### Remark 5

Modulo the expression of the connection in terms of frame coefficients one has the well known symmetry:

$$R_{kh,0j}=R_{0j,kh}$$

One deduces from the identity 7.8

$$R_{h0}\stackrel{~}{}_jX_h^\mathrm{}j\stackrel{~}{}_hX_j^jY^j(X_{jh}X_{hj})$$

## 8 Vacuum case.

In the vacuum case there are no a-priori given time lines. We can give arbitrarily on the spacetime the projection $`Y_i`$ of $`_{e_{(0)}}e_{(0)}`$ on $`e_{(i)}.`$ The quantities $`U`$ and $`f_i^j`$ being also arbitrarily given the Bianchi equations together with 6-1, 6-2; 7-3; 7-5 constitute a quasi FOSH system for the unknowns $`E_{ij},`$ $`D_{ij},`$ $`H_{ij},`$ $`B_{ij},`$ $`a_i^j,`$ $`b_i,`$ $`\omega _{ij}^h,`$ $`X_{ij}.`$ Its solution determines the spacetime metric.

## 9 Perfect fluid.

### 9.1 Fluid equations.

The stress energy tensor of a perfect fluid is

$$T_{\alpha \beta }=(\mu +p)u_\alpha u_\beta +pg_{\alpha \beta }.$$

Then

$$\rho _{\alpha \beta }=(\mu +p)u_\alpha u_\beta +\frac{1}{2}(\mu p)g_{\alpha \beta }.$$

One supposes that the matter energy density $`\mu `$ is a given function of the pressure $`p`$ and entropy $`S;`$ this last function is conserved along the flow lines

$$u^\alpha _\alpha S=0.$$
(9.1)
The Euler equations of the fluid express the generalized conservation law

$$_\alpha T^{\alpha \beta }=0,$$

and are equivalent to the equations

$$(\mu +p)u^\alpha _\alpha u^\beta +(u^\alpha u^\beta +g^{\alpha \beta })_\alpha p=0,\text{ with }u^\alpha u_\alpha =1$$

and

$$(\mu +p)_\alpha u^\alpha +u^\alpha _\alpha \mu =0$$

In our coframe they reduce to:

$$_0S=0$$
(9.2)
$$(\mu +p)Y_i+_ip=0,\text{ }Y_i\omega _{00}^i$$
(9.3)
$$_0\mu +(\mu +p)X_i^i=0.$$
(9.4)
Using the index $`F`$ of the fluid defined by

$$F(p,S)=\frac{dp}{\mu (p,S)+p},$$

the equation 9.3 reads

$$Y_i=_iF,$$
(9.5)
while 9.4 gives , modulo 9.2,

$$\mu _p^{}_0F+X_i^i=0.$$
(9.6)
The commutation relation between Pfaff derivatives

$$(_0_i_i_0)F=c_{0i}^\alpha _\alpha F$$

implies therefore

$$\mu _p^{}[_0Y_i+Y_i_0F+(f_i^jX_i^j)_jF]_i\mu _p^{}_0F_iX_h^h=0.$$
(9.7)
The use of 9.5 and 9.6 replaces $`_\alpha F`$ by functions of $`Y,X`$ and $`p,S.`$ The derivatives $`_i\mu _p^{}`$ are functions of $`Y,p,S`$ and $`_kS,`$ since

$$_i\mu _p^{}=\mu \mathrm{"}_{p^2}_ip+\mu \mathrm{"}_{pS}_iS.$$

We introduce $`S_k=_kS`$ as new unknowns. They are dragged along the flow lines by the following equation deduced from 9.2:

$$_0S_k=c_{0k}^jS_j(f_k^jX_j^k)S_j.$$
(9.8)
Following H. Friedrich we replace in 9.7 $`_iX_h^h`$ by its expression deduced from the equation

$$R_{i0}\stackrel{~}{}_hX_i^\mathrm{}h\stackrel{~}{}_iX_h^hY^h(X_{hi}X_{ih})=0,$$

and we obtain, changing names of indices

$$\mu _p^{}_0Y_h\stackrel{~}{}_jX_h^jY^h(X_{hi}X_{ih})+\mu _p^{}[Y_h_0F+$$

$$(f_h^jX_h^j)_jF]+_h\mu _p^{}_0F=0.$$
(9.9)
In 7.4 we replace $`\stackrel{~}{}_hY_i`$ by $`\stackrel{~}{}_iY_h+c_{hi}^0Y_0,`$ with

$$Y_0_0F(\mu _p^{})^1X_i^i,\text{ }c_{hi}^0X_{ih}X_{hi}$$

and we obtain

$$_0X_h^{.i}\stackrel{~}{}^iY_hc_{hi}^0Y_0Y_hY^i+(X_h^{..j}f_h^{..j})X_j^{..i}+f_j^{..i}X_h^{..j}=R_{h\mathrm{0..0}}^{\mathrm{}..i}.$$
(9.10)
The principal operator on the unknowns $`Y`$ and $`X`$ in the equations 9.9, 9.10 is diagonal by blocks. Each block, corresponding to a given index $`h,`$ is symmetric, it reads:

(
μp0
-1

-2

-3
-1000-2000-3000)fragments
μp0
-1

-2

-3
fragments-1fragments000fragments-20fragments00fragments-300fragments0\left(\begin{tabular}[]{cccc}$\mu_{p}^{\prime}\partial_{0}$&-$\partial_{1}$&-$\partial_{2}$&-$\partial_{3}$\\
-$\partial_{1}$&$\partial_{0}$&0&0\\
-$\partial_{2}$&0&$\partial_{0}$&0\\
-$\partial_{3}$&0&0&$\partial_{0}$\end{tabular}\right)

If $`\mu _p^{}>0`$ the matrix $`M^0`$ is positive definite in the C.F. frame. The system 9.9, 9.10 is a quasi FOSH system for the pairs $`Y_h,X_h^j.`$

###### Remark 6

The characteristic determinant associated with the system 9.9, 9.10, obtained by replacing $`_\alpha `$ with a covariant vector $`\xi `$ is

$$\{\xi _0^2(\mu _p^{}\xi _0^2\underset{i=1,2,3}{}\xi _i^2)\}^3$$

The roots of $`\mu _p^{}\xi _0^2_{i=1,2,3}\xi _i^2=0`$ correspond to sound waves. Their speed is at most 1 (speed of light) if and only if $`\mu _p^{}1.`$

The full system of fluid equations is quasi diagonal by blocks. The matrices $`M^\alpha _\alpha `$ corresponding to 9.2, 9.4, 9.8 reduce to $`_0.`$

### 9.2 Sources of the Bianchi equations.

In our frame the source tensor $`\rho _{\alpha \beta }`$ reduces to

$$\rho _{00}=\frac{1}{2}(\mu +3p),\text{ }\rho _{0i}=0,\text{ }\rho _{ij}=\frac{1}{2}(\mu p)\delta _{ij}$$

We have seen that $`_\alpha p`$ and $`_\alpha \mu `$ are smooth functions of $`p,S,S_i,Y`$ and $`X.`$ The same property holds for $`J_{i,0j}`$ and $`J_{i,hk}.`$

### 9.3 Conclusion.

Assembling the results of the previous sections we find the following theorem.

###### Theorem 7

The EEF (Einstein-Euler-Friedrich) system for a gravitating perfect fluid, with flow lines taken as time lines, $`U`$ and $`f_j^i\omega _{0j}^i`$ given arbitrarily, is a quasi FOSH system for the Riemann curvature tensor, the frame and connection coefficients, the density of matter, the entropy and its space derivatives.

As remarked in section 3 the system is not a usual FOSH system: the time lines are by choice timelike but the hypersurfaces $`x^0=cons\mathrm{tan}t`$ are not necessarily spacelike for the characteristic cone.

## 10 Case of dust.

### 10.1 Matter equations.

The stress energy tensor of a dust source is

$$T_{\alpha \beta }\mu u_\alpha u_\beta .$$

The flow lines, tangent to the unit vector $`u^\alpha ,`$ are then geodesics, it holds that

$$u^\alpha _\alpha u^\beta =0,$$
(10.1)
while the conservation of matter reads:

$$u^\alpha _\alpha \mu +\mu _\alpha u^\alpha =0.$$
(10.2)
We take the flow lines as time lines, i.e. $`u^\alpha =\delta _0^\alpha .`$ Then 10.2 reads

$$_0\mu +\mu X_i^i=0,$$
(10.3)
while 10-1 gives

$$\omega _{00}^iY_i=0.$$
(10.4)
Using this equation we see that, given arbitrarily on the spacetime the gauge variables $`a_0^0`$ the equations (6-1), (6-2) reduce to a quasi- diagonal system with principal operator $`_0`$ for the frame coefficients when the connection (which appears undifferentiated) is known .

### 10.2 Bianchi equations.

For a dust stress energy tensor it holds that

$$\rho _{\alpha \beta }=\mu (u_\alpha u_\beta +\frac{1}{2}g_{\alpha \beta }),$$

hence in the chosen frame

$$\rho _{00}=\frac{1}{2}\mu ,\text{ }\rho _{0i}=0,\text{ }\rho _{ij}=\frac{1}{2}\mu \delta _{ij}$$

A straightforward computation gives then, choosing the Fermi transported frame for simplicity

$$J_{i,0j}_0\rho _{ji}_j\rho _{0i}=\mu (\frac{1}{2}\delta _{ij}X_h^h+X_{ji})$$

$$J_{i,hj}_h\rho _{ji}_j\rho _{hi}=\frac{1}{2}\{_h\mu \delta _{ji}_j\mu \delta _{hi}+\mu (\omega _{hj}^i\omega _{jh}^i)\}$$

There appears a difficulty which is the presence of the space derivatives $`_h\mu `$ in $`J_{i,hj}.`$ We get rid of this problem, inspired (though with a different choice) by H. Friedrich’s treatment, by replacing the unknown $`R_{hk..j}^{\mathrm{}..i}`$ by another unknown with the same symmetries, namely

$$F_{ij,hk.}=R_{ij,hk}+\delta _{ik}R_{jh}\delta _{ih}R_{jk}\delta _{jk}R_{ih}+\delta _{ih}R_{ik}+D_{ij,hk}R$$

with

$$D_{ij,hk}\frac{1}{2}(\delta _{ik}\delta _{jh}\delta _{jk}\delta _{ih})$$

Using this unknown, the contracted Bianchi identities become

$$_0R_{:::i,hj}^0+_kF_{:::i,hj}^k0,$$
(10.5)
while the original ones become

$$_0F_{ij,hk}+_jR_{0i,hk}_iR_{0j,hk}$$

$$_0\{\delta _{ik}R_{jh}\delta _{ih}R_{jk}\delta _{jk}R_{ih}+\delta _{ih}R_{ik}+X_{ijhk}R\}$$
(10.6)
One obtains a quasi FOSH system for $`R_{0i,hk}`$ and $`F_{ij,hk}`$ by replacing on the right hand side of this equation $`R_{\alpha \beta }`$ by $`\rho _{\alpha \beta }`$: now there appears only the time derivative $`_0\mu `$ which, satisfying 10-3, can be eliminated in favor of undifferentiated terms.

Note that $`R_{ij,hk}`$ does not appear in other equations than the Bianchi equations.

### 10.3 Conclusion.

We have proved the following theorem

###### Theorem 8

The EDFF (Einstein-Dust-Ferrarese-Friedrich) system, , with flow lines taken as time lines, $`U`$ and $`f_j^i\omega _{0j}^i`$ given arbitrarily, is a quasi FOSH system for the Riemann curvature tensor, the frame and connection coefficients, the density of matter and its space derivatives.

## 11 Cauchy problem.

### 11.1 Hyperbolicity.

The fact that for the systems that we have obtained the hypersurfaces $`x^0=cons\mathrm{tan}t`$ are not necessarily spacelike for the characteristic cones poses a difficulty for the Cauchy problem, since energy estimates used in proving existence theorems rely on this property.

We have, using standard definitions

###### Proposition 9

A quasi FOSH system, with principal operator the matrix $`M^\alpha _\alpha ,`$ is a FOSH system with $`x^0`$ as a time variable if the matrix of coefficients of $`/x^0`$ is positive definite, namely if

$$\stackrel{~}{M}^0M^\alpha A_\alpha ^0U^1(M^0A_j^ib_iM^j)$$

is positive definite.

#### 11.1.1 Case of dust

###### Lemma 10

The EDFF (Einstein-Dust-Ferrarese-Friedrich) system is a FOSH system relative to $`x^0=cons\mathrm{tan}t`$ slices as long as the quadratic form

$$\overline{g}_{jh}=\underset{i=1,2,3}{}a_j^ia_h^ib_jb_h$$
(11.1)
is positive definite and $`U>0.`$

Proof.

The matrix $`\stackrel{~}{M}_{(bian)}^0`$ of coefficients of $`/x^0`$ corresponding to the Bianchi equations is after multiplication by $`U`$, setting

$$B_iA_i^jb_i,$$

$`\left(\begin{array}{cccccc}1& 0& 0& B_2& B_1& 0\\ 0& 1& 0& 0& B_3& B_2\\ 0& 0& 1& A_3& 0& B_1\\ B_2& 0& B_3& 1& 0& 0\\ B_1& B_3& 0& 0& 1& 0\\ 0& B_2& B_1& 0& 0& 1\end{array}\right)`$

Its eigenvalues are, each with multiplicity 2:

$$1,1+\sqrt{\left(B_2^2+B_3^2+B_1^2\right)},1\sqrt{\left(B_2^2+B_3^2+B_1^2\right)}$$

They are positive, because the given condition implies $`B_2^2+B_3^2+B_1^2<1.`$

The other matrices $`\stackrel{~}{M}_{(conn)}^0`$ and $`\stackrel{~}{M}_{(matter)}^0`$ are unit matrices.

#### 11.1.2 Fluid case.

###### Lemma 11

The EEF (Einstein-Euler-Friedrich) system is a FOSH system relatively to $`x^0=cons^t`$ slices as long as the quadratic form $`\overline{g}_{jh}`$ given in 11.1 is positive definite, $`U>0`$ and $`\mu _p^{}1.`$

Proof. The principal operator is symmetric and composed of blocks around the diagonal.

The matrix $`\stackrel{~}{M}_{(bian)}^0`$ is the same as in the case of dust.

The matrix $`\stackrel{~}{M}_{con}^0`$ corresponding to the connection evolution is, after multiplication by $`U`$

$`\left(\begin{array}{cccc}\mu _p^{}& B_1& B_2& B_3\\ B_1& 1& 0& 0\\ B_2& 0& 1& 0\\ B_3& 0& 0& 1\end{array}\right).`$

Its eigenvalues are $`1,`$ with multiplicity two, and:

$$\frac{1}{2}\{1+\mu _p^{}\pm 2\sqrt{\frac{1}{4}(\mu _p^{}1)^2+B_i^2)}\}$$
(11.2)
These eigenvalues are positive under the hypothesis 11-1 if $`\mu _p^{}1.`$

###### Remark 12

The condition $`\mu _p^{}1`$ expresses that the sound speed is at most equal to the speed of light.

### 11.2 Cauchy data.

The Cauchy problem is less natural in the C.F. formulation than in the usual 3+1 decomposition.

#### 11.2.1 Initial data.

An initial data set on a coordinate patch of a 3 dimensional submanifold is composed of the following elements:

1. A field of coframes $`\overline{a}_j^i^j`$and a field of covariant vectors $`\overline{a}_i^0`$

The quadratic form on $`M`$ with coefficients

$$\overline{g}_{jh}=\underset{i=1,2,3}{}\overline{a}_j^i\overline{a}_h^i\overline{b}_j\overline{b}_h$$
(11.3)
is supposed to be positive definite. The submanifold $`M_0M\times \{0\}`$ of $`M\times R`$ is then spacelike for any lorentzian metric on $`M\times R`$ which reduces on $`M`$ to

$$(\overline{U}dx^0+\overline{b}_idx^i)^2+\underset{i=1,2,3}{}(\overline{\theta }^i)^2,\text{ }\overline{U}>0$$

2. Quantities $`\overline{\omega }_{ij}^h,`$ $`\overline{X}_i^j,\overline{Y}_i.`$

3. Pairs of 2-tensors $`\overline{E}^{ij},\overline{H}_{ij},\overline{D}_{ij},\overline{B}_{ij}.`$

4.a (case of dust) A scalar function $`\overline{\mu }0.`$

4.b (perfect fluid case). Two scalar functions $`\overline{p}0`$ and $`\overline{S}0`$, and equation of state $`\mu (p,S)`$ such that $`\mu _p^{}1`$ for all $`p`$ in a neighbourhood of $`(\overline{p},\overline{S})`$.

#### 11.2.2 Constraints.

The initial data on $`M`$ do not determine on $`M`$ the derivatives $`_i\overline{}_iA_i^j\overline{b}_j_0.`$ To satisfy the quasi-constraints on $`M`$ we must first deduce from the given data and the considered evolution system the values on $`M`$ of the derivatives with respect to $`x^0`$ of the relevant unknowns. The quasiconstraints give then relations on the initial data which we call constraints and suppose to be satisfied on $`M.`$

### 11.3 Existence theorem for the reduced system.

Local existence theorems for a solution of the Cauchy problem are a direct consequence of known results for quasilinear FOSH systems. We enunciate the result in the perfect fluid case. The spaces $`H_s^{u.loc}(M)`$ are the usual uniformly local Sobolev spaces on $`M.`$

###### Theorem 13

Let there be given an initial data set in $`H_s^{u.loc}(M),s>\frac{3}{2}+1.`$ Then for any choice on $`M\times I`$ of the gauge variables in $`C^0(I,H_s^{u.loc}(M)),`$with $`I`$ an interval of $`R`$ containing $`0,`$ there exists a neighbourhood $`I^{}`$ of $`0`$ in $`R`$ such that the EEF equations admit a solution in $`C^0(I^{},H_s^{u.loc}(M))C^1(I^{},H_{s1}^{u.loc}(M)).`$

To prove that the solution satisfies, modulo the initial constraints, the original Einstein- Euler equations would be a messy task, which we will not endeavour. The uniqueness of the obtained solution gives an indirect proof, modulo the existence theorem already known \[C.B1\], \[C.B2\].

### 11.4 Irrotational flows.

We recall the equation 6-2:

$$_0b_i=a_i^h(Y_h+_h\mathrm{log}U),$$
(11.4)
with, in the case of a perfect fluid,

$$Y_h=_hF.$$

We choose then (gauge choice, the geometrical result is independent of it)

$$U=\mathrm{exp}(F).$$

The equation reduces then to

$$_0b_i=0\text{.}$$
(11.5)
Therefore it holds that $`b_i=0`$ througout the motion if it is so initially. It is the well known property of conservation of zero vorticity (irrotationality) of a perfect fluid flow.

For an irrotational flow it holds that $`X_{ij}`$ is a symmetric 2- tensor, equal up to the sign to the extrinsic curvature of the space slices $`x^0=cons\mathrm{tan}t`$ (components in the frame $`\theta ^i)`$, and $`Y_i=U^1_iU,`$ while $`\omega _{ij}^h=\overline{\omega }_{ij}^h`$ are the connection coefficients of the space metric $`\overline{g}=_{i=1,2,3}(\theta ^i)^2,`$ $`\theta ^ia_h^idx^h.`$ The pseudo covariant derivatives $`\stackrel{~}{}_i`$ are covariant derivatives $`\overline{}_i`$ in the space metric.

For an irrotational flow the EFF system is a FOSH system.

The same result holds in the case of dust the appropriate gauge choice is then $`U=cons\mathrm{tan}t.`$

###### Remark 14

The relations between frame and connection, as well as connection and curvature, for hypersurface orthogonal time lines are a special case of those found in the 3+1 decomposition \[C.B3\]. In particular, in the present notations, the Einstein constraints take the familiar form:

$$R_{i0}\overline{}_jX_i^j_iX_j^j=0$$
(11.6)
$$S_{00}\frac{1}{2}(R_{00}+\delta ^{ij}R_{ij})\frac{1}{2}(\overline{R}X_i^{.j}X_j^{.i}+(X_i^i)^2)=\rho _{00}=\frac{1}{2}(\mu +3p).$$
(11.7)
Bibliography.

\[A\] Anile M. ”Relativistic fluids and magneto- fluids” (1989) Cambridge Univ. Press.

\[AC.BY\] Anderson A., Choquet-Bruhat Y., York, J.W. (1997) ”Einstein- Bianchi hyperbolic system for General Relativity” Top. Meth. non lin. anal. 10, 353-373.

\[BY\] Brown D., York J. (1994) Physical origins of time asymmetry, Halliwell, Perez- Mercador, Zurek eds. Cambridge Univ. Press 465-474.

\[C\] Cattaneo C. (1959) Proiezioni naturali e derivate trasversa in una varieta a metrica iperbolica normale Ann. Math. pura e app. (IV) XLVIII , 361-386.

\[C.B1\] (Choquet) Foures -Bruhat Y. (1958) ”Théorèmes d’existence en mécanique des fluides relativistes” Bull. Soc. Math. France, 86 155-175.

\[C.B2\] (Choquet) -Bruhat Y. (1966) ”Etude des équations des fluides chargés relativistes inductifs et conducteurs” Comm. Math. Phys. 3, 334-357.

\[C.B3\] (Choquet) Foures -Bruhat Y. (1956) ”Sur l’intégration des équations de la Relativité Générale” J. Rat. Mech. and Anal. 5, 951-966.

\[F1\] Ferrarese G. (1963) ”Contributi a la tecnica delle proiezioni in una varietà a metrica iperbolica normale” Rendic. matem. 22, 147-168.

\[F2\] Ferrarese G. (1979) Introduzione a la dinamica riemanniana dei sistemi continui, Ed. pitagora, Bologna.

\[F3\] Ferrarese G. and Catani C. (1994) ”Generalized frames of references and intrinsic Cauchy problem in General Relativity” in ”Physics on manifolds”, Flato M., Kerner R., Lichnerowicz A. ed. Reidel, 93-108.

\[HF\] Friedrich H. (1998) ”On the evolution equations for gravitating ideal fluid bodies in General Relativity” Phys. Rev. D 57, 2317-2322.

\[KOF\] Friedrichs K. O., (1974) ”On the laws of relativistic electro- magneto- fluid dynamics Comm. pure app. math. XXVII, 749-808.

\[R\] Rendall A. (1992) ”The initial value problem for a class of general relativistic fluid bodies” J. Math. Phys. 33, 1047-1053.

Yvonne Choquet-Bruhat, Université Paris 6, tour 22-12, 4 place Jussieu 75230 Paris, France.

James W. York, Department of Physics, UNC, Chapel Hill, N.C, U.S.A.


text:
# Transitions to equilibrium state in classical ϕ⁴ lattice

## I introduction

The investigation of statistical behavior of Hamiltonian coupled oscillators from microscopic dynamics has attracted much attention since the original Fermi-Pasta-Ulam (FPU) problem \[1-4\]. The important issue in this aspect is the approach to equilibrium of nonlinear system with many degrees of freedom. It was found in the FPU study that with sufficiently small energy level various modes of the system do not go to equipartition state as conventionally expected for equilibrium . This observation has stimulated many works studying the transitions from nonequilibrium to equilibrium, based on numerical simulations of microscopic dynamics \[5-9\]. Various models of Hamiltonian system with many degrees of freedom have been considered for this purpose, among which the $`\varphi ^4`$ lattice is one of the most extensively investigated.

The model Hamiltonian reads

$`\widehat{H}`$ $`=`$ $`{\displaystyle \underset{i=1,\mathrm{}.N}{}}({\displaystyle \frac{m\widehat{p}_i^2}{2}}+{\displaystyle \frac{\mu (\widehat{q}_i\widehat{q}_{i+1})^2}{2}}+{\displaystyle \frac{\beta \widehat{q}_i^4}{4}})`$ (1)
$`\widehat{q}_{N+1}`$ $`=`$ $`\widehat{q}_1\text{ , }\widehat{p}_{N+1}=\widehat{p}_1`$ (2)
With the rescalings $`q_i=\sqrt{\mu /\beta }\widehat{q}_i,`$ $`p_i=\sqrt{m/\beta }\widehat{p}_i`$ and $`t=\sqrt{m/\mu }\widehat{t}`$, we can reduce Eq. (1) to $`H=\underset{i=1,\mathrm{}.N}{}(\frac{p_i^2}{2}+\frac{(q_iq_{i+1})^2}{2}+\frac{q_i^4}{4}).`$ So far most of papers considering the statistical behavior of the lattice $`\varphi ^4`$ model have used this form and taken the average initial energy as the control parameter. In this paper, we use an alternative while equivalent method to fix the initial average energy $`\overline{E}=\frac{1}{N}\underset{i=1,\mathrm{}N}{}E_i`$ and vary the coupling strength $`\mu `$ as the control parameter, i.e., we consider the lattice $`\varphi ^4`$ Hamiltonian

$$H=\underset{i=1,\mathrm{}.N}{}(\frac{p_i^2}{2}+\frac{\mu (q_iq_{i+1})^2}{2}+\frac{q_i^4}{4})$$
(3)
which leads to the set of canonical equations

$`\underset{i}{\overset{.}{p}}`$ $`=`$ $`{\displaystyle \frac{H}{q_i}}=q_i^3+\mu (q_{i+1}+q_{i1}2q_i)`$ (4)
$`\underset{i}{\overset{.}{q}}`$ $`=`$ $`{\displaystyle \frac{H}{p_i}}=p_i`$ (5)
We take $`\overline{E}=1.25\times 10^3`$ throughout the paper without lossing any generality.

Before going to the detail of the dynamics of Eqs. (3), we can first briefly discuss the general feature of the system for different coupling regimes. For $`\mu 0,`$ the system consists of practically independent oscillators, each has one degree of freedom, and thus has nearly periodic motion. The dynamics of the system may be quasiperiodic and the approach to equilibrium is not possible. For $`\mu \mathrm{},`$ the nonlinear terms $`q_i^4/4`$ become not important in comparison with the coupling. Thus, different modes are practically uncoupled to each other \[10-12\], and the system dynamics is quasiperiodic again. Complicated dynamics approaching equilibrium may be expected for not too small and not too large couplings (i.e., not too large and not too small energies if we fix $`\mu =1`$ and use the average energy as our control parameter). Up to date, all papers discussing the statistical behaviors of the $`\varphi ^4`$ model consider only small energy (i.e., large coupling $`\mu `$) regime. To our knowledge, the behaviors of large energy (or say, small coupling $`\mu `$) regime has never been discussed. In this paper we will study the microscopic dynamics of the system and the related statistical features in the whole coupling range by paying particular attention to the new area, the small coupling regime.

In the following, we will make a Lyapunov exponent analysis of the system dynamics in Sec. II; and study the statistical quantities, the entropy and the probability distribution of system variables in Sec. III. In Sec. IV, we will consider the scaling behaviors observed in the transitions to equilibrium, and some conclusions and possible future developments will be discussed in the last section.

## II Lyapunov exponent analysis

Given the average energy $`\overline{E}`$ and the coupling $`\mu ,`$ the dynamics of Eqs. (3) is determined by the initial conditions of the $`(q_i,p_i)`$ variables. In this paper we will test and compare two types of typical initial conditions. First, the initial energy is distributed to a single $`k=1`$ mode

$`q_i`$ $`=`$ $`0\text{ },`$ (6)
$`p_i`$ $`=`$ $`A\mathrm{sin}(2\pi i/N)+\sigma _i\text{ },\text{ }|\sigma _i|1\text{ , }A=\sqrt{4\overline{E}}`$ (7)
$`i`$ $`=`$ $`1,2,\mathrm{},N`$ (8)
i.e.,

$`Q_k`$ $`=`$ $`0\text{ , }k=0,1,2,\mathrm{}.N/2`$
$`P_1`$ $`=`$ $`A\sqrt{N/2}\text{}P_k=0\text{ , }k=0,2,3,\mathrm{}.,N/2`$
$`Q_k(t)`$ $`=`$ $`\sqrt{2/N}{\displaystyle \stackrel{N}{\underset{i=1}{}}}q_i(t)\mathrm{cos}(2\pi ki/N),\text{ }P_k(t)=\underset{k}{\overset{.}{Q}}(t)`$
This initial condition has been popularly used in studying the statistical behavior of Hamiltonian coupled oscillators. Second, the initial excitation is restricted on a single site

$`q_i`$ $`=`$ $`0\text{ },\text{ }i=1,2,\mathrm{}.,N`$ (9)
$`p_1`$ $`=`$ $`\sqrt{2N\overline{E}}\text{ and }p_i=\sigma _i\text{ , }i=2,3,\mathrm{}N`$ (10)
which is considered for the first time, to our knowledge, in this paper. Here the equipartition problem will be analysed for both mode distribution and site distribution. For avoiding the strict symmetry restriction preventing equipartitions we add a small random numbers $`\sigma _i`$ in Eqs. (4) and (5) with $`\left|\sigma _i\right|10^7.`$

Lyapunov exponent (LE) analysis is a very useful dynamical tool for revealing the statistical behavior of Hamiltonian systems . If the dynamics is quasiperiodic, the largest LE is zero, and the system will never reach the equal-probability-distribution in the energy surface for $`N2`$. Because for a Hamiltonian system with $`N`$ degrees of freedom the quasiperiodic trajectory can move only in a $`N`$-dimensional torus, while the energy surface has $`2N1`$ dimensions. Therefore, the positivity of the largest LE is a strong indication (of course, not a sufficient condition) of possible equilibrium state.

It is well known that there exist some difficulties for correctly computing LEs in Eqs. (3) when $`N`$ is large . Numerically, LE is calculated approximately with a very long but still finite time $`T`$. However, if $`N`$ is large, it may happen that the system wonders in a certain $`N`$-dimensional torus for an extremely long time and then bursts to a chaotic status. In Figs. 1(a) and (b) we take $`N=64`$ and use the $`k=1`$ mode initial condition (4), and plot the largest LE $`\lambda _1`$ against the average time $`T`$ for $`\mu =6.9\times 10^7`$ and $`2.5`$, respectively. In both case we find that $`\lambda _1`$ is maintained at zero for long time, and jumps to be positive at a certain large $`T`$. If we decrease $`\mu `$ from $`\mu =6.85\times 10^7`$ and increase $`\mu `$ from $`\mu =2.9,`$ the periods for $`\lambda _1`$ to change from zero to positive increase. It seems that there exist two critical couplings $`\mu _1`$ and $`\mu _2,`$ and zero $`\lambda _1`$ is observed for $`\mu <\mu _1`$ or $`\mu >\mu _2,`$ and chaos is identified between these two couplings $`\mu _1<\mu <\mu _2.`$ Transitions from quansiperiodicity to chaos occur at both large $`(\mu _2)`$ and small $`(\mu _1)`$ critical couplings. For approximately fixing $`\mu _1`$ and $`\mu _2,`$ we take average time $`T_0=1\times 10^7,`$ and zero $`\lambda _1`$ is assumed if no burst to positive LE accurs within the period $`T_0.`$

In Fig. 1(c) we plot $`\lambda _1`$ vs $`\mu `$ of Eqs. (3) at $`N=32`$ and for both initial conditions (4) (circles) and (5) (triangles). In each curve we find indeed critical couplings $`\mu _1(\mu _1^{^{}})`$ and $`\mu _2(\mu _2^{^{}})`$, and $`\lambda _1`$ is zero for $`\mu <\mu _1`$ $`(\mu <\mu _1^{^{}})`$ and $`\mu >\mu _2`$ $`(\mu >\mu _2^{^{}}).`$ Positive $`\lambda _1`$, indicating chaos, is observed between these two critical couplings. Another interesting point is that the circle and triangle lines have practical identical positive $`\lambda _1`$ in the region $`\mu _1^{^{}}<\mu <0.25,`$ this verifies approach to equilibrium, which is independent of the initial conditions.

In Fig. 1(d) we do the same as (c) by taking larger system $`N=64`$. For the initial condition (4), the lower (higher) critical coupling $`\mu _1(\mu _2)`$ becomes much smaller (much larger) than that in (c). This feature is intuitively expected since it is familiar to us that for large $`N`$ one needs only small interactions between the subsystems to bring the system to equilibrium states. \[note, for Eqs. (3) large coupling $`\mu `$ corresponds to small interactions between different modes\]. However, for the initial condition (5), the critical coupling $`\mu _1^{^{}}`$ increases as the system size increases. This phenomenon, which seems to be against the intuition, can be easily understood, based on the measue synchronization (MS) of Hamiltonian systems. In order to reach the equilibrium state, all the oscillators of Eqs. (3) should have identical statistical behavior, independing on their initial condition, because they are identical in microscopic dynamics. The necessary condition realizing this identity is that the trajectories of all oscillators should be identical in measure, so-called MS. The critical coupling for the MS between any two neighbor oscillators depends on their energy difference. Large coupling is needed for the MS if a pair of neighbor oscillators have large energy difference. Since for condition (5) the energy difference between the site $`i=1`$ and its neighbors increases by enlarging the system size with fixed average energy \[in (5) $`E_1E_2=N\overline{E}`$\], the critcal coupling for MS, which is nothing but $`\mu _1^{^{}}`$, increases as well for increasing the system size (actually $`\mu _1^{^{}}\sqrt{N}`$). This feature is essentially different from the system behavior with the initial condition (4).

## III Entropies and probability distributions

The dynamic behavior of LE is closely related to the statistical behavior of the system. It is interesting to further investigate the energy distribution of subsystems, because the equipartition behaviors of energies of different modes and different sites are some necessary conditions of the equilibrium state.

The energies for different modes can be computed by

$$E_m(k,t)=\frac{1}{2}P_k^2(t)+\frac{1}{2}\omega _k^2Q_k^2(t)\text{}k=0,2,\mathrm{}.\frac{N}{2}$$
(11)
where $`\omega _k=2\mathrm{sin}(\frac{k\pi }{N})`$. At the same time, the energy for different sites are represented as

$$E_s(i,t)=\frac{1}{2}p_i^2(t)+\frac{1}{4}q_i^4(t)\text{}i=1,2,\mathrm{}..N$$
(12)
For the long time behavior, we are interested in the averages

$`\overline{E}_m(k,T)`$ $`=`$ $`{\displaystyle \frac{1}{T}}{\displaystyle _0^T}E_m(k,t)𝑑t`$ (13)
$`\overline{E}_s(i,T)`$ $`=`$ $`{\displaystyle \frac{1}{T}}{\displaystyle _0^T}E_s(i,t)𝑑t`$ (14)
In actual simulations we take large while finite $`T=10^7`$.

Based on the definition of Eqs. (8), we can investigate the system statistics by defining entropies in the mode space

$`S_m`$ $`=`$ $`{\displaystyle \stackrel{N/2}{\underset{k=0}{}}}\mathrm{\Omega }_k\mathrm{ln}\mathrm{\Omega }_k`$ (15)
$`\mathrm{\Omega }_k`$ $`=`$ $`\overline{E}_m(k)/{\displaystyle \stackrel{N/2}{\underset{k^{^{}}=0}{}}}\overline{E}_m(k^{^{}})`$ (16)
and in the site space

$`S_s`$ $`=`$ $`{\displaystyle \stackrel{N}{\underset{i=1}{}}}\mathrm{\Omega }_i^{^{}}\mathrm{ln}\mathrm{\Omega }_i^{^{}}`$ (17)
$`\mathrm{\Omega }_i^{^{}}`$ $`=`$ $`\overline{E}_s(i)/{\displaystyle \stackrel{N}{\underset{i^{^{}}=1}{}}}\overline{E}_s(i^{^{}})`$ (18)
In numerical computation we measure the scaled entropy differences

$`\eta _m`$ $`=`$ $`(S_m(\mathrm{max})S_m)/S_m(\mathrm{max})`$ (19)
$`\eta _s`$ $`=`$ $`(S_s(\mathrm{max})S_s)/S_s(\mathrm{max})`$ (20)
with $`S_m(\mathrm{max})=\mathrm{ln}(N/2+1)`$, $`S_s(\mathrm{max})=\mathrm{ln}(N).`$ In Fig. 2(a) \[(b)\] we plot $`\eta _m`$ ($`\eta _s`$) vs $`\mu `$ from the initial condition (4) \[the initial condition (5)\]. We find that in (a) $`\eta _m`$ has discontinuous jumps at the two critical couplings $`\mu _1`$ and $`\mu _2`$ and it takes very small value for $`\mu _1<\mu <\mu _2`$, indicating large $`S_m`$ entropy in this coupling regime. In (b) $`\eta _s`$ has a plateau of $`\eta _s1`$ for $`\mu <\mu _1^{^{}},`$ and it jumps down to about $`\eta _s<1\times 10^6,`$ and keeps this small value in the entire regime $`\mu >\mu _1^{^{}}.`$ It is noticed that the second transition at $`\mu _2^{^{}}`$ is not manifested by $`\eta _s`$, because the MS between sites is not affected by this transition.

In order to show that in the regime $`\mu _1^{^{}}<\mu <1.65`$, the system actually reach the statistical equilibrium state for both initial conditions, we take $`\mu =0.1`$ and plot in Figs. 2(c) and (d) the probability distributions of $`x`$ and $`p`$ variables, respectively, for all the $`64`$ particles and for both initial conditions (4) and (5). Each probability distribution curve is drawn with $`1\times 10^6`$ data. From the figures it is verified that all the sites have identical distributions, showing perfect equipartion behaviors.

## IV scaling behaviors in the transitions to equilibrium

In the above numerical analysis we found two transitions of Eqs. (3) for the initial condition (4) from quansiperiodicity to chaos, one at small coupling $`\mu _1`$, transition to equipartion between different sites in real space, the other at large coupling $`\mu _2`$, transition to equipartition between various modes in Fourier space. Both $`\mu _1`$ and $`\mu _2`$ can be considerably changed by varying the system size. It is interesting to investigate the possible scaling behavior between the critical $`\mu `$ and the system size $`N.`$

In Fig. 3(a) we plot both $`\mathrm{ln}\mu _1`$ and $`\mathrm{ln}\mu _2`$ vs $`\mathrm{ln}N.`$ Between the two curves the system motion is chaotic, where statistical equipartition behavior is expected. And quasiperiodicity exists in the regions below the $`\mu _1`$ curve and above the $`\mu _2`$ curve. By quasiperiodicity we mean that quasiperiodic motion is observed in the time period $`tT=1\times 10^7.`$ And we assume, (without proof) that the actual critical couplings for the transitions are very near the plotted curves of Fig. 3(a).

Both $`\mu _1(N)`$ and $`\mu _2(N)`$ curves have well behaved power law scalings

$$\mu _1(N)\frac{1}{N^\alpha },\text{ }\alpha 3.82$$
(21)
$$\mu _2(N)N^\beta ,\text{ }\beta 1.92$$
(22)
The solid lines $`\mu _1(N)`$ and $`\mu _2(N)`$ are given by $`A/N^\alpha ,`$ $`A1\times 10^3,`$ $`\alpha 3.82`$ and $`\mu _2(N)=BN^\beta ,`$ $`B2.4,`$ $`\beta 1.92,`$ respectively. Equations (12) and (13) are numerical results. We do not know how to derive the exponents $`\alpha `$ and $`\beta `$ analytically, and even have no idea at present to explain them. We suggest these scaling may be model dependent and do not have universal behavior. For instance, the scaling (13) is different from the scaling Casetti and coworkers obtained for the FPU model.

With the initial condition (5) the critical point $`\mu _1^{^{}}`$ and $`\mu _2^{^{}}`$ are also related to $`N`$. The scalings are shown in Fig. 3(b), which read

$`\mu _1^{^{}}(N)`$ $``$ $`N^\gamma \text{ },\text{ }\gamma 0.57`$ (23)
$`\mu _2^{^{}}(N)`$ $``$ $`N^\beta \text{ },\text{ }\beta 1.95`$ (24)
## V conclusion and discussion

In conclusion, we have investigated the statistical behavior of the classical $`\varphi ^4`$ lattice from its microscopic dynamics. The particular attention focuses on the transitions to chaos, equipartition behavior and equilibrium state. By varying the coupling strength, two transitions to equipartition from both small and large couplings are observed. The scaling behaviors of these couplings against the system size are revealed numerically.

However, there are still some problems unsolved. First, it is difficult to precisely fix the transition points, because near the transitions the system may show zero LE for extremely long period, and suddenly jump to positive LE after a certain time. Thus, by a finite time simulation, one cannot guarantee zero LE, and consequently cannot definitely determine the critical conditions. In this paper, we apply an time truncation in identifying $`\mu _1`$ and $`\mu _2`$, i.e., we regard the largest LE to be zero if its zero value is kept for a given very long while still finite time \[see Ref.\]. We guess that the actual transition points are not far from the critical couplings shown in Figs. 1-3. We will try a systematical method to convincingly identify the transitions in the future work. Moreover, the mechanism underlying the scaling relations (12), (13) and (15) are still unclear, and worthwhile further investigating.

ACKNOWLEDGEMENT

This research was supported by the National Natural Science Foundation of China, the Nonlinear Science Project of China, and the Foundation of Doctoral training of Educational Bureau of China.

Captions of figures

Fig. 1 The largest LE of Eqs. (3) $`\lambda _1`$ for the fixed average energy $`\overline{E}=\frac{1}{N}\stackrel{N}{\underset{i=1}{}}E_i=1.25\times 10^3.`$ (a) The largest LE vs the time $`T`$ with $`T`$ being the time for computing $`\lambda _1.`$ $`\mu =6.9\times 10^7,`$ $`N=64,`$ and the initial condition $`(4)`$ is taken. It is shown that $`\lambda _1(T)`$ remains zero for a long period, and it bursts to be positive after a certain time. (b) the same as (a) with $`\mu =2.5`$. (c) $`\lambda _1`$ vs $`\mu `$ for $`N=32.`$ Circles are for the initial condition (4) and triangles are for the initial condition (5). The transition points $`\mu _1`$, $`\mu _2`$ and $`\mu _1^{^{}}`$, $`\mu _2^{^{}}`$ are determined by computations of $`T=1\times 10^7.`$ (d) The same as (c) with $`N=64.`$

Fig. 2 $`N=64.`$ $`\eta _m`$ \[(a)\] and $`\eta _s`$ \[(b)\] defined in Eq. (11) plotted vs $`\mu .`$ In (a) the initial condition (4) is used and in (b) Eq. (5) is used. $`\eta _m`$ has discontinuities at two transition points $`\mu _1`$ and $`\mu _2`$ of Fig. 1(d), and $`\eta _s`$ jumps down at the transition point $`\mu _1^{^{}}`$ of Fig. 1(d). (c), (d) $`\mu =0.1.`$ The $`x_i`$ \[(c)\] and $`p_i`$ \[(d)\] probability distributions of the sites $`i=1,2,\mathrm{}.N`$ for both initial conditions (4) and (5).

Fig. 3 (a) The transition couplings $`\mu _1`$ and $`\mu _2`$ plotted for different systems size with the initial condition (4) being taken. In the region between $`\mu _1`$ and $`\mu _2`$ lines, the systems motion is chaotic, and equipartition behavior can be identified. $`\mu _1`$ and $`\mu _2`$ have power law relations with the size $`N`$ as $`\mu _1`$ $`1/N^{3.82},`$ $`\mu _2N^{1.92}.`$ The exponents are measured by the slopes of the two solid lines. (b) The same as (a) with the initial condition (5) used. $`\mu _1^{^{}}`$ $`N^{0.57},`$ $`\mu _2^{^{}}N^{1.95}.`$


text:
# Dynamics of a Bose-Einstein condensate at finite temperature in an atomoptical coherence filter

## Abstract

The macroscopic coherent tunneling through the barriers of a periodic potential is used as an atomoptical filter to separate the condensate and the thermal components of a 87Rb mixed cloud. We condense in the combined potential of a laser standing-wave superimposed on the axis of a cigar-shape magnetic trap and induce condensate dipole oscillation in the presence of a static thermal component. The oscillation is damped due to interaction with the thermal fraction and we investigate the role played by the periodic potential in the damping process.

PACS numbers: 03.75.Fi, 32.80.Pj, 03.75.-b

The combination of Bose-Einstein condensates (BECs) with periodic optical potentials provides a versatile tool for getting a deeper understanding of macroscopic quantum phenomena, directly addressing the issue of phase coherence between degenerate coupled Bose gases  recently culminating in the observation of the Mott insulator phase transition . Properties of condensates confined in an optical lattice have opened new perspectives in the investigation of different dynamical phenomena such as superfluid flow , atomic Bloch oscillations , the Josephson effect  and have stimulated speculations on a wide range of phenomena including solitons , dynamical instability  and transport behavior .

The dynamics of BECs at finite temperatures, is a challenging question even in the pure harmonic case. Many efforts have been devoted to the extension of theoretical models at nonzero temperature and to the study of the low-energy collective excitations in presence of a significant thermal component .

From the experimental point of view the effects of the interactions between the condensate and the thermal fraction on the modes of excitation have been investigated only in harmonic potentials. Temperature dependence of damping and frequency shifts in quadrupole and scissors modes have been investigated in a TOP trap  while also a damped out-of-phase dipolar oscillation of the thermal cloud and the condensate has been observed in a cloverleaf trap .

In this Letter we address the problem of the dynamics of a condensate in a periodic potential at finite temperature. This is done by using an optical lattice as an atomoptical filter to induce, through a selective manipulation, a relative motion between the condensate and thermal components of a magnetically trapped mixed cloud of rubidium atoms. In the presence of the periodic potential these two components respond differently to a sudden displacement of the magnetic trap. Atoms with sufficiently long coherence lengths are allowed to tunnel through the potential barriers of a one dimensional optical lattice while atoms with a higher kinetic energy remain confined in the potential wells. While the condensate coherently flows through the optical barriers, thus performing a dipole oscillation, the thermal component is blocked in the presence of the periodic potential and provides a damping mechanism. In particular, we observe a sharp change of the damping rate in the BEC dipole oscillations by increasing the optical potential depth.

Our procedure for creating Bose condensates in an optical lattice, already described in details elsewhere , can be summarized as follows. Starting with a sample of ultracold atoms of 87Rb confined in a Ioffe-type magnetic trap in the $`(F=1,m_F=1)`$ state, we perform rf-evaporative cooling until we reach a temperature slightly above the condensation threshold. Next, we superimpose a one-dimensional optical lattice to the longitudinal $`x`$-axis of the magnetic trap by means of a far detuned, retroreflected laser beam. We then continue the evaporation process through the phase-transition temperature . In this way, the atomic cloud is prepared in an equilibrium state of the combined magnetic and optical dipole potentials. In this experiment we typically stop the evaporation when the fraction of condensed atoms is $`20\%`$ corresponding to a temperature $`T120`$ nK. In the atomic cloud region the resulting potential has the form

$`V`$ $`=`$ $`V_{mag}+V_{opt}`$ (1)
$`=`$ $`{\displaystyle \frac{1}{2}}m[\omega _x^2x^2+\omega _{}^2(x^2+y^2)]+sE_R\mathrm{cos}^2\left({\displaystyle \frac{2\pi x}{\lambda }}\right)`$ (2)
with $`m`$ the atomic mass, $`\omega _x=2\pi \times 9`$ Hz and $`\omega _{}=2\pi \times 92`$ Hz the axial and radial frequency of the magnetic potential respectively. In (2) the optical potential depth is expressed, via the dimensionless factor $`s`$, in units of the recoil energy of an atom absorbing one laser photon, $`E_R=h^2/2m\lambda ^2`$, where $`\lambda `$ is the wavelength of the laser creating the standing wave.

By varying the intensity of the laser beam, (detuned typically $`\mathrm{\Delta }\lambda 3`$ nm to the blue of the $`D_1`$ line at $`\lambda 795`$ nm), we change the value of $`s`$ up to 2.5. The optical potential is then calibrated by measuring the Rabi frequency of the Bragg transition between the momentum states $`h/\lambda `$ and $`+h/\lambda `$ induced by the standing wave . Due to the large detuning of the optical lattice, the maximum spontaneous scattering rate is $`\mathrm{\Gamma }_{sp}=0.25`$ Hz and spontaneous scattering can be neglected on the timescale of our experiment.

In order to induce the relative center-of-mass motion of the condensate and the thermal fractions in the combined trap, we first displace instantaneously ($`t_{dis}2\pi /\omega _x`$) the magnetic trapping potential in the $`x`$-direction by a distance $`\mathrm{\Delta }x=18\mu `$m sufficiently small to remain in the superfluid regime . After the displacement, the atomic cloud finds itself out of equilibrium, nevertheless due to the presence of the optical lattice, only the condensate is forced into motion by the magnetic potential gradient. Indeed, when the potential energy given by the displacement is smaller than the optical potential depth, as in our experiment, a motion along the $`x`$-direction can be established only by means of tunneling through the barriers of the standing wave. We have previously shown  that a pure condensate and a thermal cloud above Tc behave in a dramatically different way. The coherence properties of the condensate imply a macroscopic tunneling and a collective dipole motion can occur. On the contrary, the incoherent nature of the thermal cloud drastically reduces the collective tunneling probability, thus preventing the establishing of a collective oscillation. In the present work we produce mixed clouds and the optical lattice acts as a filter to separate the coherent and incoherent components. After the displacement, we wait half a period of the BEC oscillation in the harmonic trap ($`\pi /`$ $`\omega _x`$) then we move back the magnetic potential by the distance $`\mathrm{\Delta }x`$. At the end of this procedure, the thermal cloud finds itself at rest in the center of the trap, while the condensate oscillates inside the thermal cloud with an initial amplitude $`2\mathrm{\Delta }x`$. After a variable evolution time, the trapping potential is switched off and the cloud is imaged subsequently to an additional ballistic expansion of 28 ms. From the observed column density we deduce the condensate fraction, the temperature and the center-of-mass position of both the condensate and the thermal components in the trap. In Fig. 1 we report three absorption images corresponding to successive times of the evolution in the combined trap.

In order to show the effects of the optical lattice on the dynamics of a mixed cloud, we study the center-of-mass motion of the condensate and the thermal components in two different configurations. In one (harmonic case) we switch off the periodic potential after the preparatory procedure described above while in the other (combined case) we leave the lattice on also during the subsequent evolution.

Harmonic case \- In absence of the lattice both the condensate and the thermal cloud are free to oscillate. Starting from an off-equilibrium configuration in which the thermal cloud is at the center of the trap, the condensate mean-field pushes the thermal fraction out of the potential minimum. During the time evolution the oscillation amplitude of the thermal cloud rises until an equilibrium is reached (Fig. 2a). Notice that, in a pure harmonic trap, when no relative motion is induced, the BEC dipolar mode can not be affected by the presence of the thermal component . We expect that the kinetic energy lost by the condensate is essentially converted into kinetic energy of the thermal cloud. From the fit of the relative center-of-mass position (Fig. 2b) with a damped oscillation, we obtain a frequency of $`8.42\pm 0.04`$ Hz which is smaller of roughly $`5\%`$ from the measured trapping frequency. This is a consequence of the interaction between the oscillating condensate and the initially stationary thermal cloud and is similar to the observation reported in  for an oscillating thermal cloud in presence of an initially stationary condensate.

Combined case \- In presence of the optical lattice the thermal cloud stays locked at the center of the trap due to its incoherent nature, while the condensate oscillates as shown in Fig. 3. Again the coupling between the condensate ad the thermal cloud causes damping in the center-of-mass oscillations of the condensate. Here, the kinetic energy lost by the BEC during the evolution cannot be converted into the kinetic energy of the thermal cloud which is kept fixed by the periodic potential and must be therefore partly converted into internal energy of the cloud, partly absorbed by the optical lattice. Note that in this case the pure condensate would perform undamped oscillations .

We investigate the damping rate of the condensate center-of-mass motion as a function of the optical potential depth $`sE_R`$ for a fixed trap displacement. The experiment is performed in the collisionless regime where collisions between thermal atoms are negligible. In this regime Landau damping, present when a thermal bath of elementary excitations absorbs quanta of the condensate collective excitations, represents an important mechanism to explain the dynamical behaviour of trapped Bose gases .

For temperatures larger than the chemical potential, the Landau damping rate of low-energy excitations in a harmonically trapped Bose-condensed gas can be estimated using the form

$$\mathrm{\Gamma }_L=\frac{3\pi }{8}\frac{k_BTa\omega }{\mathrm{}c}$$
(3)
where $`c=\mathrm{}\sqrt{4\pi an_0}/m`$ is the sound velocity in the center of the trap, $`n_0`$ being the central density of the condensate, $`a`$ the s-wave scattering length and $`\omega `$ the frequency of the collective excitation. In our case, for small values of $`s`$, we can estimate $`\mathrm{\Gamma }_L`$ calculating the central density of the condensate in presence of the periodic potential using a variational approach with a gaussian ansatz for the condensate wave-function. For $`s=0.7`$ and assuming a fixed temperature $`T=120`$ nK, we have $`\mathrm{\Gamma }_L=1.6`$ Hz in fair agreement with the experimental results.

In Fig. 4 we show the damping rate measured as a function of the optical potential depth. Close to $`s=1`$ we observe a sudden increase of the damping rate. This occurs in coincidence with the formation of the first bound state of the condensate in the lattice . Actually the occurrence of a band structure strongly modifies the energy spectrum. In particular the Landau damping, being a resonant energy transfer, strongly depends on the form of energy spectrum  and an adequate treatment would be needed. In addition a theoretical analysis of finite temperature damping of excitations in presence of a periodic potential should also include other damping mechanisms such as inter-component damping recently suggested in . However, at present, a comprehensive theoretical analysis which accounts for all the effects introduced by the presence of the lattice is lacking.

In conclusion, we have demonstrated the operation of an atomoptical filter capable of separating the condensate and the thermal fraction in a mixed atomic cloud across Tc. This filtering technique could, e. g., be applied to spatially separate the ground state of a finite temperature BEC fully from the thermal component. In this work, we are able to control the relative motion between the center-of-masses of the condensate and of the thermal fraction in a mixed cloud, this enable us to investigate quantitative effects both in the amplitude and in the frequency of the dipole oscillations. Damping is caused by interactions between condensate and the thermal fractions. Novel features are expected in presence of the periodic lattice potential which modifies the energy spectrum of the system. Indeed by changing the optical lattice parameters a sharp rise in the damping rate has been observed when a bound state is likely to be formed.

We are indebted to Francesco Minardi for many useful discussions and for his help in the initial stage of the experiment. We acknowledge Gabriele Ferrari, Sandro Stringari, Stefano Giorgini and Lev Pitaevskii for fruitful discussions, and Maurizio Artoni for careful reading of the manuscript.

This work has been supported by the EU under Contracts No. HPRI-CT 1999-00111 and No. HPRN-CT-2000-00125, by the MURST through the PRIN 1999 and 2001 Initiatives and by the INFM Progetto di Ricerca Avanzata “Photon Matter”.