# Finite-temperature phase transitions in $S=1/2$ three-dimensional Heisenberg
magnets from high-temperature series expansions

M. G. Gonzalez  Sorbonne Université, CNRS, Laboratoire de Physique Théorique
de la Matière Condensée, LPTMC, F-75005 Paris, France  Helmholtz-Zentrum
Berlin für Materialien und Energie, Hahn-Meitner Platz 1, 14109 Berlin,
Germany  Dahlem Center for Complex Quantum Systems and Fachbereich Physik,
Freie Universität Berlin, 14195 Berlin, Germany  B. Bernu  Sorbonne
Université, CNRS, Laboratoire de Physique Théorique de la Matière Condensée,
LPTMC, F-75005 Paris, France  L. Pierre  Université Paris X, UFR SGMI
(Sciences gestion, mathématique et informatique), F-92000 Nanterre, France  L.
Messio  Sorbonne Université, CNRS, Laboratoire de Physique Théorique de la
Matière Condensée, LPTMC, F-75005 Paris, France  Institut Universitaire de
France (IUF), F-75005 Paris, France

###### Abstract

By means of high-temperature series expansions we study the singularities of
the specific heat ($c_{v}$) and the magnetic susceptibility ($\chi$) of
$S=1/2$ models presenting a phase transition belonging to the three
dimensional Heisenberg universality class. We first calculate the critical
temperature $T_{c}$ and the critical exponent $\gamma$ using the standard Dlog
Padé method on $\overline{\chi}(\beta)=\chi(\beta)/\beta$ for the
ferromagnetic Heisenberg model on the face-centered cubic, body-centered
cubic, simple cubic, pyrochlore and semi-simple cubic lattices ($\beta=1/T$).
We also explore the possibility of using this method on $\overline{\chi}(e)$
and $c_{v}(\beta)$ to calculate the critical energy $e_{c}$ and the critical
exponent $\alpha$. Finally, we adapt a method initially developed for
logarithmic singularities [Phys. Rev. B 104, 165113 (2021)] to cusp
($-1<\alpha<0$) and divergent singularities ($\gamma>0$), and propose
different interpolation methods for $c_{v}$ and $\chi$. We apply our method to
several of the previously mentioned lattices and present for each the
reconstructed $c_{v}$ and $\chi$ down to $T_{c}$.

##  I Introduction

The quantum ferromagnetic Heisenberg model was firstly introduced to explain
why certain compounds developed a spontaneous magnetization when cooled below
a given temperature (called the critical temperature), even in the absence of
an applied magnetic field. Its success in explaining these so-called
ferromagnets placed the model at a central spot in the study of quantum
magnetism. Even though the Mermin-Wagner theorem[1] precludes the existence of
magnetic order at finite temperatures to two-dimensional (2D) systems, real
quasi-2D materials often present finite-temperature phase transitions due to
remaining 3D correlations[2, 3, 4, 5]. All of these phase transitions belong
to the 3D Heisenberg universality class, defined by the symmetry breaking
$SU(2)\to U(1)$.

Phase transitions are characterized by the critical temperature $T_{c}$ and
exponents of the singularities present in the thermodynamic functions.
Obtaining them is thus an important task, and to do so many different methods
have been developed. In particular, high-temperature series expansions (HTSE)
methods such as the Dlog Padé or ratio methods are known to obtain numerically
accurate critical temperatures for simple lattices such as the simple cubic
(sc), body-centered cubic (bcc) and face-centered cubic (fcc)[6, 7, 8, 9, 10].
More recently the list was expanded with the pyrochlore[11], the diamond[12,
10], and the semi-simple cubic (ssc)[13] lattices. However, for the critical
exponents, the results are not so clear.

For example, the critical exponent $\gamma$ of the magnetic susceptibility
$\chi(T)$ has been calculated using field theory’s renormalization group on
the $N$-vector model (that is in the same universality class for $N=3$),
yielding $\gamma=1.3895(50)$[14, 15]. On the other hand, HTSE calculations in
the quantum spin-$\frac{1}{2}$ case show higher values, $\gamma=1.42(1)$ for
the ferro- and $\gamma=1.43(1)$ for antiferromagnetic cases[7, 8, 9]. Another
example of critical exponent is the less studied $\alpha$ from the specific
heat $c_{v}(T)$. This exponent is negative, what implies a non-divergent
singularity. Instead, $c_{v}(T)$ presents a cusp-like behavior that reaches a
maximum value with an infinite slope. To the best of our knowledge, the
standard HTSE Dlog Padé and ratio methods have never been used on $c_{v}$ in
the literature. However, indirect HTSE calculations through scaling relations
give $\alpha=-0.200(15)$[10]. On the other hand, the field theory result is
$\alpha=-0.122(10)$[15], showing a larger discrepancy than in the case of
$\gamma$.

Finally, it would be desirable to have reliable results not only on the
critical quantities (critical temperature $T_{c}$ and exponents) but on the
thermodynamic functions at all temperatures ($c_{v}(T)$, $\chi(T)$). In this
sense, the quantum Monte Carlo (QMC) calculations obtain reliable results, but
only on finite lattices, and only in the absence of frustration due to the
sign problem. Methods based on exact diagonalization and tensor network
algorithms can be used on frustrated systems, but only in dimensions 1 and
2[16, 17, 18, 19]. Other methods work directly in the thermodynamic limit,
like the rotationally invariant Green’s function method which obtains
qualitatively good results in 3D[20, 21]. The pseudo-Majorana functional
renormalization group provides quantitatively good results down to moderate
temperatures, but deviates from exact results at low temperatures[22, 23].
Finally, the HTSE are quasi-exact at high temperatures but fail close to the
transition temperature, even when using Padé approximants. An interpolation
scheme using HTSE was proposed in our previous article in cases where $c_{v}$
presents a logarithmic divergency, such as in the 2D-Ising or 2D-XXZ
models[24], resulting both in an evaluation of critical quantities and of
$c_{v}(T)$ for temperatures from infinite down to $T_{c}$.

In this article, we firstly revisit the standard HTSE Dlog Padé method for
Heisenberg ferromagnets on several 3D lattices such as the fcc, bcc, sc,
pyrochlore and ssc lattices. For most, we calculated higher orders in the HTSE
than previous articles. Also, we extend the Dlog Padé method to obtain
quantities such as the critical energy $e_{c}$, the $c_{v}$ critical exponent
$\alpha$, and non-universal quantities $A$ and $B$ ($c_{v}(\beta)\sim
B-A(\beta_{c}-\beta)^{-\alpha}$). Finally, we extend the previously mentioned
interpolation method[24] to cases where $c_{v}$ presents a cusp-like behavior
with a negative exponent $\alpha>-1$ and to cases where $\chi$ presents a
divergent singularity with positive exponent $\gamma>0$. For $c_{v}$, two
different interpolation methods for characterizing the singularity are
proposed, that depend on $\alpha$, $\beta_{c}$, and $A$ or $B$. For $\chi$, we
only use one of these two methods, which presents the advantage of a reduced
parameter space. For some 3D lattices, we are able to extrapolate $c_{v}$ and
$\chi$ at all temperatures down to $T_{c}$.

The remaining of the article is organized in the following way. In Sec. II, we
present the HTSE methods to study finite-temperature phase transitions, both
standard and new. In Sec. III we show our results, first for the standard,
next for the new methods. Finally, conclusions and perspectives are given in
Sec. IV.

##  II Model and Methods

The Heisenberg model is defined as

$\mathcal{H}=J\sum_{\langle ij\rangle}\mathbf{S}_{i}\cdot\mathbf{S}_{j},$ (1)

