Topological susceptibility of SU(3) pure-gauge theory

from out-of-equilibrium simulations

arXiv:2411.00620v1 [hep-lat] 1 Nov 2024

Claudio Bonanno,𝑎 Alessandro Nada𝑏 and Davide Vadacchino𝑐,∗

𝑎 Institutode Física Téorica UAM-CSIC, c/ Nicolás Cabrera 13-15, Universidad Autónoma de Madrid,
Cantoblanco, E-28049 Madrid, Spain
𝑏 Dipartimento di Fisica, Universitá degli Studi di Torino and INFN, Sezione di Torino, Via Pietro Giuria 1,

I-10125 Turin, Italy

𝑐 Centre for Mathematical Sciences, University of Plymouth, Plymouth, PL4 8AA, United Kingdom

E-mail: claudio.bonanno@csic.es, alessandro.nada@unito.it,

davide.vadacchino@plymouth.ac.uk

In JHEP 04 (2024) 126 [arXiv:2402.06561] we recently proposed an out-of-equilibrium setup to

reduce the large auto-correlations of the topological charge in two-dimensional CP 𝑁 −1 models.
Our proposal consists of performing open-boundaries simulations at equilibrium, and gradually
switching on periodic boundary conditions out-of-equilibrium. Our setup allows to exploit the
reduced auto-correlations achieved with open boundaries, avoiding at the same time unphysical
boundary effects thanks to a Jarzynski-inspired reweighting-like procedure. We present prelim-
inary results obtained applying this setup to the 4𝑑 SU(3) pure-gauge theory and we outline a
computational strategy to mitigate topological freezing in this theory.

The 41st International Symposium on Lattice Field Theory (LATTICE2024)

28 July – 3 August 2024
Liverpool, UK

∗ Speaker

© Copyright owned by the author(s) under the terms of the Creative Commons
Attribution-NonCommercial-NoDerivatives 4.0 International License (CC BY-NC-ND 4.0). https://pos.sissa.it/
The SU(3) topological susceptibility from out-of-equilibrium simulations Davide Vadacchino

1. Introduction

Over the last 50 years, lattice field theory has proven to be a successful first-principle approach
to the study of the non-perturbative regime of strongly-interacting gauge theories and in particular
Quantum Chromo–Dynamics (QCD). Crucial to this success have been the striking advances made
in the study and design of algorithms employed to explore their phase space.
Despite these advances, the simulation of theories characterized by topologically non-trivial
phase spaces close to the continuum limit remains difficult. This state of affairs is caused by
topological freezing, a severe divergence of the integrated correlation time 𝜏int of the the topological
charge as a function of the correlation length in the approach towards the continuum limit of said
theories. In contrast to other cases of critical slowing down for non-topological quantities (e.g.,
the Wilson loop), where the divergence is polynomial in the lattice correlation length 𝜉L with a
small exponent (typically 𝜏int ∼ 𝜉L𝑧 with 𝑧 ≃ 2 in the continuum limit 𝜉L → ∞), topological freezing
is characterized by a much more dramatic growth. This has been known to affect both SU(𝑁)
and Sp(2𝑁) Yang–Mills gauge theories in four dimensions, as well as other lower-dimensional
Quantum Field Theories whose vacuum state possess non-trivial topological features, see Refs. [1–
4] and Fig. 1. In practice, this causes the Markov Chain Monte Carlo (MCMC) trajectory of the
system to remain trapped in a sector of the phase space with fixed topological charge, and this
loss of ergodicity can potentially introduce unwanted biases in the estimation of topological and
non-topological physical observables.
Given the theoretical and phenomenological relevance of observables like the topological
susceptibility both for QCD hadron phenomenology [5–12] and for Beyond Standard Model
physics [4, 13–17], several proposal have been put forward over the years in the attempt of ad-
dressing the issue of topological freezing, see Refs. [3, 9, 18–36], and Refs. [37–39] for recent
reviews.
A popular strategy in QCD simulations with dynamical fermions consists in performing simu-
lations with Open Boundary Conditions on the temporal side of the lattice [40, 41]. This approach
reduces the severity of topological critical slowing down, but at the cost of enhancing finite size
effects, which manifest in the form of boundary effects. Another recent strategy that has been
extensively applied both in pure-gauge theories [9, 21, 42–48] and in the presence of dynamical
fermions [36] combines Periodic and Open Boundary Conditions in a parallel tempering frame-
work (Parallel Tempering on Boundary Conditions—PTBC) in order to avoid systematic boundary
effects. In Ref. [49] we have recently proposed a novel approach designed to mitigate topological
freezing which shares a few similarities with PTBC. This new setup is rooted on Jarzynski’s equal-
ity [50] and its application to non-equilibrium evolutions in lattice field theory (see Refs. [51–54])
to combine Open and Periodic Boundary Conditions.
In this contribution, we apply our out-of-equilibrium approach to the four dimensional SU(3)
quenched lattice gauge theory. Our aim is to perform a preliminary assessment of the viability of
our strategy by measuring the topological susceptibility 𝜒 and the integrated autocorrelation time
of the topological charge 𝜏int (𝑄 2 ) at a moderately small value of the lattice spacing. This enables
us to perform a comparison with traditional approaches, that are still viable in this regime, and
with other recent state-of-the-art algorithms, like the Parallel Tempering on Boundary Conditions.
Moreover, this affords us a clear understanding of the computational cost of our approach. This is

