New loop expansion for the Random Magnetic Field Ising Ferromagnets at zero temperature

TL;DR

This paper introduces a novel loop expansion method for the zero-temperature Random Field Ising Model, providing a more suitable framework for strongly disordered systems and revealing new cubic vertex contributions.

Contribution

The paper develops a new perturbative loop expansion around the Bethe solution tailored for T=0 disordered systems, differing from the epsilon-expansion and maintaining dimensional reduction validity.

Findings

New cubic vertices in the effective theory

One-loop corrections reveal additional terms

Dimensional reduction remains valid at this order

Abstract

We apply to the Random Field Ising Model at zero temperature (T= 0) the perturbative loop expansion around the Bethe solution. A comparison with the standard epsilon-expansion is made, highlighting the key differences that make the new expansion much more appropriate to correctly describe strongly disordered systems, especially those controlled by a T = 0 RG fixed point. This new loop expansion produces an effective theory with cubic vertices. We compute the one-loop corrections due to cubic vertices, finding new terms that are absent in the epsilon-expansion. However, these new terms are subdominant with respect to the standard, supersymmetric ones, therefore dimensional reduction is still valid at this order of the loop expansion.

Tables

Table 1. Table 1: Rules for evolving cavity fields and effective coupling in the computation of correlations at T = 0 𝑇 0 T=0 .

	$u_1^{\prime }$	$u_2^{\prime }$	$J_{12}$
$\left|h\right|>J_2+J_1$	$sign\left(h\right)J_1$	$sign\left(h\right)J_2$	0
	$h_{-}$	$h_{+}$	$\frac{J_1+J_2-\left|h\right|}{2}$
	0	$h$	$J_1$

Table 2. Table 2: Rules for computing the ground state configuration given the triplet of cavity messages ( h 1 , h 2 , J 12 ) subscript ℎ 1 subscript ℎ 2 subscript 𝐽 12 (h_{1},h_{2},J_{12}) .

	${\sigma }_1^{\ast }$, ${\sigma }_2^{\ast }$
	${\sigma }_1^{\ast }={\sigma }_2^{\ast }=sign\left(h_2\right)$
	${\sigma }_1^{\ast }={\sigma }_2^{\ast }=sign\left(h_1\right)$
	${\sigma }_1^{\ast }=sign\left(h_1\right)$, ${\sigma }_2^{\ast }=sign\left(h_2\right)$

Equations

$$\int d^{D}x\left(\frac{1}{2}(\partial \varphi (x){)}^2+\frac{1}{2}\tau \varphi (x{)}^2+\frac{1}{4}g\varphi (x{)}^4-h(x)\varphi (x)\right),$$

$$-\Delta \varphi (x)+\tau \varphi (x)+g\varphi (x{)}^3=h(x)$$

$$H[\sigma ,h]=-\frac{1}{2}\sum _{(i,j)}{\sigma }_i{\sigma }_j-\sum _i{h}_i^R{\sigma }_i,$$

$$C_{ij}\equiv \overline{{\sigma }_i^{\ast }{\sigma }_j^{\ast }}{R}_{ij}\equiv \frac{1}{2}\overline{1-{\sigma }_i^{\ast }⟨{\sigma }_i{⟩}_j}$$

$$G(x)=\sum _{L=1,\infty }\mathcal{N}(x,L)G^{\text{BL}}(L)$$

$$\mathcal{N}(x,L)\propto (2D-1{)}^L\exp \left(-x^2/(4L)\right)L^{-D/2},$$

$$\tilde{G}(p)\propto \sum _{L=1,\infty }(2D-1{)}^L\exp (-L{p}^2)G^{\text{BL}}(L).$$

$$\mathcal{H}[{\sigma }_1,{\sigma }_2]=-h_1{\sigma }_1-h_2{\sigma }_2-J{\sigma }_1{\sigma }_2+E.$$

$$P_L(h_1,h_2,J)=Q_L(h_1,h_2)\delta (J)+aL{\lambda }^L\hat{g}(h_1)\hat{g}(h_2)F_L(J)$$

$$P_L(J)=(1-aL{\lambda }^L)\delta (J)+a\rho L^2{\lambda }^L\exp (-\rho LJ)$$

$${\overline{J}}^L\propto {\lambda }^L$$

$$R^{\text{BL}}(L)\propto {\lambda }^L$$

$$P_L(h_1,h_2)=P(h_1)P(h_2)+L{\lambda }^{L}f(h_1)f(h_2),$$

$$C^{\text{BL}}(L)\propto L{\lambda }^L$$

$${\lambda }_c=\frac{1}{2D-1}.$$

$$\tilde{G}(p)\propto {\int }_0^{\infty }dL\exp \left(-L(p^2+\tau )\right){\mathcal{G}}^{\text{BL}}(L)$$

$$\tau \equiv -\log (\lambda (2D-1))$$

$$\begin{array}{rl}\tilde{C}(p)\propto & {\int }_0^{\infty }dL\exp \left(-L(p^2+\tau )\right)L=\frac{1}{(p^2+\tau {)}^2}\\ \tilde{R}(p)\propto & {\int }_0^{\infty }dL\exp \left(-L(p^2+\tau )\right)=\frac{1}{p^2+\tau }\end{array}$$

$$A(y,x;h)\equiv ⟨\sigma (y){⟩}_{h^{\ast }+h}-⟨\sigma (y){⟩}_{h^{\ast }-h}$$

$$\hat{R}(y-x)=\underset{h\to 0}{\lim }\overline{\frac{A(y,x;h)}{2h}}.$$

$$S(x;h)=\frac{1}{2}\int d^{D}yA(y,x;h).$$

$$\overline{S(x;h)}\approx h\chi \equiv h\int dy\hat{R}(y-x),$$

$${\tilde{G}}^{loop}(p)={\int }_0^{\infty }d\vec{L}\tilde{\mathcal{N}}(p,\vec{L}){\lambda }^{\Sigma (\vec{L})}{\mathcal{G}}^{\text{BL}}(\vec{L})$$

$$\tilde{\mathcal{N}}(p,\vec{L})$$

$$\mathcal{D}(\vec{L})$$

$$\tilde{\mathcal{N}}(p,\vec{L})$$

$$\mathcal{D}(\vec{L})$$

$$\begin{array}{rl}{\tilde{G}}_{{\varphi }^4}^{loop}(p)& =\frac{1}{(p^2+\tau {)}^2}\int d^{D}q\frac{1}{q^2+\tau }\\ {\tilde{G}}_{{\varphi }^3}^{loop}(p)& =\frac{1}{(p^2+\tau {)}^2}\int d^{D}q\frac{1}{q^2+\tau }\frac{1}{(p-q{)}^2+\tau }\end{array}$$

$$\frac{1}{p^2+\tau }={\int }_0^{\infty }dL\exp \left(-L(p^2+\tau )\right)$$

$${\mathcal{G}}_C^{\text{LG}}(\vec{L})=L_A(L_I+L_O){\mathcal{G}}_R^{\text{LG}}(\vec{L})=L_A$$

$$\frac{1}{(p^2+\tau {)}^2}={\int }_0^{\infty }dLL\exp \left(-L(p^2+\tau )\right),$$

Full text

\dates

This manuscript was compiled on

\verticaladjustment

-2pt

\dropcap

The critical behavior of the Random magnetic Field Ising ferromagnetic Model (RFIM) has been the subject of intense scrutiny in the last forty years. It is one of the most studied problems in statistical mechanics; however, at present, there are some basic questions to which we are unable to answer. The physics of the problem is quite clear: naively one expects that in presence of a random magnetic field a ferromagnetic transition is still possible but its critical properties should be different from that of an ordinary ferromagnet (1). These systems can be realized experimentally as diluted antiferromagnets in a field, assuming that these two classes of systems belong to the same universality class (2, 3, 4): a ferromagnetic transition is observed, however some of the measurements are quite difficult because of an exceptionally strong critical slowing down that forbids the system to thermalize well near and below the critical temperature. As far as the dynamical properties are concerned, the RFIM model shares many characteristics with other glassy systems (5).

The first study of the RFIM criticality used the perturbative renormalization group based on the same diagrammatic expansion that has been crucial to computing the critical exponents of standard ferromagnets in dimensions $D=4-\epsilon$ (6, 7). One way to treat the RFIM is to introduce an effective replicated $\phi_{4}$ model once the disorder has been integrated out. The diagrammatical rules are simple and we recall them in the Supporting Information (SI). Surprisingly, considering as usual only the most divergent diagrams and under some assumptions (see below) Dimensional Reduction (DR) holds, i.e. the critical properties (e.g. the critical exponents of this random system) are the same of a pure system in $D-2$ dimensions at each order of perturbation theory (8, 9, 10). DR is an astonishing relation as far as it connects the properties of different systems defined in different dimensions. In certain cases, DR has been rigorously established (11) 111Dimensional reduction has been proved between lattice animals in $D$ dimensions and ferromagnets in constant imaginary magnetic field in $D-2$ dimensions and between the Langevin stochastic equations in $D+1$ dimensions (one time dimension and $D$ space dimensions) and Boltzmann statistical mechanics in $D$ dimensions..