where $J$ is the exchange interaction, the sum $\langle ij\rangle$ runs over
nearest neighbors on a 3D lattice, and $\mathbf{S}_{i}$ are the quantum
spin-$\frac{1}{2}$ operators. The classical approximation consists in
replacing the operators $\mathbf{S}_{i}$ by 3D vectors. In the ferromagnetic
case ($J<0$), the quantum ground-state energy per site $e_{0}$ is exactly the
same as the classical one, namely

$e_{0}=-\frac{Z}{2}S^{2}$ (2)

where $Z$ is the coordination number of the lattice. Even though $e_{0}$ does
not change when taking into account quantum fluctuations, the critical
temperature $T_{c}$ does[8]. On the other hand, in the antiferromagnetic case
($J>0$) on a bipartite lattice, $e_{0}$ and $T_{c}$ are the same as in the
ferromagnetic case at the classical limit, but they both change in the quantum
model and $e_{0}$ is no longer known exactly.

To study these kind of finite temperature phase transitions we use high
temperature series expansions (HTSE). HTSE allow to perform a series expansion
of certain thermodynamic functions around $\beta=0$, where $\beta$ is the
inverse temperature ($\beta=1/T$). Two important functions are the free energy
per site $f$ and the ferromagnetic zero-field susceptibility per site $\chi$
[7]. Their HTSE are written:

$\beta f&=-\ln
2-\frac{1}{n_{u}}\sum_{i=1}^{n}\frac{a_{i}}{4^{i}i!}\beta^{i}+O(\beta^{n+1\textasciiacute})\\\
\overline{\chi}=T\chi=\frac{1}{4}+\frac{1}{2n_{u}}\sum_{i=1}^{m}\frac{b_{i}}{4^{i}i!}\beta^{i}+O(\beta^{m+1}),$
(3)

where $a_{i}$ and $b_{i}$ are integers times $J^{i}$, and $n_{u}$ is the
number of spins in the unit cell. These HTSE are typically known up to orders
13 to 15 for 3D lattices (see Table 1 for the order depending on the lattice).

The thermodynamic functions present singularities at the critical temperature.
However, there are several methods that can be used to extract information
about the critical point from the first coefficients of the series. We present
in Sec. II.1 the most commonly used: the Dlog Padé method, and pursue in Sec.
II.2 with the description of a new interpolation method.

###  II.1 Dlog Padé method

We assume a thermodynamic function $f(x)$ that has a power law singularity at
$x_{c}$, of type:

$f^{s}(x)=A\left(x_{c}-x\right)^{-\theta},$ (4)

such that $f(x)-f^{s}(x)$ is analytic from $x=0$ to some $x>x_{c}$. $x_{c}$ is
the critical point and $\theta$ is the critical exponent. Then, the critical
point and exponent can be obtained from the Dlog Padé method: the logarithmic
derivative

$D\log f^{s}(x)=\frac{{f^{s}}^{\prime}(x)}{f^{s}(x)}=\frac{\theta}{x_{c}-x}$
(5)

has a simple pole given by $x_{c}$, whose residue is the critical exponent
$\theta$. In practice, the critical point and exponent are determined from the
poles and residues of the Padé approximants of the HTSE of $D\log f(x)$.

This method has been used with $f(x)=\overline{\chi}(\beta)$ to obtain results
for $\beta_{c}$ and the critical exponent $\gamma$ on most of the typical 3D
lattices[8, 10, 7]. The Dlog Padé method presents a fast convergence of
$\beta_{c}$ with the HTSE order, giving several significant digits. Notably,
to the best of our knowledge, this method has never been used with other
thermodynamic functions such as $c_{v}$ to determine $\alpha$, or to obtain
the critical values of the energy $e_{c}$ and of the entropy $s_{c}$, what is
now done in Sec. III.2 and III.3.

###  II.2 Interpolation method for cusp singularities

Now we propose an alternative method to extract information on the critical
point, using the specific heat. This is an extension from our previously
introduced interpolation method for the case of logarithmic singularities[24].
In the present case (3D Heisenberg universality class), the singular behavior
of the specific heat is expressed as

$c_{v}^{s}(\beta)=-A\left(\beta_{c}-\beta\right)^{-\alpha}$ (6)

where $A$ is positive. Also, $-1<\alpha<0$ so that there is no divergency at
the critical point. Instead, $c_{v}$ reaches a maximum value $B$ with an
infinite slope (from higher temperatures).

We build a regular function $R(\beta)$ by removing the singular behavior from
the specific heat. We explore two different ways of doing this. The first one,
called the interpolation method 1 (IM1), is analogous to that of Ref. [24]:

$R(\beta)=c_{v}(\beta)-c_{v}^{s}(\beta)$ (7)

With this definition, the value $B$ of $c_{v}$ at the singularity is
$R(\beta_{c})$. The $c_{v}$-HTSE coefficients are calculated for a specific
model up to an order $n$, and the series of $c_{v}^{s}$ are known at all
orders supposing that $\beta_{c}$, $A$ and $\alpha$ are known. Thus, the
$R$-HTSE are obtained at order $n$, whose coefficients depend on $\beta_{c}$,
$A$ and $\alpha$. Compared to the case with a logarithmic divergency[24], the
parameter space has one more dimension.

The other alternative to build the regular function $R$, called interpolation
method 2 (IM2) is:

$R(\beta)=\frac{1}{A}\frac{c_{v}^{s}(\beta)}{c_{v}(\beta)-B}$ (8)

Again, the $R$-HTSE can be obtained up to order $n$ and this time, it depends
on $\beta_{c}$, $B$ (instead of $A$), and $\alpha$. Defined this way,
$R(\beta_{c})=1/A$. It is worth mentioning that other similar methods could be
developed with other regular functions $R$. But the two proposed here are
sufficiently different in nature so that if both methods give similar results,
we consider the results as trust-worthy.

The idea behind these kind of methods is that the parameters to obtain
$R(\beta)$ have to be well-chosen in order for $R$ to be a truly regular
function. This means that the singularity has to be canceled exactly. When
this is done, the Padé approximants of $R(\beta)$ will coincide down to the
critical temperature (and a little further below). The quality of a given set
of parameters $\\{\beta_{c},A,\alpha\\}$ (for IM1) or
$\\{\beta_{c},B,\alpha\\}$ (for IM2) is measured by the quality function
already introduced in Ref. [24],

$Q^{2}=\frac{2}{(n-1)n}\sum_{i=1}^{N_{\mathcal{P}}}\sum_{j=1}^{i-1}M_{\epsilon}\left(\frac{\mathcal{P}_{i}(\beta_{m})-\mathcal{P}_{j}(\beta_{m})}{\overline{F}(\beta_{c})}\right)$
(9)

where $N_{\mathcal{P}}$ is the number of Padé approximants $\mathcal{P}_{i}$
without singularities in the range $[0,\beta_{m}]$, $\beta_{m}$ is chosen
larger that $\beta_{c}$ to check the regular character of $R$ beyond the
critical point. We take $\beta_{m}=(1+\delta)\beta_{c}$ with $\delta=0.05$.
$M_{\epsilon}(x)$ is a smooth function whose value is $\lesssim 1$ when
$x\ll\epsilon$, and $\gtrsim 0$ when $x\gg\epsilon$. We use
$M_{\epsilon}(x)=1/(1+(x/\epsilon)^{8})$ with $\epsilon=0.005$. Finally,
$\overline{F}(\beta_{c})=\frac{1}{2}(\mathcal{P}_{i}(\beta_{c})+\mathcal{P}_{j}(\beta_{c}))$
is the average of the two Padé approsximants at $\beta_{c}$. This $Q$ function
represents roughly the proportion of coinciding Padé approximants down to the
critical temperature. Parameters with $Q>0.5$ are considered as good.

Once a high quality set of parameters is found, $c_{v}(\beta)$ can be
reconstructed by replacing the regular function by any of its coinciding Padé
approximants $\mathcal{P}_{i}$,

