We establish a polynomial-time approximation algorithm for partition functions of quantum spin models at high temperature. Our algorithm is based on the quantum cluster expansion of Netočný and Redig and the cluster expansion approach to designing algorithms due to Helmuth, Perkins, and Regts. Similar results have previously been obtained by related methods, and our main contribution is a simple and slightly sharper analysis for the case of pairwise interactions on bounded-degree graphs.

Classical algorithms for approximating partition functions of quantum models that make use of cluster expansions have occurred in two recent papers.1,2 In this paper, we provide a simple and concise exposition of how to construct such algorithms, with the intent of making the technique accessible to a wide audience.

A quantum spin system is modelled by a hypergraph G = (X, E). At each vertex x of G, there is a d-dimensional Hilbert space Hx with d < . The Hilbert space on the hypergraph is given by HGxXHx. An interaction Φ assigns a self-adjoint operator Φ(e) on HexeHx to each hyperedge e of G. The Hamiltonian on G is defined by HGeE(G)Φ(e). We are interested in the quantum partition functionZG(β) at inverse temperature β, defined by ZG(β)TreβHG.

In what follows, we shall focus our attention on quantum spin systems modeled by bounded-degree graphs; however, generalizations to bounded-degree bounded-rank hypergraphs are also possible. We shall assume that Φ(e)1 for every eE, where denotes the operator norm. Note that this is always possible by a rescaling of β. To state our main result, recall that a fully polynomial-time approximation scheme for a sequence of complex numbers (zn)nN is a deterministic algorithm that, for any n and ϵ > 0, produces a complex number ẑn such that znẑnϵzn in time polynomial in n and 1/ϵ.

Theorem 1.

FixΔZ+. There is a fully polynomial-time approximation scheme for the partition functionZG(β) for all graphsGof maximum degree at most Δ and all complex numbersβsuch thatβ1e4Δ.

Our algorithm is based on combining the abstract cluster expansion for quantum spin systems of Netočný and Redig3 with the algorithmic framework of Helmuth, Perkins, and Regts.4 The condition β=O1Δ is optimal under the complexity-theoretic assumption that RP (randomized polynomial time) is not equal to NP (non-deterministic polynomial time) due to the results on the hardness of approximate counting.5,6 We remark that these results concern real values of β; however, similar computational complexity transitions from P (polynomial time) to BQP-hard (bounded-error quantum polynomial time) and P to #P-hard can also be observed for complex values of β by using the methods of Refs. 7–9. For a formal definition of these complexity-theoretic notions, we refer the reader to Ref. 10.

Previous work on polynomial-time approximate counting algorithms for classical models have typically followed one of the three approaches: the correlation decay method,11,12 Markov-chain Monte Carlo,13 or interpolation-type methods.4,14 The latter two of these methods have also been used to design classical algorithms for quantum models.1,2,9,15–17 The goal of this paper is to convey the simplicity and flexibility of the third method that results from using the cluster expansion formalism. We emphasize that ideas of this type have previously been used to establish similar algorithms: for β(10e2Δ)1 with quasi-polynomial runtime1 and for β(16e3Δ)1 with polynomial runtime.2 Both the runtime of our algorithm and that of Ref. 2 are polynomials of a relatively high degree; examining our proof gives an upper bound of O(log(dΔ)) for the degree. While our results represent a modest improvement in the bound for β, we view our main contribution as being the simplicity of our analysis.

We note also that a priori information on the location of zeros of the partition function can be combined with the methods of this paper to develop polynomial-time algorithms. As noted in Ref. 4, this is an alternate route to the results of Patel and Regts18 using Barvinok’s method.14 For quasi-polynomial time results of this flavour in the quantum setting, see Ref. 1.

This paper is structured as follows: In Sec. II, we introduce the abstract cluster expansion. Then, in Sec. III, we show how the partition function of quantum spins systems admits such a cluster expansion. In Sec. IV, we use this framework to establish our approximation algorithm for the quantum partition function at high temperature. Finally, we conclude in Sec. V with some remarks and open problems.