2
The SU(3) topological susceptibility from out-of-equilibrium simulations Davide Vadacchino

103
Heat-Bath + Over-Relaxation

102

τint(Q2)
101

100
4 6 8 10
r0/a

Figure 1: Growth of the integrated auto-correlation time 𝜏int of 𝑄 2 as a function of the inverse lattice spacing.
This quantity expresses the number of updating steps required to generate two decorrelated measures of 𝑄 2 .
In this case 𝜏int (𝑄 2 ) is given in units of standard updating steps, defined as 4 lattice sweeps of over-relaxation
and 1 lattice sweep of heat-bath. Data for 𝜏int (𝑄 2 ) are taken from dataset generated in Ref. [12], while the
scale setting in units of the Sommer parameter 𝑟 0 ≃ 0.5 fm is taken from the interpolation of the results of
Ref. [58].

the first step towards undertaking the same analysis closer to the continuum limit, where traditional
approaches are doomed to fail in sampling non-trivial topology. This also provides a first solid
test-bed to implement improvements based on recent machine-learning techniques like Stochastic
Normalizing Flows, see for example Refs. [55–57], that fit into this approach quite naturally.
This proceeding is structured as follows. In Sec. 2 we briefly describe our approach and provide
details about our setup. In Sec. 3 we report on our numerical results and discuss them. We then
conclude and provide a roadmap for the future steps of this analysis in Sec. 4.

2. Setup

Our approach it based on the use of non-equilibrium evolutions and Jarzynski’s inequality to
resample ensembles of configurations produced with Open Boundary Conditions (on a defect). It
is well known that the use of Open Boundary Conditions (OBC) strongly mitigates topological
freezing. The crucial difference with Periodic Boundary Conditions (PBC) is the absence, in
the former case, of those potential barriers that separate different topological sectors, causing the
freezing of topological charge in the first place. While effective in overcoming topological freezing,
OBC come with the important drawback of introducing non-physical effects that need to be removed
by computing correlation functions sufficiently far from the boundaries. Our strategy is designed to
exploit the advantages of OBC without suffering from its pitfalls (similarly to the PTBC algorithm).
First, the so-called defect is defined: a three-dimensional subset of the lattice on which the
coupling can be tuned with the purpose of interpolating between OBC, corresponding to coupling
0, to PBC, corresponding to coupling 𝛽, like on the rest of the lattice. This is achieved by defining

3
The SU(3) topological susceptibility from out-of-equilibrium simulations Davide Vadacchino

the system on a hypercubic lattice of length 𝐿 lattice spacings, with the action
𝛽 ∑︁ (𝑛) h i
(𝑛)
𝑆 𝑛 [𝑈] = − 𝐾 𝜇𝜈 (𝑥)ℜTr 𝑃 𝜇𝜈 (𝑥) (1)
𝑁 𝑐 𝑥, 𝜇≠𝜈

(𝑛) † †
where 𝑃 𝜇𝜈 (𝑥) = 𝑈 𝜇(𝑛) (𝑥)𝑈𝜈(𝑛) (𝑥 +𝑎 𝜇)𝑈
ˆ 𝜇(𝑛) (𝑥 +𝑎 𝜈)𝑈
ˆ 𝜈(𝑛) (𝑥) is the elementary plaquette operator at
(𝑛)
site 𝑥 on the plane (𝜇, 𝜈), and 𝐾 𝜇𝜈 (𝑥) is a numerical factor used to modify the boundary conditions.
It is defined as:
(𝑛)
𝐾 𝜇𝜈 (𝑥) = 𝐾 𝜇(𝑛) (𝑥)𝐾 𝜈(𝑛) (𝑥 + 𝑎 𝜇)𝐾
ˆ 𝜇(𝑛) (𝑥 + 𝑎 𝜈)𝐾
ˆ 𝜈(𝑛) (𝑥), (2)


 𝑐(𝑛),
 𝜇 = 1, 𝑥1 = 𝐿 − 𝑎, 0 ≤ 𝑥0 , 𝑥2 , 𝑥3 < 𝐿 d ,