$c_{v}(\beta)&=\mathcal{P}_{i}(\beta)-A\left(\beta_{c}-\beta\right)^{-\alpha}\qquad\textrm{for
IM1},\\\
c_{v}(\beta)=B-\frac{\left(\beta_{c}-\beta\right)^{-\alpha}}{\mathcal{P}_{i}(\beta)}\qquad\qquad\,\,\textrm{for
IM2}.$ (10)

###  II.3 Interpolation method for divergent singularities

The two methods presented in the previous subsection can be extended to
quantities with divergent singularities, as the magnetic susceptibility, whose
singular part writes:

$\overline{\chi}^{s}(\beta)=C\left(\beta_{c}-\beta\right)^{-\gamma}$ (11)

where $C$ and $\gamma$ are positive. From this point, the two methods IM1 and
IM2 can be applied as in the previous subsection with the simplification that
no constant term has to be taken into account (the term $B$ of the previous
section can be discarded as the divergency dominates it). This leads to an
important difference between the extensions of IM1 and IM2 to divergent
singularities. For IM1, $C$ has to be taken into account and the parameter
space consists in $\\{\beta_{c},C,\gamma\\}$. But for IM2,

$R(\beta)=\frac{1}{C}\frac{\overline{\chi}^{s}(\beta)}{\overline{\chi}(\beta)}=\frac{\left(\beta_{c}-\beta\right)^{-\gamma}}{\overline{\chi}(\beta)},$
(12)

the parameter space is reduced to $\\{\beta_{c},\gamma\\}$. Thus, we will only
use IM2 to interpolate $\chi$ in the following. From this regular function,
the rest of the method is the same as described in the previous subsection,
and the susceptibility can be reconstructed from any of its coinciding Padé
approximants,

$\chi(\beta)=\beta\frac{\left(\beta_{c}-\beta\right)^{-\gamma}}{\mathcal{P}_{i}(\beta)}.$
(13)

##  III Results

We have numerically but exactly calculated the HTSE of $\beta f(\beta)$ and
${\overline{\chi}}(\beta)$ for several 3D lattices: the fcc, bcc, sc, ssc, and
pyrochlore lattices. The maximum order $n$ depends on the lattice according to
Table 1, where we get the same order for both $\beta f(\beta)$ and
${\overline{\chi}}(\beta)$. Note that we have several orders more than
previous works for several lattices[7, 10, 13, 25]. The new terms in the HTSE
are provided in Appendix A. In what follows, we study the ferromagnetic
Heisenberg model ($J=-1$ in Eq. (1)) on said lattices. However, all methods
are still valid for the antiferromagnetic case, given that the susceptibility
is accordingly defined for each lattice. Except for the pyrochlore lattice,
all the lattices mentioned above are bipartite and thus the corresponding
susceptibility is the staggered susceptibility[12]. For the pyrochlore, and
any other non-bipartite lattice, the existence and nature of a phase
transition in the antiferromagnetic case is not trivial[26, 27, 28]; and
therefore the definition of the susceptibility related to the order parameter
is more complicated.

Lattice | this article | $n_{\beta f}$ | $n_{\overline{\chi}}$
---|---|---|---
fcc | 13 | 12[7] | 14[10]
bcc | 15 | 13[7] | 14[10]
sc | 17 | 13[7] | 14[10]
ssc | 20 | - | 14[13]
pyrochlore | 16 | 13[25] | 12[25]
sc-ssc, Eq. (18) | 13 | - | -
Table 1: The order $n$ of the HTSE used for $\beta f$ and ${\overline{\chi}}$
in this article (we use the same order for both quantities), compared to
orders from other articles. The lattices are fcc (face centered cubic), bcc
(bond centered cubic), sc (simple cublic), ssc (semi simple cubic) and the
pyrochlore lattice. sc-ssc is a model interpolating between sc and ssc,
defined in the main text. The new orders are provided in Appendix A

###  III.1 Dlog Padé method applied to $\overline{\chi}(\beta)$

We use first the Dlog Padé method on $\overline{\chi}(\beta)$, which is the
standard method to obtain the values of $T_{c}$ and the critical exponent
$\gamma$ from HTSE[7, 9]. The results are shown on Fig. 1 for the fcc, bcc,
sc, pyrochlore and ssc lattices. Taking into account all the poles $\beta_{i}$
of the Padé approximants of the logarithmic derivative of
$\overline{\chi}(\beta)$, we define the density of poles as a sum of Gaussian
distributions:

$N(\beta)=\sum_{i}e^{-\frac{1}{2}\left(\frac{\beta_{i}-\beta}{\sigma}\right)^{2}}$
(14)

where $\sigma=0.0002$ for the first three lattices and $\sigma=0.005$ for the
latter. For each lattice we use the poles from the four highest orders in the
corresponding HTSE.

Figure 1: Density of poles as defined in Eq. 14 from the Dlog Padé method on
$\overline{\chi}(\beta)$ for the fcc, bcc, sc, ssc and pyrochlore lattices as
a function of $\beta-\beta_{c}$. For each curve, the poles from the 4 highest
orders are taken.

The fcc lattice has the highest coordination number $Z=12$, and the highest
critical temperature (smallest $\beta_{c}$) of all lattices studied. As a
consequence, the HTSE exploitation leads to high quality results, even though
the large $Z$ limits the highest order $n$ that can be reached. We can see
from Fig. 1 that the values of $\beta_{c}$ are concentrated around a well-
defined value, $\beta_{c}=0.4982(2)$, in agreement with the previous
calculations with the same method[7, 10]. For the bcc lattice ($Z=8$) we get
$\beta_{c}=0.7937(2)$, for the sc lattice $\beta_{c}=1.1926(2)$ and for the
pyrochlore lattice $\beta_{c}=1.39(1)$. This last value is larger than in the
sc lattice even though they both have the same coordination number. It has
been argued that this is caused by the contribution of antiferromagnetic
states to the partition function which is more important in the pyrochlore
than in the sc lattice[29, 11]. Due to frustration, the energy difference
between ferro- and antiferromagnetic states in the pyrochlore lattice is less
that in the sc lattice. Finally, for the ssc lattice, the poles are too scarce
and scattered to extract accurate values of $\beta_{c}$ (lower peaks in Fig. 1
while using a higher $\sigma$). This is not surprising in the case of the ssc,
even at orders as high as 20: because of the low coordination number $Z=3$,
$\beta_{c}$ is large and the system is close to the limit $Z=2$ (where the
system can be mapped into a 1D chain, with no singularity at finite
temperatures).

From the residues we can calculate the value of the critical exponent
$\gamma$. Since all five lattices belong to the same universality class,
$\gamma$ is the same and thus we gather all the results in Fig. 2. In this
case, we use Eq. (14) with $\beta_{i}\to\gamma_{i}$ and $\sigma=0.004$. As can
be seen in the figure, the residues from the pyrochlore and ssc lattices do
not contribute significantly to the final result. In total we get
$\gamma=1.428(10)$, in agreement with previous results[7, 10]. This is
different from the renormalization group value, $\gamma=1.3895(50)$ [14, 15],
see dashed line in Fig. 2. It was proposed that this discrepancy comes from
the low order of the HTSE, and that higher orders might bring the numbers
closer together, as it happens for the Ising model[30, 7, 9]. However, the
present inclusion of higher orders does not seem to point in that direction.
This is indicating that a lot more orders would be needed in order to see an
appreciable shift towards the renormalization group value.

Figure 2: Density of residues ($\gamma$) from the Dlog Padé method on
$\overline{\chi}(\beta)$ for the fcc, bcc, sc, ssc and pyrochlore lattices as
a function of $\gamma$. For each lattice we use the 4 highest orders. The
dashed line indicates the result from field theory renormalization group[14,
15].

