A comparative study between Monte Carlo entropic sampling method and local mean field investigations of thermal properties of spin-crossover nanoparticles based on Ising-like model

INTRODUCTION

Spin-crossover (SCO) materials^[1–6] are composed of transition metals with an electronic configuration ranging between $$ \mathrm{3d^4} $$ and $$ \mathrm{3d^7} $$ (chromium, manganese, iron, and cobalt). These metals coordinate with (Sulphur, Oxygen, or Nitrogen) atoms in octahedral symmetry. Iron (II)-based SCO solids constitute interesting materials due to their bistable electronic and vibronic nature that allows them to exhibit, as a function of various stimuli, a spin switching between a diamagnetic low-spin state (LS, $$ e_g^0 $$$$ t_{2g}^6 $$) with total spin S = 0 and a paramagnetic high-spin (HS, $$ e_g^2 $$$$ t_{2g}^4 $$) with a spin value, S = 2. Common external stimuli are temperature, pressure, electric fields, high-magnetic fields, and light^[7–20].

In Fe-based SCO materials, the spin transition between LS and HS is accompanied by an increase of $$ \sim10 $$% of the Fe-ligand bond lengths due to the weakening of the metal-ligand bonds accompanying the population of the anti-bonding orbitals $$ e_{g}^{2} $$. This results in the change of the unit cell volume by $$ \sim $$3%-5%. The SCO molecular solids are excellent candidates as sensors of temperature, pressure^[21–24], and gases, notably of hazardous volatile organic compounds, as recently demonstrated in several reports^[25,26]. They can also serve as switches of light emission in the electroluminescent devices or photoluminescent^[27,28] containing SCO complexes or to enhance plasmonic resonances^[29].

To try to explain and reproduce these switching phenomena from the LS to HS states, researchers have proposed different models, some of which are based on Ising-like^[30–35], mechano-elastic^[36], and electroelastic^[37] models, on which switchable core-shell SCO nanoparticles have been investigated in details^[38,39].

In the present study, we focus on the modeling of thermal effects in 2D square SCO nanoparticles. The Ising-like model is used here and solved in the framework of local mean field approximation (LMFA) and within the Monte Carlo entropic sampling (MCES) technique in an attempt to reveal the role of interactions with the external environment (matrix). In this work, we discuss and compare the results obtained by these two methods. Experimental results on the effect of the surface or the substrate have been reported^[35,40,41].

Although the MCES technique is a good choice to solve the Hamiltonian associated with the extended Ising-like model, it is limited in terms of computing capacity for samples containing a small number of molecules. As we would like to analyze the cases of nanoparticles with large numbers of molecules and even thin films in contact with a substrate, which cannot be studied by MCES, we plan to use an alternative formalism based on LMFA. With this contribution, the purpose of this contribution is to determine the reliability, limitations, and drawbacks of the LMFA method when applied to very small systems by comparing it with the MCES technique.

The manuscript is organized as follows. Section 2 introduces the Model and some details of the parameters used in the model. In section 3, the MCES technique is presented, while the LMFA is introduced in section 4. In section 5, we present the simulated thermal behavior of order parameters and the discussion of the obtained results. Finally, we conclude and outline the possible extensions of this work in section 6.

THE MODEL

In condensed matter physics, magnetic phenomena observed in a regular lattice are described with the well-known Ising Model based on interactions between quantum spin variables that were first introduced phenomenologically by Pauli^[42]. Later on, in 1925, Ising^[43] used the concept of spins to describe the thermal behavior of the magnetization of magnetic materials using a 1D model whose statistics were treated in the frame of a transfer matrix method. Pointing quickly the exact resolution of this model (with zero applied field) performed with a remarkable mathematical tour de force of Onsager in 1944^[44], it is also important to notice the considerable efforts of several authors in the realization of approximative methods based on variational treatments of the Ising model, which led to mean-field theory^[45], Bethe-Peierls^[46], etc., which are today widely used to analytically describe a large panel of cooperative phenomena in solid state physics and beyond.

In the field of molecular switchable materials, the Ising model^[30–35] has been adapted to study SCO materials in which two different magnetic states, HS and LS, interact to display at the macroscopic scale various physical properties that give rise to a wide variety of transition behaviors. Among them, one can quote: (ⅰ) thermally-induced hysteretic first-order phase transitions^[15], gradual or continuous^[47] conversions, incomplete^[48], and multi-step transformations^[49].

To model SCO nanoparticles, a spin fictitious operator with two eigen-values +1 and -1 is associated with each molecule describing the HS and LS states, respectively. Then, the Ising-like Hamiltonian is used to model the energy operator connected to SCO. In this respect, the Hamiltonian $$ H $$ can be written as the sum of three energetic contributions related to the isolated molecule, the interaction between molecules, and finally, between the molecules at the boundary with their immediate external environment. Thus, $$ H $$ is expressed as:

with

(ⅰ) The single molecule contribution:

where $$ N $$ is the number of molecules, and $$ \Delta $$ is the energy gap at zero Kelvin between the fundamental LS state and the excited HS state.

(ⅱ) The interaction term:

where $$ J_{i j} $$ is the coupling constant between the spins of the sites $$ i $$ and $$ j $$.

(ⅲ) The interaction between the surface molecules and their environment:

where $$ L $$ is the energy state of each molecule at the surface in addition to the previous ligand $$ \frac{\Delta}{2} $$. The variable, $$ M $$, is the number of molecules $$ \sigma_k $$ at the surface.

In the following, we give a detailed description of each of these terms.

Molecular contribution $$ \boldsymbol{H_1} $$

The two states, HS and LS, have different spin and orbital degeneracies, and the spectrum of intramolecular vibrations of the molecule strongly depends on the spin state of the central metal.

And, as it has been reported^[1,2], the density of states becomes higher in the HS state than in the LS state, which is then equivalent to considering these two states as having different effective degeneracies, denoted here by $$ g_{HS} $$ and $$ g_{LS} $$. For simplicity, we consider $$ g_{HS} $$ and $$ g_{LS} $$ as independent of temperature. From the energetic point of view, these effective degeneracies enter the Hamiltonian as entropic terms, of which only their ratio plays a relevant role. When these specific characteristics of the two states, HS and LS, are considered, the expression of the molecular term $$ H_1 $$ is written as:

(5) $$ H_1=\frac{\Delta-k_B T \ln g}{2} \sum\limits_{i=1}^N \sigma_i $$

Where $$ g=\frac{g_{H S}}{g_{L S}} $$, $$ T $$ is the temperature, and $$ k_B $$ is the usual Boltzmann constant.