The cluster expansion is a powerful tool from mathematical physics that allows one to express, via power series expansions, perturbations of a well-understood reference model. When the perturbations are sufficiently small, the power series expansions are convergent and allow one to draw many conclusions regarding correlation decay, zero-freeness, and other related properties. This method was originally introduced by Mayer in the study of imperfect gases19 but has since been greatly abstracted and simplified. The formulation we use is due to Kotecký and Preiss.20 

An abstract polymer model is a triple (C,w,), where C is a countable set whose objects are called polymers, w:CC is a function that assigns to each polymer γC a complex number wγ called the weight of the polymer, and ∼ is a symmetric compatibility relation such that each polymer is incompatible with itself. Equivalently, the incompatibility relation ≁ is a symmetric and reflexive relation. A set of polymers is called admissible if all the polymers in the set are all pairwise compatible. Note that the empty set is admissible. Let G denote the collection of all admissible sets of polymers from C. Then, the abstract polymer partition function is defined by

Our algorithm is based on reformulating the partition function of a quantum spin system in the abstract polymer model language (see Sec. III). The utility of this is due to the following fact about log(Z(C,w)).

Let Γ be a non-empty ordered tuple of polymers. The incompatibility graphHΓ of Γ is the graph with vertex set Γ and edges between any two polymers if and only if they are incompatible. Γ is called a cluster if its incompatibility graph HΓ is connected. Let GC denote the set of all clusters of polymers from C. The abstract cluster expansion20,21 is a formal power series for logZ(C,w) in the variables wγ, defined by

where φ(H) denotes the Ursell function of a graph H,

In this section, we shall show how the partition function of a quantum spin system admits an abstract polymer representation and hence an abstract cluster expansion. We return to the more general setting of hypergraphs for the remainder of this section.

Consider a quantum spin system modeled by the hypergraph G = (X, E) with interaction Φ, where at each vertex x of G, there is a d-dimensional Hilbert space Hx with d < . Recall that Φ assigns a self-adjoint operator Φ(e) on HexeHx to each hyperedge e of G. Define a polymer γ in this model to be a multiset (Eγ, mγ) of hyperedges EγE with the multiplicity function mγ:EγZ+ whose support Eγ induces a connected subgraph. Say that two polymers are compatible if and only if their supporting subgraphs are vertex disjoint. For a polymer γ, let γeEγmγ(e) denote its size and let γEγ denote the cardinality of its support; by a slight abuse of notation, we will write γ={γi}i=1γ. With these definitions, the partition function ZG(β) admits an abstract polymer model representation3,22,23 as formalized by the following lemma:

Lemma 2.
The partition functionZG(β) admits the following abstract polymer model representation:
where

Note that Sγ is the symmetric group of degree γ, and this assigns an ordering to the product of interactions. We prove Lemma 9 in  Appendix A. Note that the abstract polymer model representation holds as a formal power series in β. As an immediate corollary, we obtain a cluster expansion for log(ZG(β)).

Corollary 3.
The partition functionZG(β) admits the following cluster expansion:

For algorithms, an important quantity is the truncated cluster expansion for log(ZG(β)),

where ΓγΓγ.

We now establish our approximation algorithm. First, we show that the truncated cluster expansion provides a good approximation to log(ZG(β)). Netočný and Redig3 provided a sufficient condition for the convergence of the quantum cluster expansion based on the formalism of Kotecký and Preiss.20 In the following lemma, we follow their analysis in the setting of bounded-degree graphs. In particular, we obtain convergence criteria based on the maximum degree alone.

Lemma 4.
FixΔZ+. LetG = (V, E) be a graph of maximum degree at most Δ, and letβbe a complex number such thatβ1e4Δ. The cluster expansion for log(ZG(β)) converges absolutely,ZG(β) ≠ 0, and formZ+,

We prove Lemma 10 in  Appendix B. This lemma implies that to obtain an multiplicative ϵ-approximation to ZG(β), it is sufficient to compute Tm(ZG(β)) to order m=logV(G)/ϵ. We shall proceed by establishing an algorithm for computing Tm(ZG(β)) in time exp(O(m))V(G)O(1). Helmuth, Perkins, and Regts4 showed that such an algorithm exists, given the following three lemmas:

Lemma 5.