So far, the search for $\beta_{c}$ and $\gamma$ has been standard.
Alternatively, we can use our knowledge of $\gamma$ to get $\beta_{c}$ (see
Fig. 3). For all lattices, the residues from Fig. 1 plotted versus their poles
$\beta$ from Fig. 2 collapse on a line with a large positive slope, whose
intersection with $\gamma\simeq 1.4$ gives an approximation of $\beta_{c}$.
Thus, with a correct choice of $\beta_{c}$, the residues versus
$(\beta-\beta_{c})/\beta_{c}$ should collapse on lines crossing at the
universal $\gamma$ for all lattices. Surprisingly, the lines of the fcc, bcc,
sc and pyrochlore lattices present similar slopes, whereas it is smaller for
ssc. The ssc line gives residues at $\gamma$ for $\beta_{c}=4.20(5)$ (see Fig.
3). The behavior of these lines could help to determine critical values when
the points do not clearly accumulate near a single point $(\beta_{c},\gamma)$,
as in the ssc lattice.

Figure 3: Poles and residues from the Dlog Padé method on
$\overline{\chi}(\beta)$ for the fcc, bcc, sc, ssc and pyrochlore lattices as
a function of $(\beta-\beta_{c})/\beta_{c}$. The inset shows the poles and
residues for the pyrochlore lattice.

Another alternative is to use the Diagonal Dlog Padé method as presented in
Ref. [13]. In this method, the Padé approximants of the inverse logarithmic
derivative of $\chi$ are calculated, and only those with the same order in the
denominator and numerator (the diagonal ones) are taken into account.
Estimates of $\beta_{c}$ are then obtained from the least positive root of the
numerator (unless it is also a root of the denominator). These values converge
when the order $n$ increases, as illustrated for the pyrochlore and ssc
lattices in Table 2. In the case of the pyrochlore, we obtain
$\beta_{c}=1.39(1)$ in agreement with our previous results. However, for the
ssc lattice we get a more precise estimate $\beta_{c}=4.20(1)$.

$n$ | 8 | 10 | 12 | 14 | 16 | 18 | 20
---|---|---|---|---|---|---|---
pyrochlore | 1.314 | 1.408 | 1.380 | 1.394 | 1.394 |  |
ssc | 5.043 | 4.343 | 4.359 | 4.351 | 4.206 | 4.209 | 4.202
Table 2: Diagonal Dlog Padé results for $\beta_{c}$ for the ssc and pyrochlore
lattices.

###  III.2 Dlog Padé method applied to $c_{v}(\beta)$

So far, $\beta_{c}$ and the critical exponent $\gamma$ have been determined
using the ferromagnetic susceptibility. The problem on relying on $\chi$ is
that it depends on the order parameter, which is not generally known. Even in
cases where it is known, like for antiferromagnetic models ($J>0$) presenting
a phase transition (bipartite lattices), new HTSE of the antiferromagnetic
staggered susceptibility $\chi_{\rm AF}$ have to be calculated[8], which is
computationally expensive. On the other hand, other thermodynamic functions
such as the specific heat $c_{v}(\beta)$ and the entropy $s(\beta)$ are always
indicative of a phase transition, and ferro- and antiferromagnetic models are
connected by the transformation $\beta\rightarrow-\beta$. So the advantage is
that no new HTSE have to be calculated. However, there is also a disadvantage.
In these functions, the ferro- and antiferromagnetic singularities coexist on
the HTSE, and are always present on the positive and negative $\beta$ axis.
Keeping this in mind, we now try to characterize a phase transition in the
universality class of the Heisenberg 3D model, but without knowing the order
parameter (i.e. $\chi$): for this, we now focus on the $c_{v}(\beta)$
function.

Since $c_{v}$ behaves as $B-A\left(\beta_{c}-\beta\right)^{-\alpha}$ with
$\\-1<\alpha<0$ close to $\beta_{c}$, the Dlog Padé method cannot be used
directly (the logarithmic derivative of $c_{v}$ does not have a simple pole at
$\beta_{c}$). However, the Dlog Padé method can be used on $c_{v}-B$. Doing so
provides a good number of poles near the accepted $\beta_{c}$. As $B$ is a
priori unknown, we can select its value such that we get the highest quality
of results and then deduce $\alpha$. The issue is that there is a wide range
of $B$ values that give high quality results. However, this method allows to
obtain a well-defined dependency between the height of the peak $B$ and the
critical exponent $\alpha$. The resulting $B(\alpha)$ is displayed in Figs. 6
and 7, together with our interpolation method results.

###  III.3 Dlog Padé method applied to $\overline{\chi}(e)$

Finally, we can also study the singularity in $\overline{\chi}(e)$, where $e$
is the energy per site. To determine the type of singularity of this function
at the transition, occurring at the critical energy $e_{c}$, we start from the
singularity in $c_{v}$, which can be re-written as

$c_{v}^{s}(\beta)=B-\tilde{A}\Delta T^{-\alpha}$ (15)

where $\tilde{A}=A/T_{c}^{-2\alpha}$ and $\Delta T=T-T_{c}$. By integration we
get that close to the critical point

$\Delta e=B\Delta T+\frac{\tilde{A}}{1-\alpha}\Delta T^{1-\alpha}+o(\Delta
T^{1-\alpha})$ (16)

where $\Delta e=e-e_{c}$. In order to get information on the singularity in
$\overline{\chi}(e)$ we need the inverse $\Delta T(\Delta e)$. Since $\alpha$
is negative, the leading order is in $\Delta T$. Then

$\Delta T=\frac{\Delta e}{B}-\frac{\tilde{A}}{(1-\alpha)B^{2-\alpha}}\Delta
e^{1-\alpha}+o(\Delta e^{1-\alpha})$ (17)

Keeping only the leading order and knowing that
$\overline{\chi}^{s}(T)\propto\Delta T^{-\gamma}$ leads to the simple result
that $\overline{\chi}^{s}(e)\propto\Delta e^{-\gamma}$, and the Dlog Padé
method should give $e_{c}$ as pole and $\gamma$ as residue. However, taking
into account that $\alpha$ is between $-0.1$ and $-0.2$, the second leading
term has similar order compared to the leading term. The quotient between both
terms depends on $\Delta e^{-\alpha}$. This corresponds to a cusp-like
singularity that reaches 0 only at the critical energy $\Delta e=0$. For
example, using typical values of $A$, $B$, and $T_{c}$ on the fcc lattice,
this quotient is about 1 when $\Delta e=1$ and about 0.25 when $\Delta
e=0.0001$. In conclusion, only singularities at $\Delta e=e-e_{c}=0$ exist,
but the simple pole assumption is valid only infinitesimally close to the
critical point and therefore the HTSE are not able to represent it accurately.
Thus, the Dlog Padé method can be used to obtain values for $e_{c}$ from the
poles, but the residues cannot capture the values of $\gamma$ or $\alpha$.

Figure 4: Density of poles from the Dlog Padé method on $\overline{\chi}(e)$
for the ferromagnetic model. Results are shown as a function of $e-e_{c}$,
where $e_{c}$ is the value at which there is a peak. For each lattice, poles
from HTSE at orders from $n-3$ to $n$ are used, where $n$ is given in Tab. 1.

