Two-term expansion of the ground state one-body density matrix of a mean-field Bose gas

We consider the homogeneous Bose gas on a unit torus in the mean-field regime when the interaction strength is proportional to the inverse of the particle number. In the limit when the number of particles becomes large, we derive a two-term expansion of the one-body density matrix of the ground state. The proof is based on a cubic correction to Bogoliubov's approximation of the ground state energy and the ground state.


Introduction
We consider a homogeneous system of N bosons on the unit torus T d , for any dimension d ≥ 1. The system is governed by the mean-field Hamiltonian which acts on the bosonic Hilbert space . Here the kinetic operator −∆ is the usual Laplacian (with periodic boundary conditions). The interaction potential w is a real-valued, even function. We assume that its Fourier transform is non-negative and integrable, namely w(x) = p∈2πZ d w(p)e ip·x with 0 ≤ w ∈ ℓ 1 (2πZ d ).
In particular, w is bounded. Since w is even, w is also even.
Under the above conditions, H N is well defined on the core domain of smooth functions. Moreover, it is well-known that H N is bounded from below and can be extended to be a self-adjoint operator by Friedrichs' method. The self-adjoint extension, still denoted by H N , has a unique ground state Ψ N (up to a complex phase) which solves the variational problem Here ·, · is the inner product in H N . 1 In the present paper, we are interested in the asymptotic behavior of the ground state Ψ N ∈ H N of H N in the limit when N → ∞. More precisely, we will focus on the one-body Date: March 23, 2021. 1 We use the convention that ·, · is linear in the second argument and anti-linear in the first. density matrix γ (1) Ψ N which is a trace class operator on L 2 (T d ) with kernel γ (1) Note that γ (1) Ψ N ≥ 0 and Trγ (1) Ψ N = N . 1.1. Main result. Our main theorem is Theorem 1 (Ground state density matrix). Assume that 0 ≤ w ∈ ℓ 1 ((2πZ) d ). Then the ground state Ψ N of the Hamiltonian H N in (1) satisfies Here |u u| is the orthogonal projection on u. We use the bra-ket notation, where |u = u is a vector in the Hilbert space H and u| is an element in the dual space of H which maps any vector v ∈ H to the inner product u, v H .
To the leading order, our result implies Bose-Einstein condensation, namely lim N →∞ 1 N γ (1) Ψ N = |u 0 u 0 | in the trace norm. This result is well-known and it follows easily from Onsager's inequality 1 (see [18]). The significance of Theorem 1 is that it gives the next order correction to γ (1) Ψ N , thus justifying Bogoliubov's approximation in a rather strong sense as we will explain.