However, we know that DR cannot be always valid in the case of the RFIM. Indeed the lower critical dimension (i.e. the dimension where the ferromagnetic transition disappears) of the RFIM should be 2 plus the lower critical dimension of the standard Ising ferromagnetic model, which is 1. Therefore, according to DR, there should be no transition in the RFIM in $D=3$: this is at variance with a simple argument (12), (13), showing that the lower critical dimension for the RFIM is $D=2$. The existence of a transition in $D=3$ clearly shows that DR fails for the RFIM at low dimensions. The DR conjecture for the upper critical dimension instead (i.e. the dimension at and above which the critical exponents are those predicted by the mean field theory) is $D=6$ and is expected to be correct.

The origin of DR (and the possible cause of its failure) was identified when it was shown that the sum of the leading diagrams is related to the stationary points of the Landau-Ginzburg (LG) Hamiltonian

[TABLE]

where $h(x)$ is the random magnetic field. The stationary points of the LG effective Hamiltonian are the solutions of the equation

[TABLE]

In the case where the LG effective Hamiltonian has only one stationary point (i.e. only one local minimum) the sum of the leading diagrams gives the correct results: this is true in the region $\tau\geq 0$ due to the convexity of the LG effective Hamiltonian. When the free energy has multiple local minima, one should take the weighted sum of all the stationary points: in this case the results obtained just summing the leading diagrams are not correct. It can be shown that the LG effective Hamiltonian has many stationary points below the critical point ($\tau=0,g>0$) and the derivation of DR is not sound (14, 15). However, these multiple solutions are not seen in the standard perturbative (in $g$) expansion, and one can formally derive the DR result in any dimension.

In which dimension does DR break down? One can conceive at least three scenarios.

Dimensional reduction gives a grossly wrong result also near $D=6$. In this scenario, the difference between the correct results and those coming from the $6-\epsilon$ expansion is already non-zero at order $\epsilon^{k}$ for some $k$. In other words, the results of the $\epsilon$ expansion are not correct even as a Taylor expansion in $\epsilon$: this would be not so surprising because some assumptions in the derivation are wrong. 2. 2.

A less extreme suggestion is that this difference is exponentially small when $\epsilon\to 0$, e.g. it is of order $\exp(-A/\epsilon)$. This was argued in a non-fully conclusive way by Parisi, Dotsenko (16). The rationale for the suggestion is that since the difference is invisible in perturbation theory, it is reasonable to suppose that perturbation theory is still valid and that the difference is of order $\exp(-A/g_{R})$, where $g_{R}$ is the renormalized dimensionless coupling constant. We know that at the critical point $g_{R}$ is proportional to $\epsilon$ and this leads to the aforementioned conclusion. In this scenario the perturbation theory is correct but there are non-perturbative exponentially small terms that have to be added to get the correct result. 3. 3.

There is a critical dimension $D_{DR}$ such that for $D>D_{DR}$ dimensional reduction is correct and for $D<D_{DR}$ dimensional reduction fails. We can interpret this scenario by saying that there are two fixed points: the DR one and the one with broken DR. At $D=D_{DR}$ the DR fixed point becomes unstable and the stable one is the one with broken DR. This scenario has been strongly advocated by Tarjus and Tissier: using the Non-Perturbative Functional Renormalization group, they found that DR is no more valid for $D<D_{DR}\simeq 5.1$ (17)

Arguments in favor of each of three scenarios have been going on for a few decades. Some doubts on the correctness of LG based approaches come from the suggestion that the dominant fixed point of the RFIM is a $T=0$ one (1, 18). This means that the RG flux starting at or below the critical temperature will evolve towards zero temperature. The analytic arguments are complemented by very large scale numerical studies confirming that the critical exponents at $T=0$ and at finite $T$ are the same (19). Nowadays most of the more accurate simulations are done at $T=0$ in the Ising case where there is only a quite minor effect of critical slowing down near the critical point (20, 21). In the LG approach, and at odds with the BL one, we cannot set $T=0$ from the onset but have to take the zero-temperature limit.

In this paper, we argue that the first scenario is likely incorrect by showing that a different perturbative expansion leads to the same prediction as the $\epsilon$-expansion. As in functional renormalization group approaches, we don’t treat perturbatively the non-free field term of the Hamiltonian. The start of our computation is not the LG Hamiltonian but the Hamiltonian of the the original RFIM model. We study the RFIM at $T=0$ using a new loop-expansion around the mean field Bethe solution, recently proposed in ref. (22). This expansion has some very interesting features that make it very promising:

• •

The mean field Bethe solution is obtained solving the model on a random regular graph, usually called Bethe lattice (BL), having the same coordination number of the physical finite dimensional lattice we are interested into. On the contrary, the standard LG effective Hamiltonian can be derived introducing an auxiliary lattice whose coordination number goes to infinity (23). The finiteness of the coordination number introduces new features that produce crucial differences, as we will see later.

• •

The RFIM on a BL can be studied in great detail and many of its properties can be analytically computed. In standard mean field theory (that is valid on a fully connected lattice) the global magnetization is the only important variable and it satisfies a simple equation. On the contrary on the BL, the relevant quantity is the probability distribution of the local effective fields.

• •

Computations can be done setting the temperature straight to 0. The leading order of the perturbation theory (i.e. the theory on the BL) at $T=0$ can be easily derived. Successive orders can be can be expressed as diagrams corresponding to physical generalized loops on the lattice and their contribution for large distances can be identified and computed. Their physical interpretation is clear.

• •

Most importantly. we do not have to do assumptions on the uniqueness of the solution of some equations as it was implicitly done in the DR derivation.

In this paper, we shall see how this method is implemented. We stress that the case of the RFIM is just a particular application of the method that is much more general: a recent application shows that the mean-field hybrid transition of bootstrap percolation does not survive in finite dimension (24) and it could also be used for spin glasses (with zero or non-zero magnetic field) and for the Anderson transition.

Let us consider for definiteness an explicit realization of the RFIM. The spins are located on the $D$-dimensional lattice $\mathbb{Z}^{D}$ with nearest neighbors interactions. The Hamiltonian is

[TABLE]

where the sum is over nearest neighbors, $\sigma_{i}\in\{-1,+1\}$, and the random magnetic fields $h_{i}^{R}$ are independent and identically distributed according to a Gaussian with 0 mean and variance $\sigma^{2}_{h}$. At $T=0$, the only relevant configuration is the ground state of the system, that we call $\sigma^{*}$. In order to describe the critical behavior, we introduce the $T=0$ disconnected correlation $C_{ij}$ and response function $R_{ij}$:

[TABLE]

where the overline denotes the average over the random realizations of the fields $h^{R}$ and $\langle\cdot\rangle_{j}$ denotes evaluation on the minimum energy configuration conditioned to having spin $j$ flipped with respect to the ground state, i.e. $\sigma_{j}=-\sigma^{*}_{j}$. We note that $R_{ij}$ is the probability that site $i$ belongs to the lowest energy excitation of site $j$ (a droplet), therefore the sum over all sites $i$ yields the average droplet size, $\sum_{i}R_{ij}$. We considered the above definition of $R_{ij}$ because it is convenient for numerical experiments, while alternative choices are possible and will be discussed below. Both $C_{ij}$ and $R_{ij}$ are translational invariant, since they contain disorder averaged quantities.

Expanding around the Bethe solution

The expansion proposed in ref. (22) is an expansion around the Bethe solution: by replicating the model $M$ times, and rewiring the $M$ copies, one can show that the limit $M\to\infty$ gives the Bethe approximation, while the original model is recovered for $M=1$. A diagrammatic loop expansion can then be constructed expanding in powers of $1/M$, similarly to the standard perturbative expansion. We call this framework the $M$-layer or the BL approach. At leading order, the correlation functions on the $M$-layer are strictly related to the ones found on a BL. The latter are easy to compute since in the thermodynamic limit the BL contains no loops of finite length.

Using the $M$-layer expansion, one finds that a generic correlation $G$ (either $C$ or $R$) between the lattice origin and a point $x\in Z^{D}$ on the lattice can be written, at the leading order, as

[TABLE]

where ${\cal N}(x,L)$ is the number of non-backtracking paths going from the origin to $x$ of length $L$ on the original lattice and $G^{\text{\tiny{BL}}}(L)$ is the correlation function computed on the BL between two spins at distance $L$. The BL has the same connectivity $z=2D$ of the finite dimensional lattice and the same probability distribution for the random fields.

In the region of large $x$ and $L$ we find that in $D$ dimensions (25)