From the HTSEs of $\beta f(\beta)$ and of $\overline{\chi}(\beta)$ at order
$n$, we obtain the series of $\overline{\chi}(e)$ at order $n-1$ (because the
series of $e(\beta)$ are of order $n-1$). Then we use the Dlog Padé method on
$\overline{\chi}(e)$ and obtain the critical energies for all lattices in the
ferromagnetic case (see Fig. 4). We use Eq. (14) with $e$ (instead of $\beta$)
and $\sigma=0.01$ for the fcc, bcc, sc and pyrochlore lattices, and
$\sigma=0.001$ for the ssc lattice. We find $e_{c}=-0.87(1)$ for the fcc
($\sim 58\%$ of the ground-state energy), $e_{c}=-0.61(1)$ for the bcc ($\sim
62\%$), $e_{c}=-0.52(1)$ for the sc ($\sim 70\%$), $e_{c}=-0.57(1)$ for the
pyrochlore ($\sim 76\%$), and $e_{c}=-0.302(1)$ for the ssc ($\sim 81\%$).
Contrary to the $\overline{\chi}(\beta)$ case, the method on
$\overline{\chi}(e)$ works notably better for the ssc lattice than for the
rest. The pyrochlore lattice has a noticeably lesser amount of poles around
$e_{c}$. The reason is again that this lattice is not bipartite, and the
antiferromagnetic solutions on the positive $e$ axis are frustrated. This
leads to a large amount of poles appearing in the positive $e$ axis at values
$e^{*}<|e_{c}|$. Finally, the residues are different for all lattices,
indicating a dependency on non-universal quantities such as $A$, $B$ and
$T_{c}$. As it was expected from our previous analysis, $\alpha$ or $\gamma$
cannot be extracted.

Figure 5: Critical temperature $T_{c}$ (top panel) and difference between the
ground-state and critical energies $e_{0}-e_{c}$ (bottom panel) as a function
of the coordination number $Z$ for the fcc, bcc, sc, ssc and pyrochlore
lattices obtained from the Dlog Pade method on $\overline{\chi}$. We also show
the interpolation between the sc and ssc lattices by using an effective $Z$
(see main text).

###  III.4 $Z=2$ limit

We summarize our results obtained with the Dlog Padé method for fcc, bcc, sc,
pyrochlore and ssc lattices in Fig. 5. We plot the critical temperature
$T_{c}$ extracted from $\overline{\chi}(\beta)$, together with the difference
between the critical energy $e_{c}$ extracted from $\overline{\chi}(e)$ and
the ground state energy $e_{0}$ (known exactly for the ferromagnetic case), as
a function of the coordination number $Z$. Both $T_{c}$ and $e_{c}$ show a
linear behavior with respect to $Z$, and $Z=2$ is a critical point for the
finite temperature transitions in 3D ferromagnets, corresponding to a one-
dimensional chain, characterized by $T_{c}=0$ and $e_{c}=e_{0}$. To get more
points, we define the sc-ssc model, interpolating between sc and ssc lattices
(for which we have the HTSE up to order $n=13$), with two types of links on
the cubic lattice:

$\mathcal{H}=J_{1}\sum_{\langle
ij\rangle}\mathbf{S}_{i}\cdot\mathbf{S}_{j}+J_{1}^{\prime}\sum_{\langle
ij\rangle^{\prime}}\mathbf{S}_{i}\cdot\mathbf{S}_{j}$ (18)

in such a way that when $J_{1}=J_{1}^{\prime}=1$ we have the sc lattice. When
$J_{1}=1$ and $J_{1}^{\prime}=0$ (or vice versa), the Hamiltonian becomes that
of the ssc lattice. We also know that the coordination number goes from $Z=3$
at $J_{1}^{\prime}=0$ to $Z=6$ at $J_{1}^{\prime}=1$, so we can define an
effective coordination number
$Z_{\text{eff}}(J_{1}^{\prime})=3(J_{1}+J_{1}^{\prime})$ such that $e_{0}$ is
proportional to $Z_{\text{eff}}$. The discrepancies between pentagons and
circles of Fig. 5 at $Z=3$ and 6 (more visible for the energies) are due to
different HTSE orders (13 for the sc-ssc lattice, and 17 and 20 for the sc and
ssc, respectively). In addition, we also continued our calculations for the
frustrated case $J_{1}^{\prime}/J_{1}<0$ and found that $T_{c}$ vanishes for
$\left(J_{1}^{\prime}\right)_{c}=-0.15(3)J_{1}$.

###  III.5 Interpolation methods for $c_{v}(\beta)$

Figure 6: Top panel: values of the singularity parameters $A$ and $B$ as a
function of $\alpha$ for the fcc lattice obtained from the two interpolation
methods IM1 and IM2. We show the results for the three highest orders of the
HTSE of $c_{v}$. The pink dots are the results obtained with the Dlog Padé
method on $c_{v}(\beta)-B$. Bottom panel: the quality $Q$ (see Eq. (9)).

The interpolation methods IM1 and IM2 presented in Sec. II.2 for
$c_{v}(\beta)$ have a three dimensional parameter space:
$(\beta_{c},A,\alpha)$ and $(\beta_{c},B,\alpha)$ for IM1 and IM2,
respectively. It is one more than the similar method used for logarithmic
divergencies[24], which does not have to determine a critical exponent.
Exploring the whole parameter space is thus time-consuming, so it is
convenient to rely on other methods to narrow down some of the dimensions. In
this sense, the Dlog Padé method on $\overline{\chi}$ studied in Sec. III.1
provides accurate values for the inverse critical temperature, $\beta_{c}$.
Thus, we leave this parameter fixed. Regarding the value of $\alpha$, we know
that the renormalization group value is $-0.122(10)$ [15], while indirect
estimations from HTSE throw out $-0.200(15)$[10]. So, for this parameter, we
will search in the range $[0.05,0.30]$ at 0.01 intervals. For each value of
$\alpha$ we search for the best value of $A$ or $B$ (for IM1 or IM2) in a
range from 0.1 to 9, at 0.002 intervals. It is important that this step is
small since the peaks in $Q(A)$ for a given $\alpha$ and $\beta_{c}$ tend to
be very narrow.

Figure 7: Values of the non-universal parameters $A$ and $B$ as a function of
$\alpha$ for the fcc, bcc and sc lattices; using only the highest order in the
HTSE. Results obtained with IM1 and IM2 are shown with full and dashed lines,
while symbols correspond to the Dlog Padé method on $c_{v}(\beta)-B$.

Fig. 6 shows the $A$ and $B$ values depending on the choice of $\alpha$,
obtained through IM1 and IM2, for the fcc lattice, together with the Dlog Padé
method results on $c_{v}-B$, discussed in Sec. III.2. The results for $A$ and
$B$ show a good convergence with the HTSE order (specially at higher values of
$\alpha$). Furthermore there is a good agreement between all three methods.
Let us recall that IM1 and IM2 remove the singularity by subtracting and
dividing, respectively, such that for IM1, $A$ is a fitting parameter and $B$
is a byproduct. For IM2 it is the other way around. The quality $Q$ takes high
values throughout all the $\alpha$ range: over $80\%$ of Padé approximants
coincide past $\beta_{c}$ (up to $(1+\delta)\beta_{c}$). Even though there is
a tendency to higher $Q$ values for $\alpha$ closer to 0, it is not possible
to pick one good value for $\alpha$, even choosing more restricting values for
the $Q$-parameters, $\delta$ and $\epsilon$.

For the bcc lattice, $Q$ is between 0.7 to 0.9, whereas on the sc lattice,
goes from 0.5 to 0.7 as $\alpha$ gets closer to zero. However, having half of
the Padé approximants down to $T_{c}$ is still a very good solution, since
none of them are the same at that point when taking the raw HTSE. Also, these
lattices show a convergence of the $A$ and $B$ values with the HTSE order $n$
that is similar to the fcc lattice, using IM1 and IM2. $A$ and $B$ values for
the highest HSE order are given in Fig. 7 for the fcc,bcc and sc lattices. For
all lattices, IM1 and IM2 give similar results, specially in the case of $B$.
Differences only show up for $A$ at values of $\alpha$ far from 0.