Bogoliubov's approximation.
It is convenient to turn to the grand canonical setting. Let us introduce the Fock space For any Fock space vector Ψ = (Ψ n ) ∞ n=0 ∈ F with Ψ n ∈ H n , we define its norm by and define the particle number expectation by In particular, the vacuum state |0 = (1, 0, 0, ...) is a normalized vector on Fock space which has the particle number expectation 0|N |0 = 0.
We will denote by a * p and a p the creation and annihilation operators with momentum p ∈ 2πZ d , namely a * p = a * (u p ), a p = a(u p ), u p (x) = e ip·x . They satisfy the canonical commutation relation (CCR) The creation and annihilation operators can be used to express several operators on Fock space. For example, the number operator can be written as Similarly, the Hamiltonian H N in (1) can be rewritten as The right side of (4) is an operator on Fock space, which coincides with (1) when being restricted to H N . In the following we will only use the grand-canonical formula (4).
In 1947, Bogoliubov [4] suggested a heuristic argument to compute the low-lying spectrum of the operator H N by using a perturbation around the condensation. Roughly speaking, he proposed to first substitute all operators a 0 and a * 0 in (4) by the scalar number √ N (c-number substitution 2 ), and then ignore all interaction terms which are coupled with coefficients of order o(1) N →∞ . All this leads to the formal expression where Note that the expression (5) is formal since H N acts on the N -body Hilbert space H N while the Bogoliubov Hamiltonian H Bog acts on the excited Fock space Strictly speaking, for a * 0 a * 0 a0a0 we should rewrite it as (a * 0 a0) 2 − a * 0 a0 before doing the substitution where we have introduced the projections where the coefficients β p > 0 are determined by In fact, by using the CCR (3) it is straightforward to check that Consequently, where Note that the assumption 0 ≤ w ∈ ℓ 1 (2πZ d ) ensures that E Bog is finite. Moreover we have the uniform bounds Thus Bogoliubov's approximation predicts that the ground state energy of H N is In 2011, Seiringer [18] gave the first rigorous proof of (11). He also proved that the lowlying spectrum of H N is given approximately by the elementary excitation e p . These results have been extended to inhomogeneous trapped systems in [11], to more general interaction potentials in [12], to a large volume limit in [9], and to situations of multiple-condensation in [14,17].
Let us recall the approach in [12] which also provides the convergence of the ground state of the mean-field Hamiltonian H N in (1). Mathematically, the formal expression (5) can be made rigorous using the unitary operator introduced in [12] which is defined by Recall from [12,Proposition 4.2] that where N + is the number operator on the excited Fock space F + , Thus U N implements the c-number substitution in Bogoliubov's argument because it replaces a 0 by √ N − N + ≈ √ N (we have N + ≪ N due to the condensation). Then the formal expression (5) can be reformulated as which is rigorous since the operators on both sides act on the same excited Fock space. By justifying (14), the authors of [12] recovered the convergence of eigenvalues of H N first obtained in [18], and also obtained the convergence of eigenfunctions of H N to those of H Bog . In particular, for the ground state, we have from [ where |0 is the vacuum in Fock space. The convergence (15) holds strongly in norm of F + , and also strongly in the norm induced by the quadratic form of H Bog in F + . In particular, this implies the convergence of one-body density matrix in trace class (see (68) for a detailed explanation). Since Trγ Recall that P = |u 0 u 0 | = 1 − Q. The formula (17) looks similar to the result in Theorem 1, except that the cross term P γ Ψ N P is missing. Putting differently, to get the result in Theorem 1 we have to show that lim N →∞ Tr P γ (1) As explained in [12,Eq. (2.19)], from (16) and the Cauchy-Schwarz inequality one only obtains that the left side of (18) is of order O( √ N ). Moreover, (18) implies that thus answering an open question in [13]. As explained in [13,Section 5], (19) would follow if we could replace U N Ψ N by U B |0 (which is a quasi-free state, and thus satisfies Wick's Theorem [19,Chapter 10]). However, the norm convergence (15) is not strong enough to justify (19).

