format: text

text:
# The susceptibility of the square lattice Ising model: New developments

## Keywords

Ising susceptibility, high-temperature series, low-temperature series, scaling function, irrelevant variables, differentiably finite functions, scaling fields.

## 1 Introduction

Since Onsager’s  celebrated solution of the Ising model free energy in 1944, followed by Yang’s  proof of Onsager’s result for the spontaneous magnetization in 1952, almost half a century has passed during which time many, if not most of the world’s most able mathematical physicists have devoted themselves to the problem of elucidating the susceptibility. Attempting to list all these contributions would produce a bibliography of prohibitive length, and one that would inevitably commit many sins of omission. Therefore rather than attempt this, we will only make mention of those papers that have directly motivated our work here, and crave the forgiveness of those who we have inadvertently offended.

While much of the notation for describing the square lattice Ising model is standard, we begin by defining our notation here both for the benefit of the more casual reader and to emphasize those cases where we deviate from convention. The interactions in the two perpendicular directions are taken to be

$$K=\beta J,K^{}=\beta J^{}.$$
(1.1)
but we also often set $`K^{}=K`$ to discuss the isotropic lattice. For high temperatures, $`s=\mathrm{sinh}2K`$ and $`s^{}=\mathrm{sinh}2K^{}`$ are appropriate variables for series expansions , while for low temperatures, we use $`1/s`$ and $`1/s^{}`$ instead. Thus,

$$s^{}=\mathrm{sinh}2K^{}=1/\mathrm{sinh}2K=1/s,s_{}^{}{}_{}{}^{}=\mathrm{sinh}2K_{}^{}{}_{}{}^{}=1/\mathrm{sinh}2K^{}=1/s^{}.$$
(1.2)
In many cases, high-temperature and low-temperature formulas can be obtained from each other by Kramers-Wannier duality with a simple interchange of primes and stars. The critical temperature is defined by the condition $`s^{}=s^{}`$.

A conventional high-temperature variable is $`v=\mathrm{tanh}K`$, while an often-used low-temperature variable is $`u=\mathrm{exp}(4K)`$. The translations between these and our variables are

$$s=\mathrm{sinh}2K=2v/(1v^2)$$
(1.3)
and

$$s^{}=2u^{1/2}/(1u),$$
(1.4)
and similarly for the primed variables.

In studying the critical behavior, we will use both the variable $`t=1T_c/T`$ and more frequently

$$\tau =(1/ss)/2\text{ (isotropic)}$$
(1.5)
to parameterize deviations from the critical temperature. To leading order, $`\tau =2K_c\sqrt{2}t.`$

An elliptic parameterization will be useful, and to that end we define the elliptic modulus,