𝐾 𝜇(𝑛) (𝑥)

= (3)
 1,

elsewhere,

where the size of the three-dimensional defect is given by 𝐿 d × 𝐿 d × 𝐿 d , while 𝑐(𝑛) is an a priori
arbitrary function. The rôle of the latter is to interpolate between OBC when 𝑐 = 0 to PBC when
𝑐 = 1.
Second, an ensemble is obtained from equilibrium MCMC simulation with OBC and the
standard combination of local HB+OR updates. We call 𝑛between the number of full lattice sweeps
between successive configurations along the Markov Chain with OBC. Each configuration in this
ensemble is used as starting point for a non-equilibrium evolution towards a configuration with
PBC. In particular, iterations of the same local updating algorithm are run while gradually changing
the value of 𝑐 along a so-called protocol,
𝑛
𝑐(𝑛) = 1 − (4)
𝑛step − 1
where 𝑛step is the number of steps separating the two different boundary conditions. We remark that
in this framework the system does not have to reach equilibrium between successive steps.
Third, any desired observable is calculated on the configurations of the newly obtained ensemble
with PBC by performing an appropriate reweighting with the statistical weight 𝑒 −𝑊 ,
⟨O𝑒 −𝑊 ⟩ f
⟨O⟩ NE = , (5)
⟨𝑒 −𝑊 ⟩ f
where 𝑊 is the work spent along a non-equilibrium evolution, that can be calculated explicitly as
follows,
∑︁−1
𝑛step
𝑊= {𝑆 𝑛+1 [𝑈𝑛 ] − 𝑆 𝑛 [𝑈𝑛 ]} , (6)
𝑛=0

with 𝑈𝑛 the gauge configuration at the 𝑛th step of the out-of-equilibrium transformation.
This strategy behind Eq. (5) can be summarized as follows: we aim at sampling a target
distribution 𝑝 with PBC by reweighting a sample obtained from the prior distribution 𝑞 0 with OBC
and driven towards PBC through a sequence of non-equilibrium forward transitions with probability
𝑛step
Ö
Pf [𝑈0 , . . . , 𝑈] = P𝑛 (𝑈𝑛−1 → 𝑈𝑛 ), (7)
𝑛=1

4
The SU(3) topological susceptibility from out-of-equilibrium simulations Davide Vadacchino

eq MC
nstep non-eq MC

nstep

nstep

nstep
nbetween nbetween nbetween nbetween
Figure 2: Sketch of our out-of-equilibrium Monte Carlo setup.

where each transition probability P𝑛 is defined by the intermediate action 𝑆 𝑛 . The expectation
value of any observable O with respect to the target distribution
1 −𝑆 [𝑈 ]
𝑝[𝑈] = 𝑒 , 𝑆 ≡ 𝑆 𝑛step , (8)
𝑍
can be computed using Eq. (5), where the average over the evolutions is formally defined as
∫
1 −𝑆0 [𝑈0 ]
⟨. . . ⟩f = [d𝑈0 . . . d𝑈]𝑞 0 [𝑈0 ] Pf [𝑈0 , . . . , 𝑈] . . . , 𝑞 0 [𝑈0 ] = 𝑒 . (9)
𝑍0
We refer to Ref. [49] for a more in-depth discussion of these non-equilibrium evolutions.
A sketch of our out-of-equilibrium Monte Carlo (MC) setup is shown in Fig. 2. The horizontal
axis represents equilibrium evolutions with OBCs, used to sample the prior distribution 𝑞 0 (OBC).
The vertical axis represents out-of-equilibrium evolutions with protocol 𝑐(𝑛), used to gradually
reach the target distribution 𝑝 (PBC). The out-of-equilibrium trajectories are 𝑛step update steps long,
and are separated by 𝑛between steps, the latter performed at equilibrium.
The ratio between the partition functions of the prior 𝑞 0 and target 𝑝, 𝑍0 and 𝑍 respectively,
defines the free energy difference Δ𝐹
𝑍
𝑒 −Δ𝐹 = , (10)
𝑍0
which can also be computed from the work along the non-equilibrium evolutions using the celebrated
Jarzynski’s equality [59],

⟨exp(−𝑊)⟩ f = exp (−Δ𝐹) . (11)