[TABLE]

where $x^{2}\equiv\lVert x\rVert^{2}$, and $\propto$ denote equality up to a constant. and we obtain for the Fourier transform of eq. [5] in the small momentum region

[TABLE]

So we just need to compute how correlations behave on a BL at large distances. As shown in the SI, at $T=0$ the crucial quantity that encodes all the information about two spins $\sigma_{1}$ and $\sigma_{2}$ at positions $x_{1}$ and $x_{2}$ is the dependence of the ground state energy on their values, after that we minimize over all the other spins. For Ising spins the energy must be of the form

[TABLE]

On the BL, we can derive a recursive equation for the joint probability distribution $P_{L}(h_{1},h_{2},J)$ for two spins at distance $L$, and obtain an explicit expression for large $L$ by imposing some consistency condition. The computation is presented in the SI, and from here on we discuss the large $L$ behaviour, which is the relevant one at criticality. The result can be written in the form:

[TABLE]

where $\lambda,\hat{g},F_{L}$ and $Q_{L}$ depend implicitly on the external random field distribution and $F_{L}(J)=\rho L\exp(-\rho LJ)$ at the leading order in $L$. In spin-glass jargon, $\lambda$ is called the anomalous eigenvalue and governs the decay of ferromagnetic susceptibilities along a chain in the BL (32, 30). The expression for the marginal probability $P_{L}(J)\equiv\int dh_{1}dh_{2}\,P_{L}(h_{1},h_{2},J)$ is thus given by:

[TABLE]

The coefficient $\rho$ can be computed exactly, as shown in the SI. The result [10] is quite surprising: the quantity $J$ is either exactly 0 or is of order $1/L$ with a probability of order $L\lambda^{L}$. Denoting with $\overline{\bullet}^{L}$ averages over $P_{L}$, we obtain

[TABLE]

It can be shown that the response function $R_{12}$ receives contributions only from the event $J>|h_{1}|$, and can be easily computed. We obtain that the average response function on a BL behaves as the average effective coupling:

[TABLE]

In a similar way, if we look to $P_{L}(h_{1},h_{2})=\int dJP_{L}(h_{1},h_{2},J)$ we find that for large $L$

[TABLE]

where the function $f(h)$ is even and $P(h)$ is the Bethe distribution of cavity fields. One can see that the dominant contribution to the disconnected correlation on the BL at distance $L$ comes from the correlated fields $h_{1}$ and $h_{2}$. For this reason we find:

[TABLE]

The behaviors for $R^{\text{\tiny{BL}}}$ and $C^{\text{\tiny{BL}}}$ correspond to the known ones on a line (29, 30).

We have found that for large $L$, $G^{\text{\tiny{BL}}}(L)\approx\mathcal{G}^{\text{\tiny{BL}}}(L)\lambda^{L}$, where $\cal{G}^{\text{\tiny{BL}}}$ is a polynomial in $L$. On the BL, the correlation functions decrease exponentially and the critical point is located where their exponential decrease matches the exponential increase of the number of paths (resulting in a diverging susceptibility). In both the BL and in finite dimensional model at the zeroth-order of the $M$-layer construction, the critical point is located at

[TABLE]

Near the critical point, starting from eq. [7] we can write

[TABLE]

where the sum over $L$ has been replaced by an integral and $\tau$ is the inverse of the correlation length $\xi$ and it is given by

[TABLE]

As usual, $\tau=0$ at the critical point. The representation in eq. [16] is the equivalent of the proper time representation in a field theory context.

If we now put eqs. [12,14] in eq. [16], we obtain the leading order of the expansion of the correlation functions of a $D$-dimensional model around the BL:

[TABLE]

Comparison with the Dimensional Reduction

The formulae in eq. [18] are the same of those coming from the LG effective Hamiltonian approach, with a few crucial differences though.

• •

In the LG approach, near the critical point and in the zero-loop approximation, the equations for the stationary points are linear: the unique solution is $\phi(x)\!=\!\int\!dyR(x-y)h(y)$, where $R(x)$ is the Fourier transform of $\widetilde{R}(p)$. In the BL approach, near the critical point and in the zero-loop approximation, there is an infinite number of local minima (i.e. configurations whose energy does not decrease if a finite number of spins are flipped), but the only thermodynamically relevant configuration is the global minimum (33).

• •

In the LG approach, the response function $R(x)$ does not depend on the field and it does not fluctuate: more precisely, all possible paths give the same contribution. In the case of the BL approach, only an exponentially small number of paths gives a contribution to the response function: for large $L$ the probability of a given path to have a non zero $J$ and consequently to contribute to the response is exponentially small (29).

The above differences become strikingly evident if we consider avalanches, a well-studied phenomenon in the RFIM. We are interested in seeing the change in the magnetizations when we change by a finite amount the magnetic field in a given point. Denoting the original field at position $x$ by $h^{*}$, we define

[TABLE]

where the label denotes the field at $x$. The quantity $A(y,x;h)$ is the variation of the magnetization at $y$ when we change the magnetic field at $x$ adding or subtracting a term $h$. In the limit of small field avalanches are related to an alternative choice of the response function that we call $\hat{R}(y-x)$. It amounts to consider only (the density of) excitations with strictly zero energy cost. More precisely, we have:

[TABLE]

One can show that $\hat{R}(y-x)$ is proportional to the zero temperature limit of the two-point connected correlation function, with a factor proportional to the inverse temperature.

The avalanche size $S(x;h)$ is given by

[TABLE]

For small $h$ we have

[TABLE]

with $\chi$ the susceptibility associated to the response function $\hat{R}(y-x)$. The susceptibility diverges as $1/\tau$ at the critical point both in the LG approach and on the Bethe lattice, but new features arise when we consider the probability distribution $P(S)$ of $S(x;h)$. Let us compare what happens in the two approaches at the zeroth-order of the loop expansion.

• •

In the LG approach, since $\hat{R}$ is related to connected correlation functions we easily find that $S$ does not fluctuates: $P(S)=\delta\left(S-\frac{c}{\tau}\right)\,$ for some constant $c$. Therefore, we find: $\overline{S}\propto\tau^{-1},\quad\overline{S^{2}}\propto\tau^{-2}$. The median value of $S$ is divergent.

• •

In the BL approach, following (26), one can argue that: $P(S)\propto S^{-3/2}\exp(-S\tau^{2}).$ We thus find: $\overline{S}\propto\tau^{-1},\quad\overline{S^{2}}\propto\tau^{-3}$. The median value of $S$ is finite and the divergence of $\overline{S}$ and $\overline{S^{2}}$ stems from rare events in the tail of the distribution.

The power law divergence of $\overline{S^{2}}$ is quite different in the two approaches: $\tau^{-2}$ in LG and $\tau^{-3}$ on the BL (more details in the SI).

1-st order in the Loop-expansion in the BL approach

We have seen that at the zeroth-order the physical behavior is quite different in the BL and in the LG approaches, although the critical behavior of the correlation function is superficially similar: in the LG approach, anomalous large fluctuations do not exist, while on the BL everything is dominated by rare large fluctuations. The superficial similarity for the average two-point correlations disappears if we look to high-order correlation functions (responsible for avalanches).

At this point, it is not clear what happens when we consider the loop expansion in the BL approach. The natural question is whether this loop expansion produces the same results as in LG. We have two alternative scenarios:

• •

The difference in the high-order correlations that we have seen at the tree level (zeroth-order BL) contaminates the two-points correlations when the leading contributions coming from the loop are considered. In this case, we would have additional terms at $T=0$ that are ignored in the LG approach. This would lead to the appearance of extra terms in the $\epsilon$ expansion in $6-\epsilon$ and DR should fail already in the $\epsilon$ expansion.

• •

The difference in the high-order correlations do not produce leading discrepancies on the two-points correlation functions and the contribution of the loops is the same as in the LG approach. As a consequence, in an unexpected way, we would recover perturbative dimensional reduction in $6-\epsilon$.

One can present many hand-waving arguments in favor of the first or the second scenario. However, the proof of the pudding is in the eating. In the following we prove that at one loop the results of the BL and LG approaches are the same. This will be done presenting a computation (down to the metal) of the one loop correction in the case of the BL.

Roughly speaking the idea at the basis of the loop expansion around the BL is to start approximating, at least locally, the $D$-dimensional lattice with loopless (acyclical) graphs: these are Caley trees with self-consistent conditions at the boundary or Bethe lattices. Of course, loops are present in the $D$-dimensional lattice and their effect is introduced perturbatively, by considering a sequence of BL with a finite number of loops. The expansion is similar to the virial expansion, where the complex interaction among infinitely many particles is decomposed in terms of simpler interactions between a finite number of particles.

