ARS MATHEMATICA CONTEMPORANEA Also available at http://amc-journal.eu ISSN 1855-3966 (printed edn.), ISSN 1855-3974 (electronic edn.) ARS MATHEMATICA CONTEMPORANEA 12 (2017) 67-76 Euler's divergent series and an elementary model in Statistical Physics Bill Allombert Institut de Mathématiques de Bordeaux, UMR 5251, Université de Bordeaux, 351 cours de la Liberation, F-33405 Talence Cedex, France Jean-Paul Allouche CNRS, Institut de Mathématiques de Jussieu-PRG, Université Pierre et Marie Curie, Case 247, 4 Place Jussieu, F-75252 Paris Cedex 05, France Michel Mendes France Institut de Mathematiques de Bordeaux UMR 5251, Université de Bordeaux, 351 cours de la Liberation, F-33405 Talence Cedex, France Received 3 April 2015, accepted 3 April 2016, published online 2 May 2016 We discuss the multiple integral of a multivariate exponential taken with respect either to the Lebesgue measure or to the discrete uniform Bernoulli measure. In the first case the integral is linked to Euler's everywhere divergent power series and its generalizations, while in the second case the integral is linked to a one-dimensional model of spin systems as encountered in physics. Keywords: Euler divergent series, Abel-Plana Formula, Stirling numbers, spin system, Ising chain. Math. Subj. Class.: 33E20, 28A35, 82B44, 82D30 1 Introduction Consider the integral (N > 1) E-mail addresses: biU.aUombert@math.u-bordeaux.fr (Bill Allombert), jean-paul.allouche@imj-prg.fr (Jean-Paul Allouche), michel.mendes-france@math.u-bordeaux.fr (Michel Mendes France) ©® This work is licensed under http://creativecommons.org/licenses/by/3.0/ Abstract (1.1) 68 Ars Math. Contemp. 12 (2017) 37-50 If ^ is the Lebesgue measure on (0, to)n and H = 1, the integral is linked to the series ]T( —1)"(n!)N-1x" n>0 which, for N = 2, is attributed to Euler. If N =1 the series reduces to (1 + x) 1 (convergent for |x| < 1) and for N > 2 it diverges for all x = 0. If on the other hand, ^ is the Bernoulli measure on the set {-1,1}N then the integral reads zn (x)=2n ^ exp i-h i- x n «£{±1}N \ V=1 / j = 1 and could represent a certain spin system described in Section 6. 2 Euler's divergent series If ^ is the Lebesgue measure on (0, to)n, we suppose that H > 0 and x > 0. There is no loss of generality in the choice H =1 in Formula 1.1; take new variables vj = Huj. Integrate ' N \ N \ N Zn(x) = / exp i — | Uj | — x J^J Uj | JJ^ duj J(0,)N \ \ I I with respect to duN to obtain Z ( ) i exp l—pN-1 Ujjj -1 d ZN(x) = W-i 1 + x Uj j0dUj * Suppose N > 2 since the case N = 1 is trivial. ZN converges for all complex x outside tN-1 W Uj the negative real axis (—to, 0). Expand + x nN=11 Uj) into a formal power series / n-1 \ N-1 N-1 Zn(x) = exp ( — £ Uj ) £ ( —1)"x" ^ Ujn n duj • A0,~)N 1 \ j=1 / n>0 j=1 j=1 If we accept to permute the summation with the integral, then N-1 „ Zn(x) = ^ ( —1)"x" n u" exp(—Uj)duj = ^ ( —1)"(n!)N-1x". n>0 j=1 A0,^ n>0 What happens if N =1 or 2? The case N =1 is trivial yet interesting, 1 Z1(x) = / exp(—u — xu)du = 0 1 + x Expanding the integral with respect to x we obtain Z1 (x) = ^n>0 ( —1)"x" and for x =1 we rediscover the well-known "equality" J2n>0 ( — 1)" = 1 • B. Allombert et al.: Euler's divergent series and an elementary model in Statistical Physics 69 exp(-u) du = 0.5963 ... and therefore concluded 0 1 + u ^ ( —1)n n! = 0.5963 ... n>0 a most astonishing equality! In his beautiful book G. H. Hardy [2] discusses in detail this case N = 2. Remark 2.1. The constant du = 0.5963... is called the Euler or the Euler- Gompertz constant (see [3], [1, Section 6.2], and in particular [1, Section 6.2.4] for the name "Gompertz"). Among the numerous results related to this constant we do not resist to write the following continued fraction expansion: exp(-u) 1 + u du 1 2 - 4 - 6 8 - '.. This continued fraction expansion is sometimes attributed to Stieltjes, but in [8] Stieltjes indicated that it was studied by Laguerre. We found indeed in [5, p. 154] that Laguerre considered e times the Prym function1 eQ(a) = ff" e1-xxa-1dx and obtained as consecutive approximations of eQ(0) the sequence 4 20 124 920 7940 7' 34 209' 2546 13327' which are exactly the values of the first few truncatures of the above continued fraction (also see Laguerre [4, p. 77]). Of course eQ(0) = fc ™ exp(-u) 1+u du = Z2(1): it would thus be interesting to obtain such nice continued fraction expansions for the quantities ZN(1). More generally, a formula given by Tannery in [9, p. 1699] or an easy rewriting of a formula given by Laguerre in [4, end of Page 75] reads -dt 1 x + 1 - x + 3 — x+5- But Z2(x) exp(—u) 1 r™ --u -du = — 1 + xu x + 7 — .. e u i r™ e-t e du = 1 e1/M —dt. x J o x + u h/a DO CXJ 1 0 4 9 x e t 1 4 9 oo t x 0 1Note that there seems to be a misprint in the formula given by Laguerre, where e1 x is replaced by e see the original definition by Prym [7, p. 169]. 70 Ars Math. Contemp. 12 (2017) 37-50 Hence Z2(x) 1 + x - 2 x2 4x 1 + 3x-- 9x2 1 + 5x-- 1 + 7x - '.. 3 The Borel operator The sequence can be defined recursively by means of the so-called Borel operator B : f ^ exp(-u)f(ux)du Jo -m)/(ux)du. o The Borel operator applies the series J2n>o f (n)(0) n onto J2n>0 f (n)(0)xn. Using the relation Z0 (x) = exp(-x) and ZN+1 = BZN, we see that the integral ZN is therefore the Nth iterate BN of x ^ exp(-x), or equivalently the (N - 1)st iterate BN-1 of x ^ (1 + x)-1. 4 The Abel-Plana summation and the r function In this section we study the behavior of ZN when N goes to infinity. Note that for real x > 0, the sequence N ^ ZN (x) is bounded from above by 1 and furthermore it is increasing. Indeed let AN(x) = ZN +1(x) - ZN (x) and nN(x) = x nN=1 M. Then An(x) = i exP I - ¿ u I ( * ( , - exP( nN (x))) n . J(0^)N \ J V1 + nN(x) Since ^ - exp(-t) > 0, ZN+1(x) > ZN(x) as claimed. Therefore ZN(x) tends to a limit which we now compute. Theorem 4.1. For all real x > 0, we have lim ZN(x) = 1. Proof. Since the result is trivial for x = 0, we may assume x > 0. We note that ZN (x) can be written as a diverging series Zn(x) = £ (-1)n/N(n) n>0 where fN : s ^ r(1 + s)N-1xs is an analytic function on the half-plane K(s) > -1. By blindly applying the Abel-Plana Formula (see [6, III, formula X]) to this series, we get f-1/2+iTO r(1 + z)N-1 xz ZN(x) = - -77—,-;- J-\/2-iw 2« sin(nz) _ fr(1/2 + it)N-1x-1/2+_^dt 2 cosh(nt) 1 2 B. Allombert et al.: Euler's divergent series and an elementary model in Statistical Physics 71 or by displacing the integration contour, I/2+íto p/1 + -1Xz Zn(x) = 1 - / + Z\ , dz (4.1) Ji/2-icx, 2, sm(nz) = 1 - f r(3/2 + ,t)N-V/2+" d( (4.2) ' 2 cosh(nt) v y ' —oo The convergence of the integrals is provided by the fact the r function decreases like exp(-72 |z|) as z goes to — l/2±iro (resp. 1/2±i), and sin(nz) increases like exp(n|z|). Strictly speaking, the Abel-Plana Theorem only applies for N = 0. However, by applying the Borel operator to the right-hand side and interverting the summations by Fubini's Theorem, we find that ír(i/2 + it)N-1x-1/2+¿í , Fn(x) = —' \ , ,-dt N ( ) J-n 2 cosh(nt) satisfies the same recursion as ZN(x). Indeed, . . r , ^ r(l/2 + it)N-1(lu)-1/2+it , , BFn (x) = exp(-u) ——----dtdu N ( ) Jo P( -J-to 2cosh(nt) rr(i/2 + it)N-1x-1/2+« r~ 1/2+it . w , = —-- u-1/2+it exp(-u)dudt. J-^ 2 cosh(nt) ./o From the identity r(1/2 + it) = J'0° u-1/2+lt exp(-u)du, it follows fr(1 /2 + it)Nx 1/2+it BFn (x) = i(1/ o+ w ^-dt = Fn+1(x) ./-rc 2 cosh(nt) therefore FN = ZN. Now since |r(3/2 + it) | < ^ < 1 for all t G R, the integral (4.2) converges to 0 when N goes to infinity for all real x > 0, thus we have proved: lim ZN(x) = 1. □ To conclude this section, we note that this formula for ZN involves a single integral which is much more suitable for numerical computations than the original formula involving a multiple integral. Note also that N need not be an integer. . . 5 A differential equation It might be worthwhile to mention that the function N \ N \ N ,JN (X) = / 72 Ars Math. Contemp. 12 (2017) 111-126 is a solution of a differential equation of order N - 1 with polynomial coefficients. Indeed, the shortest way to establish this is to introduce the linear operator U defined by U(z) = (xz)'. Clearly U(xn) = (n + 1)xn so that Uk(xn) = (n + 1)kxn. Then rN-1 ry /„,\ _ TTN- 1 V"^ / i\n/ i\N-In UN-1ZN(x) = UN(-1)n(n!)N-1x n>0 = £ (-1)n(n!)N-1(n +1)N- n>0 = £ ( —1)"((n +1)!)N-1xn ; n>0 xUN-1Zn(x) = £ ( —1)n((n +1)!)n-1xn+1 n>0 = 1 - Zn(x). The function ZN (x) is thus solution of the (N - 1)-st order differential equation xUN-1y + y = 1 with initial conditions y(0) = 1, y'(0) = -1,..., y(N-2)(0) = (-1)n-2 ((N - 2)!)N-1. The reader may well criticize the above proof since it involves divergent series. There is however no problem in justifying the result by applying the operator U to the integral representation of ZN (x); the calculations are just slightly more cumbersome. Example 5.1. Z2 (x), Z3 (x), Z4 (x) are respectively solution of the equations x2y' + (x + 1)y = 1 x3y" + 3x2y' + (x + 1)y = 1 y'" + 6x3y" + 7x2y' + (x + 1)y = 1 x y +oi y + 7x y + (x The reader will recognize the numbers above as the Stirling numbers of the second kind. This can be proved by noting that both families of numbers obey the formula «n+1,k = k&n,k + «n,k-1. 6 An unconventional spin system We now assume that p is the Bernoulli measure on {-1, +1}N: zn (x)=2n e exp i -h i ^ i- x n «£{±1}^ \ V=1 / j = 1 We interpret ZN as the partition function of a certain spin system which we describe below. Conventional spin systems are discussed for example in C. J. Thompson [10]. Imagine an N-component particle, each component of which has a spin uj = ±1, and which are instantaneously influenced by the N - 1 others. The "total" spin of the particle, B. Allombert et al.: Euler's divergent series and an elementary model in Statistical Physics 73 i.e., its sign is j uj. A real external field H acts on the spins. The Hamiltonian attached to the spin system in state u = (u^ u2,..., uN) with external field -H is then given by N N ^n uj+HJ2uj. j=1 i=1 The behavior of the spin system is controlled by the partition function, in particular by its thermodynamical limit log Zn (x) lim -—-- n^TO N Theorem 6.1. For all real x > 0, ZN(x) = cosh(x) cosh(H)N - (-1)N sinh(x) sinh(H)N. Proof. By using the relation exp(-1) = cosh(t) — sinh(t), we write zn(x)=2n ^ exp ( -h (^ uj illcosh(x n uj) - sinh(x n uj)). ue{±i}^ V \j=1 J J \ j=1 j=1 J Since f]f=1 uj = ±1, cosh is even and sinh is odd, it follows that 1 ( ( N W( N \ Zn (x) = 2N exp l -H l ^ uj l l l cosh(x) - sinh(x) JJ uj l . (6.1) m£{±1}^ \ \j=1 // V j=1 / The following two formulas are easily proved by recursion on N: £ exp (h uj )) = (2cosh(H))N ue{±1}( V \j=1 JJ £ exp lH uj im uj = (2sinh(H))N. ue{±1}N V \j=W / j=1 From Equation (6.1) it follows: Zn(x) = cosh(x) cosh(H)N - (-1)N sinh(x) sinh(H)N. □ Remark 6.2. Theorem 6.1 above implies that Zn(x) ~ cosh(H)N cosh(x), so that N ^-TO log Zn (x) lim -—-= logcosh(H) n^TO N z K J which happens to be independent of x and which is continuous with respect to H. The system has no critical value of the external field and therefore presents no phase transition. 76 Ars Math. Contemp. 12 (2017) 111-126 7 A disturbed Ising chain In the preceding section we described an unconventional spin system. We now turn to the most familiar one, namely the one dimensional Ising chain (see [10]) with Hamiltonian N N j=i j=i where J is a "coupling constant". Actually this Hamiltonian corresponds to the parameters —H and - J but that makes no essential difference for our computation. We consider in fact a perturbed Ising chain with the additional term x\\\N=i uj. The Hamiltonian is therefore N N N H(u) = H ^^ Uj + Ujuj+i + x J^J j , uj+i + x|juj j=1 j=1 j=1 and the partition function is now Yn = 22N ^ exp(—H(u)) «e{±i}N which we propose to compute where we need to specify uN +1. Following most textbooks, we simplify the model by assuming that the chain is cyclic: un +i — ui. Theorem 7.1. Define A± — exp( —J) cosh(H) ± (exp(—2 J) cosh(H)2 + 2 sinh(2 J))2 , A± — exp( —J) sinh(H) ± (exp(—2 J) sinh(H)2 — 2 sinh(2 J)) 2 . Then Yn = 2N cosh(x)(AN + AN) — -tN^ sinh(x)(A\ + A-). Proof. Observe as in Section 6 that cosh x , sinh x „ yn = 2N yn 2N" yn where (N N \ — H ^ uj — J ^ ujuj+i ), j=i j=i / (N N \ N —Hj2uj—Jj2uj uj+i) n uj. j=i j=i j j=i The classical way to compute Y\ is to introduce the 2 x 2 transfer matrix L = i ¿i(l, 1) ¿i(l, —1) Li = 1 Li( —1,1) Li( —1, —1) B. Allombert et al.: Euler's divergent series and an elementary model in Statistical Physics 75 where Li(ui, «2) = exp ^ - §(« + «2) - Juiu^j . In other words j _ ( exp(-H — J) exp(J) Ll = ^ exp( J) exp(H — J) Then Yn = Li(MI,M2)LI(M2,U3) . ..LI(Un ,ui) «e{±i}N = ^ Lf (ui,ui) = Trace(Lf ) = Af + Aw «i£{±1} where A+ and A_ are the eigenvalues of L1, i.e., the solutions of A2 - 2Aexp(-J) cosh(H) + exp(-2J) - exp(2J) = 0. Therefore A± = exp(-J)cosh(H) ± (exp(-2J) cosh(H)2 + 2sinh(2J))2 The computation of Y^ is quite similar. Let L2(1, 1) L2(1, -1) where so that then ¿2 = l L2(-1,1) L2(-1,-1) H ¿2(^1,^2) = "i exp ( --^(ui + "2) - Jmi«2 L I exp(-H - J) exp( J) 2 - exp(J) - exp(H - J) Yx = L2(ui,M2)L2(u2,M3) . ..L2(u f,ui) «E{±1}W = ^ Lf (ui,ui) = Trace (L f )= Af + A- «i£{±i} where A+ and A_ are the eigenvalues of L2, i.e., the solutions of A2 + 2Aexp(-J) sinh(H) - exp(-2J) + exp(2J) = 0. Therefore A± = - exp(-J)sinh(H) ± (exp(-2J) sinh(H)2 - 2sinh(2J))2 . Finally Yf = cosh(x)(A+ + A f ) - sinh(x)(A+ + Af ). □ 68 Ars Math. Contemp. 12 (2017) 111-126 Remark 7.2. The reader will easily verify that for J = 0 we obtain the value of ZN computed in Section 6. Remark 7.3. It is not difficult to see that max{|A_|, |A_|, |A+|} < A+. Hence YN ~ 2N cosh(x)(A+) when N goes to infinity. This implies that the following limit exists, is continuous in both variables J and H, and is independent of x (as in Remark 6.2); therefore the system has no phase transition: y log YN A+ lim -zrz.- = log- N 2 = log f(eXp(—J)COSh(H) + (exp(-2J)cosh(H>2 + 2(sinM2J))^ . 8 Conclusion and acknowledgements This article illustrates a classical fact, namely that one formula may well lead to far distant and unexpected developments. Unifying themes is probably one of the most exciting aspects of mathematics. We warmly thank H. Cohen and A. Lasjaunias for their very kind help. We are very grateful to J.-Y. Yao and to the referees for their precise and constructive remarks and comments. J.-P. A. was partially supported by the ANR project ANR-12-IS01-0002 "FAN" (Fractals et Numeration). References [1] S. R. Finch, Mathematical constants, volume 94 of Encyclopedia of Mathematics and its Applications, Cambridge University Press, Cambridge, 2003, doi:10.1017/CBO9780511550447, http://dx.doi.org/10.1017/CBO97 805115504 47. [2] G. H. Hardy, Divergent Series, Oxford, at the Clarendon Press, 1949. [3] J. C. Lagarias, Euler's constant: Euler's work and modern developments, Bull. Amer. Math. Soc. (N.S.) 50 (2013), 527-628, doi:10.1090/S0273-0979-2013-01423-X, http://dx.doi. org/10.1090/S027 3-097 9-2013-01423-X. [4] Laguerre, Sur l'integrale f° ^-x^, Bull. Soc. Math. France 7 (1879), 72-81, http://www. numdam.org/item?id=BSMF_18 7 9_7_72_1. [5] E. Laguerre, Sur la reduction en fractions continues d'une fraction qui satisfait a une equation differentielle lineaire du premier ordre dont les coefficients sont rationnels, J. Math. Pures Appl. 1 (1885), 135-166, http://eudml.org/doc/2 34 408. [6] E. Lindelof, Le calcul des résidus et ses applications a la théorie des fonctions, Gauthier-Villars, Imprimeur-Libraire, 1905, http://www.gutenberg.org/ebooks/2 97 81. [7] F. E. Prym, Zur Theorie der Gammafunction, J. Reine Angew. Math. 82 (1877), 165-172, doi: 10.1515/crll.1877.82.165, http://dx.doi.org/10.1515/crll.1877.82.165. [8] T.-J. Stieltjes, Recherches sur les fractions continues, Ann. Fac. Sci. Toulouse Sci. Math. Sci. Phys. 8 (1894), J1-J122, http://www.numdam.org/item?id=AFST_18 9 4_1_8_4_ J1_0 . [9] J. Tannery, Sur les integrales euleriennes., C. R. Acad. Sci., Paris 94 (1882), 1698-1701. [10] C. J. Thompson, Mathematical statistical mechanics, The Macmillan Co., New York; Collier-Macmillan Ltd., London, 1972, a Series of Books in Applied Mathematics.