Of crucial importance is clearly the quality of the sampling of 𝑝. This can be evaluated using
two figures of merit that estimate how far from equilibrium the evolution is. The first is the reverse

5
The SU(3) topological susceptibility from out-of-equilibrium simulations Davide Vadacchino

Kullback–Leibler divergence,

eKL (𝑞 0 Pf ∥ 𝑝Pr ) = ⟨𝑊⟩ f − Δ𝐹 ≥ 0,

𝐷 (12)

which measures the reversibility of the evolution defined by a given protocol 𝑐(𝑛), and the second
is the Effective Sample Size (ESS), estimated with
−𝑊 2
ˆ ≡ ⟨𝑒 ⟩f =
ESS
1
, (13)
⟨𝑒 −2𝑊 ⟩f ⟨𝑒 −2(𝑊 −Δ𝐹 ) ⟩
f

The ESS has an interesting intuitive interpretation. It is easy to show that the variance of ⟨O⟩
obtained from 𝑛 non-equilibrium measurements is related to the corresponding variance sampled
from the target distribution 𝑝 as follows,

Var(O)NE Var(O) 𝑝
= . (14)
𝑛 𝑛 ESS
Hence, the former is larger by a factor of 1/ESS with respect to the latter, obtained by sampling 𝑝
at equilibrium.
In the above formula, the autocorrelations between evolutions have been ignored for simplicity.
However, they must be accounted for whenever an estimate of the total cost of using the above
strategy is needed. A possible way of estimating the latter was provided in Ref. [49]. We report it
here,
2𝜏int
Var(O)NE × (𝑛step + 𝑛between ) ≃ Var(O) 𝑝 × (𝑛step + 𝑛between ). (15)
ˆ
ESS
Since this quantity can be computed for any other algorithm, it can be used as a quantitative
statistics-independent metric for the efficiency of the algorithm in computing a given observable
with a target precision.

3. Numerical results

In this section we report on the preliminary results obtained for the two figures of merit ESS ˆ
and 𝐷eKL and for the topological susceptibility. All our results refer to simulations conducted for
𝛽 = 6.40 on a symmetric (𝐿/𝑎) 4 = 304 lattice. Assuming 𝑟 0 ≃ 0.5 fm for the Sommer scale,
this point corresponds to a lattice spacing 𝑎/𝑟 0 ≃ 0.1027(5), i.e., 𝑎 ≃ 0.05 fm and 𝐿 ≃ 1.5 fm,
which is a large enough box to be insensitive to finite-size effects in the topological susceptibility
with periodic boundary conditions. The integrated auto-correlation time with PBC of 𝑄 2 was
estimated to be 𝜏int (𝑄 2 ) ≃ 520(20) in units of standard Monte Carlo updating steps using standard
algorithms [7], see also the caption of Fig. 1.
The values of the figures of ESSˆ and 𝐷 eKL can be calculated as follows. The former is obtained
using Eq. (13) on the values of the work defined by Eq. (6). The values of 𝐷 eKL are also obtained
from the work by direct computation using Eqs. (11) and (12). These quantities were computed for
several values of 𝐿 d and 𝑛step and are displayed in Fig. 3 as a function of 𝑛step /(𝐿 d /𝑎) 3 . In Ref. [49],
the latter was identified as the proper scaling variable for these figures of merit in the case of the 2𝑑
CP 𝑁 −1 . Fig. 3 shows that this is also the case for the 4𝑑 SU(3) pure lattice gauge theory. Indeed,

6
The SU(3) topological susceptibility from out-of-equilibrium simulations Davide Vadacchino

β = 6.4 β = 6.4
1.2 1.0
Ld /a = 3 Ld /a = 3, nbetween = 10 Ld /a = 4, nbetween = 50
1.0 Ld /a = 4
0.8
Ld /a = 3, nbetween = 50 Ld /a = 5, nbetween = 10
Ld /a = 5 Ld /a = 3, nbetween = 100 Ld /a = 5, nbetween = 50
0.8 Ld /a = 6 Ld /a = 4, nbetween = 10 Ld /a = 6, nbetween = 10
0.6
D̃KL

ESS
0.6

ˆ
0.4
0.4

0.2 0.2

0.0 0.0
2 3 4 5 6 7 2 3 4 5 6 7
nstep /(Ld /a)3 nstep /(Ld /a)3

Figure 3: Figures of merit to quantify how far from equilibrium the evolution is as a function of the length
of the out-of-equilibrium trajectory 𝑛step and of the defect size in lattice units (𝐿 d /𝑎) 3 .