The loop expansion around the BL is an expansion in topological diagrams. The contribution of a given topological diagram can be written as the probability of finding such a topological diagram embedded in the $D$-dimensional lattice times the averaged value that the observable takes on that given structure when inserted in a loop-less and infinite Bethe lattice. Only the topologically connected part of the average of the observable has to be consider. As it is shown in ref. (22), this connectivization procedure practically corresponds to adding the value of the observable evaluated on each of the subgraphs that are obtained from the original structure by sequentially removing its lines times a factor $-1$ for each line removed. The one loop contribution comes from the two diagrams shown in Fig. 1. They look similar to standard Feynman diagrams, however their physical interpretation is quite different. In standard Feynman diagrams the loops do not have a special meaning, here instead they have a geometrical meaning. Generalizing eq. [16], the one loop contribution, in Fourier space and as a function of the incoming moment, can be written as

[TABLE]

where $\vec{L}$ is a vector containing the lengths of each line in the topological diagram, the factor $\tilde{\cal{N}}(p,\vec{L})$ accounts for the number of such topological diagrams, while $\Sigma(\vec{L})$ is the sum of all $L$’s, and $\lambda$ is the same eigenvalue on the BL as in the previous discussion. The term ${\cal G}^{\text{\tiny{BL}}}(\vec{L})$ is the generalization of ${\cal G}^{\text{\tiny{BL}}}(L)$ at the zero-th order: it is the only term depending on the model and has to be carefully computed on the BL. In the case of the two diagrams in Fig. 1, we find:

• •

for the left diagram

[TABLE]

• •

for the right diagram

[TABLE]

Setting $\mathcal{G}(\vec{L})=1$ we recover the conventional diagrams of the field theory approach in the cases of a $\phi^{4}$ interaction (left diagram) or $\phi^{3}$ interaction (right diagram) written in the Feynman proper time representation:

[TABLE]

where $\tau=-\log(\lambda(2D-1))$ as usual. In fact, we can go backward from last expression containing integrals in momentum space to the previous one: for each line we have to use the representation

[TABLE]

In this way, the integral over the loop momentum $q$ becomes Gaussian: it can be readily done and we obtain the previous results, eqs. [24,25].

It is clear that the crucial point is the computation of the function ${\cal G}^{\text{\tiny{BL}}}(\vec{L})$, since it contains all the information related to the theory we are considering. In the case of the standard LG approach, only the left diagram is present: a standard computation gives for the disconnected and the connected correlation functions

[TABLE]

and, using the representation

[TABLE]

we obtain the standard result where some lines have a single pole, $(p^{2}+\tau)^{-1}$, while others have a double pole $(p^{2}+\tau)^{-2}$. This representation can be derived also for higher orders of the perturbative expansion. The first perturbative proof (9) of DR was based on the use of the identity for the diagrams: ${\cal G}_{R}^{\text{\tiny{LG}}}(\vec{L})={\cal D}(\vec{L}).$ In this way, the denominator in $\tilde{\cal N}$ becomes ${\cal D}({\vec{L}})^{D/2-1}$ and the final expression for the diagrams is the same of a vanilla $\phi^{4}$ theory in dimensions $D-2$.

How to compute the factors ${\cal G}^{\text{\tiny{BL}}}(\vec{L})$ in the BL approach? We have to compute the average connected and disconnected correlations on a BL where we have the same local geometry ($z=2D$) plus a manually injected topological diagram. The final results can be summarized as follows

• •

The left diagram gives the same type of contribution of the diagrams of LG reproducing DR.

• •

In the region where either $L_{A}$ or $L_{B}$ is small the right diagram has a behavior quite similar to the left diagram. Nothing new comes from this diagram in this region.

• •

The real interesting region for the right diagram is when all the $L$’s are large: this gives the relevant contribution at large distances (small momentum) discussed below.

Computation of the new diagram on the BL

In order to compute factors ${\cal G}^{\text{\tiny{BL}}}(\vec{L})$ on a BL, we have to go through a sequence of simple steps. Some of them are rather lengthy yet straightforward.

We consider a BL where we add a loop of the type of (Fig. 1, right). Apart from the loop, the rest of the lattice is a standard BL with fixed connectivity $z=2D$: this means that variables $\sigma_{I}$ and $\sigma_{O}$ are the root of $z-1$ infinite tree-like branches, variables $\tau_{I}$ and $\tau_{O}$ are the root of $z-3$ tree-like branches, while the spins along the topological lines are the root of $z-2$ tree-like branches.

We are interested in computing the probability distribution of the random restricted Hamiltonian $\mathcal{H}[\sigma_{I},\sigma_{O}]$, i.e. the one we obtain after minimizing with respect to all the other variables. This 2-spins Hamiltonian is obtained from the 4-spins Hamiltonian $\mathcal{H}[\sigma_{I},\sigma_{O},\tau_{I},\tau_{O}]$, where the distances among are fixed to the values $L_{I},L_{O},L_{A},L_{B}$ shown in Fig. 1, by

[TABLE]

The 4-spins Hamiltonian can be computed by summing four statistically independent 2-spins Hamiltonian and the cavity fields on $\sigma_{I},\sigma_{O},\tau_{I},\tau_{O}$ coming from the infinite trees. A two-spin Hamiltonian for a line of length $L$ is described by two fields and a coupling, $(u_{1},u_{2},J)$, whose joint law for large $L$ can be written in the form

[TABLE]

Last equations differs from eq. [9] only in the fact that here we do not include the contribution by the external fields and the cavity fields for the spin at the extremities of the line. When we compute the probability distribution of the quantity $\mathcal{H}[\sigma_{I},\sigma_{O},\tau_{I},\tau_{O}]$ it factorizes into the product of four terms coming from each of the lines. Connectivization of the diagram, as prescribed by ref. (22), corresponds to removing the term $P(u_{1})P(u_{2})\delta(J)$ on each line. Therefore, on each line we can decide if we take the contribution $Q^{D}$ or $Q^{C}$: in the first case we have a disconnected term that can be represented with a line bearing a cross, in the second case we have a connected term that can be represented with a line without a cross. In this way the diagrammatics becomes graphically equivalent to the one of the LG approach, with the addiction of extra diagmas containing cubic vertices. We note that in this computation it is not obvious that the most divergent diagrams will be the ones with the maximal number of possible crosses, as in the standard LG expansion that leads to DR. In fact, this is not the case for diagrams with cubic vertices. We obtain the following results:

• •

the connected correlation (response) function is not renormalized at one loop, since ${\cal G}^{\text{\tiny{BL}}}_{R}(\vec{L})=0$;

• •

a factor ${\cal G}^{\text{\tiny{BL}}}_{C}(\vec{L})=\text{const}\neq 0$ appears for the disconnected correlation function when $L_{I}=L_{O}=L_{A}=L_{B}=L$. Numerically, when the four lengths are all different, the result is consistent with the behaviour ${\cal G}^{\text{\tiny{BL}}}_{C}(\vec{L})=a+b\left(\frac{L_{A}}{L_{B}}+\frac{L_{B}}{L_{A}}\right)$ .

The detailed derivation of this result is presented in the SI, together with a numerical consistency check.

At this point one should compare this new contribution to the one obtained from the diagrams coming from the standard expansion around the LG theory, looking at the power-law divergence when $\tau\to 0$ in the limit $p\to 0$.

Let’s focus on $C$. Within the LG approach, the divergence of the diagram is of order $\tau^{-5}$. Noticing that

[TABLE]

and that the Incomplete Euler Gamma function behaves as $\Gamma\left[0,(p^{2}+\tau)\right]\simeq-\log[p^{2}+\tau]$ for ${(p^{2}+\tau)\to 0}$, and using eqs. [27,29], we find that the divergence of the new diagram is at most $-\tau^{-4}\log(\tau)<\tau^{-5}$.

The new diagrams coming from cubic vertices are thus sub-dominant with respect to the standard ones in the one loop expansion of two point correlation functions.

Discussion, Conclusions and Perspectives

In this work, we have applied the new topological expansion around the Bethe solution proposed in ref. (22) to the RFIM at $T=0$, numerically and semi-analytically, obtaining consistent results. It is crucial that we expand around the Bethe solution because it is deeply different from the one found in a standard Landau-Ginzburg approach: while in the latter fluctuations do not play any role, the Bethe solution is dominated by rare fluctuations, especially at $T=0$; this is of primary importance given that the RFIM critical behavior is controlled by a $T=0$ fixed point.

A direct consequence of the fluctuations-dominated behavior at $T=0$ is that higher order correlations do not decay faster than the average correlation, $\overline{G(x)^{p}}\approx\overline{G(x)}$, and this produces an effective theory with vertices of all degrees, including cubic vertices, essentially because diagrams with multiple lines between the same vertices are allowed.

We have analyzed the first two one-loop corrections to the correlation functions due to cubic vertices, finding that they give a contribution that is divergent at $D<6$, as also happens for the standard quartic diagrams.

