Free energy and metastable states in the square-lattice $J_1$-$J_2$ Ising model

Free energy as a function of polarization is calculated for the square-lattice $J_1$-$J_2$ Ising model for $J_2<|J_1|/2$ using the random local field approximation (RLFA) and Monte Carlo (MC) simulations. Within RLFA, it reveals a metastable state with zero polarization in the ordered phase. In addition, the Landau free energy calculated within RLFA indicates a geometric slab-droplet phase transition at low temperature, which cannot be predicted by the mean field approximation. In turn, restricted free energy calculations for finite-size samples, exact and using MC simulations, reveal metastable states with a wide range of polarization values, but with only two domains. Taking into account the dependence of the restricted free energy on the nearest-neighbor correlations allows us to identify several more metastable states. The calculations also reveal additional slab-droplet transitions at $J_2>|J_1|/4$. These findings enrich our knowledge of the $J_1$-$J_2$ Ising model and the RLFA as a useful theoretical tool to study phase transitions in spin systems.


Introduction
The square-lattice J 1 -J 2 Ising model is one of the minimal extension of the standard Ising model, in which the coupling J 1 between nearest neighbors is complemented by the coupling J 2 between diagonally next-nearest neighbors.The properties of this model are of both fundamental and practical interest, in particular, since its quantum Heisenberg counterpart is relevant to the antiferromagnetism in the parent compounds of the cuprate and pnictide families of high-temperature superconductors [1][2][3].Indeed, recent state-of-the-art numerical calculations [4][5][6][7][8][9][10][11][12][13][14] confirm earlier findings [15][16][17][18][19][20][21][22] that diagonal interactions are important in describing the available experimental data for these compounds.Spin frustrations due to the J 2 coupling affect the temperature and order of the phase transition and lead to a valencebond-solid state in the transverse field Ising model [23,24] and to a quantum spin liquid state in the Heisenberg model when J 2 is close to |J 1 |/2 [25][26][27][28].
We recently highlighted the existence of metastable states with arbitrary polarization in the square-lattice J 1 -J 2 Ising model in the interval J 2 ∈ (0, |J 1 |) using Monte Carlo (MC) simulations, which was further supported by simple microscopic energy considerations [29].For the ferromagnetic ground state, i.e. for J 1 < 0 and J 2 < |J 1 |/2, these states are rectangles with polarization opposite to the surrounding, briefly considered much earlier in [30,31].Note that these states differ from the metastable states of the standard Ising model, consisting of straight stripes, into which a system with zero polarization, when applying the single-spin flip MC algorithm, relaxes after quenching only in about 1/3 of cases and only in the absence of an external field [32][33][34].Significantly, the random local field approximation (RLFA) [35], also applied in [29], points to a metastable state with zero polarization in the same J 2 coupling range, thus reflecting the appearance of microscopic metastable states, which seems impossible for the mean field approximation (MFA).
The polarization-dependent Landau free energy F (m), considered in the framework of Landau's phenomenological theory of phase transitions [36], is an effective tool for studying metastable states (including those with polarization opposite to the external field) and can be used to calculate the relaxation rate of a system to the ground state using the Landau-Khalatnikov equation [37] (see, e.g., [38] for such calculations in ferroelectrics).It should be noted that for short-range interactions, the Landau free energy obtained within MFA differs qualitatively from the restricted free energy calculated for a fixed value of polarization exactly or using the MC method for finite-size samples [39][40][41].In the former case, below the phase transition temperature, the Landau free energy has a double-well shape.In the latter case, at a temperature close to zero, the restricted free energy is determined by states with minimal energy, and the barrier between its two minima with opposite polarization corresponds to states with two large domains, the interface energy of which is proportional to the sample size L. Thus, relative to the total energy, proportional to the number of spins N = L 2 , the barrier vanishes in the thermodynamic limit [42,43].It was shown that, despite the loss of detailed information about microscopic spin configurations, the restricted free energy can be harnessed to well reproduce the MC polarization dynamics of the Ising model in good agreement with the droplet theory [44].These ideas were further developed in [45][46][47][48].The temperature dependence of the restricted free energy barrier in the three-dimensional J 1 -J 2 Ising model was analytically estimated in [31] in connection with domain growth and shrinking after low-temperature quenching.
Here, for the square-lattice J 1 -J 2 Ising model, we calculate the Landau free energy in the RLFA framework and the restricted free energy exactly for a square sample of size L = 6 and using the MC method for L = 10 and L = 14.We pay special attention to the metastable states, which appear in this model and were studied earlier in [29], and explore how they are reflected in the free energy.The features of the geometric slab-droplet phase transition in the free energy calculated by both methods are also briefly discussed.

