Surface Effects Leading to Unusual Size Dependence of the Thermal Hysteresis Behavior in Spin-Crossover Nanoparticles

Abstract

We analyze the size effect on spin-crossover transition nanoparticles in a 2D Ising-like model subject to a specific ligand-field at the surface. By anisotropic sampling method applied to the finite 2D square Ising lattices with various sizes, we determined the density of macro states by scanning the spin configurations. This information, which is independent on the system parameters, is used to exactly calculate the thermal behavior of spin-crossover nanoparticles whose ligand-field of the atoms at the surface is lower than those of the bulk. We found that decreasing the size of the nanoparticles leads to a global increase of the effective interaction, which has the consequence to enhance the width of the thermal hysteresis. This unusual behavior opens a new avenue in controlling the bistability characteristics at small scale, one of the important conditions of applicability of these materials at the nanometric scale.

1. Introduction

In recent years, there has been a growing interest in the field of spin-crossover (SCO) solids. Under various constraints, such as temperature [1] or pressure [2,3,4] variations, irradiation by visible light [5,6,7,8,9] or magnetic fields [10], SCO solids can be switched from a low-spin (LS) state to a high-spin (HS) state, and conversely [11,12]. For example, Fe(II) SCO compounds [13] are diamagnetic (S = 0, LS) and paramagnetic (S = 2, HS) in the low- and high-temperature phases, respectively. The spin-crossover transition involves both electronic transformation (spin and orbital), and structural modifications. In the case of Fe(II) compounds the metal-ligand bond lengths change by about 0.2 Å, (~10%) and the ligand-metal-ligand bond angles by 0.5°–8° [6,14]. Optical properties are also changed [5,15,16], with usually a colorless HS state and a strongly-colored LS state. Consequently, magnetic and optical measurements [17,18,19,20,21] are the major experimental techniques used for quantitative investigations of the spin transitions. In many cases, elastic interactions between the SCO units are strong enough to induce hysteresis at the thermal spin transition [22] which occurs as a first-order phase transition. Such “switchable” molecular solids are promising in terms of optical data storage [23].

For a few years, the synthesis and the design of spin-crossover nanoparticles have become possible and several interesting behaviors have been reported [24,25,26,27,28,29,30,31,32,33,34,35]. The experimental results, for most of these studies, confirmed the expected narrowing of the hysteresis loop and, in some cases, the collapse of the loop below some “critical size”, but very few of them provided a coherent set of data over the relevant size range, that is, above and below the critical size. In addition, the samples usually present a rather wide distribution of sizes, various shapes, and particles were often aggregated in various manners. All of these features obviously impact the switching behavior of the SCO nanoparticles and made the investigation on the properties of a unique nanoparticle a very challenging goal. From the theoretical point of view, there were a few studies of the thermal and photo-switching properties of SCO nanoparticles [36,37,38,39], which aimed at reproducing the collapse of the hysteresis effects at the critical size of the particles. These theoretical studies were based on Metropolis simulations on finite size lattices with free boundary conditions. Here, we introduce an additional aspect, which is the specific electronic state of the atoms located at the surface of the nanoparticle. This was suggested by a recent experimental investigation [26], on SCO nano-particles of Fe(pyrazine)[Pt(CN)₄], with well-controlled sizes above and below the critical size. In addition to the expected narrowing of the hysteresis loop upon decreasing size, this system displayed a sizable lowering of the transition temperature and an increase of the residual HS fraction. The latter features are explained by the different coordination of the surface atoms, including water molecules instead of the organic ligands. These SCO surface atoms experience a weaker ligand-field capable of trapping them in the HS state, which consequently lowers the equilibrium temperature of the particle through a negative pressure effect. For practical reasons the present investigation was performed on 2-D networks where this specific effect was located at the edges. It will be consequently referred to in terms of “edge effect”. Although the present 2-D model might be considered as a preliminary approach to future 3-D models, it is also relevant due to the particular geometry of the nano-crystals which are in shape of square platelets with a rather constant thickness. The core of the square lattice is already submitted to a quasi-uniform effect from both surfaces, irrespective of the size of the network, and the size effect effectively reduces to a 2-D problem. This viewpoint will be supported by the present simulations. We performed Monte Carlo (MC) simulations using the well-known Ising-like Hamiltonian (interacting two-level systems), widely used in literature to describe the equilibrium [40,41] and non-equilibrium [42] properties of SCO solids. This choice was dictated by the simplicity of this description, despite of its phenomenological character. The paper is organized as follows: Section 2 is devoted to the model; Section 3 to the results and their discussion; Section 4 is dedicated to a detailed comparison to experimental data, including the impact of size distributions, and Section 5, to our conclusion.