We also found that they give an extra contribution. However, this contribution is sub-dominant with respect to the one given by the usual one-loop diagram coming from the standard LG theory. This means that, within our framework and at the 1-loop order, Dimensional Reduction is still valid at $6-\epsilon$ dimensions because the most divergent diagrams remain the super-symmetric ones

Let us finally remark that the analyzed cubic vertices are really important already at the mean-field level (zero-th order). A peculiarity of the RFIM at $T=0$ on finite connectivity lattices is the existence of avalanches: this collective $T=0$ phenomenon cannot be described within the standard field theoretical treatment, while it appears naturally if $\phi^{3}$ vertices are introduced. At the critical point, the avalanches size distribution follows a power law with a nontrivial exponent $\tau$. In our framework, we easily enough find the correct mean-field value $\tau=3/2$ that cannot be computed within the standard LG approach. Avalanches have a fractal dimension $d_{f}$ that is connected to fluctuations in the $T=0$ integrated response $\chi$ via $\overline{\chi^{2}}/\overline{\chi}\propto L^{d_{f}}$. We plan to compute the one-loop correction for $\chi^{2}$, i.e. for three-point functions, and obtain in this way the $\epsilon$-expansion for $d_{f}$. \acknowThis research has been supported by the European Research Council under the European Unions Horizon 2020 research and innovation programme (grant No. 694925 – Lotglassy, G Parisi) and by the Simons Foundation (grant No. 454949, G Parisi). \showacknow

1 Standard diagrammatic rules for the random field Ising model

The standard way to compute the loop expansion for the Random Field Ising Model (RFIM) is to introduce an effective replicated $\phi_{4}$ model once the disorder has been integrated out (31). In practice one is left with few operating rules to construct Feynmann diagrams, that we briefly recall here.

The main difference with a standard $\phi_{4}$ theory is that the bare propagator is composed of two parts: a connected part, that is commonly indicated with a line, going as $\widetilde{R}(p)\propto\frac{1}{p^{2}+\tau}$, that will contribute to the connected correlation function $\overline{\langle\sigma_{i}\sigma_{j}\rangle_{c}}=\overline{\langle\sigma_{i}\sigma_{j}\rangle-\langle\sigma_{i}\rangle\langle\sigma_{j}\rangle}$ and a disconnected part, indicated with a line plus a cross, $\widetilde{C}(p)\propto\frac{1}{(p^{2}+\tau)^{2}}$, that will be the dominant contribution to the disconnected correlation function $\overline{\langle\sigma_{i}\rangle\langle\sigma_{j}\rangle}$.

In practice, to build Feynmann diagrams, one should put vertices with 4 lines, that could be connected or disconnected. The only rule in the construction of the diagrams for the expansion of the two-point connected correlation function is that at least one connected path between the two points should be present. Selecting only the most divergent diagrams at each order, one discovers that they correspond to the diagrams with the highest number of allowed crosses at each order in the perturbative expansion. They are shown in Fig. 2 for the connected correlation function up to second order.

2 Solution of the Random Field Ising Model on the Bethe Lattice at $T=0$

The solution presented in this Section is the $0^{th}$ order of the loop expansion presented in the main text. It has been presented in some detail in Refs. (29, 30), but we find useful to reported here again for completeness.

We consider a model with Hamiltonian

[TABLE]

where $J>0$ and $h_{i}^{R}$ are i.i.d. random variables extracted from a Gaussian probability distribution with zero mean and standard deviation $\sigma_{h}$. The edge set $E$ defines a Bethe lattice (BL) of finite connectivity $z$ (mathematically speaking it is a random regular graph of constant degree $z$).

Following the standard cavity method, we consider cavity fields $h_{i\to j}$ and $u_{i\to j}$ defined on each edge of the graph. They parametrize, respectively, the marginal probability distribution on $\sigma_{i}$ in the cavity graph where edge $(ij)$ has been removed, and the marginal probability distribution on $\sigma_{j}$ in the cavity graph where all edges involving vertex $j$, but $(ij)$, have been removed. At $T=0$ the self-consistency equations among cavity fields read

[TABLE]

where $\partial i$ is the set of neighbours of $i$, i.e. spins linked to $i$ via an edge of the graph.

Within the cavity method one is interested, rather than in the specific solution on a given graph, in the solution averaged over the ensemble of random graphs and random fields. To this end it is enough to solve Eqs. [33,34] in distribution sense and compute the probability distributions of cavity fields $h$ and $u$. We call $P(u)$ the latter.

Willing to compute the correlations between spins $\sigma_{0}$ and $\sigma_{L}$ that are connected by a line of length $L$ (the path linking $\sigma_{0}$ and $\sigma_{L}$ is unique on a BL in the thermodynamic limit) we need to integrate out all the spins along the line and compute the triplet $(u_{0},u_{L},J_{L})$, where $J_{L}$ is the effective coupling between $\sigma_{0}$ and $\sigma_{L}$, while $u_{0}$ and $u_{L}$ are the effective fields on $\sigma_{0}$ and $\sigma_{L}$ coming from the line. Such a triplet can be computed in a recursive way. Let us join two chains, the first one between $\sigma_{1}$ and $\tau$, characterized by the triplet $(u_{1},u_{\tau,1},J_{1})$, and the second one between $\tau$ and $\sigma_{2}$ identified by $(u_{\tau,2},u_{2},J_{2})$. In order to compute the triplet describing the effective Hamiltonian between $\sigma_{1}$ and $\sigma_{2}$ we need to sum over $\tau$ and keep only the lowest energy term (we are working at $T=0$)

[TABLE]

with $h=h_{\tau}^{R}+u_{\tau,1}+u_{\tau,2}+\sum_{k\in\partial\tau\setminus 1,2}u_{k\to\tau}$, and $u_{k\to\tau}$ are independent random variables extracted from $P(u)$. Explicit expressions for $u_{1}^{\prime}$, $u_{2}^{\prime}$ and $J_{12}$, assuming $J_{2}\geq J_{1}\geq 0$, are given in Table 1, where h_{\pm}=\big{(}h\pm\operatorname*{sign}(h)(J_{2}-J_{1})\big{)}/2. We went from the two initial triplets $(u_{1},u_{\tau,1},J_{1})$, and $(u_{\tau,2},u_{2},J_{2})$ to the new one $(u_{1}+u_{1}^{\prime},u_{2}+u_{2}^{\prime},J_{12})$, with the insertion of $z-2$ cavity fields acting on the central spin $\tau$.

In practice we start from a population $P_{L=1}(u_{0},u_{1},J_{1})$ of triplets all equal to $(0,0,J)$. To evolve the population $P_{L-1}$ into population $P_{L}$ we follow the rules summarized in Table 1, where each triplet of the population $P_{L-1}$ is joined to a triplet $(0,0,J)$ and $z-2$ cavity fields $u_{k\to\tau}$ extracted from $P(u)$ are added on the central spin.

Unfortunately this procedure is very ineffective, because at each step a constant fraction of the population (the one satisfying the condition $|h|>J_{L-1}+J$) produces a new triplet with $J_{L}=0$. Given that $J_{L}=0$ is a fixed point of the iteration, the part of the population keeping information about branches with non-zero effective couplings shrinks exponentially fast during the iteration.

To amplify this signal, we evolve two populations of the same size: one population keeps the pairs $(u_{0},u_{L})$ along branches with $J_{L}=0$, while the second one stores the triplets along branches with $J_{L}\neq 0$. At the same time we measure the probability $p_{L}=\mathbb{P}[J_{L}\neq 0]$, that is the relative weight of the second population to the first one, which is found to decay exponentially fast with $L$: $p_{L}=aL\lambda^{L}+b\lambda^{L}+o(\lambda^{L})$ as shown in Fig. 3. $\lambda$ is the largest eigenvalue associated to the linearization of the BP eqs. [33,34] around the fixed point. At the critical point, $\sigma_{h}=\sigma_{h,c}$, $\lambda(\sigma_{h,c})=1/(z-1)$ holds. The average coupling on the second population decays as $\overline{J_{L}}=\frac{a}{L}+\frac{b}{L^{2}}+o(L^{-2})$ as shown in Fig. 4. We see that $\overline{u_{0}u_{L}}\propto\frac{1}{L^{2}}$ on the population with $J\neq 0$ while $\overline{u_{0}u_{L}}\propto L\lambda^{L}$ on the population with $J=0$. Moreover $\overline{u_{0}^{2}J}-\overline{u_{0}^{2}}\,\overline{J}\propto\frac{1}{L^{2}}$ on the population with $J\neq 0$.

Once we have $p_{L}$ and the two populations at each length $L$, it is quite simple to compute correlation functions. Indeed, given a triplet $(u_{1},u_{2},J_{12})$ associated to a path, where the internal spins have been integrated out, the effective two spin Hamiltonian reads

[TABLE]