FixΔZ+, and letG = (V, E) be a graph of maximum degree at most Δ. The clusters of size at mostmcan be listed in timeexp(O(m))VO(1).

Proof.

Our proof follows that of Ref. 4, Theorem 6. First, we enumerate all connected subgraphs in G of size at most m. This can be achieved in time exp(O(m))VO(1) by Ref. 18, Lemma 3.6. Then, for each subgraph H, we enumerate all polymers (multisets) of size at most m whose corresponding subgraph in G is H. If H has size n, then there are precisely m1n1 of these, and they can be enumerated in time exp(O(m)). The enumeration of clusters in the claimed time then follows as in the proof of Ref. 4, Theorem 6.■

Lemma 6.

The Ursell functionφ(H) can be computed in timeexp(O(V(H))).

Proof.

This is a result of Ref. 24; see Ref. 4, Lemma 5.■

Lemma 7.

The weightwγof a polymerγcan be computed in timeexp(O(γ)).

We prove Lemma 11 in  Appendix B.

Lemma 8.

FixΔZ+, and letG = (V, E) be a graph of maximum degree at most Δ. The truncated cluster expansionTm(ZG(β)) can be computed in timeexp(O(m))V(G)O(1).

Proof.

We can list all clusters in G of size at most m in time exp(O(m))VO(1) by Lemma 5. For each of these clusters, we can compute the Ursell function in time exp(O(m)) by Lemma 6, and the polymer weights in time exp(O(m)) by Lemma 11. Hence, the truncated cluster expansion for log(ZG(β)) can be computed in time exp(O(m))V(G)O(1).■

Combining Lemmas 10 and 8 gives a fully polynomial-time approximation scheme for the partition function ZG(β) when G has maximum degree at most Δ and β is at most 1e4Δ. This proves Theorem 1.

We have discussed how classical algorithms based on cluster expansion methods apply to quantum spin systems at high temperature. Our focus has been on conveying the simplicity of the method, which has appeared previously in other forms.1,2 We note that it may be possible to use the Markov chain polymer approach of Ref. 25 to obtain an algorithm with an improved runtime.

For discrete classical spin systems, expansion methods have also been used at low temperatures, i.e., when β ≫ 1.4,26–29 It would be interesting to adapt these methods to quantum systems, e.g., by developing algorithms based on Pirogov–Sinai methods for quantum perturbations of classical systems.30 We remark, however, that it seems difficult to use this approach for low-temperature quantum systems with an infinite degeneracy of ground states, for example, when the set of ground states possesses a continuous symmetry.

We thank Michael Bremner, Adrian Chapman, Ashley Montanaro, and Will Perkins for helpful discussions. R.L.M. was supported by the QuantERA ERA-NET Cofund in Quantum Technologies implemented within the European Union’s Horizon 2020 Programme (QuantAlgo project) and EPSRC Grant Nos. EP/L021005/1 and EP/R043957/1. T.H. was supported by EPSRC Grant No. EP/P003656/1.

No new data were created during this study.

Lemma 2
(restatement). The partition functionZG(β) admits the following abstract polymer model representation:
where

Proof.
Our proof follows the work of Netočný and Redig.3 By definition,
We shall now rewrite the inner sum as a product over disjoint objects. Let S=(ei)i=1n denote any sequence of hyperedges, and let GS be the graph with vertex set [n] and edges between any two vertices i and j if and only if eiej. Define a sequential polymer to be a subsequence of S that corresponds to a maximally connected component of GS. Say that two sequential polymers are compatible if and only if their corresponding subgraphs in G are vertex disjoint. Let ΓS denote the set of all sequential polymers in S. It follows that
Let γ denote the length of a sequential polymer. Furthermore, let ΓG ≔ ∪SΓS denote the set of all sequential polymers in G, and let GG denote the collection of all admissible sets of sequential polymers in G. Observe that for any admissible set of sequential polymers {γi}i=1k, there are precisely i=1kγi!i=1kγi! sequences S that give rise to it. Thus, we may write
By interchanging the summations over n and k, we obtain
Now, by transforming the sum over admissible sets of sequential polymers into a sum over admissible sets of polymers and summing the weights of their permutations, we obtain
Since there are equivalent sequential polymers being distinguished in the sum over permutations, i.e., permutations of repeated hyperedges, we had to introduce a factor of 1eEγmγ(e)! to avoid overcounting. This completes the proof.■

