# A Faster Algorithm for the Free Energy in One-Dimensional Quantum Systems

Samuel O. Scalet Department of Applied Mathematics and Theoretical Physics,
University of Cambridge, United Kingdom<EMAIL_ADDRESS>
###### Abstract

We consider the problem of approximating the free energy density of a
translation-invariant, one-dimensional quantum spin system with finite range.
While the complexity of this problem is nontrivial due to its close connection
to problems with known hardness results, a classical subpolynomial-time
algorithm has recently been proposed [Fawzi et al., 2022]. Combining several
algorithmic techniques previously used for related problems, we propose an
algorithm outperforming this result asymptotically and give rigorous bounds on
its runtime. Our main techniques are the use of Araki expansionals, known from
results on the nonexistence of phase transitions, and a matrix product
operator construction. We also review a related approach using the Quantum
Belief Propagation [Kuwahara et al., 2018], which in combination with our
findings yields an equivalent result.

##  1 Introduction

Making predictions of physical quantities in thermal equilibrium for quantum
many-body systems is a central goal of condensed matter physics. Despite a
plethora of succesful algorithms succeeding this task in many settings, the
field of Hamiltonian complexity has also shown obstacles to this task that
cannot be overcome. Thereby, algorithms are usually either heuristic or
specialized to settings in which no complexity theoretic barrier rules out the
existence of efficient algorithms.

The free energy $F=E-TS$ is an example of such a quantity and has been
addressed in many previous works on both, constructive results designing
algorithms and complexity theoretic reductions proving hardness. The interest
in the free energy stems from both its immediate physical interest as the
value of a functional which variationally characterizes the Gibbs state and
also its proximal use in calculations of observables since their computation
can be reduced to the computation of derivatives of the free energy. Two
examples of regimes where provably efficient algorithms exist are lattice
systems at high temperature [17] and one dimensional systems at arbitrary
nonzero temperature [15, 10]. Remarkably, the computational efficiency is
intimately linked to the nonexistence of thermal phase transitions, a
connection that has also formally been established [12].

In this work, we focus on the one-dimensional setting and further specialize
to the case of translation-invariant systems in the thermodynamic limit, i.e.,
infinite system size, which is arguably the physically most relevant. Despite
the reduction of degrees of freedom in the problem input, the existence of
efficient algorithms is nontrivial. The hardness of a related problem, the
translation-invariant 1D ground-state problem, which arises as the zero
temperature limit of the free energy has been proven in [11, 5], see also the
discussion in [10] on the relation to hardness of the free energy density
problem.

Many algorithms in the 1D case are formulated for finite systems. Since the
concept of thermal equilibrium is derived from an infinte system limit, these
algorithms will usually be used for large system sizes to approximate the
thermodynamic limit. The error in approximating a system in the thermodynamic
limit by a finite one is standardly bounded by the inverse of the system size
and consequently this approach does not allow to perform better than linearly
in the inverse error. In the recent work [10], the authors overcome this
barrier and introduce a new algorithmic approach that directly treats an
infinite system size. The resulting algorithm runs in subpolynomial time in
the inverse error, i.e., in time less than the inverse error to any constant
power greater zero.

#### Setup and main result

Let us introduce the setup of one-dimensional translation-invariant spin
systems in the thermodynamic limit. The local Hilbert space is
$\mathbb{C}^{d}$ and the system is defined via a nearest neighbour interaction
$h\in\mathcal{B}(\mathbb{C}^{d}\otimes\mathbb{C}^{d})$ with $\|h\|\leq 1$. Let
us define the Hamiltonian for a finite chain of length $N$ as

$H_{[1,N]}=\sum_{i=1}^{N-1}h_{i,i+1},$

where $h_{i,i+1}$ has support on sites $i,i+1$. The case of higher range
interactions can be treated by appropriately blocking the system.

The free energy density for a given inverse temperature $\beta$ is defined as

$f(h)=-\lim_{N\to\infty}\frac{1}{\beta
N}\log(Z_{N})=-\lim_{N\to\infty}\frac{1}{\beta N}\log(\Tr\left[e^{-\beta
H_{[1,N]}}\right]).$ (1)

Our main result is the following subpolynomial-time algorithm for
approximating the free energy density. It improves on the asymptotics in the
prior work [10], which only gives a runtime of
$\exp(\mathcal{O}(\log(1/\varepsilon)/\log(\log(1/\varepsilon))))$.

######  Theorem 1.

There exists a classical algorithm taking as input a hermitian matrix
$h\in\mathbb{C}^{d^{2}\times d^{2}}$ and a precision $1>\varepsilon>0$ and
running in time

$\exp(\sqrt{A\log(\frac{B}{\varepsilon})}\log(A\log(\frac{B}{\varepsilon}))),$

where $A,B$ are constants depending on the inverse temperature and local
dimension only, outputting an approximation $\tilde{f}$ to the one-dimensional
translation-invariant free energy density $|f-\tilde{f}|\leq\varepsilon$.

#### Overview of Approach

We give an overview of our approach to compute the free energy density. It
breaks down into two steps: We first relate the global free energy density to
an increment of partition functions, an idea inspired by [15]. Secondly, we
use a matrix product operator (MPO) construction from [14] to compute these
partition functions individually. Combining these error estimates gives our
main result.

Compared to [15], we give a strongly simplified proof of the first step based
on Araki expansionals [4, 18].

######  Lemma 1.

There exist constants $A,\xi>0$ depending on $\beta$ only such that

$\left|e^{-\beta f(h)}-\frac{Z_{l+1}}{Z_{l}}\right|\leq Ae^{-l/\xi}$

A non-quantitative version has been proven in [10]. Our approach is closely
related but provides an explicit convergence rate. While the locality of Araki
expansionals holds up to _superexponentially_ small errors, our proof involves
a result on local indistinguishability that only allows for exponential
convergence. Despite this superexponential decay being essential to obtain a
_subpolynomial_ runtime in [10], we still achieve this property due to the
subexponential runtime in system size and inverse error of the MPO algorithm.

