Variational Autoregressive Networks with probability priors.
arXiv:2605.16020v1 [cs.LG] 15 May 2026
Piotr Białasa,∗, Piotr Korcylb , Tomasz Stebelb , Dawid Zapolskia,c a Institute of Applied Computer Science, Jagiellonian University, ul. Łojasiewicza 11, 30-348, Kraków, Poland b Institute of Theoretical Physics, Jagiellonian University, ul. Łojasiewicza 11, 30-348, Kraków, Poland c Doctoral School of Exact and Natural Sciences, Jagiellonian University, ul. Łojasiewicza 11, 30-348, Kraków, Poland
Abstract Monte Carlo methods are essential across diverse scientific fields, yet their efficiency is frequently hampered by critical slowing down—a sharp increase in autocorrelation times near phase transitions. Although deep learning approaches, such as neural-network-based samplers, have been proposed to alleviate this issue, they face another serious problem: the difficulty of training the models. This difficulty partially stems from the overly general nature of original machinelearning architectures, which often ignore underlying physical symmetries and force networks to relearn them from scratch. In this paper, we demonstrate that incorporating physical priors into the model significantly enhances performance. Building upon existing strategies that integrate spin-spin interactions, we propose a framework that utilizes a prior probability distribution as a starting point for training. Our results for the Ising model, as well as for the EdwardsAnderson spin glass model, suggest that moving away from "blank slate" models in favor of physics-informed priors reduces the training burden and facilitates the simulation of larger system sizes in discrete spin models. Keywords: Monte-Carlo, neural networks
1. Introduction Monte-Carlo methods play a very important role in many areas, ranging from social sciences, through quantitative finance, to physics. The generation of samples is often done by some variant of the Metropolis-Hastings algorithm, ∗ Corresponding author.
Email addresses: [email protected] (Piotr Białas), [email protected] (Piotr Korcyl), [email protected] (Tomasz Stebel), [email protected] (Dawid Zapolski)
which resorts to the construction of a Markov chain of consecutive configurations [1]. However, by the nature of this construction, which builds the new configuration based on the previous one, the configurations are correlated. Thus, the efficiency of the simulation depends not only on how fast the configurations can be generated but also on the number of configurations that we have to generate before they are independent. This is quantified by the autocorrelation time. Unfortunately, in many cases, this autocorrelation time grows in the interesting region, i.e. near the phase transition. This undesired phenomenon is called critical slowing down and is the bane of lattice QCD. With the advent of deep learning, there have been renewed efforts to apply machine-learning techniques to alleviate this problem. In a seminal work [2], the authors proposed using neural networks to generate samples from a target distribution for discrete spin systems1 . Although this approach significantly reduces autocorrelation times [4], the system sizes that could be simulated remain limited compared to those accessible with specialized Monte-Carlo methods. Although the generation of samples is relatively fast and the samples are not correlated, the training of the neural networks has emerged as the main computational bottleneck. We believe that the difficulty of training comes, at least in part, from the very general nature of the original proposal in [2]. If we restrict our consideration to two-state spin systems, then the only input to the training algorithm is the energy of a configuration. Moreover, the factorization of the Boltzmann probability distribution into a product of conditional probabilities required by this approach destroys the existing symmetries. Those symmetries must be relearned by the networks, which increases the difficulty of the training. Many of the improvements that allowed simulation of larger system sizes were due to the incorporation of at least some part of physics [5, 6, 7, 8, 9]. In this contribution, instead of starting from a blank slate, we propose to add a prior probability distribution as a starting point for the training. This can be thought of as an extension of the approach presented in [5, 6] and [9], where some of the spin-spin interactions were incorporated into the model, leading to faster training times. This work is organized as follows. In the next chapter, we introduce autoregressive neural generators. Next, we review some of the methods that aim to improve a simple approach presented in [2]. In the following chapter, we show how to augment this approach when we have some approximation of the desired probability distribution. The rest of the work presents the results for the ferromagnetic Ising model and the Edwards-Anderson spin glass model [10]. 1 For continuous models, like e.g. ϕ4 theory, an approach based on normalizing flows was proposed in [3]; however, in this contribution, we restrict ourselves solely to discrete spin systems.
2
2. Autoregressive neural generators Given some set of N spins, s = (s0 , s1 , . . . , sN −1 )
(1)
we will be concerned with the Boltzmann distribution2 p(s) = −1
where β = (kB T ) ration s and
1 −βE(s) e Z(β)
(2)
is the inverse temperature, E(s) is the energy of configuX Z(β) = e−βE(s) (3) s
is the so called partition function. In practice, we seldom know the partition function Z(β), so we only have access to unnormalized probability P (s) = e−βE(s) .
(4)
The idea of neural generators is to train a model to approximate the probability distribution p(s) by some other probability distribution q(s). For discrete systems, we represent the probability q(s) as a product of conditional probabilities N −1 Y q(s0 , s1 , . . . , sN −1 ) = q(s0 ) q(si |s<i ) (5) i=1
where N is the number of spins and s<i ≡ (s0 , s1 , . . . , si−1 )
(6)
sop i ≡ {sj : j op i = true}.
(7)
or in general where op is any operator used to compare indices i and j. In this way, we can use some ML model for q(si |s<i ). Any chosen model must have the autoregressive property: the output i can only depend on the inputs j < i. In the simplest example of a system of spins that can take only two values of ±1 we need a model with N inputs and N outputs. Each input corresponds to one spin and each output to the conditional probability q(si = 1|s<i ). In the following, we will use the shorthand qi (s<i ) ≡ q(si = 1|s<i ). (8) and similarly for p(si |s<i ). The authors of [2] call such models based on neural networks Variational Autoregressive Networks (VAN). Given that, we can sample the configurations s from q(s) using ancestral sampling. Once we know the 2 Please note that most distributions can be put in this form.
3
values of the spins s1 , s2 , . . . , si−1 , we sample the spin si from the distribution q(si |s<i ) and so on. The model is trained by minimizing the quantity that, with some abuse of notation (as it differs from free energy by a factor of β −1 ), we will call Fq X X q(s) q(s) q(s) log − log Z. (9) Fq = q(s) log −βE(s) = p(s) e s s Defining F = − log Z we obtain Fq = DKL (q||p) + F.
(10)
So, minimizing Fq is equivalent to minimizing the (reverse) Kullback-Leibler divergence3 DKL (q||p). This can be done using the REINFORCE algorithm as described in [2]. The Fq can be estimated by sampling from q(s): M 1 X log q(si ) + βE(si ) , Fq ≈ M i=1
si ∼ q(si ).
(11)
The choice of architecture that would embody the probabilities q(si |s<i ), while crucial, is arbitrary. In the original reference [2], the authors considered both fully connected [12] and convolutional networks with masks that ensure the autoregressive property. In practice, the distribution q(s) will never be equal to p(s), but if it is close enough, we can still use these samples to produce samples from p(s) [3, 13]. One way is to use the samples generated from q(s) as the proposal in the MetropolisHastings algorithm, accepting them with probability p(si+1 ) q(si ) P (si+1 ) q(si ) min 1, = min 1, . (12) p(si ) q(si+1 ) P (si ) q(si+1 ) This step introduces correlations, but for q(s) close to p(s) those correlations will be substantially smaller [4]. Another way is to use importance sampling. Given some observable O, we can calculate its average as X X p(s) ⟨O⟩p ≡ p(s)O(s) = q(s) O(s) q(s) s s (13) M 1 X i i i ≈ w(s)O(s ), s ∼ q(s ) M i=1 where si ∼ q(si ) denotes that the samples si are sampled with probability q(si ). The w(s) are the so called importance weights w(s) =
p(s) . q(s)
(14)
3 F is related to the so-called evidence lower bound (ELBO): ELBO=−F , see for example q q Ref. [11].
4
As we do not know p(s), but only P (s), the final approximation of ⟨O⟩p is PM 1 w̃(si )O(si ) ⟨O⟩p ≈ M 1 i=1 PM i i=1 w̃(s ) M where w̃(s) =
(15)
P (s) = Zw(s). q(s)
(16)
An interesting corollary is that ⟨w̃(s)⟩q(s) = Z ⟨w(s)⟩q(s) = Z which makes
X
q(s)
s
p(s) =Z q(s)
(17)
N
Znis =
1 X e−βE(s) , N i=1 q(s)
s ∼ q(s)
(18)
an unbiased estimator of Z [13]. This makes this approach attractive even for relatively small systems, as it allows for a precise estimation of Z and related quantities such as mutual information and entanglement entropy [14, 15, 16] and even the whole reduced density matrix [17]. Those quantities are usually hard to estimate using Monte-Carlo methods. 3. Physics aware neural generators As already pointed out in the Introduction, the approach outlined above is very general. This can be an advantage, as we only need to provide an energy function definition to train an arbitrary model. However, when simulating more restricted models, e.g. with only nearest neighbor interaction and some symmetries, the model has to unnecessarily relearn those properties. For example, the probability q(si |s<i ) depends on all spins sj with j < i. However, for models with nearest neighbor interactions this dependence can be significantly reduced. For a two dimensional system, assuming that we order spins row by row, only two rows contribute due to the Markov property. For example, in Figure 1 q(s27 |s<27 ) depends only on the spins s0 to s7 and s19 to s26 . This restricted dependency is not taken into account by simple VANs as proposed in [2]. Hierarchical autoregressive networks introduced in [7] alleviate this problem to some extent, resulting in superior training results. This approach also restores some of the translational symmetry broken by factorization (5). Other contributions also exploit the multilevel renormalization group structure of nearest neighbor statistical models [8]. Another issue is the symmetries, which play a very important role with regard to the properties of the system, especially an internal symmetry such as Z2 . On the level of the conditional probabilities, the Z2 symmetry requires qi (−s<i ) = 1 − qi (s<i ) 5
(19)
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
s0
s1
s2
s3
s4
s5
s6
s7
s0
s1
s2
0
0
s8
s9
s10
s11
s12
s13
s14
s15
s8
s9
s10
0
0
s16
s17
s18
s19
s20
s21
s22
s23
s16
s17
s18
0
0
s24
s25
s26
s27
s28
s29
s30
s31
s24
s25
s26
0
0
s32
s33
s34
s35
s36
s37
s38
s39
s32
s33
s34
0
0
s40
s41
s42
s43
s44
s45
s46
s47
s40
s41
s42
0
0
s48
s49
s50
s51
s52
s53
s54
s55
s48
s49
s50
0
0
s56
s57
s58
s59
s60
s61
s62
s63
s56
s57
s58
Figure 1: Calculating the p(s27 |s<27 ). Solid black lines indicate the spins that are already fixed. The dashed lines denote the spins that are not yet fixed and should be summed over to obtain the true conditional probability p27 (s<27 ). The dotted lines indicate padding spins that are required for efficient t4 order calculations. The inner contour represents the convolution kernel for t4 order approximation. Other approximations require smaller kernels and respectively smaller padding. This approach does not take into account the periodic boundary conditions in the vertical direction.
6
which is not enforced by the LeakyReLU activation functions typically used in VANs. This symmetry can be partially enforced by averaging probabilities in symmetric configurations [2, 18]. This approach can be extended to other symmetries, such as translational symmetry, which is also broken by factorization (5) [4, 19]. Again, this leads to better training results. The authors of [5, 6] aim to include in the model some of the interactions of the Hamiltonian H(s). In this contribution, we take the same route, but our approach is less rigorous and more general, as explained in the next section. 4. Approximate conditional probabilities Probability distribution p(s) can be factorized in the same way as q(s) p(s0 , s1 , . . . , sN −1 ) = p(s0 )
N −1 Y i=1
p(si |s<i ).
(20)
Our idea is to start with a approximation for each conditional probability p(si |s<i ) ≈ p̃(si |s<i )
(21)
and use this expression as a starting point for training the neural network 1 , n−1 −h − logit(p̃(si |s<i )) 1+e i
q(si |s<i ) =
(22)
where hn−1 is the output of the last layer of the neural network modeling q(s) i and logit is the inverse of the logistic σ function (see Figure A.9) σ(x) ≡
1 , 1 + e−x
y . 1−y
logit(y) ≡ σ −1 (y) = log
(23)
In this way, the model has only to learn the difference between distributions p̃(s) and p(s). If the approximation is good, we can expect reduced training times as well as better final quality of the model. By definition P e−βE(s) p(s≤i ) s p(si |s<i ) ≡ = P >i −βE(s) . (24) p(s<i ) s≥i e The denominator is just a normalization factor so that X p(si |s<i ) = 1.
(25)
si =±1
Combining this with the definition of logit, we obtain that P
−βE(s)
logit pi (s<i ) = log P
e−βE(s)
s>i e
s>i
7
si =1 si =−1
.
(26)
The sums in the above expression can be further decomposed as X X e−βE(s) = e−βI(s<i ) e−βE({s≥i }) s>i
(27)
s>i
where I({s}) denotes the internal energy of the set of spins {s}, i.e. , the energy of all the bounds directly connecting two spins in {s}. The E({s}) denotes the energy of all the bounds originating from any spin in the set {s}. The first factor does not depend on si , so it will cancel out in (26), leaving P logit pi (s<i ) = log P
s>i e
s>i
−βE(s≥i ) si =1
e−βE(s≥i )
.
(28)
si =−1
The formula (28) will be our starting point for approximations. In this contribution, we consider nearest neighbor spin glass systems with energy given by X E(s) = − Ji,j si sj (29) <i,j>
where Ji,j are some arbitrary link variables. For Ji,j = 1 this reduces to the ferromagnetic Ising model as a special case. We will also use another representation for the Boltzmann factor (4). As spins can take only values ±1, we have ( eβJij si sj = 1 βJij si sj e = (30) e−βJij si sj = −1 which allows us to write 1 1 (1 + si sj )eβJij + (1 − si sj )e−βJij 2 2 1 βJij 1 −βJij +e ) + si sj (eβJij − e−βJij ) = (e 2 2 = cosh(βJij ) + si sj sinh(βJij )
eβJij si sj =
(31)
= cosh(βJij ) (1 + si sj tanh(βJij )) . When Ji,j = 1 this reduces to cosh(β) (1 + si sj tanh(β)) .
(32)
We will approximate p(si |s<i ) taking only a subset of spins from s<i and sum only over a subset of spins from s>i , neglecting the dependence on any other spins. As an example, we will use the probability p27 (s<27 ) (see Figure 1). In this section, we will concentrate on the Ising model setting Ji,j = 1. The formulas for the general spin glass case will be given in Section 6.
8
The simplest approximation would be to consider only the nearest neighbor interactions with spins s19 and s26 (black links in Figure 1) which gives the expression eβ(s19 +s26 ) (1) logit p̃27 (s<27 ) = log −β(s +s ) = 2β(s19 + s26 ). 19 26 e
(33)
This corresponds to the approach taken in [5, 6]. Please note, however, that those references use different neural network architectures, so the results are not directly comparable. In the next approximation, we add the dependence on spin s20 and sum over spin s28 . Using the representation (32), we obtain X eβ(s27 s28 +s28 s20 ) s28
= (cosh β)
2
(34) X 2 (1 + s27 s28 t)(1 + s28 s20 t) = (cosh β) (1 + s27 s20 t2 ), s28
where t = tanh β and
(2)
logit p̃27 (s19 , s20 , s26 ) = 2β(s19 + s26 ) + log
1 + s20 t2 1 − s20 t2
(35)
2
≈ 2β(s19 + s26 ) + 2s20 t .
For the order t3 calculations, we add the dependence on spin s21 and additionally sum over the spins s29 , s34 and s35 . This requires the calculation of the sum X eβ(s26 s34 +s34 s35 +s35 s27 +s27 s28 +s28 s20 +s28 s29 +s29 s21 ) s28 ,s29 ,s34 ,s35
=
X
eβ(s26 s34 +s34 s35 +s35 s27 )
s34 ,s35
X
eβ(s28 s20 +s28 s29 +s29 s21 ) .
(36)
s28 ,s29
Going through the same steps as before, we obtain (3)
logit p̃27 (s19 , s20 , s21 , s26 ) = 2β(s19 + s26 ) + 2s20 t2 + 2(s21 + s26 )t3 .
(37)
Finally, we present the results for the approximation t4 , where we add dependence on spins s22 and s25 and additionally sum over spins s30 , s33 and s36 (4)
logit p̃27 (s19 , s20 , s21 , s22 , s25 , s26 ) = 2β(s19 + s26 ) + 2s20 t2 + 2(s21 + s26 )t3
(38)
+ 2(s20 + s22 + s25 )t4 . Those calculations were performed using the Mathematica software. Formulae (33), (35), (37), (38) can be easily translated to match any other spin, but they do not take into account the periodic boundary conditions in the vertical 9
β 0.40 t0 t1 t2 t3 t4 true
−874.490(7) −886.705(5) −891.059(4) −892.393(4) −900.478
βc 0.50 -709.783 −891.824(7) −910.720(8) −970.19(1) −906.461(6) −928.152(6) −997.307(8) −912.173(5) −935.544(5) −1011.365(8) −914.114(4) −938.308(5) −1017.941(7) −924.4135 −952.648 −1051.105 0.42
Table 1: Values of Fq obtained with different conditional approximations for the Ising model on the 32 × 32 lattice, t0 denotes random distribution where each configuration has same probability which is a reasonable approximation for an untrained model.
direction. That means that the formulas for the spins in the last rows are, strictly speaking, not correct. However, as those formulas are approximations anyway, we just assume that they are steps in a good direction and favor the ease of implementation over correctness. One way of testing the quality of these approximations is to calculate the Fq for each of them and compare it with the true value that can be obtained using formulas from [20]. Those numbers are presented in Table 1. Fq was estimated using the formulas (11) and (22) setting hn−1 = 0. It is worth noting i in Table 1 that already t1 provides a noticeable improvement with respect to an untrained neural network for all the temperatures considered. The value for t0 was calculated analytically, assuming a uniform distribution q(s) = 2−N , where N is the number of spins. This was done using another representation for Fq : X X Fq = q(s) log q(s) + β q(s)E(s) = β ⟨E(s)⟩q − Sq , (39) s
s
where Sq is the entropy of the distribution q(s). For a uniform distribution, this equals N log(2). If one colors the spins in a checkerboard fashion and flips every white spin, the energy of the configuration changes sign. As both flipped and non-flipped configurations are equally probable in this case, we obtain a zero average energy, leading to the final result, Fq = −N log(2), presented in Table 1. 5. Results – Ising model We started by testing our approach on a ferromagnetic Ising model on a 32 × 32 lattice. In the absence of the magnetic field, this model has a Z2 symmetry associated with the change of sign of each spin in the configuration, as such a transformation obviously leaves the energy (29) unchanged. In the thermodynamical limit, this model exhibits a continuous phase transition involving spontaneous breaking of this symmetry. In the high temperature symmetric phase, the net magnetization 1 X si (40) ⟨m(s)⟩p(s) , where m(s) = 2 L i 10
= 0.40 = 0.42 = c = 0.50
8 7 6 5 4 3 2 1 0
1.00
0.75
0.50
0.25 0.00
m
0.25
0.50
0.75
1.00
Figure 2: Magnetization distribution for the values of β considered in this reference. Obtained from Monte-Carlo simulations.
is zero. Below the critical temperature, in the broken phase, a non-zero spontaneous magnetization appears. For the finite systems considered here, there is formally no phase transition and the average magnetization should always be zero. The transition, however, manifests itself through the presence of two peaks in the magnetization histogram (see Figure 2) and for low enough temperatures those peaks become totally separated. Monte Carlo simulation may become trapped in one of those regions with non-zero magnetization. The same problem can appear in neural samplers and is known as mode collapse. We have chosen four different values of β: (0.4, 0.42, βc , 0.5), where βc =
√ 1 log(1 + 2) ≈ 0.4406868, 2
(41)
is the inverse critical temperature in the thermodynamic limit. We have chosen them to cover the symmetric phase (β = 0.40), the transition region (β = 0.42) and the broken phase (β = βc and β = 0.5) (see Figure 2). We measure the quality of training using the estimated sample size (ESS) [21, 22], 2
ESS =
⟨w̃⟩q
⟨w̃2 ⟩q
(42)
and by comparing different estimates of F . We also track magnetization to see if the Z2 symmetry is respected. For our investigations, we have chosen a very simple architecture consisting of two layers of fully connected neural networks with a leaky ReLU activation layer between them. As we wanted to assess only the influence of the probability priors, we did not add any enhancements related to symmetries or β annealing. 11
t0 t1 t2 t3 t4
0.75
0.4
0.50 0.3
t0 t1 t2 t3 t4
ESS
0.2
0.25
m 0.00 0.25 0.50
0.1
0.75 0.0
0
500
1000
era
1500
2000
1.00
2500
0
500
1000
era
1500
2000
2500
Figure 3: History of the training for Ising model on 32 lattice at critical β. ESS (left) and m (right). The horizontal axis marks the eras where an era is a 100 weights update using the batch of 4096 samples. For clarity we present a moving average over 100 weights updates.
All the models were implemented using the PyTorch library. More details are presented in Appendix A. The summary of the results is presented in Appendix B in Tables B.4 to B.6. The Monte Carlo simulations were performed using the Wolff cluster algorithm [23]. In Figure 3 we present the history of the training for the model at critical temperature. We can see that the inclusion of approximate probabilities has a significant impact on the efficiency of the training. In particular, there is a notable change between the t1 and t2 approximations. The dip in ESS for the t2 approximation corresponds to the restoration of the Z2 symmetry, as can be seen in the plot showing magnetization in the right panel of Figure 3. To test the quality of the trained models, we have generated 220 configurations for each model and calculated Fq and the estimate derived from (18) Fnis = − log Znis .
(43)
Using the configurations obtained from Monte-Carlo simulations, we have also calculated another estimator of F [24] Fmc = − log Zmc . where −1 Zmc =
q(s) e−βE(s)
(44)
M
1 X q(si ) ≈ , M i=1 e−βE(si ) p
si ∼ p(si ).
(45)
It can be proven [24] that Fnis ≥ F ≥ Fmc .
(46)
The discrepancy between those two estimators gives us an estimate of the quality of training and is an indicator of mode collapse. This can be quantified by 12
Fq Fnis Fmc Ftrue
950.5
Fq Fnis Fmc Ftrue
1049.0
951.0
1049.5
951.5
1050.0
952.0
1050.5
952.5
1051.0
t0
t1
t4
t2 t3 Approximation
t0
t1
t2 t3 Approximation
t4
Figure 4: Various estimates of F for the Ising model at critical β (left) and β = 0.5 (right). Uncertainties are much smaller than the points size.
introducing a parameter [24] w̄ ≡
Znis = e(Fmc −Fnis ) . Zmc
(47)
The values of this parameter can be found in Table B.4 in Appendix B. The various estimators of F for β = βc are compared in Figure 4 (left). We can see that starting with the approximation t2 , Fnis and Fmc coincide. The expansion in t that we use to derive the priors formally is divergent above the critical β. However, looking at the results for β = 0.5 in Figure 5, we again see a significant improvement with the increasing order of the approximation. Actually, only the t4 approximation was able to train properly and restore the Z2 symmetry. Looking closely at the estimators of F in Figure 4 (right), we observe that there remains a slight discrepancy between Fnis and Fmc resulting in w̄ slightly below one and a small bias of |M | obtained using NIS compared to Monte Carlo simulations (see Tables B.4 and B.6). 6. Edwards-Anderson model In this section, we will consider the nearest neighbor spin glass model [10] with energy given by (29), where now the link variables take on the values J = ±1 with equal probability. Spin glass models generally do not exhibit spontaneous symmetry breaking, and the phase transition, if present, is rather associated with replica symmetry breaking [10, 25]. The two dimensional model that we simulate does not exhibit a phase transition at finite β [26], but it is increasingly difficult to simulate with growing β. By denoting ti,j ≡ tanh(βJi,j ) and following the same derivation we have used in Section 4, we obtain the approximations for a general spin glass model 13
1.00
t0 t1 t2 t3 t4
0.8 0.7 0.6
0.75 0.50
0.5
0.25
ESS0.4
m 0.00
0.3
0.25
0.2
0.50
0.1
0.75
0.0
0
500
1000
era
1500
2000
1.00
2500
t0 t1 t2 t3 t4
0
500
1000
era
1500
2000
2500
Figure 5: History of the training for the Ising model at β = 0.5.
with nearest-neighbor interactions: (1)
logit p̃27 (s19 , s26 ) = 2β(s19 J19,27 + s26 J26,27 )
(48)
(2) (1) logit p̃27 (s19 , s20 , s26 ) = logit p̃27 (s19 , s26 ) + 2s20 t20,28 t28,27 (3) (2) logit p̃27 (s19 , s20 , s21 , s26 ) = logit p̃27 (s19 , s20 , s26 )+
(49)
2(s21 t21,29 t29,28 t28,27 + s26 t26,34 t34,35 t35,27 ) (50) (4) (3) logit p̃27 (s19 , s20 , s21 , s22 , s25 , s26 ) = logit p̃27 (s19 , s20 , s21 , s26 )+
2s22 t22,30 t30,29 t29,28 t28,27 + 2s25 t25,33 t33,34 t34,35 t35,27 + 2s20 t20,28 t28,36 t36,35 t35,27 .
(51)
In the case of link variables Ji,j taking only the values ±1 ti,j = tanh(βJi,j ) = Ji,j tanh(β) = Ji,j t.
(52)
and t remains the expansion parameter. In Table 2 we present the values of Fq obtained without any training, using only approximate probabilities. All the values are calculated for the same specific choice of link variables J. One should note that for β = 0.9 the value of Fq for the t4 approximation is higher then the value for the t3 approximation, indicating the possible divergence of the series in t. Those values were obtained in the same way as in the Ising model case, with the notable exception that the true value is unknown. 7. Results - Edwards-Anderson model We have trained models for one particular set of link variables Ji,j and three different values of β = (0.3, 0.6, 0.9). Because there are no analytical results for this model, the results were compared with Monte-Carlo simulations performed 14
β 0.30 t0 t1 t2 t3 t4
−794.507(3) −797.584(3) −798.043(2) −798.087(2)
0.6 -709.783 −1007.76(1) −1034.374(7) −1043.466(6) −1045.100(5)
0.9 −1282.50(2) −1354.14(1) −1379.10(1) −1378.11(1)
Table 2: Values of Fq obtained with different conditional approximations for Edward-Anderson model on 32 × 32 lattice for one specific choice of J, t0 denotes random distribution where each configuration has same probability which is a reasonable approximation for an untrained model.
t0 t1 t2 t3 t4
0.5 0.4 0.3
0.01
0.00
m 0.01
ESS
0.2 0.02 0.1 0.03
0.0
0
1000
2000
era
3000
4000
5000
0
1000
2000
era
3000
4000
t0 t1 t2 t3 t4 5000
Figure 6: Training history for the Edwards-Anderson model on 32 lattice at β = 0.6. ESS (left) and m (right). The horizontal axis marks the eras where an era is a 100 weights update using the batch of 4096 samples. For clarity we present a moving average over 100 weights updates.
using the parallel tempering technique [27]. In this section we present the results only for β = 0.6 and β = 0.9, the summary of all the results is presented in Appendix B in Tables B.7 to B.9. Looking at Figure 6, where the training history for the β = 0.6 is presented, we observe similar behavior as for the Ising model. We see a big increase in the efficiency of training from t2 to t3 . Looking at the magnetization in the same Figure (right) and estimates of F in the Figure 7 (left), we notice that the model is able to train starting with the t2 approximation, albeit with significantly lower ESS (see also Tables B.7 to B.9). The situation is different for β = 0.9, as presented in Figure 8. While we still see a fourfold increase in ESS from t1 to t3 , it is still at the level of only 5%. Looking at Figure 7 (right) and into Tables B.7 and B.9, we can see that none of the models trained properly, each exhibiting a mode collapse as indicated by the value of parameter w̄. This mode collapse seems to be associated with broken Z2 symmetry, as we have a small but non-zero average magnetization.
15
1055.0
Fq Fnis Fmc
1055.5
Fq Fnis Fmc
1416 1417 1418
1056.0
1419
1056.5
1420 1057.0
1421
1057.5
1422
1058.0
1423
t0
t1
t1
t0
t4
t2 t3 Approximation
t4
t2 t3 Approximation
Figure 7: Different estimates of F for the Edwards-Anderson model at β = 0.6 (left) and β = 0.9 (right).
0.00
t0 t1 t2 t3 t4
0.08
0.06
0.01
0.02
m
ESS 0.04
0.03
0.02
t0 t1 t2 t3 t4
0.04
0.00 0
2000
4000
era
6000
8000
10000
0.05
0
2000
4000
era
6000
8000
10000
Figure 8: Results for the Edwards-Anderson model on 32 lattice at β = 0.9. ESS (left) and m (right). The horizontal axis marks the eras where an era is a 100 weights update using the batch of 4096 samples. For clarity we present a moving average over 200 weights updates.
16
We should also mention that the best results are obtained for t3 , which is in agreement with the results in Table 2. 8. Conclusions and outlook In this contribution, we have advocated for the use of probability priors to enhance the autoregressive neural generators for spin systems. The idea is to approximate the conditional probabilities by using only a subset of neighboring spins. This can be done in a systematic way using a weak expansion in tanh β. We have derived formulas up to order t4 for the general case of nearest neighbor interactions. We have tested our approach on the two-dimensional Ising model and J = ±1 Edwards-Anderson spin glass. In both cases, the results are very encouraging, the addition of priors significantly improves the quality of training. While the references [5] and [6] also introduce the priors, albeit using a different architecture, they only use the nearest neighbors in the approximation. Our results prove that including further spins can significantly improve the training, sometimes making the difference between an architecture that can be trained and one that cannot. We have derived the approximate probabilities using an expansion in powers of tanh β and this can be done for any nearest neighbors model. However, those approximations do not have to be derived in such a systematic way. Any approximation that is closer to the target distribution than the "blank slate" approach can potentially be beneficial. The search for such approximations for different physical models is the subject of ongoing work. We consider this paper as a proof of concept, and to this end, we have used only a very simple neural network architecture. We have also omitted by design any possible enhancements, like incorporating the Z2 symmetry or β annealing, which have been proven to increase the quality of the training. All those things are orthogonal to the approach proposed here and can be used jointly. In particular, the priors can be incorporated into any autoregressive architecture. We have already incorporated the priors into a transformer based architecture with similarly promising results [28]. Acknowledgments We gratefully acknowledge Polish high-performance computing infrastructure PLGrid (HPC Center: ACK Cyfronet AGH) for providing computer facilities and support within computational grant no. PLG/2025/018811. T.S. and D.Z. kindly acknowledge the support of the Polish National Science Center (NCN) Grant No. 2021/43/D/ST2/03375. P.K. acknowledges the support of the Polish National Science Center (NCN) grant No. 2022/46/E/ST2/00346. References [1] N. Metropolis, A. W. Rosenbluth, M. N. Rosenbluth, A. H. Teller, E. Teller, Equation of state calculations by fast computing machines, The Journal 17
of Chemical Physics 21 (1953) 1087–1092. doi:10.1063/1.1699114. W. K. Hastings, Monte Carlo sampling methods using Markov chains and their applications, Biometrika 57 (1970) 97–109. doi:10.1093/biomet/57.1.97. [2] D. Wu, L. Wang, P. Zhang, Solving statistical mechanics using variational autoregressive networks, Phys. Rev. Lett. 122 (2019) 080602. [3] M. S. Albergo, G. Kanwar, P. E. Shanahan, Flow-based generative models for markov chain monte carlo in lattice field theory, Phys. Rev. D 100 (2019) 034515. doi:10.1103/PhysRevD.100.034515. [4] P. Białas, P. Korcyl, T. Stebel, Analysis of autocorrelation times in neural Markov chain Monte Carlo simulations, Phys. Rev. E 107 (2023) 015303. doi:10.1103/PhysRevE.107.015303. arXiv:2111.10189. [5] I. Biazzo, The autoregressive neural network architecture of the boltzmann distribution of pairwise interacting spins systems, Communications Physics 6 (2023) 1–10. doi:10.1038/S42005-023-01416-5;. [6] I. Biazzo, D. Wu, G. Carleo, Sparse autoregressive neural networks for classical spin systems, Machine Learning: Science and Technology 5 (2024) 025074. doi:10.1088/2632-2153/ad5783. [7] P. Białas, P. Korcyl, T. Stebel, Hierarchical autoregressive neural networks for statistical systems, Comput. Phys. Commun. 281 (2022) 108502. doi:10.1016/j.cpc.2022.108502. arXiv:2203.10989. [8] A. Singha, E. Cellini, K. A. Nicoli, K. Jansen, S. Kühn, S. Nakajima, Multilevel generative samplers for investigating critical phenomena, 2025. arXiv:2503.08918. [9] L. M. Del Bono, F. Ricci-Tersenghi, F. Zamponi, Nearest-neighbors neural network architecture for efficient sampling of statistical physics models, Machine Learning: Science and Technology 6 (2025) 025029. doi:10.1088/26322153/adcdc1. arXiv:2407.19483. [10] S. F. Edwards, P. W. Anderson, Theory of spin glasses, Journal of Physics F: Metal Physics 5 (1975) 965. doi:10.1088/0305-4608/5/5/017. [11] D. Blessing, X. Jia, J. Esslinger, F. Vargas, G. Neumann, Beyond ELBOs: A Large-Scale Evaluation of Variational Methods for Sampling, arXiv eprints (2024). arXiv:2406.07423. [12] M. Germain, K. Gregor, I. Murray, H. Larochelle, Made: Masked autoencoder for distribution estimation, 2015. arXiv:1502.03509. [13] K. A. Nicoli, S. Nakajima, N. Strodthoff, W. Samek, K.-R. Müller, P. Kessel, Asymptotically unbiased estimation of physical observables with neural samplers, Phys. Rev. E 101 (2020) 023304.
18
[14] P. Białas, P. Korcyl, T. Stebel, D. Zapolski, Rényi entanglement entropy of a spin chain with generative neural networks, Phys. Rev. E 110 (2024) 044116. doi:10.1103/PhysRevE.110.044116. [15] P. Białas, P. Korcyl, T. Stebel, Mutual information of spin systems from autoregressive neural networks, Phys. Rev. E 108 (2023) 044140. doi:10.1103/PhysRevE.108.044140. [16] A. Bulgarelli, E. Cellini, K. Jansen, S. Kühn, A. Nada, S. Nakajima, K. A. Nicoli, M. Panero, Flow-based sampling for entanglement entropy and the machine learning of defects, Phys. Rev. Lett. 134 (2025) 151601. doi:10.1103/PhysRevLett.134.151601. [17] P. Białas, P. Korcyl, T. Stebel, D. Zapolski, Estimation of the reduced density matrix and entanglement entropies using autoregressive networks, 2025. arXiv:2506.04170. [18] S.-H. Li, L. Wang, Neural network renormalization group, Phys. Rev. Lett. 121 (2018) 260601. doi:10.1103/PhysRevLett.121.260601. [19] Simulating first-order phase transition with hierarchical autoregressive networks, Phys. Rev. E 107 (2023) 054127. doi:10.1103/PhysRevE.107.054127. arXiv:2212.04955. [20] A. E. Ferdinand, M. E. Fisher, Bounded and inhomogeneous ising models. i. specific-heat anomaly of a finite lattice, Phys. Rev. 185 (1969) 832–846. doi:10.1103/PhysRev.185.832. [21] A. Kong, A note on importance sampling using standarized weights, University of Chicago Technical Reports (1992). [22] J. S. Liu, Metropolized independent sampling with comparisons to rejection sampling and importance sampling, Statistics and Computing 6 (1996) 113–119. doi:10.1007/BF00162521/METRICS. [23] U. Wolff, Collective monte carlo updating for spin systems, Phys. Rev. Lett. 62 (1989) 361–364. doi:10.1103/PhysRevLett.62.361. [24] K. A. Nicoli, C. J. Anders, T. Hartung, K. Jansen, P. Kessel, S. Nakajima, Detecting and mitigating mode-collapse for flow-based sampling of lattice field theories, Phys. Rev. D 108 (2023) 114501. doi:10.1103/PhysRevD.108.114501. [25] G. Parisi, Order parameter for spin-glasses, Physical Review Letters 50 (1983) 1946. doi:10.1103/PhysRevLett.50.1946. [26] I. Morgenstern, K. Binder, Magnetic correlations in two-dimensional spinglasses, Phys. Rev. B 22 (1980) 288–303. doi:10.1103/PhysRevB.22.288.
19
q0
q1
q2
q3
σ
σ
σ
σ
+
+
+
+
h1 0
h1 1
h1 2
h1 3
LReLU
LReLU
LReLU
LReLU
h0 0
h0 1
h0 2
h0 3
s0
s1
s2
s3
l0
l1
l2
l3
Figure A.9: Schema of the model architecture. The neural network block (light blue) can be exchanged for any autoregressive model. The li are the approximate logits.
[27] K. Hukushima, K. Nemoto, Exchange monte carlo method and application to spin glass simulations, Journal of the Physical Society of Japan 65 (1996) 1604–1608. doi:10.1143/JPSJ.65.1604. [28] P. Białas, P. Korcyl, T. Stebel, A. Stefański, D. Zapolski, Sampling twodimensional spin systems with transformers, 2026. arXiv:2604.27738. Appendix A. Implementation As our aim was only to provide the concept of proof for the probability priors, we used the simple architecture, consisting of two dense layers with LeakyReLu activation function in between. The logits of approximate probabilities are added to the output of the last layer and then passed through the logistic function (see Figure A.9). The most computationally demanding part of the training is sampling the configurations, as we need to calculate each q(si |s<i ) separately in order. The naive approach would be to propagate a whole vector of spins through the network each time (with the first i spins set to already defined values) and take only the output i, disregarding the rest. That, however, would be very wasteful. The better approach is to cache all the previous outputs of each layer and use them to calculate only the output i at each stage. This amounts to effectively multiplying the input by one row of the weight matrix at each step. That greatly reduces the amount of time needed for sampling. When computing the probability of a configuration, which is needed for loss calculations, we just pass the whole input vector through the network at once. 20
Ising EA
t0 1.0 1.0
t1 1.5 1.5
t2 1.5 1.8
t3 1.5 2.0
t4 1.5 2.2
Table A.3: Relative increase in time needed for one epoch of training (batch size = 4096).
.. . · · · J10,18 f [0] = 2β · · · J18,26 · · · J26,34 .. .
.. .
.. .
J11,19
J12,20
J19,27
J20,28
J27,35
J28,36
.. .
.. .
.. . ···
· · · J17,18
···
f [1] = 2β · · · J25,26
···
· · · J33,34 .. .
.. .
.. .
J18,19
J19,20
J26,27
J27,28
J34,35
J35,36
.. .
.. .
··· ··· ···
Figure A.10: Factors needed for t1 approximation in Edwards-Anderson model.
This optimization implies that we also need two ways of calculating the approximate probabilities: one that does it one at a time and the other all at the same time. For the Ising model, this is relatively easy. The formulas (33),(35),(37) and (38) can be written in terms of a kernel multiplication (see Figure 1). So we can use the convolution if we need to calculate all logits at once. Otherwise, we just apply the kernel at one position in the spins lattice. This calculation introduces a 50% overhead (see Table A.3). In the case of the Edwards-Anderson model, the situation is more involved: we cannot use the convolutions, as the kernel would depend on the position of the spin due to different link variables at each location. We start by precalculating all the products of link variables required in equations (48) to (51). For example, for the t2 approximation, spin i depends on three previous spins and so we need three tables of precalculated factors. Portions of those tables are presented in Figures A.10 and A.11. With those tables, we can next calculate the required logit by simple multiplication. For example, the logit for the s27 spin would be calculated as f [0 ,3 ,3]* s [: ,3 ,5] + f [1 ,3 ,3]* s [: ,4 ,4] + f [2 ,3 ,3]* s [: ,4 ,6]
assuming that s represents an array of spins as presented in Figure 1. The first dimension is the batch dimension. All the digits can be calculated at the same time using the following line of code f [0]* s [: ,0:8 ,2:10] + f [1]* s [: ,1:9 ,1:9] + f [2]* s [1:9 ,3:11]
Those calculations introduce an overhead of up to 220% (see Table A.3).
21
.. . · · · J11,19 J19,18 f [2] = 2t2 · · · J19,27 J27,26 · · · J27,35 J35,34 .. .
.. .
.. .
J12,20 J20,19
J13,21 J21,20
J20,28 J28,27
J21,29 J29,28
J28,36 J36,35
J29,37 J37,36
.. .
.. .
··· ··· ···
Figure A.11: Factors needed for t2 approximation in Edwards-Anderson model.
Appendix B. Results In this section we present the results for estimators of F and the average energy, magnetization and absolute value of magnetization. Each value was obtained using 220 samples from neural models and approximately 4 × 105 samples from Monte Carlo. The errors were calculated using the bootstrap method.
22
Fq Fnis Fmc β = 0.40000 F = −900.4783078 0.003(4) −897.999(2) −899.85(4) −900.48(1) 0.152(4) −899.547(1) −900.480(2) −900.476(4) 0.505(2) −900.1351(8) −900.4777(9) −900.475(2) 0.8139(3) −900.3745(5) −900.4777(4) −900.4767(9) 0.8303(3) −900.3844(4) −900.4777(4) −900.4770(8) β = 0.42000 F = −924.4135102 0.005(5) −921.899(2) −923.66(2) −924.42(1) 0.298(2) −923.794(1) −924.413(2) −924.414(2) 0.459(2) −924.0173(8) −924.414(1) −924.413(2) 0.6676(7) −924.2108(7) −924.4137(7) −924.414(1) 0.6798(6) −924.2207(6) −924.4138(6) −924.414(1) β = 0.44069 F = −952.6480795 0.004(4) −950.298(2) −951.88(3) −952.65(1) 0.014(8) −950.877(1) −951.954(9) −952.644(6) 0.31(1) −952.155(1) −952.647(2) −952.649(2) 0.35(1) −952.2659(8) −952.647(1) −952.648(2) 0.381(6) −952.2667(9) −952.647(1) −952.647(2) β = 0.50000 F = −1051.1049876 0.032(4) −1048.722(2) −1050.381(5) −1051.12(1) 0.13(6) −1049.273(2) −1050.418(3) −1051.103(7) 0.31(1) −1049.334(2) −1050.413(1) −1051.109(6) 0.32(2) −1049.340(2) −1050.412(1) −1051.115(7) 0.770(8) −1050.9674(4) −1051.0586(6) −1051.1042(8) ESS
t0 t1 t2 t3 t4 t0 t1 t2 t3 t4 t0 t1 t2 t3 t4 t0 t1 t2 t3 t4
w̄ 0.53(3) 1.003(4) 1.002(2) 1.0010(9) 1.0007(9) 0.47(1) 0.998(3) 1.000(2) 1.000(1) 1.000(1) 0.46(2) 0.501(5) 0.998(3) 0.999(2) 1.000(2) 0.479(5) 0.504(4) 0.499(3) 0.495(4) 0.955(1)
Table B.4: Different estimates of F for the Ising model. The value beside the β is the true value of free energy.
23
none
t1
t2
t3
t4
none
t1
t2
t3
t4
E M |M | E M |M | E M |M | E M |M | E M |M | E M |M | E M |M | E M |M | E M |M | E M |M |
variational NIS β = 0.40 −1142.12(8) −1134.7(8) −146.2(2) −94(2) 228.5(1) 205(2) −1129.71(7) −1133.8(1) 4.8(3) 0.2(6) 206.9(1) 206.2(3) −1129.28(8) −1133.81(9) −1.8(2) −0.5(4) 201.7(1) 206.0(2) −1131.83(6) −1134.00(9) −1.3(3) 0.0(3) 202.9(1) 206.3(2) −1131.76(8) −1133.92(8) 0.8(3) 0.0(3) 203.4(1) 205.9(2) β = 0.42 −1295.2(1) −1271.9(9) −411.5(3) −302(4) 437.2(2) 368(3) −1262.4(1) −1267.9(1) 5.8(4) −0.7(8) 346.2(2) 359.5(3) −1260.2(1) −1268.2(1) −3.2(4) 1.0(5) 345.1(2) 360.1(3) −1260.89(9) −1268.03(9) 0.9(4) −0.4(5) 341.2(2) 359.7(2) −1260.6(1) −1267.97(8) −0.1(3) 0.4(5) 339.9(2) 359.8(2)
Table B.5: Ising model. Part I.
24
MC −1133.9(2) 0.3(4) 206.2(2) −1133.9(2) 0.3(4) 206.2(2) −1133.9(2) 0.3(4) 206.2(2) −1133.9(2) 0.3(4) 206.2(2) −1133.9(2) 0.3(4) 206.2(2) −1267.9(2) 1.4(7) 359.9(3) −1267.9(2) 1.4(7) 359.9(3) −1267.9(2) 1.4(7) 359.9(3) −1267.9(2) 1.4(7) 359.9(3) −1267.9(2) 1.4(7) 359.9(3)
none
t1
t2
t3
t4
none
t1
t2
t3
t4
E M |M | E M |M | E M |M | E M |M | E M |M | E M |M | E M |M | E M |M | E M |M | E M |M |
variational NIS β = 0.44 −1494.26(9) −1484(1) −739.28(9) −705(1) 739.26(9) 703(1) −1479.0(1) −1468.9(9) 713.5(1) 673(3) 713.7(1) 674(3) −1468.5(1) −1467.9(2) 3.8(8) 1(1) 674.3(2) 669.6(5) −1477.8(1) −1468.2(2) −0.5(8) 2.0(8) 693.8(2) 670.4(5) −1472.7(1) −1468.1(1) −6.1(7) −1(1) 685.4(2) 670.1(4) β = 0.50 −1795.07(5) −1783(5) −939.07(2) −923(1) 939.07(2) 924(9) −1784.78(5) −1787.6(2) −934.24(2) −933.2(2) 934.24(2) 933.2(2) −1781.68(5) −1787.4(1) 932.07(2) 933.25(8) 932.07(2) 933.26(8) −1785.59(5) −1787.6(1) 931.29(3) 933.2(1) 931.29(3) 933.2(1) −1791.21(6) −1788.12(6) 0(1) 0(1) 935.39(2) 933.58(5)
Table B.6: Ising model. Part II
25
MC −1468.2(1) 1(1) 670.3(3) −1468.2(1) 1(1) 670.3(3) −1468.2(1) 1(1) 670.3(3) −1468.2(1) 1(1) 670.3(3) −1468.2(1) 1(1) 670.3(3) −1787.41(9) 0(2) 933.15(4) −1787.41(9) 0(2) 933.15(4) −1787.41(9) 0(2) 933.15(4) −1787.41(9) 0(2) 933.15(4) −1787.41(9) 0(2) 933.15(4)
Fnis Fmc β = 0.30000 0.006(1) −798.082(2) −800.86(1) −800.86(3) 0.5686(9) −800.5702(7) −800.8526(9) −800.851(2) 0.9165(1) −800.8074(3) −800.8509(3) −800.8515(6) 0.95606(6) −800.8290(2) −800.8515(2) −800.8515(4) 0.96558(5) −800.8337(2) −800.8512(2) −800.8511(4) β = 0.60000 0.006(3) −1055.045(2) −1057.43(2) −1058.18(7) 0.065(7) −1056.198(2) −1057.508(4) −1058.21(1) 0.173(3) −1057.317(1) −1058.198(2) −1058.198(5) 0.463(2) −1057.8139(7) −1058.198(1) −1058.197(2) 0.506(1) −1057.8580(9) −1058.1972(9) −1058.199(2) β = 0.90000 0.0004(2) −1415.451(2) −1418.24(5) −1423.5(2) 0.012(9) −1420.943(2) −1422.53(1) −1423.45(2) 0.048(4) −1421.178(2) −1422.517(5) −1423.44(1) 0.052(7) −1421.464(2) −1422.749(4) −1423.43(3) 0.040(4) −1421.259(2) −1422.645(4) −1423.47(1) ESS
t0 t1 t2 t3 t4 t0 t1 t2 t3 t4 t0 t1 t2 t3 t4
Fq
Table B.7: Different estimates of F for Edwards-Anderson model.
26
w̄ 1.00(3) 1.001(2) 0.9994(6) 1.0000(4) 1.0001(4) 0.47(3) 0.495(5) 1.000(5) 1.001(2) 0.998(2) 0.006(1) 0.400(8) 0.396(5) 0.51(2) 0.439(6)
none
t1
t2
t3
t4
none
t1
t2
t3
t4
E M |M | E M |M | E M |M | E M |M | E M |M | E M |M | E M |M | E M |M | E M |M | E M |M |
variational NIS β = 0.30 −592.83(5) −599.8(5) −0.67(4) 0.2(5) 25.00(2) 25.3(2) −598.58(4) −599.62(6) 0.08(3) −0.03(4) 24.98(2) 25.02(2) −599.37(4) −599.66(4) −0.10(3) −0.05(3) 24.99(2) 24.98(2) −599.60(4) −599.64(5) −0.08(3) 0.03(3) 24.99(2) 24.98(2) −599.80(4) −599.74(4) 0.06(3) 0.05(3) 24.99(2) 24.99(2) β = 0.60 −1080.63(3) −1080.4(6) 3.23(3) 4.7(3) 25.00(2) 25.0(3) −1079.02(3) −1080.6(1) 4.45(3) 4.3(1) 25.08(2) 25.04(6) −1078.66(3) −1080.75(8) −0.15(3) −0.11(7) 25.20(2) 25.24(4) −1079.74(3) −1080.78(5) −0.53(3) −0.01(5) 25.16(2) 25.15(3) −1079.94(4) −1080.70(5) −0.45(3) 0.00(4) 25.17(2) 25.24(3)
MC −599.64(5) 0.00(3) 25.00(2) −599.64(5) 0.00(3) 25.00(2) −599.64(5) 0.00(3) 25.00(2) −599.64(5) 0.00(3) 25.00(2) −599.64(5) 0.00(3) 25.00(2) −1080.78(4) −0.04(4) 25.19(2) −1080.78(4) −0.04(4) 25.19(2) −1080.78(4) −0.04(4) 25.19(2) −1080.78(4) −0.04(4) 25.19(2) −1080.78(4) −0.04(4) 25.19(2)
Table B.8: Edwards-Anderson model. Part I.
27
none
t1
t2
t3
t4
E M |M | E M |M | E M |M | E M |M | E M |M |
variational NIS β = 0.90 −1319.91(2) −1322(3) −9.60(3) −9(2) 24.36(2) 23(3) −1317.72(2) −1317.6(1) −17.13(3) −16.4(2) 28.56(2) 28.5(1) −1317.03(2) −1317.4(1) −17.26(3) −16.4(3) 28.98(2) 29.0(1) −1316.32(2) −1316.2(5) −15.27(3) −14.6(3) 28.48(2) 28.2(2) −1317.08(2) −1317.5(1) −9.39(4) −8.5(1) 28.31(2) 28.1(1)
MC −1316.86(3) 0.07(6) 28.16(3) −1316.86(3) 0.07(6) 28.16(3) −1316.86(3) 0.07(6) 28.16(3) −1316.86(3) 0.07(6) 28.16(3) −1316.86(3) 0.07(6) 28.16(3)
Table B.9: Edwards-Anderson model. Part II.
28