2. Results

In order to reproduce first-order phase transition in 1D-SCO, Linares et al. [43] have introduced a short- and a long-range interaction to the Ising-like model [40,41]. Later on, and in order to take into account the matrix effect, Chiruta et al. [44] have added an additional interaction term which accounts for the SCO molecules at the surface with the environment (matrix effect), this is termed, the matrix interaction L.

Thus, the system’s Hamiltonian can be written as follows:

$$H=\frac{\Delta -k_{B}T\ln g}{2}{\displaystyle {\sum }_{i=1}^{N_t}{\mathsf{\sigma }}_i-J{\displaystyle {\sum }_{<i,j>}{\mathsf{\sigma }}_i{\mathsf{\sigma }}_j}}-G<\mathsf{\sigma }>{\displaystyle {\sum }_{i=1}^{N_t}{\mathsf{\sigma }}_i-L{\displaystyle {\sum }_{i=1}^{N_s}{\mathsf{\sigma }}_i}}$$

(1)

In Equation (1), the HS and LS states are associated, respectively, with the eigenvalues +1 and −1 of the fictitious spin operator, σ, with their respective different degeneracies g+ and g−, the ratio of which, is denoted, $g=\frac{g_{+}}{g_{-}}$. These degeneracies account for both internal (spin, orbit, and intermolecular vibrations) and external (lattice phonons) degrees of freedom. Here, Δ is the energy difference E(HS) ─ E(LS) between the spin states of isolated molecules, and J is the coupling parameter accounting for short-range ferroelastic interactions between spin states, here limited to first-neighboring molecules; and G is the long-range part of the interaction, which is identified as due to lattice phonons and volume change accompanying the transition. The energetic contribution L (>0), limited to surface atoms, describes the additional “negative” ligand-field felt by the surface atoms. This contribution may arise from coordinated water molecules which substitute to one or several nitrogen atoms linked to the Fe center, or may also be due to the elastic effect of the surrounding polymeric matrix (with bigger lattice parameter than that of the embedded SCO nanoparticle). Let us denote by, $N_t=N_x^2$, the total number of atoms, by $N_c={\left(N_x-2\right)}^2$ the number of atoms in the bulk or the core, and by $N_s=4\left(N_x-1\right),$ the number of atoms at the "surface" (edges), $N_x \times N_x$ square lattice.

Hamiltonian Equation (1) can be re-expressed as a function of the dimensionless macroscopic variables, m, s and c, given by:

$$m={\displaystyle {\sum }_{i=1}^{N_t}{\mathsf{\sigma }}_i}$$

(2)

$$s={\displaystyle {\sum }_{<i,j>}{\mathsf{\sigma }}_i{\mathsf{\sigma }}_j}$$

(3)

and:

$$c={\displaystyle {\sum }_{i=1}^{N_s}{\mathsf{\sigma }}_i}$$

(4)

which are, respectively, associated with the total fictitious magnetization, two spin correlations, and the surface fictitious magnetization. Within these parameters, the total Hamiltonian of Equation (1) is re-written as:

$$H=\left(\frac{\Delta -k_{B}T\ln g}{2}-G<\mathsf{\sigma }>\right)m-Js-Lc$$

(5)

The associated bulk and surface order parameters of the system, denoted here $m_V$ and m_s, are then simply defined as:

$$m_V=\langle {\mathsf{\sigma }}_V\rangle =\frac{{\displaystyle {\sum }_{i=1}^{N_t}{\mathsf{\sigma }}_i}}{N_t}$$

(6)

and:

$$m_s=\langle {\mathsf{\sigma }}_s\rangle =\frac{{\displaystyle {\sum }_{i=1}^{N_s}{\mathsf{\sigma }}_i}}{N_t}$$

(7)

The HS fractions, n_s (surface), and $n_V$ (volume), connect to these two quantities through the general relation:

$$n_i=\frac{1+m_i}{2}(i=s,\text{ v})$$

(8)

The equilibrium temperature T_eq is defined as the temperature for which the average spin value <σ> is equal to zero ($m_V$ = 0), which leads to a HS fraction equal to 1/2. The present model involves two local situations, the features of which are listed in

Table 1. The ligand field as function of the molecules' positions in the lattice.