The proof in [15] instead uses a tool called quantum belief propagation [13,
8]. At first sight its range of applicability is much wider allowing for non-
translation-invariant and higher dimensional systems and also has a better
dependence on the temperature. However, the proof is conditional on an
additional ingredient, a result on the decay of correlations. While decay of
correlations has been proven for one-dimensional translation-invariant
infinite [4] and, only recently, finite [6] systems, the validity of the
condition in [15] is not entirely clear as it is needed for a _perturbed_ and
finite Gibbs state. To clarify these connections, we review this approach in
Appendix A.1 giving a slightly modified proof that replaces the decay of
correlations entirely by a local indistinguishability result that has been
proven rigorously [18]. While not the focus of this manuscript, we note that
for higher dimensional systems, the approach based on quantum belief
propagation still works and, when further restricting to the high-temperature
regime, the decay of correlations has been proven as well, even for non-
translation-invariant Gibbs states. It should be noted that in this regime,
the algorithm yields only a quasipolynomial runtime, which falls behind the
state of the art cluster expansion techniques [17], running in only polynomial
time.

We also want to comment on an alternative way of reducing an algorithm for
local observables, such as the MPO subroutine we use, to the free energy
problem. The equivalence of free energy and local observable expectation
values is folklore knowledge due to a relation known as the derivative trick,
giving expectation values of local observables as derivatives of the free
energy. A precise result can be found in [7], which inverts this relation
using a discretized integration of the derivative of the free energy.
Suffering a polynomial overhead in the inverse error, this approach does not
recover a subpolynomial runtime. We note, however, the interesting analogy
that an iteration over a discretized temperature interval in [7] and the
iteration over sites in [15], both yield a product formula for the partition
function.

#### Practical aspects

Many algorithms with provable runtime bounds, including ours, come with a
number of parameters and explicit bounds on those that guarantee the desired
accuracy. These bounds are often impractically large to actually run the
provably correct algorithm, but choosing the parameters within the limits of
the computational resources available provides an immediate heurestic approach
that often gives good approximations. To this end, it is helpful to limit the
number of these parameters, so that this heuristic choice of parameters
becomes more feasible.

In our case, the involved quantities are the length scale in Lemma 1, as well
as two parameters hidden inside the MPO part of the algorithm, a length scale
and a degree of polynomial approximation. It is a major theoretical
contribution in [14], that these parameters yield a polynomial bond dimension
_without the need for any truncation_. Such truncations are common in the
literature on tensor network algorithms, but their impact on trace norm errors
is uncontrolled111However, the 2-norm error can be bounded, which renders
these truncations useful also theoretically for pure states.. We note however,
that they have been demonstrated to be succesful in practice [9] and could be
integrated into the algorithm proposed here for a heuristic variant. We note
that in [10], an alternative algorithm for the problem we consider, the
practical implementation is much simpler since it contains only a single
parameter which can be simply chosen as large as allowed by the computational
resources and this has been demonstrated to yield good convergence results in
practice. Similarly, the original approach in [15], calculating the partition
functions by diagonalization, would yield a straightforward practical
implementation but comes at the cost of not having a theoretical subpolynomial
runtime guarantee.

We finally remark, that while our algorithm could be reduced to the local
observable problem (see [10, Appendix B]), this can already be achieved using
the MPO subroutine on its own directly. An extension to marginals of larger
subsytems has been put forward in [2], which could also be extended to the
thermodynamic limit.

##  2 Notation

We consider a multipartite spin system on a chain $\mathbb{Z}$, where to each
site we associate a local Hilbert space $\mathcal{H}_{i}\cong\mathbb{C}^{d}$
and the Hilbert space on an interval $[k,l]$ is given by
$\mathcal{H}_{[k,l]}=\bigotimes_{i=k}^{l}\mathcal{H}_{i}$. As introduced
earlier, we consider a two-site interaction
$h\in\mathcal{B}(\mathbb{C}^{d}\otimes\mathbb{C}^{d})$ and denote by
$h_{i,i+1}$ this interaction with support on sites $i,i+1$. The Hamiltonian
for a given interval is then $H_{[k,l]}=\sum_{i=k}^{l-1}h_{i,i+1}$. We denote
by $\Tr[\cdot]$, $\tr[\cdot]=\tr[\cdot]/d^{N}$ the unnormalized and normalized
respectively. For notational convenience the trace will always be considered
to be defined for the space on which the operators act. We denote the
partition function by $Z_{N}=\Tr[\exp(-\beta H_{N})]$ keeping the inverse
temperature implicit. As before the free energy density is given as
$f(h)=-1/\beta\lim_{N\to\infty}\log(Z_{N})/N$. We introduce the notation

$\rho_{[k,l]}=\frac{e^{-\beta H_{[k,l]}}}{\Tr[e^{-\beta H_{[k,l]}}]}$

for the _Gibbs state_ on a finite interval .

We introduce a matrix product operator (MPO) on a subsystem
$\mathcal{H}_{[k,l]}$ as a set of $(k-l+1)d^{2}$ $D\times D$ matrices
$M_{r,s}^{i}$ with $1\leq r,s\leq d$, $k\leq i\leq l$ defining a state
$\rho[M]$ as

$\Tr[\rho[M]\bigotimes_{i=k}^{l}\ket{r_{i}}\bra{s_{i}}]=\Tr[\prod_{i=k}^{l}M_{r_{i},s_{i}}^{i}].$

While MPOs are not used directly in our work, note that the above formula
implies directly that the trace of an MPO $\Tr[\rho[M]]$ can be computed in
time polynomial in $(k-l,D,d)$.

##  3 Reduction to local partition functions

In this section we prove Lemma 1. We start by introducing some definitions and
results regarding Araki expansionals [4]. The main result of this work is a
proof of the analyticity of the free energy, for translation-invariant one-
dimensional finite-range systems. This turns out to be closely related to
various results on the decay of correlations as well as the computational
complexity of those systems. The main object of interest is the non-hermitian
operator

$E_{n}=e^{-\beta H_{[1,n]}}e^{\beta H_{[2,n]}}.$

Due the noncommutativity of the terms in the Hamiltonian, this operator is in
general supported on the whole interval $[1,n]$. The following result shows,
however, that it can be well approximated by a local operator and is uniformly
bounded in $n$. Note, that the result holds at arbitrary temperature in the
case of finite-range interactions. The version of the result we cite has been
given for systems with exponentially decaying interactions above a threshold
temperature, but extend to all temperatures in the finite-range setting.

######  Lemma 2 ([18, Proposition 4.2]).

There is a constant $\mathcal{G}$ and a superexponentially decaying function
$\delta(\cdot)$ only depending on the temperature, but not on $n$ such that
for all $a\geq n$

  1. 1.

$\|E_{n}\|,\|E_{n}^{-1}\|\leq\mathcal{G}$

  2. 2.