Lemma 4.
(restatement). FixΔZ+. LetG = (V, E) be a graph of maximum degree at most Δ, and letβbe a complex number such thatβ1e4Δ. The cluster expansion for log(ZG(β)) converges absolutely,ZG(β) ≠ 0, and formZ+,

Proof.
For convenience, we normalize the trace so that Tr(I)=1. Note that this is equivalent to a rescaling of the partition function by a multiplicative factor. We introduce a polymer γx to every vertex x in G consisting of only that vertex. We define γx to be incompatible with every polymer that contains x. Then, we have
For a vertex x, there are at most (eΔ)n2 connected subgraphs with n edges that contain x (Ref. 31, Lemma 2.1). Furthermore, for such a subgraph, there are precisely k1n1 polymers (multisets) γ with γ=k that correspond to it. Thus, we may write
By interchanging the summations over n and k, we obtain
By taking β1e4Δ, we have
Fix a polymer γ. By summing over all vertices in γ, of which there are at most γ+1, we obtain
Now, by applying the main theorem of Ref. 20 with a(γ)=γ+1 and d(γ)=γ, we obtain that the cluster expansion converges absolutely in β, ZG(β) ≠ 0, and
We complete the proof by summing over all vertices in G,

Lemma 7

(restatement). The weightwγof a polymerγcan be computed in timeexp(O(γ)).

Proof.
For convenience, let n=γ and fix an enumeration γ1, …, γn of the multiset of edges in γ. By an inclusion–exclusion argument as in the derivation of Ryser’s formula for the permanent,32 we have
Thus, we may write
The first sum is over all subsets of [n], of which there are 2n. For each of these subsets A, we diagonalize the sum of the interactions ∑iAΦ(γi) to obtain the eigenvalues. This can be achieved in time exp(O(n)); here, we are using our assumption that the single-spin Hilbert spaces Hx are d-dimensional with d < . The trace may then be evaluated in time exp(O(n)) by evaluating the sum of the nth powers of the eigenvalues. Hence, the weight of a polymer can be computed in time exp(O(n))=exp(O(γ)).■

