Unveiling Phase Transitions In 1D Systems With Short-Range Interactions

do The statement that any phase transition is related to the appearance or disappearance of long-range spatial correlations precludes a ﬁnite transition temperature in one-dimensional (1D) systems. In this paper we demonstrate that the 1D Ising model with short-range exchange interactions exhibits a second-order phase transition at a ﬁnite temperature relying on the proper choice of the order parameter. To accomplish this, we combined analytical calculations and high-precision entropic sampling simulations and chose a slightly diﬀerent order parameter, namely the module of the magnetization. Notably, we detected a phase transition with a corresponding critical temperature around 15 K, which is in excellent agreement with experimental results. Our study indicates that an inappropriate choice of the order parameter may mask phase transitions in one-dimensional systems.


I. INTRODUCTION
The impossibility of a phase transition in (d ≤ 2)systems with continuous symmetries and short-range interactions has been established after the publication of the Mermin-Wagner theorem [1]. Particularly, almost all statistical mechanics textbooks, supported by the Perron-Frobenius theorem, prohibits a finite temperature phase transition in one-dimensional (1D) systems with discrete symmetries [2]. Nonetheless, many theoretical works have pursuit such transition under different arguments. The development of heteroepitaxy techniques ushered in the manufacturing of continuous, extremely low dimensional metal matrices on non magnetic substrates [3][4][5], allowing experimental characterization and theoretical studies of 1D magnetic systems [6]. For instance, the experimental observation of an 1D chain of Co atoms deposited on a Pt surface enduring a ferromagnetic phase at a temperature around 15 K [7] seems to contradict classical theoretical predictions whether there are true thermodynamic phase transitions in 1D systems with short range interactions [8]. Some attempts to explain this apparent discrepancy have been made. Curilef et al. investigated the critical temperature of the Ising model including a long-range term through a power law that decays over large interparticle distances [9]. However they made no comparisons to the results of Ref. [7]. In its turn Li and Liu [10] proposed a Heisenberg model with a large magnetic anisotropy and an additional external magnetic field to describe a monatomic chain of Co. Afterwards, Vindigni et al. [6] developed a combined experimental and theoretical study to explain the behavior of the magnetization of 1D magnetic systems, using the Heisenberg model with an anisotropic term, solved by the transfer matrix technique. Recently a study was performed on the magnetic properties of an 1D Au-Co chain on a copper surface using the density functional theory and kinetic Monte Carlo (KMC) [11]. They performed a KMC simulation of a 1D chain of Heisenberg anisotropic magnets taking the module of the magnetization as the order parameter. They found a phase transition at non zero temperature. The authors argued that this does not violate the Mermin-Wagner theorem due to the nonzero anisotropic term.
The foregoing works treated an 1D magnetic chain using a Heisenberg Hamiltonian with an additional anisotropic or a long-range term. Notwithstanding the first approximation for an 1D magnetic system should be the Ising model, for which we have now a new scenario. This model was proposed by Wilhelm Lenz to the doctoral work of his student Ernst Ising to describe the behavior of a magnetic material in one dimension. Their solution showed that at zero magnetic field a second order phase transition occurs only at T = 0. However, Lars Onsager solved analytically the two-dimensional version of this model and encountered a phase transition at a finite temperature [12]. Since then, the Ising model became a cornerstone to condensed matter physics and science in general [13][14][15][16][17].
In this paper, in contrast with the common understanding, we reveal a phase transition at finite temperature in the 1D Ising model. To drive the reasoning, well-established theories should be revisited. In 1966, Kadanoff determined the scaling laws for a second order phase transition of the Ising model [18] based on a Widom's argument [19] on the consequences of considering the free energy as a homogeneous function of parameters describing the proximity to the critical point, such as ϵ = (T − T c )/T c . He assumed that a large lattice can be divided into smaller ones of still large sizes, but much smaller than the system correlation length. The zeroth order approximation is based on the assumption that the free energy of an isolated system of side L is an analytical function for ϵ, but singular for the system size. Kadanoff considered two possibilities for the specific heat divergence which is given by ϵ −α . The first one assumes the exponent α ̸ = 0 and the second α = 0, the latter representing the 2D Ising model since the results obtained with α = 0 corroborate the experiments [18]. Hence, we can state that in a second-order phase transition the specific heat is an analytical function of temperature and singularities occur as the system size increases.
Another important issue is the Landau's theory of second-order phase transitions [8], which is based on a continuous transition from an ordered state to a disordered one. A second-order phase transition may be related to some symmetry property, such as the Curie point where a ferromagnetic substance passes into the paramagnetic phase. To describe this event, Landau defined an order parameter as a quantity that assumes nonzero values (positive or negative) for an ordered phase and zero for the disordered one. Furthermore, he proposed that the singularity is observed when the thermodynamic potential is expanded in terms of the order parameter.
This paper is organized as follows. Section II describes the chose of the order parameter. Section III calculates the partition function using the counting for small lattices and serie analysis. Section IV describes the entropic sampling simulations. Section V combines the analytical and simulational results. Finally, section VI presents our conclusions.