to a very good degree of approximation, the values of 𝐷 eKL and ESSˆ are seen to collapse on a single
3
curve when represented as a function of 𝑛step /(𝐿 d /𝑎) . Verifying this scaling is key to understanding
how to tune the simulation parameters as the lattice spacing is reduced, and to estimate the cost of
a simulation as a function of 𝛽.
We now move on to our evaluation of the topological susceptibility. The latter was defined as
follows,

1
𝑎4 𝜒 ≡ ⟨𝑄 2 ⟩ , (16)
𝑉

with 𝑉 = (𝐿/𝑎) 4 and where 𝑄 was obtained using a standard clover discretization and computed
after 30 cooling steps. We performed several runs with different values of 𝑛step and 𝐿 d . In the
left-hand side panel of Fig. 4, we display the non-equilibrium values thus obtained as a function
of 𝐷 eKL . The value of the topological susceptibility obtained with traditional algorithms, see
Ref. [12], is represented as a black horizontal line. In all cases, the non-equilibrium and traditional
determinations are in agreement within 1.5 standard deviations at most. This shows that, for the
values of 𝐷 eKL attained in our settings, the susceptibility obtained with non-equilibrium methods is
not significantly affected by the details of the simulation and it reproduces the equilibrium estimates.
Finally, in order to assess the efficiency of the non-equilibrium method we have estimated
the cost of computing the topological susceptibility with this approach, using the left-hand side of
Eq. (15). The latter quantity was computed for each of the available values of 𝑛step and 𝐿 d and is
displayed in the right-hand side panel of Fig. 4 as a function of the integrated autocorrelation time
𝜏int (𝑄 2 ). Smaller autocorrelation times are obtained for larger values of either 𝐿 d or 𝑛between , or both,
as expected. In the regime of small autocorrelations, 𝜏int (𝑄 2 ) ∼ 0.5–1.5, the values of the cost gather
around 1, which appears to be the optimal regime for our setup. Given that at the currently explored
value of 𝛽 this was achieved for values of 𝑛between and 𝑛step of the order of a few hundreds, this means
that the cost of 1 decorrelated out-of-equilibrium trajectory in terms of number of updates only
improves by a small factor on the results obtained with the standard algorithm, 𝜏int (𝑄 2 ) = 520(20).
However, the gain obtained with our method is expected to increase at larger values of 𝛽, owing to
the improved scaling of 𝜏int (𝑄 2 ) obtained with OBC.

7
The SU(3) topological susceptibility from out-of-equilibrium simulations Davide Vadacchino

×10−5 β = 6.4 ×10−7

1.0 4
Ld /a = 3, nbetween = 50 Ld /a = 5, nbetween = 10

Var(Q2 /V ) ×(nstep + nbetween )

Standard MCMC Ld /a = 5, nbetween = 10
Ld /a = 3, nbetween = 50 Ld /a = 5, nbetween = 50 Ld /a = 3, nbetween = 100 Ld /a = 5, nbetween = 50
0.9
Ld /a = 3, nbetween = 100 Ld /a = 6, nbetween = 10 Ld /a = 4, nbetween = 50 Ld /a = 6, nbetween = 10
3
0.8 Ld /a = 4, nbetween = 50
hQ2 i/V

0.7 2

0.6
1
0.5

0.5 0.6 0.7 0.8 0.9 1.0 1.1 0.5 1.0 1.5 2.0 2.5 3.0
D̃KL τint (Q2 )

Figure 4: Left: results for the topological susceptibility. Right: results for the variable defined in Eq. (15),
which expresses the efficiency of the out-of-equilibrium algorithm.

4. Conclusions and future outlooks

In this proceeding we presented a preliminary study of the topological susceptibility of the pure
SU(3) lattice gauge theory using out-of-equilibrium simulations. We implemented the same setup
already proposed in and tested in Ref. [49], and our results are encouraging. The figures of merit
ˆ and 𝐷
ESS eKL are found to only depend on the scaling variable 𝑛step /(𝐿 d /𝑎) 3 rather than on 𝑛step and
𝐿 d separately. Moreover, our estimate of the topological susceptibility is statistically compatible
with the estimates present in the literature, with an integrated auto-correlation time for 𝑄 2 that is
slightly smaller than the one attained with traditional algorithms.
In the near future we plan to extend the present preliminary study in two different directions.
Firstly, by probing larger values of 𝛽 and by studying systematically how the cost of the out-of-
equilibrium behaves as the lattice spacing is reduced. Secondly and most importantly, by combining
the setup discussed and described above with the discrete neural network layers composing Normal-
izing Flows. Stochastic Normalizing Flows fit naturally in the out-of-equilibrium setup, and have the
potential of improving the efficiency of our approach even further, see for example the improvement
obtained for SU(3) gauge theory in Ref. [57]. This would be achieved by reducing the amount of
updating steps necessary to obtain a fixed ESS, at the price of performing a once-and-for-all training
procedure with moderate costs.