The J 1 -J Ising model
The square-lattice J 1 -J 2 Ising model Hamiltonian reads where each spin s i takes the value +1 or −1.The sums are over nearest 〈i, j〉 and diagonal next-nearest 〈〈i, j〉〉 neighbors, as well as over each spin coupled to the external field h i at its position.In what follows, we set the values of the coupling constants J 1 = −1 and J 2 < 1/2, which correspond to the ferromagnetic ground state (the case J 2 > 1/2 with a striped antiferromagnetic ground state is similar in many aspects, but has a more complex spin topology and will be considered separately).Note that the model is invariant with respect to the simultaneous change of the sign of J 1 and the replacement of homogeneous polarization with Néel checkerboard one, corresponding to the antiferromagnetic order of the parent compounds of cuprate superconductors [49].

Landau free energy in the random local field approximation
RLFA is based on the exact formula for the average spin [35,50]: where β = 1/T is the inverse temperature in energy units.The local field, h s i = − j J i j s j , acting on the spin s i is caused by all spins s j coupled with it.
The brackets in Eq. ( 2) correspond to thermal averaging, which is performed with probability distribution [35,51] where the product is taken over all spins s j coupled to s i , and m j = 〈s j 〉 = me iqr j is the thermally averaged polarization at position r j , the variation of which in space is determined by the propagation vector q.The model applies to both ferromagnetic and ferroelectric materials, so we call m polarization, which can be magnetic (magnetization) or electric.Here we consider only homogeneous polarization m and external field h, which corresponds to q = (0, 0).Note that Eq. ( 3) implies that, within RLFA, the fluctuations of each spin are considered independent.
Eq. ( 2) follows from equating to zero the derivative of the Landau free energy F (m) [44], which corresponds to thermodynamic equilibrium at a fixed value of polarization m [52].In Figure 1: (a) The RLFA solution for J 2 = 0.3 (see also Fig. 2 in [29] for different values of J 2 ).(b) The Landau free energy F 0 (m) as a function of temperature within RLFA for J 2 = 0.3.Red dots correspond to the local maximum of F 0 (m) at each temperature (i.e.barrier), dark blue dots correspond to its global minimum (stable states), and light blue dots correspond to its local minimum (metastable states).A metastable state with zero polarization exists at temperatures from zero to T 0 and from T 1 to T c ≈ 1.26 < T 2 , see also (a), at which a first order phase transition occurs according to RLFA.order to reconstruct F (m) and obtain its correct dependence on the external field h, we rewrite Eq. (2) in the form ∂ F/∂m= f (m) − h, integration of which over m yields F (m) − F (0)= F 0 (m), where F 0 (m) = m 0 f (m)d m − hm is the Landau free energy minus the integration constant, which depends only on temperature.
The RLFA solution and the Landau free energy F 0 (m) calculated in this way for J 2 = 0.3 are shown in Fig. 1.Here and in what follows, only values m > 0 are considered due to the symmetry F (m) = F (−m) in the absence of an external field.The calculated free energy indeed points to the metastable state with m = 0, discussed in [29].The dependence of the metastable state barrier on temperature and coupling constant J 2 is discussed in more detail in Appendix A.1.It should be noted that within RLFA the transition turns out to be first order for 0.25 ≲ J 2 ≲ 1.25 [29], which is comparable to the cluster MFA results [53][54][55][56], while recent more accurate calculations narrow this interval to a small region around J 2 = 1/2 [53,[57][58][59][60]. RLFA also strongly underestimates the critical temperature T c just below J 2 = 1/2 [29], which could be attributed to increased frustration in this region.On the other hand, just above this value, where the ground state is striped antiferromagnetic, the accuracy of RLFA is quite high.It is also worth noting that the RLFA Landau free energy becomes flat in the limit J 2 → 1/2 (see Fig. 2), which indicates the absence of a phase transition at this point [29].
At low temperature, the RLFA Landau free energy shows a kink at a polarization value m c ≈ 0.5 for J 2 = 0.3, see Fig. 1b and Fig. 2. A similar kink observed in the restricted free energy [44,61,62] is attributed to the geometrical phase transition [63,64] and corresponds to the transition from a slab (for m < m c ), as the most probable configuration, to the droplet (for m > m c ), which will be discussed in the next section.For J 2 = 0, the RLFA prediction for m c ≈ 0.53 (Fig. 2) is close to the exact value of the critical polarization m c = 0.5 [63].We emphasize that RLFA is able to predict the geometrical phase transition, in contrast to MFA.We also checked that the mean field cluster approximation (even for cluster sizes up to 4 × 4 spins), formulated as in [65], does not predict this transition, despite the good accuracy of the approximation in describing ferroelectric phase transitions [66][67][68][69].We can explain this as T = 0.1 follows.An essential part of obtaining the Landau free energy, as discussed above, is finding the dependence of the external field h on the polarization m using Eq. ( 2).Within RLFA, the right-hand side of this equation is the sum of tanh functions with all possible values of the local field in their argument, while within MFA it contains only one term with the local field corresponding to the mean field.For different values of m, the solution of the equation corresponds to different local fields and, as a consequence, different spin configurations, as shown in Appendix A.2.

Restricted free energy calculated exactly for small samples
Now we want to look at possible metastable states in the free energy of a sample of N spins restricted to a certain value of the total spin M and defined as [39,44] where the sum is over all possible energies E and n(M , E) is the number of configurations for a given set of values, which would be the density of states in the continuous approximation.While the correspondence between polarization and total spin is obvious, m = M /N , here we prefer to consider the variable M , since it takes only integer values that are easy to relate with the corresponding spin configurations.For small samples, the sum in (4) can be computed exactly.For a square sample of size L = 6, yielding the total number of spins N = 36, the result for J 2 = 0.3 is shown in Fig. 3.In all calculations we apply periodic boundary conditions (free boundary conditions are considered in Appendix B.2).The critical temperature for L = 6 is equal to T c = 1.67, while for L = 100 (which practically corresponds to an infinite sample size) it is T c = 1.26 [29].
Configurations that contribute to the free energy F (M ) at zero temperature (i.e. have the lowest energies) for different values of M are shown in Fig. 4. At zero temperature, at M < 16, the restricted free energy is flat with a couple of pits, also mentioned in [44].These pits arise from spin configurations with a completely flat boundary between two slabs with opposite spin directions, see configurations with M = 0 (M 0 ) and M = 12 (M 12 ) in Fig. 4. When a spin flips on this boundary (transitions M 0 → M 2 and M 12 → M 14 in Fig. 4) the energy increases by 4J 1 .The same restricted free energy per spin F 0 /N as a function of total spin M and temperature T .The free energy is defined only for integer even values of M and linearly interpolated between them.For ease of comparison with the RLFA results in Fig. 1, we set here F = 0 at M = 0 at each temperature.Dark blue dots correspond to the global minimum of F 0 /N at each temperature.
0 12  2 144 4 144  6 144  8 144  10 144 12 12  14 144  16 216  18 36  20 72  22 576   24 72  26 288  28 36  30 72  32 72  34 36   36 When the last spin in the row flips, the energy decreases by this value (transition M 10 → M 12 in Fig. 4).The distance between two neighbor pits of F (M ) is equal to 2L, since any single-spin flip changes the total moment by ∆M = ±2.Metastable states inherent in the J 1 -J 2 Ising model can be observed as local minima at M = 18, 20, and 24.It is interesting to note that while the transitions M 20 → M 22 and M 24 → M 26 are accompanied by an increase in energy, states M 18 and M 20 have the same energy, but are not connected to each other by a single spin flip.Thus, each of these two states are metastable in the sense of single-spin flip Monte Carlo dynamics.The same applies to state M 28 in Fig. 4, which cannot be transformed to state M 30 with the same energy by a single-spin flip and is therefore metastable.
At M ≳ N /2, the configurations that contribute to the free energy F (M ) at zero temperature correspond mainly to the droplet (starting from configuration M 18 in Fig. 4), but this depends on J 2 , as indicated in the bottom row of Fig. 4. At J 2 > 1/4 there is a reverse transition to the slab phase for M 24 , and then back to the droplet phase for larger M .The same thing happens for M 22 when J 2 > 1/3.These threshold values of J 2 can be justified by explicitly calculating the energies of these states.We show the restricted free energy for various values of J 2 and the dependence of the metastable state barriers on J 2 and temperature in Appendix B.1.
Since the lowest energies that determine F (M ) at zero temperature imply a minimum number of domains in spin configurations (Fig. 4), various metastable states with a large number of small domains of size at least 2 × 2, observed after low-temperature quenching in [29], do not cause local minima of F (M ) (see Fig. 3 and Fig. 4).However, we can try to reveal these states as local minima of the restricted free energy depending on the nearest spin correlation, e.g.such a function F (M , c 1 ) with c 1 = 〈i, j〉 s i s j was calculated for the standard Ising model in [47].This correlation can be easily and unambiguously calculated and is convenient for determining metastable states in the sense of Monte Carlo dynamics, since after a single-spin flip it changes only by ∆c 1 equal to 0, ±4, ±8.Here we calculate the restricted free energy F (M , c 1 , c 2 ), also taking into account the correlation with the next-nearest neighbor, c 2 = 〈〈i, j〉〉 s i s j .The result for the energy E(M , c 1 , c 2 ) and entropy, S(M , c 1 , c 2 ) = log n, which together constitute the free energy F = E − T S, is shown in Fig. 5 as a color map.Note that in this case the energy is simply a linear function, since Eq. ( 1) for the Hamiltonian becomes H = J 1 c 1 + J 2 c 2 − hM .In the absence of an external field (h = 0), it does not depend on M , but possible correlations c 1,2 do.
Metastable states, the large dark blue dots in Fig. 5, were determined under the condition that any single-spin flip (with the only possible coordinate changes being ∆M = ±2 and ∆c 1,2 = 0, ±4, ±8) will result in equal or greater energy.Thus, these states have the minimum possible energies for a given M and turn out to be the same for all J 2 ∈ (0, 0.5).Their configurations are shown in Fig. 6, which also includes configuration M 28 as discussed above.In addition to the metastable states for F (M ), we reveal for F (M , c 1 , c 2 ) only four higher-energy metastable states with M = 0 and one with M = 12.The metastability criterion used was too stringent to identify higher-energy metastable states that might not have all possible transitions ∆c 1,2 that were tested for energy changes ∆E ≥ 0. Indeed, only a subset of states with a given set of values (M , c 1 , c 2 ) could be metastable, having different connectivity via a single-spin flip with other states, which seems impossible to take into account., with total spin M , nearest and next-nearest neighbor correlations c 1 and c 2 , and number of configurations n for a given set of values.

Restricted free energy calculated for larger samples by Monte Carlo method
For larger square samples, with L = 7 and L = 8, the free energy (4) can only be calculated using powerful supercomputers, given the large number of 2 N configurations for N spins.Alternatively, it can be calculated approximately with sufficiently high accuracy using the Monte Carlo method.We use the Wang-Landau algorithm [70][71][72], which has proven to be very efficient for this purpose at low temperature.It consists in performing a random walk in polarization and energy space to extract an estimate for the number of configurations n(M , E) from Eq. ( 4) that produces a flat histogram.Using the Wang-Landau algorithm, we reproduce the exact results for L = 6 with high accuracy and obtain similar results for L = 10, see Fig. 7, where a much larger number of metastable states are clearly visible at low temperature for J 2 = 0.3.Note that the restricted free energy for larger samples could also be calculated, as was done, for example, in [44] for J 2 = 0.However, for our purposes, namely to show how metastable states are reflected in the free energy, the size L = 10 seems optimal, and all metastable states are clearly visible and convenient for analysis.For larger samples, the total spin M and free energy F will scale as the number of spins N , but the height of the energy barriers will remain the same, i.e. (or half of it for free boundary conditions) for M ≲ N /2 and 4J 2 for M ≳ N /2 (see Fig. 15 in App.B.2 for L = 14 and both periodic and free boundary conditions).The restricted free energy at L = 10 for various values of J 2 and the dependence of the metastable states barrier on J 2 and temperature are shown in Appendix B.1 and compared with the case L = 6.The properties of the restricted free energy depending on nearest-neighbor correlations F (M , c 1 , c 2 ) is similar to the case L = 6.Using the same procedure described above for the case L = 6, several more metastable states can be identified.The energy of possible states and their entropy are shown in Appendix B.3.

Discussion
The primary goal of this work was to consider the metastable states in the J 1 -J 2 Ising model, recently found in the RLFA solution and MC simulations of low-temperature quenching [29], from a free energy perspective.We found that the Landau free energy F (m) as a function of polarization m calculated within RLFA has a local minima at zero polarization below the ferromagnetic phase transition for J 2 ∈ (0, 1/2) (Fig. 1), which indicates a metastable state.
At the same time, the restricted free energy F (M ) as a function of the total spin M calculated exactly for a small sample size L = 6 (Fig. 3) and using MC simulations for L = 10 (Fig. 7) and larger (Fig. 15 in App.B.2) has local minima corresponding to metastable states with various values of polarization.Some of them, namely at M ≲ N /2, are due to long stripes with an activation energy of 4J 1 of a spin flip on a flat domain boundary (Fig. 4).In the standard Ising model, the system can become stuck in these states with a final polarization following a Gaussian distribution after zero-temperature quenching from an initially random configuration with zero polarization [33].
Metastable states at M ≳ N /2 are caused by droplet-shaped domains with an activation energy of 4J 2 of a spin flip in their corner, at least at J 2 < 1/4 for both sample sizes L = 6 and L = 10 (Fig. 13 in App.B.1).At J 2 > 1/4, the dependence of the barrier height on J 2 changes (it is interesting whether this threshold value of J 2 holds true for any sample sizes).This corresponds to a change in the sequence of minimal energy configurations for increasing total spin M (see Fig. 4) with a return to the slab phase and then back again, which may be important for some applications.Indeed, the importance of the geometrical slab-droplet transition for various physical situations, including the dewetting transition between hydrophobic surfaces, was highlighted in [64].
As J 2 → 1/2, some of the barriers disappear, while the rest tend to the value 4|J 1 | at a distance of 2L from each other (Fig. 12 and Fig. 13 in App.B.1).The ground state with two slabs of different widths becomes extensively degenerate for large L [23,58] with T c = 0 [57,73].However, when free boundary conditions are applied, the fully polarized state always represents the global energy minimum at zero temperature, even for J 2 = 1/2, when the energy gap is 2|J 1 | (see Fig. 15 in App.B.2).This is due to the difference in the energy contributions of spins in the corners of domains compared to bulk spins.Note that the position of the slab-droplet phase transition is the same for both boundary conditions at J 2 = 0.
It should be noted that the metastable states into which the system relaxes after lowtemperature quenching in [29] are not exactly the same as shown in Fig. 4 and which determines the free energy at zero temperature.Even considering nearest-neighbor correlations, we were only able to additionally reveal a few highly correlated metastable states in the restricted free energy (Fig. 6).However, combining the metastable states shown in Fig. 6 with the proper domain spacing will also result in a metastable state, but with higher energy, bringing them closer to those observed in [29].Since the energy of the metastable states in [29] is much higher, the system is more likely to get stuck in them, relaxing during quenching on the way to thermal equilibrium.Metastable states like in Fig. 6 can in principle be reached after quenching at non-zero temperature after a sufficiently long relaxation time and domain coarsening, with a higher probability for those closer to the equilibrium polarization.At the same time, any of these states will be reached inevitably if the total spin is conserved during quenching, as in the Kawasaki [74] two-spin exchange algorithm, which is relevant for models describing transport phenomena caused by spatial inhomogenity such as diffusion, heat conduction, etc.
Although the zero polarization of the metastable state and the low height of the barrier proportional to the temperature near zero (see Fig. 8 in App.A.1) is not exactly what follows from MC calculations, where the barrier heights are much higher and decrease with temperature (see Fig. 14 in App.B.1), the fact of even a rough indication of the metastable state by RLFA is very valuable.Another valuable RLFA prediction that turns out to be quite accurate is the geometric slab-droplet phase transition at zero temperature (Fig. 1b and Fig. 2).The reason why RLFA is so effective in this situation, in our opinion, is that by definition it takes into account the local field due to all possible configurations of spins interacting with the central spin, not just the mean field.The probability of these configurations, in turn, is determined by the polarization (see Appendix A.2).
Finally, we will mention some recent advances in the experimental observation of metastable states using sub-picosecond optical pulses, which we believe could be applied to reveal metastable states discussed here and in [29].For instance, in the quasi-two-dimensional antiferromagnet Sr 2 IrO 4 , a long-range magnetic correlation along one direction was converted into a glassy condition by a single 100-fs-laser pulse [75].Atomic-scale PbTiO 3 /SrTiO 3 superlattices, counterpoising strain and polarization states in alternate layers, were converted by sub-picosecond optical pulses to a supercrystal phase in [76].In a layered dichalcogenide crystal of 1T -TaS 2 , a hidden low-resistance electronic state with polaron reordering was reached as a result of a quench caused by a single 35-femtosecond laser pulse [77].See also the references to relevant superconducting and magnetic materials with next-nearest-neighbor interactions mentioned in Introduction and [29].

Conclusion
The Landau free energy calculated within RLFA for the square-lattice J 1 -J 2 Ising model has a local minimum at zero polarization for J 2 ∈ (0, |J 1 |/2) at low temperature in the ferromagnetic state, indicating a metastable state.This reflects the appearance of local minima (metastable states) with various polarization values in the restricted free energy, which was calculated exactly and using MC simulations.The restricted free energy as a function of nearest-neighbor correlations shows several more metastable states, although all of them represent only the lowest energy part (with a minimum number of domains) of the metastable states observed after low-temperature quenching in [29].We also show that RLFA predicts the slab-droplet phase transition as a kink in the polarization dependence of the Landau free energy at low temperature, and explain this by averaging over different spin configurations within RLFA.At the same time, the restricted free energy reveals additional slab-droplet transitions at J 2 > |J 1 |/4, which may be important for applications.Finally, we believe that metastable states should be taken into account when considering the low-temperature behavior of materials described by the J 1 -J 2 Ising model.In addition, easy-to-use RLFA can help reveal the presence of metastable states and geometrical phase transitions in more complex systems, e.g., with site or bond disorder and spin tunneling in a transverse field.for the metastable state at m = 0, which appears in the restricted free energy F (m) calculated within RLFA, as a function of temperature T .Only J 2 < 0.31 are considered, which implies T 0 < T 1 (for temperatures definition see Fig. 1a), since these two temperatures become undefined at larger J 2 [29].For polarization m from zero to about 1/2, the external field h is close to zero (Fig. 9(a)) and is determined by configurations that provide a zero local field h s (see Fig. 10(a) and (b)).This corresponds to a central spin flip without a change in energy, and as a result, the Landau free energy is almost independent of m.This is similar to the behavior of the restricted free energy in the slab phase (see Sec. 4), where a spin flip at the slab boundary does not change the energy.
At m ≳ 1/2, the most probable configurations of neighboring spins change (see Fig. 10(c)) and yield the local field value of h s = 2 (Fig. 9(b)).This implies that a further increase in m will be accompanied by a minimization of domain boundaries and a decrease in F (m), similar to the behavior of the restricted free energy in the droplet phase (Sec.4).This is how the slab-droplet phase transition is reflected within RLFA.
Note that the same mechanism is responsible for the predicting of the zero-polarization metastable state at J 2 ∈ (0, 1/2).The number of local field values is squared in this case (see Fig. 11), since for each local field contribution from the J 1 interaction (Fig. 10) there are the same number of contributions from the J 2 coupling.Within RLFA, the tanh functions in Fig. 11 (and Fig. 9) do not move along the h axis (see Eq. 2), but only their relative weights change with m according to Eq. 3. Thus, for small polarization m in Fig. 11(a), h > 0 (and hence the derivative of the Landau free energy ∂ F/∂m > 0), but at some value of m the weight of configurations with negative local fields increases to such an extent that h becomes zero again (which corresponds to the Landau free energy barrier) and then negative, as in Fig. 11(b).Thus, the entire effect is due to the redistribution of weights of different local fields in Eq. 2 (and the corresponding configurations of the neighboring spins) as m changes.

B Restricted free energy properties B.1 Metastable state barrier dependence on J 2 and T
In Fig. 12 we show side by side the restricted free energy at zero temperature for several values of J 2 , from −0.1 up to 0.5, calculated exactly for L = 6 and using the MC method for L = 10.It can be observed that the height of the barrier of metastable states, ∆F M = F M − F M −2 , depends on J 2 .This is explicitly shown in Fig. 13.At L = 6, for M = 22 at J 2 < 1/3 and for M = 26 at J 2 < 1/4, the barrier height is 4J 2 and is determined by the spin flip at the corner of the droplet.The thresholds for J 2 follow from the energy ratio of the competing configurations for M 22 , M 24 , and M 26 in Fig. 4. At L = 10, the barrier heights for metastable states in the droplet phase (i.e. at M ≳ N /2) also correspond to the spin flip at the corner of the droplet and are equal to 4J 2 for J 2 < 1/4 (Fig. 13b).For J 2 > 1/4, this dependence changes due to slab-droplet and vice versa transitions, as in the case of L = 6.Some barriers, starting from a certain point, gradually disappear, others approach the barrier of the slab metastable state, equal to 4|J 1 |.At higher temperatures, other higher energy configurations in addition to those shown in Fig. 4 contribute to the partition function in Eq. ( 4) for each value of M .This affects the dependence of the above discussed energy barriers on temperature, which for J 2 = 0.3 is shown in Fig. 14.At L = 6, the metastable state barrier at M = 22 disappears at T ≈ 0.65, which is close to the corresponding temperature T 0 ≈ 0.6 from the RLFA solution (see Fig. 1).The case of L = 10 is very similar.We verified that the linear dependence of barriers height on temperature and their disappearance at a temperature close to T 0 , obtained within RLFA, are also valid for other values of J 2 .Note that the linear dependence on T at low temperature follows directly from the definition of the free energy F = E − T S, where E is the energy and S is the entropy.4) for several values of J 2 (listed in the legend and valid for both plots) using the MC method for L = 14, (a) for periodic boundary conditions, (b) for free boundary conditions.The curves for J 2 = 0.3 and 0.5 are shifted down by 75 and 120, respectively.The solid lines provide guides to the eye.