At this stage, we see that from the point of view of a single site Hamiltonian, as far as $$ \Delta- $$$$ k_B T \ln g>0 $$, LS states, for which $$ \sigma=-1 $$, are favored. When $$ \Delta-k_B T \ln g<0 $$, the HS states, for which $$ \sigma=+1 $$ are favored. At the equilibrium temperature $$ T_{e q}=\frac{\Delta}{k_B \ln g} $$, there is a "crossover" situation between these two configurations. The physical quantity that allows the thermal behavior of the system is the HS fraction, denoted here as $$ N h s $$, which represents the probability of the population being in the HS state. In terms of spin values, it represents the probability of occupation of the spin value $$ \sigma=+1 $$. It is straightforward to establish that $$ N h s $$ relates to the average value of the fictitious magnetization $$ <\sigma> $$ as follows:

(6) $$ Nhs=\frac{1+<\sigma>}{2} $$

Simple Boltzmann statistics performed on Hamiltonian (1) leads to the following relations,

(7) $$ <\sigma>=\tanh \beta\left(\frac{\Delta-k_B T \ln g}{2}\right) $$

and

(8) $$ { Nhs }=\frac{1}{2}\left[1+\tanh \beta\left(\frac{\Delta-k_B T \ln g}{2}\right)\right] $$

from which the curves of Figure 1 showing thermal dependence of the HS magnetization and fraction are easily derived.

Figure 1. Typical thermal evolution (A) of the average magnetization $$ <\sigma> $$; (B) of the high spin fraction $$ Nhs $$. The calculation parameters are: $$ \Delta / k_B=1, 300\ \mathrm{\sim K} $$ and $$ \ln (g)=6.01 $$.

Interaction term $$ \boldsymbol{H_2} $$

The interaction term can be split into two contributions: (ⅰ) those arising from nearest neighbors, which are parametrized using a short-range interaction term $$ J $$ (assumed here to be independent of the site); and (ⅱ) those including all the other molecules, which are parametrized with a long-range interaction term $$ G $$, which is written in the frame of the meanfield approach. In this respect, $$ H_2 $$ is written as:

(9) $$ H_2=-J \sum\limits_{<i, j>} \sigma_i \sigma_j-G<\sigma>\sum\limits_{i=1}^N \sigma_i $$

Considering solely this interaction term $$ \mathrm{H}_2 $$, its resolution leads to the thermal evolution of $$ <\sigma> $$, as given in Figure 2, which is typical of second-order phase transition. Indeed, at 0 K, the system is either in the ordered "HS" state with all the spins with $$ \sigma = 1 $$ or in the "LS" state with $$ \sigma = -1 $$. When temperature increases, the average value of the spin either decreases or increases, as shown in the upper and lower curves. The temperature at which the average value is zero is called the order-disorder temperature $$ T_{O/D} $$ or critical temperature $$ T_C $$.

Figure 2. Thermal evolution of the average magnetization $$ <\sigma> $$ obtained in a 6 $$ \times $$ 6 square lattice with the following parameters: $$ J/k_B $$ = 60 K and $$ G/k_B $$ = 172 K.

The region on the left-hand side below the upper curve (respectively above the lower curve) with respect to this $$ \mathrm{T}_{\mathrm{OD}} $$ temperature corresponds to an ordered phase, while the right-hand side corresponds to the disordered phase. For an infinite square lattice with $$ G = 0 $$, the critical temperature $$ T_C = 2.269\ \mathrm{\sim J} $$ is obtained in line with the "Onsager solution", where $$ T_C $$ is given by the equation $$ \tanh \frac{2 J}{k_B T_C}=\frac{1}{\sqrt{2}} $$.

However, for $$ G \neq 0 $$, there is no exact analytical solution since the mean-field term acts as an "external field" linearly dependent on the magnetization. On the other hand, for a finite system, $$ T_C $$ depends on the lattice size. It is worth noting the symmetric character of the thermal dependence of the net magnetization around $$ <\sigma> $$ = 0 for the ordered phases, obtained by starting from +1 or -1, as shown in Figure 2, calculated using J and G in a 6 $$ \times $$ 6 square lattice. It is important to mention that this symmetric character remains valid when including the long-range interaction term in the Hamiltonian $$ H_2 $$, although the expression of $$ T_C $$ now becomes dependent on $$ J $$ and $$ G $$^[14].

RESOLUTION BY MONTE CARLO ENTROPIC SAMPLING METHOD

In this section, we will distinguish three kinds of short-range interaction terms, namely, $$ J_{bb} $$, the interaction between molecules in the bulk, $$ J_{ss} $$, the interaction between the molecules at the surface, and $$ J_{bs} $$, the interaction between the molecules at the surface and in the bulk. The corresponding Hamiltonian is written as:

with

where the respective total and surface magnetizations, $$ m_s $$ and $$ m_t $$, respectively, are given by:

$$ H $$ can be re-written as:

Where $$ s_b $$, $$ s_s $$, and $$ s_{bs} $$ are the short-range correlation functions related to the bulk, surface, and bulk-surface neighboring sites, respectively. Their expressions are given by the following relations:

The thermal average value of the fictitious magnetization, $$ <\sigma> $$, is calculated by the following expression:

Where $$ Z $$ is the partition function $$ \ E_i $$ is the energy of the macrostate $$ i $$, and $$ d\left(m_{t_i}, s_{b_i}, s_{s_i}, s_{b s_i}, m_{s_i}\right) $$ is the density of this macrostate. Here, $$ Z $$ is given by:

and

In Equations (24) and (25), $$ N_L $$ is the number of different configurations with the same five values $$ m_t, s_b, s_s, s_{b s} $$, and $$ m_s $$. The density of states $$ d\left(m_t, s_b, s_s, s_{b s}, m_s\right) $$ is calculated by entropic sampling^[48–50], and Equation (24) is solved by numerical techniques such as bisection. From the thermal average values $$ <\sigma> $$, the HS fraction $$ Nhs $$ is calculated using the relation (6).

Let us recall the basic principles of MCES^[50–52]. It consists in introducing the detailed balance equation of the Monte Carlo (MC) procedure, as shown in:

This suited biased distribution $$ P $$ is designed to favor configurations belonging to weakly degenerated macrostates while dampening those within the highly degenerated macrostates. It can be expressed as:

In this case, the balance equation can be written as:

Using $$ d\left(m_{t_i}, s_{b_i}, s_{s_i}, s_{b s_i}, m_{s_i}\right) $$ as a bias, a MC sampling is run; it is termed a "Monte Carlo step" and yields a histogram of the frequency of the macrostates $$ H\left(m_{t_i}, s_{b_i}, s_{s_i}, s_{b s_i}, m_{s_i}\right) $$, also written as $$ H_i\left(m_{t_i}, s_{b_i}, s_{s_i}, s_{b s_i}, m_{s_i}\right) $$.