Acknowledgments

We thank L. Giusti, A. Patella and F. Sanfilippo for insightful comments and discussions.
Numerical simulations have been performed on the DiRAC Data Intensive service at Cambridge.
The work of C. Bonanno is supported by the Spanish Research Agency (Agencia Estatal de In-
vestigación) through the grant IFT Centro de Excelencia Severo Ochoa CEX2020-001007-S and,
partially, by grant PID2021-127526NB-I00, both funded by MCIN/AEI/10.13039/ 501100011033.
A. Nada acknowledges support by the Simons Foundation grant 994300 (Simons Collaboration
on Confinement and QCD Strings), from the European Union - Next Generation EU, Mission 4
Component 1, CUP D53D23002970006, under the Italian PRIN “Progetti di Ricerca di Rilevante
Interesse Nazionale – Bando 2022” prot. 2022TJFCYB and from the SFT Scientific Initiative of
INFN.

8
The SU(3) topological susceptibility from out-of-equilibrium simulations Davide Vadacchino

References

[1] B. Alles, G. Boyd, M. D’Elia, A. Di Giacomo and E. Vicari, Hybrid Monte Carlo and
topological modes of full QCD, Phys. Lett. B 389 (1996) 107 [hep-lat/9607049].

[2] L. Del Debbio, G. M. Manca and E. Vicari, Critical slowing down of topological modes,
Phys. Lett. B 594 (2004) 315 [hep-lat/0403001].

[3] C. Bonati and M. D’Elia, Topological critical slowing down: variations on a toy model,
Phys. Rev. E 98 (2018) 013308 [1709.10034].

[4] Bennett, D. K. Hong, J.-W. Lee, C. J. D. Lin, B. Lucini, M. Piai et al., Sp(2N) Yang-Mills
theories on the lattice: Scale setting and topology, Phys. Rev. D 106 (2022) 094503
[2205.09364].

[5] M. Cè, C. Consonni, G. P. Engel and L. Giusti, Non-Gaussianities in the topological charge
distribution of the 𝑆𝑈 (3) Yang–Mills theory, Phys. Rev. D 92 (2015) 074502 [1506.06052].

[6] M. Cè, M. Garcia Vera, L. Giusti and S. Schaefer, The topological susceptibility in the
large-𝑁 limit of SU(𝑁) Yang-Mills theory, Phys. Lett. B 762 (2016) 232 [1607.05939].

[7] C. Bonati, M. D’Elia and A. Scapellato, 𝜃 dependence in 𝑆𝑈 (3) Yang-Mills theory from
analytic continuation, Phys. Rev. D 93 (2016) 025028 [1512.01544].

[8] C. Bonati, M. D’Elia, P. Rossi and E. Vicari, 𝜃 dependence of 4D 𝑆𝑈 (𝑁) gauge theories in
the large-𝑁 limit, Phys. Rev. D 94 (2016) 085017 [1607.06360].

[9] C. Bonanno, C. Bonati and M. D’Elia, Large-𝑁 𝑆𝑈 (𝑁) Yang-Mills theories with milder
topological freezing, JHEP 03 (2021) 111 [2012.14000].

[10] A. Athenodorou and M. Teper, The glueball spectrum of SU(3) gauge theory in 3 + 1
dimensions, JHEP 11 (2020) 172 [2007.06422].

[11] A. Athenodorou and M. Teper, SU(N) gauge theories in 3+1 dimensions: glueball spectrum,
string tensions and topology, JHEP 12 (2021) 082 [2106.00364].

[12] C. Bonanno, The topological susceptibility slope 𝜒 ′ of the pure-gauge SU(3) Yang–Mills
theory, JHEP 01 (2024) 116 [2311.06646].

[13] C. Bonati, M. D’Elia, M. Mariti, G. Martinelli, M. Mesiti, F. Negro et al., Axion

phenomenology and 𝜃-dependence from 𝑁 𝑓 = 2 + 1 lattice QCD, JHEP 03 (2016) 155
[1512.06746].

[14] P. Petreczky, H.-P. Schadler and S. Sharma, The topological susceptibility in finite
temperature QCD and axion cosmology, Phys. Lett. B 762 (2016) 498 [1606.03145].

[15] S. Borsanyi et al., Calculation of the axion mass based on high-temperature lattice quantum
chromodynamics, Nature 539 (2016) 69 [1606.07494].

