Nonequivalent Antiferromagnetically Coupled Sublattices Induce Two-Step Spin-Crossover Transitions: Equilibrium and Nonequilibrium Aspects

: As a continuation to the previously published work (Yalçın et al. (2022)), we investigate the equilibrium and nonequilibrium properties of the spin-crossover systems, with a specific focus on the nonequivalent sublattice, and compare these properties with those of the equivalent sublattices. We used the lowest approximation of the cluster variation method (LACVM) to derive the static equations for the order parameters of the two sublattices and determine high-spin fraction in relation to temperature and external magnetic field in a spin-crossover system. At a low temperature, the transition from stable high-spin (HS) state where n HS = 1 occurs in the plateau region, where n HS = 0.5 for nonequivalent sublattices. The order parameters for non-equivalent sublattices exhibit different states at the transition temperature. Also, we study the nonequilibrium properties of the order parameters and high-spin fraction using the path probability method (PPM). With the current model, we obtain and analyze the relaxation curves for the order parameters S a , S b , and high-spin fraction. These curves demonstrate the existence of bistability at low temperatures. At the end of this study, we present the flow diagram that shows the order parameters for different temperature values. The diagram exhibits states that are stable, metastable, and unstable.


Introduction
Spin-crossover (SCO) systems have attracted significant attention because of their potential application to memory devices, molecular sensors, and displays, as well as to recent developments related to the barocaloric effect and refrigerant effect [1][2][3][4][5].The SCO phenomenon refers to the transition between two molecular spin states: the low-spin state (LS), which is diamagnetic, and the high-spin state (HS), which is paramagnetic.The LS state is energetically favorable and dominates at low temperatures, whereas the HS state dominates at high temperatures due to its favorable energy.This transition specifically occurs in first-row transition metal complexes with an electronic configuration ranging from d 4 to d 7 .
In many investigations, it has been observed that, when spin-crossover systems are exposed to light in the green to red wavelength range, they complete a spin transition from a low spin state to a high spin state, even at low temperatures.The process that has been described is referred to as light-induced excited spin state trapping [21][22][23][24][25][26][27].Furthermore, both experimental and theoretical studies have documented the presence of light-induced hidden hysteresis and bistability in spin transitions, as reported in the references [28].In addition, the application of pressure plays an important role as an external stimulus [29,30].Some SCO systems change from first-order transitions to continued transitions or incomplete transitions [31].
In this work, we extend our previously published work [49] with a focus on nonequivalent sublattices.We study the equilibrium and nonequilibrium properties of the system using the lowest approximation of the cluster variation method (LACVM) [50] and the path probability method (PPM) [51].The Ising-like Hamiltonian model that we have analyzed for SCO systems considers both equivalent and nonequivalent sublattices.We first investigate the thermodynamic properties for different values of the intersublattice and intrasublattice interaction parameters, ligand field energy, and degeneracy.Additionally, we present the magnetic evolution of the order parameters and the high-spin fraction at various fixed temperatures.For nonequivalent sublattices in the warming process, we found that the high-spin fraction n HS displayed a stable HS state in the low temperature region and jumped to the intermediate metastable/stable state HS-LS at n HS = 0.5.The relaxation curves of the order parameters and high-spin fraction with bistable states were observed.For the same initial values of the order parameters, two-step relaxation curves are present.The flow diagram of the order parameters at T = 60 for nonequivalent sublattices shows two stable states, two unstable states, and metastable states.Our results are in agreement with other experimental and theoretical findings.The self-consistent equation of the system was performed using the Newton-Raphson and Runge-Kutta methods.
The present work is organized as follows: Section 2 describes the model briefly and presents the derivations of static and dynamic equations.In Section 3, we obtain the numerical results and offer a discussion.In Section 4, we end with a conclusion.

Description of the Model and Derivations of Static and Dynamic Equations
In order to describe the two-step spin-crossover systems, an Ising-like Hamiltonian with ferromagnetic intrasublattice interactions (J a > 0, J b > 0) and antiferromagnetic intersublattice interactions (J ab < 0) was introduced in [32,52] as follows: where S a i = S a j = ±1, S b i = S b j = ±1, h = ∆ 2 − k B T 2 ln g (∆ is the ligand field energy, g = g HS /g LS is the degeneracy, T is the absolute temperature, and k B (= 1) the Boltzmann constant).All exchange coupling parameters appearing in Equation (1) are in units of k B T. In the presence of an external magnetic field H, the Zeeman effect is considered (H = gµB, µ is the Bohr magneton and B the applied magnetic field [39]).Generally, g depends on temperature through the temperature dependence of the phonon density of the lattice [39,40].Sublattice order parameters are defined as: where X a 1 and X b 1 are average fractional values of each spin state with value +1 on a and b sublattices, respectively, while X a 2 and X b 2 are average fractional value of each spin state with value −1 on a and b sublattices, respectively.By solving Equation (1), Helmholtz free energy can be expressed in terms of sublattice order parameters defined as follows: where β = 1/k B T, F F = E − TS entropy ) is the free energy and N (N = N a + N b ) presents the total number of points on sublattices.After using the LACVM [50], one can obtain the following self-consistent equations for two sublattices to investigate the static properties: After some calculations from Equation ( 4), and after applying the two possible values of the spin states [53], we can find the following: with The order parameters S a and S b for each sublattice are directly related to the high-spin fraction n HS of a spin-crossover system defined as: where S = S a +S b 2 .In order to study the dynamic aspects of the above spin system, one can obtain dynamic equations by using the path probability method introduced by Kikuchi [51].In the PPM, the rate of change of a state variable is written as: where X ij is the path probability rate for the system to go from state i to state j, with For the rates X ij , we use the following recipe, called recipe II, introduced by Kikuchi: where k ij is the rate constants, with k ij = k ji , Z is the partition function, and E is the internal energy of the system.The partition function Z is calculated for the Hamiltonian model (1), and, for the two sublattices, it can be written as: where, after many calculations, we found that e a 1 = e β(−2J a S a −J ab S b +h+H) and e a 2 = e β(2J a S a +J ab S b −h−H) for sublattice a and e b 1 = e β(−2J b S b −J ab S a +h+H) and e b 2 = e β(2J b S b +J ab S a −h−H) for sublattice b.The dynamic equations for two sublattices can be written as: or From Equation ( 12), after many calculations, and after applying the two possible values of the spin states as in the static part, we can obtain the dynamic equation for the order parameters: where k 1 is the rate constant.If we take dS a dt = 0 and dS b dt = 0, the above equations return Equation (4) that present static equations.