$\|E_{n}-E_{a}\|,\|E_{n}^{-1}-E_{a}^{-1}\|\leq\delta(n)$

The same techniques are crucial in proving a _local indistinguishability_
result, which we will also require in our analysis.

######  Lemma 3 ([18, Lemma 4.15]).

There exist constants $\mathcal{K},\alpha>0$ only depending on $\beta$ such
that for any $n^{\prime}\geq n$ and any observable $A$ supported on $[1,k]$,
we have

$\left|\Tr[A\rho_{[1,k+n]}]-\Tr[A\rho_{[1,k+n^{\prime}]}]\right|\leq\|A\|\mathcal{K}e^{-\alpha
n}.$

We now proceed to the proof of Lemma 1, a quantitative version of [10, Lemma
4.1] in an analogous fashion.

###### Proof of Lemma 1.

We decompose the ratio of partition functions using Araki expansionals

$\displaystyle\frac{Z_{m+1}}{Z_{m}}=\frac{\Tr[e^{-\beta
H_{[1,m+1]}}]}{\Tr[e^{-\beta H_{[1,m]}}]}=\frac{\Tr[E_{m+1}e^{-\beta
H_{[2,m+1]}}]}{\Tr[e^{-\beta H_{[1,m]}}]}=\frac{\Tr[\tr_{1}[E_{m+1}]e^{-\beta
H_{[1,m]}}]}{\Tr[e^{-\beta H_{[1,m]}}]}$

Subtracting the same decomposition for $Z_{l+1}/Z_{l}$ choosing $m\geq l$, we
find

$\displaystyle\left|\frac{Z_{m+1}}{Z_{m}}-\frac{Z_{l+1}}{Z_{l}}\right|$
$\displaystyle\leq\left|\Tr[\tr_{1}[E_{l/2}-E_{l+1}]\frac{e^{-\beta
H_{[1,l]}}}{\Tr[e^{-\beta H_{[1,l]}}]}]\right|$ (2)
$\displaystyle\quad+\left|\Tr[\tr_{1}[E_{l/2}]\left(\frac{e^{-\beta
H_{[1,m]}}}{\Tr[e^{-\beta H_{[1,m]}}]}-\frac{e^{-\beta
H_{[1,l]}}}{\Tr[e^{-\beta H_{[1,l]}}]}\right)]\right|$
$\displaystyle\quad+\left|\Tr[\tr_{1}[E_{m+1}-E_{l/2}]\frac{e^{-\beta
H_{[1,m]}}}{\Tr[e^{-\beta H_{[1,m]}}]}]\right|$
$\displaystyle\leq\delta(l/2+1)+\mathcal{G}\mathcal{K}e^{-\alpha
l}+\delta(l/2+1)$ $\displaystyle\leq Ae^{-l/\xi},$

where we used Lemma 2 and Lemma 3. The last line follows from upper bounding
the superexponential decay by an exponentially decaying function and combining
all constants. Note that

$f(h)=-\lim_{m\to\infty}\frac{1}{\beta
m}\log(Z_{m})=-\lim_{m\to\infty}\frac{1}{\beta
m}\sum_{i=1}^{m}\log(\frac{Z_{i}}{Z_{i-1}}),$

where $Z_{0}=1$, which implies

$f(h)=-\frac{1}{\beta}\log(\lim_{m\to\infty}\frac{Z_{m+1}}{Z_{m}})$

if the limit exists and is positive. The existence of the limit follows from
its Cauchy convergence that holds due to Equation (2) and its positivity
follows from [16, Lemma 3.5]

$\left|\log(\frac{Z_{m+1}}{Z_{m}})\right|\leq\beta\|h\|.$

Taking the limit of 2 concludes the proof. ∎

##  4 MPO algorithm

With this expression at hand we can proceed to the proof of our main result.
Deviating from the approach in [15] where the analogue expression from the
previous lemma is simply evaluated using exact diagonalization, we use the
following algorithm from [14] to construct a matrix product operator
approximation of the unnormalized Gibbs states and compute their traces.

######  Theorem 2 ([14, Theorem 3]).

There exists an algorithm taking as input the description of an input
hermitian matrix $h\in\mathbb{C}^{d^{2}\times d^{2}}$ with $\|h\|\leq 1$ and
an inverse temperature $\beta$ and outputting a matrix product operator $M$ on
$n$ sites such that

$\left\|M-e^{-\beta H_{[1,n]}}\right\|_{1}\leq\varepsilon\left\|e^{-\beta
H_{[1,n]}}\right\|_{1}.$

The bond dimension is given by

$D\leq\exp(C\sqrt{\beta\log(n/\varepsilon)}\log(\beta\log(n/\varepsilon))),$

and the algorithm runs in time

$t\leq
n\beta\exp(C\sqrt{\beta\log(n/\varepsilon)}\log(\beta\log(n/\varepsilon))),$

where $C$ is a numerical constant.

######  Remark 1.

We note that the cited work contains a case distinction for small and large
$\varepsilon$ compared to $\beta$. Since we always consider $\beta$ as a
parameter and only write the dependence on $\varepsilon$, we choose the case
of small $\varepsilon$. This threshold can be absorbed into the constant $B$
in Theorem 1.

The Theorem comes with an explicit temperature dependence, which we included
in the above formulation. This dependence will not carry over to our result
since the bottleneck is given by the approximation step using Araki
expansionals.

Combining Lemma 1 with this MPO algorithm we prove our main result.

###### Proof of Theorem 1.

The algorithm is given in Algorithm 1.

Parameter: inverse temperature $\beta$

Require: subroutine ${MPO}_{\textrm{Gibbs}}$ constructing an MPO approximation
of the unnormalized thermal state from Theorem 2, polynomial time routine to
compute the trace of a given MPO $\Tr_{MPO}$ (standard contraction), explicit
constants $D,E$ (see below)

Input: local dimension $d$, Hamiltonian term $h\in\mathbb{C}^{d^{2}\times
d^{2}}$ such that $\|h\|\leq 1$, error $\varepsilon$

Output: $\tilde{f}$ approximation to the free energy $f(h)$

1 $l\leftarrow C\log(1/\varepsilon)$;

2
$\tilde{Z}_{l}\leftarrow\Tr_{MPO}[{MPO}_{\textrm{Gibbs}}(l,h,\varepsilon\beta\exp(-\beta(2l+1))/12)]$;