B.2 Periodic vs free boundary conditions and larger samples
Changing the boundary conditions from periodic to free does not affect the restricted free energy at a qualitative level, but there are some differences, see Fig. 15.In the latter case, the spin flip energy on a flat domain interface near the edge of the sample is 2|J 1 |, which is half as much.This determines the energy barrier for metastable states in the slab phase and the decrease in energy when the entire side of the droplet dissolves in the droplet phase.At J 2 = 0, the sequence of configurations that determines the restricted free energy at zero temperature remains the same (Fig. 4), but their degeneracy is less due to the lack of translational symmetry (two domains are always at the edge of the sample).The slab-droplet transition is present at the same value of M , but the dependence on J 2 is different.As a consequence, there are more local minima at J 2 = 1/2.The fully polarized state always represents the global energy minimum at zero temperature, even for J 2 = 1/2, when the energy gap is 2|J 1 |.This is due to the difference in the energy contributions of spins in the corners of domains compared to bulk spins.

B.3 Dependence on nearest-neighbor correlations for larger samples
As in the case of L = 6, we plot the energy E(M , c 1 , c 2 ) and entropy S(M , c 1 , c 2 ) as functions of the total spin M and the nearest and next-nearest correlations (c 1 and c 2 ), which are cal- They form a quasi-hypersurface in the states space, including only states with the lowest energies.