Another interesting feature is that the values of $B$ (the peak height) are
very similar for the three lattices, while the values of $A$ seem to change
slowly as the coordination number $Z$ changes. A similar thing happens with
the parameters of the singularities on Ising models on 2D, where the values
are very similar but not universal in different lattices[24]. For the
pyrochlore lattice, $Q$ takes lower values, between 0.3 and 0.4. However, we
can still extract values of $A$ and $B$. They are smaller than in the sc
lattice, even though both have the same coordination number $Z$. Taking these
four lattices, $A$ and $B$ decrease as $T_{c}$ decreases. When $Z=2$ there is
no finite-temperature phase transition, so it might be interesting to see how
this limit is reached in terms of $A$ and $B$. Finally, for the ssc lattice,
no clear peak can be determined.

Even though it is not possible to determine the critical exponent $\alpha$
with these methods, we obtain well-defined functions for $A(\alpha)$ and
$B(\alpha)$. Thus, we can reconstruct $c_{v}$ above $T_{c}$ for a supposed
value of $\alpha$. We tried a workaround to determine $\alpha$, using the
reconstructed $c_{v}$ to calculate the critical energy $e_{c}$ by integration,
which depends on the parameters. For the fcc lattice, the values of $e_{c}$ go
from $-0.860$($-0.857$) to $-0.867$ ($-0.867$) for IM1 (IM2) as $\alpha$
changes from $-0.3$ to $-0.05$. Again, higher values of $\alpha$ show a better
agreement between methods. These values are in agreement with the Dlog Padé
estimation from the previous section ($e_{c}=-0.87(1)$). Unluckily, the Dlog
Padé method does not offer sufficiently precise values of the energy and the
function $e_{c}(\alpha)$ obtained by integration changes very little. So that
in the end, it is not possible to use this extra information to determine the
value of $\alpha$, which remains elusive. For the bcc lattice these values go
from $-0.605$ ($-0.604$) to $-0.610$ ($-0.610$) for IM1 (IM2), again in
agreement with the Dlog Padé result $e_{c}=-0.61(1)$. For the sc lattice these
values go from $-0.510$ ($-0.509$) to $-0.513$ ($-0.513$) for IM1 (IM2). These
values are in agreement with the Dlog Padé result $e_{c}=-0.52(1)$.

Figure 8: Reconstructed $c_{v}$ from the two interpolation methods for the
fcc, bcc, sc, and pyrochlore lattices. We show the results for two different
values of $\alpha$. The inset shows a zoom for the fcc lattice close to the
critical temperature.

We show in Fig. 8 the reconstructed $c_{v}$ for the four lattices (fcc, bcc,
sc and pyrochlore) using the best values of $A$ and $B$ at the accepted
$\beta_{c}$ for two limiting values of $\alpha$, $-0.2$ and $-0.1$. Both
methods give the same curves, and the differences between the two values of
$\alpha$ can only be seen very close to the corresponding critical points (see
inset for fcc) through a very different value of the peak height $B$, as can
be seen in the previous figures of $B(\alpha)$. However, this issue only
exists at exactly the critical temperature, so that it does not affect the
comparison with experimental results. The sharp theoretical peaks cannot be
captured by experiments in real compounds. To sum up, we have a good precision
for every temperature above the critical temperature $T_{c}$ obtained from
finite high-temperature series expansions. Thus, this method extrapolates the
specific heat from HTSE down to almost the critical temperature for the phase
transitions of several ferromagnetic Heisenberg models.

###  III.6 Interpolation method for $\chi(\beta)$

Finally, we apply the interpolation method IM2 to $\overline{\chi}(\beta)$ (as
explained in Sec. II.3 for divergent singularities). The parameter space is
two-dimensional $\\{\beta_{c},\gamma\\}$, and the region of high quality
values is narrow so that both parameters have to be calculated using a fine
mesh.

Fig. 9 shows the quality $Q$ as a function of $\beta_{c}$ and $\gamma$ for the
fcc, bcc, sc, and pyrochlore lattices. We also show in white circles the
results from the poles and residues obtained from the Dlog Padé method. For
all lattices, the poles and residues are concentrated around the large $Q$
region from IM2 and reciprocally, the higher $Q$ values are obtained close to
the line of poles and residues from the Dlog Padé method. This illustrates a
close connection between both methods. Specifically, the fcc lattice presents
$Q=1.00$ around $\beta_{c}=0.4981(2)$ and $\gamma=1.422(2)$, the bcc lattice
presents $Q=0.93$ for $\beta_{c}=0.7938(2)$ and $\gamma=1.420(3)$, the sc
lattice presents $Q=0.64$ for $\beta_{c}=1.1935(10)$ and $\gamma=1.44(1)$. All
of these are in agreement with the Dlog Padé results. For the pyrochlore and
ssc (not shown) lattices, the results are not so clear. The pyrochlore lattice
has a large cloud of values $Q=0.41(1)$ along a well-defined line around
$\beta_{c}=1.382(5)$ and $\gamma=1.36(3)$. However, the diagonal Dlog Padé
results lie closer to the end point of this cloud. For the ssc lattice, there
are just a few points around $Q=0.23(2)$ with $\beta_{c}=4.20(2)$ and
$\gamma=1.35(2)$, but the quality is too low to consider them reliable. With
all this data, we reconstruct $\chi(\beta)$ for all lattices (Fig. 10).

Figure 9: Interpolation method IM2 results for $\overline{\chi}(\beta)$ for
the fcc, bcc, sc and pyrochlore lattices. The quality of results $Q$ is shown
in color scale in a region of the parameter space $\\{\beta_{c},\gamma\\}$
close to the best values. White circles indicate the poles and residues from
the standard Dlog Padé method. Figure 10: Reconstructed $\chi$ from the best
values ${\beta_{c},\gamma}$ from the interpolation method IM2. The values of
$\beta_{c}$ are shown by the dashed lines in the color of the corresponding
lattice.

##  IV Conclusions and perspectives

We have studied the finite-temperature phase transition that occurs in
ferromagnetic quantum Heisenberg models on 3D lattices by using several
methods derived from the HTSE. We used the standard Dlog Padé method and
estimated $\beta_{c}$ and $\gamma$ of the fcc, bcc, sc, pyrochlore and ssc
lattices. For some of them, results are given for larger orders of HTSE than
in the previous works. We have also explored possible extensions of these
methods. While standard calculations involve $\chi(\beta)$, we have obtained
the critical energy $e_{c}$ using $\chi(e)$. We also used the Dlog Padé method
on $c_{v}(\beta)-B$ (with $B=c_{v}(\beta_{c})$) to obtain $B(\alpha)$.

Then we presented new interpolation methods to obtain $c_{v}(T)$ and $\chi(T)$
for $T>T_{c}$. These methods are efficient for the fcc, bcc, and sc lattices,
but less for the pyrochlore and ssc lattices. For $c_{v}(T)$, we are not able
to estimate a precise value of the critical exponent $\alpha$, but the methods
provide accurate relationships between the three important parameters at the
singularity, $A$, $B$, and $\alpha$. Thus, if any of them is known, the other
two can be deduced. We have also shown that the interpolated $c_{v}(T)$ has a
very small dependence on $\alpha$ as soon as $T$ is slightly above $T_{c}$,
the main difference being in the value of the peak at $T_{c}$. The
interpolation method also allows to obtain $\chi(T)$ above $T_{c}$ with
reliable values of $T_{c}$ and $\gamma$.

In conclusion, we have probed several different methods based on HTSE to study
finite-temperature phase transitions. These methods allowed us to obtain
accurately several quantities related to the critical points, such as critical
exponents, critical temperatures, and parameters related to the singularities.
As a summary, we present the main numerical results in Table 3, where DLP
stands for the Dlog Padé results and IM stands for the interpolation method.
It is important to note that even if we have approached only the ferromagnetic
Heisenberg model on the most common lattices, these methods are suitable for
studying any kind of systems with the same type of phase transitions.