3
$\tilde{Z}_{l+1}\leftarrow\Tr_{MPO}[{MPO}_{\textrm{Gibbs}}(l+1,h,\varepsilon\beta\exp(-\beta(2l+1))/8)]$;

4 $\tilde{f}_{\beta}\leftarrow-(1/\beta)\log(\tilde{Z}_{l+1}/\tilde{Z}_{l})$;

Algorithm 1 Algorithm for computing the free energy density.

Since all expressions here are ratios of traces, we can switch to normalized
traces $\tr[\cdot]=\Tr[\cdot]/d^{l}$. Recall, that we denote the partition
functions for the system of size $l$ by $Z_{l}$. We further denote by
$\tilde{Z}_{l}$, $\tilde{Z}_{l+1}$ the corresponding approximation from the
MPO algorithm, i.e., we run the the algorithm from Theorem 2 at an error level
$\varepsilon_{MPO}$ to be chosen later and contract the resulting tensor
network to obtain a trace estimate. The latter step can be done in time
polynomial in the bond dimension, physical dimension and $l$.

Our task is to provide explicit bounds on the parameters $l$ and
$\varepsilon_{MPO}$ to obtain the desired accuracy.

The approximation error can be split up into a term caused by the truncation
in $l$ and one term due to the approximation errors from the MPO

$\left|f(h)-\tilde{f}\right|\leq\left|f(h)-\left(-\frac{1}{\beta}\right)\log(\frac{Z_{l+1}}{Z_{l}})\right|+\left|-\frac{1}{\beta}\log(\frac{Z_{l+1}}{Z_{l}})-\left(-\frac{1}{\beta}\right)\log(\frac{\tilde{Z}_{l+1}}{\tilde{Z}_{l}})\right|$
(3)

In order to control the error propagation for the logarithm and fractions, we
need analytic bounds on the partition functions and their ratios. The absolute
values of the traces can be bounded by

$e^{-\beta l}\leq Z_{l}\leq e^{\beta l}.$

Furthermore, from [16, Lemma 3.6] we have

$\left|\log(\frac{Z_{l+1}}{Z_{l}})\right|\leq\beta.$

For the first term in Equation (3) by the mean value theorem we estimate

$\left|-\frac{1}{\beta}\log(e^{-\beta
f})-\left(-\frac{1}{\beta}\right)\log(\frac{Z_{l+1}}{Z_{l}})\right|\leq\frac{1}{\beta}e^{\beta}\left|\frac{Z_{l+1}}{Z_{l}}-e^{-\beta
f}\right|$

using the derivative of the logarithm and a lower bound of $e^{-\beta}$ on its
argument for all points in the intervall spanned by $e^{-\beta f}$ and
$Z_{l+1}/Z_{l}$. Invoking Lemma 1, we can bound the first term by
$\varepsilon/2$ by choosing

$l\geq\log(\frac{2e^{\beta}A}{\beta\varepsilon})\xi.$

For the second term in Equation (3) we proceed similarly. To lower bound the
derivative of the logarithm we require the error bound

$\left|\frac{Z_{l+1}}{Z_{l}}-\frac{\tilde{Z}_{l+1}}{\tilde{Z}_{l}}\right|\leq\frac{1}{2}e^{-\beta}$
(4)

which will follow from the input accuracy to the MPO subroutine. Analogously
to before, this upper bounds the derivative of the logarithm by $2e^{\beta}$
in the interval with end points $Z_{l+1}/Z_{l}$,
$\tilde{Z}_{l+1}/\tilde{Z}_{l}$. Under this condition a bound on the second
term in Equation (3) reads

$\left|-\frac{1}{\beta}\log(\frac{Z_{l+1}}{Z_{l}})-\left(-\frac{1}{\beta}\right)\log(\frac{\tilde{Z}_{l+1}}{\tilde{Z}_{l}})\right|\leq\frac{1}{\beta}2e^{\beta}\left|\frac{Z_{l+1}}{Z_{l}}-\frac{\tilde{Z}_{l+1}}{\tilde{Z}_{l}}\right|\leq\frac{\varepsilon}{2},$

i.e., we need an error of

$\left|\frac{Z_{l+1}}{Z_{l}}-\frac{\tilde{Z}_{l+1}}{\tilde{Z}_{l}}\right|\leq\frac{\varepsilon\beta
e^{-\beta}}{4}$ (5)

for the approximate ratio of partition functions. Note that this is a stronger
requirement than Equation (4) if $\varepsilon<\beta/2$ so, since we consider
only the asymptotics of small $\varepsilon$ for fixed $\beta$, we consider the
latter condition included.

Using the inequality

$\left|\frac{Z_{l+1}}{Z_{l}}-\frac{\tilde{Z}_{l+1}}{\tilde{Z}_{l}}\right|\leq\frac{1}{Z_{l}}\left|Z_{l+1}-\tilde{Z}_{l+1}\right|+\frac{\tilde{Z}_{l+1}}{Z_{l}\tilde{Z}_{l}}\left|\tilde{Z}_{l}-Z_{l}\right|,$

and recalling

$\displaystyle Z_{l}^{-1}$ $\displaystyle\leq e^{\beta l}$
$\displaystyle\frac{\tilde{Z}_{l+1}}{\tilde{Z}_{l}}$
$\displaystyle\leq\frac{Z_{l+1}}{Z_{l}}-\left|\frac{Z_{l+1}}{Z_{l}}-\frac{\tilde{Z}_{l+1}}{\tilde{Z}_{l}}\right|\leq\frac{3}{2}e^{\beta}$

we see that by running the algorithm in Theorem 2 at error levels

$\displaystyle\left|Z_{l}-\tilde{Z}_{l}\right|$
$\displaystyle\leq\varepsilon_{MPO,l}e^{\beta(l-1)}:=\frac{\varepsilon\beta}{12}e^{-\beta(l+2)}$
$\displaystyle\left|Z_{l+1}-\tilde{Z}_{l+1}\right|$
$\displaystyle\leq\varepsilon_{MPO,l+1}e^{\beta
l}:=\frac{\varepsilon\beta}{8}e^{-\beta(l+1)}$

equation (5) follows. This proves correctness of the algorithm. Note that that
the inverse input errors for the MPO algorithm are polynomial in
$\varepsilon^{-1}$ as $l$ has only logarithmic dependence.

Let us now turn our attention to the runtime of the algorithm. The relevant
part here is the runtime of the algorithm constructing the MPO descriptions of
the truncated thermal states and its resulting bond dimension, which
polynomially enters the runtime of the contraction. The remaining steps of the
algorithm are evaluations of elementary functions and their contribution to
the runtime can be neglected.