with $h_{1}=h^{R}_{1}+u_{1}+\sum_{k1=1}^{z-1}u_{k1}$, $h_{2}=h^{R}_{2}+u_{2}+\sum_{k2=1}^{z-1}u_{k2}$, and $u_{k1}$, $u_{k2}$ extracted from $P(u)$. At zero temperature the Gibbs measure is concentrated on the ground state $(\sigma^{*}_{1},\sigma^{*}_{2})$ of the Hamiltonian, that can be easily computed using the rules listed in Table 2.

Since we are at $T=0$, the disconnected correlation function is given by

[TABLE]

where $\sigma^{*}$ is the ground state configuration, computed according to the rules listed in Table 2. The connected correlation function $C^{con}_{ij}=\overline{\langle\sigma_{i}\sigma_{j}\rangle_{c}}$ is ill defined since it is identically equal to zero at $T=0$, therefore we work with the response

[TABLE]

where $\langle\cdot\rangle_{j}$ denotes the expectation over the ground state of the system conditioned to the flipping of the spin $\sigma_{j}$ , i.e. $\langle\sigma_{j}\rangle_{j}=-\sigma^{*}_{j}$. This can be achieved adding a field $h_{j}=-\sigma^{*}_{j}\cdot\infty$ on the spin $\sigma_{j}$. An alternative and more general definition of response would be $R_{ij}=\frac{1}{2}\ \overline{\sigma^{*}_{j}(\sigma^{*}_{i}-\,\langle\sigma_{i}\rangle_{j})}$. Since in the RFIM the couplings are ferromagnetic, the spin $\sigma_{i}$ can only flip in the same direction of the flip of $\sigma_{j}$, therefore the two definitions are equivalent and $0\leq R_{ij}\leq 1$. It can be shown that the response can be expressed as the zero temperature limit of an opportunely normalized connected correlation function, that is

[TABLE]

In terms of the probability law of random triplets $(h_{1},h_{2},J_{12})$ the correlation functions can be written as

[TABLE]

On chains, $h_{1}$ and $h_{2}$ are positively correlated,therefore $C_{L}\geq R_{L}$. Notice that only events with a non-zero effective couplings contribute to the response function: this is the reason why amplifying the population of cavity messages with $J_{L}\neq 0$ is mandatory to have a precise measurement of correlations in the $T=0$ limit.

In Fig. 5 we show the connected and disconnected correlation functions at distance $L$, averaged over the population of the triplets generated as explained before, in a BL with fixed connectivity $z=3$, at zero temperature and critical standard deviation $\sigma_{h,c}=1.037$ for the external field. We find the the connected correlation function decays as $R_{L}\propto\lambda^{L}$, with $\lambda=\frac{1}{z-1}$, while the disconnected correlation function is larger and decays as $C_{L}=(aL+b)\lambda^{L}$, as already found analytically in Refs. (29, 30).

The corresponding susceptibilities can be computed by summing over all the vertices of the graph

[TABLE]

where $n_{L}$ is the number of spins at distance $L$ that in a BL is $n_{L}=\frac{z}{z-1}\cdot(z-1)^{L}$. Substituting $n_{L}$ and $C_{L}$ in the equation for $\chi$ one gets:

[TABLE]

At the critical point $\lambda(\sigma_{h,c})=\frac{1}{z-1}$ and the susceptibility diverges. How can we relate this computation of the susceptibility on the Bethe lattice to the perturbative expansion for a finite dimensional model in dimension $D$ around the Bethe theory? Following ref. (22), the zeroth order expansion for the susceptibility is just eq. [42] where $n_{L}$ is replaced with the number of non-backtracking paths of length $L$ starting from a point in a $D$ dimensional lattice: $n_{L}\propto(2D-1)^{L}$. Therefore at zeroth order the expansion predicts a divergence located at the same critical point $\sigma_{h,c}$ of a Bethe lattice with connectivity $z=2D$.

2.1 One loop BL correction

The first order in the BL expansion considers the presence of structures with one spatial loop, as reported in Fig. 6. In this section, we compute the one-loop correction to the connected and disconnected correlation function, coming from the topological structure in the right part of Fig. 6, that is the one that gives an additional term with respect to the standard Landau-Ginsburg (LG) expansion. Following the prescription of ref. (22), we should compute the correlation on a BL in which such structure has been manually injected and subtracting the values of the correlation computed on the two paths $L_{I}+L_{A}+L_{O}$, $L_{I}+L_{B}+L_{O}$, supposed as independent.

Operatively we build the loop putting together four triplets or couples extracted independently from the populations of single Bethe lines obtained as explained in the previous section: two of length $L_{I}$ and $L_{O}$ for the external lines, and two of length $L_{A}$ and $L_{B}$ for the internal lines of the loop. The internal lines of the loop, with triplets $(u_{\tau_{I}}^{A},u_{\tau_{O}}^{A},J_{A})$ and $(u_{\tau_{I}}^{B},u_{\tau_{O}}^{B},J_{B})$, will just result in a new triplet whose coupling is the sum of the couplings: $J_{T}=J_{A}+J_{B}$ and whose fields are the sum of the fields $u_{\tau_{I}}^{T}=u_{\tau_{I}}^{A}+u_{\tau_{I}}^{B}$. Then the new triplet is attached to the external legs, as in eq. (35), with the only difference being that in $\tau_{I}$ and $\tau_{O}$ there are just $z-3$ additionally cavity fields (instead of $z-2$ ones) extracted from $P(u)$. We then end up with a new triplet describing the loop. We compute the correlations implied by this triplet and subtract the correlation implied by the paths $L_{I}+L_{A}+L_{O}$ and $L_{I}+L_{B}+L_{O}$ considered as independent: This is the one-loop contribution to the correlation function, that we will indicate with $G_{loop}(\vec{L})=G_{loop}({L_{I},L_{O},L_{A},L_{B}})$. In the following we report the one-loop results for $G=C,R$, that we obtain looking numerically to the behaviors when $L_{A}$ or $L_{B}$ are small or large.

Let us first analyze the response. We know that the contribution to $R$ given by the external legs should, in any case, be proportional to $\lambda^{(L_{I}+L_{O})}$, because the diagram should be connected to contribute to the connected correlation function. Thus we concentrate on the internal legs. First of all, we fix also $L_{B}$ to a finite value $L_{B}=3$ and we measure the contribution of the loop to the response function as a function of $L_{A}=L$. We measure the behavior

[TABLE]

as shown in Fig. 7. As expected, the behaviour is the same of the one-loop diagram coming from the standard theory, that is the left diagram of Fig. 6. In fact, two $\phi_{3}$ vertices reduces to a tadpole $\phi_{4}$ vertex once one internal line is fixed to a finite length. For a tadpole $\phi_{4}$ diagram, LG theory predicts that the maximal divergent contribution comes from the second diagram of Fig. 2, that has indeed the same behavior of as in eq. [43].

Next, we fix $L_{A}=L_{B}$ and we measure $R_{loop}(L_{I},L_{O},L,L)$. It receives contributions from two different diagrams: we call “Contribution A” the one coming from the loop with both $J_{A}$ and $J_{B}$ different from zero, that is, in the language of eq. [36] of the main text, taking the contribution coming from $Q^{C}_{L}(u_{1},u_{2},J)$ on both lines; “Contribution B” instead, is the one from a loop with just one coupling different from zero, that is, taking a contribution coming from $Q^{C}_{L}(u_{1},u_{2},J)$ on one line, and a contribution coming from $Q^{D}_{L}(u_{1},u_{2})\delta(J)$ on the other line. (Please note that it is not possible to take contribution from $Q^{D}_{L}(u_{1},u_{2})\delta(J)$ on both internal lines because it will result in a disconnected loop that gives zero contribution to connected correlation functions). Separately, contribution A and B, multiplied by their occurrence probabilities, have a dominant behaviour in $L$ of the type $(a+bL)\lambda^{L_{I}+L_{O}+2L}$, but they have opposite sign. When summing the two contributions, the dominant term in $L$ is exactly cancelled, and the total contribution is left with the subdominant part

[TABLE]

as shown in fig. 8. As described in the main text, we can write $R_{loop}(\vec{L})=\lambda^{\Sigma(\vec{L})}{\cal G}_{R}(\vec{L})+o(\lambda^{\Sigma(\vec{L})})$, with $\Sigma(\vec{L})=L_{I}+L_{O}+L_{A}+L_{B}$. Eq. [44] corresponds to ${\cal G}_{R}(\vec{L})=0$.

Things are different for the disconnected correlation function. Also in this case, first of all we fix $L_{B}$ to a finite value and we measure the behaviour of the loop as a function of $L_{A}=L$. We measure the behavior

[TABLE]

that again is the same contribution of the tadpole diagram from the usual $\phi_{4}$ LG theory. Then we put $L_{A}=L_{B}=L_{I}=L_{O}=L$, and we observe

[TABLE]