By construction,

The resulting restricted density of states is calculated after the correction for the bias was applied:

The flat character of the histogram $$ H\left(m_t, s_b, s_s, s_{b s}, m_s\right) $$ has a convenient convergence criterion to obtain a very good estimation of the density of the states. In our case, seven iterations of the $$ 10^6 $$ Carlo Steps were necessary. After the table of the values of $$ \left(m_t, s_b, s_s, s_{b s}, m_s\right) $$ and $$ \mathrm{d}\left(m_t, s_b, s_s, s_{b s}, m_s\right) $$ is obtained, the quasi-exact partition function can be calculated using the expression (25).

LOCAL MEAN-FIELD APPROACH

In the present work, three types of sites related to three lattice regions are considered: the atoms located in the bulk ($$ \mathrm{N_b} $$), which are surrounded by four first-neighbors and those atoms located on the surface and, more precisely, on the edge ($$ \mathrm{N_e} $$) and on the corner ($$ \mathrm{N_c} $$), which interact with three and two first-neighbors, respectively. The molecules located at the surface (edge and corner) interact with an external environment (matrix effect), which gives them specific properties.

The case of a system comprising 36 molecules is represented in Figure 3.

Figure 3. Schematic representation of a 2D 6 $$ \times $$ 6 square-shaped lattice: filled blue circles represent bulk sites ($$ \mathrm{N_b} $$); filled green and red circles represent edge ($$ \mathrm{N_e} $$) and corner ($$ \mathrm{N_c} $$) sites, respectively, interacting with their immediate environment (matrix).

Note that in the framework of a classical mean-field approximation^[53] and the case of a 2D square-shaped lattice, the coordination number is q = 4 [Figure 4].

Figure 4. Mean-field approximation in a 2D square-shaped lattice. The $$ \sigma_i $$ spin is surrounded by four nearest neighbors with the average magnetization $$ <\sigma> $$.

The consideration of three types of sites leads us to distinguish three local average order parameters $$ <\sigma_\alpha> $$ and three coordination numbers $$ \mathrm{q}_\alpha $$. In the LMFA, which amounts to taking into account a local Hamiltonian for each type of site $$ \alpha $$ with $$ \alpha=b, e, c $$, corresponding to bulk, edge, and corner, respectively.

For each region, $$ \alpha $$, we apply the LMFA. In the case of regions connecting each other only via the long-range interaction, each molecule $$ \sigma_{i} $$ is surrounded by $$ q_\alpha<\sigma_\alpha> $$, and the local Hamiltonian, given by Equation (14), can be expressed as Equation (32) below:

with

where $$ q_\alpha $$ is the number of interactions between a molecule and its first-neighbors, $$ z_\alpha $$ is the number of interactions $$ L $$ between a molecule and the external environment, and $$ J_\alpha $$ is the interaction between the nearest neighbor $$ (\mathrm{nn}) $$ molecules. For molecules located in the edge and in the corner, the average value of the number of interactions with the external environment is set equal to 1. Table 1 summarizes the characteristics of each site.

Here, for simplicity of calculations, we have taken z = 1 (the same coordination number) everywhere in the surface (edge and corner), meaning that we are considering the same connectivity for all molecules at the surface.

The short-range interaction term $$ J_\alpha $$ depends on the location of the molecules. Two interaction constants are considered in our calculations:

$$ J_{bulk}^{LMFA} $$ for the molecules located in the bulk; $$ J_{surf}^{LMFA} $$ for the molecules located at the surface.

These two interaction terms, $$ J_{bulk}^{LMFA} $$ and $$ J_{surf}^{LMFA} $$, are related to the interaction terms used in the MCES method by weighting the various interactions of each molecule in each region. Let us consider the case of the bulk region. We denote by $$ N_{bb} $$, $$ N_{bs} $$, the respective bulk-bulk and the bulk-surface bonds, leading to a total number of bonds including at least one bulk site ($$ N_{bb} $$ + $$ N_{bs} $$).

Let $$ J_{bb}^{MCES} $$ and $$ J_{bs}^{MCES} $$ be the bulk-bulk and bulk-surface short-range interactions used in MCES method. The total interaction energy term in the bulk region is given by ($$ N_{bb} $$×$$ J_{bb}^{MCES} $$) + ($$ N_{bs} $$×$$ J_{bs}^{MCES} $$). Therefore, the bulk value $$ J_{bulk}^{LMFA} $$ used in Equation (32) is obtained by averaging over $$ J_{bb}^{MCES} $$, $$ J_{bs}^{MCES} $$ weighted by the numbers of their corresponding bonds, as follows:

In a similar way, for the surface sites, Equation (35) is written as:

where $$ N_{ss} $$ and $$ J_{ss}^{MCES} $$ are the number surface-surface bonds and the surface-surface shortrange interaction term used in the MCES method, respectively.

The average spin state can be written as follows:

So, for each region $$ b $$ (bulk), $$ e $$ (edge), and $$ c $$ (corner), the average spin state is given by:

which is equivalent to:

The weighted average of the three local order parameters associated with the three regions, bulk, edge, and corner, leads to the average fictitious spin of the system expressed as:

and the high-spin fraction, Nhs, which is the probability that the HS state is occupied and which is written as Equation (6): $$ Nhs = \frac{1+<\sigma>}{2} $$.

The purpose is to solve a system of three equations $$ \xi_1 $$, $$ \xi_2 $$, and $$ \xi_3 $$, whose three order parameters $$ <\sigma_b> $$, $$ <\sigma_e> $$, and $$ <\sigma_c> $$, connected by Equation (43), are simultaneous solutions. The Newton method, which is appropriate for rapid convergences, has been used. Calculations are based on three points $$ x_1 $$, $$ x_2 $$, and $$ x_3 $$ close to the solution and which represent $$ <\sigma_b> $$, $$ <\sigma_e> $$, and $$ <\sigma_c> $$, respectively.

The calculation of the three points "close to the solution" is performed in two steps: Let X be the space in which values of ($$ x_1 $$, $$ x_2 $$, $$ x_3 $$) are scanned. Let $$ \xi $$ be the space from which ($$ \xi_1 $$, $$ \xi_1 $$, $$ \xi_1 $$) are calculated as a function of ($$ x_1 $$, $$ x_2 $$, $$ x_3 $$). Then, the first step consists of mapping space X on space $$ \xi $$. Those values of ($$ x_1 $$, $$ x_2 $$, $$ x_3 $$) that map onto ($$ \xi_1 $$, $$ \xi_2 $$, $$ \xi_3 $$) close to zero within $$ \mathrm{10_{-6}} $$, are then selected to perform the second step: application of the Newton-Raphson technique.