From Theorem 2 runtime and bond dimension are given by

$t\leq
l\beta\exp(C\sqrt{\beta\log(l/\varepsilon_{MPO})}\log(\beta\log(l/\varepsilon_{MPO})))$

and

$D=\exp(C\sqrt{\beta\log(l/\varepsilon_{MPO})}\log(\beta\log(l/\varepsilon_{MPO})))$

$\varepsilon_{MPO}$ is the input error for the MPO construction algorithm. The
runtime of the contraction is $\operatorname{poly}(l,D,d)$. The temperature
dependence of the runtimes can be simplified by combining all constants into
constants $\tilde{A},\tilde{B}$ and we obtain an overall bound on the runtime

$t_{\textrm{total}}\leq\exp(\sqrt{\tilde{A}\log(\frac{\tilde{B}}{\varepsilon})}\log(\tilde{A}\log(\frac{\tilde{B}}{\varepsilon}))),$

which closes the proof.

∎

#### Acknowledgements

I would like to thank Álvaro Alhambra and Hamza Fawzi for helpful discussions.
I acknowledge support from the UK Engineering and Physical Sciences Research
Council (EPSRC) under grant number EP/W524141/1.

## References

  * [1] Álvaro M. Alhambra  “Quantum Many-Body Systems in Thermal Equilibrium”  In _PRX Quantum_ 4 American Physical Society, 2023, pp. 040201  DOI: 10.1103/PRXQuantum.4.040201
  * [2] Álvaro M. Alhambra and J. Cirac  “Locally accurate tensor networks for thermal states and time evolution”  In _PRX Quantum_ 2, 2021, pp. 040331  DOI: 10.1103/PRXQuantum.2.040331
  * [3] Anurag Anshu, Srinivasan Arunachalam, Tomotaka Kuwahara and Mehdi Soleimanifar  “Sample-efficient learning of interacting quantum systems”  In _Nature Physics_ 17.8 Springer ScienceBusiness Media LLC, 2021, pp. 931–935  DOI: 10.1038/s41567-021-01232-0
  * [4] Huzihiro Araki  “Gibbs states of a one dimensional quantum lattice”  In _Communications in Mathematical Physics_ 14.2 Springer ScienceBusiness Media LLC, 1969, pp. 120–157  DOI: 10.1007/bf01645134
  * [5] Johannes Bausch, Toby Cubitt and Maris Ozols  “The Complexity of Translationally Invariant Spin Chains with Low Local Dimension” Publisher: Springer Science and Business Media LLC  In _Annales Henri Poincaré_ 18.11, 2017, pp. 3449–3513  DOI: 10.1007/s00023-017-0609-7
  * [6] Andreas Bluhm, Ángela Capel and Antonio Pérez-Hernández  “Exponential decay of mutual information for Gibbs states of local Hamiltonians”  In _Quantum_ 6 Verein zur Forderung des Open Access Publizierens in den Quantenwissenschaften, 2022, pp. 650  DOI: 10.22331/q-2022-02-10-650
  * [7] Sergey Bravyi, Anirban Chowdhury, David Gosset and Pawel Wocjan  “On the complexity of quantum partition functions”  In _Nature Physics_ 18.11, 2022, pp. 1367–1370  DOI: 10.1038/s41567-022-01742-5
  * [8] Ángela Capel, Massimo Moscolari, Stefan Teufel and Tom Wessel  “From decay of correlations to locality and stability of the Gibbs state”  arXiv, 2023  arXiv: http://arxiv.org/abs/2310.09182
  * [9] Bin-Bin Chen, Lei Chen, Ziyu Chen, Wei Li and Andreas Weichselbaum  “Exponential Thermal Tensor Network Approach for Quantum Lattice Models”  In _Phys. Rev. X_ 8 American Physical Society, 2018, pp. 031082  DOI: 10.1103/PhysRevX.8.031082
  * [10] Hamza Fawzi, Omar Fawzi and Samuel O. Scalet  “A subpolynomial-time algorithm for the free energy of one-dimensional quantum systems in the thermodynamic limit”  In _Quantum_ 7, 2023, pp. 1011  DOI: 10.22331/q-2023-05-22-1011
  * [11] Daniel Gottesman and Sandy Irani  “The Quantum and Classical Complexity of Translationally Invariant Tiling and Hamiltonian Problems” Publisher: Theory of Computing  In _Theory of Computing_ 9.2, 2013, pp. 31–116  DOI: 10.4086/toc.2013.v009a002
  * [12] Aram W. Harrow, Saeed Mehraban and Mehdi Soleimanifar  “Classical algorithms, correlation decay, and complex zeros of partition functions of Quantum many-body systems”  In _Proceedings of the 52nd Annual ACM SIGACT Symposium on Theory of Computing_ ACM, 2020  DOI: 10.1145/3357713.3384322
  * [13] M.. Hastings  “Quantum belief propagation: An algorithm for thermal quantum systems”  In _Physical Review B_ 76.20 American Physical Society (APS), 2007  DOI: 10.1103/physrevb.76.201102
  * [14] Tomotaka Kuwahara, Álvaro M. Alhambra and Anurag Anshu  “Improved Thermal Area Law and Quasilinear Time Algorithm for Quantum Gibbs States”  In _Physical Review X_ 11.1 American Physical Society (APS), 2021  DOI: 10.1103/physrevx.11.011047
  * [15] Tomotaka Kuwahara and Keiji Saito  “Polynomial-time Classical Simulation for One-dimensional Quantum Gibbs States”, 2018  arXiv:1807.08424 [quant-ph]
  * [16] Marco Lenci and Luc Rey-Bellet  “Large deviations in quantum lattice systems: One-phase region” Publisher: Springer Science and Business Media LLC  In _Journal of Statistical Physics_ 119.3-4, 2005, pp. 715–746  DOI: 10.1007/s10955-023-03193-5
  * [17] Ryan L. Mann and Tyler Helmuth  “Efficient algorithms for approximating quantum partition functions”  In _Journal of Mathematical Physics_ 62.2, 2021, pp. 022201  DOI: 10.1063/5.0013689
  * [18] David Pérez-García and Antonio Pérez-Hernández  “Locality Estimates for Complex Time Evolution in 1D”  In _Communications in Mathematical Physics_ 399.2, 2023, pp. 929–970  DOI: 10.1007/s00220-022-04573-w