1.3.
Outline of the proof. To prove Theorem 1 we have to extract some information going beyond Bogoliubov's approximation. Roughly speaking, we will refine (14) by computing exactly the term of order O(N −1 ). Our proof consists of three main steps.
Step 1 (Excitation Hamiltonian). After implementing the c-number substitution, instead of ignoring all terms with coefficients of order o(1) N →∞ , we will keep all terms of order O(N −1 ). More precisely, in Lemma 7 below we show that in an appropriate sense, where The formula (20) is obtained by a direct computation using the actions of U N as in [12], plus an expansion of √ N − N + and (N − N + )(N − N + − 1) in the regime N + ≪ N . The advantage of using G N is that it is well-defined on the full Fock space F + . This idea has been used to study the norm approximation for the many-body quantum dynamics in [8].
Step 2 (Quadratic transformation). Then we conjugate the operator on the right side of (20) by the Bogoliubov transformation U B in (7). In Lemma 8 we prove that where and R 2 is an error term whose expectation against the ground state is of order O(N −3/2 ). Note that in C N we keep only cubic terms with three creation operators or three annihilation operators. These are the most problematic terms. All other cubic terms, as well as all quartic terms, are of lower order and can be estimated by the Cauchy-Schwarz inequality (the quartic terms always come with a factor N −1 instead of N −1/2 and this helps).
As we will see, the energy contribution of the cubic term C N is of order O(N −1 ). Thus (21) implies that which improves (11). Moreover, for the ground state we have which in turn implies the norm approximation (up to an appropriate choice of the phase factor for Ψ N ) and the following bound on the one-body density matrix Unfortunately the latter bound is still weaker than (18). Thus the desired result (18) cannot be obtained within Bogoliubov's theory.
Step 3 (Cubic transformation). To factor out the energy contribution of the cubic term C N in (21), we will use a cubic transformation. It is given by where From the assumption w ∈ ℓ 1 (2πZ d ) and the bounds (10) we have the summability Here we insert the cut-off 1 ≤N in the definition of U S to make sure that it does not change the particle number operator N + too much; see Lemma 5 for details. The choice of the cubic transformation above can be deduced on an abstract level. Consider an operator of the form A = A 0 + X where X stands for some perturbation. Then, in principle, we can remove X by conjugating A with e S provided that In our situation, A 0 = p =0 e(p)a * p a p and X = C N , allowing to find S explicitly in (25). In Lemma 9 we prove that with an error term R 3 whose expectation against the ground state is of order O(N −3/2 ). This allows us to obtain the following improvements of (22), (23) and (24).
Then the ground state energy of the Hamiltonian H N in (1) satisfies Consequently, we have the norm approximation (up to an appropriate choice of the phase factor for Ψ N ) As we will explain, Theorem 2 implies (18) and thus justifies Theorem 1. The idea of using cubic transformations has been developed to handle dilute Bose gases in [20,2,3,16,1], where the interaction potential has a much shorter range but the interaction strength is much larger in its range. In this case, the contribution of the cubic terms is much bigger, and Bogoliubov's approximation has to be modified appropriately to capture the short-range scattering effect. Results similar to (11) have been proved recently for the Gross-Pitaevskii limit [2] and for the thermodynamic limit [20,10]. It is unclear to us how to extend Theorem 1 to the dilute regime.
Our work shows that in the mean-field regime, in contrast to the dilute regime, the cubic terms are smaller, and they actually contribute only to the next order correction to Bogoliubov's approximation (there are also some quadratic and quartic terms which contribute to the same order of the cubic term). On the other hand, it is interesting that the contribution of the cubic terms is not visible in the expansion of the one-body density matrix in Theorem 1; putting differently the approximation in Theorem 1 can be guessed using only Bogoliubov's theory (although its proof requires more information).
There have been also remarkable works concerning higher order expansions in powers of N −1 in the mean-field regime; see [15] for a study of the ground state, [7] for the low-energy spectrum, and [6,5] for the quantum dynamics. These works are based on perturbative approaches which are very different from ours. Note that the method of Bossmann, Petrat and Seiringer in [7] also gives access to the higher order expansion of the reduced density matrices (see [7,Eq. (3.15)] for a comparison). We hope that our rather explicit strategy complements the previous analysis in [15,6,5,7] concerning the correction to Bogoliubov's theory in the mean-field regime.
Organization of the paper. In Section 2 we will derive some useful estimates for the particle number operator N + . Then we analyze the actions of the transformations U N , U B , U S in Sections 3, 4, 5, respectively. Finally, we prove Theorem 2 in Section 7 and conclude Theorem 1 in Section 7.

Moment estimates for the particle number operator
In this section we justify the Bose-Einstein condensation by showing that the ground state has a bounded number of excited particles. As explained in [18], the uniform bound on the expectation of N + follows easily from Onsager's inequality (2). For our purpose, we will need uniform bounds for higher moments of N + . The following lemma is an extension of [13, Lemma 5].
Proof. As in [13,Lemma 5], from the operator inequality we obtain for s = 1, 2, 3. Let us assume that s ∈ N is even. We will show that Ψ N , N s+1 + Ψ N ≤ C. Since Ψ N is a ground state of H N , it solves the Schrödinger equation Consequently, we get the identity The left side of (30) can be estimated using (28) and (29) as For the right side of (30), since using (4) and the CCR (3) we write Now we take the expectation against Ψ N and estimate. For the first term on the right side of (32), by the Cauchy-Schwarz inequality, we get for a given j Here we have used that a 0 a 0 commutes with Similarly, for the second term, we have as before. For the third term, we can bound and proceed similarly for other terms. Thus in summary, from (32) we get Inserting (31) and (33) into (30), we obtain By the Cauchy-Schwarz inequality we get We can now use and obtain Telescoping this inequality and using [13,Lemma 5] we arrive at a bound on Ψ N , N s+1 + Ψ N that is uniform in N . This gives the desired result for odd powers of N + . Finally, using (34), we obtain the bound for any s ∈ N and this ends the proof.
In order to put Lemma 3 in a good use, we will also need the fact that the moments of N + are essentially stable under the actions of the Bogoliubov transformation and the cubic transformation.
Lemma 4. Let U B be given in (7). Then Lemma 5. Let U S = e S be given in (25). Then for all t ∈ [−1, 1] and k ∈ N, The results in Lemma 4 and Lemma 5 are well-known. For the completeness, let us quickly explain the proof of Lemma 5, following the strategy in [2, Proposition 4.2] (the proof of Lemma 4 is similar and simpler).
Proof of Lemma 5. Take a normalized vector Φ ∈ F + and define Then Here we have used It is obvious that |Θ k (N + )| ≤ C k (N + + 1) k−1 . Combining with the summability (27) and the Cauchy-Schwarz inequality we obtain Thanks to the cut-off, we can bound Since the latter bound is uniform in Φ, we get the desired operator inequality.
We will also need the following refinement of Lemma 5.
Proof. Take a normalized vector Φ ∈ F + and define g(t) = Φ, e tS N k + e −tS Φ , ∀t ∈ [−1, 1]. Then proceeding similarly to (37), we have with f (t) being defined in the proof of Lemma 5. Using (38) and the Cauchy-Schwarz inequality we obtain From Grönwall's lemma, it follows that The latter bound is uniform in Φ and it implies the desired conclusion.