Results and Discussion
In this section, we computed results using the static and dynamic equations that were introduced in the previous section.First, we are going to investigate how the intersublattice parameter, intrasublattice parameter, degeneracy, ligand field energy, and external magnetic field influence the phase transition properties of a spin-crossover system.In the analysis of equilibrium properties, we examine both the equivalent and nonequivalent sublattices.Second, we compute results for the nonequilibrium properties of the system.

Static Behaviors
Figure 1 shows the temperature evolution of the order parameters (S a , S b ) and highspin fraction (n HS ) for different selected values of the intersublattice interaction parameter (J ab ) in the case of equivalent (J a = J b ) and nonequivalent (J a ̸ = J b ) sublattices.For equivalent sublattices in Figure 1a, the order parameter transition was gradual when the intersublattice interaction parameter (J ab ) was weak.As the value of J ab increased, the interaction between the two sublattices became stronger, leading to a first-order transition characterized by a hysteresis of the order parameters.This property was observed when the reduced temperature was increased and when it was decreased.While increasing the coupling parameter (J ab ) values, the width of the hysteresis became enlarged.The thermal behavior of the high-spin fraction n HS depends on the J ab , plotted in Figure 1b.For small values of the coupling parameter J ab = −400 (orange curve) and −300 (green curve), the SCO transition was gradual.The thermal properties of the high-spin fraction n HS underwent a change in shape as J ab increased, i.e., J ab = −200 (magenta curve), −100 (red curve), or −50 (blue curve), and the transition of the spin-crossover (SCO) exhibited thermal hysteresis.By increasing the coupling parameter J ab , the width of the hysteresis became expanded, as demonstrated by the simulations performed for the order parameters (S a , S b ).Next, we studied the thermal properties of the order parameters  ,  , and high spin fraction  without an external magnetic field for various chosen values of in trasublattice parameters (  ,  ).In Figure 2a,b,  =  are equivalent sublattices whereas, in Figure 2c,d,  ≠  are non-equivalent sublattices.When the coupling pa rameters were small, such as  =  = 15 and 50 for equivalent sublattices and  = 15 and  = 14 for nonequivalent sublattices, the transition was gradual and hysteresi was not observed.When the intrasublattice parameter increased to values such as  =  = 80 and 100, a first-order transition occurred, leading to hysteresis during the warm ing and cooling processes.This was observed for the order parameters  ,  , and for th high-spin fraction  .When the sublattices were nonequivalent, with  = 50 and  = 49 (magent curve), the order parameter transitions were first-order, resulting in hysteresis between the warming and cooling processes.The high-spin fraction  transition occurred in two steps and exhibited hysteresis between the warming and cooling processes.Increas ing the  and  values resulted in a first-order transition that shifted to lower temper atures, as depicted in Figure 2. The intrasublattice coupling parameter has been addressed in different theoretical models [52,55].The additional parameters utilized in our theoreti cal model are:  = −15, ∆= 400, g = 100 and  = 0.0.Figure 1c,d display the phase transitions of the order parameters and the high-spin fraction n HS for nonequivalent sublattices.The sublattices had different values of J a = 200 and J b = 199.These phase transitions are shown for specific values of the intersublattice interaction parameter (J ab ).We analyzed the phase transition of order parameters and the high-spin fraction for J ab = −50 (indicated by the blue curve) in order to explain the intermediate plateau region for the high-spin fraction n HS .At low temperatures, the order parameters of two sublattices, S a and S b , were both equal to −1, indicating a lowspin state.As the temperature increased, the order parameters remained the same until a bifurcation occurred at T ↑ ∼ 177.9.At this point, the order parameter of sublattice b underwent a first-order phase transition from S b = −0.78(low-spin LS state) to S b = 1 (high-spin HS state).The order parameter of sublattice a remained in a low-spin state until the temperature reached approximately T ↑ ∼ 202.7, at which point it underwent a first-order phase transition from S a = −0.72 (LS state) to S a = 1 (HS state).At temperatures above T ↑ ∼ 202.7, both sublattices became equivalent again.During the cooling process, a comparable situation arose where the breaking of the symmetry between the two sublattices occurred at a temperature of T ↓ ∼ 90.2.This was when the order parameter of sublattice b underwent a first-order phase transition from S b = 0.9 (HS state) to S b = −1 (LS state).The order parameter of sublattice a remained in a high-spin state until the temperature reached around T ↓ ∼ 10.6.At this temperature, it underwent a first-order phase transition from S a = 1 (HS state) to S a = −1 (LS state).At temperatures below T ↓ ∼ 10.6, both sublattices exhibited the same state, which is referred to as an LS state.When the J ab values decreased, we observed that the order parameters S a and S b exhibited similar properties for the warming process.It is important to emphasize that the value of S a (=1) remained constant during the cooling process, and no phase transition was observed for J ab values of −200 (magenta curve) and −400 (orange curve) from S a = 1 to S a = −1 as it appeared for S b .
The high-spin fraction n HS for the nonequivalent sublattices (Figure 1d) exhibited an abrupt (first-order phase transition) behavior, with two steps and separate branches in the plateau region between the hysteresis loop during both the warming and cooling processes.This phenomenon is also supported by the experimental and theoretical study of Fe II (5 − NO 2 − sal − N(1, 4, 7, 10)) [54].The properties of n HS for J ab = −50 (blue curve) arose due to the presence of various combinations of spin states on the plateau (n HS = 1 2 ).The order parameters (S a , S b ) exhibited different spin states during the heating process, from T 1/2 ∼177.9 to ∼202.7, as discussed in Figure 1c.During the cooling process, the high-spin fraction n HS exhibited similar properties.There were different combinations of spin states on the plateau (n HS = 1 2 ), with values ranging from T 1/2 ∼ 90.2 to ∼10.6 .While J ab values decreased, the transition became two-stepped with the presence of hysteresis.The width of the hysteresis loop in the second step was slightly reduced with decreasing J ab (−200 and −400), and the whole cycle was shifted to relatively high temperatures.At very low temperatures in cooling mode, an incomplete spin transition was present on the plateau (n HS = 1 2 ) for both J ab values.The selected value parameters are presented in the figure caption.
Next, we studied the thermal properties of the order parameters S a , S b , and high-spin fraction n HS without an external magnetic field for various chosen values of intrasublattice parameters (J a , J b ).In Figure 2a,b, J a = J b are equivalent sublattices, whereas, in Figure 2c,d, J a ̸ = J b are non-equivalent sublattices.When the coupling parameters were small, such as J a = J b = 15 and 50 for equivalent sublattices and J a = 15 and J b = 14 for nonequivalent sublattices, the transition was gradual and hysteresis was not observed.When the intrasublattice parameter increased to values such as J a = J b = 80 and 100, a firstorder transition occurred, leading to hysteresis during the warming and cooling processes.This was observed for the order parameters S a , S b , and for the high-spin fraction n HS .
When the sublattices were nonequivalent, with J a = 50 and J b = 49 (magenta curve), the order parameter transitions were first-order, resulting in hysteresis between the warming and cooling processes.The high-spin fraction n HS transition occurred in two steps and exhibited hysteresis between the warming and cooling processes.Increasing the J a and J b values resulted in a first-order transition that shifted to lower temperatures, as depicted in Figure 2. The intrasublattice coupling parameter has been addressed in different theoretical models [52,55].The additional parameters utilized in our theoretical model are: J ab = −15, ∆ = 400, g = 100 and H = 0.0.
Figure 3 shows the temperature dependence of the S a , S b , and n HS for four different values of degeneracy in the absence of an external magnetic field.As it is shown in Figure 3a,b, in the case of equivalent sublattices (J a = J b = 50), the transition was first-order for large values of g (1000/magenta curve and 5000/green curve) and the whole hysteresis shifted to a lower temperature region.When the degeneracy became g (100/red curve and 10/blue curve), the hysteresis disappeared and the transition became gradual and shifted to the higher temperature side.We can conclude that, for small values of g, the transition is gradual and shifted to the high temperature side.
In the case of nonequivalent sublattices, when the intrasublattice parameters were J a = 51 and J b = 50, the transition was first-order for S a , S b , while, for n HS , it was firstorder with two steps, as depicted in Figure 3c,d.The high-spin fraction n HS shifted to the left (lower temperature) as the degeneracy g increased.For non-equivalent sublattices, we found three typical two-step transitions during the warming and cooling processes for the high-spin fraction for the selected values: (I) first-order transition with the presence of hysteresis at the lower transition temperature and gradual transition at the higher transition temperature when g = 10; (II) first-order transition with the presence of hysteresis at the lower and higher transition temperature when g = 100; (III) first-order transition with a hidden, stable low-spin state that remained stable until T = 30 during the warming, and an incomplete spin transition on the plateau during the cooling processes at the lower transition temperature.The first-order transition with the presence of hysteresis during the warming and cooling processes at the higher transition temperature was present in the case of g = 1000.Similar properties were found for n HS in the case of g = 5000.Parameter values are given in the text of the figure caption.
Figure 4 displays the temperature dependence of the order parameters S a , S b , and high-spin fraction n HS for different values of ligand field energy in the absence of an external magnetic field.It is interesting to emphasize that, for equivalent sublattices in Figure 4a,b, the order parameters S a , S b , and the high-spin fraction n HS displayed a stable high-spin (HS) state in the low temperature region and jumped to relax the low-spin state, where ∆ = 85/blue curve, when the temperature increased.During the transition from a low-spin (LS) state to a high-spin (HS) state, a thermal hysteresis loop appeared in the warming and cooling phases.In the case of nonequivalent sublattices, when the intrasublattice parameters were  = 51 and  = 50, the transition was first-order for  ,  , while, for  , it was first order with two steps, as depicted in Figure 3c,d.The high-spin fraction  shifted to the left (lower temperature) as the degeneracy g increased.For non-equivalent sublattices, we found three typical two-step transitions during the warming and cooling processes for the high-spin fraction for the selected values: (I) first-order transition with the presence o hysteresis at the lower transition temperature and gradual transition at the higher transi tion temperature when g = 10; (II) first-order transition with the presence of hysteresis a the lower and higher transition temperature when g = 100; (III) first-order transition with a hidden, stable low-spin state that remained stable until T = 30 during the warming, and an incomplete spin transition on the plateau during the cooling processes at the lower transition temperature.The first-order transition with the presence of hysteresis during the warming and cooling processes at the higher transition temperature was present in the case of g = 1000.Similar properties were found for  in the case of g = 5000.Pa rameter values are given in the text of the figure caption.
Figure 4 displays the temperature dependence of the order parameters  ,  , and high-spin fraction  for different values of ligand field energy in the absence of an ex ternal magnetic field.It is interesting to emphasize that, for equivalent sublattices in Fig ure 4a,b, the order parameters  ,  , and the high-spin fraction  displayed a stable high-spin (HS) state in the low temperature region and jumped to relax the low-spin state The existence of this phenomenon was theoretically predicted by Miyashita et al. [56] and confirmed also by an experimental study on Rb 0.64 Mn[Fe(CN) 6 ] 0.88 •1.7H 2 O [57].Other theoretical studies have also observed this phenomenon when investigating the spincrossover transition of Co-Fe Prussian blue analogs [58].A study of pressure-induced phase transitions in SCO materials also confirmed the presence of a metastable state at low temperatures [30].The stable high-spin state was not observed in the low temperature region for ∆ = 105/red curve, the transition was abrupt, and the hysteresis was present.When the ligand field energy became larger, e.g., ∆ = 125/green curve, the transition was gradual, and the stable state and hysteresis loop were not observed.
The nonequivalent sublattices for J a = 51 and J b = 50 and different values of ligand field energy are shown in Figure 4c,d.The transition of the order parameter S a , for ∆ = 85/blue curve, in the heating and cooling processes was located between 0.8 and 1.0, with the presence of hysteresis.The hysteresis value was ∆T = T ↑ − T ↓ = 26.Sublattice b exhibited a first-order transition of the order parameter S b with the presence of hysteresis and also exhibited the same value (∆T = 26) as that of S a .The transition for S b occurred between −1.0 and 0.79.For ∆ = 105/red curve, the transition became abrupt, and the hysteresis loops narrowed and shifted to a larger temperature region.When ∆ became larger (= 125/green curve), at low temperatures, the order parameters of the two sublattices, S a and S b , became equal to −1, indicating a low-spin state.As the temperature increased, the order parameters remained the same until T ↑ ∼ 43.At this point, the order parameter of sublattice b underwent a first-order phase transition from S b = −0.8(LS state) to S b = 1 (HS state).The order parameter of sublattice a remained in a low-spin state until the temperature reached approximately T ↑ ∼ 95.4, at which point it underwent a first-order phase transition from S a = −0.36(LS state) to S a = 0.68 (HS state).At higher temperatures, both sublattices became equivalent.The value of S a during the cooling process remained in a high-spin state at low temperatures, and no phase transition was observed from    In Figure 4d, the temperature dependences of n HS are shown.When ∆ = 85/blue curve and other parameter values matched those given in the figure caption, in the warming process, we found that the high-spin fraction n HS displayed a stable HS state in the low temperature region and jumped to the intermediate metastable/stable state HS-LS at n HS = 0.5.The study published by Gindulescu et al. [59] shows that compounds exhibiting a two-step thermal transition display a metastable state at low temperatures which cor-responds to a high-spin fraction (n HS = 0.5).At higher region temperatures, the thermal hysteresis was present separately from the stable HS state at low temperatures.For higher values of ∆(= 125), the stable state disappeared and a hidden LS state occurred.This state remained stable until T = 43.8,when the system transitioned from the LS state to the intermediate HS-LS state, leading to a two-step transition.The hidden spin state results align well with experimental findings obtained through the application of pressure or the light-induced excited spin state trapping process [28,[60][61][62][63].
Figure 5  The final simulations for the static properties in this study concern the case where the difference between the intrasublattice parameters is very large,  = 200 and  = 100 We obtained the order parameters and high-spin fraction, which appear in Figure 6a,b The other parameter values were:  = −50, ∆= 900, g = 400 and  = 0.0.Simila properties are depicted in Figure 1c,d.We obtained this separately to study the dynami The final simulations for the static properties in this study concern the case where the difference between the intrasublattice parameters is very large, J a = 200 and J b = 100.We obtained the order parameters and high-spin fraction, which appear in Figure 6a

Time-Dependent Behaviors
Finally, we studied the dynamic properties of the system by solving Equation which represent the equation of motion.Figure 7 presents the time evolution of the o parameters  ,  , and  at three different fixed temperatures:  = 75, 150, and The equivalent sublattice parameters for this case were as follows:  =  = 200,  −50, ∆= 900, g = 400,  = 1, and  = 0.0.These parameters are identical to t shown in Figure 1a,b.Figure 7a,b display the relaxation of  ,  , and  for diffe initial values of the order parameters (−1 to +1 with a step of 0.25) below the hyste temperature ( = 75).The order parameters  ,  , and  evolved from the initial ues to the stable state of the system.In the stable state,  and  had a value of −1 while  had a value of 0 (LS).For the fixed temperature  = 150, the relaxa curves from different initial values of the order parameters were more complex.As sh in Figure 7c,d, the relaxation is characteristic of a typical bistable behavior.The   evolved in time from the initial values to one of the stable states −1 (LS) or +1 ( Similarly, the  changed over time from the initial values to one of the stable stat (LS) or +1 (HS).At  = 200, the  ,  , and  saturated HS (+1) from the initial ues, indicating that the system was in a stable state.

Time-Dependent Behaviors
Finally, we studied the dynamic properties of the system by solving Equation ( 13), which represent the equation of motion.Figure 7 presents the time evolution of the order parameters S a , S b , and n HS at three different fixed temperatures: T = 75, 150, and 200.The equivalent sublattice parameters for this case were as follows: J a = J b = 200, J ab = −50, ∆ = 900, g = 400, k 1 = 1, and H = 0.0.These parameters are identical to those shown in Figure 1a,b.Figure 7a,b display the relaxation of S a , S b , and n HS for different initial values of the order parameters (−1 to +1 with a step of 0.25) below the hysteretic temperature (T = 75).The order parameters S a , S b , and n HS evolved from the initial values to the stable state of the system.In the stable state, S a and S b had a value of −1 (LS), while n HS had a value of 0 (LS).For the fixed temperature T = 150, the relaxation curves from different initial values of the order parameters were more complex.As shown in Figure 7c,d, the relaxation is characteristic of a typical bistable behavior.The S a and S b evolved in time from the initial values to one of the stable states −1 (LS) or +1 (HS).Similarly, the n HS changed over time from the initial values to one of the stable states 0 (LS) or +1 (HS).At T = 200, the S a , S b , and n HS saturated HS (+1) from the initial values, indicating that the system was in a stable state.
Figures 8 and 9 present the time dependence of the order parameters S a , S b , and n HS at four different fixed temperatures: T = 5, 15, 35, and 125.Figures 8a and 9a show S a , S b , and n HS at a fixed temperature of T = 5, which corresponds to the existence of the hidden stable state, and the system staying in the HS state for initial values of order parameters 0.9 and 1.The curves starting from the initial values of the order parameters S a = S b ≤ 0.8 relaxed faster at the LS stable state for S a , S b , and n HS .Similar bistability at low temperatures was observed previously by Enachescu et al. [67].At a constant temperature of T = 15, which corresponds to a stable LS (low spin) state of the S a , S b , and n HS order parameters, the relaxation curve exhibited a sigmoidal shape for higher initial values from the HS to the LS state.This can be observed in Figures 8b and 9b.As the initial values of the order parameters decreased, the sigmoidal shape of the relaxation curve diminished.These properties are confirmed by the experimental results previously obtained for the [Fe 0.5 Zn 0.5 (btr) 2 (NCS) 2 ]H 2 O compound [68], as well as by the experimental and theoretical results observed previously by Chakraboty et al. [69].That is strong evidence of bistable properties, a finding which has been confirmed b several theoretical and experimental studies [15,70,71].At a constant temperature of  125, the system exhibits a stable state, which can be confirmed by the flow diagram.The relaxation curves for different initial values of the order parameters at a temperature of T = 35, within the region where bistability occurs, are shown in Figures 8c and 9c.The specified system parameters (J a = J b = 50, J ab = −50, ∆ = 85/blue curve, g = 10, k 1 = 1, and H = 0.0) in Figure 4a,b confirm the presence of thermal hysteresis around the transition temperature for the equivalent sublattices.The relaxation curves in this region converged toward either one of the two metastable/stable states, depending on the initial values of the order parameters.
That is strong evidence of bistable properties, a finding which has been confirmed by several theoretical and experimental studies [15,70,71].At a constant temperature of T = 125, the system exhibits a stable state, which can be confirmed by the flow diagram.In the previous figures, we assumed that the sublattices a and b were equivale However, we investigated the dynamic properties of nonequivalent sublattices, as show in Figures 10 and 11.Different relaxation curves of the  and  for different initial v ues of the order parameters  =  = 0.75, 0.9, 1 were found at a fixed temperature  5, as shown in Figure 10a.The sublattice a represented by the order parameters  f these three selected values stayed in the HS (=1) steady state, whereas the sublattice represented by the order parameters  relaxed in the LS (=1) steady state.For all t initial values of the order parameters  =  < 0.75, all the curves relaxed into the steady state.With respect to the high-spin fraction  for the same conditions and the initial values of the order parameters  =  = 0.75, 0.9, 1, the relaxation from H converted and reached a plateau at about  = 0.5, depicted in Figure 11a.At this pl eau, as shown, the sublattices were in different spin states (HS-LS).For other initial valu of the order parameters, the  relaxation curve reached LS.At the constant tempe ture of  = 35 (Figure 10b), for the order parameters  and  , the relaxation curv transitioned  ↔  from different initial values of the order parameters.This typi bistability property of the order parameters corresponds to the nonequivalent intrasubl tice parameters.The properties of  and  reported here are comparable with expe mental results for nonequivalent sublattices of the [Fe(etz) ](BF ) spin-crossover sy tem observed previously by Hinek et al. [72,73].While all the relaxation curves of the hig In the previous figures, we assumed that the sublattices a and b were equivalent.However, we investigated the dynamic properties of nonequivalent sublattices, as shown in Figures 10 and 11.Different relaxation curves of the S a and S b for different initial values of the order parameters S a = S b = 0.75, 0.9, 1 were found at a fixed temperature T = 5, as shown in Figure 10a.The sublattice a represented by the order parameters S a for these three selected values stayed in the HS (=1) steady state, whereas the sublattice b represented by the order parameters S b relaxed in the LS (=1) steady state.For all the initial values of the order parameters S a = S b < 0.75, all the curves relaxed into the LS steady state.With respect to the high-spin fraction n HS for the same conditions and for the initial values of the order parameters S a = S b = 0.75, 0.9, 1, the relaxation from HS converted and reached a plateau at about n HS = 0.5, depicted in Figure 11a.At this plateau, as shown, the sublattices were in different spin states (HS-LS).For other initial values of the order parameters, the n HS relaxation curve reached LS.At the constant temperature of T = 35 (Figure 10b), for the order parameters S a and S b , the relaxation curves transitioned HS ↔ LS from different initial values of the order parameters.This typical bistability property of the order parameters corresponds to the nonequivalent intrasublattice parameters.The properties of S a and S b reported here are comparable with experimental results for nonequivalent sublattices of the [Fe(etz) 6 ](BF 4 ) 2 spin-crossover system observed previously by Hinek et al. [72,73].While all the relaxation curves of the high-spin fraction n HS starting from different initial values (Figure 11b) converged to a stable state with a high-spin fraction n HS = 0.5, the dynamic properties of the high-spin fraction matched the steady-state static properties shown in Figure 4c at the temperature of T = 35.
We obtained the relaxation curve in the bistable regions of S a , S b , and n HS at a temperature of T = 60.The relaxation curves started with different initial values of the order parameters and relaxed into different, separated steady states, a phenomenon which is clear evidence of self-organization in the bistable region.Varret et al., during the experimental study of the HS fraction of Fe(ptz) 6 (BF 4 ) 2 , presented the results for n HS (starting values of 0.4 and 0.6) at a constant temperature which remained invariant during more than 5 h of irradiation [74].The bistable state can be confirmed even for the high-spin fraction at this fixed temperature.In the high-temperature (T = 100) region, the relaxation curves from the initial values converged into the stable HS state.The selected value parameters for Figures 10 and 11     We analyzed the dynamic properties of the system for nonequivalent sublattices cases where the difference between the intrasublattice parameters was higher ( = 20 and  = 100), as shown in Figure 6 at the equilibrium property.Other parameters wer  = −50, ∆= 900, g = 400,  = 1, and  = 0.0.As we reported before,  ,  , an  evolved from the initial values to the stable state of the system at a fixed temperatu of  = 5.The values of  and  were both −1 (LS), while the value of  was 0 (LS At a temperature of  = 125, the order parameters  ,  , and  in the bistable regio reached the following values during relaxation from different initial values:  ≃ ±1 fo sublattice a,  ≃ −0.92 for sublattice b, and the high-spin fraction at  = 0.5 an  = 0.02.It is important to emphasize that the different initial values of the order p rameters, e.g.,  = −0.5 and  = 0.5 (curves with a violet color), exhibited a distin relaxation pattern.The  relaxed more rapidly from its initial value to  ≃ −1 whereas the  relaxation process exhibited a two-step behavior with a sigmoidal natu until it reached  ≃ −0.92.The relaxation process of the high-spin fraction exhibited two-step behavior with a sigmoidal nature, which is typical for nonlinear relaxation cooperative origin [43], as presented in Figure 13b.The experiments [33] observed th form of relaxation.
At a temperature of  = 150, a bistable state was observed, as depicted in Figure We analyzed the dynamic properties of the system for nonequivalent sublattices in cases where the difference between the intrasublattice parameters was higher (J a = 200 and J b = 100), as shown in Figure 6 at the equilibrium property.Other parameters were: J ab = −50, ∆ = 900, g = 400, k 1 = 1, and H = 0.0.As we reported before, S a , S b , and n HS evolved from the initial values to the stable state of the system at a fixed temperature of T = 5.The values of S a and S b were both −1 (LS), while the value of n HS was 0 (LS).At a temperature of T = 125, the order parameters S a , S b , and n HS in the bistable region reached the following values during relaxation from different initial values: S a ≃ ±1 for sublattice a, S b ≃ −0.92 for sublattice b, and the high-spin fraction at n HS = 0.5 and n HS = 0.02.It is important to emphasize that the different initial values of the order parameters, e.g., S a = −0.5 and S b = 0.5 (curves with a violet color), exhibited a distinct relaxation pattern.The S a relaxed more rapidly from its initial value to S a ≃ −1, whereas the S b relaxation process exhibited a two-step behavior with a sigmoidal nature until it reached S b ≃ −0.92.The relaxation process of the high-spin fraction exhibited a two-step behavior with a sigmoidal nature, which is typical for nonlinear relaxation of cooperative origin [43], as presented in Figure 13b.The experiments [33] observed this form of relaxation.
At a temperature of T = 150, a bistable state was observed, as depicted in Figure 12c.The order parameters of sublattice a and b both reached a steady state with values of S a ≃ ±1 and S b ≃ ±0.91, respectively.The relaxation of the high-spin fraction took place at two separate steady states, with n HS = 0.52 and n HS = 0.48, depending on the initial values of the order parameters.The nonequilibrium properties corresponded with the static properties observed at the separate branches in the plateau region between the hysteresis loop during both the warming and cooling processes, as depicted in Figure 6b.The order parameters S a , S b , and n HS , depicted in Figures 12d and 13d, were obtained at a constant temperature of T = 225.These order parameters were relaxed from various initial values and reached the HS state.
Magnetochemistry 2024, 10, x FOR PEER REVIEW 18 o two separate steady states, with  = 0.52 and  = 0.48, depending on the initial v ues of the order parameters.The nonequilibrium properties corresponded with the st properties observed at the separate branches in the plateau region between the hyster loop during both the warming and cooling processes, as depicted in Figure 6b.The or parameters  ,  , and  , depicted in Figures 12d and 13d, were obtained at a const temperature of  = 225.These order parameters were relaxed from various initial val and reached the HS state.The above relaxation results can also be observed using the flow diagrams in a t dimensional  −  phase space describing the solutions of the dynamic equation.flow diagram of  and  are presented for the rate constant  = 1 and for diffe temperatures in Figures 14 and 15, respectively.The flow diagrams of the order para ters  and  for equivalent sublattices are shown in Figure 14a-d, under a cons temperature (T = 5, 15, 35, and 125); all other parameters correspond to the situation Figures 4a and 8.The open circle represents the stable equilibrium solution, which co sponds to the lowest values of the free energy or the deepest minimum.On the other h the filled circle represents the unstable solution or state, which corresponds to local m ima (the peaks) or a saddle point in these figures.
The flow diagrams of the order parameters  and  for nonequivalent sublat are shown in Figure 15a-d, under a constant temperature (T = 5, 35, 60, and 100); all o parameters correspond to the situations in Figures 4c and 10.In addition, in Figure 15 open circle represents the stable equilibrium state, whereas the filled circle represents unstable state.Moreover, the filled square depicted in Figure 15c represents the meta ble state.This state is reached when the system relaxes into it, but it does not corresp to the deepest minimum.Instead, it corresponds to a secondary minimum.The nonli differential equation (Equation ( 13)) can be numerically solved using the Runge-K method.4a and 8.The open circle represents the stable equilibrium solution, which corresponds to the lowest values of the free energy or the deepest minimum.On the other hand, the filled circle represents the unstable solution or state, which corresponds to local maxima (the peaks) or a saddle point in these figures.
Consequently, with these dynamic studies and the flow diagrams of  and  , w verified the stable, metastable, and unstable branches of the order parameters and high spin fraction that we had found in the equilibrium study.The flow diagrams of the order parameters S a and S b for nonequivalent sublattices are shown in Figure 15a-d, under a constant temperature (T = 5, 35, 60, and 100); all other parameters correspond to the situations in Figures 4c and 10.In addition, in Figure 15, the open circle represents the stable equilibrium state, whereas the filled circle represents the unstable state.Moreover, the filled square depicted in Figure 15c represents the metastable state.This state is reached when the system relaxes into it, but it does not correspond to the deepest minimum.Instead, it corresponds to a secondary minimum.The nonlinear differential equation (Equation ( 13)) can be numerically solved using the Runge-Kutta method.
Consequently, with these dynamic studies and the flow diagrams of S a and S b , we verified the stable, metastable, and unstable branches of the order parameters and high-spin fraction that we had found in the equilibrium study.

Conclusions
We used the lowest approximation of the cluster variation method and the path pro ability method to analyze the equilibrium and nonequilibrium properties of the order p rameters  ,  , and high-spin fraction  in the spin-crossover system of the two-su lattice model.Equivalent and nonequivalent sublattices were considered for the presente model.The main focus of the study is on how the intersublattice parameter, intrasublatti parameter, degeneracy, ligand field energy, and external magnetic field affect the pha transition properties of the spin-crossover system.The SCO system exhibited gradual an first-order phase transitions with varying values of the model parameters.For nonequi alent sublattices, the order parameters had an intermediate temperature, where subla tices have different states for specific values of the antiferromagnetic ( < 0) intersubla tice interaction.The high-spin fraction  exhibited a first-order phase transition wi two-step properties and a plateau region between the hysteresis loop during both th warming and cooling processes.In addition, we show the presence of a stable high-sp state of  at a low temperature for both equivalent and nonequivalent sublattices.It noteworthy that, for the first time for nonequivalent sublattices, we observed a transitio from the stable high-spin state to the intermediate stable HS-LS state at a value of  0.5 in the low temperature region.The dynamic properties were observed.The phenom enon of bistability occurred at low temperatures.The relaxation curves started with d ferent initial values of the order parameters and relaxed into different, separated stead

Conclusions
We used the lowest approximation of the cluster variation method and the path probability method to analyze the equilibrium and nonequilibrium properties of the order parameters S a , S b , and high-spin fraction n HS in the spin-crossover system of the two-sublattice model.Equivalent and nonequivalent sublattices were considered for the presented model.The main focus of the study is on how the intersublattice parameter, intrasublattice parameter, degeneracy, ligand field energy, and external magnetic field affect the phase transition properties of the spin-crossover system.The SCO system exhibited gradual and first-order phase transitions with varying values of the model parameters.For nonequivalent sublattices, the order parameters had an intermediate temperature, where sublattices have different states for specific values of the antiferromagnetic (J ab < 0) intersublattice interaction.The high-spin fraction n HS exhibited a first-order phase transition with two-step properties and a plateau region between the hysteresis loop during both the warming and cooling processes.In addition, we show the presence of a stable high-spin state of n HS at a low temperature for both equivalent and nonequivalent sublattices.It is noteworthy that, for the first time for nonequivalent sublattices, we observed a transition from the stable high-spin state to the intermediate stable HS-LS state at a value of n HS = 0.5 in the low temperature region.The dynamic properties were observed.The phenomenon of bistability occurred at low temperatures.The relaxation curves started with different initial values of the order parameters and relaxed into different, separated steady states, a phenomenon which is clear evidence of self-organization in the bistable region.It is important to emphasize that, for the different initial values of the order parameters, the relaxation process of the order parameters and the high-spin fraction exhibited a two-step behavior with a sigmoidal nature.With these dynamic studies and the flow diagrams of S a and S b , we also verified the existence of the stable, metastable, and unstable branches of the order parameters and high-spin fraction that have been found in some earlier equilibrium studies [75,76].
plateau ( = ) for both  values.The selected value parameters are presented in th figure caption.

Figure 1 .
Figure 1.Thermal behaviors of the order parameters  ,  , and high-spin fraction  in the ab sence of the external magnetic field for different selected values of the intersublattice paramete ( ).In panels (a,b),  =  = 200, that is, the sublattices are equivalent, while, in panels (c,d)  = 200 and  = 199, that is, the sublattices are nonequivalent.Other parameters are the same ∆= 900, g = 400, and  = 0.0.

Figure 1 .
Figure 1.Thermal behaviors of the order parameters S a , S b , and high-spin fraction n HS in the absence of the external magnetic field for different selected values of the intersublattice parameter (J ab ).In panels (a,b), J a = J b = 200, that is, the sublattices are equivalent, while, in panels (c,d), J a = 200 and J b = 199, that is, the sublattices are nonequivalent.Other parameters are the same: ∆ = 900, g = 400, and H = 0.0.

Figure 2 .
Figure 2. Thermal behaviors of the order parameters  ,  , and high-spin fraction  in the absence of the external magnetic field for different selected values of the intrasublattice parameter ( ,  ).In panels (a,b),  =  , that is, the sublattices are equivalent, while, in panels (c,d),  ≠  that is, the sublattices are nonequivalent.Other parameters are the same:  = −15, ∆= 400, g = 100, and  = 0.0.

Figure 3
Figure3shows the temperature dependence of the  ,  , and  for four different values of degeneracy in the absence of an external magnetic field.As it is shown in Figure3a,b, in the case of equivalent sublattices ( =  = 50), the transition was first-order for large values of g (1000/magenta curve and 5000/green curve) and the whole hysteresis

Figure 2 .
Figure 2. Thermal behaviors of the order parameters S a , S b , and high-spin fraction n HS in the absence of the external magnetic field for different selected values of the intrasublattice parameter (J a , J b ).In panels (a,b), J a = J b , that is, the sublattices are equivalent, while, in panels (c,d), J a ̸ = J b , that is, the sublattices are nonequivalent.Other parameters are the same: J ab = −15, ∆ = 400, g = 100, and H = 0.0.

Figure 3 .
Figure 3. Thermal behaviors of the order parameters  ,  , and high-spin fraction  in the ab sence of the external magnetic field for different selected values of the degeneracy parameter (g).In panels (a,b),  =  = 50, that is, the sublattices are equivalent, while, in panels (c,d),  = 51 and  = 50, that is, the sublattices are nonequivalent.Other parameters are the same:  = −50, ∆= 300, and  = 0.0.

Figure 3 .
Figure 3. Thermal behaviors of the order parameters S a , S b , and high-spin fraction n HS in the absence of the external magnetic field for different selected values of the degeneracy parameter (g).In panels (a,b), J a = J b = 50, that is, the sublattices are equivalent, while, in panels (c,d), J a = 51 and J b = 50, that is, the sublattices are nonequivalent.Other parameters are the same: J ab = −50, ∆ = 300, and H = 0.0.
S a = 1 to S a = −1.The sublattice b, during the cooling process, underwent a first-order phase transition from S b = 0.56 (HS state) to S b = −0.79(LS state).At low temperatures, S b reached a value of −1.Magnetochemistry 2024, 10, x FOR PEER REVIEW 10 of 25

Figure 4 .
Figure 4. Thermal behaviors of the order parameters  ,  , and high-spin fraction  in the ab sence of the external magnetic field for different selected values of the ligand field parameter (∆).In panels (a,b),  =  = 50, that is, the sublattices are equivalent, while, in panels (c,d),  = 51 and  = 50, that is, the sublattices are nonequivalent.Other parameters are the same:  = −50, g = 10, and  = 0.0.

Figure 5
Figure 5 displays the magnetic field dependencies of  ,  , and  for various tem perature values and specific values of the intrasublattice parameter.When  =  = 50 the sublattices are equivalent, while, when  = 51 and  = 50, the sublattices are non equivalent.The other parameters were:  = −50, g = 10, and ∆= 85.In the case o equivalent sublattices,  ,  , and  exhibited a first-order transition at  = 40 when the external magnetic field was increased and decreased.The transition became gradua (smooth) at  = 50 and  = 55 (Figure 5a,b).The order parameters of the system exhib ited a first-order transition with hysteresis and  exhibited a two-step transition with hysteresis at temperatures of 40, 50, and 55 for nonequivalent sublattices (Figure 5c,d) Theoretical research has recently shown similarities in properties when examining the im

Figure 4 .
Figure 4. Thermal behaviors of the order parameters S a , S b , and high-spin fraction n HS in the absence of the external magnetic field for different selected values of the ligand field parameter (∆).In panels (a,b), J a = J b = 50, that is, the sublattices are equivalent, while, in panels (c,d), J a = 51 and J b = 50, that is, the sublattices are nonequivalent.Other parameters are the same: J ab = −50, g = 10, and H = 0.0.

Figure 5 .
Figure5.The external magnetic field evolution of the order parameters S a , S b , and high-spin fraction n HS for different selected values of temperature (T).In panels (a,b), J a = J b = 50, that is, the sublattices are equivalent, while, in panels (c,d), J a = 51 and J b = 50, that is, the sublattices are nonequivalent.Other parameters are the same: J ab = −50, g = 10, and ∆ = 85.

12 Figure 6 .
Figure 6.Thermal behaviors of the order parameters  and  (a) and high-spin fraction  in the absence of the external magnetic field.In panels (a,b):  = 200 and  = 100, that is sublattices are nonequivalent.Other parameters are:  = −50, ∆ = 900, g = 400, and  = 0

Figure 6 .
Figure 6.Thermal behaviors of the order parameters S a and S b (a) and high-spin fraction n HS (b) in the absence of the external magnetic field.In panels (a,b): J a = 200 and J b = 100, that is, the sublattices are nonequivalent.Other parameters are: J ab = −50, ∆ = 900, g = 400, and H = 0.0.

Magnetochemistry 2024 , 2 Figure 7 .
Figure 7. Relaxation curves of the order parameters  and  in panels (a,c,e) and of the high spin fraction  in panels (b,d,f) at various temperatures.These temperatures include those below the hysteretic temperature ( = 75), above the hysteretic temperature ( = 200), and within th bistability region ( = 150).The curves are plotted for different initial values of the order parame ters  and  .Using the following parameter values:  =  = 200 (for equivalent sublattices, a shown in Figure 1a,b,  = −50, ∆= 900, g = 400,  = 1, and  = 0.0).

Figures 8
Figures 8 and 9 present the time dependence of the order parameters  ,  , and  at four different fixed temperatures:  = 5, 15, 35, and 125.Figures 8a and 9a show  ,  , and  at a fixed temperature of  = 5, which corresponds to the existence o the hidden stable state, and the system staying in the HS state for initial values of orde parameters 0.9 and 1.The curves starting from the initial values of the order parameter

Figure 7 .
Figure 7. Relaxation curves of the order parameters S a and S b in panels (a,c,e) and of the high-spin fraction n HS in panels (b,d,f) at various temperatures.These temperatures include those below the hysteretic temperature (T = 75), above the hysteretic temperature (T = 200), and within the bistability region (T = 150).The curves are plotted for different initial values of the order parameters S a and S b .Using the following parameter values: J a = J b = 200 (for equivalent sublattices, as shown in Figure 1a,b, J ab = −50, ∆ = 900, g = 400, k 1 = 1, and H = 0.0).

Figure 8 .
Figure 8. Relaxation curves of the order parameters S a and S b at various temperatures.These temperatures include those below the hysteretic low temperature stable state (T = 5) (a), below the hysteretic temperature (T = 15) (b), within the bistability region (T = 35) (c), and above the hysteretic temperature (T = 125) (d).The curves are plotted for different initial values of the order parameters S a and S b .Using the following parameter values: J a = 50, J b = 50 (representing equivalent sublattices, as shown in Figure 4a), J ab = −50, ∆ = 85, g = 10, k 1 = 1, and H = 0.0.

Figure 9 .
Figure 9. Relaxation curves of the high-spin fraction n HS at various temperatures.These temperatures include those below the hysteretic low temperature stable state (T = 5) (a), below the hysteretic temperature (T = 15) (b), within the bistability region (T = 35) (c), and above the hysteretic temperature (T = 125) (d).The curves are plotted for different initial values of the order parameters S a and S b .Using the following parameter values: J a = 50, J b = 50 (representing equivalent sublattices, as shown in Figure 4b), J ab = −50, ∆ = 85, g = 10, k 1 = 1, and H = 0.0.
are the same as in Figure 4c for ligand field energy ∆ = 85.Magnetochemistry 2024, 10, x FOR PEER REVIEW 16 of state with a high-spin fraction  = 0.5, the dynamic properties of the high-spin fracti matched the steady-state static properties shown in Figure 4c at the temperature of  35.

Figure 10 .
Figure 10.Relaxation curves of the order parameters  and  in the case of  = 51,  = (representing nonequivalent sublattices, as shown in Figure 4c) at various temperatures.These te peratures include those below the hysteretic low temperature ( = 5) (a), below the hysteretic te perature ( = 35) (b), within the bistability region ( = 60) (c), and above the hysteretic temperat ( = 100) (d).Different initial values of the order parameters  and  are considered.Other rameters are the same values as in Figure 4c:  = −50, ∆ = 85, g = 10,  = 1, and  = 0.0.We obtained the relaxation curve in the bistable regions of  ,  , and  at a te perature of  = 60.The relaxation curves started with different initial values of the ord parameters and relaxed into different, separated steady states, a phenomenon which clear evidence of self-organization in the bistable region.Varret et al., during the expe

Figure 10 .
Figure 10.Relaxation curves of the order parameters S a and S b in the case of J a = 51, J b = 50 (representing nonequivalent sublattices, as shown in Figure 4c) at various temperatures.These temperatures include those below the hysteretic low temperature (T = 5) (a), below the hysteretic temperature (T = 35) (b), within the bistability region (T = 60) (c), and above the hysteretic temperature (T = 100) (d).Different initial values of the order parameters S a and S b are considered.Other parameters are the same values as in Figure 4c: J ab = −50, ∆ = 85, g = 10, k 1 = 1, and H = 0.0.

Figure 11 .
Figure 11.Relaxation curves of the high-spin fraction n HS in the case of J a = 51, J b = 50 (representing non-equivalent sublattices, as shown in Figure 4c) at various temperatures.These temperatures include those below the hysteretic low temperature (T = 5) (a), below the hysteretic temperature (T = 35) (b), within the bistability region (T = 60) (c), and above the hysteretic temperature (T = 100) (d).Different initial values of the order parameters S a and S b are considered.Other parameters are the same values as in Figure 4c: J ab = −50, ∆ = 85, g = 10, k 1 = 1, and H = 0.0.

Figure 13 .
Figure 13.Relaxation curves of the high-spin fraction n HS (a) below the hysteretic temperature (T = 5), (b) in the bistability (T = 125) region, (c) in the flat (T = 150) region of the two-step spin-crossover system, and (d) above the hysteretic temperature (T = 225) for different initial values of the order parameters S a , S b .With the parameter values: J a = 200, J b = 100 (nonequivalent sublattices, Figure 6b), J ab = −50, ∆ = 900, g = 400, k 1 = 1, and H = 0.0.The above relaxation results can also be observed using the flow diagrams in a twodimensional S a − S b phase space describing the solutions of the dynamic equation.The flow diagram of S a and S b are presented for the rate constant k 1 = 1 and for different temperatures in Figures 14 and 15, respectively.The flow diagrams of the order parameters S a and S b for equivalent sublattices are shown in Figure 14a-d, under a constant temperature (T = 5, 15, 35, and 125); all other parameters correspond to the situations in Figures4a and 8.The open circle represents the stable equilibrium solution, which corresponds to the lowest values of the free energy or the deepest minimum.On the other hand, the filled circle represents the unstable solution or state, which corresponds to local maxima (the peaks) or a saddle point in these figures.

Figure 14 .
Figure 14.Flow diagram of the system for order parameters S a and S b with different temperature values.(a) for T = 5, (b) for T = 15, (c) for T = 35, and (d) for T = 125, with other parameter values: J a = 50, J b = 50 (equivalent sublattices, Figures 4a and 8), J ab = −50, ∆ = 85, g = 10, k 1 = 1, and H = 0.0.The filled circle corresponds to the unstable state and the open circle corresponds to the stable state.

Figure 15 .
Figure 15.Flow diagram of the system for order parameters S a and S b with different temperature values.(a) for T = 5, (b) for T = 35, (c) for T = 60, and (d) for T = 100, with other parameter values: J a = 51, J b = 50 (nonequivalent sublattices, Figures 4c and 10), J ab = −50, ∆ = 85, g = 10, k 1 = 1, and H = 0.0.The filled circle corresponds to the unstable state, the filled square corresponds to the metastable state, and the open circle corresponds to the stable state.