Two-fluid hydrodynamics of cold atomic bosons under influence of the quantum fluctuations at non-zero temperatures
aa r X i v : . [ c ond - m a t . qu a n t - g a s ] J a n Two-fluid hydrodynamics of cold atomic bosons under influence of the quantumfluctuations at non-zero temperatures
Pavel A. Andreev ∗ Faculty of physics, Lomonosov Moscow State University, Moscow, Russian Federation, 119991. andPeoples Friendship University of Russia (RUDN University),6 Miklukho-Maklaya Street, Moscow, 117198, Russian Federation (Dated: January 29, 2021)Ultracold Bose atoms is the physical system, where the quantum and nonlinear phenomena playcrucial role. Ultracold bosons are considered at the small finite temperatures. Bosons are consideredas two different fluids: Bose-Einstein condensate and normal fluid (the thermal component). Anextended hydrodynamic model is obtained for both fluids, where the pressure evolution equationsand the pressure flux third rank tensor evolution equations are considered along with the continuityand Euler equations. It is found that the pressure evolution equation contains zero contributionof the short-range interaction. The pressure flux evolution equation contains the interaction whichgives the quantum fluctuations in the zero temperature limit. Here, we obtain its generalization forthe finite temperature. The contribution of interaction in the pressure flux evolution equation whichgoes to zero in the zero temperature limit is found. The model is obtained via the straightforwardderivation from the microscopic many-particle Schrodinger equation in the coordinate representa-tion.
I. INTRODUCTION
Small temperature bosons are studied in terms of two-fluid hydrodynamics consisting of the Bose-Einstein con-densate (BEC) and normal fluid [1]. Each fluid is consid-ered in terms of two hydrodynamic equations: the conti-nuity and Euler equations. It is assumed that the BECcan be completely described by the concentration andvelocity field, or, in other terms, by the Gross-Pitaevskiiequation [1], since BEC is the collection of particles in thesingle quantum state. However, the normal fluid modelrequires an truncation of the set of hydrodynamic equa-tions. The pressure of normal fluid existing in the Eulerequation for the normal fluid is an independent function.Equation for the pressure evolution provides an expres-sion for the pressure perturbations via the perturbationsof other functions. Application of the equation of statefor pressure makes the model more simple, but equationof state for pressure leads to the less accurate model.Moreover, the kinetic pressure in the Euler equationfor the BEC is usually chosen to be equal to zero [2],[3]. Since the kinetic pressure is related to the occupa-tion of the excited states However, there is the nonzeropart caused by the quantum fluctuations [4]. The pres-sure evolution equation of the weakly interacting bosonscontains no interaction, but next equation in the chainof the quantum hydrodynamic equations (the equationfor the pressure flux third rank tensor) contains the in-teraction causing the depletion of the BECs at the zerotemperature. Therefore, it is necessary to consider thepressure flux evolution equation both for the BEC and forthe normal fluid at the analysis of the small temperatureinfluence. ∗ Electronic address: [email protected]
The quantum depletion of the BECs is the appear-ance of the bosons in the excited states while system iskept at the zero temperature. So, some energy of thecollective motion is transferred to the individual motionof a portion of particles. It is caused by the quantumfluctuations related to the interparticle interaction. Thequantum fluctuations are considered in literature for along time. Mostly, their theoretical analysis is basedon the Bogoliubov-de Gennes approach [5], [6], [7], [8].The quantum fluctuations in BECs are studied experi-mentally as well [9], [10], [11]. This method is general-ized for the dipolar BECs, where the quantum fluctua-tions plays crucial role at the description of the dipolarBECs of lantanoids. The dipolar lantanoid BECs re-veal the large scale instability causing the splitting ofthe cloud of atoms on the number of macroscopic drops.This highly nonlinear phenomena is called the quan-tum droplet formation [12–31]. Therefore, reassemble ofbosons in smaller compact groups causes the stabiliza-tion of the system. These studies give the motivation forthe study of quantum fluctuations. However, here we re-strict our analysis with the short-range interaction onlyand no dipole-dipole interaction is discussed. Moreover,the influence of the finite temperature is the necessarypart of complete model of these phenomena.Fundamental feature of the collective dynamics in thespectrum of the sound waves. Two distinct sound veloci-ties exist in finite temperature ultracold Bose gas [1], [32].The two-fluid model shows that the slower mode (secondsound) is associated with the BEC component, while thefaster mode (first sound) is associated with the thermalcomponent. Generalized expressions for the speeds ofsounds are obtained within developed model.Derivation of two-fluid hydrodynamics for the finitetemperature bosons in the limit of small temperature,where the large fraction of the bosons is located in theBEC state is given from the microscopic motion in accor-dance with the quantum hydrodynamic method [2], [3],[33], [34]. The microscopic dynamics is described by theSchrodinger equation in the coordinate representation.Collection of the macroscopic functions is presented todescribe the collective effects in ultracold bosons. Thelast includes the concentration of particles, the velocityfield and the pressure tensor. The derivation of basicequations is made for all bosons distributed on the lowerenergy level and the excited levels as the single fluid. Thedecomposition on two fluids is made on the microscopicscale. After general structure of equations is obtainedfor the arbitrary temperature and arbitrary strength ofinteraction, an approximate calculation of functions pre-senting the interaction is made for the regime of short-range interaction. Hence, the small parameter related tothe small area of interaction potential is used. It givesa specification for general model, but also the first ordercontribution on the small parameter is applied. The fur-ther truncation is made at calculation of the interactionterms for weak interaction and small temperatures.This paper is organized as follows. In Sec. II majorsteps of derivation of hydrodynamic equations from theSchrodinger equation are demonstrated, where the pres-sure evolution equation (the quantum Bohm potentialevolution equation) and the third rank tensor evolutionequation are obtained along with the continuity and Eu-ler equations. In Sec. III calculation of interaction inthe Euler equation, the pressure evolution equation, andthe pressure flux evolution equation is demonstrated. InSec. IV presents the suggested version of the extendedtwo-fluid quantum hydrodynamic model for the ultracoldfinite temperature bosons. In Sec. V the limiting regimeof derived model for the BEC is obtained under influenceof the quantum fluctuations at the zero temperature. InSec. VI a brief summary of obtained results is presented.
II. ON DERIVATION OF HYDRODYNAMICEQUATIONS FROM MICROSCOPIC QUANTUMDYNAMICSA. Basic definitions of quantum hydrodynamicsand the Euler equation derivation
On the microscopic level we do not have notion of tem-perature. Hence we consider system of interacting bosonsgoverned by the Schrodinger equation ı ¯ h∂ t Ψ = ˆ H Ψ withthe following Hamiltonianˆ H = N X i =1 (cid:18) ˆ p i m i + V ext ( r i , t ) (cid:19) + 12 X i,j = i U ( r i − r j ) , (1)where m i is the mass of i-th particle, ˆ p i = − ı ¯ h ∇ i is themomentum of i-th particle. The last term in the Hamil-tonian (1) is the boson-boson interaction U ij . We do notspecify the form of interaction. However, the derivationpresented below employs that the interaction has finitevalue on the small distances between particles and shows fast decay at the increase of the interparticle distance.Definitely, no distinguishing between bosons in the BECstate and bosons in other states is made at this stage.Separation of all bosons on two subsystems is made interms of collective variables.Hydrodynamic model usually made for each species ofparticles. If we consider a single species then all massesequal to each other.Distribution of particles in a trap, waves, solitons, os-cillations of form of trapped particles are described by theconcentration of particles. Concentration is an essentialmacroscopic function both for classical and for quantumfluids. The module of the macroscopic wave function inthe Gross-Pitaevskii equation gives the square root ofconcentration of bosons in the BEC state. Therefore, westart the derivation of quantum hydrodynamic equationsfrom the definition of concentration. The quantum me-chanics is based on notion of point-like objects in spitethe wave nature of quantum objects. So, the eigenfunc-tion of the coordinate operator in the coordinate repre-sentation is the delta function ˆ xψ x ′ ( x ) = x ′ ψ x ′ ( x ), wherenormalized wave function is ψ x ′ ( x ) = δ ( x − x ′ ). Obvi-ously, the operation of concentration in the coordinaterepresentation of quantum mechanics is the sum of deltafunctions ˆ n = P Ni =1 δ ( r − r i ). Moreover, it is supportedby the general rule for quantization. We need to takecorresponding classical function.Transition to description of the collective motion ofbosons is made via introduction of the concentration [2],[33]: n = Z dR N X i =1 δ ( r − r i )Ψ ∗ ( R, t )Ψ(
R, t ) , (2)which is the first collective variable in our model. Othercollective variables appear during the derivation. Equa-tion (2) contains the following notations dR = Q Ni =1 d r i isthe element of volume in 3 N dimensional configurationalspace, with N is the number of bosons. Concentration(2) is the sum of partial concentrations n = n n + n b de-scribing the distribution of BEC n b and normal fluid n n in the coordinate space.The equation for evolution of concentration (2) can beobtained by acting by time derivative on function (2).The time derivative acts on the wave functions under theintegral while the time derivatives of the wave functionare taken from the Schrodinger equation. Obtain thecontinuity equation for concentration (2) after straight-forward calculations ∂ t n + ∇ · j = 0 , (3)where the new collective function called the current ap-pears as the following integral of the wave function j ( r , t ) = Z dR N X i =1 δ ( r − r i ) × × m i (Ψ ∗ ( R, t )ˆ p i Ψ( R, t ) + c.c. ) , (4)with c.c. is the complex conjugation.Both introduced collective functions n ( r , t ) and j ( r , t )are quadratic forms of the wave function. Each of themcan be splitted on two parts related to the BEC andnormal fluid. Hence we have n = n n + n b and j = j n + j b .No microscopic definitions are introduced for the partialfunctions n n , n b , j n , and j b . Therefore, the continuityequation (3) splits on two partial continuity equations ∂ t n a + ∇ · j a = 0 , (5)where subindex a stands for b and n .Continue the derivation of hydrodynamic equationsand consider the time evolution of the particle current(4). Act by time derivative on function j (4) and use theSchrodinger equation with Hamiltonian (1). It leads tothe general form of the current evolution equation ∂ t j α + ∂ β Π αβ = − m n∂ α V ext + 1 m F αint , (6)whereΠ αβ = Z dR N X i =1 δ ( r − r i ) 14 m [Ψ ∗ ( R, t )ˆ p αi ˆ p βi Ψ( R, t )+ˆ p α ∗ i Ψ ∗ ( R, t )ˆ p βi Ψ( R, t ) + c.c. ] (7)is the momentum flux, and F αint = − Z ( ∂ α U ( r − r ′ )) n ( r , r ′ , t ) d r ′ , (8)with the two-particle concentration n ( r , r ′ , t )= Z dR N X i,j =1 ,j = i δ ( r − r i ) δ ( r ′ − r j )Ψ ∗ ( R, t )Ψ(
R, t ) . (9)It is necessary to split equation (6) on two equations foreach subsystem of bosons. In current form equation (6)consist of superposition of functions which are quadraticforms of the wave function. Hence, each term can besplitted on two parts and we find two similar equationsfor the currents ∂ t j αa + ∂ β Π αβa = − m n a ∂ α V ext + 1 m F αa,int . (10)The first and third terms are proportional to the concen-tration and the current. therefore, they require no com-ments. Nontrivial difference between two current evolu-tion equation appears at further analysis of the momen-tum flux Π αβ and the interaction F αint . However, we pointout some difference which appear for the momentum flux Π αβ . Its structure is obtained in many papers (see forinstance [2] after equation (52), [35] equation (24))Π αβ = nv α v β + p αβ + T αβ , (11)where p αβ is the pressure tensor, and T αβ is the ten-sor giving the quantum Bohm potential, its approximateform can be written in the following form T αβ = − ¯ h m (cid:20) ∂ α ∂ β n − ∂ α n · ∂ β nn (cid:21) . (12)Tensor T αβ (12) is obtained for noninteracting particleslocated in the single quantum state.Basically, the pressure tensor p αβ is defined via thewave function Ψ( R, t ). However, it requires some ma-nipulations with the wave function and introduction ofa number of intermediate function. Hence we do notpresent its explicit form. Nevertheless, the pressure ten-sor is related to the distribution of bosons on quantumstates with energies above E min . Therefore, for bosonsin the BEC state we have p αβB = 0 if no quantum fluctu-ations are considered andΠ αβB = n B v αB v βB + T αβB + p αβqf , (13)where T αβB is the function of n B (12) if there is no in-teraction. Distribution of particles on different quantumstates does not allow to get full expression (12), but thefirst term. However, it can be used as an equation of statefor noninteracting limit. The normal fluid bosons havenonzero pressure p αβn = 0. Hence, all terms in presenta-tion (11) exists in this regime. Expression (12) appearsfor bosons in the single state in the absence of interaction.Hence, it is an approximate expression for the weakly in-teracting bosons being in the BEC state. It is even lessaccurate for normal fluid bosons, but we use it as anequation of state for the quantum part of the momentumflux. B. The pressure evolution equation
Extending the set of hydrodynamic equations we canderive the equation for the momentum flux evolution. Itcan be expected that this equation brings extra infor-mation for the normal fluid bosons only. However, thequantum fluctuations give contribution in the evolutionof the kinetic pressure of BECs in the limit of zero tem-perature via the divergence of the third rank tensor. Ifthe temperature is nonzero we have two partial kineticpressures for the BEC and for the normal fluid. Considerthe time evolution of the momentum flux (7) using theSchrodinger equation with Hamiltonian (1).It gives to the following expression: ∂ t Π αβ = ı ¯ h Z dR N X i =1 δ ( r − r i ) 14 m [ ˆ H ∗ Ψ ∗ ( R, t )ˆ p αi ˆ p βi Ψ( R, t ) − Ψ ∗ ( R, t )ˆ p αi ˆ p βi ˆ H Ψ( R, t ) + ˆ p α ∗ i ˆ H ∗ Ψ ∗ ( R, t )ˆ p βi Ψ( R, t ) − ˆ p α ∗ i Ψ ∗ ( R, t )ˆ p βi ˆ H Ψ( R, t ) − c.c. ] (14)The part of the presented terms contains the Hamilto-nian ˆ H under action of the momentum operators. Wepermute the Hamiltonian ˆ H and the operators acting onit. Hence, the result of permutation presented by termswhere no operators act on the Hamiltonian ˆ H . However,the terms containing the corresponding commutators ap-pear. Therefore, all terms are combined in two groups: ∂ t Π αβ = ı ¯ h Z dR N X i =1 δ ( r − r i ) 14 m [ ˆ H ∗ Ψ ∗ · ˆ p αi ˆ p βi Ψ − Ψ ∗ ˆ H ˆ p αi ˆ p βi Ψ + ˆ H ∗ ˆ p α ∗ i Ψ ∗ · ˆ p βi Ψ − ˆ p α ∗ i Ψ ∗ · ˆ H ˆ p βi Ψ − c.c. ]+ ı ¯ h Z dR N X i =1 δ ( r − r i ) 14 m [ − Ψ ∗ [ˆ p αi ˆ p βi , ˆ H ]Ψ+[ˆ p α ∗ i , ˆ H ∗ ]Ψ ∗ · ˆ p βi Ψ − ˆ p α ∗ i Ψ ∗ · [ˆ p βi , ˆ H ]Ψ − c.c. ] (15)The first group of terms in expression (15) gives the di-vergence of flux of tensor Π αβ . The second group of termscontains the commutators. This group leads to the con-tribution of interaction in the momentum flux evolution.It gives the momentum flux evolution equation ∂ t Π αβ + ∂ γ M αβγ = − m j β ∂ α V ext − m j α ∂ β V ext + 1 m ( F αβ + F βα ) , (16)where the momentum flux is the full flux of all bosonsΠ αβ = Π αβn + Π αβb , the splitting on two subspecies is tobe made later, F αβ = − Z [ ∂ α U ( r − r ′ )] j β ( r , r ′ , t ) d r ′ (17)is the force tensor field, M αβγ = Z dR N X i =1 δ ( r − r i ) 18 m i (cid:20) Ψ ∗ ( R, t )ˆ p αi ˆ p βi ˆ p γi Ψ( R, t )+ˆ p α ∗ i Ψ ∗ ( R, t )ˆ p βi ˆ p γi Ψ( R, t ) + ˆ p α ∗ i ˆ p γ ∗ i Ψ ∗ ( R, t )ˆ p βi Ψ( R, t )+ˆ p γ ∗ i Ψ ∗ ( R, t )ˆ p αi ˆ p βi Ψ( R, t ) + c.c. (cid:21) (18)is the current (flux) of the momentum flux, and j ( r , r ′ , t ) = Z dR X i,j = i δ ( r − r i ) δ ( r ′ − r j ) × × m i [Ψ ∗ ( R, t )ˆ p i Ψ( R, t ) + c.c. ] (19)is the two-particle current-concentration function.If quantum correlations are dropped function j α ( r , r ′ , t ) splits on product of the current j α ( r , t ) andthe concentration n ( r ′ , t ). Interaction in the momentumflux evolution equation (16) is presented by symmetrizedcombinations of tensors F αβ , which is the flux or currentof force.Partial momentum flux equations appear as ∂ t Π αβa + ∂ γ M αβγa = − m j βa ∂ α V ext − m j αa ∂ β V ext + 1 m ( F αβa + F βαa ) , (20)where M αβγ = M αβγB + M αβγn , with M αβγa = n a v αa v βa v γa + v αa ( p βγa + T βγa ) + v βa ( p αγa + T αγa )+ v γa ( p αβa + T αβa ) + Q αβγa + T αβγa + L αβγa . (21)The pressure is the average of the square of the thermalvelocity, when tensor Q αβγa is the average of the prod-uct of three projections of the thermal velocity. Function L αβγa presents quantum-thermal terms. For the BEC wehave p αβB = 0, Q αβγB = 0, L αβγB = 0, since it has no contri-bution of the excited states. For symmetric equilibriumdistributions we have Q αβγn = 0, L αβγn = 0. We general-ize this conclusion for nonequilibrium states as the trivialequations of state for these functions. Tensor T αβγa is T αβγa = − ¯ h m n a ( ∂ α ∂ β v γa + ∂ α ∂ γ v βa + ∂ β ∂ γ v αa ) . (22)This definition of tensor T αβγ differs from equation (27)in Ref. [35] by extraction of the quantum Bohm poten-tials written together with pressure tensors in equation(21). Equation (27) in Ref. [35] contains approximateform of the quantum Bohm potential T αβ . Equation (21)includes the quantum Bohm potential in its general form.Moreover, expression (22) is an exact formula obtainedwith no assumption about structure of the many-particlewave function like the first term in equation (23) in Ref.[35].Equations (3)-(20) are obtained in general form. Theshort-range nature of the inter-particle interaction is notused. Moreover, the traditional hydrodynamic equationsare presented via the velocity field and the pressure ten-sor while equations (3)-(20) are written via the currentand the momentum flux.The method of the introduction of the velocity fieldin the equations of quantum hydrodynamics of spinlessparticles is presented in Refs. [2], [35]. The method ofcalculation of the terms containing interaction for theshort-range interaction limit is also described in Refs.[2], [35]. Let us present results of application of thesemethods for finite temperature bosons. Moreover, weconsider the short-range interaction in the first order bythe interaction radius. C. Appearance of the quantum fluctuations in thethird rank tensor evolution equation
Derivation of the quantum fluctuations requires thecalculation of the time evolution of the current of themomentum flux M αβγ (18). The method of derivationis similar to the equations obtained above. The timederivative of tensor M αβγ acts on the wave function inits definition. The time derivative of the wave func-tion is replaced by the Hamiltonian (1) in accordancewith the many-particle microscopic Schrodinger equation ı ¯ h∂ t Ψ = ˆ H Ψ. It leads to the following expression: ∂ t M αβγ = ı ¯ h Z dR N X i =1 δ ( r − r i ) 18 m i (cid:20) ˆ H ∗ Ψ ∗ · ˆ p αi ˆ p βi ˆ p γi Ψ − Ψ ∗ ˆ p αi ˆ p βi ˆ p γi ˆ H Ψ + ˆ p α ∗ i ˆ H ∗ Ψ ∗ · ˆ p βi ˆ p γi Ψ + ˆ p α ∗ i ˆ p γ ∗ i ˆ H ∗ Ψ ∗ · ˆ p βi Ψ − ˆ p α ∗ i Ψ ∗ · ˆ p βi ˆ p γi ˆ H Ψ − ˆ p α ∗ i ˆ p γ ∗ i Ψ ∗ · ˆ p βi ˆ H Ψ+ˆ p γ ∗ i ˆ H ∗ Ψ ∗ · ˆ p αi ˆ p βi Ψ − ˆ p γ ∗ i Ψ ∗ · ˆ p αi ˆ p βi ˆ H Ψ − c.c. (cid:21) (23)The part of the presented terms contains the Hamilto-nian ˆ H under action of the momentum operators. Wepermute the Hamiltonian ˆ H and the operators acting onit. Hence, the result of permutation presented by termswhere no operators act on the Hamiltonian ˆ H . However,the terms containing the corresponding commutators ap-pear. Therefore, all terms are combined in two groups: ∂ t M αβγ = ı ¯ h Z dR N X i =1 δ ( r − r i ) 18 m i (cid:20) ˆ H ∗ Ψ ∗ · ˆ p αi ˆ p βi ˆ p γi Ψ − Ψ ∗ ˆ H ˆ p αi ˆ p βi ˆ p γi Ψ + ˆ H ∗ ˆ p α ∗ i Ψ ∗ · ˆ p βi ˆ p γi Ψ + ˆ H ∗ ˆ p α ∗ i ˆ p γ ∗ i Ψ ∗ · ˆ p βi Ψ − ˆ p α ∗ i Ψ ∗ · ˆ H ˆ p βi ˆ p γi Ψ − ˆ p α ∗ i ˆ p γ ∗ i Ψ ∗ · ˆ H ˆ p βi Ψ+ ˆ H ∗ ˆ p γ ∗ i Ψ ∗ · ˆ p αi ˆ p βi Ψ − ˆ p γ ∗ i Ψ ∗ · ˆ H ˆ p αi ˆ p βi Ψ − c.c. (cid:21) + ı ¯ h Z dR N X i =1 δ ( r − r i ) 18 m i (cid:20) − Ψ ∗ [ˆ p αi ˆ p βi ˆ p γi , ˆ H ]Ψ+[ˆ p α ∗ i , ˆ H ∗ ]Ψ ∗ · ˆ p βi ˆ p γi Ψ + [ˆ p α ∗ i ˆ p γ ∗ i , ˆ H ∗ ]Ψ ∗ · ˆ p βi Ψ − ˆ p α ∗ i Ψ ∗ · [ˆ p βi ˆ p γi , ˆ H ]Ψ − ˆ p α ∗ i ˆ p γ ∗ i Ψ ∗ · [ˆ p βi , ˆ H ]Ψ+[ˆ p γ ∗ i , ˆ H ∗ ]Ψ ∗ · ˆ p αi ˆ p βi Ψ − ˆ p γ ∗ i Ψ ∗ · [ˆ p αi ˆ p βi , ˆ H ]Ψ − c.c. (cid:21) . (24) The first group of terms leads to the divergence of the fluxof tensor M αβγ . The second group of terms containingthe commutators presents the interactions.Final form of tensor M αβγ evolution equation can beexpressed in the following terms: ∂ t M αβγ + ∂ δ R αβγδ = ¯ h m n∂ α ∂ β ∂ γ V ext − m Π βγ ∂ α V ext − m Π αγ ∂ β V ext − m Π αβ ∂ γ V ext + 1 m F αβγqf + 1 m ( F αβγ + F βγα + F γαβ ) , (25)where F αβγqf = ¯ h m Z [ ∂ α ∂ β ∂ γ U ( r − r ′ )] n ( r , r ′ , t ) d r ′ (26)is the quantum force contribution leading to the quantumfluctuations, and F αβγ = − Z [ ∂ α U ( r − r ′ )]Π βγ ( r , r ′ , t ) d r ′ (27)is the interaction contribution containing nonzero limitin the classical regime, withΠ αβ ( r , r ′ , t ) = Z dR X i,j = i m i δ ( r − r i ) δ ( r ′ − r j ) ×× [Ψ ∗ ˆ p αi ˆ p βi Ψ + (ˆ p βi ) ∗ Ψ ∗ ˆ p αi Ψ + c.c. ] . (28)Tensor Π αβ ( r , r ′ , t ) can be simplified in the correlation-less regime to the following form Π αβ ( r , r ′ , t ) = Π αβ ( r , t ) · n ( r ′ , t ). However, the correlations caused by the sym-metrization of the bosonic many-particle wave functionare used below.Terms F αβγ and F αβγqf are the third rank force tensorsdescribing the interparticle interaction. However, equa-tion (25) contains the flux of tensor M αβγ which is thefourth rank tensor appearing in the following form: R αβγδ = Z dR N X i =1 δ ( r − r i ) 116 m i (cid:20) Ψ ∗ ˆ p αi ˆ p βi ˆ p γi ˆ p δi Ψ+ˆ p α ∗ i Ψ ∗ ˆ p βi ˆ p γi ˆ p δi Ψ + ˆ p β ∗ i Ψ ∗ ˆ p αi ˆ p γi ˆ p δi Ψ + ˆ p γ ∗ i Ψ ∗ ˆ p αi ˆ p βi ˆ p γi Ψ+ˆ p δ ∗ i Ψ ∗ ˆ p αi ˆ p βi ˆ p γi Ψ + ˆ p α ∗ i ˆ p δ ∗ i Ψ ∗ ˆ p βi ˆ p γi Ψ+ˆ p α ∗ i ˆ p γ ∗ i Ψ ∗ ˆ p βi ˆ p δi Ψ + ˆ p γ ∗ i ˆ p δ ∗ i Ψ ∗ ˆ p αi ˆ p βi Ψ + c.c. (cid:21) . (29)Equation (25) is obtained for bosons with the arbitrarytemperature. It can be separated on two equations fortwo following subsystems: the BEC and the normal fluid.All terms in equation (25) are additive on the particles.Therefore, they are additive on the subsystems. Hence,the structure of the partial equations is identical to thestructure of equation (25): ∂ t M αβγa + ∂ δ R αβγδa = − m n a ∂ α ∂ β ∂ γ V ext − m Π βγa ∂ α V ext − m Π αγa ∂ β V ext − m Π αβa ∂ γ V ext + 1 m F αβγa,qf + 1 m ( F αβγa + F βγαa + F γαβa ) , (30)where subindex a equal B for the BEC and n for thenormal fluid.The fourth rank kinematic tensor R αβγδa (29) has thefollowing form after the introduction of the velocity fieldvia the Madelung transformation of the many-particlewave function: R αβγδa = n a v αa v βa v γa v δa + v αa v δa ( p βγa + T βγa ) + v βa v δa ( p αγa + T αγa ) + v γa v δa ( p αβa + T αβa )+ v αa v γa ( p βδa + T βδa ) + v βa v γa ( p αδa + T αδa ) + v αa v βa ( p γδa + T γδa )+ v αa Q βγδa + v βa Q αγδa + v γa Q αβδa + v δa Q αβγa + Q αβγδa + T αβγδa + L αβγδa . (31)This structure shows some similarity to the representa-tions for the second rank tensor momentum flux (11) andfor the third rank tensor (21), where the higher rank ten-sors are partially transformed via the concentration, ve-locity field and, if possible, via tensors of smaller rank.However, this transformation is partial since there is thetensor of the equal rank, but defined in the comovingframe. Moreover, this final tensor is splitted on few parts.It is two parts for the second rank tensor momentumflux, where we have the kinetic pressure (quasi-classicalpart of thermal nature) and the quantum Bohm potential(the quantum part). There are three parts for the thirdrank tensor M αβγ . They are the quasi-classical part ofthermal nature, the quantum part, and the combinedthermal-quantum part. For the fourth rank tensor wealso have three parts: the quasi-classical part of thermalnature Q αβγδa , the quantum part T αβγδa , and the com-bined thermal-quantum part L αβγδa .Developed model shows that arbitrary quantum sys-tem can be modeled via the hydrodynamic equationswhich are traditionally associated with the fluid dynam-ics. Quantum systems demonstrates that each parti-cle shows the properties of the wave and this wave-likebehavior is incorporated in the quantum hydrodynamic model. This conclusion follows from the fact that thequantum hydrodynamic is derived from the Schrodingerequations which contains these information. These gen-eral concept is illustrated for the ultracold bosons, butthe quantum hydrodynamic method can be applied toother physical systems. This similarity between quan-tum behavior and the dynamics of fluids recently foundunusual realization. It is experimentally found that clas-sic fluid objects demonstrate the quantum-like behavior[36], [37], [38]. It is observed as the millimetric dropletwalking on the surface of vibrating fluids, where the mo-tion of droplets is affected by the resonant interactionwith their own wave field [37], [38]. Systems walkingdroplets demonstrate various quantum effects [39], [40],[41]. III. CONTRIBUTION OF INTERACTION INTHE QUANTUM HYDRODYNAMICEQUATIONS
Equations (6), (16), and (25) contain terms describinginteraction. Approximate forms of these force fields ofdifferent tensor ranks are necessary to get a truncatedset of equations. In our case, it is necessary to includethe short-range nature of the potential of the interparti-cle interaction. Moreover, the weak interaction limit isconsidered. These two assumptions are used to get sim-plified form of F α , F αβ , F αβγ and F αβγqf in this section. A. Interaction terms in the Euler equation
The short-range interaction in the Euler for the singlespecies of quantum particles can be written as the di-vergence of the symmetric quantum stress tensor F α = − ∂ β σ αβ .The first order by the interaction radius approximationgives the following expression for the quantum stress ten-sor (see also [2]) σ αβ ( r , t ) = − Z dR X i,j.i = j δ ( r − R ij ) ×× r αij r βij | r ij | ∂U ( r ij ) ∂ | r ij | Ψ ∗ ( R ′ , t )Ψ( R ′ , t ) , (32)where R ′ = { ..., R ij , ..., R ij , ... } with vector R ij locatedat i -th and j -th places.Expression (32) can be rewritten in terms of two-particle concentration σ αβ ( r , t ) = − T r ( n ( r , r ′ , t )) Z d r r α r β r ∂U ( r ) ∂r , (33)where the notion of trace is used T rf ( r , r ′ ) = f ( r , r ) . (34)Consideration of the short-range interaction leads tothe separation of integral containing the potential of in-teraction. So, the characteristic of interaction does notdepend on the motion or position of particles. This inte-gral simplifies in the following way Z r α r β r ∂U∂r d r = 13 δ αβ Z rU ′ d r = − δ αβ Z U d r . (35)The last integral in this expression is denoted as g = R U d r .The two-particle concentration can be calculated in theweakly interacting limit (see [2]) n ( r , r ′ , t ) = n ( r , t ) n ( r ′ , t )+ | ρ ( r , r ′ , t ) | + ℘ ( r , r ′ , t ) , (36)where n ( r , t ) = X f n f ϕ ∗ f ( r , t ) ϕ f ( r , t ) (37)is the expression of concentration (2) in terms of the sin-gle particle wave functions ϕ f ( r , t ), ρ ( r ′ , r , t ) = X f n f ϕ ∗ f ( r , t ) ϕ f ( r ′ , t ) (38)is the density matrix, and ℘ ( r , r ′ , t ) = X f n f ( n f − | ϕ f ( r , t ) | | ϕ f ( r ′ , t ) | , (39)The last term in equation describes interaction of pairsof particles being in the same quantum state. It can beseen from the existence of single quantum number g inall wave are single-particle wave functions.Expression (36) can be substituted in the general ex-pression of the force field (8). However, equation (8) doesnot contain information about the short-range nature ofconsidered interaction. The first and second terms arerelated to particles located in different quantum states.It cannot be seen from equation (36), but it follows fromintermediate terms which can be found in Ref. [2].The trace of the two-particle concentration enteringthe quantum stress tensor has the following form T rn ( r , r ′ , t ) ≈ n ) ′ + n B , (40)where the first term on the right-hand side symbol ’means that the product of concentrations is related tothe particles in different quantum states. Therefore, thefirst term has no n B contribution from selfaction of BEC.The dropped terms are described in Ref. [42].Present explicit contribution of the BEC concentration n B and the concentration of normal fluid n n in the firstterm on the right-hand side of equation (40):( n ) ′ = (( n B + n n )( n B + n n )) ′ = ( n n + 2 n B n n ) . (41)The last term in equation (40) appears for particlesbeing in BEC. While the first term on the right-hand side in equation (40) related to interaction of particles beingin different quantum states. Hence, it gives contributionfor the interaction between BEC and normal fluid and for the interaction between bosons belonging to normalfluid.Full expression of the quantum stress tensor for thebosons at finite temperature can be written in terms ofthe concentration of BEC and the concentration of nor-mal fluid: σ αβ = 12 gδ αβ (2 n n + 4 n B n n + n B ) . (42)If we consider dynamics of BEC or normal fluid wecannot use the notion of the quantum stress tensor σ αβ for the interaction of subspecies as it is for the interactionof different species.The first (last) term in equation (42) contains the self-action of the normal fluid (of the BEC). The secondterm in equation (42) presents the interaction betweenthe BEC and normal fluid.If we consider dynamics of BEC we need to extractforce caused by the BEC and normal fluid acting on theBEC. This force is the superposition of a part of the sec-ond term in equation (42) and the last term in equation(42): F αB = − gn B ∂ α (2 n n + n B ) . (43)The second term in equation (42) can be rewritten asfollows F α = − g ( n B ∂ α n n + n n ∂ α n B ). The first part ofthis expression is used in equation (43).If we consider dynamics of normal fluid it means thatthe source of field in the first term of n can be the normalfluid and the BEC, hence the last term gives no contri-bution in this case in equation (42): F αn = − gn n ∂ α ( n n + n B ) . (44)The nonsymmetric decomposition allows to use the no-tion of the NLSE. It is necessary condition to have theGP equation at finite temperatures. Moreover, the non-symmetric form is traditionally used in literature [32].Same chose is made at analysis j β below.
1. Nonlinear Schrodinger equations
Dropping the pressure of normal fluid and using thequantum Bohm potential in form (12) we find a closed setof hydrodynamic equations. Introducing the macroscopicwave function for both the BEC and the normal fluid forthe potential velocity fields as Φ a = √ n a e ımφ a / ¯ h , where φ a is the potential of the velocity field v a = −∇ φ a . ı ¯ h∂ t Φ B = − ¯ h ∇ m + V ext + g ( n B + 2 n n ) ! Φ B , (45)and ı ¯ h∂ t Φ n = − ¯ h ∇ m + V ext + 2 g ( n B + n n ) ! Φ n . (46)The kinetic energy (the first term on the right-handside of equations (45) and (46)) corresponds to the appli-cation of the noninteracting limit for the quantum Bohmpotential for the BEC and for the normal fluid.The pressure of the normal fluid is dropped in equation(46).Equations (45), (46) correspond to equations 127-129given in Ref. [1] while there is a difference in the form ofpresentation.Therefore, the account of the pressure evolution to-gether with the pressure flux evolution gives the gen-eralization of the model presented with the nonlinearSchrodinger equations (45), (46). Necessity of additionalequations is demonstrated in Refs. [4], [43], [44] if onewants to include the quantum fluctuations. B. Interaction terms in the pressure evolutionequation
General form of the pressure evolution equation con-tains the interaction via the force second rank tensorfield. Its main contribution is obtained in the first or-der by the interaction radius. The result appears in thefollowing form F αβ ( r , t ) = 18 m ∂ γ Z dR X i,j ; i = j δ ( r − R ij ) ×× r βij r γij | r ij | ∂U ( r ij ) ∂ | r ij | (cid:20) Ψ ∗ ( R ′ , t )(ˆ p α (1) + ˆ p α (2) )Ψ( R ′ , t ) + c.c. (cid:21) − m Z dR X i,j ; i = j δ ( r − R ij ) r αij r γij | r ij | ∂U ( r ij ) ∂ | r ij | ×× (cid:20) ( ∂ γ (1) − ∂ γ (2) )Ψ ∗ ( R ′ , t )(ˆ p α (1) − ˆ p α (2) )Ψ( R ′ , t )+Ψ ∗ ( R ′ , t )( ∂ γ (1) − ∂ γ (2) )(ˆ p α (1) − ˆ p α (2) )Ψ( R ′ , t ) + c.c. (cid:21) . (47)Form (47) appears at the expansion of the force tensorfield (17) using the short-range nature of interaction (see[2] for the method described for the force field, or [35] forapplication of this method to fermions).For the force tensor field F αβ we can present the inter-mediate expressions like equations (33) and (36) obtainedfor the force field F α = − ∂ β σ αβ . However, similar ex-pressions obtained for F αβ are rather large. Hence, westart the presentation with equation similar to equation(40) obtained after taking trace of the intermediate ex-pressions.Therefore, we obtain the following simplification ofequation (47) for the force tensor field F αβ : F αβ = − g m ∂ β [2( n Λ α ) ′ + n B Λ αB ] + ı ¯ h g m [2( nr αβ ) ′ − α Λ β ) ′ + n B r αβB − Λ αB Λ βB ] + c.c., (48)where we use the intermediate functions Λ α and r αβ withthe following definitions:Λ α = X f n f ϕ ∗ f ˆ p α ϕ f = mj α − ı ¯ h ∂ α n, (49)and r αβ = X f n f ϕ ∗ f ˆ p α ˆ p β ϕ f = m (cid:18) nv α v β + p αβ − ¯ h m X f n f a f ∂ α ∂ β a f (cid:19) − ı m ¯ h ∂ α ( nv β ) + ∂ β ( nv α )] . (50)The calculation of functions Λ α and r αβ includes theMadelung transformation of the single-particle wavefunctions ϕ f ( r, t ) = √ a f e ıS f . Next, we use the fol-lowing definitions of the velocity field and the pres-sure tensor in terms of the single-particle wave functions nv α = P f n f a f (¯ h∂ α S f /m ), and p αβ = P f n f a f u αf u βf ,where u αf = (¯ h∂ α S f /m ) − v α .Let us represent terms like ( n Λ α ) ′ in the explicit form: F αβ = − g m ∂ β [2 n n Λ αn +2 n n Λ αB + 2 n B Λ αn + n B Λ αB ]+ ı ¯ h g m [2 n n r αβn − αn Λ βn + 2 n n r αβB − αn Λ βB n B r αβn − αB Λ βn + n B r αβB − Λ αB Λ βB ] + c.c.. (51)Further calculation gives the representation of tensor F αβ in term of hydrodynamic functions: F αβ = − g m ∂ β [2 n n j αn + 2 n n j αB + 2 n B j αn + n B j αB ]+ g m (cid:20) n n ( ∂ j n + ∂ j n ) − j n ∂ n n − j n ∂ n n +2 n n ( ∂ α j βB + ∂ β j αB ) − j αB ∂ β n n − j βB ∂ α n n +2 n B ( ∂ α j βn + ∂ β j αn ) − j αn ∂ β n B − j βn ∂ α n B + n B ( ∂ α j βB + ∂ β j αB ) − j αB ∂ β n B − j βB ∂ α n B (cid:21) . (52)The momentum flux evolution equation contains thesymmetric combination of the force tensor fields F αβ : F αβ + F βα = − gm (cid:20) j αn ∂ β n n + j βn ∂ α n n )+2( j αn ∂ β n B + j βn ∂ α n B )+2( j αB ∂ β n n + j βB ∂ α n n ) + j αB ∂ β n B + j βB ∂ α n B (cid:21) . (53)The zero temperature analysis demonstrates that thereis nonzero pressure for the BECs, caused by the quantumfluctuations entering the set of hydrodynamic equationsvia the evolution of the pressure flux [4], [43], [44]. Thepressure also exists for the normal fluid. So, we makedecomposition of the momentum flux evolution equationon two partial equations for Π αβn and Π αβB . Formallythis decomposition is presented with equation (20). Tocomplete this procedure we need to split the force tensorfield F αβ + F βα = F αβB + F βαB + F αβn + F βαn , where F αβB + F βαB = − gm (cid:20) j αB ∂ β n n + j βB ∂ α n n )+ j αB ∂ β n B + j βB ∂ α n B (cid:21) , (54)and F αβn + F βαn = − gm (cid:20) j αn ∂ β n n + j βn ∂ α n n ) +2( j αn ∂ β n B + j βn ∂ α n B ) (cid:21) . (55)After extraction of the pressure tensor p αβ from themomentum flux evolution Π αβ we have extra contribu-tion of the interaction in the pressure evolution equationin compare with equations (20). It contains the followingcontribution F αβ − v β F α .Using equations (43), (44), (55), (54) find F αβ + F αβ − v α F β − v β F α = 0 for the BECs and for the normal fluid.A pressure evolution equation is used in [45] for bosonsabove the critical temperature. Equation 4 of Ref. [45]contains the force in the following form n v · F whichgenerally differs from F αβ − v β F α obtained above. C. The short-range interaction in the third ranktensor evolution equation
The third rank tensor M αβγ (25) evolution equationcontains two kinds of the third rank force tensors F αβγ (27) and F αβγqf (26). Consider them separately.
1. Quasi-classical third rank force tensor
Tensor F αβγqf (26) is proportional to the Planck con-stant, so it goes to zero in the classical limit. The thirdrank force tensor F αβγ (27) is different, it has nonzerolimit in the quasiclassical regime. However, we are inter-ested in the value of tensor F αβγ (27) in one of quantumregimes for the degenerate bosons.We calculate the third rank force tensor F αβγ (27) inthe first order by the interaction radius appears in thefollowing form F αβγ ( r , t ) = 18 m ∂ µ Z dR X i,j ; i = j δ ( r − R ij ) r µij r γij | r ij | ∂U ( r ij ) ∂ | r ij | (cid:20) Ψ ∗ ( R ′ , t )ˆ p α (1) ˆ p β (1) Ψ( R ′ , t ) + ˆ p β ∗ (1) Ψ ∗ ( R ′ , t )ˆ p α (1) Ψ( R ′ , t ) + c.c. (cid:21) − m Z dR X i,j ; i = j δ ( r − R ij ) r µij r γij | r ij | ∂U ( r ij ) ∂ | r ij | (cid:20) ( ∂ µ (1) − ∂ µ (2) )Ψ ∗ ( R ′ , t )ˆ p α (1) ˆ p β (1) Ψ( R ′ , t )+Ψ ∗ ( R ′ , t )( ∂ µ (1) − ∂ µ (2) )ˆ p α (1) ˆ p β (1) Ψ( R ′ , t )+ ˆ p β ∗ (1) ( ∂ µ (1) − ∂ µ (2) )Ψ ∗ ( R ′ , t )ˆ p α (1) Ψ( R ′ , t )+ ˆ p β ∗ (1) Ψ ∗ ( R ′ , t )( ∂ µ (1) − ∂ µ (2) )ˆ p α (1) Ψ( R ′ , t )+ c.c. (cid:21) . (56)Here, the part of expression for F αβγ containing the interaction potential appears as the independent multiplier. Ithas same form as the integral in the Euler equation (35). Hence, tensor F αβγ is proportional to the Groos-Pitaevskiiinteraction constant.Further calculation in the weakly interacting limit, following the method described in Ref. [2], gives an intermediaterepresentation of the third rank force tensor: F αβγ ( r , t ) = − g m " Π αβB ∂ γ n B + Π αβ ∂ γ n + (cid:18) n X f n f ∂ γ ϕ ∗ f ˆ p α ˆ p β ϕ f + ı ¯ h Λ γ r αβ − ı ¯ h Λ α ∗ κ γβ + ı ¯ h Λ β κ αγ + c.c. (cid:19) , (57)0where κ αβ = X f n f p α ϕ ∗ f · ˆ p β ϕ f . (58)Function κ αβ is the nonsymmetric tensor. It has symmetric real part and the antisymmetric imaginary part: κ αβ = m (cid:18) nv α v β + p αβ + ¯ h m X f n f ∂ α · a f ∂ β a f (cid:19) − ım ¯ h [ v α ∂ β n − v β ∂ α n + X f n f a f ( u α ∂ β a f − u β ∂ α a f )] . (59)No specific notation is introduced for the third rank tensor P f n f ∂ γ ϕ ∗ f ˆ p α ˆ p β ϕ f . In our calculations we need itsimaginary part multiplied by 2: X f n f ∂ γ ϕ ∗ f ˆ p α ˆ p β ϕ f + c.c. = 2 m " ∂ γ n · v α v β − ∂ α n · v β v γ − ∂ β n · v α v γ − nv γ ( ∂ β v α + ∂ α v β )+ X f n f a f ( ∂ γ a f ) u αf u βf − ∂ β p αγ − ∂ α p βγ + 12 X f n f a f ( u βf ∂ α u γf + u αf ∂ β u γf )+ v α X f n f a f ( ∂ γ a f · u βf − ∂ β a f · u γf ) + v β X f n f a f ( ∂ γ a f · u αf − ∂ α a f · u γf ) − ¯ h m X f n f ∂ γ a f · ∂ α ∂ β a f (cid:19) . (60)Equation (57) includes term Π αβB ∂ γ n B which describes full contribution of the BEC in F αβγ .The second term in equation (57) is an analog of the first term in equation (36). All following terms in equation(57) are the analog of the second term in equation (36). It can be interpreted as the exchange interaction.Further calculation of the (57) gives the partially truncated expression mainly presented via the macroscopichydrodynamic functions: F αβγ ( r , t ) = − g (cid:20) αβB ∂ γ n B +4Π αβ ∂ γ n + ∂ α n (cid:18) Π βγ + ¯ h m ∂ β ∂ γ n (cid:19) + ∂ β n (cid:18) Π αγ + ¯ h m ∂ α ∂ γ n (cid:19) + ∂ γ n (cid:18) Π αβ − ¯ h m ∂ α ∂ β n (cid:19) + n (cid:18) ∂ γ n · v α v β − ∂ α n · v β v γ − ∂ β n · v α v γ − nv γ ( ∂ α v β + ∂ β v α ) − ∂ β p αγ − ∂ α p βγ + X f n f ( ∂ γ a f ) u αf u βf + X f n f a f ( u βf ∂ α u γf + u αf ∂ β u γf ) + 32 v α X f n f ( ∂ γ a f · u βf − ∂ β a f · u γf )+ 32 v β X f n f ( ∂ γ a f · u αf − ∂ α a f · u γf ) − h m X f n f ∂ γ a f · ∂ α ∂ β a f (cid:19)(cid:21) . (61)Equation for the third rank tensor M αβγ evolution contains the symmetric combination of the third rank forcetensors (61) which cam be presented in the following form: F αβγ + F βγα + F γαβ = − g (cid:20) βγB ∂ α n B + Π αγB ∂ β n B + Π αβB ∂ γ n B ) + 4(Π βγ ∂ α n + Π αγ ∂ β n + Π αβ ∂ γ n )+ ∂ α n (cid:18) βγ + ¯ h m ∂ β ∂ γ n (cid:19) + ∂ β n (cid:18) αγ + ¯ h m ∂ α ∂ γ n (cid:19) + ∂ γ n (cid:18) αβ + ¯ h m ∂ α ∂ β n (cid:19) − n ( ∂ α Π βγ + ∂ β Π αγ + ∂ γ Π αβ ) − h m n∂ α ∂ β ∂ γ n (cid:21) . (62)The pressure flux evolution equation obtained as the reduction of the third rank tensor M αβγ evolution equationcontains the following combination of the force fields: F αβγ + F βγα + F γαβ − mn ( F α Π βγ + F β Π αγ + F γ Π αβ ) = g ∂ α ( n Π βγ ) + ∂ β ( n Π αγ ) + ∂ γ ( n Π αβ )] ′ + g ¯ h m [3 n∂ α ∂ β ∂ γ n − ∂ α n · ∂ β ∂ γ n − ∂ β n · ∂ α ∂ γ n − ∂ γ n · ∂ α ∂ β n ] ′ , (63)where symbol [] ′ specifies that product of functions describing the BEC is excluded similarly equations (40) and (41).The quantum part can be represented F αβγ + F βγα + F γαβ − mn ( F α Π βγ + F β Π αγ + F γ Π αβ ) = g ∂ α ( n Π βγ ) + ∂ β ( n Π αγ ) + ∂ γ ( n Π αβ )] ′ + g ¯ h m [ ∂ α ( n∂ β ∂ γ n − ∂ β n · ∂ γ n ) + ∂ β ( n∂ α ∂ γ n − ∂ α n · ∂ γ n ) + ∂ γ ( n∂ α ∂ β n − ∂ α n · ∂ β n )] ′ . (64)It can be considered as F αβγ + F βγα + F γαβ − mn ( F α Π βγ + F β Π αγ + F γ Π αβ ) = ˜ F αβγ + ˜ F βγα + ˜ F γαβ , where˜ F αβγ = g ∂ α ( n Π βγ ) + g ¯ h m [ ∂ α ( n∂ β ∂ γ n − ∂ β n · ∂ γ n )] ′ . (65)It is the derivative of the second rank tensor.If we consider the zero temperature limit we find F αβγ = ( − g/ m )Π αβB ∂ γ n B . The transition from the equation ofevolution for tensor M αβγ to the equation of evolution for tensor Q αβγ , which is the sibling of M αβγ , but the pressureflux Q αβγ is defined in the comoving frame leads to the canceling of such term. Hence, the nonzero contribution comesfrom the quantum part of the third rank force tensor F αβγqf . The nonzero temperature gives a nonzero contributionof F αβγ at the transition to the pressure flux evolution equation.Equation (65) is obtained for all bosons. We need to separate it on the force acting on the BEC and the forceacting on the normal fluid.Finally, we obtain ˜ F αβγn + ˜ F βγαn + ˜ F γαβn = g ∂ α ( n n Π βγn ) + ∂ β ( n n Π αγn ) + ∂ γ ( n n Π αβn )+ n n ∂ α (Π βγB + ∂ β Π αγB + ∂ γ Π αβB ) + Π βγn ∂ α n B + Π αγn ∂ β n B + Π αβn ∂ γ n B ]+ g ¯ h m (cid:20) n n ∂ α ∂ β ∂ γ n n + 3 n n ∂ α ∂ β ∂ γ n B − ∂ α n n · ∂ β ∂ γ n n − ∂ β n n · ∂ α ∂ γ n n − ∂ γ n n · ∂ α ∂ β n n − ∂ α n n · ∂ β ∂ γ n B − ∂ β n n · ∂ α ∂ γ n B − ∂ γ n n · ∂ α ∂ β n B − ∂ α n B · ∂ β ∂ γ n n − ∂ β n B · ∂ α ∂ γ n n − ∂ γ n B · ∂ α ∂ β n n (cid:21) , (66)and ˜ F αβγB + ˜ F βγαB + ˜ F γαβB = g n B ∂ α (Π βγn + ∂ β Π αγn + ∂ γ Π αβn )+Π βγB ∂ α n n + Π αγB ∂ β n n + Π αβB ∂ γ n n ] + g ¯ h m (cid:20) n B ∂ α ∂ β ∂ γ n n − ∂ α n B · ∂ β ∂ γ n n − ∂ β n B · ∂ α ∂ γ n n − ∂ γ n B · ∂ α ∂ β n n − ∂ α n n · ∂ β ∂ γ n B − ∂ β n n · ∂ α ∂ γ n B − ∂ γ n n · ∂ α ∂ β n B (cid:21) . (67)Equations (66) and (67) gives final expressions for the quasi-classic force fields in the pressure flux evolution equationsfor two fluid model.
2. The third rank force tensor describing the quantumfluctuation
The quantum fluctuations in the zero temperatureBECs is caused by tensor F αβγqf (26). Its major con- tribution can be found in the first order by the interac-2tion radius approximation. Here, we consider the smallnonzero temperature regime of F αβγqf for the bosons. So,we obtain its generalization for the two-fluid model. Thequantum third rank force tensor F αβγqf (26) is calculatedin the first order by the interaction radius F αβγqf = − ¯ h m ∂ δ Z dR X i,j.i = j δ ( r − R ij ) ×× r δij ∂ αi ∂ βi ∂ γi U ( r ij )Ψ ∗ ( R ′ , t )Ψ( R ′ , t ) . (68)In formula (68) for tensor F αβγqf we have separation ofthe integral containing the interaction potential, as wehave at the calculation of other force fields above. How-ever, here we obtain different integral R r α ∂ β ∂ γ ∂ δ . Cal-culation of this integral leads to the second interactionconstant given below in the simplified expression for thequantum third rank force tensor F αβγqf = ¯ h m g I αβγδ ∂ δ T rn ( r , r ′ , t ) , (69)where g = 23 Z d r U ′′ ( r ) , (70)and I αβγδ = δ αβ δ γδ + δ αγ δ βδ + δ αδ δ βγ . (71)Calculation of the two-particle concentration leads tothe following expression: F αβγqf = ¯ h m g I αβγδ ∂ δ (2 n n + 4 n B n n + n B ) . (72)Here we have ∂ δ (2 n n + 4 n B n n + n B ). Next, we openbrackets and find (4 n n ∂ δ n n + 4 n n ∂ δ n B + 4 n B ∂ δ n n +2 n B ∂ δ n B ). The source of field is under derivative, whilethe multiplier in front corresponds to the species underthe action of field. Hence the first two terms correspondto the force acting on the normal fluid, while the thirdand fourth terms correspond to the force acting on theBEC.Therefore, we can separate expression (72) on twoparts corresponding to the BEC and to the normal fluid: F αβγqf,B = ¯ h m g I αβγδ (2 n B ∂ δ n n + n B ∂ δ n B ) , (73)and F αβγqf,n = ¯ h m g I αβγδ ∂ δ ( n n ∂ δ n n + n B ∂ δ n n ) . (74) IV. HYDRODYNAMIC EQUATIONS FOR TWOFLUID MODEL OF BOSONS WITH NONZEROTEMPERATURE
This section provides the final set of equations obtainedin this paper. The method of the introduction of the ve-locity field and corresponding representation of the hy-drodynamic equations is not described in this paper. Itcan be found in number other papers, majority of detailsare given in Refs. [2], [35].In this regime we have two continuity equations: ∂ t n B + ∇ · ( n B v B ) = 0 , (75)and ∂ t n n + ∇ · ( n n v n ) = 0 . (76)The Euler equation for bosons in the BEC state mn B ( ∂ t + v B · ∇ ) v αB + ∂ β T αβB + ∂ β p αβqf + gn B ∂ α n B = − n B ∂ α V ext − gn B ∂ α n n , (77)where the quantum Bohm potential is given by equation(12) in the noninteracting limit.The Euler equation for bosons in the excited statescorresponding to the nonzero temperature mn n ( ∂ t + v n · ∇ ) v αn + ∂ β p αβn,eff +2 gn n ∂ α n n = − n n ∂ α V ext − gn n ∂ α n B , (78)where the effective pressure tensor p αβn,eff = p αβn + T αβn .The effective pressure evolution equation for normalboson fluid is also a part of developed and applied hy-drodynamic model ∂ t p αβn,eff + v γn ∂ γ p αβn,eff + p αγn,eff ∂ γ v βn + p βγn,eff ∂ γ v αn + p αβn,eff ∂ γ v γn + ∂ γ T αβγn + ∂ γ Q αβγn = 0 . (79)Moreover, we have the pressure (the quantum Bohmpotential) p αβB,eff = T αβB + p αβqf evolution equation for theBEC ∂ t p αβB,eff + v γB ∂ γ p αβB,eff + p αγB,eff ∂ γ v βB + p βγB,eff ∂ γ v αB + p αβB,eff ∂ γ v γB + ∂ γ T αβγB + ∂ γ Q αβγqf = 0 . (80)Let us to point out the following property of the quan-tum Bohm potential that it satisfies the following equa-tion for the arbitrary species a∂ t T αβa + v γa ∂ γ T αβa + T αγa ∂ γ v βa + T βγa ∂ γ v αa T αβa ∂ γ v γa + ∂ γ T αβγa = 0 . (81)It is expected that approximate form of the quantumBohm potential (12) satisfies equation (81) existing atthe zero interaction. Hence, substitute (12) in equation(81) with the zero right-hand side: ∂ β ∂ γ n a · ( ∂ γ v αa − ∂ α v γa ) + ∂ α ∂ γ n a · ( ∂ γ v βa − ∂ β v γa )+ 13 ∂ γ n a · ( ∂ β ∂ γ v αa + ∂ α ∂ γ v βa − ∂ α ∂ β v γa )+ 13 n a [ △ ( ∂ β v αa + ∂ α v βa ) − ∂ α ∂ β ( ∇ · v a )] = 0 , (82)where the continuity equation is used for the time deriva-tives of concentration. Make the condition of the po-tentiality of the velocity field v a = ∇ φ a . Include it inequation (82) and find that this equation is satisfied.Hence, we obtain the simplified form of the pres-sure evolution equations, where the traditional quantumBohm potential is extracted: ∂ t p αβqf,B + v γB ∂ γ p αβqf,B + p αγqf,B ∂ γ v βB + p βγqf,B ∂ γ v αB + p αβqf,B ∂ γ v γB + ∂ γ Q αβγqf,B = 0 , (83)and ∂ t p αβn + v γn ∂ γ p αβn + p αγn ∂ γ v βn + p βγn ∂ γ v αn + p αβn ∂ γ v γn + ∂ γ Q αβγn = 0 . (84)Equation for the evolution of quantum-thermal part ofthe third rank tensor is [4], [43]: ∂ t Q αβγqf + ∂ δ ( v δB Q αβγqf )+ Q αγδqf ∂ δ v βB + Q βγδqf ∂ δ v αB + Q αβδqf ∂ δ v γB = ¯ h m g I αβγδ (cid:18) n B ∂ δ n B + 2 n B ∂ δ n n ) (cid:19) − m n B ∂ α ∂ β ∂ γ V ext + ˜ F αβγB + ˜ F βγαB + ˜ F γαβB + 1 mn ( p αβqf,eff ∂ δ p γδqf,eff + p αγqf,eff ∂ δ p βδqf,eff + p βγqf,eff ∂ δ p αδqf,eff ) , (85)and ∂ t Q αβγn + ∂ δ ( v δn Q αβγn )+ Q αγδn ∂ δ v βn + Q βγδn ∂ δ v αn + Q αβδn ∂ δ v γn = ¯ h m g I αβγδ (cid:18) n n ∂ δ n n + n n ∂ δ n B (cid:19) − m n n ∂ α ∂ β ∂ γ V ext + ˜ F αβγn + ˜ F βγαn + ˜ F γαβn + 1 mn ( p αβn,eff ∂ δ p γδn,eff + p αγn,eff ∂ δ p βδn,eff + p βγn,eff ∂ δ p αδn,eff ) , (86)where ˜ F αβγa is not presented explicitly since equations(66) and (67) show that required expressions are ratherlarge.Refs. [2], [3]. Hydrodynamic model for fermions withpressure evolution is derived in Refs. [35], [46], [47].Terms proportional to p αβa,eff ∂ δ p γδa,eff appears in thepressure flux evolution equation, but it leads to the con-tribution beyond the chosen approximation [48], [49].Term containing the external potential − m n a ∂ α ∂ β ∂ γ V ext goes to zero for the parabolictrap. However, it can give some nontrivial contributionfor other form of potentials. V. BEC DYNAMICS UNDER THE INFLUENCEOF THE QUANTUM FLUCTUATIONS
Developed model shows that there is nontrivial evolu-tion equation for the pressure and the pressure flux ofthe BEC. Therefore, the well-known model of BEC is ex-tended in spite the fact that the kinetic pressure tensoris expected to be equal to zero due to zero temperature.However, the quantum fluctuations lead to the nonzerooccupation numbers for the excited states.If we need to consider pure BEC we need to drop thecontribution of the normal fluid in the model presentedabove. Therefore, let us summarize the BEC model inparabolic traps: ∂ t n B + ∇ · ( n B v B ) = 0 , (87) mn B ( ∂ t + v B · ∇ ) v αB + ∂ β ( p αβqf + T αβB )+ gn B ∂ α n B + n B ∂ α V ext = 0 , (88) ∂ t p αβqf + v γB ∂ γ p αβqf + p αγqf ∂ γ v βB + p βγqf ∂ γ v αB + p αβqf ∂ γ v γB + ∂ γ Q αβγqf = 0 , (89)and ∂ t Q αβγqf + ∂ δ ( v δ Q αβγqf ) + Q αγδqf ∂ δ v β + Q βγδqf ∂ δ v α + Q αβδqf ∂ δ v γ = ¯ h m g I αβγδ n∂ δ n. (90)This simplified model is reported earlier in Refs. [4], [43],[44], where the dipole-dipole interaction is also consid-ered. It has been demonstrated that the quantum fluc-tuations gives mechanisms for the instability of the small4amplitude perturbations [4]. Moreover, the dipolar partof the quantum fluctuations creates conditions for thebright soliton in the repulsive BECs [44]. The developedin previous section model provides proper generalizationof earlier model giving the small temperature contribu-tion. VI. CONCLUSION
Revision of the two-fluid model for the finite temper-ature ultracold bosons has been presented through thederivation of the number of particles, the momentum,the momentum flux, and the third rank tensor balanceequations. The derivation has been based on the traceof the microscopic dynamics of quantum particles via theapplication of the many-particle microscopic Schrodingerequation. Hence, the microscopic Schrodinger equationdetermines the time evolution for the macroscopic func-tions describing the collective motion of bosons.General equations have been derived for the arbitrarystrength of interaction and the arbitrary temperatures.The set of equations has been restricted by the thirdrank kinematic tensor (the flux of pressure). The trun-cation is made for the low temperatures weakly interact-ing bosons after the derivation of the general structure ofhydrodynamic equations. Therefore, the thermal part ofthe fourth rank kinematic tensor has been taken equal tozero. Next, the terms containing the interaction poten-tial have been considered for the short-range interaction.The small radius of interaction provides the small pa-rameter for the expansion. The expansion is made in the force field in the Euler equation, the force tensor field inthe momentum flux equation, and the third rank forcetensor in the pressure flux evolution equation. The firstterm in expansion on the small interparticle distance hasbeen considered in each expansion, which corresponds tothe first order by the interaction radius.This model allows to gain the quantum fluctuationswhich is the essential property of the BECs. More-over, the interaction causing the quantum fluctuationshas been consider at the finite temperature.The functions obtained in the first order by the inter-action radius have been expressed via the trace of two-particle functions. The two-particle functions have beencalculated for the weakly interacting bosons.The single species of 0-spin bosons has been consid-ered. Therefore, the single fluid hydrodynamics has beenderived. Next, it has been included that the concentra-tion of particles, the current of particles (the momentumdensity), the momentum flux, and the current of the mo-mentum flux are additive functions. Consequently, theycan be easily splitted on two parts: the BEC and thenormal fluid of bosons (the non BEC part). Hence, thetwo fluid model of single species of bosons is obtained.This separation on two fluids has been made in generalform of equations.
VII. ACKNOWLEDGEMENTS
Work is supported by the Russian Foundation for BasicResearch (grant no. 20-02-00476). This paper has beensupported by the RUDN University Strategic AcademicLeadership Program. [1] F. Dalfovo, S. Giorgini, L. P. Pitaevskii, and S. Stringari,Rev. Mod. Phys. , 463 (1999).[2] P. A. Andreev, L. S. Kuz’menkov, Phys. Rev. A ,053624 (2008).[3] P. A. Andreev, Laser Phys. , 035502 (2019).[4] P. A. Andreev, arXiv:2005.13503 (accepted to Chaos).[5] T. D. Lee, K. Huang, and C. N. Yang, Phys. Rev. ,1135 (1957).[6] L. Pitaevskii and S. Stringari, Phys. Rev. Lett. , 4541(1998).[7] E. Braaten and J. Pearson, Phys. Rev. Lett. , 255(1999).[8] G. E. Astrakharchik, R. Combescot, X. Leyronas, and S.Stringari, Phys. Rev. Lett. , 030404 (2005).[9] K. Xu, Y. Liu, D. E. Miller, J. K. Chin, W. Setiawan,and W. Ketterle, Phys. Rev. Lett. , 180405 (2006).[10] A. Altmeyer, S. Riedl, C. Kohstall, M. J. Wright, R.Geursen, M. Bartenstein, C. Chin, J. Hecker Denschlag,and R. Grimm, Phys. Rev. Lett. , 040401 (2007).[11] S. B. Papp, J. M. Pino, R. J. Wild, S. Ronen, C. E.Wieman, D. S. Jin, and E. A. Cornell, Phys. Rev. Lett. , 135301 (2008).[12] H. Kadau, M. Schmitt, M. Wenzel, C. Wink, T. Maier,I. Ferrier-Barbut, T. Pfau, Nature , 194 (2016). [13] I. Ferrier-Barbut, H. Kadau, M. Schmitt, M. Wenzel, andT. Pfau, Phys. Rev. Lett. , 215301 (2016).[14] D. Baillie, R. M. Wilson, R. N. Bisset, and P. B. Blakie,Phys. Rev. A , 021602(R) (2016).[15] R. N. Bisset, R. M. Wilson, D. Baillie, P. B. Blakie, Phys.Rev. A , 033619 (2016).[16] F. Wachtler and L. Santos, Phys. Rev. A , 061603R(2016).[17] F. Wachtler and L. Santos, Phys. Rev. A , 043618(2016).[18] P. B. Blakie, Phys. Rev. A , 033644 (2016).[19] A. Boudjemaa and N. Guebli, Phys. Rev. A , 023302(2020).[20] V. Heinonen, K. J. Burns, and J. Dunkel, Phys. Rev. A , 063621 (2019).[21] B. A. Malomed, Physica D , 108 (2019).[22] E. Shamriz, Z. Chen, and B. A. Malomed, Phys. Rev. A , 063628 (2020).[23] Z. Li, J.-S. Pan, and W. Vincent Liu, Phys. Rev. A ,053620 (2019).[24] E. Aybar and M. O. Oktel, Phys. Rev. A , 013620(2019).[25] P. Examilioti, and G. M. Kavoulakis, J. Phys. B: At. Mol.Opt. Phys. , 175301 (2020). [26] T. Miyakawa, S. Nakamura, H. Yabu, Phys. Rev. A ,033613 (2020).[27] F. Bottcher, Jan-Niklas Schmidt, J. Hertkorn, Kevin S.H. Ng, Sean D. Graham, M. Guo, T. Langen, and T.Pfau, arXiv:2007.06391.[28] R. N. Bisset, L. A. P. Ardila, and L. Santos, Phys. Rev.Lett. , 025301 (2021).[29] Y. Wang, L. Guo, S. Yi, and T. Shi, Phys. Rev. Research , 043074 (2020).[30] M. J. Edmonds, T. Bland, and N. G. Parker, J. Phys.Commun. , 125008 (2020)[31] D. Baillie and P. B. Blakie, Phys. Rev. A , 043606(2020).[32] A. Griffin, Phys. Rev. B , 9341 (1996).[33] L. S. Kuz’menkov, S. G. Maksimov, and V. V. Fe-doseev, Theor. Math. Fiz. , 136 (2001) [Theoreticaland Mathematical Physics , 110 (2001)].[34] P. A. Andreev, L. S. Kuz’menkov, Prog. Theor. Exp.Phys. , 053J01 (2019).[35] P. A. Andreev, arXiv:2001.02764.[36] J. W. M. Bush, Y. Couder, T. Gilet, P. A. Milewski, andA. Nachbin, Chaos , 096001 (2018). [37] Y. Couder, S. Protiere, E. Fort, and A. Boudaoud, Na-ture , 208 (2005).[38] J. W. M. Bush, Annu. Rev. Fluid Mech. , 269–292(2015).[39] T. Cristea-Platon, P. J. Saenz, and J. W. M. Bush, Chaos , 096116 (2018).[40] A. Chowdury, A. Ankiewicz, N. Akhmediev, and W.Chang, Chaos , 123116 (2018).[41] N. B. Budanur, and M. Fleury, Chaos , 013122 (2019).[42] P. A. Andreev, Int. J. Mod. Phys. B , 1350017 (2013).[43] P. A. Andreev, arXiv:2007.15045.[44] P. A. Andreev, arXiv:2009.12720.[45] G. M. Kavoulakis, C. J. Pethick, and H. Smith, Phys.Rev. A , 2938 (1998).[46] P. A. Andreev, arXiv:1912.00843.[47] P. A. Andreev, K. V. Antipin, M. I. Trukhanova, LaserPhys. , 015501 (2021).[48] I. Tokatly, O. Pankratov, Phys. Rev. B , 15550 (1999).[49] I. V. Tokatly, O. Pankratov, Phys. Rev. B62