Excitation Hamiltonian
In this section, we study the action of the transformation U N in (12). By conjugating H N with U N , we can factor out the contribution of the condensation. More precisely, we have Lemma 7. We have the operator identity on F ≤N p+ℓ a * k−ℓ a p a k and the error term R 1 satisfies the quadratic form estimate Moreover, we have the operator inequality on F + Proof. A straightforward computation using the relations (12) shows that This operator identity holds on F ≤N + . For further analysis, we will expand √ N − N + and (N − N + )(N − N + − 1), making the effective expressions well-defined on the whole Fock space F + . This idea has been used before in [8]. Here it suffices to use and The operator inequalities (40) and (41) with G N given in the statement of Lemma 7 and with the error term By the Cauchy-Schwarz inequality, we have the quadratic form estimates This completes the first part of Lemma 7. Now let us turn to the operator inequality on the Fock space F + . We have proved that Let us compare the right side of (42) with the corresponding version without the cut-off 1 ≤N . First, consider the terms commuting with N + . Since this operator is not smaller than its product with the cut-off 1 ≤N . Moreover, using Finally, consider Moreover, since X changes the number of particles by at most 2, we have Hence, combining with the above bound on ±X we find that This completes the proof of the operator inequality on F + in Lemma 7.

Quadratic transformation
Recall that the Bogoliubov transformation U B in (7) diagonalizes H Bog as in (9). In this section, we will study the action of U N on the operator G N . We have Lemma 8. Let G N be given in Lemma 7. Then we have the operator identity on F + and the error term R 2 satisfies Proof. Let us decompose Non-cubic terms. Let us prove that (44) Let us start with the quadratic terms involving a * p a * −p . Using (8)

and the CCR (3) we have
Obviously the constant in (45) satisfies Moreover, the other terms in (45) can be rewritten as All of the sums in (46), (47), (48) are of the general form (43)-(44), thanks to the uniform bounds (10). The quadratic terms involving a * p a p can be treated similarly.
Next, consider It is straightforward to see that, except the constant all other terms in (49) can be written as in (43), with the corresponding bound (44) following from (10). By the same argument, we can show that the terms involving a * p a p N + , a * p+ℓ a * q−ℓ a p a q and N + (N + − 1) are of the general form (43)-(44). Next, let us bound the terms of the general form (43)-(44). We consider the case s ≥ t (the other case is treated similarly). By the Cauchy-Schwarz inequality ...,ms,n 1 ,...,nt   A m 1 ,...,ms,n 1 ,...,nt a * m 1 . . . a * ms a n 1 . . . a nt + h.c. . . a * n 1 (N + + 5) s−1 a n 1 . . . a nt for all ε > 0. Note that if min(t, s) ≥ 1, then on the right side of (50) we can replace (N + + 1) t+s−1 by N t+s−1 + . In particular, for the non-cubic term D N , using (50) with ε = N −1/2 and t + s ≤ 4 we get Cubic terms. By using (8) we have By using (10), we can write the last sum of (52) as p,q,r A p,q,r a * p a q a r with sup p q,r Using (50) with ε = 1, t = 1, s = 2, we get ± p,q,r A p,q,r a * p a q a r + h.c. ≤ Here N + ≤ N 2 + since the spectrum of N + is {0, 1, 2, ...}. The second sum on the right side of (52) can be treated by the same way. Thus from (52) and its adjoint, we have In particular, (54) implies that 0|U B C N U * B |0 = 0. Therefore, from (9), (51) and (54) we obtain the desired conclusion of Lemma 8.