| $\beta_{c}$ | $\gamma$ | $e_{c}$ | $\alpha=-0.1$ | $\alpha=-0.2$
---|---|---|---|---|---
Lattice | DLP | IM | DLP | IM | DLP | IM |  $A$ (IM) |  $B$ (IM) |  $B$ (DLP) |  $A$ (IM) |  $B$ IM |  $B$ (DLP)
fcc | 0.4982(2) | 0.4981(2) | 1.426(2) | 1.422(2) | -0.87(1) | -0.862(5) | 3.25(5) | 2.89(3) | 2.90(5) | 2.41(5) | 1.75(2) | 1.73(2)
bcc | 0.7937(2) | 0.7938(2) | 1.419(2) | 1.420(3) | -0.61(1) | -0.607(3) | 3.10(5) | 2.85(2) | 2.85(5) | 2.17(6) | 1.71(2) | 1.71(2)
sc | 1.1926(2) | 1.1935(10) | 1.433(3) | 1.44(1) | -0.52(1) | -0.511(2) | 2.71(3) | 2.63(3) | 2.62(5) | 1.80(3) | 1.57(2) | 1.56(3)
pyrochlore | 1.39(1) | 1.382(5) |  | 1.36(3) | -0.57(1) | -0.578(3) | 2.15(15) | 2.1(1) |  | 1.33(3) | 1.22(3) |
ssc | 4.20(1) | 4.20(2) |  | 1.35(2) | -0.302(1) |  |  |  |  |  |  |
Table 3: Summary of results obtained in this article for the fcc, bcc, sc,
pyrochlore, and ssc lattices. DLP stands for Dlog Padé method and IM indicates
the corresponding interpolation method for each quantity.

## Acknowledgments

This work was supported by the French Agence Nationale de la Recherche under
Grant No. ANR-18-CE30-0022-04 LINK

##  Appendix A HTSE for the $S=1/2$ models

In this section we present the new terms added to the known HTSE of $\beta f$
and $\overline{\chi}$. These are written in terms of Eq. 3 and, unless
otherwise stated, $n_{u}=1$.

###  A.1 fcc lattice

For this lattice we calculated one more order in the $\beta f$-HTSE with
respect to Ref. [7]:

$a_{13}=-421817449494804480\ J^{13}$ (19)

###  A.2 bcc lattice

For this lattice we calculated two more orders in the $\beta f$-HTSE with
respect to Ref. [7] and one more order for $\overline{\chi}$-HTSE with respect
to Ref. [10]:

$\displaystyle a_{14}$ $\displaystyle=$ $\displaystyle 246102905022713856\
J^{14}$ $\displaystyle a_{15}$ $\displaystyle=$ $\displaystyle
9001661201883684864\ J^{15}$ (20) $\displaystyle b_{15}$ $\displaystyle=$
$\displaystyle-88805626440393148956672\ J^{15}$ (21)

###  A.3 sc lattice

For this lattice we calculated four more orders in the $\beta f$-HTSE with
respect to Ref. [7] and three more orders for $\overline{\chi}$-HTSE with
respect to Ref. [10]:

$\displaystyle a_{14}$ $\displaystyle=$ $\displaystyle-27667884260938752\
J^{14}$ $\displaystyle a_{15}$ $\displaystyle=$ $\displaystyle
2908030732698175488\ J^{15}$ $\displaystyle a_{16}$ $\displaystyle=$
$\displaystyle 122264703581556307968\ J^{16}$ $\displaystyle a_{17}$
$\displaystyle=$ $\displaystyle-7238339805811283361792\ J^{17}$ (22)

$\displaystyle b_{15}$ $\displaystyle=$ $\displaystyle-220236945885669801984\
J^{15}$ $\displaystyle b_{16}$ $\displaystyle=$ $\displaystyle
12562562473105938481152\ J^{16}$ $\displaystyle b_{17}$ $\displaystyle=$
$\displaystyle-722105535259151290073088\ J^{17}$ (23)

###  A.4 ssc lattice

For this lattice we calculated the complete $\beta f$-HTSE and six more orders
for $\overline{\chi}$-HTSE with respect to Ref. [13]. In this case, $n_{u}=4$.

$\displaystyle a_{1}$ $\displaystyle=$ $\displaystyle 6\ J^{1}$ $\displaystyle
a_{2}$ $\displaystyle=$ $\displaystyle 18\ J^{2}$ $\displaystyle a_{3}$
$\displaystyle=$ $\displaystyle 36\ J^{3}$ $\displaystyle a_{4}$
$\displaystyle=$ $\displaystyle-324\ J^{4}$ $\displaystyle a_{5}$
$\displaystyle=$ $\displaystyle-3600\ J^{5}$ $\displaystyle a_{6}$
$\displaystyle=$ $\displaystyle 20592\ J^{6}$ $\displaystyle a_{7}$
$\displaystyle=$ $\displaystyle 788256\ J^{7}$ $\displaystyle a_{8}$
$\displaystyle=$ $\displaystyle-267552\ J^{8}$ $\displaystyle a_{9}$
$\displaystyle=$ $\displaystyle-292582656\ J^{9}$ $\displaystyle a_{10}$
$\displaystyle=$ $\displaystyle-2338428672\ J^{10}$ $\displaystyle a_{11}$
$\displaystyle=$ $\displaystyle 158857763328\ J^{11}$ $\displaystyle a_{12}$
$\displaystyle=$ $\displaystyle 3398565523968\ J^{12}$ $\displaystyle a_{13}$
$\displaystyle=$ $\displaystyle-111431579830272\ J^{13}$ $\displaystyle
a_{14}$ $\displaystyle=$ $\displaystyle-5116416515020800\ J^{14}$
$\displaystyle a_{15}$ $\displaystyle=$ $\displaystyle 83476825611595776\
J^{15}$ $\displaystyle a_{16}$ $\displaystyle=$ $\displaystyle
9038092645962510336\ J^{16}$ $\displaystyle a_{17}$ $\displaystyle=$
$\displaystyle-20724045060052942848\ J^{17}$ $\displaystyle a_{18}$
$\displaystyle=$ $\displaystyle-18839898190998133604352\ J^{18}$
$\displaystyle a_{19}$ $\displaystyle=$
$\displaystyle-253691725481243238334464\ J^{19}$ $\displaystyle a_{20}$
$\displaystyle=$ $\displaystyle 45155117370822689756676096\ J^{20}$ (24)

$\displaystyle b_{15}$ $\displaystyle=$ $\displaystyle-572943639086174208\
J^{15}$ $\displaystyle b_{16}$ $\displaystyle=$ $\displaystyle
20821759189681766400\ J^{16}$ $\displaystyle b_{17}$ $\displaystyle=$
$\displaystyle 1524898473896350777344\ J^{17}$ $\displaystyle b_{18}$
$\displaystyle=$ $\displaystyle-22745675831893506785280\ J^{18}$
$\displaystyle b_{19}$ $\displaystyle=$
$\displaystyle-4446640583932914089459712\ J^{19}$ $\displaystyle b_{20}$
$\displaystyle=$ $\displaystyle-17616386456676250248806400\ J^{20}$ (25)

###  A.5 pyrochlore lattice

For this lattice we calculated three more orders in the $\beta f$-HTSE and
four more orders for $\overline{\chi}$-HTSE with respect to Ref. [25]. In this
case, $n_{u}=4$.

$\displaystyle a_{14}$ $\displaystyle=$ $\displaystyle-1532961802218000384\
J^{14}$ $\displaystyle a_{15}$ $\displaystyle=$
$\displaystyle-44591351194841260032\ J^{15}$ $\displaystyle a_{16}$
$\displaystyle=$ $\displaystyle 6653104879154138357760\ J^{16}$ (26)