##  Appendix A Quantum Belief Propagation

In this Appendix we give an alternative proof of Lemma 1 based on [15]. We
deviate in a few technical details, in particular, we do not need a result on
the decay of correlations. While the decay of correlations has been proven for
translation-invariant infinite [4], and translation-invariant finite [6] Gibbs
states in one dimension, the requirement in [15] is for a finite but perturbed
Gibbs state. We believe that this gap could be filled, however, our approach
circumvents this issue by making use of a known local indistinguishability
result instead.

The significantly simpler proof in the main text is restricted to the
translation-invariant one-dimensional setting. The Quantum Belief Propagation
(QBP) used in this chapter is a general perturbation formula that applies to
arbitrary quantum systems satisfying Lieb-Robinson bounds. However, to
conclude a result equivalent to Lemma 1, we require the additional property of
local indistinguishability that is satisfied, i.e., for systems at high
temperature, but usually also comes with alternative algorithmic approaches
such as cluster expansion [17].

This section is structured as follows. We start by giving the exact QBP
formula for a perturbation equivalent to adding sites in the partition
function. The formula gives a perturbation operator that relates a perturbed
unnormalized Gibbs state to the original one. The second step is to prove
locality properties of these perturbation operators. Here we slightly deviate
from previous formulations and give an approximation in terms of truncated
Hamiltonians rather than partial traces of the perturbation operators. We
believe this formulation could also be convenient for further work using the
QBP. We conclude by combining this formula with the local indistinguishability
in Lemma 3.

###  A.1 Exact Quantum Belief Propagation

The Quantum Belief Propagation is originally due to [13]. It gives an analytic
expression for a perturbed (unnormalized) Gibbs state. We will use this result
to express the ratio of partition functions with an additional site as of
Lemma 1.

We make use of the time evolution operator for a given Hamiltonian

$\Gamma^{t}_{H_{0}}(A)=e^{itH_{0}}Ae^{-itH_{0}}$

and start by giving the main definitions necessary for the QBP.

######  Definition 1 (Exact Quantum Belief Propagation).

The following equations define the Quantum Belief Propagation Operators for
any perturbed Hamiltonian $H(s)=H_{0}+sV$, $0\leq s\leq 1$.

$\displaystyle\Phi_{\beta}^{H(s)}(V)$
$\displaystyle=\int_{-\infty}^{\infty}f_{\beta}(t)e^{iH(s)t}Ve^{-iH(s)t}dt$
$\displaystyle\tilde{f}_{\beta}(\omega)$
$\displaystyle=\begin{cases}\frac{\tanh(\beta\omega/2)}{\beta\omega/2}&\textrm{
for }\omega\neq 0\\\ 1&\textrm{ for }\omega=0\end{cases}$

The inverse Fourier transform is given by [3, Appendix B]

$f_{\beta}(t)=\frac{1}{2\pi}\int_{-\infty}^{\infty}e^{i\omega
t}\tilde{f}_{\beta}(\omega)dt=\frac{2}{\beta\pi}\log\frac{e^{\pi|t|/\beta}+1}{e^{\pi|t|/\beta}-1}.$

Now this defines the QBP operators as follows

$\displaystyle\eta$ $\displaystyle=\eta(1)$ $\displaystyle\eta(s)$
$\displaystyle=\mathcal{T}\exp(-\frac{\beta}{2}\int_{0}^{s}\Phi^{H(s)}_{\beta}(V)ds)$
$\displaystyle=\mathds{1}+\sum_{k\geq
1}\left(-\frac{\beta}{2}\right)^{k}\int_{0}^{s}\int_{0}^{s_{1}}\cdots\int_{0}^{s_{k-1}}\Phi^{H(s_{1})}_{\beta}(V)\ldots\Phi^{H(s_{k})}_{\beta}(V)ds_{1}\ldots
ds_{k}.$

Note that $\|f_{\beta}\|_{L^{1}}=1$ and therefore a straightforward
application of the triangle inequality yields

$\left\|\Phi^{H(s)}_{\beta}(V)\right\|\leq\|V\|.$

In the following, we will specialize our statements to the choice of
Hamiltonian and perturbation relevant for us. To reduce the amount indices we
fix a length scale $L$ and introduce

$\displaystyle H(s)$ $\displaystyle=H_{[1,L]}+sh_{L,L+1}$ $\displaystyle
H_{l}(s)$ $\displaystyle=H_{[L-l,L]}+sh_{L,L+1}.$

By fixing the energy units we keep the assumption $\|h\|\leq 1$. For the sake
of brevity we refrain from giving the explicit dependence of constants on
$\beta$ in the following. We start from the following formulation of Quantum
Belief Propagation for $H_{0}=H(0)$ and $V=h_{L,L+1}$. Note that in the
following, all traces are understood to act on the system $[1,L+1]$, which is
compensated for by a factor $d$.

######  Lemma 4 ([8, Proposition 6]).

We have

$\frac{Z_{L+1}}{Z_{L}}=d\frac{\Tr[\exp(-\beta H(1))]}{\Tr[\exp(-\beta
H(0))]}=d\Tr\left[\eta\frac{\exp(-\beta H(0))}{\Tr[\exp(-\beta
H(0))]}\eta^{\dagger}\right]$

Our aim in the following subsections is to quantify the locality of the QBP
operators in order to approximate this expression by an expectation value of a
local observable in the Gibbs state of a smaller system.

###  A.2 Locality of the QBP Operators

The first step is to approximate the operators $\Phi_{\beta}^{H(s)}$ using
time evolution operators of truncated Hamiltonians. Although the proof is
closely related to the known approach [8], the following Lemma provides a
slightly different formulation. Rather than projecting $\eta$ to a local
operator, we redefine it using an approximation of the time evolution.

######  Lemma 5.

There are some constants $A,B$ depending only on $\beta$ such that for
$0<\varepsilon<1$ and for

$l\geq A\log(\frac{1}{\varepsilon})+B$

we have

$\left\|\Phi_{\beta}^{H(s)}(h_{L,L+1})-\Phi_{\beta}^{H_{l}(s)}(h_{L,L+1})\right\|\leq\varepsilon.$

###### Proof.

The defintion of $\Phi^{H(s)}_{\beta}(h_{L,L+1})$ contains time-evolved
interactions for all times with a time-depdent weight $f_{\beta}(t)$. We split
the definition into short time scales, for which we use a Lieb-Robinson bound
(see Appendix B), and long time-scales, whose contribution can be bounded
using the decay of $f_{\beta}(t)$.