Figure 2 :
Figure 2: The Landau free energy F 0 (m) calculated within RLFA (solid lines) for several values of J 2 at temperature T = 0.1, which shows a kink at polarization around m = 0.5.The dashed lines are the Landau free energies according to MFA.

Figure 3 :
Figure 3: (a) Restricted free energy F calculated exactly by Eq. (4) at J 2 = 0.3 for sample size L = 6 as a function of total spin M for a set of temperatures T below and above the phase transition.(b)The same restricted free energy per spin F 0 /N as a function of total spin M and temperature T .The free energy is defined only for integer even values of M and linearly interpolated between them.For ease of comparison with the RLFA results in Fig.1, we set here F = 0 at M = 0 at each temperature.Dark blue dots correspond to the global minimum of F 0 /N at each temperature.

Figure 4 :
Figure 4: Configurations, denoted M n , that contribute to the restricted free energy F (M ) at zero temperature, i.e. have minimal energies for a given M , for all possible values of the total spin M at L = 6; n is their degeneracy due to translational and rotational symmetry.The configurations that change for J 2 > 1/4 and J 2 > 1/3 are shown separately in the bottom row (they affect the dependence of the energy barrier ∆F M on J 2 as shown in Appendix B).

Figure 5 :
Figure 5: (a) Energy E as a color map for all possible states with total spin M and nearest-neighbor correlations c 1 and c 2 , (b) entropy, S = log n, of these states.The restricted free energy is F = E − T S, where T is the temperature.Large dark blue dots in (a) correspond to metastable states and one stable state with M = 36.