Solving the system amounts to calculating $$ h_1 $$, $$ h_2 $$, and $$ h_3 $$ such that:

which allows to write in the vicinity of points $$ x_1 $$, $$ x_2 $$, and $$ x_3 $$:

The following system of three equations with three unknowns is thus obtained as:

And where $$ h_1, h_2 $$, and $$ h_3 $$ are calculated as follows:

In a general way, the $$ \xi_i $$ functions with $$ \mathrm{i} = 1 $$ (bulk), 2 (corner), and 3 (edge) and their derivatives can be expressed as follows:

With the $$ h_i $$ obtained in this way, the solutions $$ \left(x_i\right)_{s o l} $$ are deduced by the relation:

RESULTS AND DISCUSSION

We are interested in the effect of the coupling between surface and bulk molecules, and for that, we perform numerical simulations by MCES and by LMFA. To proceed, we choose a set of realistic SCO thermodynamic parameters, and for this purpose, we chose the data derived from experimental data of one of these SCO compounds: [$$ \mathrm{Fe(btr)_{2}(NCS)_{2}} $$], btr = 4, 4′‐bis‐1, 2, 4‐triazole^[54], whose molecular structure is shown in Figure 5. The energy gap $$ \Delta $$ and ln($$ g $$) are derived from the enthalpy change ($$ \Delta H \approx 11 $$ kJ/mol) and the entropy change ($$ \Delta S \approx 50 $$ J/mol/K), respectively: $$ \Delta = \Delta H/R \approx 1, 300 $$ K), ln⁡($$ g $$) = $$ \Delta S/R \approx 6.01 $$, and $$ R $$ is the perfect gas constant. The equilibrium temperature of the system $$ T_{eq} $$ is deduced by $$ T_{eq}=\frac{\Delta/{k_B}}{ln⁡(g)} $$, which leads to $$ \approx 216.3 $$ K. In the first step, in the case of a 6 $$ \times $$ 6 configuration, the values of the short-range interaction parameters $$ J_{bb}^{MCES} $$ and $$ J_{ss}^{MCES} $$ are set to 60 K, and the intensity of the short-range interaction parameter $$ J_{bs}^{MCES} $$ between bulk and surface is gradually increased from 0 to 60 K. In the second step, the size effect is studied with each method. It should be added that the dimensions of the square lattice used in these simulations do not correspond, for the moment, to synthesized nanoparticles.

Figure 5. Schematic structure of the SCO compound iron (II) [$$ \mathrm{Fe(btr)_{2}(NCS)_{2}} $$], comprising a bis-triazole (btr) whose molecular formula are $$ \mathrm{C_{4}N_{6}H_{4}} $$ and a cyano group NCS.

The case $$ \boldsymbol{J_{bs}^{MCES} = 0} $$

Figure 6 displays the results obtained in the case of a 6 $$ \times $$ 6 square-shaped lattice, and the values of the interaction parameters $$ J $$ used by MCES and by LMFA, calculated from Equations (34) and (35), are gathered in Table 2. The curves highlight a two-step behavior in the form of two well-separated hysteresis loops. The hysteresis width is defined as $$ \Delta T = T_{up} - T_{down} $$, where $$ T_{up} $$ is the ascending thermal transition temperature and $$ T_{down} $$ is the descending thermal transition temperature. The average temperature between $$ T_{up} $$ and $$ T_{down} $$ is the equilibrium temperature $$ T_{eq} $$, which corresponds to a HS fraction equal to 1/2($$ Nhs = 1/2 $$).

Figure 6. Comparison of the thermal evolution of the HS fraction in a SCO system for a 6 $$ \times $$ 6 square lattice: MCES (black squares), LMFA (blue up triangles). The computational parameters are: $$ J_{bb}^{MCES})/k_{B} $$ = 60 K, $$ J_{ss}^{MCES}/k_{B} $$ = 60 K, $$ J_{bs}^{MCES}/k_{B} = 0 \mathrm{K} $$, $$ J_{bb}^{LMFA}/k_{B} $$ = 36 K, $$ J_{ss}^{LMFA}/k_{B} $$ = 33.33 K, $$ \Delta/k_{B} $$ = 1, 300 K, $$ G/k_{B} $$ = 172 K, $$ L/k_{B} $$ = 290 K, and ln($$ g $$) = 6.01.

Table 2

Correspondence between the values of short-range interactions $$ \boldsymbol{J/k_B} $$ used by LMFA and by MCES method in the case of a 6 $$ \times $$ 6 square lattice when $$ \boldsymbol{J_{bs}^{MCES}/k_B} $$= 0 K

$$ \boldsymbol{J_{bb}^{MCES}/k_{B} \mathrm{[K]}} $$	$$ \boldsymbol{J_{ss}^{MCES}/k_{B} \mathrm{[K]}} $$	$$ \boldsymbol{J_{bs}^{MCES}/k_{B} \mathrm{[K]}} $$	$$ \boldsymbol{J_{bulk}^{LMFA}/k_{B} \mathrm{[K]}} $$	$$ \boldsymbol{J_{surf}^{LMFA}/k_{B} \mathrm{[K]}} $$
60	60	0	36	33.33

Figure 6 illustrates a qualitative and quantitative agreement between the calculations obtained by MCES and LMFA. The hysteresis loops, observable for each method, at lower temperatures, $$ \sim144 $$ K, are related to the behavior of surface molecules (edge and corner). These surface molecules commute from LS to HS state before the bulk molecules, and they drive the thermal transition. The hysteresis loops, positioned for each method at $$ 185 $$ K, are related to the bulk molecules. The two methods of simulation lead to identical values of the transition temperature in the surface ($$ T_{eq surf}^{LMFA} = T_{eq surf}^{MCES}\sim144 $$ K) and in the bulk ($$ T_{eq bulk}^{LMFA} = T_{eq bulk}^{MCES} \sim185 $$ K). The bulk and the surface appear as two systems that evolve almost independently of each other. Between 150 and 182 K, the curves reveal, both by MCES and by LMFA, a long intermediate plateau (mixture of LS and HS configurations) around $$ Nhs \approx $$ 0.55.

It can be noticed that the hysteresis widths $$ \Delta T $$ are slightly larger by LMFA than by MCES. The values of the thermal hysteresis widths are summarized in Table 3.

Table 3