that corresponds to ${\cal G}_{C}(L,L,L,L)=const$.

We now want to compute the one-loop correction to the susceptibilities. In an analogous way to eq. [42], we should account for all the subgraphs of the type of the right dyagram in Fig. 6 that are presents in a finite-dimensional lattice. The computation can be done exactly using the number of non-backtracking paths, however, the large $L$ behavior of this counting factor is also captured if we assume that the number of paths from 0 to $x$ of length $L$ is given by the random walk probability of reaching $x$ in time $L$ in $D$ dimensions multiplied by the number of generic non-backtracking paths of length $L$ starting from 0: $n_{L}(x)\propto\frac{(z-1)^{L}}{L^{D}}e^{-x^{2}/(2L)}$. In the same way we can compute $n_{L_{A},L_{B}}(x,y)$, defined as the number of paths of length $L_{A}$ and $L_{B}$ that have the same starting and ending point $x$ and $y$: $n_{L}(x,y)\propto\frac{(z-1)^{L_{A}+L_{B}}}{(L_{A}L_{B})^{D}}\left(e^{-\frac{(x-y)^{2}}{2}(\frac{1}{L_{A}}+\frac{1}{L_{B}})}\right)$. The one loop correction to the susceptibility associated to a generic correlation function $G$ is thus:

[TABLE]

Replacing the sums over $x,y,d$ with integrals, and performing the integrals we obtain:

[TABLE]

It is now clear that the contribution to the connected susceptibility is always finite at the critical point, because $R_{loop}(\vec{L})$ decays more rapidly than $(z-1)^{-\Sigma(\vec{L})}$. Things are different for the disconnected susceptibility. In this case, we can apply the Ginzburg criterion to identify the upper critical dimension: we look at $\frac{d\chi}{d\lambda}$, that is divergent at the critical point in $D\leq 6$. At this point, one could think that Dimensional Reduction is broken at $D=6$. In fact for $D<6$, we have cubic diagrams, not of the type of the super-symmetric ones, that are important. However, as we explained in the main text, once we compare their divergence with the divergence of the standard one-loop $\phi_{4}$ diagrams, we discover that their contribution is sub-dominant with respect to the usual $\phi_{4}$ term. This implies that DR is still valid at $6-\epsilon$ dimensions.

To conclude, we just mention that until now we do not know the exact behavior of $C_{loop}$ as a function of the length of the legs. Having measured $C_{loop}(L,L,L,L)\propto\lambda^{4L}$ we can think to different cases:

• •

A) $C_{loop}(\vec{L})\propto\lambda^{L_{I}+L_{A}+L_{B}+L_{O}}$

• •

B) $C_{loop}(\vec{L})\propto\lambda^{L_{I}+L_{A}+L_{B}+L_{O}}(\frac{L_{A}}{L_{B}}+\frac{L_{B}}{L_{A}})$

• •

C) $C_{loop}(\vec{L})\propto\lambda^{L_{I}+L_{A}+L_{B}+L_{O}}(\frac{L_{I}}{L_{B}}+\frac{L_{I}}{L_{A}}+\frac{L_{O}}{L_{B}}+\frac{L_{O}}{L_{A}})$

We expect that the presence of terms of the type $\lambda^{L}/L$ should signal the presence of squared disconnected correlation function 222We somehow expect the presence of important square correlations at zero temperature, see the Conclusions. In fact we measured numerically that $\overline{(u_{0}u_{L})^{2}}\propto\lambda^{L}/L$, and we expect that $\overline{(\langle\sigma_{0}\rangle\langle\sigma_{L}\rangle)^{2}}\propto\overline{(u_{0}u_{L})^{2}}$, given that $\overline{\langle\sigma_{0}\rangle\langle\sigma_{L}\rangle}\propto\overline{u_{0}u_{L}}$. To understand which terms are present, we measure $C_{0}\equiv C_{loop}(L_{I},L_{O},L,L)$ and $C_{1}\equiv C_{loop}(L_{I},L_{O},2L,L)$, at fixed, finite values of $L_{I}$ and $L_{O}$, and we look at them as a function of $L$. We obtain the behaviours: $C_{0}\propto\lambda^{2L}$, $C_{1}\propto\lambda^{3L}$. The ratio $Q(L_{I},L_{O})\equiv\frac{C_{0}/\lambda^{2L}}{C_{1}/\lambda^{3L}}$ is independent from $L_{I}$ and $L_{O}$. This result tells us that the case C) is not present. Indeed this is what we expected: the case C) corresponds to a connected loop, but we know that the connected correlation, that can receive contribution only by a connected loop, is not renormalized at one loop. We thus expect that the connected loop gives no contribution to $C$, as found. We numerically find that $Q=0.96$: if the situation A) were the only present, $Q=1$, while in the case B) $Q=0.8$: to recover the measured $Q=0.96$ we need a linear combination of the two cases. From the numerical computation, we thus expect the one-loop contribution to the disconnected correlation function to have the form

[TABLE]

with $a=1$, $b=0.1$.

3 BL results for the distribution of Avalanches

The distribution of the size of the avalanches $s$ at the critical point is expected to be

[TABLE]

with $\rho$ the critical exponent for the avalanches whose value on the BL is $\rho_{MF}=\frac{3}{2}$. This distribution can be obtained in the framework of percolation on the BL (see (26) and refs. therein).

We explained in the main text that $s$ is proportional to the susceptibility associated to the connected correlation function: $\chi=\sum_{x}\langle\sigma_{0}\sigma_{x}\rangle_{c}$. Given that $\overline{\chi}\propto\tau^{-1}$, the distribution [47] in the MF region implies that $\overline{\chi^{2}}\propto\tau^{-3}\neq\overline{\chi}^{2}$. This result cannot be recovered from the LG theory. In this case, in fact, the global magnetization is the only important variable, there are no fluctuations in the magnetization nor in the susceptibility, for which therefore we can write $\overline{\chi^{2}}=\overline{\chi}^{2}\propto\tau^{-2}$. Let us now look in detail to what are the field-theoretical predictions on $\chi^{2}=\sum_{x,y}\langle\sigma_{0}\sigma_{x}\rangle_{c}\langle\sigma_{0}\sigma_{y}\rangle_{c}$, for which we have to look to three point functions. If we admit that there are only $\phi_{4}$ vertices, as in the MF FC model, the diagram with no loop is the left one in Fig. 10. Giving that each line corresponds to a connected propagator and thus bring a factor $\tau^{-1}$, the left diagram will be associated at the critical point to a divergence of the type $\overline{\chi^{2}}\propto\tau^{-2}$, recovering the FC MF result. If now we imagine that the associated field theory includes also $\phi_{3}$ vertices, the situation will change: the right diagram in Fig. 10 is possible, leading to a critical behaviour: $\overline{\chi^{2}}\propto\tau^{-3}$.

We have announced that, following ref. (22), diagrams with $\phi_{3}$ vertices should be present in the field-theoretical description of the RFIM at $T=0$ when expanding around the finite connectivity Bethe solution: in this section, we have shown that their presence is perfectly compatible with the MF description of the avalanches, for which we can recover the critical exponent $\rho^{MF}=\frac{3}{2}$, in contrast to the standard FC $\phi_{4}$ theory that cannot justify the probability distribution of the avalanches.

In ref. (27), the following connection between avalanches and DR is stated: DR breaks down due to avalanches if they are “big enough”, more precisely if the fractal dimension $d_{f}$ of the largest typical critical avalanches satisfies the condition $d_{f}=D-d_{\phi}$, with $D$ the spatial dimension and $d_{\phi}$ the scaling dimension of the field near the relevant zero-temperature fixed point. In a way analogous to Ref. (27), we have seen that the $\phi_{3}$ diagrams will not automatically destroy DR: in particular, at one loop they are sub-dominant with respect to the standard $\phi_{4}$ ones implying that at $D=6-\epsilon$, DR is preserved.

4 Ansatz for coupling and fields at distance L

In this section, we verify our previous numerical results using a different method. We introduce an Ansatz $P_{L}(u_{0},u_{L},J)$ for the joint distribution of the effective coupling and fields between two spins at distance $L$ in a BL, which should capture the leading behviour at large $L$. We assume the form

[TABLE]

and check it’s consistency. $P(u)$ is the already mentioned Bethe distribution of cavity fields, while $g(u)$ is the eigenfunction associated to the largest eigenvalue $\lambda$ with respect to a perturbation of $P(u)$ (34). $g(u)$ is symmetric, therefore $g^{\prime}(u)$ is anti-symmetric. We impose $\int_{-\infty}^{\infty}g(u)du=1$.

We impose normalization:

[TABLE]

obtaining the relation $c_{0}=a$.

The functional form 48 has to reproduce itself when attaching two chains to create a new one of length $L_{1}+L_{2}=L$, according do the rules of Table 1. Using the symbol $\otimes$ to denote the iteration in distribution of two chains according to these rules, or the addition of a field to an extremity of a chain, we have to check that