$\displaystyle\left|\Phi_{\beta}^{H(s)}(h_{L,L+1})-\Phi_{\beta}^{H_{l}(s)}(h_{L,L+1})\right|$
$\displaystyle\leq\sup_{|t|\leq
a}\left\|\Gamma^{t}_{H(s)}(h_{L,L+1})-\Gamma^{t}_{H_{l}(s)}(h_{L,L+1})\right\|$
$\displaystyle\phantom{=}+4\left\|h_{L,L+1}\right\|\int_{a}^{\infty}f_{\beta}(t)dt$

The integral can be upper bounded by (see [1, Appendix B])

$\int_{a}^{\infty}f_{\beta}(t)dt\leq\frac{4}{\pi^{2}(e^{\pi
a/\beta}-1)}\leq\frac{\varepsilon}{8}$

if we set

$a=\frac{\beta}{\pi}\log(\frac{8}{\pi^{2}\varepsilon}+1).$

Now, by Theorem 3 setting

$\displaystyle l$
$\displaystyle\geq\frac{E\beta}{\pi}e^{2D}\log(\frac{8}{\pi^{2}\varepsilon}+1)$
$\displaystyle l$
$\displaystyle\geq\frac{C\beta}{D\pi}\log(\frac{8}{\pi^{2}\varepsilon}+1)+\frac{1}{D}\log(\frac{4}{\varepsilon})$

and choosing $D=1$, we have that the norm bound on the first term is also
$\leq\varepsilon/2$. Combining the constants in the conditions on $l$ closes
the proof. ∎

Replacing the $\Phi$ by these local versions in the definition of QBP
operators immediately yields local approximations for those.

######  Lemma 6.

For $l\geq A\log(1/\varepsilon)+B$, where $A,B$ are $\beta$-dependent
constants, we have

$\displaystyle\|\eta-\eta_{l}\|$ $\displaystyle\leq\varepsilon$
$\displaystyle\|\eta_{l}\|,\|\eta\|$ $\displaystyle\leq\exp(\beta/2).$

Here, $\eta_{l}$ is an operator with support on sites $L-l\ldots L+1$:

$\eta_{l}=\mathcal{T}\exp(-\frac{\beta}{2}\int_{0}^{1}\Phi^{H_{l}(s)}_{\beta}(h_{L,L+1})ds)$

###### Proof.

Recall that

$\left\|\Phi^{H(s)}_{\beta}(h_{L,L+1})\right\|,\left\|\Phi^{H_{l}(s)}_{\beta}(h_{L,L+1})\right\|\leq\|h_{L,L+1}\|\leq
1.$

By using the triangle inequality we find

$\displaystyle\|\eta-\eta_{l}\|$ $\displaystyle\leq\sum_{k\geq
1}\left(\frac{\beta}{2}\right)^{k}\int_{0}^{1}\int_{0}^{s_{1}}\cdots\int_{0}^{s_{k-1}}\Big{\|}\Phi^{H(s_{1})}_{\beta}(h_{L,L+1})\ldots\Phi^{H(s_{k})}_{\beta}(h_{L,L+1})$
$\displaystyle\qquad-\Phi^{H_{l}(s_{1})}_{\beta}(h_{L,L+1})\ldots\Phi^{H_{l}(s_{k})}_{\beta}(h_{L,L+1})\Big{\|}ds_{1}\ldots
ds_{k}$ $\displaystyle\leq\sum_{k\geq
1}\left(\frac{\beta}{2}\right)^{k}\int_{0}^{1}\int_{0}^{s_{1}}\cdots\int_{0}^{s_{k-1}}\|h_{L,L+1}\|^{k-1}$
$\displaystyle\qquad\sum_{m=1}^{k}\|\Phi^{H(s_{m})}_{\beta}(h_{L,L+1})-\Phi^{H_{l}(s_{m})}_{\beta}(h_{L,L+1})\|ds_{1}\ldots
ds_{k}$ $\displaystyle\leq\varepsilon_{QBP}\sum_{k\geq
1}\frac{k(\beta/2)^{k}}{k!}$
$\displaystyle=\varepsilon_{QBP}\frac{\beta}{2}e^{\beta/2}$
$\displaystyle\leq\varepsilon,$

where we assumed

$l\geq A\log(\frac{\beta e^{\beta/2}}{2\varepsilon})+B$

for $A,B$ such that
$\varepsilon_{QBP}\leq\varepsilon\frac{2}{\beta}e^{-\beta/2}$ from Lemma 5.
Redefining the constants we conclude the proof. ∎

###  A.3 Local Indisinguishability and Alternative Proof of Lemma 1

In the previous sections we showed that the ratio of partition functions can
be written in terms of QBP operators and that those have approximations by
local operators. Connecting this with the local indistinguishability that we
already saw in Lemma 3, we can reprove Lemma 1, which we restate here.

######  Lemma 7 (Restatement of Lemma 1).

For some $\beta$-dependent constants $A,\xi$, we have

$\left|e^{-\beta f(h)}-\frac{\Tr[\exp(-\beta H_{[1,l+1]})]}{\Tr[\exp(-\beta
H_{[1,l]})]}\right|\leq\varepsilon$

if $l\geq\xi\log(A/\varepsilon)$

###### Proof.

We will start from the ratio of partition functions at finite $L\geq l$. The
following bounds will hold uniformly in $L$ and thereby give the desired
result in the limit. Note that in the following the traces in each fraction
are on equally large subsystems which causes the factors of $d$. By writing
$l=k_{1}+k_{2}$ as a sum of the decay lengths for the QBP and the local
indistinguishability, we estimate