Values of the thermal hysteresis width for molecules located in the bulk and in the surface obtained by MCES, $$ \boldsymbol{\Delta T_{bulk}^{MCES}} $$ and $$ \boldsymbol{\Delta T_{surf}^{MCES}} $$ and by LMFA, $$ \boldsymbol{\Delta T_{bulk}^{LMFA}} $$ and $$ \boldsymbol{\Delta T_{surf}^{LMFA}} $$

$$ \boldsymbol{\Delta T_{bulk}^{MCES} \mathrm{[K]}} $$	$$ \boldsymbol{\Delta T_{bulk}^{LMFA} \mathrm{[K]}} $$	$$ \boldsymbol{\Delta T_{surf}^{MCES} \mathrm{[K]}} $$	$$ \boldsymbol{\Delta T_{surf}^{LMFA} \mathrm{[K]}} $$
The simulation parameters are: $$J_{bb}^{MCES}/k_B = 60$$ K, $$J_{ss}^{MCES}/k_B =$$ 60 K, $$J_{bs}^{MCES}/k_B = 0$$, $$J_{bb}^{LMFA}/k_B = 36\ \mathrm{K}$$, $$J_{ss}^{LMFA}/k_B = 33.33$$ K, $$\Delta/k_B = 1, 300$$ K, $$G/k_B = 172$$ K, $$L/k_B = 290$$ K, and ln($$g$$) = 6.01.
4.8	7	7.6	10.5

The case $$ \boldsymbol{J_{bs}^{MCES} \neq} $$ 0

The value of the short-range interaction parameter $$ J_{bs}^{MCES}/k_B $$ is gradually increased up to 60 K. Figure 7, Figure 8 and Figure 9, show the results obtained respectively for $$ J_{bs}^{MCES}/k_B = 20 $$, $$ J_{bs}^{MCES}/k_B = 40 $$ and $$ J_{bs}^{MCES}/k_B = 60 $$ K in the case of a 6 $$ \times $$ 6 lattice. The values of the interaction parameters $$ J $$ used by MCES and by LMFA, obtained from Equations (34) and (35), are gathered in Table 4.

Figure 7. Comparison of the thermal evolution of the HS fraction in a SCO system for a 6 $$ \times $$ 6 square lattice: MCES (black squares), LMFA (blue up triangles). The computational parameters are: $$ J_{bb}^{MCES}/k_B = 60\ \mathrm{K} $$, $$ J_{ss}^{MCES}/k_B = 60\ \mathrm{K} $$, $$ J_{bs}^{MCES}/k_B = 20\ \mathrm{K} $$, $$ J_{bb}^{LMFA}/k_B = 44\ \mathrm{K} $$, $$ J_{ss}^{LMFA}/k_B = 42.22\ \mathrm{K} $$, $$ \Delta/k_B = 1, 300\ \mathrm{K} $$, $$ G/k_B = 172\ \mathrm{K} $$, $$ L/k_B = 290\ \mathrm{K} $$, and ln($$ g $$) = 6.01.

Figure 8. Comparison of the thermal evolution of the HS fraction in a SCO system for a 6 $$ \times $$ 6 square lattice: MCES (black squares), LMFA (blue up triangles). The computational parameters are: $$ J_{bb}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{ss}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{bs}^{MCES}/k_{B} = 40\ \mathrm{K} $$, $$ J_{bb}^{LMFA}/k_{B} = 52\ \mathrm{K} $$, $$ J_{ss}^{LMFA}/k_{B} = 51.11\ \mathrm{K} $$, $$ \Delta/k_{B} = 1, 300\ \mathrm{K} $$, $$ G/k_{B} = 172\ \mathrm{K} $$, $$ L/k_{B} = 290\ \mathrm{K} $$, and ln($$ g $$) = 6.01.

Figure 9. Comparison of the thermal evolution of the HS fraction in a SCO system for a 6 $$ \times $$ 6 square lattice: MCES (black squares), LMFA (blue up triangles). The computational parameters are: $$ J_{bb}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{ss}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{bs}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{bb}^{LMFA}/k_{B} = 60\ \mathrm{K} $$, $$ J_{ss}^{LMFA}/k_{B} = 60\ \mathrm{K} $$, $$ \Delta/k_{B} = 1, 300\ \mathrm{K} $$, $$ G/k_{B} = 172\ \mathrm{K} $$, $$ L/k_{B} = 290\ \mathrm{K} $$, and ln($$ g $$) = 6.01.

Table 4

Correspondence between the values of short-range interactions $$ \boldsymbol{J/k_B} $$ used by MCES and by LMFA for a 6 $$ \times $$ 6 lattice when $$ \boldsymbol{J_{bs}^{MCES}/k_B} $$= 20, 40, and 60 K

$$ \boldsymbol{J_{bb}^{MCES}/k_{B} \mathrm{[K]}} $$	$$ \boldsymbol{J_{ss}^{MCES}/k_{B} \mathrm{[K]}} $$	$$ \boldsymbol{J_{bs}^{MCES}/k_{B} \mathrm{[K]}} $$	$$ \boldsymbol{J_{bulk}^{LMFA}/k_{B} \mathrm{[K]}} $$	$$ \boldsymbol{J_{surf}^{LMFA}/k_{B} \mathrm{[K]}} $$
60	60	20	44	42.22
60	60	40	52	51.11
60	60	60	60	60

The phase diagrams of the 6 $$ \times $$ 6 system, obtained by LMFA and by MCES, are plotted in Figure 10 and Figure 11, in the space coordinates, temperature vs. $$ J_{bs}^{MCES}/k_{B} $$, for the bulk and for the surface molecules, respectively.

Figure 10. Phase diagram T vs. $$ J_{bs}^{MCES}/k_{B} $$ for a 6 $$ \times $$ 6 system. The red, blue and black squares correspond, respectively, to the upper and lower transitions ($$ T_{up} $$) and ($$ T_{down} $$) temperatures of the thermal HS fraction and to the equilibrium temperature of the bulk by LMFA. The red, blue and black stars correspond, respectively, to the upper and lower transitions ($$ T_{up} $$) and ($$ T_{down} $$) temperatures of the thermal HS fraction and to the equilibrium temperature of the bulk by MCES. The computational parameters are: $$ J_{bb}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{ss}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ \Delta/k_{B} = 1, 300\ \mathrm{K} $$, $$ G/k_{B} = 172\ \mathrm{K} $$, $$ L/k_{B} = 290\ \mathrm{K} $$, and ln($$ g $$) = 6.01.