9
The SU(3) topological susceptibility from out-of-equilibrium simulations Davide Vadacchino

[16] M. P. Lombardo and A. Trunin, Topology and axions in QCD, Int. J. Mod. Phys. A 35 (2020)
2030010 [2005.06547].

[17] A. Athenodorou, C. Bonanno, C. Bonati, G. Clemente, F. D’Angelo, M. D’Elia et al.,

Topological susceptibility of N 𝑓 = 2 + 1 QCD from staggered fermions spectral projectors at
high temperatures, JHEP 10 (2022) 197 [2208.08921].

[18] W. Bietenholz, P. de Forcrand and U. Gerber, Topological Susceptibility from Slabs, JHEP 12
(2015) 070 [1509.06433].

[19] A. Laio, G. Martinelli and F. Sanfilippo, Metadynamics surfing on topology barriers: the
𝐶𝑃 𝑁 −1 case, JHEP 07 (2016) 089 [1508.07270].

[20] M. Lüscher, Stochastic locality and master-field simulations of very large lattices, EPJ Web
Conf. 175 (2018) 01002 [1707.09758].

[21] M. Hasenbusch, Fighting topological freezing in the two-dimensional 𝐶𝑃 𝑁 −1 model, Phys.

Rev. D 96 (2017) 054504 [1706.04443].

[22] L. Giusti and M. Lüscher, Topological susceptibility at 𝑇 > 𝑇c from master-field simulations
of the SU(3) gauge theory, Eur. Phys. J. C 79 (2019) 207 [1812.02062].

[23] A. Florio, O. Kaczmarek and L. Mazur, Open-Boundary Conditions in the Deconfined

Phase, Eur. Phys. J. C 79 (2019) 1039 [1903.02894].

[24] L. Funcke, K. Jansen and S. Kühn, Topological vacuum structure of the Schwinger model
with matrix product states, Phys. Rev. D 101 (2020) 054507 [1908.00551].

[25] G. Kanwar, M. S. Albergo, D. Boyda, K. Cranmer, D. C. Hackett, S. Racanière et al.,

Equivariant flow-based sampling for lattice gauge theory, Phys. Rev. Lett. 125 (2020)
121601 [2003.06413].

[26] K. A. Nicoli, C. J. Anders, L. Funcke, T. Hartung, K. Jansen, P. Kessel et al., Estimation of

Thermodynamic Observables in Lattice Field Theories with Deep Generative Models, Phys.
Rev. Lett. 126 (2021) 032001 [2007.07115].

[27] D. Albandea, P. Hernández, A. Ramos and F. Romero-López, Topological sampling through

windings, Eur. Phys. J. C 81 (2021) 873 [2106.14234].

[28] G. Cossu, D. Lancastera, B. Lucini, R. Pellegrini and A. Rago, Ergodic sampling of the
topological charge using the density of states, Eur. Phys. J. C 81 (2021) 375 [2102.03630].

[29] S. Borsanyi and D. Sexty, Topological susceptibility of pure gauge theory using Density of
States, Phys. Lett. B 815 (2021) 136148 [2101.03383].

[30] G. Papamakarios, E. T. Nalisnick, D. J. Rezende, S. Mohamed and B. Lakshminarayanan,

Normalizing flows for probabilistic modeling and inference, J. Mach. Learn. Res. 22 (2021) 1.

10
The SU(3) topological susceptibility from out-of-equilibrium simulations Davide Vadacchino

[31] P. Fritzsch, J. Bulava, M. Cè, A. Francis, M. Lüscher and A. Rago, Master-field simulations
of QCD, PoS LATTICE2021 (2022) 465 [2111.11544].

[32] R. Abbott et al., Normalizing flows for lattice gauge theory in arbitrary space-time
dimension, 2305.02402.

[33] T. Eichhorn, G. Fuwa, C. Hoelbling and L. Varnhorst, Parallel tempered metadynamics:

Overcoming potential barriers without surfing or tunneling, Phys. Rev. D 109 (2024) 114504
[2307.04742].

[34] D. Howarth and A. J. Peterson, Topological charge unfreezing with AMReX, 2312.11599.

[35] D. Albandea, G. Catumba and A. Ramos, The Strong CP Problem in the Quantum Rotor,
2402.17518.

[36] C. Bonanno, G. Clemente, M. D’Elia, L. Maio and L. Parente, Full QCD with milder
topological freezing, JHEP 08 (2024) 236 [2404.14151].