Description	Total	Edges	Core
Number of atoms	$N_x^2$	4(Nx − 1)	(Nx − 2)2
Ligand-field	-	$\frac{\Delta -2L-k_{B}T\ln g}{2}$	$\frac{\Delta -k_{B}T\ln g}{2}$

. The key point is that the specific weaker effective ligand fields (Δ − 2L) at the surface act as a negative (but non uniform) pressure on the system through the interaction of surface atoms with their nearest neighboring bulky atoms. As a matter of fact, these contributions should be included in the total ligand-field.

Hamiltonian Equation (1) was exactly solved in the canonical approach. The parameter values were chosen from typical data in spin-crossover literature, such as the molar entropy change ΔS ≈ 50 J/K/mol leading to ln (g) = ΔS/R ≈ 6.01 (where R is the gas constant) and energy gap Δ/k_B = 1300 K leading to a typical value of the equilibrium temperature T_eq = Δ/(k_B ln(g)) ≈ 216.3 K. The short-range interaction parameter, J, is taken equal to J/k_B = 15 K and the long-range parameter, G, is taken equal to G/k_B = 172.7 K. The interaction between the edge molecules and the environment, L, is taken equal to L/k_B = 120 K.

The exact calculations were performed on square lattices of size $N_x\times N_x$ up to $N_x = 12$, using open boundary conditions in order to fairly take into account for the experimental situation of SCO nanoparticles having a specific ligand field at the surface. Next, the thermal dependence of the HS fraction, n_HS(T), for several lattice sizes were determined. The results are summarized in Figure 1. The system was warmed from T = 175 K to 225 K and cooled down to the initial temperature 175 K by 1 K steps. The total HS fraction, n_HS(T), was obtained from Equation (8) where:

$$m=<\mathsf{\sigma }>=\frac{{\displaystyle \sum _{i=1,NL}\frac{m_i}{N}d(m_i{s}_i{c}_i)\exp \left(-\frac{1}{k_{B}T}(-h_f{m}_i-J{s}_i-L{c}_i)\right)}}{{\displaystyle \sum _{i=1,NL}d(m_i{s}_i{c}_i)\exp \left(-\frac{1}{k_{B}T}(-h_f{m}_i-J{s}_i-L{c}_i)\right)}}$$

(9)

and d(m,s,c) is the number of configurations for a given set of values, calculated using Monte Carlo Entropic Sampling (MCES) method, described in the next section. NL is the number of distinct configurations of states m, s, c (m, s, and c defined by the Equations (2)–(4)) where the “field” h_f has the expression:

$$h_f=-\left(\frac{\Delta -k_{B}T\ln g}{2}-G\langle \mathsf{\sigma }\rangle \right)$$

(10)

Similarly, the average of the edge molecules is given by the following equation:

$$<c>=\frac{{\displaystyle \sum _{i=1,NL}\frac{c_i}{N_c}d(m_i{s}_i{c}_i)\exp \left(-\frac{1}{k_{B}T}(-h_f{m}_i-J{s}_i-L{c}_i)\right)}}{{\displaystyle \sum _{i=1,NL}d(m_i{s}_i{c}_i)\exp \left(-\frac{1}{k_{B}T}(-h_f{m}_i-J{s}_i-L{c}_i)\right)}}$$

(11)

3. Monte Carlo Entropic Sampling

We used the Monte Carlo entropic sampling method [45,46] to generate the table that contains the dimensionless macroscopic variables, m, s, c, and their density, d(m,s,c). Therefore, we used the principle of MCES described by Shteto et al. [45] which consists in introducing in the detailed balance equation of the MC procedure, a suited biased distribution P in order to favor configurations belonging to weakly degenerate macrostates and to dampen those belonging to the highly degenerate macrostates.

The balance equation is given by:

$$P_{i}W(i\to j)=P_{j}W(j\to i)$$

(12)

and the biasing probability was chosen as the inverse of the desired restricted density of states.

$$P_i\propto \frac{1}{d(m,s,c)}$$

(13)

In this case, the balance equation can be written:

$$\frac{W(i\to j)}{W(j\to i)}=\frac{P_j}{P_i}=\frac{d(m_i,s_i,c_i)}{d(m_j,s_j,c_j)}$$

(14)

Since, in the first Monte Carlo step, the density of the state d(m,s,c) is unknown we put all d(m,s,c) equal to 1. So, after iteration i the density will be d_i(m,s,c). Then, using d_i(m,s,c) as a bias, a MC sampling is run; it is termed a “Monte Carlo stage” and yields a histogram of the frequency of the macrostates [46], H_i(m,s,c):