Figure 11. Phase diagram T vs. $$ J_{bs}^{MCES}/k_{B} $$ for a 6 $$ \times $$ 6 system. The red, blue and black squares correspond, respectively, to the upper and lower transitions ($$ T_{up} $$) and ($$ T_{down} $$) temperatures of the thermal HS fraction and to the equilibrium temperature of the surface by LMFA. The red, blue and black stars correspond, respectively, to the upper and lower transitions ($$ T_{up} $$) and ($$ T_{down} $$) temperatures of the thermal HS fraction and to the equilibrium temperature of the surface by MCES. The computational parameters are: $$ J_{bb}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{ss}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ \Delta/k_{B} = 1, 300\ \mathrm{K} $$, $$ G/k_{B} = 172\ \mathrm{K} $$, $$ L/k_{B} = 290\ \mathrm{K} $$, and ln($$ g $$) = 6.01.

By MCES, as can be seen in Figure 10 and Figure 11, when the interaction parameter $$ J_{bs}^{MCES}/k_{B} $$ increases, the equilibrium temperature of the surface molecules, $$ T_{{eq}\ {surf}}^{MCES} $$, is shifted to higher temperatures, whereas the equilibrium temperature of the bulk molecules, $$ T_{{eq}\ {bulk}}^{MCES} $$, is shifted to lower temperatures. The hysteresis loops, quite distinct for $$ J_{bs}^{MCES}/k_{B} = 0 $$, therefore, tend to overlap. The thermal hysteresis loop $$ \Delta T_{bulk}^{MCES} $$ increases ($$ \sim4.8\ \mathrm{K} $$ wide when $$ J_{bs}^{MCES}/k_{B} = 0 $$ and $$ \sim8\ \mathrm{K} $$ wide when $$ J_{bs}^{MCES}/k_{B} = 60 $$). Simultaneously, the thermal hysteresis loop $$ \Delta T_{surf}^{MCES} $$ decreases ($$ \sim7.6\ \mathrm{K} $$ wide for $$ J_{bs}^{MCES}/k_{B} = 0 $$ and $$ \sim4\ \mathrm{K} $$ wide for $$ J_{bs}^{MCES}/k_{B} = 60 $$).

In LMFA, as shown in Figure 10 and Figure 11, an opposite behavior is observed for surface molecules. First of all, when the interaction parameter $$ J_{bs}^{MCES}/k_{B} $$ increases, the equilibrium temperature of the surface molecules, $$ T_{{eq}\ {surf}}^{LMFA} $$, is shifted to lower temperatures. Then, the thermal hysteresis loop $$ \Delta T_{surf}^{LMFA} $$ increases ($$ \sim10.5\ \mathrm{K} $$ wide for $$ J_{bs}^{MCES}/k_{B} = 0 $$ and $$ \sim38.5\ \mathrm{K} $$ wide for $$ J_{bs}^{MCES}/k_{B} = 60 $$). These results are summarized in Table 5 and the results obtained when $$ J_{bs}^{MCES}/k_{B} = 0 $$ are recalled for a better understanding of the results. The LMFA drastically accentuates the hysteresis effect, especially for surface molecules. The system seems constrained, and the tendency for the two hysteresis to overlap is not observed.

Table 5

Values of the thermal hysteresis width for molecules located in the bulk and in the surface in MCES, $$ \boldsymbol{\Delta T_{bulk}^{MCES}} $$ and $$ \boldsymbol{\Delta T_{surf}^{MCES}} $$, and in LMFA, $$ \boldsymbol{\Delta T_{bulk}^{LMFA}} $$ and $$ \boldsymbol{\Delta T_{surf}^{LMFA}} $$ as a function of $$ \boldsymbol{J_{bs}^{MCES}/k_{B}} $$

$$ \boldsymbol{J_{bs}^{MCES}/k_{B} \mathrm{[K]}} $$	$$ \boldsymbol{\Delta T_{bulk}^{MCES} \mathrm{[K]}} $$	$$ \boldsymbol{\Delta T_{bulk}^{MCES} \mathrm{[K]}} $$	$$ \boldsymbol{\Delta T_{bulk}^{LMFA} \mathrm{[K]}} $$	$$ \boldsymbol{\Delta T_{surf}^{LMFA} \mathrm{[K]}} $$
The computational parameters are: $$\Delta/k_B = 1, 300 \mathrm{K}$$, $$G/k_B = 172\ \mathrm{K}$$, $$L/k_B = 290\ \mathrm{K}$$, and ln($$g$$) = 6.01.
0	4.8	7.6	7	10.5
20	6.4	6.2	17.5	20.36
40	8.0	5.0	30.5	29
60	8.8	4.0	45.0	38.5

The case $$ \boldsymbol{J_{bs}^{MCES} = 0} $$

The results obtained and presented in Figure 12, highlight three important points. First of all, as previously observed in the case of the 6 $$ \times $$ 6 lattice, when $$ J_{bs}^{MCES}/k_{B} = 0 $$, the equilibrium temperatures obtained in the bulk and in the surface are strictly identical with the two methods. The values of the equilibrium temperatures are reported in Table 7. Then, the equilibrium temperatures, in the bulk and in the surface, are slightly higher than those observed in the case of the 6 $$ \times $$ 6 lattice (See also Table 8 for a comparison). This point may be explained by the fact that for larger sizes of particles, $$ T_{eq} $$ and $$ T_{c.} $$ shift to higher temperatures due to the increase in the number of couplings.

Figure 12. Comparison of the thermal evolution of the HS fraction in a SCO system with a size 7 $$ \times $$ 7 showing an excellent agreement between the MCES method (black squares) and the LMFA (blue up triangles). The calculation parameters are: $$ J_{bb}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{ss}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{bs}^{MCES}/k_{B} = 0\ \mathrm{K} $$, $$ J_{bb}^{LMFA}/k_{B} = 40\ \mathrm{K} $$, $$ J_{ss}^{LMFA}/k_{B} = 32.72\ \mathrm{K} $$, $$ \Delta/k_{B} = 1, 300\ \mathrm{K} $$, $$ G/k_{B} = 172\ \mathrm{K} $$, $$ L/k_{B} = 290\ \mathrm{K} $$, and ln($$ g $$) = 6.01.

Table 7

Values of the equilibrium temperatures of the bulk and surface molecules obtained in MCES and in LMFA, as a function of the system size (6 $$ \times $$ 6 and 7 $$ \times $$ 7 lattices)