Figure 6 :
Figure 6: Configurations of metastable states, denoted M n c 1 ,c 2

4J 1 LFigure 7 :
Figure 7: Restricted free energy F as a function of total spin M for J 2 = 0.3 calculated by the MC method at several temperatures T below and above the phase transition.The sample size is L = 10.

Figure 8 :
Figure 8: (a) Restricted free energy F 0 (m) within RLFA for J 2 = 0.3 and polarization limited by m ∈ (0, 0.5) to show the appearance at low temperature of a barrier at m ̸ = 0 whose height first increases and then decreases as the temperature approaches zero.(b) The barrier height ∆F bar = F bar − F (0)for the metastable state at m = 0, which appears in the restricted free energy F (m) calculated within RLFA, as a function of temperature T .Only J 2 < 0.31 are considered, which implies T 0 < T 1 (for temperatures definition see Fig.1a), since these two temperatures become undefined at larger J 2[29].

6 - 1 Figure 9 :
Figure 9: The right-hand side of the RLFA equation (2) for J 2 = 0 and T = 0.1 (solid dark blue) as a function of the external field h and its left-hand size (solid purple), which is just the polarization m.The intersection of two curves for a given value of m yields the dependence h(m).Dashed lines correspond to MFA.

Figure 10 :Figure 11 :
Figure10: In the Ising model (J 1 = −1, J 2 = 0), the local field h s = − 〈 j〉 J i j s j at the central spin position (in light blue) is determined by its four nearest neighbors with spin values of +1 (white) or −1 (dark blue).The value of the local field and the number of configurations n due to rotational symmetry is indicated above each configuration as h s n .h s changes its sign when spins are flipped.Thus, the local field can have five different values when J 2 = 0.