$$H_i(m,s,c)\propto d(m,s,c)\frac{1}{d_i(m,s,c)}$$

(15)

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

$$d_{i+1}(m,s,c)\propto d_i(m,s,c)H_i(m,s,c)$$

(16)

After the table that consist from m, s, c, and d(m,s,c) is obtained, the exact partition function can be calculated using the following expression:

$$Z={\displaystyle \sum _{(m,s,c)}d(m,s,c)\exp (-\mathsf{\beta }(-hm-Js-Lc))}$$

(17)

and so all the thermodynamic properties of the systems can be derived analytically.

4. Results and Discussion

Figure 1 illustrates the thermal behavior of the HS fraction, n_HS(T), calculated for various particle sizes. Upon decreasing particle size $\left(N_x\times N_x\right)$ the transition temperature is downward shifted and the width of the thermal hysteresis loop progressively increases as a result of surface effects.

The detailed discussion of the results of the model follows. This behavior contrasts with the classical behavior of the thermal hysteresis with particle size, where indeed it monotonously vanishes at small sizes. To help the reader in identifying the stability of all solutions drawn in the thermal hysteresis, we highlighted in Figure 1b, the stable, metastable and unstable states accompanying the thermally-induced first-order spin-crossover transition. The stable and metastable states are easily identified: HS (LS) is stable for all temperatures bigger (lower) than T_eq and becomes metastable for T < T_eq (T > T_eq). For all of these branches, increasing the temperature results in the increase of the HS fraction: that is a normal behavior. However, in the middle region of the thermal hysteresis (blue dots in Figure 1b), the HS fraction becomes a decreasing function with temperature ($\frac{\partial n_{HS}}{\partial T}<0$) which obviously corresponds to an instability. The situation is somehow similar to the mechanical instability region of a van der Waals gas, where the volume increases under a positively applied pressure.

We could separate the contributions of surface atoms and core atoms in the thermal behavior of the total HS fraction of Figure 1. Figure 2 displays two chosen cases of nanoparticles with sizes (4 × 4) and (12 × 12) depicting all contributions. One can clearly see, that in both cases, the surface (blue curves) starts to transform much earlier than the core (green curves), which clearly indicates that the surface is driving the thermal LS to HS transitions. It is however obvious that for very bigger sizes (unreachable with MCES method), the contribution of the surface will be marginal and the core will dominate. This tendency is already visible in the case of 12 × 12 particle, in which we remark that the curve of the total response (red) is much more close to that of the bulky atoms.

4.1. Size Dependence of the Equilibrium Temperature

We found that the size dependence of the transition temperature could be described by analytical laws (easily extended to 3D models). Here, it worth mentioning that even for small sizes, exact analytical solutions of Hamiltonian Equation (1) in 2D are out of reach [47], contrary to our previous studied case where the surface atoms were fixed in the HS state [48].

We introduce here a simple idea that the transition temperature of the system is still the result of a null total effective ligand-field. Taking advantage of the data listed in Table 1, the transition temperature, T_eq, for a square lattice N_x × N_x is the solution of the following equation:

$$\frac{\Delta -k_B{T}_{eq}\ln g}{2}\times {\left(N_x-2\right)}^2+\frac{\Delta -2L-k_B{T}_{eq}\ln g}{2}\times 4\left(N_x-1\right)=0$$

(18)

which expresses that the total ligand-field is equal to zero at T = T_eq.

Solving Equation (18) gives the analytical expression of the lattice size-dependence of the transition temperature, as:

$$T_{eq}=\frac{N_c}{N_t}{T}_{eq}^{bulk}+\frac{N_s}{N_t}{T}_{eq}^{surf},$$

(19)

which can be written as:

$$N_x^2{T}_e{}_q={(N_x-2)}^2{T}_e{}_q^{bulk}+4(N_x-1)T_e{}_q^{surf},$$

(20)

leading to a parabolic and linear contribution of the bulk and surface contribution to $N_x^2{T}_{eq}$. These tendencies have been checked numerically in the simulations, and shown in Figure 3, which indicates that from the phenomenological view point, they are followed.

where $T_{eq}^{bulk}$ and $T_{eq}^{surf}$ are the transition temperatures of the bulk and the surface, whose expressions are:

$$T_{eq}^{bulk}=\frac{\Delta }{k_B\ln g}\approx 216.3K\text{ and }{T}_{eq}^{surf}=\frac{\Delta -2L}{k_B\ln g}\approx 176.3K$$

