{'fraction_non_alphanumeric': 0.03343991473487876, 'fraction_numerical': 0.009059419131361577, 'mean_word_length': 4.665157232704402, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
The self-organized multi-lattice Monte Carlo simulation

The self-organized Monte Carlo simulations of 2D Ising ferromagnet on the square lattice are performed. The essence of devised simulation method is the artificial dynamics consisting of the single-spin-flip algorithm of Metropolis supplemented by the random walk in the temperature space. The walk is biased to the critical region through the feedback equation utilizing the memory-based filtering recursion instantly estimating the energy cumulants. The simulations establish that the peak of the temperature probability density function is located nearly the pseudocritical temperature pertaining to canonical equilibrium. In order to eliminate the finite-size effects, the self-organized approach is extended to multi-lattice systems, where feedback is constructed from the pairs of the instantaneous running fourth-order cumulants of the magnetization. The replica-based simulations indicate that several properly chosen steady statistical distributions of the self-organized Monte Carlo systems resemble characteristics of the standard self-organized critical systems.


Introduction
The Monte Carlo (MC) simulation methods are nonperturbative tools of the statistical physics developed hand in hand with increasing power of nowadays computers. The benchmark for testing of MC algorithms represents the exactly solvable Ising spin model [1]. Between algorithms applied to different variants of this model, the using of the single-spin-flip Metropolis algorithm [2,3] prevails due to its simplicity. Nevertheless the principal problems have emerged as a consequence of the accuracy and efficiency demands especially for the critical region.
As a consequence of this several methods enhancing MC efficiency have been considered. The procedures of great significance are: finite-size-scaling relations and renormalization group based algorithms [4] - [7], cluster algorithms [8,9] lowering the critical slowing down, histogram and reweighting techniques [10] interpolating the stochastic data, and multicanonical ensemble methods [11] overcoming the tunneling between coexisting phases at 1st order transitions. Even with the mentioned modifications and related MC versions, a laborious and human assisted work is needed until a satisfactory accuracy of results is achieved. This is trivial reason why utilization of the self-organization principles have attracted recent attention of MC community.
In the present paper we deal with the combining of the self-organization principles and MC dynamics. Of course, this effort has a computational rather than physical impact. Our former aim was to design the temperature scan seeking for the position of the critical points of the lattice spin systems. But this original aim has been later affected by the general empirical idea of the self-organized criticality (SOC) [12], originally proposed as a unifying theoretical framework describing a vast class of systems evolving spontaneously to the critical state. The SOC examples are sand pile, forestfire [13] and game-of-life [14] models. The SOC property can be defined through the scale-invariance of the steady state asymptotics reached by the probability density functions (pdf's) constructed for spatial and temporal measures of dissipation events called avalanches. The dynamics of standard SOC systems is governed by the specific nonequilibrium critical exponents linked by the scaling relations [15] in analogy to standard phase transitions.
Notice that dynamical rules of standard SOC systems are functions of the microscopic parameters uncoupled to the global control parameters. On the contrary, in the spin systems governed by the single-spin-flip Metropolis dynamics, the spins are flipped according to prescriptions depending on their neighbors, but also on the global or external control parameters, like temperature or selected parameters of the Hamiltonian. The modification, by which the MC spin dynamics should be affected to mimic the SOC properties, is discussed in [16]. In agreement with [17], any such modification needs the support of nonlinear feedback mechanism ensuring the critical steady stochastic state. It is clear that the feedback should be defined in terms of MC instant estimates of the statistical averages admitting the definition of critical state. The real feedback model called probability-changing cluster algorithm [22] was appeared without any reference to the general SOC paradigm. The alternative model was presented in [18], where temperature of the ferromagnet Ising spin system was driven according to recursive formula corresponding to general statistical theory [19]. This example was based on the mean magnetization leading to series of the temperature moves approaching the magnetic transition. Despite the success in the optimized localization of critical temperature of the Ising ferromagnet, the using of term SOC seems to be not adequate for this case. The reason is the ab-sence of the analysis of spatio-temporal aspects of MC dynamics, which can be considered as a noncanonical equilibrium [20,21], due to residual autocorrelation of the sequential MC sweeps.
Regarding the above mentioned approaches, the principal question arises, if the MC supplemented by feedback, resembles or really pertains to branch of the standard SOC models which are well-known from the Bak's original paper [12] and related works.
The plan of our paper is the following. Sec.2 is intended to generalization of the averaging relevant for the implementation of the self-organization principles. The details of feedback construction based on the temperature gradient of specific heat are discussed in Sec.3. These proposals are supplemented by the simulation results carried out for 2D Ising ferromagnet. The details of the multi-lattice self-organized MC simulations, stabilizing true critical temperature via fourth-order magnetization cumulants, are presented in Sec.4. The oversimplified mean-field model of self-organized MC algorithm is discussed in Sec.5. Several universal aspects of the self-organized MC dynamics are outlined by replica simulations in Sec.6. Finally, the conclusions are presented.

The running averages
As we already mentioned in introduction, the mechanism by which many body system is attracted to the critical (or pseudocritical) point should be mediated by the feedback depending on the instantaneous estimates of the statistical averages. In this section we introduce the running averages important to construct the proper feedback rules.
Consider the MC simulation generating the sequence of configurations { X t ′ , t ′ = 1, 2, . . . , t } according to importance sampling update prescription of Metropolis [3] producing the canonical equilibrium Boltzmann distribution as a function of the constant temperature T . The estimate of the canonical average A t of some quantity A can simply calculated from the series of t sampled real values A t ′ ≡ A(X t ′ ) reweighted by w t ensuring the trivial normalization The summation given by Eq.(1) is equivalent to recurrence showing how the average changes due to terminal contribution A t . Consider the generalized averaging, where w t is replaced by the constant parameter 0 < η ≪ 1 which is independent of t: The consequence of setting w t → η is the convolution defined by the modified weights which undergo to normalization For finite initial choice t ′ = 1 it yields where (1−η) t ≪ 1 if the time of averaging is sufficiently large η ≪ 1 ≪ t. It should be remarked that the generalized averages labeled by A η,t,T are equivalent to gamma filtered [23] fluctuating inputs A t . Note that the term gamma originated from the analytic form of the weight w η,t,t ′ . The filtering represents the application of the selection principle suppressing the information older than the memory depth ∝ 1/η.

The specific heat feedback
To attain the critical region self-adaptively we construct the temperature dependent feedback changing the temperature in a way enhancing extremal fluctuations. The running estimates of averages are necessary to predict (with share of the uncertainty) the actual system position and the course in the phase diagram leading to critical point. The pseudocritical temperature T c (N ) of some finite system consisting of N degrees of freedom is defined by the maximum C(N, T c (N )) of the specific heat C(N, T ). To form an attractor nearly T c (N ), we propose the following dynamics of the temperature random walker biased by the gradient where C is MC estimate of the specific heat. From that follows that the energy fluctuations are controlled by the temperature representing the additional slowly varying degree of freedom. It is assumed here and in further that T t remains constant during N random microscopic moves (1MC step per N ). The sign function occurring in Eq.(9) is used to suppress the extremal fluctuations of F C t causing the unstable boundless behavior of T t . From several preliminary simulations it can be concluded that the replacement of sign function by some smooth differentiable function (e.g. tangent hyperbolic or arcus tangent) seems to be irrelevant for keeping of smaller dispersion of T t . The non-constant temperature steps are constrained by |T t+1 − T t | < ∆ due to action of the pseudorandom numbers r t drawn from the uniform distribution within the interval (0, 1).
Very important for further purposes is the quasiequilibrium approximation A t ≃ A η,t,Tt justified under the restrictions t eq,A η ≪ 1, t eq,A ∆ ≪ 1, where t eq,A is the equilibration time of A. With help of this approximation, the running averages from Eq.(4) can be generalized using which allow the averaging under the slowly varying temperature. Here T t is the temperature for which the last For later purposes we also introduce the zero passage time τ A pertaining to A. The time is defined as a measure of the stage during which the sign(A t ) is invariant. More formal definition of τ A requires the introducing of two auxiliary times t ′ , t ′′ . The first time t ′ defines the instant, where whereas t ′′ counts the events for which The counting is finished for To be thorough, the conditions A } k=1 . From the above definition it follows that random walk in temperature is unidirectional (the sign of T t+1 − T t is ensured) for t ′′ = 1, 2, . . . , τ F C . The arithmetic average of τ F C is related to temperature dispersion Here the specific heat is expressed in terms of the energy E cumulants E t , E 2 t . In the frame of the quasi-static approximation it is assumed: Tt . Subsequently, using properties of energy cumulants with equilibrium Boltzmann weights, the temperature derivative of C can be approximated by The Ising ferromagnet is simulated in further to study the effect of feedback defined by Eqs.(9), (10) and (16). However, it is worthwhile to note that many of the presented results are of general relevance. Given the spin of the square L × L lattice with the periodic boundary conditions, the Ising Hamiltonian can be defined in exchange coupling units where nn means that summation running over the spin nearest neighbors. In general, the dynamics of SOC systems exhibits two distinct regimes. During the transient regime the proximity of critical or pseudocritical point is reached. The second steady regime is called here the noncanonical equilibrium in analogy with [21]. In this regime the attraction to critical point is affected by the critical noise. This general classification is confirmed by our results shown in Figs.1-3. In Fig.1 we see the stochastic paths of T t pertaining to different initial values of energy cumulants and spin configurations. The paths are attracted by T c (N ) with some uncertainty in noncanonical equilibrium. For the sufficiently narrow steady pdf's, T c (N ) can be approximated by where N av is the number of inputs. The stationary pdf of T t walk is shown in Fig.2a and pdf of F C t with non-Gaussian flat tails is depicted in Fig.2b.
The alternative quantity capable for the characterization of the noncanonical equilibrium is the autocorrelation The simulations results depicted in Fig.3a evidenced that minimum time for which the anticorrelation (K τ < 0) occur is of the order O( τ F C ). As we see from Fig.3b, the power-law dependence can be identified within the region of vanishing τ F C . More profound discussion of this fact is presented in Sec.6.

Multi-lattice simulations
In this section we try to avoid the problem of finitesize-scaling related to true equilibrium critical temperature T c and critical exponents. The problem is solved by the multi-lattice self-organized simulations based on the dynamical rules treating the information from running averages of magnetization. The considerations are addressed to models, where 2nd order phase transitions take place. The proposal is again applied to 2D Ising ferromagnet on the square lattice. The quantity indicating deviations of magnetization order parameter m L = (1/L 2 ) ij s i from the gaussianity is the fourth-order cumulant The standard way leading to true T c = lim N →∞ T c (N ) is the construction of the temperature dependences U L l ,T , U Ls,T for two lattices L s = L l . Then T c follows from the condition of the scale-invariance According it the self-organized multi-lattice MC simulation method consists of the following three main points repeated in the canonical order for the counter t = 1, 2, . . . : 3. The temperature shift biased to eliminate difference If L l > L s the ordering of cumulants in F U t,ls is chosen subject to assumption Any modification of the Eqs. (22) and (23) is possible when the preliminary recognition of the critical point neighborhood is performed. The parametric tuning recovers that stabilization of the noncanonical equilibrium via feedback F U t,ls requires smaller η and ∆ than single-lattice simulations based on action of F C t . The Eq.(22) can be generalized for n r lattices, i.e. for n rp = nr 2 (n r − 1) competing lattice pairs labeled by l, s: where 1/n rp term rescales additive contributions. The presented method also offers the continuous checking of estimated critical exponents. It comes from the standard assumption that canonical equilibrium magnetization exhibits critical scaling |m is the scaling function and β/ν is the ratio of magnetization (β) and the correlation length (ν) critical exponents, respectively. If the temperature fluctuates nearly T c , the equilibrium finite-size-scaling relation changes to |m Lj | η,t,Tt=Tc ≃ L −β/ν j f m (0). For n r ≥ 2 lattices and sufficiently small η, ∆, the following arithmetic average can be defined Similar to treatment of the steady temperature fluctuations, the quantity β/ν ≃ N −1 av Nav t=1 (β/ν) η,t can be defined. The simulations carried out for cases n r = 2, 3 are compared in Fig.4. In agreement with expectation, the localization of T c for n r = 3 with recursion taken from Eq.(25) is much subtle than for n r = 2. In addition, the statistics of (β/ν) η,t is weakly depending on η and ∆, which seems to be logical due to universality of exponents in the canonical equilibrium limit (T c = const).

The mean-field analysis of algorithm
In this section we present calculations aimed to understand how the attractivity of critical point arises and how the noncanonical equilibrium is attained by means of feedback. Only a rough approximation of the complex simulation process is considered, where spin degrees of freedom are replaced by the unique magnetization (per site) term m(t). Furthermore, it assumes that the selected central spin s(t) flips in a mean field created by its z neighbors. Let π(t) denotes the probability of the occurrence of s(t) = 1/2 state, then the probability of s(t) = −1/2 is 1 − π(t). The master equation for π(t) can be written in the form The Glauber's [24] heat bath dynamics with the transition probabilities W [−+] and W [+−] between states s(t) = ±1/2 is preferred in comparison to nondifferentiable Metropolis form due to analyticity arguments relevant for formulation by means of differential equation. Within the mean-field approximation it can be assumed In the above expression τ f is the time associated with the spin flip process. The expression takes into account ±zm(t) variations of energy belonging to flips from s = ±1/2 to s = ∓1/2 given by the effective single-site Hamiltonian −z s(t) m(t). Assuming that π = m + 1/2 and using Eqs.(27), (28) we obtain Subsequently, the feedback differential equation of the temperature variable is suggested in the form where m c is the "nucleation" parameter of the ferromagnetic phase, and α > 0 is the constant parameter. Unlike the works [18,19], where feedback consisting of |m| term is considered, the m 2 is absorbed to the feedback F m (t) proposed here to ensure the analyticity. For  6 The comparison of MC and SOC dynamics.
In the section we discuss the universal aspects of the non-equilibrium self-organized MC dynamics applied to the canonical Ising model. As is already mentioned, the attributes of the SOC systems are avalanches reflected by the power-law pdf distributions. We follow with the construction of certain temporal characteristics by supposing their uncertain links to avalanches. By using Eqs.(12)-(14) the evolution of any quantity can be mapped to the sequence of passage times. The example of this view represents pdf of τ F C depicted in Fig.2. Because of the substantial difference in the exponents of pdf's belonging to τ F C and τ F U , no universality attributes are indicated. More encouraging should be to find of pdf's independent of the feedback type. The natural way toward this aim seems to be the investigation of the passage time sequences linked to the order parameter of the canonical equilibrium of given system. In the case of the Ising model the ordering is described by the magnetization, or, eventually by isolated spin value. Therefore, it seems to be logical to define the passage times τ m , τ s given by Eqs.(12)-(14) (corresponding to A = m L , s i , where i is the arbitrary but fixed site position). The simulation results are depicted in Fig.6. In their structure, the following attributes relevant for interpretation in terms of SOC can be identified: (I) the power-law behavior with the unique exponent φ pertaining to different feedbacks F C t , F U t , i.e. to single-lattice and twolattice systems; (II) the interval of dependence from (I) broaden with the size of lattices  Below T c pdf splits into separable contributions fitted here by the bimodal distribution P τs = b 0 P τs (p 0 ) + b 1 P τs (p 1 ) with parameters b 0 = 0.003, b 1 = 0.004, p 0 = 0.61, p 1 = 0.02. The supplementary analysis of statistics of the successive time differences of m L recovers that b 1 P τs (p 1 ) term originates from the mechanism of the long-time tunneling among nearly saturated states of the opposite polarity. From the figure it also follows that power-law short-time regime described by Eq.(33) is formed only if the feedback mechanism is activated. This conditional occurrence of universality can be considered as an additional (IV)th attribute relevant for identification of SOC. The noncanonical equilibrium attained by the self-organized MC dynamics for L = 10 leads to P τs dependence, which can be approximated by the fit with parameters b 2 = 0.00264, τ 2 = 43.521. Let us to note that for the sand pile model [12] the spatial measure of avalanche is associated with the energy integral taken during the stage following distur- bance. In the case of MC the analog of spatial measure can be the extremal magnetization where t In agreement with SOC attribute labeled I, the independence of feedback type is identified.
It should be noticed that there remains the principal problem of the link between SOC and presented selforganized MC method focused to the critical region. The problem is to identify the specific algorithmic segment which should be interpreted as an analogue of a disturbance initializing the avalanche. Fortunately, the advanced MC approach exists through which a disturbance can be absorbed into MC algorithm with minor violation of the original dynamics. In general, the approach of interest based on the coevolution of given system and its replica is known under the term damage spreading technique [25]. To apply it, let us consider the self-organized MC referential system labeled here as {1} t , which incorporates instant cumulants (single with F C t or multi-lattice based on F U t ), spin configurations and instant temperature T . Consider also the replica counterpart {2} t of the system {1} t . As is known, the parallel simulation of {1} t and {2} t should be applied with the identical pseudorandom sequences. Canonically, the measure of damage effect is then defined through the single-time two-replica difference where m 2. the replica temperature is modified according where constant ǫ T causes the small disturbance of the coincidence of the contents of {1} t (k) Evidently, the temperature disturbance plays role similar to adding of the grain to the sand pile. The sand pile is stabilized if the rest state occurring for t = t Using this, the simulated path is mapped onto the sequence {D (k) max } k=1 . Consequently, the pdf's can be extracted which are depicted in Fig.8. They show that the effective exponents centred nearly −2.3 fit simulated pdf's fairly well. As in the case labeled I, pdf's related to τ (k) D and D (k) max are weakly susceptible to the feedback choice. Since both temporal and spatial power-law attributes of universality are indicated, the standard SOC paradigm can be considered as a framework adaptable for the analysis of the suggested self-organized MC dynamics.

Conclusions
Several versions of MC algorithm combining the selforganization principles with the MC simulations have been designed. The substantial feature of method is the establishing of running averages coinciding with gamma filtering of the noisy MC signal. The simulations are combined with the mean-field analysis describing the motion of temperature near to magnetic transition point. The replica-based simulations indicate that pdf distributions of passage times in a noncanonical equilibrium attain the interval of the power-law behavior typical for the standard SOC pdf distributions. We hope that the present contribution will stimulate further self-organized studies of diverse lattice models, e.g. those related to percolation problem.