$$k=\{\begin{array}{cc}s^{}/s^{}=ss^{}\hfill & \text{for }T>T_c\text{,}\hfill \\ s^{}/s^{}=1/ss^{}\hfill & \text{for }T<T_c\text{.}\hfill \end{array}$$
(1.6)
Let $`\sigma _{i,j}`$ be the spin at lattice site $`(i,j)`$ and define the two-point function

$$C(M,N)=\sigma _{0,0}\sigma _{M,N}.$$
(1.7)
In terms of the correlation functions the susceptibility is

$$\beta ^1\chi =(C(M,N)^2)$$
(1.8)
where $``$ is the magnetization.

In 1956 Syozi and Naya  presented an approximation to the anisotropic high-tempera-ture susceptibility which gave the correct critical point, correct critical exponent, an amplitude estimate that was wrong by less than 1%, and reproduced the first 8 terms of the series expansion. It was also exact along the disorder line111For the fully anisotropic triangular lattice Ising model, with coupling constants $`v_i=\mathrm{tanh}J_i/kT,`$ $`i=1,2,3,`$ the condition for the disorder line is that $`v_1v_2+v_3=0.`$ Along this line the correlations decay exponentially, and the partition function factorizes . The disorder condition is of no particular interest in the case of the nearest neighbor square lattice Ising model. In 1976 a celebrated paper by Wu, McCoy, Tracy and Barouch  showed how the high- and low-temperature expansions of the susceptibility could be understood in terms of a multi-particle expansion, with an odd number of particles being appropriate at high temperatures and an even number at low temperatures. With this interpretation it became clear that the result of  was just the lowest order, or one-particle approximation to the susceptibility.

In terms of the elliptic modulus $`k`$ (1.6) the high- and low-temperature susceptibilities can be written

$$\beta ^1\chi _+=k^{\frac{1}{2}}(1k^2)^{\frac{1}{4}}\underset{l=0}{\overset{\mathrm{}}{}}\widehat{\chi }^{(2l+1)}$$
(1.9)
and

$$\beta ^1\chi _{}=(1k^2)^{\frac{1}{4}}\underset{l=1}{\overset{\mathrm{}}{}}\widehat{\chi }^{(2l)}$$
(1.10)
respectively, where $`\widehat{\chi }^{(j)}`$ is the sum over all lattice separations of the $`j`$-particle contribution to the two-point function, and was first given  as a $`2j`$-fold multiple integral. It was subsequently shown that this integral can be reduced to a $`j`$-fold integral of the form

$$\widehat{\chi }^{(j)}=\frac{k^{j/2}}{(2\pi )^jj!}𝑑u_1\mathrm{}𝑑u_j(G^{(j)})^2f^{(j)},$$
(1.11)
where $`G^{(j)}`$ is a fermionic determinant and $`f^{(j)}`$ is an algebraic function. This reduction has been achieved by various routes . The factor, $`G^{(j)}`$, appearing in the integrand, which can be expressed in terms of Pfaffians , has been found to have a product form, first by Palmer and Tracy in the low temperature regime , and then independently by Yamada in both the high- and low-temperature regimes . From the product form it readily follows that the first non-zero term in (1.11) is $`2^{j(1j)}k^{j^2/2}.`$ From the original expression one could only conclude that each integral entered at order $`k^{j/2},`$ so clearly massive cancellations occur. While this comes as a surprise if handling the integral  directly, it follows straightforwardly  from the product form.

In terms of the elliptic modulus, the first terms in the high- and low-temperature expansions for the isotropic ($`K=K^{}`$) susceptibility are

$$\widehat{\chi }_{\text{iso}}^{(1)}=\frac{k^{\frac{1}{2}}}{(1k^{\frac{1}{2}})^2}$$
(1.12)
and

$$\widehat{\chi }_{\text{iso}}^{(2)}=\frac{(1+k^2)\mathrm{E}(1k^2)\mathrm{K}}{3\pi (1k)(1k^2)},$$
(1.13)
where E and K are the complete elliptic integrals222We trust there is no confusion with the coupling constant $`K`$. of the second and first kind respectively. Anisotropic versions of these formulae are given in section 3. Simple forms for higher terms in the expansion are not known.

In  one of us gave compelling evidence that unlike the free-energy and spontaneous magnetization, the anisotropic susceptibility $`\chi (K,K^{})`$ is not differentiably finite. A series in $`n`$ variables, $`f(𝐳)`$, is said to be differentiably finite or D-finite if and only if it satisfies a system of $`n`$ partial differential equations of the form

$$P_{i,0}(𝐳)f(𝐳)+P_{i,1}(𝐳)\frac{}{z_i}f(𝐳)+\mathrm{}+P_{i,k_i}(𝐳)\frac{^{k_i}}{z_i^{k_i}}f(𝐳)=0,$$
(1.14)
where the $`P_{i,j}(𝐳)`$ are polynomials and for each $`i=1,\mathrm{},n`$, $`P_{i,k_i}(𝐳)`$ is not the null polynomial, see e.g. Proposition 2.2 in . Thus the expression for the susceptibility was shown to be in a different—and less tractable—class of function than other known properties of the Ising model. The evidence for this was based on the observation (not proved) that the anisotropic susceptibility $`\chi (v,v^{}),`$ as a function of $`v`$ with $`v^{}`$ fixed has a natural boundary on the unit circle $`|v|=1.`$

For the isotropic susceptibility, another of us  provided strong confirmation (though again, not a proof) of this observation by showing that the circle $`|s|=1`$ in the complex $`s=\mathrm{sinh}2K`$ plane is a natural boundary. These two observations are discussed further in section 3.2, following a discussion of the general anisotropic case in section 3.1.

In section 3.3 we also prove the important result that while $`\chi (k)`$ is not D-finite, $`\widehat{\chi }^{(j)}(k)`$ is D-finite for all $`j.`$

Two other directions in which we have achieved substantial progress are in the generation of series coefficients for the individual series $`\widehat{\chi }^{(j)},`$ continuing work initiated in , and even greater progress in obtaining the coefficients of the total series $`\chi `$ using nonlinear partial difference equations for the correlation functions .

In order to generate the series for the total susceptibility $`\chi _+`$ or $`\chi _{}`$ without computing separately the $`j`$-particle contributions, a more efficient method of series generation is obtained by first returning to the expression (1.8) of the susceptibility as the sum over all lattice separations of the two-point correlation functions.

In the scaling limit, the two-point functions were found to satisfy a nonlinear differential equation of Painlevé type  which was then solved to give the leading scaling terms $`|\tau |^{7/4}`$ and $`|\tau |^{3/4}`$ exactly . In 1980, a discrete analogue of this equation was discovered by McCoy and Wu  which holds for the two-point functions of the lattice Ising model at arbitrary temperature, and which reduces to the Painlevé equation in the scaling limit. In the same year, a simple set of partial difference equations was derived by Perk , which reproduced the equation of McCoy and Wu and provided an additional equation. Also in the same year, an unrelated set of difference equations, which can be used to compute the correlation functions on the diagonal, $`M=N`$, was obtained by Jimbo and Miwa . Many of these developments are described in some detail in .

The difference equations are valid for arbitrary temperature and were used by Kong, et al.  to obtain exactly the leading “short distance” constant terms in the susceptibility both at the ferromagnetic and anti-ferromagnetic points at $`T=\pm T_c.`$ This work was later extended to give the amplitudes of the term $`\tau \mathrm{log}|\tau |`$ . Here we dramatically extend that work by obtaining all terms in the “short-distance” part333When we speak of the “short-distance” part, we include the analytic background term, often ignored in scaling discussions of the critical region. $`B_{\mathrm{f}/\mathrm{af}}`$ (see (1.16) and (1.17)) of the susceptibility to O$`(\tau ^{14}).`$ Details necessary for the generation of $`C(M,N)`$ appear in section 4 while in section 6.1 we describe the numerical analysis of the $`C(M,N)`$ that leads us to conclude that the “short-distance” terms have the form

$$B_{\mathrm{f}/\mathrm{af}}=\underset{q=0}{\overset{\mathrm{}}{}}\underset{p=0}{\overset{\sqrt{q}}{}}b_{\mathrm{f}/\mathrm{af}}^{(p,q)}\tau ^q(\mathrm{log}|\tau |)^p\text{ (isotropic)}.$$
(1.15)
Although clearly nonanalytic at $`\tau =0`$ we have denoted these “short-distance” terms in (1.15) by $`B,`$ as a further reminder that they also include the analytical “background.” The actual coefficients in (1.15) can be found in the Appendix.

The quadratic difference equations of Perk  can also be used to generate high- and low-temperature series for $`\chi `$ and as shown in section 4 the series coefficients can be obtained in polynomial time!

While some people have expressed the view that a polynomial time algorithm for the computation of the series coefficients constitutes a solution, it is clearly preferable to have a closed form expression. Nevertheless, a polynomial time algorithm is equivalent to a complete solution if one seeks only the series coefficients, as to expand any closed form expression also takes polynomial time. As the history of the development of these key nonlinear recurrences described above shows, the ingredients for such an algorithm have existed unexploited in the literature for many years.

Analysis of the resulting series of hitherto unimaginable length combined with the “short-distance” knowledge contained in (1.15) leads to a solidly based conjecture specifying completely the remaining “scaling” part of the susceptibility $`\chi `$ of the isotropic Ising model. Near the anti-ferromagnetic point the “short-distance” terms are complete and we simply have for $`T>T_c`$

$$\beta ^1\chi _{\mathrm{af}}=B_{\mathrm{af}}\text{ (isotropic)}.$$
(1.16)
Near the ferromagnetic point we conjecture (for $`T>T_c`$ or $`T<T_c`$)

$$\beta ^1\chi _\pm =C_{0_\pm }(2K_c\sqrt{2})^{7/4}|\tau |^{7/4}F_\pm +B_\mathrm{f}\text{ (isotropic)}$$
(1.17)
where the scaling-amplitude functions $`F_\pm `$ have (possibly asymptotic) expansions in integer powers of $`\tau `$ without any of the powers of $`\mathrm{log}|\tau |`$ present in $`B_{\mathrm{f}/\mathrm{af}}.`$ The leading terms are

$$F_\pm =1+\tau /2+5\tau ^2/8+3\tau ^3/1623\tau ^4/38435\tau ^5/768+f_\pm ^{(6)}\tau ^6+\mathrm{O}(\tau ^7)\text{ (isotropic)}$$
(1.18)
where $`f_+^{(6)}f_{}^{(6)}.`$ In fact the breakdown in equality is dramatic, and we estimate that $`f_+^{(6)}=0.1329693327\mathrm{},`$ and $`f_{}^{(6)}=6.330746944\mathrm{},`$ where more accurate values of these and further terms in (1.18) through order $`\tau ^{15}`$ are given in the Appendix. In section 6.2 we describe the “short-distance” subtraction and analysis on which the assumed form of the expansion of the scaling-amplitude function $`F_\pm `$ (1.18) is based, while the analysis leading to the numerical values of the coefficients in (1.18) is sketched in section 6.3.

Aharony and Fisher  have predicted a scaling-amplitude function $`F(\mathrm{A}\&\mathrm{F})`$ that is equal above and below $`T_c`$ on the assumption that the Ising model critical region can be described entirely by two nonlinear scaling fields. Our exact result (1.18) is clearly different and furthermore the explicit expansion (cf. eqs. (22-24) in )

$$F(\mathrm{A}\&\mathrm{F})=1+\tau /2+5\tau ^2/8+3\tau ^3/1611\tau ^4/19217\tau ^5/384+97\tau ^6/3072+\mathrm{O}(\tau ^7)$$
(1.19)
differs from (1.18) at order $`\tau ^4.`$ This is unequivocal evidence for the presence of at least one, and almost surely two, irrelevant operators444Note that $`\tau `$ changes sign in (1.18) and (1.19) as we change from $`T>T_c`$ to $`T<T_c.`$. There is further possible evidence for irrelevant operators in the “short-distance” terms (1.15) which contain powers of $`\mathrm{log}|\tau |`$ beyond the first starting at $`\tau ^4(\mathrm{log}|\tau |)^2`$, and thus are not of the “energy” form given by the nonlinear field analysis.

An important numerical study investigating corrections to scaling was that of Gartenhaus and McCullough  who confirmed the $`F(\mathrm{A}\&\mathrm{F})`$ form in (1.19) through O$`(\tau ^3)`$ and provided a good estimate of the term linear in $`\tau `$ in $`B_\mathrm{f}`$ in (1.17). Estimates of terms to O$`(\tau ^2\mathrm{log}|\tau |)`$ in $`B_{\mathrm{af}}`$ in (1.16) were obtained in . An attempt  to go beyond this using longer series than available in  was inconclusive other than to indicate the necessity of terms beyond that predicted by Aharony and Fisher .

Our numerical work began in part as a modest attempt to improve on  but expanded to where it now clearly quantifies the effect of nonlinear scaling fields  as well as irrelevant operators. As anisotropy is a marginal operator, extending the present calculation to the anisotropic square lattice would, we expect, be extremely helpful in better understanding the effect of the irrelevant operator(s) we have identified. Both this and a study of the susceptibility on the triangular and hexagonal lattices are projects we hope to tackle in the near future.

The layout of the paper is as follows. In section 2 we define the model and give some useful parameterizations. In section 3 we show how the key integrals referred to above may be simplified, and provide a short proof of the assertion that $`\widehat{\chi }^{(j)}`$, so defined, is D-finite, even though, as we have seen, $`\chi `$ is presumably not. In section 4 the computation of the susceptibility series from the correlation functions by means of nonlinear partial difference equations is shown to be achievable in polynomial time. In section 5 we first discuss the isotropic series in $`k`$, and then the $`q`$ series form of the susceptibility. In this subsection some regularity features of the coefficients are discovered and the consequences partially developed. In section 6 we summarize our numerical work, state our conjecture giving the complete analytic structure of the isotropic susceptibility, and quantify the effect of irrelevant variables in the scaling fields. Finally, in section 7 we review scaling theory as it applies to the two-dimensional Ising model, and comment on the relevance of our results to this theory and to the renormalization group. The high- and low-temperature series for $`\chi _{\mathrm{iso}}`$ are given on the WWW at site www.ms.unimelb.edu.au/$``$tonyg.

## 2 Definitions and notation

In this section of the paper it will be convenient to formulate the general anisotropic model. Our later numerical work is confined to the isotropic model.

The complementary modulus of the elliptic modulus $`k`$ (1.6) is given by $`k^{}=\sqrt{1k^2}`$. Here the prime is not related to the anisotropy. The moduli, $`k`$ and $`k^{}`$, are related to the elliptic nome, $`q`$, by

$$k=\left(\frac{\theta _2(0,q)}{\theta _3(0,q)}\right)^2,k^{}=\left(\frac{\theta _4(0,q)}{\theta _3(0,q)}\right)^2.$$
(2.1)
The connection between the two moduli yields a well-known theta function identity. (See section 13.20 of Bateman  for these and other formulas.)

In terms of these variables the magnetization takes a particularly simple form

$$=(1(ss^{})^2)^{\frac{1}{8}}=k^{{}_{}{}^{}\frac{1}{4}}=\underset{n=1}{\overset{\mathrm{}}{}}\frac{1q^{2n1}}{1+q^{2n1}}$$
(2.2)
for $`T<T_c,`$ and $`0`$ otherwise.

The theta functions introduced above have useful infinite sum and product forms, which we will subsequently use. For zero argument these are

$`\theta _2(0,q)`$ $`=2{\displaystyle \underset{n=0}{\overset{\mathrm{}}{}}}q^{(n+1/2)^2}=2q^{1/4}{\displaystyle \underset{n=1}{\overset{\mathrm{}}{}}}(1q^{2n})(1+q^{2n})^2,`$ (2.3)
$`\theta _3(0,q)`$ $`=1+2{\displaystyle \underset{n=1}{\overset{\mathrm{}}{}}}q^{n^2}={\displaystyle \underset{n=1}{\overset{\mathrm{}}{}}}(1q^{2n})(1+q^{2n1})^2,`$ (2.4)
$`\theta _4(0,q)`$ $`=1+2{\displaystyle \underset{n=1}{\overset{\mathrm{}}{}}}(1)^nq^{n^2}={\displaystyle \underset{n=1}{\overset{\mathrm{}}{}}}(1q^{2n})(1q^{2n1})^2.`$ (2.5)
Thus,

$$k^{1/4}=\sqrt{2}q^{1/8}\underset{n=1}{\overset{\mathrm{}}{}}\frac{1+q^{2n}}{1+q^{2n1}},k_{}^{}{}_{}{}^{1/4}=\underset{n=1}{\overset{\mathrm{}}{}}\frac{1q^{2n1}}{1+q^{2n1}}.$$
(2.6)
In the elliptic parameterization, one variable is taken to be either $`k`$ or $`q`$. We take the other variable to be the anisotropy parameter. Onsager  used the variable $`a`$ defined by

$$\mathrm{sn}ia=\{\begin{array}{cc}is^{}\hfill & \text{for }T>T_c\text{,}\hfill \\ is^{}\hfill & \text{for }T<T_c\text{.}\hfill \end{array}$$
(2.7)
We will also need the related variable

$$\mathrm{sn}ia^{}=\{\begin{array}{cc}i/s^{}\hfill & \text{for }T>T_c\text{,}\hfill \\ i/s^{}\hfill & \text{for }T<T_c\text{.}\hfill \end{array}$$
(2.8)
In some contexts, it will be useful, following ref. , to use instead the variables $`\alpha `$ and $`\alpha ^{}`$ defined by

$`\mathrm{cot}\alpha `$ $`=i\sqrt{k}\mathrm{sn}ia=\sqrt{s^{}s^{}}=\sqrt{s^{}/s},`$ (2.9)
$`\mathrm{cot}\alpha ^{}`$ $`=i\sqrt{k}\mathrm{sn}ia^{}=\sqrt{1/(s^{}s^{})}=\sqrt{s/s^{}}.`$ (2.10)
Obviously $`\mathrm{cot}\alpha =\mathrm{tan}\alpha ^{}`$. The isotropic values of these variables are $`\mathrm{sn}ia=\mathrm{sn}ia^{}=i/\sqrt{k}`$ and $`\alpha =\alpha ^{}=\pi /4`$.

We observe that $`q^{1/2}=x`$ where $`x`$ is the variable of ref. , defined by

$$e^{2K}=x^{1/2}\underset{n=1}{\overset{\mathrm{}}{}}\frac{(1x^{8n7})(1x^{8n1})}{(1x^{8n5})(1x^{8n3})},$$
(2.11)
in which the spontaneous magnetization has a simple product form. The elliptic modulus of ref.  is related to our $`k`$ by a Landen transformation \[28, section 15.6\].

In reference  it was also noted that the first terms of the isotropic ($`K=K^{}`$) high-temperature susceptibility have a simple product expression, which breaks down at order $`q^{8/4}`$. This is explained by the expansion of the high-temperature susceptibility in multi-particle states described in the Introduction. The contribution of one-particle states has a product form, and the first three-particle state contributes at order $`q^{8/4}`$. This one-particle product form is

$$\widehat{\chi }_{\text{iso}}^{(1)}=q^{1/4}\underset{n=1}{\overset{\mathrm{}}{}}\frac{(1+q^{n/4})^2(1q^{4n})^2}{(1q^{n/4})^2(1+q^n)^2}=\frac{q^{1/4}}{\theta _3(0,q)}\underset{n=1}{\overset{\mathrm{}}{}}\frac{(1+q^{n/4})^2(1q^{2n})^3}{(1q^{n/4})^2}.$$
(2.12)
(Strictly speaking the product form above is not the product form of ref.  which actually breaks down only at order $`q^{9/4}`$. They differ by a factor which first contributes at order $`q^{8/4}`$ and the absence of this factor in ref.  exactly compensates for the addition of the first term of $`\widehat{\chi }_{\text{iso}}^{(3)}`$.)

Attempts to do the same for the low-temperature series have failed as we do not know of a product form for $`\widehat{\chi }_{\text{iso}}^{(2)}`$. The function K has one, namely

$$\mathrm{K}=\frac{1}{2}\pi \theta _3^2(0,q),$$
(2.13)
while the function E is given in terms of $`q`$ by the formula

$$\mathrm{E}=\frac{\theta _3^4(0,q)+\theta _4^4(0,q)}{3\theta _3^4(0,q)}\mathrm{K}\frac{1}{12\mathrm{K}}\frac{\theta _1^{\prime \prime \prime }(0,q)}{\theta _1^{}(0,q)}.$$
(2.14)
## 3 Integral formulae for $`\widehat{\chi }^{(j)}`$

In this section we present integral expressions for the $`\widehat{\chi }^{(j)}`$ and show that they define D-finite functions.

### 3.1 Trigonometric form of integrals

We will start with the form given by Yamada  for the $`j`$-particle contribution $`\widehat{\chi }^{(j)}`$ in terms of elliptic variables, and derive from it an expression in terms of trigonometric/hyperbolic variables. The integral expression is

$$\widehat{\chi }^{(j)}=\frac{k^{j/2}}{(2\pi )^jj!}_0^{4\mathrm{K}}𝑑u_1\mathrm{}_0^{4\mathrm{K}}𝑑u_j\left(G^{(j)}\right)^2\frac{1+_{n=1}^jx_n}{1_{n=1}^jx_n}\frac{1+_{n=1}^jz_n}{1_{n=1}^jz_n},$$
(3.1)
where the modulus of the complete elliptic integral $`\mathrm{K}`$ is $`k`$, and this same modulus is assumed in all Jacobi elliptic functions which appear below. Considered as a function of one of the $`u_n`$, the integrand has a single simple pole on the real axis which derives from the last factor in the integrand. The contour of integration is deformed in the vicinity of the pole so that only half the residue is taken. The function $`G^{(j)}`$ can be written as the product

$$G^{(j)}=\underset{1m<nj}{}h_{mn},$$
(3.2)
and $`z_n`$, $`x_n`$ and $`h_{mn}`$ are given by

$`z_n=e^{i\omega _n}`$ $`={\displaystyle \frac{\mathrm{sn}\frac{1}{2}(u_n+ia^{})}{\mathrm{sn}\frac{1}{2}(u_nia^{})}}={\displaystyle \frac{\mathrm{sn}ia^{}\mathrm{cn}u_n+\mathrm{cn}ia^{}\mathrm{sn}u_n}{\mathrm{sn}ia^{}\mathrm{dn}u_n\mathrm{dn}ia^{}\mathrm{sn}u_n}},`$
$`x_n=e^{\gamma _n}`$ $`=k\mathrm{sn}\frac{1}{2}(u_n+ia^{})\mathrm{sn}\frac{1}{2}(u_nia^{})=k{\displaystyle \frac{\mathrm{cn}ia^{}\mathrm{cn}u_n}{\mathrm{dn}ia^{}+\mathrm{dn}u_n}},`$
$`h_{mn}`$ $`=\sqrt{k}\mathrm{sn}\frac{1}{2}(u_mu_n).`$ (3.3)
Note that the trigonometric/hyperbolic variables, $`\omega _n`$ and $`\gamma _n`$ are related to the elliptic variables, $`u_n`$ and $`a^{}`$ by the functional equation

$$\mathrm{exp}(\pm \frac{1}{2}i\omega _n\frac{1}{2}\gamma _n)=\sqrt{k}\mathrm{sn}\frac{1}{2}(u_n\pm a^{}).$$
(3.4)
The choice of $`a^{}`$ rather than $`a`$ as anisotropy parameter is arbitrary because of the symmetry under interchange of horizontal and vertical lattice axes. Choosing $`a^{}`$ at this point results in expressions for $`\omega _n`$ and $`\gamma _n`$ which are equivalent to those of Onsager.

The mapping between the elliptic parameterization and the trigonometric/hyperbolic one was described in ref.  for the isotropic case. The formulas below are a generalization of this mapping. It is simpler to use the trigonometric parameterization than it is to use the elliptic parameterization for numerical series generation.

The variables $`\omega _n`$ and $`\gamma _n`$ defined above satisfy the identities given in Appendix 2 of Onsager’s paper ,

$`\mathrm{cosh}\gamma _n`$ $`=(\mathrm{cn}ia\mathrm{dn}u_nk\mathrm{dn}ia\mathrm{cn}u_n)/M_n`$ (3.5)
$`\mathrm{sinh}\gamma _n`$ $`=ik_{}^{}{}_{}{}^{2}\mathrm{sn}ia/M_n`$ (3.6)
$`\mathrm{cos}\omega _n`$ $`=(\mathrm{dn}ia\mathrm{cn}u_nk\mathrm{cn}ia\mathrm{dn}u_n)/M_n`$ (3.7)
$`\mathrm{sin}\omega _n`$ $`=k_{}^{}{}_{}{}^{2}\mathrm{sn}u_n/M_n,`$ (3.8)
with

$$M_n=\mathrm{dn}ia\mathrm{dn}u_nk\mathrm{cn}ia\mathrm{cn}u_n.$$
(3.9)
Notice that it is $`a`$ which appears in these formulas rather than $`a^{}`$. The mapping from the variables $`u_n`$ and $`a`$ to the variables $`\omega _n`$ and $`\gamma _n`$ is conformal, as is seen from the formulas, also given by Onsager

$`\omega _n/u_n`$ $`=\gamma _n/a=k_{}^{}{}_{}{}^{2}/M_n=i\mathrm{sinh}\gamma _n/\mathrm{sn}ia`$
$`\gamma _n/u_n`$ $`=\omega _n/a=k_{}^{}{}_{}{}^{2}k\mathrm{sn}ia\mathrm{sn}u_n/M_n.`$ (3.10)
Finally, Onsager gives the functional equation

$$\mathrm{cot}\frac{1}{2}(\omega _ni\gamma _n)=(1+k)\mathrm{sc}\frac{1}{2}(u_n+ia)\mathrm{nd}\frac{1}{2}(u_n+ia).$$
(3.11)
To make the change of variables, we first note that $`\varphi _n=\pi `$ corresponds to $`u_n=0`$, and $`\varphi _n=0`$ corresponds to $`u_n=2\mathrm{K}`$. Using eq. (3.10), we see that

$$_\pi ^\pi 𝑑\omega \mathrm{cot}\alpha /\mathrm{sinh}\gamma \mathrm{}=\sqrt{k}_0^{4\mathrm{K}}𝑑u\mathrm{},$$
(3.12)
which implies

$$\widehat{\chi }^{(j)}=\frac{\mathrm{cot}^j\alpha }{j!}_\pi ^\pi \frac{d\omega _1}{2\pi }\mathrm{}_\pi ^\pi \frac{d\omega _{j1}}{2\pi }\left(\underset{n=1}{\overset{j}{}}\frac{1}{\mathrm{sinh}\gamma _n}\right)\left(\underset{1i<kj}{}h_{ik}\right)^2\frac{1+_{n=1}^jx_n}{1_{n=1}^jx_n}.$$
(3.13)
The condition $`\omega _1+\mathrm{}+\omega _j=0`$ mod $`2\pi `$, which results from having performed one of the integrations, is assumed. In terms of $`\omega _n`$, $`k`$ and $`\alpha `$, the quantities $`x_n`$ and $`\mathrm{sinh}\gamma _n`$ can be expressed as

$`x_n`$ $`=\mathrm{cot}^2\alpha \left[\xi \mathrm{cos}\omega _n\sqrt{(\xi \mathrm{cos}\omega _n)^2(\mathrm{cot}\alpha )^4}\right],`$ (3.14)
$`\mathrm{sinh}\gamma _n`$ $`=\mathrm{cot}^2\alpha \sqrt{(\xi \mathrm{cos}\omega _n)^2(\mathrm{cot}\alpha )^4}`$ (3.15)
with

$$\xi =\left(1+1/(k\mathrm{cot}^2\alpha )\right)^{1/2}\left(1+k/\mathrm{cot}^2\alpha \right)^{1/2}=\left(1+(s^{})^2\right)^{1/2}\left(1+(s^{})^2\right)^{1/2}.$$
(3.16)
In the isotropic case, these reduce to

$`x_n`$ $`=s+s^1\mathrm{cos}\omega _n\sqrt{(s+s^1\mathrm{cos}\omega _n)^21},`$ (3.17)
$`\mathrm{sinh}\gamma _n`$ $`=\sqrt{(s+s^1\mathrm{cos}\omega _n)^21}.`$ (3.18)
Here $`s=\sqrt{k}`$ for high-temperature and $`s=1/\sqrt{k}`$ for low temperature, although the distinction is irrelevant as the dependence of the integrand on $`s`$ and $`1/s`$ is symmetric. Finally

$`h_{ik}`$ $`=\mathrm{cot}\alpha {\displaystyle \frac{\mathrm{sin}\frac{1}{2}(\omega _i\omega _k)}{\mathrm{sinh}\frac{1}{2}(\gamma _i+\gamma _k)}}={\displaystyle \frac{1}{\mathrm{cot}\alpha }}{\displaystyle \frac{\mathrm{sinh}\frac{1}{2}(\gamma _i\gamma _k)}{\mathrm{sin}\frac{1}{2}(\omega _i+\omega _k)}}`$ (3.19)
$`={\displaystyle \frac{2(x_ix_k)^{1/2}\mathrm{cot}\alpha \mathrm{sin}\frac{1}{2}(\omega _i\omega _k)}{1x_ix_k}}.`$ (3.20)
When $`j=1`$ or $`j=2`$, the integrals can be rewritten in terms of known functions as was noted in the introduction for the isotropic case. When $`j=3`$ Glasser  showed that the integrals can be written as an integral involving the square root of a polynomial of degree higher than 4. Such integrals are often called hyperelliptic integrals, and are a special case of Abelian integrals.

The corresponding anisotropic expressions are

$`\widehat{\chi }^{(1)}=`$ $`{\displaystyle \frac{1}{(k^{1/2}k^{1/2})^2}}\left[2\mathrm{csc}2\alpha +\sqrt{(k^{1/2}k^{1/2})^2+4\mathrm{csc}^22\alpha }\right],`$ (3.21)
$`\widehat{\chi }^{(2)}=`$ $`{\displaystyle \frac{k^{1/2}}{1+k}}\sqrt{(k^{1/2}k^{1/2})^2+4\mathrm{csc}^22\alpha }\widehat{\chi }_{\text{iso}}^{(2)}.`$ (3.22)
As we show in section 5, we can say more about the general form of the expansion, based on inspection of the long series published by Nickel , as well as more recent extensions reported here, giving hope that there is still more regularity to be found.

### 3.2 Natural boundaries

As mentioned in the introduction, there were observations, first in  and then in , that strongly suggested the susceptibility of the Ising model is a function with a natural boundary unlike the free-energy or magnetization. We expand and clarify those observations here in the light of our new knowledge of the general anisotropic case discussed in section 3.1 and the additional numerical work on the isotropic limit described in section 6.

First, we expect that as described in  , the $`\chi ^{(j)}`$ given by the integrals (3.13) multiplied by the factors outside the summation in eq. (1.9) or (1.10) will be singular at the symmetry points of the integrand and where the denominator factor $`1_{n=1}^jx_n`$ vanishes. The symmetry point condition requires all $`\omega _n`$ to be equal and given by $`\omega _n=\omega =2\pi m^{}/j,`$ $`m^{}=1,2,\mathrm{}j.`$ The vanishing of the denominator factor requires the $`x_n`$, now all equal, to be given by $`x_n=x=\mathrm{exp}(2\pi im/j),`$ $`m=1,2,\mathrm{}j.`$ Equivalently, from the explicit formula (3.14),

$$\mathrm{cot}^2(\alpha )(\xi \mathrm{cos}(2\pi m^{}/j))=\mathrm{cos}(2\pi m/j).$$
(3.23)
With $`\mathrm{cot}^2(\alpha )=s^{}/s`$ from (2.9) and $`\xi `$ given by (3.16), we find (3.23) can be reduced to

$$\mathrm{cosh}(2K)\mathrm{cosh}(2K^{})\mathrm{sinh}(2K)\mathrm{cos}\frac{2\pi m}{j}\mathrm{sinh}(2K^{})\mathrm{cos}\frac{2\pi m^{}}{j}=0$$
(3.24)
with $`m,m^{}=1,2,\mathrm{}j`$ as discussed above. It will be noted that the left-hand side of (3.24) is the denominator in the Onsager integral for the free-energy and thus we find the (to us) surprising result that the singularity of $`\chi ^{(j)}`$, a property of the Ising model in a magnetic field, is intimately connected with a property in zero field.

The full $`\chi `$, being a sum of $`\chi ^{(j)}`$, will naively be expected to be singular at a dense set of points and thus have the Onsager line (3.24) as a natural boundary. The presence of natural boundaries has implications for expansions about the physical singularity points $`s=\pm 1`$ that are necessary to understand corrections to scaling. We briefly explore some of these implications but restrict ourselves for simplicity to the ferromagnetic point $`s=1`$ in the isotropic model. We also make the plausible assumption that the singularity in each $`\chi ^{(j)}`$ closest to $`\tau =0`$ is the most important for determining, in expansions of $`\chi `$, the $`\tau ^p`$ large $`p`$ asymptotics and this considerably simplifies the discussion555Note that singularities on $`|s|=1`$ map to points on the imaginary axis in the complex $`\tau `$ plane..

Let $`\tau _j=i\mathrm{sin}\theta _j`$ be that singularity in $`\chi ^{(j)}`$ that is closest to the ferromagnetic $`\tau =0`$; $`\theta _j`$ is fixed by $`\mathrm{cos}\theta _j=(1+\mathrm{cos}\varphi _j)/2`$ with $`\varphi _j=2\pi /j.`$ Choose the branch-cut arising from this singularity to lie along the imaginary $`\tau `$ axis and directed away from $`\tau =0.`$ Take $`\tau =i𝒯`$ to be a point on the branch-cut. Provided the (positive) deviation

$$\delta \theta _j=\mathrm{arcsin}𝒯\theta _j$$
(3.25)
is not too large, the discontinuity across the cut can be estimated from the linearized singularity equations (14) in  and (12) in . For $`j>2`$ a more general result that includes the first order correction is

$`\mathrm{Disc}(\chi ^{(j)})`$ $`=`$ $`C_ji^{j^2}/\mathrm{sin}^2\varphi _j[\delta \theta _j\mathrm{sin}\theta _j/\mathrm{sin}^2\varphi _j]^{(j^23)/2}`$ (3.26)
$`\times `$ $`\left[1+\delta \theta _j\mathrm{sin}\theta _j\left\{{\displaystyle \frac{\mathrm{cos}\theta _j}{4\mathrm{sin}^2\theta _j}}+{\displaystyle \frac{j^24}{8\mathrm{cos}\theta _j(1+\mathrm{cos}\theta _j)}}{\displaystyle \frac{1}{2}}\right\}+\mathrm{O}(\delta \theta _j^2)\right]`$
where

$$C_j=2\sqrt{2}(2/\pi )^{(j3)/2}(j/2)^{(j^24)/2}(\mathrm{\Pi }_{m=1}^j\mathrm{\Gamma }(m))/\mathrm{\Gamma }((j^21)/2)3.7655j^{1/12}2^j\mathrm{exp}(j^2/4).$$
(3.27)
The last equality is valid only in the large $`j`$ limit. To obtain the discontinuity in say $`\chi _+`$, we must first sum the discontinuities in $`\chi ^{(j)}`$, $`j`$ odd, and this we can crudely estimate by integrating over $`j`$, keeping only the leading exponential factors in (3.25) and (3.26) and making a small angle approximation as well. Essentially the same formula is obtained for the discontinuity in $`\chi _{}`$ on summing over $`j`$ even; in either case we find

$$\mathrm{Disc}(\underset{j}{}\chi ^{(j)})_{\sqrt{2}\pi /𝒯}^{2\pi /𝒯}𝑑j[(𝒯j/(\sqrt{2}\pi )1)/(2\sqrt{\mathrm{e}})]^{j^2/2}$$
(3.28)
where the sum is over $`j`$ odd or $`j`$ even, and the lower limit is simply the restriction to those $`j`$ values that contribute, while the upper limit roughly defines the limit of validity of the linearized approximation. The precise value of this limit is not important since the integrand has a maximum well below the limit. For large $`j`$ only the maximum of the integrand matters and (3.28) reduces to

$$\mathrm{Disc}(\chi )\mathrm{exp}(39.76/𝒯^2)$$
(3.29)
in which the numerical coefficient of $`1/𝒯^2`$ is $`\pi ^2x^3/(2(x1))`$ with $`x`$ the smallest real solution of $`2\mathrm{log}((x1)/2)+1/(x1)=0.`$ Keeping terms such as the correction term in (3.26) and the $`2^j`$ in (3.27) in the steepest descent analysis lead to O$`(1/𝒯)`$ corrections to the exponent in (3.29) so we conjecture that the right hand side of (3.29) is exactly the leading exponential. It is interesting that the discontinuity (3.29) is similar to that found in weak coupling field theory expansions, but the mechanism producing the cut here does not seem to be related in any way to instantons.

The additional singularities that $`\chi `$ has at $`\tau =0,`$ namely terms such as the divergent “scaling” $`|\tau |^{7/4}`$ or the “short- distance” powers of $`\mathrm{log}|\tau |,`$ are a complication we do not know how to handle in any rigorous fashion666 If we knew that the “short-distance” corrections $`B_\mathrm{f}`$ in (1.17) formed a convergent sequence, then a subtraction process similar to that described in section 6.2 could be carried out here. The numerical evidence from section 6.2 is that the point $`\tau =0`$ in the scaled and pole subtracted $`\chi _\pm `$ is no longer a branch point singularity and hence the final result (3.30) is justifiable and would apply directly to the scaling-amplitude function $`F_\pm .`$ Unfortunately, we do not have any independent information to suggest that the “short-distance” sum (1.15) is convergent and thus the argument leading to (3.30) is at best suggestive that one or the other (or both) of the “short-distance” or “scaling” sequences are asymptotic.. As a consequence we will simply ignore them and make the simplest, yet reasonable, assumption that they make an additive contribution not relevant for understanding the effect of the natural boundary. Then the cut discontinuity (3.29) would imply a divergent behavior in the $`\tau `$ expansion of $`\chi `$. That is to say, the coefficient of $`\tau ^p`$ in the limit $`p\mathrm{}`$ will diverge as $`\mathrm{\Gamma }(p/2)/a^{p/2}`$ with $`a39.76`$. This follows from a contour integral around the origin distorted to run on either side of the cut imaginary axis. The contribution of the cut discontinuity to the coefficient, $`C_p`$, of $`\tau ^p`$ in the expansion is then

$$C_p𝑑𝒯/𝒯^{p+1}\mathrm{exp}(a/𝒯^2)\mathrm{\Gamma }(p/2)/a^{p/2}.$$
(3.30)
We know little about the cut discontinuity on the circle $`|s|=1`$ other than what we have deduced near $`s=1`$ as given by (3.29). However, the fact that the amplitudes of the singularities of $`\chi ^{(j)}`$ on $`|s|=1`$ vary dramatically with order $`j`$ almost certainly implies there will be no cancellation of singularities in the sum of $`\chi ^{(j)}`$ that defines $`\chi .`$ Furthermore, the variation with order means that there is no length scale at which, as one approaches the circle of singularities $`|s|=1,`$ $`\chi _{\mathrm{iso}}`$ can be smooth. It is these two points taken together that we consider overwhelming evidence that, in the isotropic case at least, $`\chi `$ has a natural boundary that is the entire $`|s|=1`$ circle. We do not imply by this that all points on the circle are equally “singular”. As argued above, the existence of an asymptotic expansion about $`s=1`$ seems likely. A very different situation arises at a point such as $`s=i.`$ While the singularity in $`\chi ^{(j)}`$ on the circle $`s=\mathrm{exp}(i\theta )`$ nearest $`s=1`$ lies at a distance $`\mathrm{\Delta }\theta =\mathrm{O}(1/j)`$ for large $`j,`$ the corresponding nearest distance from $`s=i`$ is $`\mathrm{\Delta }\theta =\mathrm{O}(1/j^2).`$ Furthermore this latter singularity is larger in its leading amplitude than the former by a factor roughly $`j^{(j^2)/2}.`$ The reduction in distance and dramatic increase in amplitude suggests that an asymptotic expansion about $`s=i`$ is not possible, but this has not been proved.

For the anisotropic case we have not analyzed (3.13) in detail so we do not have the necessary amplitude information to make the same claim directly. However, we can take the extreme anisotropic limit of $`s^{}`$ infinitesimal (but not $`0`$) and find that the Onsager line (3.24) has now come very close to the circle $`|v|=1.`$ At this point we can connect to the work in . There it was observed that if $`\chi (K,K^{})`$ was written as

$$\chi (v,v^{})=H_n(v)v_{}^{}{}_{}{}^{n}$$
(3.31)
the $`H_n`$ would be singular at a dense set of points on $`|v|=1`$ as $`n\mathrm{}.`$ Furthermore, if as in the discussion above, $`v^{}`$ is infinitesimal (but not $`0`$) the amplitudes of the singularities vary dramatically with order $`n`$ and the same conclusion as in the isotropic case is reached. Given that the Onsager line (3.24) is a natural boundary in two extremes, it seems highly probable that it is also a natural boundary at all intermediate anisotropy values.

We conclude this section by contrasting the above behavior of the susceptibility with that of the free-energy and magnetization which are only singular at an isolated set of points, not a dense set. This is precisely what one expects for a D-finite function, as in that case we have the following777Due to M. Bousquet-Mélou, private communication.

###### Theorem 1

Let $`f(x,y)=_{n0}y^nH_n(x)`$ be a D-finite series in $`y`$ with rational coefficients. For $`n0,`$ let $`S_n`$ be the set of poles of $`H_n(x)`$; let $`S=_nS_n.`$ Then $`S`$ has only a finite number of accumulation points.

This is observed in practice. The anisotropic magnetization and free-energy each have exactly one accumulation point , while the (non-D-finite) susceptibility appears to have an infinite number .

### 3.3 $`\chi ^{(j)}`$ is D-finite

This remarkable result, that $`\chi ^{(j)}`$ (or equivalently $`\widehat{\chi }^{(j)}`$) is D-finite while $`\chi `$ is not, follows from the results of Lipshitz  (see also Zeilberger ), who gives several basic definitions and theorems concerning D-finiteness of series in several variables. The integrand of the trigonometric/hyperbolic form of $`\widehat{\chi }^{(j)}`$ is an algebraic function of the variables $`c_j=\mathrm{cos}\omega _j`$ and $`k^{1/2}`$, and is thus a D-finite function in these same variables. Then by Theorem 2.7 of ref. , integrating over one or more variables preserves D-finiteness which implies the result. The term D-finite is synonymous with holonomic in much of the literature.

Kashiwara and Kawai  have shown that any Feynman diagram is holonomic, whereas an infinite sum of such diagrams may not be. This is just the phenomenon we observe here. At first glance it appears that Kashiwara’s definition of holonomic differs from that used here, but this is not so. The point is that the definition of D-finite functions of more than one variable requires that the underlying system of partial differential equations be such that only a finite number of initial conditions are needed to specify the function. Such systems are called “maximally over-determined” or “holonomic” in the analysis literature. In the single variable case, the question of a finite number of initial conditions is clearly automatically satisfied.

Motivated by the above observation, we have attempted to find linear differential equations with polynomial coefficients in $`k^{1/2}`$, or equivalently, linear recurrences for the series coefficients, of $`\widehat{\chi }_{\text{iso}}^{(3)}`$ and $`\widehat{\chi }_{\text{iso}}^{(4)}`$. With the available series of order $`k^{257/2}`$ and $`k^{62}`$ respectively, calculated using the methods of numerical integration described in , we have ruled out any such recurrences of depth 14 with coefficients of degree 15 for $`\widehat{\chi }_{\text{iso}}^{(3)}`$, and of depth 6 with coefficients of degree 7 for $`\widehat{\chi }_{\text{iso}}^{(4)}`$. Thus while these functions are provably D-finite, it is clear that the generating differential equation will be a fairly cumbersome object.

We have also attempted to fit these series as polynomials in the complete elliptic integrals K and E with polynomial coefficients in $`k`$ and obtained similar negative results.

## 4 Correlation functions and difference equations

In this section we give details of our more efficient method of series generation based on summing the correlation functions, obtained by means of nonlinear recurrences, as outlined in the Introduction.

The equations we used for generation of the very long series are the ones given by Perk . Here we present a slight generalization due to McCoy and Wu  which keeps track of the separate multi-particle components. The expansions of the two-point correlation functions in multi-particle components are analogous to the corresponding expansions (1.9) and (1.10) for the susceptibility

$$C(M,N;\lambda )=\{\begin{array}{cc}k^{1/2}(1k^2)^{1/4}_{n=0}^{\mathrm{}}\lambda ^{2n+1}\widehat{C}^{(2n+1)}(M,N)\hfill & \text{for }T>T_c\hfill \\ (1k^2)^{1/4}_{n=0}^{\mathrm{}}\lambda ^{2n}\widehat{C}^{(2n)}(M,N)\hfill & \text{for }T<T_c\text{,}\hfill \end{array}$$
(4.1)
with $`\widehat{C}^{(0)}(M,N)=1`$ and for $`j>0`$

$`\widehat{C}^{(j)}(M,N)=`$ (4.2)
$`{\displaystyle \frac{\mathrm{cot}^j\alpha }{j!}}{\displaystyle _\pi ^\pi }{\displaystyle \frac{d\omega _1}{2\pi }}\mathrm{}{\displaystyle _\pi ^\pi }{\displaystyle \frac{d\omega _j}{2\pi }}\left({\displaystyle \underset{n=1}{\overset{j}{}}}{\displaystyle \frac{1}{\mathrm{sinh}\gamma _n}}\right)\left({\displaystyle \underset{1i<kj}{}}h_{ik}\right)^2\left({\displaystyle \underset{n=1}{\overset{j}{}}}x_n\right)^M\mathrm{cos}(N{\displaystyle \underset{n=1}{\overset{j}{}}}\omega _n).`$
The fugacity $`\lambda `$ is associated with the number of particles, and $`C(M,N)=C(M,N;1)`$.

With these definitions the quadratic partial difference equations are

$`s^2[C(M,N;\lambda )^2C(M,N1;\lambda )C(M,N+1;\lambda )]`$
$`+[C^{}(M,N;\lambda )^2C^{}(M1,N;\lambda )C^{}(M+1,N;\lambda )]=0`$ (4.3)
$`s_{}^{}{}_{}{}^{2}[C(M,N;\lambda )^2C(M1,N;\lambda )C(M+1,N;\lambda )]`$
$`+[C^{}(M,N;\lambda )^2C^{}(M,N1;\lambda )C^{}(M,N+1;\lambda )]=0`$ (4.4)
$`ss^{}[C(M,N;\lambda )C(M+1,N+1;\lambda )C(M,N+1;\lambda )C(M+1,N;\lambda )]`$
$`=C^{}(M,N;\lambda )C^{}(M+1,N+1;\lambda )C^{}(M,N+1;\lambda )C^{}(M+1,N;\lambda ).`$ (4.5)
The object $`C^{}(M,N;\lambda )`$ is the correlation function on the dual lattice and is obtained from $`C(M,N;\lambda )`$ by replacing $`s^{}`$ with $`s^{}`$ and $`s`$ with $`s_{}^{}{}_{}{}^{}`$. Equation (4) holds for all $`M`$ and $`N`$. When $`M=N=0`$ equations (4) and (4) must be replaced by

$`C^{}(1,0;1)=\sqrt{1+s^2}sC(0,1;1)`$ (4.6)
$`C^{}(0,1;1)=\sqrt{1+s_{}^{}{}_{}{}^{2}}s^{}C(1,0;1).`$ (4.7)
We do not know of $`\lambda 1`$ versions of these equations.

These equations are nearly enough to determine all two-point functions completely. For the isotropic expansion, all that is lacking is either the high or the low temperature set of diagonal correlations ($`M=N`$). From either one of these the other can be obtained using equation (4) with $`M=N`$. When $`\lambda 1`$ we have used the integral formula (4.2) to compute the diagonal correlation functions. For the rest of this section we focus on the case $`\lambda =1`$ where two superior methods for obtaining the diagonal correlations are available. From the purely computational point of view the difference equations of Jimbo and Miwa  are almost certainly the most efficient and to be preferred. However, from the point of view of understanding the analytical structure of the correlations the original Toeplitz determinants  are better. We also find that for the numerical computations we have carried out so far the evaluation of the determinants is only a small fraction of the total project time so efficiency is not yet an issue.

We will restrict ourselves in the following to the isotropic lattice. In that case and for $`N>0,`$ the diagonal correlation $`C(N,N)`$ is the determinant of an $`N\times N`$ Toeplitz matrix with elements $`a_{i,j}=a_{ij}`$ that are the integrals

$$a_n=(2\pi )^1_0^\pi 𝑑\theta (s\mathrm{exp}(2i\theta )/s)\mathrm{exp}(2i\theta n)/\sqrt{\tau ^2+\mathrm{sin}^2\theta }.$$
(4.8)
These integrals apply both above and below $`T_c`$ and furthermore since $`\tau \tau `$ corresponds to $`s1/s,`$ one can establish from the integral in (4.8) the relations $`a_{n1}(\tau )=a_n(\tau ).`$ Explicit formulae for $`a_n`$ for small $`n`$ can be given in terms of elliptic integrals E and K but these expressions are not particularly enlightening and are not needed here. Rather we need the series expansions of (4.8), either in $`s,`$ $`1/s`$ or $`\tau `$ depending on the application.

The series expansions in $`s`$ and $`1/s`$ both for $`a_n`$ and $`C(N,N)`$ are completely straightforward with computer packages such as Maple that automatically handle the multiple precision arithmetic required. Furthermore, these same packages can be set to treat the $`C(M,N)`$ in the recursion formulae (4)–(4) as series and thus very little programming is necessary to generate the susceptibility. Admittedly, some steps need to be taken to conserve time and/or memory resources but this is very hardware dependent and will not be described here. What is worth noting, however, is that for a series to order $`N`$ the recursion formulae require O$`(N^2)`$ multiplications of series of length $`N`$ and thus in a naive implementation, O$`(N^4)`$ multiplications. Since the word length grows linearly with $`N`$ the algorithm has complexity of at most O$`(N^6)`$. There are more efficient ways to multiply long series and numbers with a large number of digits  but we have not found it necessary to explore these options.

Timing proportional to $`N^6`$ is observed in practice, in our implementation of the recursion in Maple. We have generated high-temperature series of order 323 and low-temperature series of order 646 for the isotropic susceptibility. The entire calculation took 123 hours on a 500MHz DEC Alpha with 21164 processor running Maple V version 5.1. We have also obtained shorter anisotropic series in this way (either the nearest off-diagonal correlation functions, or additional assumptions are needed). As lattice anisotropy is a marginal operator, we hope that an extension of this calculation will be very revealing.

The series in $`\tau `$ is most easily obtained by expressing the Toeplitz element integral (4.8) in terms of hypergeometric functions. To show this connection we start by writing

$$\pi a_n/2=\sqrt{1+\tau ^2}A_s\tau A_c$$
(4.9)
in terms of the real integrals

$`A_s={\displaystyle _0^{\pi /2}}𝑑\theta \mathrm{sin}(\theta )\mathrm{sin}(\nu \theta )/\sqrt{\tau ^2+\mathrm{sin}^2\theta },`$
$`A_c={\displaystyle _0^{\pi /2}}𝑑\theta \mathrm{cos}(\theta )\mathrm{cos}(\nu \theta )/\sqrt{\tau ^2+\mathrm{sin}^2\theta },`$ (4.10)
where $`\nu =2n+1.`$ The required symmetry $`a(\nu ,\tau )=a(\nu ,\tau )`$ is explicit in eqs. (4.9) and (4) and by standard integration by parts manipulation one can show that the cosine integral satisfies the differential equation

$$(1+\tau ^2)(d/d\tau )\tau (d/d\tau )A_c\nu ^2\tau A_c=0$$
(4.11)
while the sine integral can be expressed as the derivative

$$A_s=(\tau /\nu )(d/d\tau )A_c.$$
(4.12)
Furthermore, a direct evaluation of the integral in (4) for small $`\tau `$ yields $`A_c=\mathrm{}_\nu +\mathrm{o}(\tau )`$ with

$$\mathrm{}_\nu =\mathrm{log}(|\tau |/4)+\psi (\nu /2)/2+\psi (\nu /2)/2\psi (1/2),$$
(4.13)
and this initial condition together with (4.11) completely determines $`A_c.`$ Since (4.11) can be recognized as the hypergeometric differential equation in the variable $`z=\tau ^2,`$ we can immediately write the solution as

$`A_c=`$ $`\mathrm{}_\nu F(\nu /2,\nu /2;1;\tau ^2)+{\displaystyle \underset{k=1}{\overset{\mathrm{}}{}}}(\nu /2)_k(\nu /2)_k/(k!)^2(\tau ^2)^k`$
$`\times (\psi (k+1)\psi (1)+\psi (\nu /2)/2\psi (k+\nu /2)/2+\psi (\nu /2)/2\psi (k\nu /2)/2).`$ (4.14)
On Taylor expansion we now obtain

$`\pi a_n/2`$ $`=\nu ^1+\mathrm{}_\nu \tau +[\nu (\mathrm{}_\nu 1/2)+\nu ^1]\tau ^2/2+\nu ^2(\mathrm{}_\nu 1)\tau ^3/4`$ (4.15)
$`+[\nu ^3(\mathrm{}_\nu 5/4)\nu 2\nu ^1]\tau ^4/16+[\nu ^4(\mathrm{}_\nu 3/2)4\nu ^2(\mathrm{}_\nu 1/2)]\tau ^5/64+\mathrm{O}(\tau ^6)`$
which can be extended as required.

From the expression (4.15) one can conclude that $`C(N,N),`$ as an $`N\times N`$ determinant, will contain logarithmic terms $`\tau ^q(\mathrm{log}|\tau |)^p`$ with $`qp`$ and $`pN`$—barring cancellation. However, there is cancellation and the key qualitative features of the cancellation can be deduced simply by using an alternative representation for the determinant. In particular, by systematically subtracting rows and columns one can show that an equivalent determinant has matrix elements $`a_{i,j}`$ which are of the same integral form as the $`a_n=a_{i,j}`$ in eq. (4.8) except for the replacement of the exponential $`\mathrm{exp}(2i\theta n)`$ by the product $`\mathrm{exp}(i\theta n)(2\mathrm{sin}(\theta ))^{i+j2}.`$ The new matrix is no longer of Toeplitz form but for our purposes here it is the better representation because of the powers of $`\mathrm{sin}(\theta )`$ which have the effect of shifting the $`\mathrm{log}|\tau |`$ singularity of the integrals to higher order in $`\tau .`$ Indeed one can show from the new integrals that the leading singular behavior of each matrix element $`a_{i,j}`$ is proportional to $`\tau ^{i+j1}\mathrm{log}|\tau |`$ and this is sufficient to show that in the logarithmic terms $`\tau ^q(\mathrm{log}|\tau |)^p`$ in $`C(N,N)`$ one must have $`qp^2.`$ We have not attempted to pursue this argument to deduce $`C(N,N)`$ analytically but rather have resorted to a numerical small $`N`$ evaluation of $`C(N,N)`$ using (4.15) and then fitting to obtain formulae valid for general $`N.`$ Our results for the leading logarithm term are summarized by the expression

$`\sqrt{s}C(N,N,\tau )`$ $`=`$ $`C(N,N,\tau =0){\displaystyle \underset{p=0}{\overset{\mathrm{}}{}}}4^p(\mathrm{log}|\tau |+_N)^p(N\tau /4)^{p^2}`$
$`\{{\displaystyle \underset{k=1}{\overset{p1}{}}}(N^2k^2)^{pk}\}\{1+(1+2(N^2p^2))\tau ^2/8+\mathrm{O}(\tau ^3)\}`$
in which we have used $``$ to denote the discrete logarithm, i.e.

$$_N=\psi (N+1)/2+\psi (N)/2\psi (1)\mathrm{log}(4).$$
(4.17)
The product factor in the first braces of (4) is to be understood as unity for $`p<2`$ and in the final brace pair the coefficients of $`\tau ^q`$ are polynomials in $`N`$ of degree $`q.`$ Furthermore, with the $`\sqrt{s}`$ factor extracted explicitly as in (4) these coefficients of $`\tau ^q`$ vanish if $`q=2k+1`$ with $`k<p.`$ The critical correlation factor in (4) is

$`C(N,N,\tau =0)`$ $`=`$ $`{\displaystyle \underset{n=1}{\overset{N}{}}}\mathrm{\Gamma }^2(n)/(\mathrm{\Gamma }(n+1/2)\mathrm{\Gamma }(n1/2)),`$ (4.18)
and approaches $`A/N^{1/4}`$ as $`N\mathrm{}`$ where  $`\mathrm{log}(A)=3\zeta ^{}(1)+\mathrm{log}(2)/12`$. We have used (4), extended or truncated to some order in $`\tau `$, as input to the quadratic recursion formulae to generate all correlation products $`\sqrt{s}C(m,n)`$ as series in $`\tau `$ within an octant $`mn0,`$ $`m+n2N+1.`$ While most of our results are numerical, they are consistent with the assumption that the structure observed on the diagonal remains true on the octant. That is, if we define $`n=\mu N,`$ $`0\mu 1`$ and set $`m+n=2N`$ (even shell) or $`m+n=2N+1`$ (odd shell) then on these even/odd shells for large $`N,`$ the correlations are of the form

$$\sqrt{s}C_{\mathrm{e}/\mathrm{o}}(\mu )=N^{1/4}\underset{p=0}{\overset{\mathrm{}}{}}(\mathrm{log}|\tau |+_N)^p(N\tau /4)^{p^2}A_{\mathrm{e}/\mathrm{o}}^{(p)}(\mu )$$
(4.19)
where the $`A_{\mathrm{e}/\mathrm{o}}^{(p)}(\mu )`$ are still Taylor series in $`\tau `$ as in (4) but with coefficients that are now (possibly asymptotic) Laurent series in $`1/N`$ rather than polynomial in $`N.`$ The highest power of $`N`$ multiplying $`\tau ^q`$ remains $`N^q.`$

Considerable care must be exercised in the numerical work since the recursion formulae are unstable. A toy recursion relation of structure similar to the ones we use in the Ising study illustrates this nicely. Let

$$C_{m+1}^n=(2(C_m^n)^2C_m^{n1}C_m^{n+1})/C_{m1}^n$$
(4.20)
which is to be applied to all possible $`n`$ and iterated forward in $`m.`$ The recursion (4.20) has as a solution a constant, say $`M,`$ but is susceptible to a steady state growth of errors so that

$$C_m^nM+ϵ(1)^n\alpha ^m$$
(4.21)
is also a possible solution. On substituting (4.21) into (4.20) one finds the growth constant $`\alpha `$ must satisfy $`\alpha ^26\alpha +1=0`$ so that

$$\alpha =3+\sqrt{8},\text{ }\mathrm{log}_{10}(\alpha )=0.765\mathrm{}.$$
(4.22)
This is very close to what we find in our numerical work and implies that if we want some number of digits $`D`$ that are accurate at the outer edge of the octant on a shell specified by $`n+m=2N+1`$ we need to start with digits $`D_0D+1.53N.`$ In practice we have worked as high as $`N=146`$ with $`D_0=380`$ using the automatic multiple precision facility of Maple. In section 6.1 we use these results to calculate the “short-distance” (including analytic background) terms in the susceptibility.

## 5 Conjectured short-distance structure

In this section we state some conjectures for the short distance behavior of the $`\lambda 1`$ model introduced in section 4. We arrived at these conjectures by inspection of series obtained from a combination of the integrals (4.2) and the difference equations (4)–(4). Our interest in this model is motivated by several considerations. Firstly, it enables us to see how the analytic structure of $`\chi `$ evolves as successive contributions $`\widehat{\chi }^{(2j)}`$ or $`\widehat{\chi }^{(2j+1)}`$ are added. Secondly, the correlation function $`C(M,N;\lambda )`$ for small values of $`M`$ and $`N`$ may be required as initial conditions for certain series generation algorithms. Thirdly, we hope that the presence of an additional parameter which can be varied will provide some insight into the Ising model itself. Finally, the deformations of the elliptic functions that appear in our conjectures may be of intrinsic mathematical interest.

Following the lead of ref.  we make a change of variable from the modulus, $`k`$, to the nome $`q`$. Examination of the $`j`$-particle contributions to the isotropic susceptibility and two-point functions reveals that there is much regular structure. We arrive at exact conjectures for $`\widehat{C}^{(j)}(0,0)`$, $`\widehat{C}^{(j)}(1,0)`$, $`\widehat{C}^{(j)}(1,1)`$, $`\widehat{C}^{(j)}(2,0)`$ and $`\widehat{C}^{(j)}(2,1)`$ as functions of $`j`$ and $`q`$. These provide the first terms in the short distance expansion of the susceptibility.

### 5.1 $`q`$-series in the Ising susceptibility

Although more complete results have been obtained for the correlation functions, we will demonstrate the method by which we derived our conjectures using the susceptibility as an example. To make our observations more concrete, we reproduce tables of series coefficients for $`(1k^2)^{1/4}\widehat{\chi }_{\text{iso}}^{(j)}`$. As an example of how to interpret the table, we read $`(1k^2)^{1/4}\widehat{\chi }_{\text{iso}}^{(3)}=8(k^{9/2}/2^9+4k^6/2^{12}+16k^{13/2}/2^{13}+4k^7/2^{14}+20k^{15/2}/2^{15}+84k^8/2^{16}+\mathrm{}).`$

Inspecting Table 1, we make the conjecture that the series can be sensibly decomposed as

$$\begin{array}{c}2^j(\sqrt{k}/2)^{j^2}[(1+c_{0,1}k^2+c_{0,2}k^4+\mathrm{})\hfill \\ \hfill +(\sqrt{k}/2)^j(4+k+c_{1,2}k^2+c_{1,3}k^3+\mathrm{})\\ \hfill +(\sqrt{k}/2)^{2j}(c_{2,0}+c_{2,1}k+c_{2,2}k^2+\mathrm{})+\mathrm{}]\end{array}$$
(5.1)
In fact, we observe that as $`j`$ tends to larger and larger values, the first row of coefficients tends towards the expansion in $`k`$ of $`2^jq^{j^2/4}/\theta _3(0,q)`$.888This may be checked by reverting the product form given in (2.6) for $`k=k(q)`$ and using the identity $`\theta _3(0,q)=\sqrt{2\mathrm{K}/\pi }`$ to obtain the expansion of $`\theta _3(0,q)`$ in $`k`$. The first correction comes in at order $`q^{j(j+1)/4}`$. Hence we make the change of variable from $`k`$ to $`q`$ in $`\widehat{\chi }_{\text{iso}}^{(j)}`$ and divide the result by $`2^jq^{j^2/4}/\theta _3(0,q)`$. For all $`j`$ we obtain a series of the form $`1+4q^{j/4}+cq^{(j+1)/4}+\mathrm{}`$ The sequence of terms starting at order $`q^{j/4}`$ appears again to be fitted by a recognizable function of $`q`$, at least until contributions appear at order $`q^{2j/4}`$ or $`q^{3j/4}`$. Subtracting this assumed product form yields a series whose first correction enters at order $`q^{2j/4}`$. The coefficients of the terms of orders lying between $`q^{2j/4}`$ and $`q^{3j/4}`$ are not independent of $`j`$ as was the case previously, but vary linearly with $`j`$. The constant part has a product form, and the $`j`$-dependent part may as well be we do not have sufficiently many terms to make a firm conjecture. Likewise the terms between $`q^{3j/4}`$ and $`q^{4j/4}`$ depend quadratically on $`j`$, and the correction at $`q^{4j/4}`$ appears to vary as the fourth power of $`j`$.

From the long series we have produced, we have been able to conjecture that

$$\begin{array}{cc}& 2^j(1k^2)^{1/4}\widehat{\chi }_{\text{iso}}^{(j)}=\frac{q^{j^2/4}}{\theta _3(0,q)}[1+4q^{j/4}\underset{n=1}{\overset{\mathrm{}}{}}\frac{1+q^{n1/2}}{1+q^n}\hfill \\ & +4q^{2j/4}\left(\underset{n=1}{\overset{\mathrm{}}{}}(1q^{2n1})^4(1+q^{n1/2})^4+(j+1)4q^{1/4}/\theta _2^2(0,q^{1/2})\right)\hfill \\ & +4q^{3j/4}(f_3^{(0)}(q)+(j+1)(j+2)f_3^{(2)}(q))+O\left[q^{4j/4}\right]]\hfill \end{array}$$
(5.2)
where

$`f_3^{(0)}`$ $`=1+6q^{1/2}+26q+17q^{3/2}81q^255q^{5/2}+285q^3\mathrm{}`$ (5.3)
$`f_3^{(2)}`$ $`=1+q^{1/2}5q4q^{3/2}+15q^2+10q^{5/2}39q^3+\mathrm{}`$ (5.4)
These are consistent with the expressions for the correlation functions in the following section. That is to say, summing the correlation functions given in the following subsection gives terms that agree, as far as they should, with the above expression. Similarly, summing (5.2) over $`j`$ gives a series that agrees to the appropriate (low) order with the known expansion for $`\chi _{iso.}`$

### 5.2 $`q`$-series in the two-point functions

In this subsection we write down some conjectured results for the short-distance correlation functions which we have derived empirically. For these cases, unlike the susceptibility, the general correction term is apparent from the series and we are able to formulate exact conjectures. Let us define operators, $`\mathrm{\Phi }_0`$ and $`\mathrm{\Phi }_1`$, which convert power series in $`z`$ to power series in $`q`$ according to the rules

$`\mathrm{\Phi }_0{\displaystyle \underset{n=0}{\overset{\mathrm{}}{}}}c_nz^n`$ $`={\displaystyle \underset{n=0}{\overset{\mathrm{}}{}}}c_nq^{n^2/4},`$ (5.5)
$`\mathrm{\Phi }_1{\displaystyle \underset{n=0}{\overset{\mathrm{}}{}}}c_nz^n`$ $`={\displaystyle \underset{n=0}{\overset{\mathrm{}}{}}}c_nq^{n(n+1)/4}.`$ (5.6)
Then we conjecture that

$`2^j(1k^2)^{1/4}\widehat{C}^{(j)}(0,0)={\displaystyle \frac{1}{\theta _3(0,q)}}\mathrm{\Phi }_0{\displaystyle \frac{z^j(1z^2)}{(1+z^2)^{j+1}}},`$ (5.7)
$`2^j(1k^2)^{1/4}\widehat{C}^{(j)}(1,0)={\displaystyle \frac{1}{\theta _3(0,q)}}\left({\displaystyle \frac{2q^{1/8}\theta _3(0,q^{1/2})}{\theta _2(0,q^{1/2})}}\right)^{1/2}\mathrm{\Phi }_1{\displaystyle \frac{z^j(1z)}{(1+z^2)^{j+1}}},`$ (5.8)
$`2^j(1k^2)^{1/4}\widehat{C}^{(j)}(1,1)={\displaystyle \frac{2(j+1)}{\theta _2(0,q)\theta _3^2(0,q)}}\mathrm{\Phi }_0{\displaystyle \frac{z^{j+1}(1z^2)}{(1+z^2)^{j+2}}},`$ (5.9)
$`2^j(1k^2)^{1/4}\widehat{C}^{(j)}(2,0)={\displaystyle \frac{1}{q^{1/4}\theta _3(0,q)}}\left({\displaystyle \frac{2q^{1/8}\theta _3(0,q^{1/2})}{\theta _2(0,q^{1/2})}}\right)^{4/2}\mathrm{\Phi }_0{\displaystyle \frac{z^{j+1}}{(1+z^2)^{j+1}(1z^2)}}`$
$`{\displaystyle \frac{16}{\theta _3(0,q)\theta _2^4(0,q^{1/2})}}q^{1/4}{\displaystyle \frac{d}{dq^{1/4}}}\left[\mathrm{\Phi }_0{\displaystyle \frac{z^{j+2}}{(1+z^2)^{j+2}(1z^2)}}\right]4[{\displaystyle \frac{1}{\theta _3(0,q)}}`$
$`+{\displaystyle \frac{\theta _2(0,q)}{2\theta _3^2(0,q)\theta _2(0,q^4)}}{\displaystyle \frac{8}{\theta _3(0,q)\theta _2^4(0,q^{1/2})\theta _2(0,q^4)}}q^{1/4}{\displaystyle \frac{d}{dq^{1/4}}}\left[\mathrm{\Phi }_0{\displaystyle \frac{z^2}{(1+z^2)^2(1z^2)}}\right]]`$
$`\times \mathrm{\Phi }_0{\displaystyle \frac{z^{j+2}}{(1+z^2)^{j+2}(1z^2)}},`$ (5.10)
and
$`2^j(1k^2)^{1/4}\widehat{C}^{(j)}(2,1)={\displaystyle \frac{1}{\theta _3(0,q)}}(8q^{5/16}{\displaystyle \frac{\theta _2(0,q^{1/4})}{\theta _2^5(0,q^{1/2})}}{\displaystyle \frac{d}{dq^{1/4}}}\left[\mathrm{\Phi }_1{\displaystyle \frac{z^{j+1}(1+z)}{(1+z^2)^{j+1}(1z^2)}}\right]`$
$`+\left(\mathrm{\Phi }_1{\displaystyle \frac{z(1+z)}{1z^4}}\right)^1\mathrm{\Phi }_1{\displaystyle \frac{z^{j+1}(1+z)}{(1+z^2)^{j+1}(1z^2)}}`$
$`\times [\theta _4(0,q)8q^{5/16}{\displaystyle \frac{\theta _2(0,q^{1/4})}{\theta _2^5(0,q^{1/2})}}{\displaystyle \frac{d}{dq^{1/4}}}\left[\mathrm{\Phi }_1{\displaystyle \frac{z(1+z)}{1z^4}}\right]])`$ (5.11)
As noted above, these are necessary, but not sufficient, for the generation of the $`j`$-particle contributions to the correlation functions, being some of the initial conditions for the recurrences.

## 6 Scaling form of the susceptibility.

The main result in this section is a conjecture completely specifying the analytic structure of the susceptibility $`\chi `$ of the isotropic Ising model in the vicinity of the critical point both in the neighborhood of the ferromagnetic point $`s=1`$ and the anti-ferromagnetic point $`s=1`$. The conjecture, contained in eqs. (1.151.18) is based on what we believe is overwhelming numerical evidence that is obtained by disentangling the “short- distance” and “scaling” parts of $`\chi `$ in a manner described below.

In section 6.1 we give the assumptions and numerical procedures we use to derive the “short-distance” part of $`\chi `$ in (1.15) with the coefficients listed in the Appendix. Then, in section 6.2 we describe the “short-distance” subtraction and analysis on which the behavior of the scaling-amplitude function $`F_\pm `$ shown in (1.18) is based. Finally, in section 6.3 we outline our fitting programs to determine the coefficients in the functions $`F_\pm `$. Since there are no confluent singularities in $`F_\pm `$ whose amplitudes need to be found, the fitting procedure is very well-conditioned and the coefficients in the Appendix are as determined to an accuracy of up to 20 digits.

### 6.1 “Short-distance” term.

That the “short-distance” contribution to $`\chi `$ can be obtained from numerical values of $`C(m,n)`$ for small $`|m|`$ and $`|n|`$ relies on certain assumptions about the behavior of the expansion coefficients in eq. (4.19). In particular, we assume that (4.19) remains valid up to $`N`$ of the order $`1/\tau `$ where it can, in principle, be matched term by term to a large distance expansion that properly describes the roughly exponential $`\mathrm{exp}(N\tau )`$ decay of correlations as $`N\mathrm{}.`$ Explicit matching formed the basis of the previous calculations of terms in the “short-distance” $`\chi `$ (cf. ) but this becomes extremely cumbersome at higher order. Our ability to go to high order here rests on the fact that we dispense with such explicit matching and rely instead on power counting to uniquely identify those terms that contribute to the “scaling” and the “short-distance” parts of $`\chi `$ separately. We believe this is tantamount to the scaling argument that in the critical region there is a single length scale proportional to $`1/\tau ^\nu `$ with $`\nu =1`$ and thus, in a way to be made more precise below, we can deduce the power law of the large distance contribution of any set of terms varying as $`N^p`$ by simply replacing $`N^p`$ with $`1/\tau ^p.`$ Terms whose variation is as a fractional power of $`\tau `$ (with possibly logarithmic multipliers) are discarded as assumed contributions to the “scaling” part of $`\chi .`$ Terms whose variation is predicted to be an integer power of $`\tau `$ (with possibly logarithmic multipliers) are assumed to be part of the “short-distance” $`\chi `$ and are treated more carefully.

To make the argument and assumptions more explicit we begin with some definitions that will be useful also for the subsequent analysis. Let the two dimensional sum (1.8) defining $`\chi `$ be reduced to a one dimensional sum by combining the contributions from all sites on the even and odd squares $`|m|+|n|=2N`$, $`2N+1`$ and then further combining these into sum and difference combinations which are necessary for separating the ferromagnetic and anti-ferromagnetic contributions. That is we write

$$\sqrt{s}C_{N_\pm }=\underset{k}{}C_{N_\pm }^{(k)}\tau ^k=\sqrt{s}(C_\mathrm{e}(\mu )\pm C_\mathrm{o}(\mu ))$$
(6.1)
where the first sum simply defines the coefficients $`C_{N_\pm }^{(k)}`$ while the second is the actual lattice sum over the octant correlations defined by (4.19) and extended by symmetry over the full square. The coefficients $`A_{N_\pm }^{(l,k)}`$ in the expansion

$$C_{N_\pm }^{(k)}=\underset{l=0}{\overset{L}{}}A_{N_\pm }^{(l,k)}(\mathrm{log}|\tau |+_N)^l$$
(6.2)
have the large $`N`$ asymptotic expansions

$$A_{N,\mathrm{f}}^{(l,k)}=\underset{p=0}{}A_\mathrm{f}^{(l,k,p)}N^{3/4+kp},\text{ }A_{N,\mathrm{af}}^{(l,k)}=\underset{p=0}{}A_{\mathrm{af}}^{(l,k,p)}N^{1/4+kp},$$
(6.3)
as assumed in (4.19) based on numerical evidence. The ferromagnetic and anti-ferromagnetic cases have been treated separately in (6.3) to emphasize that shell subtraction in (6.1) reduces the leading power of $`N`$ by unity. The upper limit $`L`$ on the logarithmic powers in (6.2) depends on $`k`$ as discussed in connection with eq. (4) but its precise value is not needed in the following discussion. We also introduce the partial sums

$`S_{N_\pm }^{(k)}={\displaystyle \underset{n=0}{\overset{N}{}}}C_{n_\pm }^{(k)}={\displaystyle \underset{l=0}{\overset{L}{}}}R_{N_\pm }^{(l,k)}(\mathrm{log}|\tau |)^l`$ $`=`$ $`{\displaystyle \underset{l=0}{\overset{L}{}}}(\beta _\pm ^{(l,k)}+\delta R_{N_\pm }^{(l,k)})(\mathrm{log}|\tau |)^l`$ (6.4)
where the $`R_{N_\pm }^{(l,k)}`$ in the second equality in (6.4) are numerical coefficients that are directly generated by the quadratic recursion relations and the subsequent lattice summations. In the final equality in (6.4) these expansion coefficients have been formally split into two terms. This formal separation is to be understood in an asymptotic sense with $`\beta _\pm ^{(l,k)}`$ a constant that is the $`N`$ independent part of $`R_{N_\pm }^{(l,k)}`$ as $`N\mathrm{}`$ and which we can understand as an “integration constant”. The remainder $`\delta R_{N_\pm }^{(l,k)}`$, we will show in the case that $`l=L`$, has an asymptotic expansion like the $`A_{N_\pm }^{(l,k)}`$ in (6.2) except for an extra factor of $`N`$ from the summation. For $`lL`$ the $`\delta R_{N_\pm }^{(l,k)}`$ can be expressed as sums of such expansions with multiplying logarithmic factors $`_N.`$

We now come to our key assumption. The partial sums in (6.4) in most cases diverge as $`N\mathrm{}`$ because of the presence of large fractional powers of $`N`$ in the $`\delta R_{N_\pm }^{(l,k)}`$. However, we assume that if the formal matching of the short and large distance expansions for $`C(m,n)`$ had been carried out, these partial sums would in fact converge and furthermore the contribution of each power of $`N`$ term could be estimated by the replacement $`N1/\tau .`$ Since all the powers of $`N`$ in $`\delta R_{N_\pm }^{(l,k)}`$ are fractional, the result is that these terms contribute only to the “scaling” part of $`\chi .`$ Thus the “short-distance” part of $`\chi `$ comes entirely from the “integration constant” term $`\beta _\pm ^{(l,k)}`$ in (6.4), and on explicitly reinserting the factors of $`\sqrt{s}`$ and $`\tau ^k`$ we get

$$B_{\mathrm{f}/\mathrm{af}}=\underset{k}{}\underset{l=0}{\overset{L}{}}\beta _\pm ^{(l,k)}\tau ^k(\mathrm{log}|\tau |)^l/\sqrt{s}=\underset{p=0}{\overset{\mathrm{}}{}}\underset{q}{}b_{\mathrm{f}/\mathrm{af}}^{(p,q)}\tau ^q(\mathrm{log}|\tau |)^p$$
(6.5)
at the ferromagnetic and anti-ferromagnetic points. The last equality in (6.5) is eq. (1.15) and follows simply by expanding $`\sqrt{s}`$ in a series in $`\tau .`$ The rest of this section comprises a discussion of the numerical scheme we use to isolate the “integration constants” efficiently.

A technical problem arises in that one must somehow isolate the constant $`\beta _\pm ^{(l,k)}`$ from an $`N`$-dependent sequence $`R_{N_\pm }^{(l,k)}`$ and this can be a very unstable procedure if the $`R_{N_\pm }^{(l,k)}`$ are combinations involving logarithms. Fortunately there are no logarithms in $`R_{N_\pm }^{(L,k)}`$ and one can iteratively remove the logarithmic factors in the remaining $`R_{N_\pm }^{(l,k)}`$ by working in sequence from $`l=L,L1,\mathrm{}`$ to $`l=0`$. To understand precisely how to do the subtractions we first relate the unknown $`N`$dependent structure of the $`R_{N_\pm }^{(l,k)}`$ to the known, or rather assumed, simple structure of the $`A_{N_\pm }^{(l,k)}.`$ From the definition of $`S_{N_\pm }^{(k)}`$ as a partial sum it follows that $`S_{N_\pm }^{(k)}S_{N1_\pm }^{(k)}=C_{N_\pm }^{(k)}`$ and by comparing coefficients in eqs. (6.4) and (6.2) we get

$$R_{N_\pm }^{(l,k)}R_{N1_\pm }^{(l,k)}=\underset{m=l}{\overset{L}{}}\left(\genfrac{}{}{0pt}{}{m}{l}\right)A_{N_\pm }^{(m,k)}(_N)^{ml}.$$
(6.6)
In the case that $`l=L`$ there are no logarithmic factors on the right hand side of (6.6) and the equation is just a discrete first order differential equation in $`N`$ whose solution is an integration constant additive to an asymptotic series of the same form as in (6.3) except for an extra power of $`N.`$ The easiest way to determine the “integration constant” $`\beta _\pm ^{(L,k)}`$ is to fit it together with unknown coefficients defining the asymptotic series part to a sequence of $`R_{N_\pm }^{(L,k)}`$. This is a stable numerical procedure, although if $`k`$ is large, $`\beta _\pm ^{(L,k)}`$ is sub-dominant to many much larger terms and it is essential to calculate in high precision. For example, we might start the recursive calculation of $`R_{N_\pm }^{(L,k)}`$ as described in the previous section with $`D_0=380`$ digits, end with $`D=155`$ digits at $`N=146,`$ and continue with this number of digits in a $`71\times 71`$ matrix inversion to obtain the $`\beta _\pm ^{(L,k)}`$ of interest while discarding the remaining 70 coefficients of the asymptotic series! The results are slightly more accurate than the final answers that are given in the Appendix.

For $`lL`$ the presence of logarithmic factors on the right-hand side of (6.6) makes fitting to $`R_{N_\pm }^{(l,k)}`$ impractical. Instead we define subtracted functions

$$F_{N_\pm }^{(l,k)}=R_{N_\pm }^{(l,k)}+\underset{m=l+1}{\overset{L}{}}\left(\genfrac{}{}{0pt}{}{m}{l}\right)(R_{N_\pm }^{(m,k)}\beta _\pm ^{(m,k)})(_N)^{ml}$$
(6.7)
which we use as replacements for $`R_{N_\pm }^{(l,k)}`$ in all fitting procedures. That the fitting functions $`F_{N_\pm }^{(l,k)}`$ will give the same $`\beta _\pm ^{(l,k)}`$ is obvious since the added terms in (6.7) all contain logarithmic factors and are thus not $`N`$-independent. That the $`F_{N_\pm }^{(l,k)}`$ are logarithm free follows from

$$F_{N_\pm }^{(l,k)}F_{N1_\pm }^{(l,k)}=A_{N_\pm }^{(l,k)}+\underset{m=l+1}{\overset{L}{}}\left(\genfrac{}{}{0pt}{}{m}{l}\right)(F_{N1_\pm }^{(m,k)}\beta _\pm ^{(m,k)})(_{N1}_N)^{ml}$$
(6.8)
and an argument by induction starting from $`l=L`$. Eq. (6.8) can be verified using the definition (6.7) and the subtraction eq. (6.6).

To summarize, we obtain $`\beta _\pm ^{(l,k)}`$ numerically by using (6.7) iteratively, starting with $`l=L1,`$ to generate $`F_{N_\pm }^{(l,k)}`$ using procedures similar to that described above for the initial $`F_{N_\pm }^{(L,k)}`$ which are just the $`R_{N_\pm }^{(L,k)}.`$ The results are given in the Appendix.

Our arguments above are at best plausible and to provide a rigorous proof justifying our numerical procedure we believe one would have to do three things. First, one would have to show the $`C(m,n)`$ expansion has the form (4.19) we deduced on numerical grounds. Second, one would have to show there exists a corresponding asymptotic expansion valid at large $`N.`$ Third, one must show both expansions have a sufficiently large domain of validity that the matching we assumed could in fact be carried out to arbitrarily high order. On the other hand, we have numerical evidence for the validity of our procedure as described in the next section. The “short-distance” terms we have calculated, when used in a series subtraction process, leave a residual in high order series coefficients that is consistent with the complete elimination of all “short-distance” terms, both singular and analytic, to the O$`(\tau ^{14})`$ we have worked. While the main intent of section 6.2 is actually quite different, it is highly unlikely that the cancellations necessary to yield the expansion of the scaling-amplitude function $`F_\pm `$ of eq. (1.18) in pure integer powers of $`\tau `$ would have occurred had there been an error, numerical or otherwise, in the results of this section.

We conclude this section with a toy model example to illustrate the possible convergence properties of the short-distance expansion. The main impediment to extending the calculation described above is that we do not have simple analytical expressions for the $`\tau `$ expansion of the $`C(M,N)`$ in general. An exception is on the diagonal where we know, cf. eq. (4), that the coefficient of $`\tau ^{p^2}(\mathrm{log}|\tau |)^p`$ is

$$C_N^{(p)}=4^{pp^2}C(N,N,\tau =0)N^p\underset{k=1}{\overset{p1}{}}(1N^2/k^2)^{pk}$$
(6.9)
and a simple expression for the asymptotic expansion of $`C(N,N,\tau =0)`$ at large $`N`$ can be found in Au-Yang and Perk . We can define the diagonal partial sums $`S_N^{(p)}=_{n=0}^NC_n^{(p)}`$ and, as in the analysis described above, ask for the short-distance coefficient $`b_{\text{diag}}^{(p)}`$ which is the $`N`$ independent part of $`S_N^{(p)}`$ in the limit $`N\mathrm{}.`$ This term $`b_{\text{diag}}^{(p)}`$ could be viewed as a partial contribution to the short-distance terms of interest but since we do not know what cancellations will occur when we include all $`C(M,N),`$ we prefer to consider it only as a toy result that is suggestive for the convergence of the short-distance terms with order $`p.`$

The term $`b_{\text{diag}}^{(p)}`$ is also equal to the $`ϵ`$ independent part of $`_{n=0}^{\mathrm{}}C_n^{(p)}e^{ϵn}`$ in the limit $`ϵ0.`$ This latter expression is more convenient since we can add to it any term such as $`1/(1e^ϵ)^{p^2+3/4k}`$ for integer $`k`$ without contributing to any $`ϵ`$ independent term. If we expand such terms in series in $`e^ϵ,`$ we obtain as an equivalent sum

$$\underset{n=0}{\overset{\mathrm{}}{}}\left[C_n^{(p)}\underset{k=0}{\overset{K}{}}g_k\mathrm{\Gamma }(n+p^2+3/4k)/n!\right]e^{ϵn}.$$
(6.10)
Now choose the $`g_k`$ in (6.10) such that the divergent terms in the asymptotic $`n`$ expansion of $`C_n^{(p)}`$ cancel the divergent terms in the Gamma functions. Then the $`n`$ sum becomes convergent even with $`ϵ=0.`$ With the cancellation extended to include also some slowly decaying terms in $`n,`$ thus requiring $`K>p^2,`$ one obtains the explicit formula

$$b_{\text{diag}}^{(p)}=\underset{n=0}{\overset{\mathrm{}}{}}\left[C_n^{(p)}\underset{k=0}{\overset{K}{}}g_k\mathrm{\Gamma }(n+p^2+3/4k)/n!\right]$$
(6.11)
which is very convenient for numerical work. We find the $`b_{\text{diag}}^{(p)}`$ calculated from (6.11) for $`p=1,2,\mathrm{},20`$ are

$`1.30`$ $`\times 10^1,`$ $`4.29`$ $`\times 10^4,`$ $`1.58`$ $`\times 10^5,`$ $`9.12`$ $`\times 10^9,`$
$`6.31`$ $`\times 10^{11},`$ $`8.31`$ $`\times 10^{14},`$ $`2.38`$ $`\times 10^{15},`$ $`3.73`$ $`\times 10^{17},`$
$`2.08`$ $`\times 10^{17},`$ $`1.18`$ $`\times 10^{17},`$ $`3.61`$ $`\times 10^{16},`$ $`1.78`$ $`\times 10^{14},`$ (6.12)
$`6.67`$ $`\times 10^{11},`$ $`5.80`$ $`\times 10^7,`$ $`5.15`$ $`\times 10^1,`$ $`1.44`$ $`\times 10^{+6},`$
$`5.34`$ $`\times 10^{+14},`$ $`8.13`$ $`\times 10^{+23},`$ $`2.06`$ $`\times 10^{+35},`$ $`2.70`$ $`\times 10^{+47}.`$
The sequence in (6.1) clearly shows the asymptotic nature of the toy expansion. Furthermore the magnitude of the terms is in semi-quantitative agreement with $`\mathrm{\Gamma }(p^2/2)/a^{p^2/2}`$ from eq. (3.30). Thus while we cannot conclude that this will be how the short-distance susceptibility terms $`B_{\mathrm{f}/\mathrm{af}}`$ in (1.15) will behave, it is at least encouraging to note that there is no evidence for any behavior more singular than that predicted by our natural boundary analysis.

### 6.2  Series “proof” of $`F_\pm `$ behavior

Any numerical analysis program to deduce functions such as the scaling-amplitude functions $`F_\pm `$ from their series expansions is in essence a fitting routine and always presupposes a knowledge of the analytic structure of the result. Although we have very little exact knowledge of $`F_\pm `$ we do know that the multiplier of $`|\tau |^{7/4}`$ in $`\chi ^{(2)}`$ (cf. eq. (1.13)) contains $`\mathrm{log}|\tau |`$ terms and there is numerical evidence  that this is true of the $`\chi ^{(n)}`$ for $`n>2`$ also. On the other hand, Gartenhaus and McCullough  found that the series for $`\chi `$ were consistent with the absence of logarithmic corrections in $`F_+`$ through order $`\tau ^3`$ and this was subsequently confirmed for $`F_{}`$ as well . Scaling arguments on the question of logarithmic corrections are necessarily inconclusive because of the lack of information on amplitudes, some of which may vanish. We know of no scaling argument that either definitely requires the presence of logarithmic terms or can definitely exclude them.

The natural boundary arguments given in section 3.2 preclude the possibility that $`\chi `$ has a convergent rather than asymptotic expansion about $`\tau =0`$, but we know of no analytical argument that shows whether this applies to both the “scaling” and “short-distance” terms or just to one or the other.999A simple ratio analysis of the available terms in the expansion of $`\sqrt{s}F_\pm `$ in the variable $`\tau ^2`$ (see appendix) is indicative of a convergent series with radius of convergence $`1.2`$. This implies a conjugate pair of singularities at $`\tau 1.1i`$. If however the series is asymptotic, this observation will fail to hold. Neither could we distinguish the convergent from the asymptotic expansions from the data we have. Although there are singularities in $`\chi `$ dense on the line $`i\tau i`$, the singularities close to $`\tau =0`$ are extremely weak and we have not detected any singularity closer to the ferromagnetic $`\tau =0`$ point than that arising out of $`\widehat{\chi }^{(6)}`$ at $`\mathrm{}(\tau )=\sqrt{7}/40.66.`$ Thus for all practical (numerical) purposes we can ignore the possible asymptotic nature of the “scaling” terms.

The possibility that the $`F_\pm `$ might be expanded in a series in integer powers of $`\tau `$ can be confirmed numerically without evaluating any of the coefficients in the expansion. The trick is to generate the series for the susceptibilities divided by the factor $`(1k^2)^{1/4}`$. Any term in $`F_+`$ or $`F_{}`$ that is not a pure integer power law $`\tau ^p`$, $`p>1`$, will necessarily contribute to the high order coefficients of the series. Of course, the leading two terms $`1+\tau /2`$ in $`F_\pm `$, which occur as poles proportional to $`1/\tau ^2`$ and $`1/\tau `$ in the scaled $`\chi _\pm `$, will also contribute, but the amplitude of these terms is known and so their contribution can be subtracted. Similarly, given the accurate “short-distance” amplitudes in the Appendix, the contribution of all “short-distance” terms can be similarly eliminated through order $`\tau ^{14}.`$ We find that the remaining high order series coefficients are plausibly consistent with the “short-distance” $`\mathrm{O}(\tau ^{15})`$ that has not been subtracted and thus that there is no numerical evidence for any powers of $`\tau `$ other than pure integers in the scaling-amplitude functions $`F_\pm `$. The rest of this section gives details of the analysis that is the basis for this conclusion while the following section reports on our procedures for estimating the coefficients in the $`\tau `$ expansions of $`F_\pm `$.

As outlined above, our search for possible terms other than those with pure integer powers of $`\tau `$ in the $`F_\pm `$ involves the observation of the high order series terms in the scaled and pole subtracted susceptibility functions

$$\mathrm{\Delta }\chi _+=\beta ^1\chi _+/(1s^4)^{1/4}(2K_c\sqrt{2})^{7/4}C_0^+2^{1/2}/(1s)^2,$$
(6.13)
$$\mathrm{\Delta }\chi _{}=\beta ^1\chi _{}/(1s^4)^{1/4}(2K_c\sqrt{2})^{7/4}C_0^{}2\sqrt{2}s^2/(s^21)^2=\underset{m=1}{\overset{\mathrm{}}{}}K_{2m}^{}s^{2m}$$

applicable for $`T>T_c`$ and $`T<T_c`$ respectively. Since the procedures for $`T>T_c`$ and $`T<T_c`$ are not different in principle, we will restrict the discussion below to the $`T<T_c`$ case only.

To implement the “short-distance” subtraction we can restrict ourselves to determining $`b_n`$, the contribution to the coefficient of $`s^n`$ arising from the $`s^1=1`$ singularity in $`(1s^4)^{1/4}(\mathrm{log}(\tau ))^p/\sqrt{s}.`$101010The denominator $`\sqrt{s}`$ factor here and in the following is included for the convenience of allowing us to use the “short-distance” amplitudes exactly as tabulated in the Appendix. It also simplifies the final formula (6.14). The contribution from $`\tau ^q(1s^4)^{1/4}(\mathrm{log}(\tau ))^p/\sqrt{s}`$ follows trivially from $`b_n`$ by $`q`$ repetitions of the derivative operation $`D_\tau b_n=(b_{n+1}b_{n1})/2`$ because of the simple form $`\tau =(s^1s)/2.`$ Note also that the contribution from the $`1/s=1`$ singularity is identical except for an overall sign $`(1)^n`$ and thus simply requires that in the end we set all odd power amplitudes to zero and double the even ones. To determine $`b_n`$ is a standard exercise in complex variable contour integration; we deform the contour in an $`s^1`$ integral with integrand $`s^{n+1}/(2\pi i)(1s^4)^{1/4}(\mathrm{log}(\tau ))^p/\sqrt{s}`$ to surround the branch-cut $`1s^1<\mathrm{}`$ and in a final step set $`s^1=\mathrm{exp}x.`$ The result is

$$b_n=\pi ^1_0^{\mathrm{}}𝑑x\mathrm{exp}(nx)(2\mathrm{sinh}2x)^{1/4}\mathrm{}\{\mathrm{exp}(i\pi /4)(\mathrm{log}\mathrm{sinh}xi\pi )^p\}$$
(6.14)
and although the integral (6.14) cannot in general be done in closed form, expansions valid asymptotically for large $`n`$ are easy to generate with computer algebra packages such as Maple. Because our series are particularly long these asymptotic expansions are essentially exact and are also very convenient for the subsequent $`D_\tau `$ differentiations.

It is also necessary to subtract and/or smooth out the contributions from complex singularities on the circle $`|s^1|=1.`$ The only significant ones for the present calculation are those at $`s^1=\pm i`$ and $`\pm \mathrm{exp}(\pm i\pi /3).`$ The subtractions of the leading contributions from the latter singularities arising from $`\widehat{\chi }^{(4)}`$ can be obtained directly from eq. (28) in . The result is

$`K_{2m}^{}`$ $`K_{2m}^{}+1536/(5005\pi )(4\sqrt{2}/(2m)!)[3^{1/4}\mathrm{sin}(2m\pi /3+\pi /4)\{(13/2)_{2m}`$
$`+13/4(15/2)_{2m}+1265/1632(17/2)_{2m}5365/384(19/2)_{2m}\}`$
$`3^{1/4}\mathrm{cos}(2m\pi /3+\pi /4)\{5/2(15/2)_{2m}+75/8(17/2)_{2m}`$
$`+53321/20672(19/2)_{2m}\}]`$ (6.15)
where $`(z)_n=\mathrm{\Gamma }(z+n)/\mathrm{\Gamma }(z)`$ is a Pochhammer symbol.

We find that the further smoothing to give the remainder

$$R_n^{}=n^{1/2}D_b^4n^{1/2}(n^4D_a^4)^3n^3K_n^{}$$
(6.16)
is more than adequate. Here as in  $`D_ag_n=(g_{n1}+g_{n+1})/2`$ suppresses the contributions from $`s^1=\pm i`$ while $`D_bg_n=(g_{n2}+g_n+g_{n+2})/3`$ reduces what has not been subtracted by (6.2). Note that the smoothing in eq. (6.16) scales up the amplitude contributions from the singularities at $`s^1=\pm 1`$ by a factor $`n^{15}.`$

We find that as we include higher and higher orders of the “short-distance” terms in the subtraction process the remainder amplitudes (6.16) decrease in a smooth fashion. The amplitude $`R_n^{}`$ we obtain after having subtracted all terms through $`\mathrm{O}(\tau ^{14})`$ listed in the Appendix is, at $`n=600`$, $`4.1\times 10^7`$ compared to $`1.5\times 10^{43}`$ one gets without pole subtraction in (6.13). The residual is very reasonably the amplitude we would expect from the $`\mathrm{O}(\tau ^{15})`$ terms and this has been confirmed by extending the “short-distance” amplitude sequence in the Appendix by two terms using a crude Padé analysis. The result of this additional subtraction is to reduce $`R_n^{}`$ at $`n=600`$ to $`(1\text{ to }1)\times 10^6.`$ Analysis of the high temperature series leads to a similar conclusion and also that the anti-ferromagnetic “short-distance” terms in the Appendix represent the susceptibility at this point completely with nothing left out.

Although our analysis does not “prove” the absence of powers of $`\tau `$ other than pure integers, we can put very stringent bounds on the amplitudes of any possible singular terms. To make this comparison concrete, we suppose either of $`F_\pm `$ contains the singular term $`A_p\tau ^p\mathrm{log}|\tau |`$ and for simplicity assume $`p`$ integer. The contribution of this term to the amplitude of the coefficient of $`s^n`$ for $`T>T_c`$ or $`1/s^n`$ for $`T<T_c`$ in the high/low temperature series expansion of (6.13) will be about $`A_p\mathrm{\Gamma }(p1)/n^p`$ relative to the pole contribution and this is to be compared to the observed amplitude we obtain after subtracting or smoothing away the known singularities as best we can. Our current results are the amplitude bounds

$`|A_p|`$ $`<`$ $`10^{35}300^p/\mathrm{\Gamma }(p1),\text{ }T>T_c,`$ (6.17)
$`|A_p|`$ $`<`$ $`10^{37}600^p/\mathrm{\Gamma }(p1),\text{ }T<T_c,`$ (6.18)
and these bounds essentially exclude any singularity with reasonable amplitude, scaling as $`\tau ^p,`$ for all $`p`$ less than about 15.

### 6.3  $`F_\pm `$ coefficient analysis

The task of determining $`F_\pm `$ numerically is enormously simplified by the a priori knowledge—strictly speaking a conjecture based on the numerical work of the last section—that $`F_\pm `$ has an expansion in integral powers of $`\tau `$ near $`\tau =0.`$ The absence of any confluent terms means that many different analyses can be used efficiently and we report here on two independent calculations that give essentially identical results, thus again confirming the above conjecture. We have not seriously attempted to optimize our analysis to give the most accurate numerical values possible and thus if it ever becomes necessary one could almost certainly improve on our coefficients as given in the Appendix.

One particular feature of the $`F_\pm `$ expansion is worth noting here. When the functions are scaled by $`\sqrt{s}`$ the resulting series are numerically consistent with series in even powers of $`\tau `$ only. We have also noted a similar simplifying role played by $`\sqrt{s}`$ in the “short-distance” terms which ultimately trace back to the Toeplitz determinant giving the diagonal correlations $`C(N,N)`$ and to the quadratic recursion relations for general $`C(M,N)`$.111111See also the discussion of the quadratic recursion relations in Itzykson and Drouffe  where rescaling by $`\sqrt{s}`$ was used to simplify the scaling limit. This result is not entirely unexpected. The singular part of the free energy in zero field is an even function of $`\tau `$ and the magnetization $`=(1s^4)^{1/8}`$ is $`s^{1/4}(\tau )^{1/8}`$ times an even function of $`\tau .`$ Nonlinear scaling field analysis then predicts, in the absence of corrections, that the “scaling” part of the susceptibility is $`s^{1/2}|\tau |^{7/4}`$ times an even function of $`\tau `$. Thus although our results for $`F_\pm `$ are not consistent with the complete absence of correction terms as discussed in section 7 the prediction that $`\sqrt{s}F_\pm `$ is even in $`\tau `$ does appear to be preserved to all orders.

To determine the coefficients in $`F_\pm `$ we return to the unscaled $`\chi _\pm `$ of eqs. (1.9,1.10) so that the terms in $`F_\pm `$ are now singularities of the function and contribute to the high order coefficients in the series. The cases $`T<T_c`$ and $`T>T_c`$ are again similar; for $`T<T_c`$ the unwanted contributions from the “short-distance” part of $`\chi `$ are subtracted as in the $`s`$-plane analysis described in section 6.2. We use essentially the same smoothing except for a $`1/4`$ shift in power necessitated by the difference in the singularity structure generated by the $`\chi `$ rescaling of section 6.2. That is, we replace the remainder eq. (6.16) by

$$R_n^{}=n^{1/2}D_b^4n^{3/4}(n^4D_a^4)^3n^{13/4}K_n^{}$$
(6.19)
and reduce the remainders $`R_n^{}`$ in (6.19) by a least squares fitting to the unknown $`F_{}`$ coefficients in a procedure similar to that described in . Fitting intervals $`\mathrm{\Delta }n>128`$ are typically used, and an FFT of the residuals is very useful as a diagnostic to interpret the observed oscillations in the residuals in terms of $`\chi `$ singularities on the circle $`|s^1|=1`$. Because the highest order terms in $`F_{}`$ are not fixed (i.e. known) unlike the “short-distance” terms, they tend to float and become effective amplitudes that incorporate all the higher order effects including the “short-distance” contributions that have not been subtracted. One technical result of this is that the residuals we observe in (6.19) are some six orders of magnitude smaller than the residuals we obtained in (6.16). In part this means that whereas (6.16) was more than adequate as a smoothing operation, (6.19) is marginally so. In addition, we observe here for the first time one of the singularities from $`\widehat{\chi }^{(6)},`$ and using the cut information given in eqs. (3.26, 3.27) can subtract it as

$`K_{2m}^{}`$ $``$ $`K_{2m}^{}+7^815/(16\pi \widehat{m}^{35/2})(3\sqrt{7}/\pi ^2)^{3/4}\{\mathrm{cos}(2\widehat{m}\mathrm{arccos}(3/4)\pi /8)`$
$``$ $`215\sqrt{7}/(32\widehat{m})\mathrm{sin}(2\widehat{m}\mathrm{arccos}(3/4)\pi /8)+\mathrm{O}(1/\widehat{m}^2)\},`$
where $`\widehat{m}=m1/4`$.

A similar analysis has been carried out for $`T>T_c`$ and our estimates for the coefficients of $`F_\pm `$ are given in the Appendix. The coefficients through O$`(\tau ^5)`$ are unambiguously rational and have been fixed in the final fittings. We have also set to zero all coefficients of odd powers of $`\tau `$ in the product $`\sqrt{s}F_\pm `$ to O$`(\tau ^{15})`$ but have allowed variable coefficients of $`\tau ^{16},`$ $`\tau ^{17}`$ and $`\tau ^{18}.`$ These latter coefficients are, as expected, sensitive to whether we stop the “short-distance” subtraction at O$`(\tau ^{14})`$ or add in the additional terms we have estimated by Padé methods. The terms quoted in the Appendix on the other hand are completely stable and thus we believe reliable except possibly for the last digit.

When we relax the constraint of zero amplitude on individual odd $`\tau ^k`$ terms in $`\sqrt{s}F_\pm `$ for integer $`5k15`$ we find no significant improvement in our fits and the resulting amplitudes are consistent with zero. For example, when $`T<T_c`$, we find best fit coefficients of $`\tau ^k`$, $`k`$ odd, that in absolute magnitude are all less than $`4\times 10^{2k33}`$. For $`T>T_c`$ the corresponding bounds can be as much as $`100`$ times larger. But in both cases these bounds are of a magnitude similar to what we estimate is the uncertainty in the even order coefficients given in the Appendix. Thus we believe $`\sqrt{s}F_\pm `$ even in $`\tau `$ to be an exact symmetry of the scaling-amplitude function.

An alternative analysis using the traditional variable $`v=\mathrm{tanh}K`$ was also carried out. The natural boundary singularities at $`|s|=1`$ are mapped to two circles $`|v\pm 1|=\sqrt{2}`$. In this expansion variable, the ferromagnetic and anti-ferromagnetic critical points are at $`v=\pm (\sqrt{2}1)`$ respectively, and all other points on the two circles are farther away from the origin. Hence the amplitudes of any other singularities are exponentially damped and may be neglected in the analysis. The two analyses are in complete agreement, but the more detailed $`s`$-plane analysis provides greater precision. For contemplated future analyses on other lattices, for which less is known about the natural boundary singularities, it may be necessary to use the $`v`$-plane analysis.

## 7 Expected scaling form of the susceptibility.

The basic scaling Ansatz for the singular part of the free-energy of the two-dimensional Ising model is

$`f_{\mathrm{sing}}(g_t,g_h,\{g_{u_j}\})`$ $`=`$ $`g_t^2\mathrm{log}|g_t|\stackrel{~}{Y}_\pm (g_h/|g_t|^{y_h/y_t},\{g_{u_j}|g_t|^{y_j/y_t}\})`$ (7.1)
$`+`$ $`g_t^2Y_\pm (g_h/|g_t|^{y_h/y_t},\{g_{u_j}|g_t|^{y_j/y_t}\}).`$
Here $`g_t,g_h,g_{u_j}`$ are nonlinear scaling fields associated with the thermal field $`t,`$ the magnetic field $`h`$ and the irrelevant fields $`\{u_j\}.`$ The exponents $`y_t,y_h>0`$ are the thermal and magnetic exponents, and $`y_j<0`$ are the irrelevant exponents.121212For the two-dimensional Ising model, $`y_t=1`$ and $`y_h=\frac{15}{8}`$. This scaling Ansatz assumes only a single resonance, between the identity and the energy. That the dimensions are integers implies higher powers of $`\mathrm{log}t`$. The surprise is the power of $`t`$ at which higher powers of $`\mathrm{log}t`$ enter. We thank Andrea Pelissetto for this clarification.. The nonlinear scaling fields have expansions

$`g_t`$ $`=`$ $`{\displaystyle \underset{n0}{}}a_{2n}(t,u)h^{2n},\text{ }a_0(0,u)=0,`$
$`g_h`$ $`=`$ $`{\displaystyle \underset{n0}{}}b_{2n+1}(t,u)h^{2n+1},`$
$`g_{u_j}`$ $`=`$ $`{\displaystyle \underset{n0}{}}c_{2n}(t,u)h^{2n},`$ (7.2)
where $`a_{2n},b_{2n+1},c_{2n}`$ are smooth functions of $`t`$ and $`u\{u_j\}.`$

If irrelevant fields are neglected, then the known zero field free energy forces the equalities $`\stackrel{~}{Y}_+(0)=\stackrel{~}{Y}_{}(0)`$ and $`Y_+(0)=Y_{}(0)`$. Furthermore, the absence of logarithmic terms in the known magnetization and the divergent part of the susceptibility requires the derivatives $`\stackrel{~}{Y}_\pm ^{}(0)`$ and $`\stackrel{~}{Y}_\pm ^{\prime \prime }(0)`$ to vanish. Aharony and Fisher have conjectured , almost certainly correctly, that there are no logarithms multiplying the leading power law divergence of all higher order field derivatives, in which case the $`\stackrel{~}{Y}_\pm `$ are constants and analyticity on the critical isotherm for $`h0`$ demands $`\stackrel{~}{Y}_+=\stackrel{~}{Y}_{}.`$ With all these constraints built in, one can show the scaling Ansatz (7.1) together with the field expansions (7) lead to

$`f(t,h=0)`$ $`=`$ $`A(a_0(t))^2\mathrm{log}|a_0(t)|+A_0(t),`$
$`(t<0,h=0)`$ $`=`$ $`Bb_1(t)|a_0(t)|^\beta ,`$ (7.3)
$`\beta ^1\chi _\pm (t,h=0)`$ $`=`$ $`C_\pm (b_1(t))^2|a_0(t)|^\gamma Ea_2(t)a_0(t)\mathrm{log}|a_0(t)|+D(t),`$ (7.4)
where $`A,`$ $`B,`$ $`C_\pm `$ and $`E`$ are constants and $`\beta =1/8,`$ $`\gamma =7/4.`$ The free energy and magnetization eqs. (7) determine the scaling field coefficients $`a_0(t)`$ and $`b_1(t)`$ in (7). The presence of irrelevant scaling fields will be expected to manifest themselves as deviations from the predicted form of the susceptibility in (7.4) and/or as deviations from the unique prediction for the coefficient of $`C_\pm .`$

Working as usual in the temperature variable $`\tau =(1/ss)/2`$, we write the predicted isotropic susceptibility from (7.4) as

$$\beta ^1\chi _\pm (\tau ,h=0)=C_{0_\pm }(2K_c\sqrt{2})^{7/4}|\tau |^{7/4}F(\mathrm{A}\&\mathrm{F})E_0/(2K_c\sqrt{2})\tau \mathrm{log}|\tau |e_0(\tau )+D_0(\tau )$$
(7.5)
where $`F(\mathrm{A}\&\mathrm{F})`$ has already been given in eq. (1.19). The “short-distance” contribution to $`\chi `$ is here predicted to be given by the sum of the term containing $`e_0(\tau ),`$ which arises as the mixing of the first two terms in the expansion of $`g_t`$ in (7), and the analytic $`D_0(\tau ).`$

The clear implication of eq. (1.18), which shows the exact $`F_\pm `$ is not equal to $`F(\mathrm{A}\&\mathrm{F})`$ is that irrelevant variables do play a role, and the multipliers $`F_\pm `$ represent the contribution of a number of scaling fields. While there are suggestions in the literature for what these scaling fields might be131313 An excellent early discussion of the types of corrections that might be found, together with a search for some of them can be found in Blöte and den Nijs . Developments in our understanding of the predictions of conformal field theory  lead us to believe that a virtually complete explanation of corrections to scaling is obtainable, at least in principle. Our analysis supports the conclusion of Barma and Fisher that a correction-to-scaling term with exponent $`\theta =4/3`$ is absent for the pure $`S=1/2`$ Ising model susceptibility considered here. A mechanism in terms of generators of the energy family of the Virasoro algebra is adduced by Caselle et al.  which gives rise to corrections at order $`\tau ^4`$ as observed. it is unlikely that a unique identification could be made here since we are dealing with a single isolated model with no free parameters to vary. A corresponding analysis of the anisotropic square, and the triangular and hexagonal lattices is likely to be enlightening in this regard. In  a study of difference equations is given, which implicitly outlines what is needed to obtain the difference equations for the hexagonal and triangular lattices. Additional material in this respect can also be found in .

Note however that two different effects manifest themselves. At fourth order in $`t`$ (or, equivalently, $`\tau `$) scaling under the assumption of only two nonlinear scaling fields breaks down, as evidenced by the difference between the coefficients of $`\tau ^4`$ in eqs. (1.18) and (1.19). However the corresponding high- and low-temperature amplitudes still satisfy $`C_0^+/C_0^{}=C_j^+/C_j^{}`$ for all $`j5.`$ For $`j>5,`$ not only does simple scaling fail to hold, but this equality also breaks down.

In the vicinity of the anti-ferromagnetic point in the high-temperature phase, $`\chi `$ is given exclusively by a “short-distance” term (1.16). Both $`B_{\mathrm{af}}`$ and $`B_\mathrm{f}`$ from eqs. (1.16) and (1.17) have expansions of the same form (1.15), where the sum over $`p`$ is restricted to $`p^2q.`$ The coefficients in this expansion can be determined from the short-distance correlations, and accurate values for the expansions of $`B_{\mathrm{f}/\mathrm{af}}`$ are given in the Appendix through O$`(\tau ^{14}).`$

Again we note that there are terms in these “short-distance” functions that are not of the Aharony and Fisher  predicted form (7.5) based on the absence of irrelevant variables.

We conclude with some speculative remarks. We cannot account physically for terms of order $`t^q(\mathrm{log}|t|)^p,`$ with $`p2`$ and $`qp^2`$ in $`B_f`$, though we can see their origin mathematically, as discussed below eq. (4.15) and in footnote 12. While terms without logarithms and terms of order $`t^q(\mathrm{log}|t|)`$ are expected, it is surprising (to us) that higher powers of $`\mathrm{log}|t|`$ enter at the orders they do.

From conformal field theory, we have predictions for the irrelevant exponents $`y_j=2,4,6,\mathrm{}.`$ The fact that $`f_+^{(k)}=f_{}^{(k)}`$ for $`k<6,`$ though $`f_\pm ^{(k)}`$ is not equal to the corresponding term in (1.19) for $`k=4,5`$ suggests the presence of only a single irrelevant operator contributing at order $`\tau ^4,`$ while the breakdown of high-low temperature symmetry in $`F_\pm `$ at O$`(\tau ^6)`$ suggests that more than one scaling operator couples to the lattice magnetization at this order. A corresponding study to that reported here on the triangular and honeycomb lattices, as well as on the anisotropic square lattice is likely to be enlightening, and we hope to report on this in future.

Acknowledgments We would like to thank Larry Glasser, Mireille Bousquet-Mélou and John Cardy for many helpful discussions and calculations, Andrea Pelissetto for adding to our understanding of scaling theory and the provision of a preprint, Barry McCoy, Alan Sokal and Michael Fisher for invaluable comments on earlier versions of the manuscript, and Douglas Abraham for useful comments on the history of the Ising model. AJG and WO would like to thank the Australian Research Council for financial support, and JHHP thanks the NSF for support in part by NSF Grant no. PHY 97-22159.

## Appendix

Here we list the expansion coefficients $`b_{p,q}`$ of the “short-distance” functions $`B_{\mathrm{f}/\mathrm{af}}`$ defined in (1.15) to O$`(\tau ^{14})`$. The prefactor $`(\sqrt{1+\tau ^2}+\tau )^{1/2}=1/\sqrt{s}`$ is understood to be expanded in a series in $`\tau .`$ The leading constant and coefficient of $`\tau \mathrm{log}|\tau |`$ have been reported previously . The $`F_\pm `$ series as deduced in section 6.3 are also given.

$`B_\mathrm{f}=(\sqrt{1+\tau ^2}+\tau )^{1/2}`$
$`[0.104133245093831026452160126860473433716236727314`$
$`0.07436886975320708001995859169799500328047632028\tau `$
$`0.0081447139091195995371542858655723893266057740\tau ^2`$
$`+0.004504107712232015926355020852986970591364528\tau ^3`$
$`+0.23961879425472180967837072450742931180586109\tau ^4`$
$`0.0025399505953392329612162686121616238176205\tau ^5`$
$`0.235288909669962491804066210882350821445764\tau ^6`$
$`+0.00191570753170091141409998516033460855797\tau ^7`$
$`+0.2143400966115384518711435705343125612378\tau ^8`$
$`0.000883215706003328768611915486246075323\tau ^9`$
$`0.19422062840719623752953468278129284679\tau ^{10}`$
$`+0.0000072335097772632765778839359680528\tau ^{11}`$
$`+0.177102037555467190714704023746648559\tau ^{12}`$
$`+0.0006888110962684387331860926084517\tau ^{13}`$
$`0.16279253648974618861881216566686\tau ^{14}`$
$`+\mathrm{log}|\tau |`$
$`(0.032352268477309406090656526721221666637730948898\tau `$
$`0.0057755293796884630091487564013201013677152980\tau ^3`$
$`+0.059074961290345476578516085774495545264759330\tau ^4`$
$`+0.00305849157585622544005057759535229287938174\tau ^5`$
$`0.0591662722088409053375931018028970139567911\tau ^6`$
$`0.002067088393167114141650194740281136875636\tau ^7`$
$`+0.05424693070421409615112542698595864778919\tau ^8`$
$`+0.0010601025315498815900774416057541651837\tau ^9`$
$`0.049300253157082567741316861339709144063\tau ^{10}`$
$`0.00026830064161204203467706137637400358\tau ^{11}`$
$`+0.0450270525719569186212816103126356308\tau ^{12}`$
$`0.00034332683257234543036792535081332\tau ^{13}`$
$`0.041428586463052869356803144137620\tau ^{14})`$
$`+(\mathrm{log}|\tau |)^2`$
$`(0.0093915698711458721317953318727075770649513654\tau ^4`$
$`0.00869592546287923802156416645191752987912922\tau ^6`$
$`+0.007669481493104540876445085447422616885330\tau ^8`$
$`+0.00015428438297902275440225213783285077606\tau ^9`$
$`0.0068054076881441249098452112921129773269\tau ^{10}`$
$`0.000310520937481414524040686012525223279\tau ^{11}`$
$`+0.00611386643219454473116391019937140965\tau ^{12}`$
$`+0.000444606198235804033861443998682830\tau ^{13}`$
$`0.0055571002151161308034896964314679\tau ^{14})`$
$`+(\mathrm{log}|\tau |)^3`$
$`(0.000015771569138451840480001012621461738178\tau ^9`$
$`+0.0000344282066208887553647799856857753380\tau ^{11}`$
$`0.0000524427177487226174161583779149393\tau ^{13})],`$
$`B_{\mathrm{af}}=(\sqrt{1+\tau ^2}+\tau )^{1/2}`$
$`[0.1588665229609474882333592313690210116925239008416`$
$`+0.149566836938535905194382029433591286374711207262\tau `$
$`+0.01071222587983288033470968550659996768542030678\tau ^2`$
$`+0.0127530188399624019539552078052153609134674971\tau ^3`$
$`0.011741188869656263932121387296300743594029390\tau ^4`$
$`0.01406604087566590060620992322775625815515533\tau ^5`$
$`+0.0131064546156258402249424759665220798848681\tau ^6`$
$`+0.012239696625538370626786005459530159711716\tau ^7`$
$`0.01184019404541084813958321002995560636877\tau ^8`$
$`0.0105854093023116661362507232392645147231\tau ^9`$
$`+0.010151560037724359473553197335905170854\tau ^{10}`$
$`+0.00908000411233119371610549453140718281\tau ^{11}`$
$`0.0085420122287896456879459087815054030\tau ^{12}`$
$`0.00771702694013238358077176900242074\tau ^{13}`$
$`+0.007123677682511208149032476379667\tau ^{14}`$
$`+\mathrm{log}|\tau |`$
$`(0.1553171901580110585934133538932734529992121600305\tau `$
$`+0.03206714814586975221843437287457551882247161782\tau ^3`$
$`0.0077168875724615093064542922101962689299768599\tau ^4`$
$`0.015675211573817078943430169665269657287132910\tau ^5`$
$`0.00028554245153720354627897919087710530890677\tau ^6`$
$`+0.0096072545027321808179041903535130201897217\tau ^7`$
$`+0.004835406420625092236673413378307375358908\tau ^8`$
$`0.00606499034448050751379194815071149626812\tau ^9`$
$`0.0073400150414474023562454611060746360875\tau ^{10}`$
$`+0.003910356521403913091050321141297009252\tau ^{11}`$
$`+0.00870842744568158003036434762719697635\tau ^{12}`$
$`0.002697783010884752101384006121375890\tau ^{13}`$
$`0.0094056230380765607719474925088649\tau ^{14})`$
$`+(\mathrm{log}|\tau |)^2`$
$`(0.01153371437882328027949011442761203640684043805\tau ^4`$
$`0.011311734920691560067535056532207842716405684\tau ^6`$
$`+0.0100457687111988577404299867962466051265974\tau ^8`$
$`0.000475698571097159420906450182271928179428\tau ^9`$
$`0.00878397202228689639470985683437717938463\tau ^{10}`$
$`+0.0011571801729636538264100914359800686355\tau ^{11}`$
$`+0.007680651109512704070606639646988801296\tau ^{12}`$
$`0.0018650912616201532939412831153215046\tau ^{13}`$
$`0.00674470189451526288478200059343432\tau ^{14})`$
$`+(\mathrm{log}|\tau |)^3`$
$`(0.0000578997194764877297760067221144062249541\tau ^9`$
$`0.00016991508824012890240796446744935908812\tau ^{11}`$
$`+0.00032664884687465587957270016883093909\tau ^{13})]`$
$`F_+`$ $`=1+\tau /2+5\tau ^2/8+3\tau ^3/1623\tau ^4/38435\tau ^5/7680.1329693327418753330\tau ^6`$
$`0.05899768720427100\tau ^7+0.121586869804903\tau ^8+0.0766007994119\tau ^9`$
$`0.10751871874\tau ^{10}0.078346589\tau ^{11}+0.0960583\tau ^{12}+0.07592\tau ^{13}`$
$`0.087\tau ^{14}0.1\tau ^{15}+\mathrm{}`$
$`=(\sqrt{1+\tau ^2}+\tau )^{1/2}(1+\tau ^2/2\tau ^4/120.1235292285752086663\tau ^6`$
$`+0.136610949809095\tau ^80.13043897213\tau ^{10}+0.1215129\tau ^{12}0.113\tau ^{14}+\mathrm{})`$
$`F_{}`$ $`=`$ $`1+\tau /2+5\tau ^2/8+3\tau ^3/1623\tau ^4/38435\tau ^5/7686.330746944662603289734\tau ^6`$
$``$ $`3.1578864931646349782\tau ^7+5.46225118896595954\tau ^8+3.521655160482472\tau ^9`$
$``$ $`4.6602157191837\tau ^{10}3.40963923001\tau ^{11}+4.055875878\tau ^{12}+3.2008085\tau ^{13}`$
$``$ $`3.59746\tau ^{14}2.985\tau ^{15}+\mathrm{}`$
$`=`$ $`(\sqrt{1+\tau ^2}+\tau )^{1/2}(1+\tau ^2/2\tau ^4/126.321306840495936623067\tau ^6`$
$`+`$ $`6.25199747046024329\tau ^85.6896599756180\tau ^{10}+5.142218271\tau ^{12}`$
$``$ $`4.67472\tau ^{14}+\mathrm{})`$
The last digit in each term above may not be reliable.

For completeness we give also the leading susceptibility amplitudes evaluated to higher accuracy than reported in :

$`C_0^+`$ $`=`$ $`1.000815260440212647119476363047210236937534925597789(2K_c\sqrt{2})^{7/4}\sqrt{2}`$
$`C_0^{}`$ $`=`$ $`1.000960328725262189480934955172097320572505951770117(2K_c\sqrt{2})^{7/4}\sqrt{2}/(12\pi ).`$