$$ \boldsymbol{\mathrm{Lattice size}} $$	$$ \boldsymbol{T_{{eq}\ {bulk}}^{MCES} \mathrm{[K]}} $$	$$ \boldsymbol{T_{{eq}\ {bulk}}^{LMFA} \mathrm{[K]}} $$	$$ \boldsymbol{T_{{eq}\ {surf}}^{MCES} \mathrm{[K]}} $$	$$ \boldsymbol{T_{{eq}\ {surf}}^{LMFA} \mathrm{[K]}} $$
The computational parameters are: $$J_{bb}^{MCES}/k_B = 60\ \mathrm{K}$$, $$J_{ss}^{MCES}/k_B = 60\ \mathrm{K}$$, $$J_{bs}^{MCES}/k_B = 0\ \mathrm{K}$$, $$J_{bb}^{LMFA}/k_B = 40\ \mathrm{K}$$, $$J_{ss}^{LMFA}/k_B = 32.72\ \mathrm{K}$$, $$\Delta/k_B = 1, 300\ \mathrm{K}$$, $$G/k_B = 172\ \mathrm{K}$$, $$L/k_B = 290\ \mathrm{K}$$, and ln($$g$$) = 6.01.
6 $$ \times $$ 6	185.0	185.0	144.0	144.0
7 $$ \times $$ 7	187.75	187.75	148.6	148.6

Table 8

Values of the thermal hysteresis width for molecules located in the bulk and in the surface obtained by MCES, $$ \boldsymbol{\Delta T_{bulk}^{MCES}} $$ and $$ \boldsymbol{\Delta T_{surf}^{MCES}} $$, and by LMFA, $$ \boldsymbol{\Delta T_{bulk}^{LMFA}} $$ and $$ \boldsymbol{\Delta T_{surf}^{LMFA}} $$

$$ \boldsymbol{\mathrm{Lattice size}} $$	$$ \boldsymbol{T_{bulk}^{MCES} \mathrm{[K]}} $$	$$ \boldsymbol{T_{bulk}^{LMFA} \mathrm{[K]}} $$	$$ \boldsymbol{T_{surf}^{MCES} \mathrm{[K]}} $$	$$ \boldsymbol{T_{surf}^{LMFA} \mathrm{[K]}} $$
The computational parameters are: $$J_{bb}^{MCES}/k_B = 60\ \mathrm{K}$$, $$J_{ss}^{MCES}/k_B = 60\ \mathrm{K}$$, $$J_{bs}^{MCES}/k_B = 0\ \mathrm{K}$$, $$J_{bb}^{LMFA}/k_B = 40\ \mathrm{K}$$, $$J_{ss}^{LMFA}/k_B = 32.72\ \mathrm{K}$$, $$\Delta/k_B = 1, 300\ \mathrm{K}$$, $$G/k_B = 172\ \mathrm{K}$$, $$L/k_B = 290\ \mathrm{K}$$, and ln($$g$$) = 6.01.
6 $$ \times $$ 6	4.8	7	7.6	10.5
7 $$ \times $$ 7	10.4	14	2.8	5

Finally, in a larger particle of size 7 $$ \times $$ 7, the hysteresis associated with the bulk are wider than the hysteresis associated with the surface molecules by LMFA and by MCES. The proportion of molecules in the bulk ($$ \sim51 $$%) being greater than that of molecules in the surface ($$ \sim49 $$%), the bulk molecules drive the transition. In the 6 $$ \times $$ 6 compound, the hysteresis related to surface molecules was greatest in both MCES and LMFA. These results are summarized in Table 9.

Table 9

Values of the thermal hysteresis width for molecules located in the bulk and in the surface by MCES, $$ \boldsymbol{\Delta T_{bulk}^{MCES}} $$ and $$ \boldsymbol{\Delta T_{surf}^{MCES}} $$, and by LMFA, $$ \boldsymbol{\Delta T_{bulk}^{LMFA}} $$ and $$ \boldsymbol{\Delta T_{surf}^{LMFA}} $$, for a 6 $$ \times $$ 6 and a 7 $$ \times $$ 7 SCO particle

$$ \boldsymbol{\mathrm{Lattice size}} $$	$$ \boldsymbol{T_{bulk}^{MCES} \mathrm{[K]}} $$	$$ \boldsymbol{T_{bulk}^{LMFA} \mathrm{[K]}} $$	$$ \boldsymbol{T_{surf}^{MCES} \mathrm{[K]}} $$	$$ \boldsymbol{T_{surf}^{LMFA} \mathrm{[K]}} $$
The computational parameters are: $$J_{bb}^{MCES}/k_B = 60\ \mathrm{K}$$, $$J_{ss}^{MCES}/k_B = 60\ \mathrm{K}$$, $$J_{bs}^{MCES}/k_B = 60\ \mathrm{K}$$, $$J_{ss}^{LMFA}/k_B = 60\ \mathrm{K}$$, $$\Delta/k_B = 1, 300\ \mathrm{K}$$, $$G/k_B = 172\ \mathrm{K}$$, $$L/k_B = 290\ \mathrm{K}$$, and ln($$g$$) = 6.01.
6 $$ \times $$ 6	8.8	45	4	38.5
7 $$ \times $$ 7	23.2	48	0	32

The case $$ \boldsymbol{J_{bs}^{MCES} \neq 0} $$

The phase diagrams of the 7 $$ \times $$ 7 system obtained by LMFA and by MCES are plotted in Figure 13 and Figure 14, in the space coordinates, temperature vs. $$ J_{bs}^{MCES}/k_{B} $$, respectively, for the bulk and for the surface molecules.

Figure 13. Phase diagram T vs. ($$ J_{bs}^{MCES}/k_{B} $$) for a 7 $$ \times $$ 7 system. The red, blue and black squares correspond, respectively, to the upper and lower transitions ($$ T_{up} $$) and ($$ T_{down} $$) temperatures of the thermal HS fraction and to the equilibrium temperature of the bulk by LMFA. The red, blue and black stars correspond, respectively, to the upper and lower transitions ($$ T_{up} $$) and ($$ T_{down} $$) temperatures of the thermal HS fraction and to the equilibrium temperature of the bulk by MCES. The computational parameters are: $$ J_{bb}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{ss}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ \Delta/k_{B} = 1, 300\ \mathrm{K} $$, $$ G/k_{B} = 172\ \mathrm{K} $$, $$ L/k_{B} = 290 K $$, and ln($$ g $$) = 6.01.

Figure 14. Phase diagram T versus ($$ J_{bs}^{MCES}/k_{B} $$) for a 7 $$ \times $$ 7 system. The red, blue and black squares correspond, respectively, to the upper and lower transitions ($$ T_{up} $$) and ($$ T_{down} $$) temperatures of the thermal HS fraction and to the equilibrium temperature of the surface by LMFA. The red, blue and black stars correspond, respectively, to the upper and lower transitions ($$ T_{up} $$) and ($$ T_{down} $$) temperatures of the thermal HS fraction and to the equilibrium temperature of the surface by MCES. The computational parameters are: $$ J_{bb}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{ss}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ \Delta/k_{B} = 1, 300\ \mathrm{K} $$, $$ G/k_{B} = 172\ \mathrm{K} $$, $$ L/k_{B} = 290 K $$, and ln($$ g $$) = 6.01.

