{'fraction_non_alphanumeric': 0.03833593293039579, 'fraction_numerical': 0.017742249951257553, 'mean_word_length': 4.045247756055576, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 4, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
On the Convergence of the Linear Delta Expansion for the Shift in T_c for Bose-Einstein Condensation

The leading correction from interactions to the transition temperature T_c for Bose-Einstein condensation can be obtained from a nonperturbative calculation in the critical O(N) scalar field theory in 3 dimensions with N=2. We show that the linear delta expansion can be applied to this problem in such a way that in the large-N limit it converges to the exact analytic result. If the principal of minimal sensitivity is used to optimize the convergence rate, the errors seem to decrease exponentially with the order in the delta expansion. For N=2, we calculate the shift in T_c to fourth order in delta. The results are consistent with slow convergence to the results of recent lattice Monte Carlo calculations.


I. INTRODUCTION
The critical temperature T c for Bose-Einstein condensation in an ideal gas of identical bosons was deduced in 1924 by Albert Einstein [1]: where n is the number density of the bosons. An obvious next question is the shift ∆T c in the critical temperature due to the interaction between the bosons. If the potential between two bosons is short-ranged, their scattering amplitude at sufficiently low energy is a constant proportional to the s-wave scattering length a. The Bose gas with such an interaction has been studied extensively since the pioneering work [2] of Lee and Yang in 1957 on bosons with a hard-sphere potential. However the shift in T c from such an interaction remained an unsolved problem until recently. A major breakthrough was made by Baym et al. [3], who showed that the problem could be reduced to a nonperturbative calculation in an effective 3-dimensional statistical field theory with O(2) symmetry. This demonstrated conclusively that the leading order shift was linear in a: ∆T c T c = c n 1/3 a, where c is a numerical constant. Using a self-consistent calculation of the propagator, they estimated the coefficient to be c = 2.9 [3]. Baym, Blaizot and Zinn-Justin calculated the coefficient analytically in the large-N limit and obtained c = 2.33 [4]. The first correction in the 1/N expansion was calculated by Arnold and Tomasik and it reduces c to 1.72 [5]. Definitive calculations of the coefficient for the physical case N = 2 were carried out by Kashurnikov, Prokof'ev, and Svistunov [6] and by Arnold and Moore [7] by applying Monte Carlo methods to a lattice formulation of the effective field theory with the result c = 1.32 ± 0.02. The second order correction to ∆T c /T c proportional to (an 1/3 ) 2 has since been calculated by Arnold, Moore and Tomasik [8]. Before the recent breakthrough, there had been many previous attempts to calculate the leading shift in T c [9,10,11,12]. The results for the coefficient c in (2) ranged from −1 to +4.7. The definitive solution of this problem provides a way to assess the reliability of the various calculational methods. One particularly simple method that can be applied to this problem is the linear δ expansion [13], also known as optimized perturbation theory [14] or variational perturbation theory [15]. With this method, an arbitrary parameter m is introduced into the theory and calculations are carried out using perturbation theory in a formal expansion parameter δ whose value is 1. The rate of convergence can be improved by adjusting the parameter m at each order in δ by some criterion such as the principal of minimal sensitivity. It was first applied to the problem of calculating ∆T c by de Souza Cruz, Pinto and Ramos [16]. At second, third, and fourth orders in δ, the predictions are c = 3.06, 2.45, and 1.48, respectively [17], which seem to be converging to the lattice Monte Carlo result.
In this paper, we revisit the application of the linear δ expansion to the calculation of the shift in T c . In Section II, we review the reduction of the calculation of ∆T c to a problem in a 3-dimensional statistical field theory. In section III, we introduce the linear δ expansion and review the results of Refs. [17]. We present two alternative precriptions for applying this method to the shift in T c . In Section IV, we show that in the large-N limit, these alternative prescriptions converge to the analytic result of Ref. [4] for some range of m. We study the improvement in the convergence rate that is obtained by using the principal of minimal sensitivity to choose a value of m that depends on the order in δ. We show that there is a family of solutions to this equation that seems to improve the convergence rate from a power law in the order of the δ expansion to exponential. In Section V, we calculate the result for N = 2 at 3 rd and 4 th orders in δ. In Section VI, we compare our results with the previous calculations using the linear δ expansion of Refs. [16,17]. Our results are consistent with slow convergence to the Monte Carlo results of Refs. [6,7].

II. EFFECTIVE THEORY
We consider a system of identical bosons of mass m at a temperature T and number density n close to the critical temperature T c for Bose-Einstein condensation. If the effects of interactions are small, Einstein's formula (1) implies that the interparticle spacing n −1/3 and the thermal wavelength λ T = (2πh 2 /mk B T ) 1/2 are comparable. We assume that they are both large compared to the range R of the two-body potential. In this case, the interaction between the bosons can be characterized entirely by the s-wave scattering length a, which we assume to have magnitude of order R. Thus we require a ≪ n −1/3 , λ T . This system can be described by a second-quantized Schrödinger equation. The corresponding imaginary-time Lagrangian is where µ is the chemical potential associated with the particle number density. In this paper, we generally use units such thath = k B = 1. Treating the system at finite temperature using the imaginary time formalism, the field ψ(x, t) can be decomposed into Fourier modes ψ n (x) with Matsubara frequencies ω n = 2πnT. The zero frequency mode ψ 0 (x) is defined by If the temperature is near the critical temperature T c , the chemical potential µ is small in magnitude compared to T . For momenta small compared to (mT ) 1/2 , the nonzero frequency modes have Euclidean energies of order T , while the zero frequency modes have Euclidean energies of order |µ|. Thus at large distances much greater than λ T , the nonzero modes decouple from the dynamics, leaving an effective theory of only the zero modes [3]. A naive derivation of the effective lagrangian is obtained simply by substituting ψ(x, t) → ψ 0 (x) in (3). The action then reduces to where the effective lagrangian is It is convenient to replace the complex field ψ 0 by the 2-component real field φ = (φ 1 , φ 2 ) defined by The effective lagrangian density for the dimensionally reduced theory can then be written as where the parameters are This O(2) field theory can be generalized to an O(N) field theory by replacing φ by an N-component real field φ = (φ 1 , ..., φ N ). This naive derivation of the effective lagrangian is supported by a more sophisticated analysis within the effective field theory framework [8]. This analysis reveals that there are also higher dimension operators in the effective lagrangian (8), such as ( φ 2 ) 3 and ∇ 2 φ · ∇ 2 φ, but they are not needed to calculate the shift in T c to leading order in an 1/3 . It also reveals that there are corrections to the coefficients r and u that can be calculated perturbatively as expansions in a/λ T . Again, they are not needed to calculate ∆T c to leading order in an 1/3 . Finally, since the scale T has been removed from the problem, the effective field theory defined by the effective lagrangian (8) requires an ultraviolet cutoff. Fortunately, the leading order shift in T c is dominated by the infrared region of momentum space and is insensitive to the ultraviolet cutoff.
From Einstein's formula (1), the leading order shift ∆T c in the critical temperature at fixed number density is related to the leading order shift ∆n c in the critical density at fixed temperature by the simple relation The number density in the original theory is n = ψ * ψ , where the angular brackets indicate the expectation value in the state corresponding to a homogeneous gas. The number density in the effective theory, to the accuracy required to calculate the shift in T c to leading order in an 1/3 , is where n 0 is the short-distance contribution to the number density, which depends on the ultraviolet cutoff. The angular brackets indicate the statistical average of the operator φ 2 . It can be expressed as where Σ(p) is the self-energy of the φ field. The integral over p is divergent and requires an ultraviolet cutoff. The number density for an ideal gas is given by the same expression with a = 0, which implies Σ(p) = 0: Again, the integral over p is divergent and requires an ultraviolet cutoff. To the accuracy required to calculate the shift in T c to leading order in an 1/3 , the constant n 0 is the same as in (12). The critical point is determined by the condition that the correlation length associated with the propagator in (13) is infinite: This condition reduces to r = 0 if a = 0. The change ∆n c in the critical density at fixed temperature T due to interactions is the difference between (12) with r = −Σ(0) and (14) with r = 0. The expression (11) for the fractional shift in T c can therefore be written as where ∆ is a quantity that can be calculated within the effective field theory: The integral over p is convergent and therefore no longer requires an ultraviolet cutoff. We have replaced the prefactor 2 in (17) by N to generalize it to the O(N) case. Comparing (16) with (2) and using (10) and (1), we find that the coefficient c is The effective theory defined by the lagrangian (8) is superrenormalizable. Although the theory has ultraviolet divergences, ∆ is dominated by infrared physics and it is insensitive to the ultraviolet cutoff. The relevant length scales are therefore those set by the parameters r and u of the effective field theory. The critical point can be reached by tuning the parameter r. At the critical point, there is only one relevant length scale and it is set by the parameter u. Since ∆ has dimensions of length, it must be proportional to u by simple dimensional analysis. The determination of the coefficient of u requires a nonperturbative calculation. Since u is linear in a as a → 0, this demonstrates that the leading term in the shift in T c is linear in a.
The calculation of the first order shift ∆T c in the critical temperature has been reduced to the calculation of the quantity ∆ in (17) for N = 2. The first calculation of this quantity within a systematically improvable framework was a calculation to leading order in the 1/N expansion by Baym, Blaizot, and Zinn-Justin [4]: The next-to-leading order (NLO) correction in the 1/N expansion was calculated by Arnold and Tomasik [5]. The first two terms in the 1/N expansion are For the physical case N = 2, ∆ has been calculated nonperturbatively by Kashurnikov, Prokof'ev, and Svistunov [6] and by Arnold and Moore [7] using Monte Carlo simulations of a lattice formulation of the O(2) field theory. The two results are consistent, but Ref. [7] quotes a smaller error: If we set N = 2 in the 1/N expansion in (20), the coefficient of −u/48π 2 is 1 at leading order and 0.74 at next-to-leading order. These coefficients seem to be approaching the lattice Monte Carlo result 0.57.

III. LINEAR δ EXPANSION
Any attempt to calculate ∆T c using ordinary perturbation theory in u is doomed to failure. In the symmetric phase where r > 0, perturbation theory is an expansion in powers of u/r 1/2 . If we naively let r → 0 to find the critical behavior, the perturbative corrections are infrared divergent, with the degree of divergence becoming more severe at each successive order. It is useful to think of perturbation theory as an expansion in powers of uξ, where ξ is the correlation length. Perturbation theory breaks down in the critical region, because uξ is of order 1 and every order of perturbation theory gives a leading order contribution.
The linear δ expansion (LDE) is a particularly simple method for obtaining nonperturbative results using perturbative techniques. Calculations are carried out as in ordinary perturbation theory by evaluating a finite number of Feynman diagrams. The name "linear δ expansion" comes from the fact that it can be defined by a lagrangian whose coefficients are linear in a formal expansion parameter δ. If L is the lagrangian for the system of interest, the lagrangian that generates the LDE has the form where L 0 is the lagrangian for an exactly solvable theory, The lagrangian L δ interpolates between L 0 when δ = 0 and L when δ = 1. The linear δ expansion is much easier to apply than the nonlinear δ expansion in which the lagrangian that interpolates between L 0 and L involves an operator that is a nonlinear function of δ [23]. The lagrangian L 0 for the exactly solvable theory typically involves an arbitrary parameter m. At δ = 1, the lagrangian (22) is independent of m. If a perturbation series in δ converges at δ = 1, it must converge to a value that is independent of m. However at any finite order in the δ expansion, results depend on m. As m varies over its physical range, the prediction for the observable often extends out to ±∞. Thus some prescription for m is required to obtain a definite prediction. A simple prescription for m is the principle of minimal sensitivity (PMS) that the derivative with respect to m should vanish. The LDE supplemented by the PMS prescription is often referred to as optimized perturbation theory [14] or variational perturbation theory [15]. A prescription such as PMS may improve the convergence rate of the LDE, because it allows the variable m to change with the order n of the LDE.
The convergence properties of the LDE have been studied extensively for the quantum mechanics problem of the anharmonic oscillator. The hamiltonian is The LDE is generated by the hamiltonian H δ = (1 − δ)H 0 + δ H, which interpolates between H 0 at δ = 0 and H at δ = 1. It can be expressed as The LDE has been proven to converge for appropriate order-dependent choices of the variable m that include the PMS criterion as a special case. Duncan and Jones showed that the finite temperature partition function converges exponentially, with the errors at n th order decreasing as exp(−bT n 2/3 /g 1/3 ) where b is a numerical constant [19]. Guida, Konishi, and Suzuki proved that the energy eigenvalues E n converge uniformly in g as n → ∞ [20]. Janke and Kleinert showed that the leading coefficient in the strong-coupling expansion for the ground state energy E 0 converges exponentially, with the errors decreasing like exp(−cn 1/3 ) where c is a numerical constant [21]. Bellet, Garcia, and Neveu have made numerical studies of the LDE optimized by the PMS criterion for the anharmonic oscillator [22]. The PMS criterion is a polynomial equation in the variable m, so its solutions are in general complex-valued. The authors of Ref. [22] found that the solutions at different orders n in the LDE can be organized into distinct families. The most accurate predictions for the ground state energy E 0 and for the quantity q 2 are obtained by using families of complex solutions. The resulting values of E 0 and q 2 are of course complex-valued, but the imaginary parts are small. The authors of Ref. [22] suggested that the real parts could be taken as the predictions of the LDE and the imaginary parts could be used as an indication of the error.
The anharmonic oscillator can be regarded as a Euclidean field theory with a single real-valued field in 1 space dimension. The statistical field theory defined by the lagrangian (8) is a generalization to a multicomponent field in 3 space dimensions. The most serious obstacles to generalizing the convergence proofs for the anharmonic oscillator to this more complicated problem come from the infrared and ultraviolet regions of momentum space. It is reasonable to expect the convergence behavior to be similar if appropriate infrared and ultraviolet cutoffs are imposed on the field theory. The quantity ∆ defined in (17) is insensitive to the ultraviolet region, so we do not expect any complications from taking the ultraviolet momentum cutoff to ∞. However, ∆ is very sensitive to the infrared region, so convergence in the presence of an infrared cutoff gives no information about the behavior of the LDE in the limit when the infrared momentum cutoff goes to 0.
We proceed to apply the LDE to the O(N) statistical field theory defined by the lagrangian L eff in (8). Our choice for the solvable field theory is the free field theory with mass m. The lagrangian L δ defined by setting L = L eff in (22) can be written Calculations are carried out by using δ as a formal expansion parameter, expanding to a given order in δ, and then setting δ = 1. Thus the Feynman diagrams generated by the LDE will involve a propagator with mass parameter m 2 , a 4-point vertex proportional to δ u, and a mass counterterm proportional to δ(r − m 2 ).
In order to apply the LDE to the shift in T c , we need a prescription for generalizing the quantity ∆ defined in (17) to the field theory whose lagrangian is L δ = L 0 + L int . The prescription must have a well-defined expansion in powers of δ, and it must reduce to (17) when δ = 1. One possibility, which we will refer to below as prescription II, is to use the expression (17), where Σ(p) is the self-energy for the field theory defined by L δ . Any such prescription defines ∆(u, m, δ) as a function of three variables. The n th order approximation in the LDE is then obtained by truncating the expansion in powers of δ at n th order to obtain a function ∆ (n) (u, m, δ) and then setting δ = 1: This sequence of approximations may not converge as n → ∞, but if it does converge, we expect it to converge to the correct value ∆. At any finite order n, the approximation (28) depends on m. We will see that as m varies over its physical range from 0 to +∞, the range of the predictions for ∆ extends out to +∞ for n even and out to −∞ for n odd. Thus some prescription for m is required to obtain a definite prediction. One prescription is the principle of minimal sensitivity (PMS): After setting δ = 1, ∆ (n) is a function of u and m only. By dimensional analysis, a solution m to the PMS criterion must be proportional to u and the value of ∆ (n) at such a point must also be proportional to u. An alternative prescription is the fastest apparent convergence (FAC) criterion: By allowing the variable m to change with the order n, prescriptions such as PMS or FAC may improve the convergence rate of the LDE. The LDE was first applied to the problem of calculating the shift in T c by de Souza Cruz, Pinto and Ramos [16]. Their prescription for the quantity ∆ defined in (17) is where Σ(p) is the self-energy for the field theory whose lagrangian is L δ . This prescription agrees with (17) at δ = 1. Note that in addition to the explicit dependence on δ in (31), there is also implicit dependence on δ through Σ(p), which can be calculated as a power series in δ.
The subscript I on ∆ in (31) distinguishes this prescription from alternative prescriptions for ∆ that we will label II and III. In the Feynman diagrams for evaluating Σ(p) in the LDE, all the propagators have mass m. Thus when the prescription I for ∆ is expanded in powers of δ, all the integrals involve only propagators of mass m, with the exception of the subtraction term in (31), which involves a massless propagator. The quantity ∆ was calculated to 2 nd order in δ in Ref. [16] and to 3 rd and 4 th order in Ref. [17]. After truncating the expansion in powers of δ and setting δ = 1, they determined m by applying the PMS prescription (29). At n th order in the δ expansion, this reduces to a polynomial equation in m with n roots, most of which are complex-valued. For a particular choice of the roots at 2 nd , 3 rd , and 4 th order in δ, the real parts of ∆/(−u/48π 2 ) are 1.31, 1.05, and 0.68, respectively.
This sequence seems to be approaching the lattice Monte Carlo result 0.57. In Ref. [17], the authors also applied the LDE to the second order shift in T c with encouraging results. Since the result for ∆ is known analytically in the large N limit, one would also like to compare the predictions of the LDE to this analytic result. Unfortunately, as we shall see, the prescription I of Ref. [17] is not well-behaved in this limit. We therefore present two alternative prescriptions for ∆ that are well-behaved in the large-N limit: These both agree with (17) at δ = 1. With prescription III, the LDE for ∆ generates integrals involving only propagators of mass m. With prescription II, the integrals involve both propagators of mass m and massless propagators. Our prescription III in (33) is closely related to the prescription I in Ref. [17]. The difference between them is given by a simple integral that is evaluated in (A.22): Because of the additional term in ∆ I proportional to N √ 1 − δ, the limit δ → 1 does not commute with the large-N limit. At n th order in the δ expansion, the terms that dominate as N → ∞ are the Nm term and the highest order term in the δ expansion which is proportional to N n u n /m n−1 . The PMS criterion gives a value of m that scales like N 1−1/n u and the corresponding value of ∆ scales like N 2−1/n u. The correct result is that ∆ scales like Nu in the large N limit. Thus prescription I has the wrong behavior in the large-N limit at any order in the LDE.

IV. LARGE N LIMIT
In this section, we study the behavior of the linear δ expansion (LDE) in the large-N limit defined by N → ∞, u → 0 with Nu fixed. We calculate ∆ to all orders in δ using prescriptions II and III. We show that the LDE converges to the correct analytic result for m/Nu greater than some minimum value, with errors that seem to decrease roughly like n −1/3 . We also study the improvement in the convergence obtained by using the PMS criterion for m. If we use the family of real solutions of the PMS criterion, the errors in ∆ seem to decrease like n −1/2 . If we use a particular family of complex solutions, the errors seem to decrease exponentially with n.
The quantity ∆ defined in (17) has the diagrammatic representation shown in Fig. 1, where the line with the shaded blob represents the propagator [p 2 + Σ(p) − Σ(0)] −1 . In the large-N limit, the leading contribution to the subtracted self-energy Σ(p) − Σ(0) comes from the series of bubble-chain diagrams whose 1 st and 4 th members are shown in Figs. 2 and 3. Since these diagrams are of order 1/N, the leading contribution at large N is obtained by expanding the expression (17) for ∆ to first order in Σ(p) − Σ(0). For the prescription II given in (32), the resulting expression is The corresponding expression for the prescription III is obtained by replacing p 2 → p 2 + m 2 (1 − δ) in the denominator of (35). Although Σ(p) is ultraviolet divergent, the difference Σ(p) − Σ(0) is finite. At any finite order in the LDE, Σ(p) − Σ(0) scales as p 2 at small p and as ln p for large p, so the integral in (35) is convergent. In addition to the diagrams for Σ(p) generated by the interaction term δ u φ 4 in (27), we must also take into account insertions of δ r = −δ Σ(0) and insertions of δ m 2 . The effect of the δ r insertions is to replace each 1PI self-energy subdiagram in Σ(p) by the subdiagram with the zero-momentum limit subtracted. The effect of the δ m 2 insertions is to replace m 2 by m 2 (1 − δ). We first consider the prescription II for ∆. The expression for the large-N diagram for ∆ with n + 1 loops is given in (A.23). The integral I n is a constant multiplied by m −(n−1) . After taking into account the δ m 2 insertions, the sum of the large-N diagrams is The integrals I n are expressed as 1-dimensional integrals in (B.23). Inserting these expres- sions into (36), we obtain where µ is the dimensionless mass parameter and A(y) is the function The prediction for ∆ at n th order in the δ expansion is obtained by expanding the expression (37) as a power series in δ, truncating after order δ n , and then setting δ = 1. It is easy to show that if the LDE for (37) converges, it converges to the correct analytic result (19). After interchanging the order of the sum and the integral, the sum can be evaluated and we obtain Upon taking the limit δ → 1, all dependence on the mass parameter µ disappears: After evaluating the integral over y, our final result is (19). Note that no prescription for µ was required to get this result. Thus, if the δ expansion converges at some value of µ, it should converge to (19). The manipulations that were used to show the convergence to (19) involved several interchanges of limits. First, we summed up some terms in the δ expansion to all orders to obtain the sum over diagrams in (36). Next, we interchanged the order of the sum over diagrams and the integral over the last momentum for each diagram. Next we took the limit δ → 1. Finally we evaluated the last momentum integral. It is difficult to translate the conditions for the validity of each of these steps into a statement about the convergence of the LDE. However, the convergence can be easily studied numerically. In Fig. 4, we show ∆ as a function of the dimensionless variational parameter µ for several orders in the δ expansion: n = 3 and n = 2 j + 3, j = 0, ..., 5. The horizontal line is the analytic result (19). The results are consistent with convergence to (19) for all µ greater than a critical value µ c which is close to 0.7. Our best estimate for the endpoint of the convergence region is µ c = 0.706 ± 0.003. If µ < µ c , ∆ seems to diverge to −∞ for n even and +∞ for n odd. For any fixed value of µ > µ c , the convergence with n is very slow. In column 1 of Table I,  we give the values of ∆ at µ = 1.039, which is the location of the minimum of the n = 3 curve in Fig. 4. We define the fractional error ε n by (Taking the real part of ∆ is superfluous for real values of µ, but it will be important later when we consider complex values of µ.) In Fig. 5, we show a log-log plot of ε n as a function of n. The points lie close to a straight line, indicating that the errors decrease like a power of n. The dotted line is ε n = 0.70n −0.37 , which goes through the points for n = 19 and 35. Thus the errors seems to decrease roughly like n −1/3 . The rate of convergence can be improved by changing the value of µ at each order in the LDE. One possibility is to use the PMS criterion (29) to determine µ. At n th order in the δ expansion, this condition reduces to a polynomial in µ of order n − 2. For n even, there are no real roots. For n odd, there is always one real root. These roots correspond to the minima of the curves in Fig. 4. The values of ∆/(−Nu/96π 2 ) are given in column 2 of Table I. There is slight improvement compared to column 1 in the convergence toward the correct value 1. The fractional errors defined in (42) are shown as a function of n in Fig. 5. The points lie close to a straight line, indicating that the errors decrease like a power of n. The dashed line is ε n = 0.78n −0.46 , which goes through the points for n = 19 and 35. Thus the errors seem to decrease roughly as n −1/2 .
The PMS criterion has either zero or one real solution, depending on whether n is even or odd. However there are always n−2 complex-valued solutions. In the case of the anharmonic oscillator, there are families of complex solutions with much better convergence properties than families consisting of purely real solutions [22]. In our problem, at any odd order n, the real solution always has the largest value of Re∆, which is the value farthest from the correct result (19). Thus this family of solutions gives the slowest possible convergence rate. However there is a strong anticorrelation between the errors in Re∆ and Im∆. This is illustrated in Fig. 6, which is a scatter plot of |Im∆| vs. Re∆ for the solutions to the PMS  Table I. Note that Im∆ for these solutions shows no sign of converging to 0, and in fact iseems to be converging to a value close to 0.4 Re∆. However the improvement in the convergence of Re∆ compared to the real values of ∆ in column 2 is evident. A log-log plot of the fractional errors in Re∆ is shown in Fig. 4. The downward curvature of the points indicates that the errors decrease faster than any power of n. The solid line is ε n = 0.32(0.93) n , which goes through the points for n = 19 and 35. Up to the order to which we have calculated, the errors are decreasing faster than this exponential. We obtain similar results if we use the prescription III for ∆ given in (33). The sum of the large-N diagrams is the same as in (36), except that the integrals I n are replaced by J n . The integrals J n are expressed as 1-dimensional integrals in (B.24). Interchanging the order of the sum and the integral and then evaluating the sum, we obtain Taking the limit δ → 1, all dependence on µ disappears and after evaluating the integral we again obtain (19). Thus, if the δ expansion converges, it should converge to the correct result. The convergence can be again studied numerically. We find that it seems to converge to (19) for all µ > µ c and it seems to diverge for µ < µ c , with µ c ≈ 0.7. In column 4 of  Table I, we give the values of ∆ for various orders n at µ = 0.947, which is the location of the minimum of ∆ for n = 3. The convergence at a fixed value of µ is very slow. As in the case of prescription II, the errors seem to decrease roughly as n −1/3 . More rapid convergence can again be obtained by using the PMS criterion for µ. The PMS criterion has no real solutions for n even and 1 real solution for n odd. The values of ∆ for the real solution at various orders n are given in column 5 of Table I. This is the family of solutions with the slowest convergence rate for Re∆. The errors seem to decrease roughly as n −1/2 . A family of solutions with a nearly optimal convergence rate for Re∆ consists of those with the maximal value of |Im∆|. The values of ∆ for this family are shown in column 6 of Table  I. The errors in Re∆ are numerically larger than those with prescription II, but they also seem to decrease exponentially with n.

V. FINITE N
In this section, we use the linear δ expansion (LDE) to calculate ∆ for finite N. We calculate ∆ to 4 th order in δ using both of the prescriptions II and III. Setting N = 2 and using the PMS criterion, we obtain the predictions relevant to the shifts in T c for Bose-Einstein condensation.
The first term in the LDE for ∆ is second order in δ. It is obtained by calculating the 3-loop diagram involving the 2-loop self-energy diagram in Fig. 2. The contribution to ∆ of order δ 3 is obtained by calculating the 4-loop diagram involving the 3-loop propagator correction in Fig. 7. The contribution of order δ 4 is obtained by calculating the 5-loop diagram involving the 4-loop propagator corrections in Fig. 8 Inserting the values for the integrals from Appendix B, truncating at successively higher orders in δ and then setting δ = 1, we obtain the results at 2 nd , 3 rd and 4 th orders in the LDE: where µ is the dimensionless mass parameter In the large-N limit, this reduces to the mass parameter introduced in (38). The accuracy of the 4 th order result in (47) is limited by the accuracy to which we have evaluated the integrals I 4b and I 4c . In Fig. 9, At order δ 2 , the PMS criterion has no solutions, because (45) is a monotonic function of µ. At order δ 3 , the PMS criterion has a single real solution µ (3) = 1.03886(N + 8)/(N + 2) and the resulting value for ∆ is At 4 th order in δ, ∆ is a monotonic function for real-valued µ, so the PMS criterion has no real solutions for µ. It does however have a conjugate pair of complex solutions whose values at N = 2 are µ = 2.0772 ± 1.773i. The corresponding values of ∆ are complex, but we can take the real part as the prediction for ∆ at 4 th order in the LDE: For N = 2, the 3 rd and 4 th order approximations in (49) and (50) differ from the lattice Monte Carlo result of (21) by about 66% and 62%, respectively. These predictions are not as accurate as in the large N limit, where the errors in the 3 rd and 4 th order approximations are about 52% and 44%, respectively. Furthermore, the percentage decrease in the error in going from 3 rd order to 4 th order is only about half as large for N = 2 as it is for large N.
With the prescription III, the expression for ∆ to 4 th order in δ is the same as (44), except that the integrals I n are replaced by J n . Inserting the value for the integrals J n from Appendix B, the results at 2 nd , 3 rd and 4 th orders in the LDE are ∆ (2) The accuracy of the 4 th order result in (53) is limited by the accuracy to which we have evaluated the integrals J 4b and J 4c . At order δ 2 , the PMS criterion has no solution. At order δ 3 , it has one real solution for µ. At 4 th order in δ, it has a conjugate pair of complex solutions. The corresponding values of ∆ are given in Table II. The values of ∆ at 3 rd order In Fig. 9, we show ∆ (n) as a function of µ for n = 2, 3, and 4. The values of ∆ (n) at the solutions to the PMS criterion are given in column 1 of Table II. At orders 3 and 4, the PMS criterion has a conjugate pair of complex solutions. As the predictions of the LDE, the authors of Ref. [17] chose the positive value for ∆ (2) and the real parts of the complex values for ∆ (3) and ∆ (4) . These predictions differ from the lattice Monte Carlo results by about 131%, 85%, and 20%, respectively.

VI. CONCLUSIONS
The calculation of the shift in T c to leading order in an 1/3 can be reduced to a nonperturbative calculation of a quantity ∆ in a 3-dimensional O(2) field theory at its critical point. This quantity can be generalized to an O(N) field theory and it can be calculated systematically using the 1/N expansion. We have shown that the linear δ expansion (LDE) can be applied to the calculation of ∆ in such a way that the large N and δ → 1 limits commute. The prescriptions II and III defined in (32) and (33) have this property, but the prescription I introduced in Ref. [16] does not.
For prescription II and III, we calculated ∆ in the large-N limit to all orders in the LDE, and studied its convergence properties. We found that it seems to converge to the correct analytic result for µ > µ c , where µ c ≈ 0.7. The convergence for fixed µ is extremely slow: the errors seem to decrease roughly like n −1/3 . The convergence can be accelerated by using the PMS criterion to choose a value of µ that depends on the order in δ. The PMS criterion is a polynomial equation with multiple roots, so there are various families of predictions for ∆. The family consisting of the real solutions of the PMS criterion gives the slowest possible convergence: the errors in ∆ seem to decrease roughly like n −1/2 . The family consisting of the complex solutions of the PMS criterion with the maximal value of |Im∆| have nearly optimal convergence: the errors in Re∆ seem to decrease exponentially with n. Unfortunately, accurate predictions still require calculations to high orders in the δ expansion. To achieve 10% accuracy requires computing to about 18 th order in δ.
For the case N = 2 relevant to Bose-Einstein condensation, we calculated ∆ to 4 th order in the LDE using both prescriptions II and III. The quantity ∆ was previously calculated to 4 th order using prescription I in Ref. [17]. The PMS criterion has complex solutions beginning at 3 rd order for prescription I and at 4 th order for prescriptions II and III. We take the values of Re∆ (n) for these solutions to be the predictions for ∆. The 3 rd and 4 th order predictions from prescription II are smaller than the lattice Monte Carlo result by about 66% and 62%, respectively. The 2 nd , 3 rd and 4 th order predictions from prescription I are larger than the lattice Monte Carlo result by about 131%, 85%, and 20%, respectively. We argue that the 20% error in the 4 th order prediction from prescription I is fortuitously small. Prescription I differs from prescription III by the additive term given in (34). The expansion of the factor √ 1 − δ truncated at n th order in δ and evaluated at δ = 1 converges to 0 as n → ∞. The prediction from prescription I should therefore approach the prediction from prescription III as n increases, which suggests that it will converge to the correct answer from below. Because the n = 2 prediction from prescription I is greater than the correct result by about a factor of 2.3, there must be a crossover in n where the prediction goes from values greater than the correct result to smaller values. The high accuracy of the 4 th order result suggests that this crossover occurs at n = 5. Thus we expect higher order predictions from prescription I to first overshoot the correct value and then approach it from below.
Before the lattice Monte Carlo calculation of Ref. [6,7], there were many previous attempts to calculate the coefficient c in the expression (2) for the leading order shift in T c for Bose-Einstein condensation. For most of these methods [3,9,10,11,12], there is no known systematic way of improving the calculation to make it more accurate. Such a method is best regarded as just a model of the critical O(N) field theory whose accuracy can only be determined by comparison with a more reliable calculational method. A systematically improvable method has two great advantages over such a model. First, the accuracy of the prediction can be improved with additional effort. Second, one can make a reliable estimate of the error in the prediction. One example of a systematically improvable method is the 1/N expansion [3,5], which has been used to calculate ∆T c to leading order [3] and to next-to-leading order [5]. Comparison with the lattice Monte Carlo results show that these first two approximations have errors of 77% and 30%, respectively. Our demonstration that the LDE converges in the large N limit suggests that the LDE may also be a systematically improvable method for N = 2.
Although the convergence rate of the optimized LDE appears to be exponential in the large-N limit, it is still rather slow. One must calculate to about 18 th order in δ to achieve 10% accuracy. For the case N = 2 relevant to Bose-Einstein condensation, we have calculated ∆ to 4 th order in the LDE. It may be feasible to extend the calculation to 5 th order, but it would be very difficult to extend it to much higher order. The slow convergence in the large-N limit suggests that even if the LDE also converges for N = 2, a strict expansion in δ is not useful for quantitative calculations. It may however be possible to use order-dependent mappings to change the expansion in δ into a more rapidly converging expansion [24].
In conclusion, we have shown that the LDE for the quantity ∆ that determines the leading order shift in T c converges in the large N limit if we use an appropriate prescription. If the PMS criterion is used to optimize the convergence, the errors seem to decrease exponentially in the order of the LDE. This suggests that the LDE may also be a systematically improvable approximation scheme in the case N = 2 relevant to Bose-Einstein condensation. With the use of order-dependent mappings, it may be possible to to develop the LDE into a general and powerful tool for quantitative calculations in superrenormalizable field theories, even at a critical point.
where q denotes the 3-dimensional integral over the momentum q : There is a simple analytic expression for the function D 2 (p): If one uses the prescription III for ∆ given in (32), the corresponding expressions for the diagrams are obtained by replacing each factor of 1/p 2 in the integrand by D 1 (p). We denote the integral corresponding to I n by J n . If one uses the prescription I of Ref. [17], there is also a diagram that is 0 th order in δ: The remaining integrals must be evaluated numerically. Three of the 5-loop integrals can be reduced to integrals over a single momentum. The integral I 4e in (A.20) is already in such a form. The integrals I 4a in (A. 16) and I 4d in (A. 19) can be reduced to this form by first integrating over p : (B.4) The resulting expression for the integrals are (B.6) Reducing these integrals to 1-dimensional integrals and evaluating then numerically, we obtain The integrals J n required if ∆ is calculated using the prescription III differ from the integrals I n calculated above only in the substitution 1/p 2 → D 1 (p). The 3-loop integral can be evaluated analytically:  The integrals I n that survive in the large-N limit can all be reduced to integrals over a single momentum. After interchanging orders of integration in (A.24) and integrating over p, the remaining integral reduces to