$\displaystyle\left|\frac{Z_{L+1}}{Z_{L}}-\frac{Z_{l+1}}{Z_{l}}\right|$ (6)
$\displaystyle\qquad=\left|d\frac{\Tr[\exp(-\beta H(1))]}{\Tr[\exp(-\beta
H(0))]}-d\Tr[\eta_{k_{1}+k_{2}}\frac{\exp(-\beta
H_{k_{1}+k_{2}})}{\Tr[\exp(-\beta
H_{k_{1}+k_{2}})]}\eta_{k_{1}+k_{2}}^{\dagger}]\right|$ (7)
$\displaystyle\qquad\leq d\left|\Tr[\eta\frac{\exp(-\beta
H_{[1,L]})}{\Tr[\exp(-\beta
H_{[1,L]})]}\eta^{\dagger}]-\Tr[\eta_{k_{1}}\frac{\exp(-\beta
H_{[1,L]})}{\Tr[\exp(-\beta H_{[1,L]})]}\eta_{k_{1}}^{\dagger}]\right|$ (8)
$\displaystyle\qquad\phantom{\leq}+d\left|\Tr[\eta_{k_{1}}\frac{\exp(-\beta
H_{[1,L]})}{\Tr[\exp(-\beta
H_{[1,L]})]}\eta_{k_{1}}^{\dagger}]-\Tr[\eta_{k_{1}}\frac{\exp(-\beta
H_{k_{1}+k_{2}})}{\Tr[\exp(-\beta
H_{k_{1}+k_{2}})]}\eta_{k_{1}}^{\dagger}]\right|$ (9)
$\displaystyle\qquad\phantom{\leq}+d\left|\Tr[\eta_{k_{1}}\frac{\exp(-\beta
H_{k_{1}+k_{2}})}{\Tr[\exp(-\beta
H_{k_{1}+k_{2}})]}\eta_{k_{1}}^{\dagger}]-\Tr[\eta_{k_{1}+k_{2}}\frac{\exp(-\beta
H_{k_{1}+k_{2}})}{\Tr[\exp(-\beta
H_{k_{1}+k_{2}})]}\eta_{k_{1}+k_{2}}^{\dagger}]\right|$ (10)

Despite being defined as the _approximate_ QBP operator for $H_{[1,L]}$, note
that in the first line we used the fact that by its definition $\eta_{l}$ is
the _exact_ QBP for the Hamiltonian $H_{l}$. This feature is due to our
definition of $\eta_{l}$ in terms of truncated time-evolution operators.

For the term (8), note that it is an expectation value of an observable that
obeys the norm bound

$\left\|\Tr_{L+1}\left[\eta^{\dagger}\eta-\eta_{k_{1}}^{\dagger}\eta_{k_{1}}\right]\right\|\leq
d\left\|\eta^{\dagger}\eta-\eta_{k_{1}}^{\dagger}\eta_{k_{1}}\right\|\leq
d(\left\|\eta\right\|+\left\|\eta_{k_{1}}\right\|)\|\eta-\eta_{k_{1}}\|\leq\frac{\varepsilon}{3}$

by choosing $k_{1}$ according to Lemma 6. Analogously by applying the locality
result twice we bound the term (10).

The term (9) can be bounded using local indistinguishability as of Lemma 3
applied to the local observable

$\Tr_{L+1}\left[\eta_{k_{1}}^{\dagger}\eta_{k_{1}}\right],$

which obeys the bound from Lemma 6 and for the decay length $k_{1}+k_{2}-1$.
Hence, adding a condition to the minimum size of $k_{1}$ (or $k_{2}$), we can
bound this term by $\varepsilon/3$ as well.

Combining all the conditions on $k_{1}$ and $k_{2}$, which gives a sufficient
condition for $l$ and taking the limit $L\to\infty$, we obtain the result. ∎

##  Appendix B Lieb-Robinson Bounds for Truncated Time Evolution

We restate Theorem 2.4 from [18]. We modify the formulation in that we assume
2-local interactions instead of only exponentially decaying ones. For finite
interactions the result becomes valid for all temperatures. The result does
_not_ require translation-invariance, but we only state it for our choice of
perturbed Hamiltonian.

We recall the time evolution operator

$\Gamma_{H}^{t}(A)=e^{itH}Ae^{-itH}$

######  Theorem 3 ([18, Theorem 2.4]).

For a an operator $A$ with support on $L,L+1$ and $\beta$-dependent constants
$C,E$,

$\displaystyle\left\|\Gamma_{H(s)}^{t}(A)-\Gamma_{H_{l}(s)}^{t}(A)\right\|$
$\displaystyle\leq\|A\|2e^{8|t|\beta}\sum_{k\geq l+1}\Omega_{k}^{*}(4|t|)$
$\displaystyle\leq 2\|A\|e^{C|t|-Dl}.$

for any choice of $D>0$ as long as $l\geq E|t|e^{2D}$.

###### Proof.

We recall some definitions from [18] specialising to the case of a 2-local 1D
Hamiltonian

$\Omega_{0}=\Omega_{1}=\Omega_{2}=2\beta\qquad\qquad 0=\Omega_{3}=\ldots$
$\displaystyle\Omega_{k}^{*}(x)$
$\displaystyle=\sum_{n=1}^{\infty}\left(\sum_{\begin{subarray}{c}\alpha\in\mathbb{N}^{n}\\\
|\alpha|=k\end{subarray}}\prod_{j=1}^{n}\Omega_{\alpha_{j}}\right)\frac{x^{n}}{n!}$
$\displaystyle\leq\sum_{n\geq\lceil k/2\rceil}\frac{(6x\beta)^{n}}{n!}$

We make repeated use of the tail bound for the exponential series

$\sum_{k\geq l}\frac{x^{k}}{k!}\leq\frac{e^{x}x^{l}}{l!}$

We start from [18, Theorem 2.4]

$\displaystyle\|\Gamma_{H(s)}^{t}(A)-\Gamma_{H_{l}(s)}^{t}(A)\|$
$\displaystyle\leq\|A\||\textup{supp}(A)|e^{8|t|\beta}\sum_{k\geq
l+1}\Omega_{k}^{*}(4|t|)$ $\displaystyle\leq 2\|A\|e^{8|t|\beta}\sum_{k\geq
l+1}e^{24|t|\beta}\frac{(24|t|\beta)^{\lceil k/2\rceil}}{\lceil k/2\rceil!}$
$\displaystyle\leq
2\|A\|e^{32|t|\beta}e^{24|t|\beta}\frac{(24|t|\beta)^{\lceil
l/2\rceil}}{\lceil l/2\rceil!}$ $\displaystyle\leq
2\|A\|e^{56|t|\beta}e^{\lceil l/2\rceil\log(24|t|\beta)-\lceil
l/2\rceil\log(\lceil l/2\rceil)}$

Choosing $l\geq 48|t|\beta e^{2D}=E|t|e^{2D}$ we find

$\|\Gamma_{H(s)}^{t}(A)-\Gamma_{H_{l}(s)}^{t}(A)\|\leq 2\|A\|e^{C|t|-Dl}$

where $C=56\beta$. ∎