When the value of the interaction parameter $$ J_{bs}^{MCES}/k_{B} $$is increased up to 60 K, the three following items are especially significant.

First, for each method, the curves presented in Figure 12 and Figure 13 reveal that the equilibrium temperatures of the surface and of the bulk in the 7 $$ \times $$ 7 lattice are higher than those obtained in the 6 $$ \times $$ 6 lattice [Figure 10 and Figure 11].

In addition, as already observed in the case of a 6 $$ \times $$ 6 lattice and by the LMFA, when the short-range interaction $$ J_{bs}^{MCES}/k_{B} $$ is increased, the equilibrium temperature of the surface molecules decreases. We further note that the increase of the lattice size leads to a decrease of the width of the hysteresis in the surface. By MCES, the hysteresis transition observed for surface molecules disappears ($$ \Delta T_{surf}^{MCES} = 0 $$ K) in favor of an abrupt transition [Table 10].

Simulations using the LMFA were carried for very large lattices, 20 $$ \times $$ 20, 50 $$ \times $$ 50 and 200 $$ \times $$ 200. Figure 15 illustrates the results obtained for two sets of parameters corresponding to $$ J_{bs}^{MCES}/k_{B} = 0\ \mathrm{K} $$ and $$ J_{bs}^{MCES}/k_{B} = 60\ \mathrm{K} $$ and shows that, upon increasing the lattice size, the two-step transition vanishes. The curve closely resembles a single hysteresis loop for the 200 $$ \times $$ 200 lattice, and the equilibrium temperature tends towards that of the bulk {$$ T_{eq}= \Delta/[k_{B} \mathrm{ln}⁡(g) \approx 216\ \mathrm{K} $$}. The results obtained in the case of the 7 $$ \times $$ 7 lattice suggest that if we could do simulations on larger lattices with MCES, a single hysteresis loop with an equilibrium temperature close to 216 K would also be obtained.

Figure 15. Thermal dependence of the HS fraction studied with the LMFA for square-shaped particles of different sizes: 20 $$ \times $$ 20 (blue up triangles), 50 $$ \times $$ 50 (orange circles), and 200 $$ \times $$ 200 (green stars). The computational parameters are: (A) $$ J_{bb}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{ss}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{bs}^{MCES}/k_{B} = 0\ \mathrm{K} $$, $$ J_{bb}^{LMFA}/k_{B} = 36\ \mathrm{K} $$, $$ J_{ss}^{LMFA}/k_{B} = 33.33\ \mathrm{K} $$; (B) $$ J_{bb}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{ss}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{bs}^{MCES}/k_{B} = 60\ \mathrm{K} $$, $$ J_{bb}^{LMFA}/k_{B} = 60\ \mathrm{K} $$, $$ J_{ss}^{LMFA}/k_{B} = 60\ \mathrm{K} $$. For all calculations: $$ \Delta/k_{B} = 1, 300\ \mathrm{K} $$, $$ G/k_{B} = 172\ \mathrm{K} $$, $$ L/k_{B} = 290\ \mathrm{K} $$, and ln($$ g $$) = 6.01.

CONCLUSION

The present contribution, a simulation work, with the purpose of giving some hints to colleagues who synthesize spin-transition nanoparticles, is focused on meticulous investigations of the thermal properties of 2D-SCO nanoparticles using an adapted two-state Ising-like Hamiltonian. At this end, three types of interaction terms are considered: long and short-range interactions originating from electronic and elastic molecular couplings inside the lattice, to which we added a specific energetic contribution accounting for the interactions between surface molecules with their surrounding environment (inside and outside the lattice). The model is solved in two distinct ways, which involve different levels of sophistication: (ⅰ) MCES techniques, which are exact methods and (ⅱ) LMFA based on mean-field theory but taking into account the lattice shape. In the former case, three types of two-body short-range interaction terms are involved between bulk-bulk ($$ J _{bb} $$), surface-surface ($$ J_{ss} $$), and bulk-surface ($$ J_{bs} $$) sites. The thermal dependence of the "order parameter", $$ <\sigma> $$, that is, the average fictitious magnetization, is calculated by using the Boltzmann distribution obtained from the calculated density of the macrostates. On the other hand, in the LMFA, two interaction terms corresponding to three types of sites, corner, edge, and bulk, are included in the calculations, taking into account the dependence of the coordination number on the site location. The results of the two methods have been compared for different square lattice sizes, and a good quantitative agreement is obtained between them in the case $$ J_{bs} $$ = 0, where both equilibrium temperatures and hysteresis widths were found in very good agreement. When the bulk-surface interaction $$ J_{bs} $$ is increased, the hysteresis observed by MCES in the bulk and in the surface tend to be close to each other and finally overlap for strong Jbs values. In contrast, by LMFA, this phenomenon is not observed, and furthermore, the hysteresis effects are accentuated, particularly for surface molecules. Moreover, with the MCES method, the increase of $$ J_{bs} $$ results in a shift of the surface (resp. bulk) transition temperatures to higher (resp. lower) values, while in LMFA, an opposite effect is obtained for surface molecules.

Our comparison of the two techniques (MCES and LMFA) shows that when dealing with small samples, LMFA yields results that are in perfect agreement with MCES in the case where the interaction between the surface and the bulk is zero ($$ J_{bs} $$ = 0). However, as the interaction $$ J_{bs} $$ increases, causing important interferences between the two regions of the lattice, the results provided by LMFA increasingly differ from those obtained with the MCES method, although the general trends remain the same. While considering a very weak interaction between the surface and the bulk in the same compound may look unnatural, it may be possible in nanostructured core-shell nanoparticles, where the core and the shell are of different compositions. In addition, it is important to keep in mind that the systems considered are particularly small, and as the size of the system increases, the predominance of the surface, and then $$ J_{bs} $$ interaction, vanishes, which, in turn, means that LMFA will most likely still provide correct results given the large number of molecules.

Finally, it is worth mentioning that the surface properties have been limited here to the outer layer of the lattice, while the concept of the surface may integrate some inner layers close to the interface material/air, which also have specific properties different from those of the bulk. Extensions of the present work to include the above-cited aspects are in progress, as well as other packing clusters such as rhombohedral or faced centered cubic (FCC). The results of this simulation work can also be extended to explain the behavior of polycrystalline SCO films on the $$ \mathrm{Al_{2}O_{3}} $$ substrate^[40] and the influence of the substrate on the functionality of SCO molecular materials^[41].

DECLARATIONS