$\displaystyle b_{13}$ $\displaystyle=$ $\displaystyle 85410304429842432\
J^{13}$ $\displaystyle b_{14}$ $\displaystyle=$
$\displaystyle-1996581576629084160\ J^{14}$ $\displaystyle b_{15}$
$\displaystyle=$ $\displaystyle-507546664875986436096\ J^{15}$ $\displaystyle
b_{16}$ $\displaystyle=$ $\displaystyle 35052604281755859025920\ J^{16}$ (27)

## References

  * Mermin and Wagner [1966] N. D. Mermin and H. Wagner, Absence of Ferromagnetism or Antiferromagnetism in One- or Two-Dimensional Isotropic Heisenberg Models, Phys. Rev. Lett. 17, 1133 (1966).
  * Sengupta  _et al._ [2003] P. Sengupta, A. W. Sandvik, and R. R. P. Singh, Specific heat of quasi-two-dimensional antiferromagnetic Heisenberg models with varying interplanar couplings, Phys. Rev. B 68, 094423 (2003).
  * Yasuda  _et al._ [2005] C. Yasuda, S. Todo, K. Hukushima, F. Alet, M. Keller, M. Troyer, and H. Takayama, Néel Temperature of Quasi-Low-Dimensional Heisenberg Antiferromagnets, Phys. Rev. Lett. 94, 217201 (2005).
  * Tsirlin  _et al._ [2012] A. A. Tsirlin, A. Möller, B. Lorenz, Y. Skourski, and H. Rosner, Superposition of ferromagnetic and antiferromagnetic spin chains in the quantum magnet BaAg2Cu[VO${}_{4}{]}_{2}$, Phys. Rev. B 85, 014401 (2012).
  * Tsirlin  _et al._ [2013] A. A. Tsirlin, O. Janson, S. Lebernegg, and H. Rosner, Square-lattice magnetism of diaboleite Pb2Cu(OH)4Cl2, Phys. Rev. B 87, 064404 (2013).
  * Baker  _et al._ [1967] G. A. Baker, H. E. Gilbert, J. Eve, and G. S. Rushbrooke, High-Temperature Expansions for the Spin-½ Heisenberg Model, Phys. Rev. 164, 800 (1967).
  * Oitmaa and Bornilla [1996] J. Oitmaa and E. Bornilla, High-temperature-series study of the spin-$\frac{1}{2}$ Heisenberg ferromagnet, Phys. Rev. B 53, 14228 (1996).
  * Oitmaa and Zheng [2004] J. Oitmaa and W. Zheng, Curie and Néel temperatures of quantum magnets, Journal of Physics: Condensed Matter 16, 8653 (2004).
  * Oitmaa  _et al._ [2006] J. Oitmaa, C. Hamer, and W. Zheng,  _Series expansion methods for strongly interacting lattice models_ (Cambridge University Press, 2006).
  * Kuz’min [2019] M. D. Kuz’min, Extended high-temperature series for the Spin- Heisenberg ferromagnet, Philosophical Magazine Letters 99, 338 (2019).
  * Lohmann  _et al._ [2014] A. Lohmann, H.-J. Schmidt, and J. Richter, Tenth-order high-temperature expansion for the susceptibility and the specific heat of spin-$s$ Heisenberg models with arbitrary exchange patterns: Application to pyrochlore and kagome magnets, Phys. Rev. B 89, 014415 (2014).
  * Oitmaa [2018] J. Oitmaa, Diamond lattice Heisenberg antiferromagnet, Journal of Physics: Condensed Matter 30, 155801 (2018).
  * Kuz’min  _et al._ [2020] M. D. Kuz’min, R. O. Kuzian, and J. Richter, Ferromagnetism of the semi-simple cubic lattice, The European Physical Journal Plus 135, 750 (2020).
  * Le Guillou and Zinn-Justin [1980] J. C. Le Guillou and J. Zinn-Justin, Critical exponents from field theory, Phys. Rev. B 21, 3976 (1980).
  * Guida and Zinn-Justin [1998] R. Guida and J. Zinn-Justin, Critical exponents of theN-vector model, Journal of Physics A: Mathematical and General 31, 8103 (1998).
  * Prelovšek and Kokalj [2018] P. Prelovšek and J. Kokalj, Finite-temperature properties of the extended Heisenberg model on a triangular lattice, Phys. Rev. B 98, 035107 (2018).
  * Prelovšek and Kokalj [2020] P. Prelovšek and J. Kokalj, Similarity of thermodynamic properties of the Heisenberg model on triangular and kagome lattices, Phys. Rev. B 101, 075105 (2020).
  * Prelovšek  _et al._ [2020] P. Prelovšek, K. Morita, T. Tohyama, and J. Herbrych, Vanishing Wilson ratio as the hallmark of quantum spin-liquid models, Phys. Rev. Research 2, 023024 (2020).
  * Gauthé and Mila [2022] O. Gauthé and F. Mila, Thermal Ising transition in the spin-1/2 J1-J2 Heisenberg model, arXiv preprint 10.48550/arxiv.2201.02171 (2022).
  * Juhász Junger  _et al._ [2009] I. Juhász Junger, D. Ihle, and J. Richter, Thermodynamics of layered Heisenberg magnets with arbitrary spin, Phys. Rev. B 80, 064425 (2009).
  * Müller  _et al._ [2015] P. Müller, J. Richter, A. Hauser, and D. Ihle, Thermodynamics of the frustrated J1-J2 Heisenberg ferromagnet on the body-centered cubic lattice with arbitrary spin, The European Physical Journal B 88, 159 (2015).
  * Niggemann  _et al._ [2021] N. Niggemann, B. Sbierski, and J. Reuther, Frustrated quantum spins at finite temperature: Pseudo-Majorana functional renormalization group approach, Phys. Rev. B 103, 104431 (2021).
  * Niggemann  _et al._ [2022] N. Niggemann, J. Reuther, and B. Sbierski, Quantitative functional renormalization for three-dimensional quantum Heisenberg models, SciPost Phys. 12, 156 (2022).
  * Gonzalez  _et al._ [2021] M. G. Gonzalez, B. Bernu, L. Pierre, and L. Messio, Logarithmic divergent specific heat from high-temperature series expansions: Application to the two-dimensional XXZ Heisenberg model, Phys. Rev. B 104, 165113 (2021).
  * Derzhko  _et al._ [2020] O. Derzhko, T. Hutak, T. Krokhmalskii, J. Schnack, and J. Richter, Adapting Planck’s route to investigate the thermodynamics of the spin-half pyrochlore Heisenberg antiferromagnet, Phys. Rev. B 101, 174426 (2020).
  * Schäfer  _et al._ [2020] R. Schäfer, I. Hagymási, R. Moessner, and D. J. Luitz, Pyrochlore $S=\frac{1}{2}$ Heisenberg antiferromagnet at finite temperature, Phys. Rev. B 102, 054408 (2020).
  * Astrakhantsev  _et al._ [2021] N. Astrakhantsev, T. Westerhout, A. Tiwari, K. Choo, A. Chen, M. H. Fischer, G. Carleo, and T. Neupert, Broken-Symmetry Ground States of the Heisenberg Model on the Pyrochlore Lattice, Phys. Rev. X 11, 041021 (2021).
  * Hering  _et al._ [2022] M. Hering, V. Noculak, F. Ferrari, Y. Iqbal, and J. Reuther, Dimerization tendencies of the pyrochlore Heisenberg antiferromagnet: A functional renormalization group perspective, Phys. Rev. B 105, 054426 (2022).
  * Schmalfuß  _et al._ [2005] D. Schmalfuß, J. Richter, and D. Ihle, Green’s function theory of quasi-two-dimensional spin-half Heisenberg ferromagnets: Stacked square versus stacked kagomé lattices, Phys. Rev. B 72, 224405 (2005).
  * Guttmann [1987] A. J. Guttmann, On the critical behaviour of self-avoiding walks, Journal of Physics A: Mathematical and General 20, 1839 (1987).