L = 6 ,Figure 12 :Figure 13 :
Figure12: Restricted free energy F as a function of total spin M at zero temperature T = 0, calculated according to Eq. (4) for several values of J 2 (listed in the legend and valid for both plots), (a) exactly for the sample size L = 6, (b) using the MC method for L = 10.The solid lines provide guides to the eye.

L = 6 , J 2 =Figure 14 :L = 14 ,
Figure 14: Restricted free energy barrier height, ∆F M = F M − F M −2 , as a function of temperature for several values of the total spin M and J 2 = 0.3.(a) The sample size is L = 6.(b) L = 10.The three upper curves correspond to M = 2, 22 and 42, and in the middle is M = 62.

Figure 15 :
Figure15: Restricted free energy F as a function of total spin M at zero temperature T = 0, calculated according to Eq. (4) for several values of J 2 (listed in the legend and valid for both plots) using the MC method for L = 14, (a) for periodic boundary conditions, (b) for free boundary conditions.The curves for J 2 = 0.3 and 0.5 are shifted down by 75 and 120, respectively.The solid lines provide guides to the eye.

Figure 16 :
Figure 16: Possible states with the total spin M and nearest-neighbor correlations c 1 and c 2 , and (a) their energy E as a color map and (b) entropy S = log n for L = 10.The restricted free energy is F = E − T S, where T is the temperature.Large dark blue dots correspond to metastable states and one stable state with M = 100.