[TABLE]

where $Q_{z-2}$ is the distribution of the sum of $z-2$ cavity fields extracted from $P(u)$ plus the random external field. More explicitely, we have

[TABLE]

where the functions $f_{1}$, $f_{2}$ and $f_{3}$ can be deduced from Table 1. A careful computation of the leading order terms in the right hand side, shows that the Ansatz (48) holds, provided that

[TABLE]

where

[TABLE]

It turns out that $\rho$ and $c_{1}$ are left undetermined. The final form for the Ansatz is thus given by:

[TABLE]

In this form, the Ansatz is normalized and stable under the merging of two chains up to $o(L\lambda^{L})$ terms. The coefficient $\rho$ can be derived using a few additional arguments, see next Section. This form for the Ansatz is compatible with what we know from the previous sections: the coupling $J$ either is exactly 0 or, with a probability of order $L\lambda^{L}$ , is of order $1/L$. The two effective fields and the coupling are independently distributed when conditioning on the event $J_{L}>0$. The two fields have correlation of order $O(L\lambda^{L})$ when conditioning on $J_{L}=0$ instead. Moreover, the Ansatz reproduces the numerical behavior of the correlation functions at length $L$.

Now we want to reproduce the numerical results for the loop contribution to the correlation functions using the Ansatz. In order to obtain the joint distribution $P_{\text{loop}}(u_{I},u_{O},J_{\text{loop}})$ of the effective fields and coupling for two spin at the extremities of a loop as in the right part of Fig. 6, we convolve the two internal branches $L_{A}$ and $L_{B}$, yielding a distribution that we call $(P_{L_{A}}*P_{L_{B}})$ on the two internal spins, and attaching the external legs $L_{I}$ and $L_{O}$:

[TABLE]

We have already said that the loop contribution to the observable is given by the value of the observable computed on the loop minus the observable computed on the two paths $L_{I}+L_{A}+L_{O}$ and $L_{I}+L_{B}+L_{O}$ considered as independent. We can easily obtain this loop correction defining the “topologically connected” loop distribution $\tilde{P}_{\text{loop}}(u_{I},u_{O},J_{\text{loop}})$ as in eq. (55) but substituting to $P_{L}$ the (improper) distribution $\widetilde{P}_{L}$ given by

[TABLE]

that is, the same as $P_{L}$ but without its asymptotic term. In this way, the loop correction is just the mean value of the observable on $\tilde{P}_{\text{loop}}(u_{I},u_{O},J_{\text{loop}})$. The loop correction for both the connected and disconnected correlation function computed on $\tilde{P}_{\text{loop}}(u_{I},u_{O},J_{\text{loop}})$ gives 0. While this is in agreement with the numerical computation for $R$, we had been expecting a non-zero contribution for $C$. However, being the Ansatz consistent up to order $O(L\lambda^{L})$, it could only give a contribution to $C_{loop}(L,L,L,L)=O(L\lambda^{\Sigma(\vec{L})})$ that in fact is not present from the numerical analysis (Please notice that higher contributions are prohibited for symmetry reasons). Thus the analytical Ansatz predictions are fully compatible with the numerical results up to the chosen order.

To go to next order, we should introduce terms $O(\lambda^{L})$ in the ansatz. Unfortunately, the addition of new terms in the original Ansatz makes the computation much more involved, and we did not perform it entirely. In particular these new terms should take into account correlations between fields and coupling in the $J\neq 0$ part, as found from the numerical analysis. However terms in the $J=0$ part can be added without much effort, in particular we added the terms $-b_{0}\lambda^{L}g(u_{0})g(u_{L})-b_{1}\lambda^{L}g^{\prime}(u_{0})g^{\prime}(u_{L})-b_{2}\lambda^{L}g^{\prime\prime}(u_{0})g^{\prime\prime}(u_{L})$ and checked how they behave under iteration. Imposing normalization and self-consistency we do find that $b_{0}=b_{2}=0$ and $b_{1}=\frac{1}{(2\hat{P}(h=0))^{2}}$. The addition of this new term gives no contribution to $R$ while it gives a contribution $C_{loop}(L,L,L,L)=O(\lambda^{\Sigma(\vec{L})})$ for the disconnected correlation function, as found from the numerical computation. We stress however that we lack some terms coming from the correction of the $J\neq 0$ part of the Ansatz to order $O(\lambda^{L})$ that do not allow us to compute exactly $C_{loop}(L,L,L,L)$ at order $O(\lambda^{\Sigma(\vec{L})})$.

5 Computation of the mean coupling decay

The analytical Ansatz presented in the previous Section requires the knowledge of 2 parameters: $c_{1}$ and $\rho$. Here we show how to compute the latter in a very effective way. We follow the ideas of Ref. (29), but correcting an error made in that work.

In practice we are interested in computing the mean value of the effective coupling at distance $L$ along the branches of the BL, where the coupling is non zero

[TABLE]

Without loss of generality and to make analytical expressions more compact we fix the single link coupling to $J=1$ hereafter.

A possible numerical method has been already discussed in the previous sections and consists in evolving a population of triplets $(u_{1},u_{2},J)_{L}$ reweighted in a such a way that triplets with $J\neq 0$ do not decrease exponentially fast but remain constant in number: this trick allows to follow triplets with non-zero effective coupling for a long enough time to measure accurately the exponent $\rho$. As an example we show in Figure 11 the inverse of the mean effective coupling as a function of $L$, measured at the critical point $\sigma_{h,c}\simeq 1.037$ for $z=3$.

Although the fit shown in Figure 11 is very good and provides an estimate to $\rho=0.755(1)$ we have to remind that the reported uncertainty only represent the statistical error given the fitting function. It is much more difficult to estimate the systematic error, that would depend — among others — on the corrections to the asymptotic scaling. For this reason we would be much more confident if we could derive an analytical expression for $\rho$.

Given that the effective couplings becomes very small even on the BL branches where they are non-zero, we would like to exploit this observation to better study the asymptotic distribution of cavity messages. Let us consider the equations for updating the triplets reported in Table 1 and let us rewrite it in a more explicit form, concentrating on the messages acting on the spin at distance $L$ (we ignore the messages arriving on the spin at the root). Schematically we have that, adding one new link, the messages change according to the following rules ($*$ messages are irrelevant in the present computation)

[TABLE]

where $h_{z-2}\sim Q_{z-2}$, i.e. $h_{z-2}=h^{R}+\sum_{i=1}^{z-2}u_{i}$ in distribution, and the functions are defined as follows

[TABLE]

From the above expressions we understand that during the evolution with probability $\mathbb{P}[1+J<|h|]$ the effective coupling becomes null, but we are interested in the complementary events, when the effective coupling remains non-zero. With probability $\mathbb{P}[|h|<1-J]$ the coupling remains unaltered and with probability $\mathbb{P}[1-J<|h|<1+J]$ it decreases. We notice that the last event becomes very rare in the limit of small $J$, because the random variable $h$ has a continuous probability density function $P_{h}$ with no Dirac deltas in 1 or -1, so $\mathbb{P}[1-J<|h|<1+J]\simeq(P_{h}(1)+P_{h}(-1))2J$.

In practice in the large $L$ limit, when all effective couplings are very small, $J_{L}\ll 1$, the evolution proceeds essentially by keeping the $J_{L}$ constant until it jumps directly to $J_{L}=0$.

Let us move now to the analysis of the cavity messages $u_{L}$. Assuming we are in the large $L$ limit and all effective couplings are very small, we can work under the above hypothesis that the effective coupling stays constant in $L$ until it becomes null. So hereafter we fix $J_{L}=J$, where $J\ll 1$ is a small constant. We call $Q_{J}$ the probability distribution of the cavity messages on the branches where the effective coupling is fixed to $J$. From Eq. [64] it is easy to derive that asymptotically $Q_{J}$ has support in $(-1+J,1-J)$ and satisfies the following equation

[TABLE]

where $\mathbb{I}$ is the indicator function and the normalizing factor $\lambda(J)$ is given by

[TABLE]

where $u^{(J)}\sim Q_{J}$. In practice $\lambda(J)$ is the rate of survival of a non-zero effective coupling equal to $J$. From this we can obtain the probability distribution of couplings in the $J\ll 1$ limit

[TABLE]

Given that we are mostly interested in studying the decay at the critical point it is worth reminding that at criticality $\lambda(0)=1/(z-1)$ holds and thus we have

[TABLE]

In Figure 12 we show data for $\lambda(J)$ computed at criticality for $z=3$ together with best interpolation via the following function

[TABLE]

The curve interpolates perfectly the data within the statistical uncertainties and it returns an estimate $\rho=0.754(1)$, compatible with the numerical estimate coming from the triplets evolution described in previous sections.