Cubic transformation
To factor out the cubic term C N in Lemma 8, we will use a cubic renormalization. We will prove Lemma 9. Let C N be the cubic term in Lemma 8 and let U S be given in (25). Then we have the operator identity on Fock spacee F + Proof. Recall that from Lemma 8 we have Thanks to Lemma 5 and Lemma 6, we find that Thus this error term is part of R 3 . For the main term, we use U S = e S and the Duhamel formula Controlling C N + [S, dΓ(ξ)]. Since dΓ(ξ) commutes with N + and (e(p + q) + e(p) + e(q))η p,q a * p+q a * −p a * −q 1 ≤N + h.c.
which is equivalent to where 1 >N = 1 − 1 ≤N = 1(N + > N ). Thanks to the summability (10), we can use the Cauchy-Schwarz inequality similarly to (50) (with ε = 1) to get Combining with Lemma 5 and Lemma 6 we obtain By the Cauchy-Schwarz inequality as in (50), we get It remains to bound [S > , C N ]. From the explicit form of S > and C N , it is straightforward to check that On the other hand, we observe that and that C N does not change the number of particles more than 3. Therefore, Moreover, it is obvious that 0|[S > , C N ]|0 = 0 for N ≥ 10. Thus from (59) and (60) we obtain Combining with Lemma 5 we conclude that Conclusion. Inserting (58) and (61) in (57) we find that Combining with (55) we deduce that Taking the expectation of the latter bound again the vacuum, we find that Thus we obtain the desired conclusion This completes the proof of Lemma 9.

Proof of Theorem 2
Proof. We will prove the ground state energy estimate Upper bound. We use the following N -body trial state Then by the variational principle and the operator inequality on F + in Lemma 7 we have By Lemma 4 and Lemma 5 we know that Consequently, Combining with Lemma 9 we find that In the last estimate, we have also used the simple upper bound which will be justified below.
Lower bound. Let Ψ N be the ground state of H N and denote Φ := U S U B U N Ψ N ∈ F + . By Lemmas 3, 35 and 36, we have Then from the operator identity on F ≤N + in Lemma 7 it follows that Next, using Lemma 9 together with two simple estimates: From (63), since Φ, N + Φ ≥ 0 we obtain the desired energy lower bound ). This and the obvious upper bound E N ≤ w(0)(N/2) imply the simple estimate (62). Thus the matching energy upper bound is valid, and hence we conclude that Ground state estimates. By comparing the ground state energy expansion (64) with the lower bound (63) we deduce that Let us write Φ = (Φ j ) ∞ j=0 with Φ j ∈ H j + . We can choose a phase factor for Ψ N such that Φ 0 ≥ 0. Then Putting back the definition Φ = U S U B U N Ψ N we obtain the norm approximation This completes the proof of Theorem 2.

Proof of Theorem 1
Proof. Let Ψ N be the ground state for H N . As explained in the introduction, we will decompose γ (1) |0 strongly in the quadratic form of H Bog on F + . Moreover, it is easy to see that (see e.g. [13, Proof of Theorem 1]). Therefore, in the limit N → ∞,
Since {u p } p =0 is an orthonormal basis for H + , we have Using the excitation map U N and the relations (12) we can decompose Therefore, by the Cauchy-Schwarz inequality For the second sum in (70), using the Cauchy-Schwarz inequality, the simple bound and Lemma 3, we find that To control the first sum in (70), we will use the bound from Theorem 2: Also, from Lemma 3, Lemma 4 and Lemma 5 it follows that Φ, (N + + 1) 4 Φ ≤ C.
Using the action of Bogoliubov's transformation in (8) and the uniform bounds (10) we obtain The right side of (77) can be bounded using the Cauchy-Schwarz inequality, the summability Similarly, the terms involving I 2 (p) are bounded by Consequently, Inserting (79) and (76) in (75) and (74) we obtain Using the latter bound and (71), we deduce from (70) that This implies (69) and completes the proof of Theorem 1.