II. THE ORDER PARAMETER
The Hamiltonian of a one-dimensional Ising model with periodic boundary conditions when an external field is present is such that: where σ i = ±1, J is the coupling constant and H the external magnetic field. Taking the magnetization as the order parameter and H ̸ = 0 its value per site in the ground state configuration is conditioned to the field sign, being m + = 1 when the field is positive and m − = −1 when negative. If H = 0 there is a singularity where both phases have the same probability. This situation is the same as described by Landau for a first order transition, i.e. an equilibrium between two phases. Calculating the mean magnetization at T = 0 we have where P {↑↑↑} (P {↓↓↓} ) denotes the probability of finding a configuration with all spins up (down). Substituting the values of m + and m − we obtain ⟨m⟩ = 0. At high temperatures the most probable configurations are those where the energies match k B T and in the limit T → ∞ they are configurations where H L = 0, the maximum energy. For these configurations the magnetization is zero since the number of positive and negative spins is the same. Hence, if H = 0, the magnetization does not behave as an order parameter in the study of a second order phase transition in the 1D Ising model, since it does not distinguish an ordered configuration from a disordered one.
Since the magnetization cannot be used as the order parameter for H = 0 the corresponding correlation length is not a proper criterion for validate the existence of a phase transition. To understand the reason, let us take a look at the relation between the correlation length and the spin-spin correlation function.
The spin-spin correlation function is defined, for lattice systems, as Hence, for a configuration in the ordered phase, g(⃗ r) must be zero. Let us exam its value for the 1D Ising model when H = 0 and T = 0. The term ⟨σ i ⟩ = ⟨σ j ⟩ = m is null due to the existence of two antisymmetric ordered configurations in the ground state, matching the Eq. 2. This statement is based on the fact that the probabilities of finding the system in an ordered configuration with spins up or down are equal. Clearly, both configurations are not present simultaneously in real systems. However, in the Transfer Matrix calculation one considers the superposition of both configurations, yielding a zero average total magnetization. This result is unrealistic since in the thermodynamic limit the system phase is determined by the initial conditions and only one ordered configuration is allowed. Several textbooks showing the solution for the 1D Ising model consider its experimental realization where two antisymmetric solutions are possible, allowing the correlation function to be zero at T = 0 [2,[20][21][22]. However, they use the Transfer Matrix result for T > 0. The bottom line is simply whether or not ⟨σ i ⟩ is zero for all T . Textbooks show that the average magnetization is zero for any temperature T ̸ = 0 and, using a mathematical trick, promptly prove that it assumes the value +1 or −1 for T = 0 [2,21,22]. This clearly suggest that for an extremely low temperature, say T = 10 −23 K, the magnetization is null, however for T = 0 it will be either +1 or −1. This is certainly a non-physical situation. The entire scientific community has accepted this assertion through decades. Hence, it will be very difficult to convince anyone that this is an inconsistent strategy. Definitely, this mathematical maneuver is tailor-made to conclude that ⟨σ i ⟩⟨σ j ⟩ = 1, allowing to achieve T c = 0.
In its turn, the term ⟨σ i σ j ⟩ is always positive. Thus, g(⃗ r) = ⟨σ i σ j ⟩ is always greater than zero for the Ising model when T = 0 contravening the requirement of g(⃗ r) = 0 for an ordered configuration. This result is acceptable when we use the results from the Transfer Matrix for any value of T .
The correlation function is, therefore, not precisely defined for H = 0 and its value, obtained from the Transfer Matrix, seems to be not appropriate. Another important aspect is that the definition of the correlation length is similar to the magnetic susceptibility, hence if the latter diverges, the former also diverges. It is related to the magnetic susceptibility as where χ = lim H→0 ( ∂M ∂H ). Otherwise, if we choose as the order parameter the module of the magnetization per site such that |m + | = |m − | = 1, we will obtain at T = 0 ⟨|m|⟩ = P {↑↑↑} + P {↓↓↓} = 1, whereas for an disordered configuration we get ⟨|m|⟩ = 0. Thus, the order parameter that indeed distinguishes an ordered phase from a disordered one is the module of the magnetization. This is a well-known fact for the 2D case and was used recently to explore a phase transition in a 1D Heinsenberg chain of spins [11].

III. PARTITION FUNCTION
If we adopt the module of the magnetization as the order parameter and write the Hamiltonian in terms of this order parameter then we can calculate a partition functionZ usingH such that when H = 0 it will yield reliable results. In this case the solution via transfer matrix becomes unfeasible, because the module in the second term is not factorable.
Since an exact solution for an infinite lattice in this case is not available, what can be done is the direct calculation of the partition function for small sizes. Defining E j and M j as where the index j denotes a particular configuration, we can writeH = −JE j − H|M j | and the partition function may be written as wherẽ In Fig. 1 we depict the possible configurations for a L = 4 lattice, from which we can write the partition function for this lattice size as Z 4 =g 4 (4, 4)e 4βH e 4βJ +g 4 (0, 2)e 2βH +g 4 (0, 0) +g 4 (−4, 0)e −4βJ =2e 4βH e 4βJ + 8e 2βH + 4 + 2e −4βJ , whereg 4 (E, M ) is the number of configurations with dimensionless interaction energy E and order parameter M . One can see thatg 4 (4, 4) = 2,g 4 (0, 2) = 8, g 4 (0, 0) = 4 andg 4 (−4, 0) = 2, where the first term corresponds to the ordered states, the second one to the configurations with one unpaired spin and the remaining to configurations with two unpaired spins, in which case they may be neighbors or not, corresponding to the termsg 4 (0, 0) andg 4 (−4, 0), respectively.
Using the same counting process for L = 8 we obtaiñ

IV. SIMULATIONAL DETAILS
Entropic sampling simulations are based on the Wang-Landau method [23] where we introduced some changes to improve accuracy and help saving CPU time [24][25][26]. We halt the simulations when the sixteenth Wang-Landau level (f 15 ) becomes flat. Once the joint density of states is available one can obtain the canonical averages as k, l = 0, 1, 2, ..., which are present in the estimate of thermodynamic quantities as mean energy, magnetization, heat capacity, susceptibility and cumulants. The joint density of states for a lattice of size L can be used in Eq. (11). Also, setting k = 0, l = 1 and taking |M | before g L (E, M ), we obtain via simulations the same magnetizations for these lattice sizes for H = 0 as

V. RESULTS
In the present work we restrict ourselves to magnetic spin models, nonetheless the proper choice of the order parameter may trigger important features in a variety of systems. This was evinced recently with the observation of the coexistence of first-and second-order phase transition in a correlated oxide [27]. Another challenge is the study of the Bell-Lavis lattice model to water where the choice of the order parameter is crucial to unveil the metastable behavior of the system [28,29]. It is worth noting that the Perron-Frobenius theorem does not apply to this case, since its demonstration assumes that the partition function can be factored.
Any thermodynamical quantity can be obtained once one achieves the partition function. Here we have computed the average value of the module of the magnetiza-tionM L , since this is the appropriate quantity to observe a phase transition in our model. We can easily getM L using the Helmholtz free energyF L = −k B T lnZ L , with M L = lim H→0 (− ∂F L ∂H ), yielding and Inspecting the numerators and denominators of the two quantities above we see that they can be described by the following series where ⌊⌋ denotes the largest integer not greater than the argument and L should be a power of 2.
In order to confirm these predictions, we carried out entropic sampling simulations for L = 4, 8, 256, 512, and 1024 to construct the joint density of states. In Fig.2 we display the analytical results and those of the simulations for these lattice sizes, assuming J = 1. One can see that they match with great accuracy. It is also clear that the system displays finite-size effects, suggesting a finite-size scaling study. It is important to note that for lattice sizes greater than 8 the simulations corroborate our conjectures on the series pointed out above. The analytical expression of the partition function for L > 8 and not a power of 2 is unattainable. Nonetheless, the simulations are free to go ahead to as large and diverse lattices as necessary.
We can also obtain a susceptibility of the order parameter asχ L = lim H→0 ( ∂M L ∂H ) yielding , and .
Notice that the second term is exactly the square of the average magnetization. Then, the first term should be ⟨M 2 L ⟩. The denominator is equal to the magnetization and the numerator can be written as These series for the numerator and denominator can be found in Refs. [30][31][32].
In Fig. 3, we show the behavior of the magnetic susceptibility considering both the standard and the new order parameter m ′ = | ∑ i σ i | obtained from simulation and the Transfer Matrix technique for the 1D Ising model.
One can see that the simulational results using the magnetization as the order parameter corroborate the Transfer Matrix ones, although an FSS behavior is noticeable. Also, the susceptibility shows a peak around the region where the analytical value diverges when we use the absolute value of the magnetization as the order parameter. It is important to state that both results were taken from the same simulational procedure, the only difference being the order parameter used.
The same behavior appears in the 2D Ising model, as one can see in Fig. 4, where is shown the magnetic susceptibility for square lattices of sizes L = 16 and 32. One can see that the divergence of the susceptibility starts near the critical temperature, even though its divergence continues for T < T c , up to T = 0. This can be considered as an FSS effect similar to the 1D case. Notice also that without a maximum in the susceptibility, to determine the critical temperature is unfeasible. Choosing the order parameter as m ′ = | ∑ i σ i | the divergence of the susceptibility is associated with the limit L → ∞ as suggested by Kadanoff [18].
Since the simulation and analytical results agree perfectly, we will use the simulation procedure to study the finite-size behavior. Following Refs. [33,34] we can define a set of thermodynamic quantities related to logarithmic derivatives of the magnetization that scales as for j = 1, 2, ..., 6, where t = (T − T c )/T c is the reduced temperature. These cumulants allow to estimate the critical exponent ν before determining the critical temperature. We performed simulations for L = 256, 320, 384, 448, 512, 640, 768, 896, 1024, and 1152 with n = 24,24,20,20,20,16,16,16,12, and 12 independent runs for each size, respectively. In Fig. 5 we depict the behavior of these quantities, obtaining 1/ν = 0.34159(72) as the mean value of the slopes of the six best fitting straight lines, yielding ν = 1/(1/ν) ± 1/(1/ν) 2 ∆(1/ν) = 2.9274 (33). The exponent ν was obtained in studies of 1D Ising model with long-range interactions, as J ij = J 0 /|i − j| 1+σ , such that the larger is σ the closer the system is of the short-range interactions. In Ref. [35] they obtained ν = 2.602721 for 10 interacting spins and σ = 1.0, while in Ref. [36] the estimate using Monte Carlo histogram method and σ = 1.0 is ν = 2.42 (1). These results suggest that our ν consists in a short-range limit for the long-range interactions systems.
where a is a constant, to determine T c for an infinite lattice. This relation is valid for the quantities χ, D K2 , D K3 , D K4 , and U 4 , defined in Refs. [33,34]. Theses results are shown in the Fig. 6. For an infinite lattice (L → ∞) we obtain five estimates for the critical temperature, the mean value yielding T c = 0.171008(63). This finite value for the critical temperature contradicts the Transfer Matrix result. There is only one way to support this result is to find a one-dimensional magnetic material that bears a phase transition at any finite temperature. The findings of a finite-temperature phase transition for 1D systems with short-range interactions are not by any means a surprise in the scientific community, as can be seen in the investigation of domain walls in Fe nanostripes [37] or the study of a ferromagnetic order in a line of cobalt atoms [7]. Thereon, several theoretical works tried to explain the phenomenum of ferromagnetism in one-dimension using the Heisenberg spin to model the systems [6,10,11]. An excellent candidate to design a one-dimensional Ising chain is cobalt due to the strong anisotropic field [38]. A study about the compound Ca 3 Co 2 O 6 , using orbital image technique, showed that the high anisotropic field freezes the spin in the easy axis direction [39]. This behavior endorses the comparison between cobalt chains and the 1D Ising model, causing the use of the Heisenberg model not strictly necessary. Therefore, we must look more closely at Gambardella's [7] work.
To investigate the cobalt chains Gambardella et al. used the X-ray magnetic circular dichroism (XMCD) technique and found that the magnetic spin moment is m s ∼ 2.03µ B and the magnetic anisotropic energy per atom is 2.0 ± 0.2 meV , which is an extremely high value when compared to the bulk one. These values confirm the results obtained in Ref. [39]. Gambardella et al. also obtained a blocking temperature of T B = 15 ± 5 K, for which the ferromagnetic order disappear. Therefore, we will use 2J = 15 meV to compare our results with the one obtained by Gambardella et al. We then obtained our critical temperature as T c = 0.171008(63)J/k B = 14.8835(55) K, where k B = 8.6173324(78) × 10 −5 eV/K, in excellent agreement with the experiment.
It is worth noting that in this work we have analyzed only the discontinuity on the order parameter since quantities such as specific heat show no finite size behavior, in contrast with 2D systems. In Fig. 7 we show how the simulation result for L = 1024 lattice fits that obtained analytically by transfer matrix technique. This behavior is known as Schottky effect, first observed in two-state systems, and in general is not accepted as a phase transition signature. More complex systems also exhibit this effect and a complete understanding of its origin is still lacking [40]. We believe that the appearance of this effect in our model is due to the fact that the first term in Eq. 6 is still factorable when computing the partition function, what does not happen in higher dimensions.

VI. CONCLUSIONS
In summary, we see that the choice of the order parameter for the 1D Ising model as the module of the magnetization and not simply as the magnetization, whose average does not distinguish the ordered state from the disordered one, reveals the existence of a second-order phase transition with a finite and well-defined critical temperature, in contrast to the belief that such a transition will never be exhibited, as stated by almost all statistical mechanics textbooks and supported by the Perron-Frobenius theorem. However, as explained above, we see that it does not apply to the present case, once the factorization of the partition function becomes impossible when we take the module of the magnetization as the order parameter instead of the magnetization itself.
Our results open up a wide range of opportunities to understand phase transitions in one dimension and may have potential ramifications and applications, leading subsequently to possible advances in a much broader context.