1.
A. W.
Harrow
,
S.
Mehraban
, and
M.
Soleimanifar
, in
Proceedings of the 52nd Annual ACM SIGACT Symposium on Theory of Computing
(
ACM
,
2020
), pp.
378
386
; arXiv:1910.09071.
2.
T.
Kuwahara
,
K.
Kato
, and
F. G.
Brandão
,
Phys. Rev. Lett.
124
,
220601
(
2020
); arXiv:1910.09425.
3.
K.
Netočný
and
F.
Redig
,
J. Stat. Phys.
117
,
521
(
2004
); arXiv:math-ph/0404018.
4.
T.
Helmuth
,
W.
Perkins
, and
G.
Regts
,
Probab. Theory Relat. Fields
176
,
851
(
2020
); arXiv:1806.11548.
5.
A.
Sly
and
N.
Sun
, in
53rd Annual IEEE Symposium on Foundations of Computer Science (FOCS)
(
IEEE
,
2012
), pp.
361
369
; arXiv:1203.2602.
6.
A.
Galanis
,
D.
Štefankovič
, and
E.
Vigoda
,
Combinatorics, Probab. Comput.
25
,
500
(
2016
); arXiv:1203.2226.
7.
M. J.
Bremner
,
R.
Jozsa
, and
D. J.
Shepherd
, in
Proceedings of the Royal Society of London A: Mathematical, Physical and Engineering Sciences
(
The Royal Society
,
2010
), p.
rspa20100301
; arXiv:1005.1407.
8.
L. A.
Goldberg
and
H.
Guo
,
Comput. Complexity
26
,
765
(
2017
); arXiv:1409.5627.
9.
R. L.
Mann
and
M. J.
Bremner
,
Quantum
3
,
162
(
2019
); arXiv:1806.11282.
10.
S.
Arora
and
B.
Barak
,
Computational Complexity: A Modern Approach
(
Cambridge University Press
,
2009
).
11.
D.
Weitz
, in
Proceedings of the Thirty-Eighth Annual ACM Symposium on Theory of Computing
(
ACM
,
2006
), pp.
140
149
.
12.
J.
Liu
,
A.
Sinclair
, and
P.
Srivastava
,
J. Stat. Phys.
174
,
287
(
2019
); arXiv:1704.06493.
13.
M.
Jerrum
and
A.
Sinclair
,
SIAM J. Comput.
22
,
1087
(
1993
).
14.
A.
Barvinok
,
Combinatorics and Complexity of Partition Functions
(
Springer
,
2016
), Vol. 274.
15.
S.
Bravyi
,
Quantum Inf. Comput.
15
,
1122
(
2015
); arXiv:1402.2295.
16.
S.
Bravyi
,
D.
Gosset
, and
R.
Movassagh
,
Nat. Phys.
(published online
2020
).
17.
E.
Crosson
and
S.
Slezak
, arXiv:2002.02232 (
2020
).
18.
V.
Patel
and
G.
Regts
,
SIAM J. Comput.
46
,
1893
(
2017
); arXiv:1607.01167.
19.
J.
Mayer
and
M.
Mayer
,
Statistical Mechanics
(
John Wiley & Sons, Inc.
,
1940
).
20.
R.
Kotecký
and
D.
Preiss
,
Commun. Math. Phys.
103
,
491
(
1986
).
21.
S.
Friedli
and
Y.
Velenik
,
Statistical Mechanics of Lattice Systems: A Concrete Mathematical Introduction
(
Cambridge University Press
,
2017
).
22.
Y. M.
Park
,
J. Stat. Phys.
27
,
553
(
1982
).
23.
D.
Ueltschi
, “
Discontinuous phase transitions in quantum lattice systems
,” Ph.D. thesis,
Verlag nicht ermittelbar
,
1998
.
24.
A.
Björklund
,
T.
Husfeldt
,
P.
Kaski
, and
M.
Koivisto
, in
49th Annual IEEE Symposium on Foundations of Computer Science
(
IEEE
,
2008
), pp.
677
686
; arXiv:0711.2585.
25.
Z.
Chen
,
A.
Galanis
,
L. A.
Goldberg
,
W.
Perkins
,
J.
Stewart
, and
E.
Vigoda
, in
Approximation, Randomization, and Combinatorial Optimization. Algorithms and Techniques (APPROX/RANDOM 2019)
(
Schloss Dagstuhl-Leibniz-Zentrum fuer Informatik
,
2019
); arXiv:1901.06653.
26.
M.
Jenssen
,
P.
Keevash
, and
W.
Perkins
, in
Proceedings of the Thirtieth Annual ACM-SIAM Symposium on Discrete Algorithms
(
SIAM
,
2019
), pp.
2235
2247
; arXiv:1807.04804.
27.
C.
Liao
,
J.
Lin
,
P.
Lu
, and
Z.
Mao
, in
Approximation, Randomization, and Combinatorial Optimization. Algorithms and Techniques (APPROX/RANDOM 2019)
(
Schloss Dagstuhl-Leibniz-Zentrum fuer Informatik
,
2019
); arXiv:1903.07531.
28.
C.
Borgs
,
J.
Chayes
,
T.
Helmuth
,
W.
Perkins
, and
P.
Tetali
, “
Efficient sampling and counting algorithms for the Potts model on Zd at all temperatures
,” in
Proceedings of the 52nd Annual ACM SIGACT Symposium on Theory of Computing
(
ACM
,
2020
) pp.
738
751
, .
29.
C.
Carlson
,
E.
Davies
, and
A.
Kolla
, arXiv:2003.01154 (
2020
).
30.
C.
Borgs
,
R.
Kotecký
, and
D.
Ueltschi
,
Commun. Math. Phys.
181
,
409
(
1996
).
31.
C.
Borgs
,
J.
Chayes
,
J.
Kahn
, and
L.
Lovász
,
Random Struct. Algorithms
42
,
1
(
2013
); arXiv:1002.0115.
32.
H. J.
Ryser
,
Combinatorial Mathematics
(
Mathematical Association of America
,
1963
), Vol. 14.