(21)

Equation (19) predicts that the transition temperature has two limiting values, namely $T_{eq}=T_{eq}^{surf}=173.6$ for $N_x=2$ for the smallest nanoparticle size, and $T_{eq}=\frac{\Delta }{k_B\ln g}\approx 216\mathrm{K}$ corresponding to the bulk transition temperature, reached for an infinite lattice ($N_x\to \infty$).

The relevance of Equation (19) is supported by the comparison to the exact results derived from entropic sampling method. As shown in Figure 4, where the global $T_{eq}$ is plotted vs $1/N_x$, the analytical prediction (blue dots) given by Equation (19) is in excellent agreement with the data arising from the simulation (red squares).

We should mention that the present model involves the assumption that the effective degeneracy $g$ is constant. Thus, $g$ is assumed as independent on the strain at the surface and as well it does not depend on the particle size.

4.2. The Size Dependence of the Hysteresis Loop and the Occurrence of the First-Order Transition

To understand the size effects on the thermal behavior of the HS fraction, we performed MC simulations using the set of parameters given in Figure 1, for which the studied system presents a gradual transition for the biggest studied size (i.e., $N_x=12$). For that, the short-range interaction parameter was set to J/k_B = 15 K. It is worth mentioning that in the usual Ising-like model, the first-order phase transition takes place only when the condition T_O.D. > T_eq is satisfied. Here, T_O.D. is the order-disorder (or Curie) temperature of the corresponding pure Ising model, obtained from Hamiltonian of Equation (1) by putting Δ/k_B = 0, L/k_B = 0 and g = 1 (ln(g) = 0). For J = 0 values, the interactions between SCO sites are of pure long-range nature and so T_O.D = G/k_B~172.7 K (minimum value). For a large lattice size, the equilibrium temperature is dominated by the contribution of the bulk material and leads here to $T_{eq} ~ T_{eq}^{bulk}~216$ K, which is larger than T_OD, which allows to say that for infinite lattice, the thermal behavior of the HS fraction is that of a gradual conversion. When $J\ne 0$, the analytical expression of T_OD is out of reach. So, to determine T_OD for different lattice sizes, we performed MC simulations for Δ/k_B = 0 K, G/k_B = 172.7 K, J/k_B = 15 K, L/k_B = 0 K, and ln(g) = 0. The results are summarized in Figure 5, which shows a net increase of T_OD with the system size, following a simple power law, $T_{OD}~\sqrt{N_x-4}$, as demonstrated in Figure 5c.

The values of T_OD vary here from ~215 K for a 4 × 4 lattice to ~224 K for a 12 × 12 lattice. On the other hand, we found from Figure 1, that T_eq ~204 K (<T_O.D = 224) for the lattice 12 × 12. Consequently, according the simple criteria of the simple Ising-like model, the system is expected to show a first-order transition, while a gradual conversion is obtained in the simulations. However, for the lattice 4 × 4, we obtained a first-order transition with a hysteresis at T_eq~186 K, which satisfies the usual criteria of occurrence of first-order transitions in Ising-like model. Figure 6 showing the behavior of the thermal hysteresis width versus the system sizes summarizes the results. The system undergoes a gradual SCO transition for sizes larger than $N_x=6$, and shows first order transitions below this critical value. Interestingly, $N_x=6$ corresponds to the value at which we observe between the number of particles in the bulk ${\left(N_x-2\right)}^2$ and that at the surface, $4\left(N_x-1\right)$. This fact demonstrates that the surface is playing the major role in driving the spin transition in this problem, as observed in experiments. The competition between the core and surface molecules are driven by the long-range, G, and the surface, L, interactions, and then the critical size to obtain the first-order transition is a function of the rate of these interaction parameters.

5. Conclusions

The impact of specific conditions at the edges of SCO nanoparticles have been investigated in a 2-D model in which a smaller ligand field was affected to surface atoms, based on the experimental observation of a decreasing transition temperature and an increasing thermal hysteresis width upon particle size reduction. This specific boundary condition basically acts as a negative pressure, shifting downwards the equilibrium temperature, however, in a size-dependent way. The interplay between this equilibrium temperature variation and the expected Curie temperature, arising from a pure Ising model, explains the enhancement of the width of the thermal hysteresis loop upon size reduction.

Further aspects related to the nanoparticle problem remain to be explored: the existence of a re-entrant first-order spin transition on decreasing the particle size will be analyzed using the same ideas developed throughout the manuscript.