[37] J. Finkenrath, Review on Algorithms for dynamical fermions, PoS LATTICE2022 (2023)
227 [2402.11704].

[38] P. A. Boyle, Advances in algorithms for solvers and gauge generation, 2401.16620.

[39] J. Finkenrath, Future trends in lattice QCD simulations, PoS EuroPLEx2023 (2024) 009.

[40] M. Lüscher and S. Schaefer, Lattice QCD without topology barriers, JHEP 07 (2011) 036
[1105.4749].

[41] M. Lüscher and S. Schaefer, Lattice QCD with open boundary conditions and twisted-mass
reweighting, Comput. Phys. Commun. 184 (2013) 519 [1206.2809].

[42] M. Berni, C. Bonanno and M. D’Elia, Large-𝑁 expansion and 𝜃-dependence of 2𝑑 𝐶𝑃 𝑁 −1

models beyond the leading order, Phys. Rev. D 100 (2019) 114509 [1911.03384].

[43] C. Bonanno, Lattice determination of the topological susceptibility slope 𝜒 ′ of 2𝑑 CP 𝑁 −1

models at large 𝑁, Phys. Rev. D 107 (2023) 014514 [2212.02330].

[44] C. Bonanno, M. D’Elia, B. Lucini and D. Vadacchino, Towards glueball masses of large-N
SU(N) pure-gauge theories without topological freezing, Phys. Lett. B 833 (2022) 137281
[2205.06190].

[45] J. L. Dasilva Golán, C. Bonanno, M. D’Elia, M. García Pérez and A. Giorgieri, The twisted
gradient flow strong coupling with parallel tempering on boundary conditions, PoS
LATTICE2023 (2024) 354 [2312.09212].

[46] C. Bonanno, M. D’Elia and L. Verzichelli, The 𝜃-dependence of the 𝑆𝑈 (𝑁) critical
temperature at large 𝑁, JHEP 02 (2024) 156 [2312.12202].

[47] C. Bonanno, C. Bonati, M. Papace and D. Vadacchino, The 𝜃-dependence of the Yang-Mills
spectrum from analytic continuation, JHEP 05 (2024) 163 [2402.03096].

11
The SU(3) topological susceptibility from out-of-equilibrium simulations Davide Vadacchino

[48] C. Bonanno, J. L. Dasilva Golán, M. D’Elia, M. García Pérez and A. Giorgieri, The SU(3)
twisted gradient flow strong coupling without topological freezing, Eur. Phys. J. C 84 (2024)
916 [2403.13607].

[49] C. Bonanno, A. Nada and D. Vadacchino, Mitigating topological freezing using

out-of-equilibrium simulations, JHEP 04 (2024) 126 [2402.06561].

[50] C. Jarzynski, Nonequilibrium Equality for Free Energy Differences, Phys. Rev. Lett. 78
(1997) 2690 [cond-mat/9610209].

[51] M. Caselle, G. Costagliola, A. Nada, M. Panero and A. Toniato, Jarzynski’s theorem for
lattice gauge theory, Phys. Rev. D 94 (2016) 034503 [1604.05544].

[52] M. Caselle, A. Nada and M. Panero, QCD thermodynamics from lattice calculations with
nonequilibrium methods: The SU(3) equation of state, Phys. Rev. D 98 (2018) 054513
[1801.03110].

[53] O. Francesconi, M. Panero and D. Preti, Strong coupling from non-equilibrium Monte Carlo
simulations, JHEP 07 (2020) 233 [2003.13734].

[54] A. Bulgarelli and M. Panero, Entanglement entropy from non-equilibrium Monte Carlo
simulations, JHEP 06 (2023) 030 [2304.03311].

[55] H. Wu, J. Köhler and F. Noé, Stochastic normalizing flows, Advances in Neural Information
Processing Systems 33 (2020) 5933.

[56] M. Caselle, E. Cellini, A. Nada and M. Panero, Stochastic normalizing flows as

non-equilibrium transformations, JHEP 07 (2022) 015 [2201.08862].

[57] A. Bulgarelli, E. Cellini and A. Nada, Sampling SU(3) pure gauge theory with Stochastic
Normalizing Flows, in 41st International Symposium on Lattice Field Theory, 9, 2024,
2409.18861.

[58] S. Necco and R. Sommer, The 𝑁 𝑓 = 0 heavy quark potential from short to intermediate
distances, Nucl. Phys. B 622 (2002) 328 [hep-lat/0108008].

[59] C. Jarzynski, Equilibrium free-energy differences from nonequilibrium measurements: A

master-equation approach, Phys. Rev. E56 (1997) 5018 [cond-mat/9707325].

12