<figure><img src="image_1.png" /></figure>

HAL
open science

# Electro-Elastic Modeling of Thermal Spin Transition in Diluted Spin-Crossover Single Crystals

Karim Affes, Yogendra Singh, Kamel Boukheddaden

<figure><img src="image_2.png" /></figure>

## To cite this version:

Karim Affes, Yogendra Singh, Kamel Boukheddaden. Electro-Elastic Modeling of Thermal Spin Transition in Diluted Spin-Crossover Single Crystals. International Journal of Molecular Sciences, 2022, 23 (22), 10.3390/ijms232213854. doi:10.3888/ijms.232213854

HAL Id: hal-03888936

https://hal.science/hal-03888936

Submitted on 21 Dec 2022

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_3.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_4.png" /></figure>

Article

# Electro-Elastic Modeling of Thermal Spin Transition in Diluted Spin-Crossover Single Crystals

Karim Affes <fnref n="1" />,<fnref n="2" />, Yogendra Singh <fnref n="1" />,† and Kamel Boukheddaden <fnref n="1" />,*

**Abstract:** Spin-crossover solids have been studied for many years for their promising applications as optical switches and reversible high-density memories for information storage. This study reports the effect of random metal dilution on the thermal and structural properties of a spin-crossover single crystal. The analysis is performed on a 2D rectangular lattice using an electro-elastic model. The lattice is made of sites that can switch thermally between the low-spin and high-spin states, accompanied by local volume changes. The model is solved by Monte Carlo simulations, running on the spin states and atomic positions of this compressible 2D lattice. A detailed analysis of metal dilution on the magneto-structural properties allows us to address the following issues: (i) at low dilution rates, the transition is of the first order; (ii) increasing the concentration of dopant results in a decrease in cooperativity and leads to gradual transformations above a threshold concentration, while incomplete spin transitions are obtained for big dopant sizes. The effects of the metal dilution on the spatiotemporal aspects of the spin transition along the thermal transition and on the low-temperature relaxation of the photo-induced metastable high-spin states are also studied. Significant changes in the organization of the spin states are observed for the thermal transition, where the single-domain nucleation caused by the long-range elastic interactions is replaced by a multi-droplet nucleation. As to the issue of the relaxation curves: their shape transforms from a sigmoidal shape, characteristic of strong cooperative systems, into stretched exponentials for high dilution rates, which is the signature of a disordered system.

**Keywords:** spin transition; electro-elastic modelling; Monte Carlo simulations; diluted molecular materials; thermal hysteresis; spatiotemporal properties

<figure><img src="image_5.png" /></figure>

check for
updates

**Citation:** Affes, K.; Singh, Y.; Boukheddaden, K. Electro-Elastic Modeling of Thermal Spin Transition in Diluted Spin-Crossover Single Crystals. *Int. J. Mol. Sci.* **2022**, *23*, 13854. https://doi.org/10.3390/ijms232213854

Academic Editors: Mateusz Reczyński and Szymon Chorazy

Received: 16 October 2022

Accepted: 2 November 2022

Published: 10 November 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_6.png" /></figure>

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

## 1. Introduction

One of the most intriguing examples of molecular bistability is the ability of some metal complexes to switch between low-spin (LS) and high-spin (HS) states. Among these functional materials, the spin-crossover (SCO) molecular solids are excellent candidates for reversible memory devices [1–3], sensors of temperature, pressure, and gases [4]. They can also serve as switches of light emission in the electroluminescent devices containing spin-crossover (SCO) complexes [5,6], or to enhance plasmonic resonances [7].

The case of iron(II)-based SCO materials, with 3d⁶ configuration, is the most common in the literature. In these ferrous SCO compounds, the five d-orbitals of iron(II) are then split into three bonding, t₂g and two anti-bonding, e_g, orbitals. When the strength of the ligand field, Δ = *E*(*t*₂g) − *E*(*e*_g), is much stronger (resp. weaker) than the electrons’ pairing energy, the metal ion can be in the LS *t*₂g*e*_g⁰ diamagnetic (resp. HS *t*₂g*e*_g<fnref n="2" /> paramagnetic) state with a total spin *S* = 0 (resp. *S* = 2). In Fe-based SCO materials, the spin transition between LS and HS is accompanied by an increase of ~10% of the Fe-ligand bond lengths

<footnote n="1">¹ GEMaC, CNRS, Université de Versailles Saint Quentin en Yvelines, 45 Avenue des Etats Unis, CEDEX, 78035 Versailles, France</footnote>

<footnote n="2">² Laboratoire des Matériaux Multifonctionnels et Applications, Faculté des Sciences de Sfax, Université de Sfax, Route de la Soukra km 3.5, Sfax 3000, Tunisia</footnote>

<footnote>* Correspondence: kamel.boukheddaden@uvsq.fr</footnote>

<footnote>† In memory of our regretted friend and colleague Dr. Yogendra Singh who passed away recently at the age of 29.</footnote>
due to the weakening of the metal–ligand bonds. This causes the changes in local volume of ~3–5% [8,9], as well as ligand metal–ligand bond angles by 0.5–80°, affecting their magnetic, vibrational and optical [10] properties. Nowadays, a large panel of experimental techniques is used to detect the spin transition. Among them, one can cite magnetic [11–13] and optical [14,15], Mössbauer spectroscopy [16], optical microscopy (OM) [17–23], X-ray diffraction [24], differential scanning calorimetry [25], diffuse reflectivity [26], etc. OM, in particular, allowed to investigate the spatiotemporal features of the spin transition on a unique single crystal. This enabled the imaging of the real-time transformation of the crystals along the thermal spin transition, with remarkable spatial (~0.3 µm) and temporal (up to 0.01 s) resolutions. External stimuli, such as pressure [5], temperature [27], electric field [28], magnetic field [29] and light [30], can also cause significant effects on the spin transition mechanism.

In crystalline state, the thermal properties of SCO materials can lead to a wide range of different thermal behaviors. For example, we mention the case of (i) first-order thermal spin transition (cooperative behavior), (ii) gradual spin transition [31] and (iii) multistep or incomplete spin transition [32]. In the latter case, a spatial self-organization of the spin states [33,34] is sometimes observed in the plateau region.

The elastic interactions between the SCO molecules are recognized as an essential ingredient in the origin of their cooperative character. The interactions between the SCO units are mainly caused by the volume differences between HS and LS molecules. If this volume difference is strong enough and the SCO molecules are well linked, then thermally induced first-order spin transitions take place. The cooperative nature of the solid results from combined short- and long-range interactions. Here, the short-range interactions relate to the shape and the distance between neighboring molecules. In contrast, the long-range interactions are caused by the local volume changes, which delocalize in the whole crystal through acoustic phonons [35–37].

Understanding the source of the interactions is critical, then, to constructing SCO materials with a high degree of cooperativity. The microscopic mechanism of interaction has been extensively studied, and various theoretical models [38] have been proposed to explain the experimental findings in a contextually sound manner. Among these models, we quote for example, (i) the Ising-like model [39], also used to describe SCO nanoparticles, (ii) the atom–phonon model [40,41], and the (iii) elastic models [42–54], which have been recently introduced. While the Ising models suffer from their qualitative nature and a lack of structural information, the elastic model accounting for the volume change at the transition appears much more realistic and has attracted a lot of interest these last years.

The study of metal ion dilution in solids is a typical experimental approach for elucidating the function of these elastic interactions. The mechanism and the range of the interactions between spin crossover units may be then studied by diluting the spin crossover metal centers. Within this technique, it is possible to regulate both the number of interconnecting centers and the intensity of those interconnections. The SCO molecules become increasingly isolated as they are randomly substituted by non-active impurities, disrupting the lattice cooperativity leading to more gradual transitions.

Previous experimental studies of metal dilution on SCO solids [55] have been performed on two-step SCO compounds [56] as well as on one-step SCO compounds, including hysteretic ones [57–61]. Diluting SCO complexes with isostructural but magnetically silent impurities has led to considerable information on the link between intermolecular contacts and the degree of cooperativity in the solid state. For example, the mixed crystal series [FeₓZn₁₋ₓ(2 - pic₃)Cl₂EtOH (x = 0.007 to 1; 2-pic = 2-picolylamine) was the first series for which this approach was introduced [62], and the SCO curve first had a sharp step but grew more gradual as the Zn/Fe ratio rose. It also shows, in this material, that the dilution process energetically favors the HS state by lowering the transition temperature. Another study [63] reporting on the effects of zinc dilution concerned the mononuclear two-step spin crossover compounds [FeₓZn₁₋ₓ(bapbpy)(NCS)₂ (bapbpy = N₆, N'₆ - di(pyridin-2-yl) - 2,2' - bipyridine - 6 - 6' - diamine). It
demonstrated that increasing the zinc concentration results in the vanishing of the cooperativity at the threshold value $x = 0.24$, with a decrease in the transition temperature. The role of cooperative elastic interactions between the SCO metal centers were validated by these findings, but not yet clearly explained. Other dilution investigations have shown similar results [62,64,65].

The dopant ion size profoundly affects the SCO characteristics. One of the most important reasons for this is that ions with a size larger than that of the HS Fe(II) cause a negative internal pressure in the LS lattice, thus favoring the HS state. When diluting Fe-based SCO compounds, alternative metal cations (Ni$^{2+}$, Co$^{2+}$, or Mn$^{2+}$), that are not magnetically quiet, are also used to investigate the effects of different-sized metal ions [64,66]. Some of them, like Ni, have a neutral role because their size is located between that of Fe(II)-HS and Fe(II)-LS [55,61,65].

Another aspect of the metal dilution relates to its effect on the relaxation properties of the HS photo-induced metastable state. In this context, it is important to mention that A. Hauser and co-authors who discovered the LIESST (light-induced excited spin state trapping) effect [67] have also been the precursors in dynamical studies of diluted SCO solids. The influence of metal dilution on the light conversion [58,68–73] was studied for both excitation and relaxation processes. When iron(II) ions are diluted with a non SCO-active dopant, the HS-LS relaxation behavior changes from a self-accelerated to an exponential or stretched process owing to the decrease in cooperativity and to the inhomogeneous character of the lattice. On the other hand, from the thermodynamic point of view, the spin equilibrium temperature, $T_{eq}$, defined as the temperature for which $G_{HS} - G_{LS} = 0$, where G is Gibbs free energy, was shown to be significantly influenced by metal dilution. Correlatively, for hysteretic spin transitions, the hysteresis loop’s width declined and eventually vanished with increasing dilutions.

The additional insight gained through these studies will be a fundamental step toward the explanation of the effect of dilution on the structural and thermal properties of the SCO materials. The primary goal of the present study is then to investigate the theoretical consequences of metal dilution on SCO systems. The investigations are made using the electro-elastic model, which has the merit of combining both electronic and elastic degrees of freedom of SCO phenomena. This model is extended here to include the effect of contaminants on SCO complexes, in which the interactions are modeled using springs.

The manuscript is organized as follows: Section 2 presents an extended version electro-elastic model adapted for diluted SCO materials and describes the simulation procedure. Section 3 reports the effects of dilution on the thermal spin transition and low-temperature relaxation of photo-induced metastable HS fractions. The effect of the dilution on the spatiotemporal aspects of the spin transition is also discussed. In Section 4, we conclude and outline the possible developments of the present work.

## 2. The Model and Simulations Procedure

We studied a distortable 2D square lattice with a fixed topology and open boundary conditions. The crystal is constituted by SCO units and dopants, also hereafter called defects. Each lattice site is associated with a fictitious spin, $S_i$, whose values are $S_i = -1$, for SCO site in LS state, $S_i = +1$ for a SCO site in HS and $S_i = 0$ for a site occupied by a dopant (defect).

The total Hamiltonian of the system, taking into account for the spin states and lattice positions, is written:

$$H = \sum_{i} (\Delta - k_{B} T \ln g) S_{i} + A \sum_{(i,j)} (r_{ij} - R(S_{i}, S_{j}))^{2} + B \sum_{(i,k)} (r_{ik} - R'(S_{i}, S_{k}))^{2} \quad (1)$$

The first term of the Hamiltonian (1) contains the energetic contribution ($\Delta = E_{HS} - E_{LS}$) arising from the difference in the ligand–field energy between the HS and LS states at 0 K. The term, $k_B T \ln g$, results from the entropic contribution ($k_B \ln g$) originating from the electronic, intramolecular vibrations and phonons changes between the LS and HS states.
This degeneracy ratio, $g = \frac{g_{HS}}{g_{LS}}$, between the two states is assumed here, for simplicity, as temperature-independent and must be seen as an effective degeneracy. For the dopant, which does not show any spin transition, this energy is zero.

The second and third contributions in (1) are the elastic interactions between nearest neighbors (nn) and next-nearest neighbors (nnn) sites, respectively, whatever their nature (SCO or defect). The parameters, A and B, are the corresponding nn and nnn elastic constants, assumed here as independent from the nature (SCO, defect) of the connected sites. The quantity, $r_{ij} = \|\vec{r}_i - \vec{r}_j\|$ (resp. $r_{ik}$) is the instantaneous distance between the nn (resp. nnn) nodes i and j (resp. i and k). In the expression (1), the parameters $R(S_i, S_j)$ [R'(S_i, S_k)] are the equilibrium bond lengths between two nn nodes i and j (resp. i and k), which depend on the spin state and the nature (SCO or defect) of the connected sites. Owing to the square symmetry of the lattice (see Figure 1), we simply take, $R'(S_i, S_k) = \sqrt{2}R(S_i, S_k)$. The different spin and lattice configurations of nn and nnn are summarized in Table 1. There, $R_{HH}(R'_{HH})$, $R_{LL}(R'_{LL})$ and $R_{HL}(R'_{HL})$ correspond to the nn (nnn) equilibrium bond lengths between HS-HS, LS-LS and HS-LS sites, respectively. Similarly, $R_{0H}(R'_{0H})$, $R_{0L}(R'_{0L})$ and $R_{00}(R'_{00})$ represent the respective nn (nnn) equilibrium bond length distances of HS-defect, LS-defect and defect-defect distances, respectively.

<figure><img src="image_7.png" /><figcaption>**Figure 1.** (a) Schematic representation of the configuration of the nearest (nn) and next-nearest (nnn) neighbor elastic interactions, considered in the present two-dimensional square lattice model. The red central site, labelled (i,j), is surrounded by four nn and four nnn sites (dark blue sites). The connections between the sites represent the springs. Each site can be in LS, HS or a dopant site. (b) Schematic view of the equilibrium nn bond lengths between two SCO sites. $R_{LL}$, $R_{HL}$ and $R_{HH}$ are the equilibrium distances between LS-LS, HS-LS and HS-HS units, respectively.</figcaption></figure>

<table><caption>**Table 1.** Spin states and bond length configurations used in the present model. H (L) denotes a HS (LS) site, while 0 represents the defect.</caption><thead><tr><th>Configuration</th><th>Bond Lengths for nn</th><th>Bond Lengths for nnn</th></tr></thead><tbody><tr><td><b>H-H</b></td><td>R<sub>HH</sub> = R(+1, +1)</td><td>R'(+1, +1) = &radic;2R<sub>HH</sub></td></tr><tr><td><b>L-L</b></td><td>R<sub>LL</sub> = R(-1, -1)</td><td>R'(-1, -1) = &radic;2R<sub>LL</sub></td></tr><tr><td><b>H-L</b></td><td>R<sub>HL</sub> = R(&plusmn;1, &mpmn;1) = (R<sub>HH</sub> + R<sub>LL</sub>)/2</td><td>R'(&plusmn;1, &mpmn;1) = &radic;2R<sub>HL</sub></td></tr><tr><td><b>0-H</b></td><td>R<sub>0H</sub> = R(0, +1) = R(+1, 0) = &frac12;(R<sub>HH</sub> + R<sub>00</sub>)</td><td>R'<sub>0H</sub> = R<sub>0H</sub>&radic;2</td></tr><tr><td><b>0-L</b></td><td>R<sub>0L</sub> = R(0, -1) = R(-1, 0) = &frac12;(R<sub>LL</sub> + R<sub>00</sub>)</td><td>R'<sub>0L</sub> = R<sub>0L</sub>&radic;2</td></tr><tr><td><b>0-0</b></td><td>R<sub>00</sub> = R(0, 0)</td><td>R'<sub>00</sub> = R<sub>00</sub>&radic;2</td></tr></tbody></table>

## 2.1. The General Expression of the Equilibrium Lattice Distance and Model Parameter Values

After a careful analysis of the constraints of Table 1 on the expressions of the bond lengths for the different spin configurations, we have been able to find a general formula Equation (2) linking the equilibrium bond length, $R(S_i, S_j)$, with the spin states $(S_i, S_j)$.
This kind of expression is very useful and allows a better and more efficient Monte Carlo (MC) coding for the numerical simulations. The latter is written as:

$$R(S_i, S_j) = R_{HL} + \frac{\delta_R}{4}(S_i + S_j)S_iS_j + R_{00}(S_i^2 - 1)(S_j^2 - 1) + \frac{1}{2}(S_i^2 - S_j^2)^2 [(S_i + S_j + 1)R_{0H} - (S_i + S_j - 1)R_{0L}], \quad (2)$$

where $\delta_R = R_{HH} - R_{LL}$ is the lattice misfit between the pure (non-diluted) LS and HS phases. The reader can easily check that all configurations of Table 1 are included in this formula.

Throughout this report, we will consider the following relations between the equilibrium distances, $R_{HL} = \frac{R_{HH}+R_{LL}}{2}$, $R_{0H} = \frac{R_{00}+R_{HH}}{2}$ and $R_{0L} = \frac{R_{00}+R_{LL}}{2}$, which decreases the number of free parameters from six to three (see Table 1). By fixing $R_{HH}$ and $R_{LL}$, only one free parameter, $R_{00}$, remains which will be used as a variable to monitor the thermal properties of the diluted SCO materials.

The following parameter values, derived from the experimental literature, are used in the present model. The ligand-field energy is taken constantly throughout the study as $\Delta = 450$ K [73] and is assumed as being independent of the presence of the defect, whose effect is then purely elastic. The latter value leads to an enthalpy change $\Delta H = 3.5$ kJ·mol⁻¹, in good agreement with the experimental data available in the literature [74]. For the degeneracy ratio, we take the usual value, $g = 150$ [75] leading to a molar entropy change at the transition $\Delta S = N_A k_B \ln g = 41.63$ J·K⁻¹·mol⁻¹, $N_A$ being the Avogadro number. Similarly, we assume the degeneracy ratio as independent of the dilution effects as well as on temperature, as previously indicated. Within these values, the transition temperature is evaluated to, $T_{eq} = \frac{\Delta H}{\Delta S} = \frac{\Delta}{k_B \ln g} \approx 90$ K for the pure SCO lattice.

The values of the nn equilibrium distances in the LS and HS states are, respectively taken as equal to $R(-1, -1) = 1$ nm and $R(+1, +1) = 1.2$ nm, which leads for the HL configuration to $R(+1, -1) = \frac{R(+1,+1)+R(-1,-1)}{2} = 1.1$ nm. As already stated, the equilibrium nnn distances are simply considered as $\sqrt{2}R(+1, +1) = 1.7$ nm, $\sqrt{2}R(-1, -1) = 1.4$ nm and $\sqrt{2}R(+1, -1) = 1.55$ nm for HS-HS, LS-LS and HS-LS sites, respectively.

The nn molecules of the lattice interact elastically with the elastic constants $A = 24,000$ Knm⁻², and we take for the nnn sites, $B = 0.3A$, whose role is to maintain the stability of the 2D lattice. The above values of $A$ and $B$ lead to a bulk modulus, $E \approx \frac{(A+2B)}{R} = 6$ GPa, which is in quite fair agreement with data of the SCO literature [76,77] ranging between 5 and 20 GPa.

## 2.2. Simulations Method

The Hamiltonian (1) is solved by MC simulations. First, we generate a 2D square lattice with the bond lengths of the HS phase. The initial lattice coordinates of a site at the position $(i, j)$ are $X(i, j) = R_{HH}(i-1)$ and $Y(i, j) = R_{HH}(j-1)$.

One of the critical points, when we perform a metal dilution relating to the positions of the defects in the lattice, let us consider a lattice with $N$ nodes and fixed concentration, $p$, of defects. Then, there will be $\frac{N!}{(Np)!(N-Np)!}$ possible lattice configurations for the spatial distribution of the $Np$ defects among the $N$ lattice sites. At the level of one single crystal, the defect positions emerging from the chemical synthesis of the material are definitely fixed (we neglect diffusion effects). Therefore, one can draw randomly the defect positions and proceed to the simulations with the selected configuration.

However, it is clear that the results will depend on the spatial distribution of the dopants in the lattice. Thus, the results of the spin transition phenomenon may differ from one diluted microcrystal to another, even taken from the same batch. In the case of magnetic studies, for example, a large number of SCO powder particles or microcrystals is used. The emerging magnetic response is then that of the average responses of several particles, having more or less the same concentration of dopant but with different spatial distribution from one particle to another. In this work, we did not consider the possible situation where
the defects may prefer, for affinity reasons, to be aggregated. The interactions between the defects (except elastic ones) are then ignored.

Along the simulations, several spatial distributions of the defects are generated randomly for each fixed dopant concentration. Once the defect positions are decided, they are fixed throughout the simulations. To ensure the random character of the spatial distribution of defects, using different random seeds we generate several spatial defect distributions at fixed rates of dilution. Then, we study the thermodynamic properties of each lattice configuration. At the end, the results are averaged on the ensemble of the studied lattices.

In practice, the stochastic algorithm of MC simulations is performed as follows. For a site "i" randomly selected with spin ($S_i = 0, \pm 1$) and position $\vec{r}_i$, a new spin value $S'_i(= -S_i)$ is set without a position change. This new spin value is accepted or rejected by the usual Metropolis criterion. It is worth noticing that we do apply the usual Metropolis criterion, even for the defect for which the spin flip does not change its total energy. Whatever the result of the spin flip attempt (the new spin value is accepted or rejected), the lattice positions are relaxed mechanically by a slight motion of each node (randomly selected) with a quantity $||\delta\vec{r}_i|| \ll ||\vec{r}_{ij}||$, where $||\vec{r}_{ij}||$ is the distance between nn. This procedure is repeated 100 times so as to reach the mechanical equilibrium. Here, the reader can realize the importance of Equation (2), which allows the calculation of the equilibrium distance for any spin value (0, +1 or -1). In the present study we take the amplitude of the random displacements as $\delta r = 0.02$ nm. This value remains quite small compared to the LS equilibrium distance $R_{LL} = 1$ nm. Next, we randomly select a new site, and we repeat the procedure for all lattice sites. Once all the nodes of the lattice are inspected for the spin change, we define such a step as the unit of Monte Carlo step "MCS". In these simulations, 1000 MCS are used to reach the thermal equilibrium without storing the data. For each MCS, the lattice is visited $10 \times N$ times to ensure the mechanical equilibrium, which means that each node is moved $10^4 \times N$ times. The same number of MCS is used to average the physical quantities, like the high spin fraction, $n_{HS}$, defined as $n_{HS} = \frac{1+\langle S_i \rangle}{2}$, the average lattice bond length $\langle r \rangle = \frac{\sum_{ij}|r_{ij}|}{N_{nn}}$, where $\langle S_i \rangle = \frac{\sum_i S_i}{N}$ and $N_{nn} \approx 2N$ is the number of nn sites.

### 3. Results and Discussion

We have considered a square lattice ($N = 40 \times 40$) containing active spin crossover atoms and randomly distributed inactive defects with a fixed concentration, $p$. In addition to the concentration of defects, which is a natural relevant control parameter, we use the equilibrium bond length, $R_{00}$, between nn defects, chosen in the interval $[R_{LL}, R_{HH}]$, to tune the thermal properties of this diluted SCO lattice.

#### 3.1. Thermal and Structural Properties for $R_{00} = R_{LL}$

This section is devoted to the analysis of the effect of the concentration of dopant on the thermal and structural properties of an SCO lattice. Concentrations of defects going from 1% to 30% are considered with a defect-defect nn equilibrium distance equal to that of LS. This situation corresponds to the substitution of Fe by an inactive metal with an ionic radius equal to that of Fe in the LS state. As a consequence, it is expected that the defect favors structurally the LS lattice. To avoid elastic frustration problems, we started all simulations from the low-temperature region by preparing the lattice with the bond lengths of the LS state. Therefore, in the LS state, the lattice has zero elastic energy excess since all nn distances are equal to the equilibrium bond length, $R_{LL}$. The temperature is then varied from 1 to 200 K on heating and cooling, with a temperature step of 1 K. As already explained, the simulation is conducted with ten generated random lattices having the same concentration of dopant. Then, the physical quantities of interest are averaged over this sample. The obtained thermal dependencies are presented in Figure 2.
<figure><img src="image_8.png" /><figcaption>**Figure 2.** Thermal dependence of (a) the HS fraction, *n*<sub>HS</sub> and (b) average bond lengths ⟨r⟩, for different concentrations of defects when *R*<sub>00</sub> = *R*<sub>LL</sub>. (c) Evolution of transition temperatures *T*<sub>LS→HS</sub> and *T*<sub>HS→LS</sub> from LS to HS (resp. HS to LS) as a function of percentage of defects. (d) Variation of lattice parameter <r> ><sub>LS</sub> and <r> ><sub>HS</sub> at low and high temperature, respectively as a function of defects concentration. (e) Comparison between MC and analytic equilibrium temperature, *T*<sub>eq</sub>. See text for more explanations.</figcaption></figure>

Figure 2a,b, summarizing the respective temperature evolution of the HS fraction, *n*<sub>HS</sub> and average bond length, ⟨r⟩, show that the pure SCO lattice (0% defects) exhibits an abrupt spin transition with hysteresis width (Δ*T* = 80 K). When the concentration of defects increases, the cooperativity of the SCO reduces, as clearly seen in the reduction in the thermal hysteresis width. At the threshold value of defect concentration, *p*<sub>c</sub> = 0.225, and the first-order spin transition changes into a gradual transition, which corresponds to a simple Boltzmann population between two-degenerate levels. Figure 2c, depicting the phase diagram, *p*-*T*, shows that beyond the critical value, *p*<sub>c</sub>, the two branches of the thermal hysteresis join and overlap, while keeping a linear increase with the concentration of defects. On the other hand, the equilibrium transition temperature, defined as
$T_{eq} = \frac{T_{HS→LS} + T_{LS→HS}}{2}$, shifts towards high-temperature regions as the concentration of defects increases. This fact confirms that the small size dopant induces a positive pressure in the lattice, something which delays the LS to HS spin transition on heating and precipitates the HS to LS conversion on cooling. This behavior is in good agreement with the data of Figure 2c, showing a linear behavior of the upper, $T_{LS→HS}$, and lower, $T_{HS→LS}$, switching temperatures, recalling the behaviors of thermal hysteresis of SCO materials under-external pressure [5]. The remarkable difference in the slopes of $T_{LS→HS}$ and $T_{HS→LS}$, as function of the concentration of defect, is attributed to the kinetic effects of the MC procedure. Indeed, the switching temperature corresponds to the occurrence of a spinodal point, that is the point at which the lifetime of the metastable state falls down, which is crucially temperature-dependent. Consequently, the lower transition temperature $T_{HS→LS}$, will be more sensitive to internal pressure than $T_{LS→HS}$, which itself is located in a region where thermal fluctuations are already dominant.

To better understand these behaviors, we monitored the thermal variation of the average bond length, $\langle r \rangle$, and the dependence for different $R_{00}$ values. The results are given in Figure 2b. First of all, it is interesting to remark on the clear difference in the shape of hysteresis loops of $n_{HS}(T)$ and $\langle r \rangle(T)$. In particular, Figure 2b reveals the presence of a monotonous decrease in the average HS bond length at high-temperature conditions (200 K for example), while this effect is absent in $n_{HS}(T)$ (Figure 2a). Thus, although $\langle r \rangle_{HS}$ becomes lower than $R_{HH}$ when increasing $p$ values, the spin system still reaches the full HS state (see Figure 2a). On the other hand, the absence of an effect in the LS state is easily understandable, since $R_{00}$ is chosen here as being equal to $R_{LL}$.

One of the important consequences of the reduction in the average HS bond length is the decrease in the lattice misfit between the HS and the LS phases. Knowing that the strength of the elastic interactions depends on the quantity, $(A + 2B)(\langle r \rangle_{HS} - \langle r \rangle_{LS})^2$, one can easily understand that by decreasing $\langle r \rangle_{HS}$. This in turn becomes smaller than $R_{HH}$, the steepest character of the thermal hysteresis for relatively small concentrations of defects transforms to bent hysteresis, with lower slopes of the heating and cooling branches.

### 3.1.1. Analytical Study through a Homogenous Elastic Lattice

Our goal here is to understand the effect of the dopant's concentration on the behavior of the average lattice bond length in the low- and high-temperature limits. For that, we perform an analytical study, based on the minimization of the expression of the total energy (elastic + electronic) in the HS and LS states with respect to the lattice parameter for a fixed concentration of defect. The calculation is conducted assuming a homogeneous lattice with a unique bond length, $x$. The energy densities of the HS and LS states, derived from Hamiltonian, (2) write,

$$ \frac{E_{HS}^{tot}}{N} = (\Delta - k_B T \ln g)(1-p) + (A + 4B)[(1-p)^2(x - R_{HH})^2 + p(1-p)(x - R_{0H})^2 + p^2(x - R_{00})^2], \quad (3) $$

$$ \frac{E_{LS}^{tot}}{N} = -(\Delta - k_B T \ln g)(1-p) + (A + 4B)[(1-p)^2(x - R_{LL})^2 + p(1-p)(x - R_{0L})^2 + p^2(x - R_{00})^2]. \quad (4) $$

The mechanical equilibrium is obtained by minimizing the total energy with respect to the $x$ variable in the HS and LS phases, which leads to the following general expressions of relaxed bond lengths:

$$ x_{relax}^{HS} = \frac{p^2 R_{00} + p(1-p)R_{0H} + (1-p)^2 R_{HH}}{p + (1-p)^2}, \quad (5) $$
$$x_{relax}^{LS} = \frac{p^2 R_{00} + p(1-p)R_{0L} + (1-p)^2 R_{LL}}{p + (1-p)^2}. \qquad (6)$$

Using now $R_{00} = R_{LL}$ and $R_{0H} = R_{HL}$, we obtain:

$$x_{relax}^{LS} = R_{LL} \text{ and } x_{relax}^{HS} = R_{HH} - \frac{\delta_R}{2} \frac{p(p+1)}{p + (1-p)^2} \qquad (7)$$

Figure 2d summarizes the comparison between the analytical and MC data for the dopant concentration-dependence of the bond lengths, showing a good agreement between the two set of results.

To determine the analytical expression of the transition temperature, we first define the value of the relaxed bond length at the transition temperature, $T_{eq}$, as $x_{1/2} = \frac{x_{relax}^{HS} + x_{relax}^{LS}}{2}$. Then, we calculate $T_{eq}$ using the relation $E_{tot}^{HS}(x_{1/2}) = E_{tot}^{LS}(x_{1/2})$, which remains valid as long as we neglect the configurational entropy of spins contribution. This leads to:

$$x_{1/2} = R_{HL} - \frac{p(1+p)}{[p + (1-p)^2]} \delta R, \qquad (8)$$

$$T_{eq}(p) = T_{eq}(0) + (A + 4B) \frac{2p\delta R}{k_B l n g} \left[ \frac{(1+p)(2-p)}{[p + (1-p)^2]} \delta R + R_{LL} \right]. \qquad (9)$$

Figure 2d,e present the p-dependence of $x_{1/2}$ (that is $\langle r \rangle$) and $T_{eq}$ obtained with Equations (8) and (9) and MC simulations. It appears that, while a good agreement is obtained for the evolution of $\langle r \rangle$ (Figure 2d), that of the temperature (Figure 2e) is subtler. Indeed, a fair agreement exists between the analytical and MC data for p values ranging between 0 and 30%, meeting the assumption of the homogenous concentration of defects. However, a discrepancy emerges for strong concentrations of defects. This is attributed to the fact that the assumption of homogenous distribution of dopants, used in the analytical approach, is less valid for $p > 0.3$. Overall, the increase in the transition temperature with dopant concentration is clearly due to $R_{00} = R_{LL}$, which clearly plays the role of an internal positive pressure in the lattice.

### 3.1.2. Spatiotemporal Aspects of the Phase Transition

The spatiotemporal studies of the thermoinduced spin transition concern the analysis of the mechanisms of formation (nucleation and growth) of HS and LS domains, along with the cooling and heating branches of the thermal hysteresis. For that, the spatial configurations of the spin-lattice states are stored at selected temperatures. Figure 3 displays a set of snapshots of the spin and lattice configurations along the thermal hysteresis loop of Figure 2a. Here, red, blue and green dots represent the HS, LS and dopant sites, respectively.

Figure 3 reveals, that upon heating ($n_{HS} = 0.25$) for a crystal with 8% of defects, the HS phase grows exclusively and quite easily from the free corners (due to surface effects), and then extends over the whole crystal, by avoiding the defects. It is also remarked that the spatiotemporal aspects of the HS → LS transition (on cooling) are quite different from those observed on heating. Indeed, on cooling the LS phase mainly nucleates around the defects and deploys over the entire crystal. In contrast, on heating the HS species grows far from the defects; these maintain their surrounding sites in the LS state, switching later at the end of the process.

One can also notice a significant lattice distortion near to LS (resp. HS) domains during the heating (resp. cooling) process. On the other hand, as the number of defects increases ($p$ value), the proportion of active SCO centers reduces, and then the lattice distortion is expected to be reduced fairly, due to the reduction in the global volume change.
<figure><img src="image_9.png" /><figcaption>**Figure 3.** Selected crystal snapshots showing both electronic and macroscopic distortion configurations during the heating and cooling processes starting from LS state for a concentration, *p* = 8% of defects and *R*<sub>00</sub> = *R*<sub>LL</sub>. The red (blue) dots represent the HS (LS) sites and green dots stand for the dopant sites.</figcaption></figure>

3.2. *Thermal and Mechanical Properties for R<sub>00</sub> = R<sub>HL</sub>*

In this part, we study the case, R₀₀ = RHL = RHH + RLL / 2, where the dopant is expected
to have a “neutral” elastic role. The concentration of defects varied in the range of 1–30%.
Since R₀₀ is neither equal to RHH, nor to RLL, the initial elastic lattice (HS or LS) will be
frustrated. Consequently, it is important to make sure, during the simulation, that the
lattice is well relaxed at each temperature step. So, starting from the relaxed structure of
the crystal, we vary the temperature in the range 1–200 K on heating and cooling, with a
thermal step of 1 K. The thermoinduced spin transition for the different concentrations
of defect is illustrated in Figure 4 in terms of HS fraction, nHS, and average lattice bond
lengths, ⟨r⟩.

Similarly to panels (a) and (b) of Figure 2, Figure 4a,b show that the presence of
defects reduces the cooperativity inside the SCO lattice, which changes from an abrupt
thermal transition to gradual transition to a high concentration (≥ 24%). Here, the equi-
librium temperature, T<sub>eq</sub>, is slightly shifted towards higher temperature. To understand
this decrease in cooperativity, we monitored the thermal variation of the average lattice
parameters in Figure 4b. One can notice a difference in the shape of hysteresis loops of n<sub>HS</sub>
and ⟨r⟩ at low and high temperatures. Indeed, increasing the proportion of defects results
in an increase in ⟨r⟩<sub>LS</sub> and as decrease in ⟨r⟩<sub>HS</sub>, as shown in Figure 4c. As a result,
the global change in lattice parameter (δ⟨r⟩ = ⟨r⟩<sub>HS</sub> − ⟨r⟩<sub>LS</sub>) between the HS and LS
phases decreases linearly with p.

As already explained, the reduction in the lattice misfit between the LS and HS states
causes a reduction in the strength of the elastic interactions. This effect is responsible for
the long-range nature of the cooperative interactions in SCO solids. Thus, the decrease in
δ⟨r⟩ as function of *p* results in the vanishing of the thermal hysteresis width on increasing
the dilution rate. In other words, the defects prevent the propagation of the elastic strain
and isolate the active SCO species, leading to a cut in the interaction paths.

The analytical dependence of <r>_{HS} and <r>_{LS} can be obtained by minimizing
the total energy in the HS and LS states. Using expressions (5) and (6), we find for the
present case:
$$x_{relax}^{LS} = R_{HL} - \frac{(1-p)(2-p)}{p + (1-p)^2} \frac{\delta R}{4} \quad \text{and} \quad x_{relax}^{HS} = R_{HL} + \frac{(1-p)(2-p)}{p + (1-p)^2} \frac{\delta R}{4} \quad (10)$$

<figure><img src="image_10.png" /><figcaption>**Figure 4.** (a,b) Thermal dependence of the HS fraction, $n_{HS}$ and average bond lengths $\langle r \rangle$ for different concentration of defects when $R_{00} = R_{HL}$. (c) Variation of lattice parameter $\langle r \rangle_{LS}$ and $\langle r \rangle_{HS}$ at low and high temperature, respectively, as a function of defect concentration. Filled circles and squares are MC data, and full pink and blue lines are analytical predictions deduced from Equation (10). (d) Evolution of transition temperatures, $T_{LS \to HS}$ (resp. $T_{HS \to LS}$), from LS to HS (resp. HS to LS) as a function of the concentration of defects showing a weak p-dependence.</figcaption></figure>

For the relatively small p values, Equation (10) gives for relaxed LS and HS bondlengths, the respective expressions $x_{relax}^{LS} = R_{LL} + p\delta R/4$ and $x_{relax}^{HS} = R_{HH} - p\delta R/4$, which well fit with the MC data of Figure 4c.

On the other hand, from the equality, $E_{tot}^{HS}(x_{1/2}) = E_{tot}^{LS}(x_{1/2})$, we easily obtain the expressions of the transition temperature, $T_{eq}$, and the corresponding value of the average lattice bond length, $x_{1/2}$:

$$T_{eq} = \frac{\Delta}{k_B l n_g} \quad \text{and} \quad x_{1/2} = R_{HL}. \qquad (11)$$

Equation (11) indicates that we expect a constant equilibrium temperature, independent of the concentration of dopants. This is in quite good agreement with the MC results of HS fraction and $\langle r \rangle$ (Figure 4a,b), which show that the hysteresis loops are almost all centered around $T = 90$ K for the different p values. In addition, Figure 4d, which summarizes the behavior the p-dependence of the limiting switching temperatures of the hysteresis, has an almost symmetrical shape. It is thus deduced that the MC equilibrium temperature, $T_{eq}$, slightly changes with p, in agreement with Equation (11). This trend is attributed to the neutral elastic character of the defect, whose size $R_{00} = R_{HL}$, is in the middle between those of HS and LS bond lengths.
## Spatiotemporal Aspects

Here, we investigate the spatiotemporal properties of the thermal hysteresis of Figure 4. The results are presented in Figure 5. We notice that upon heating (for $n_{HS}$ = 0.25) for a crystal with 6% of defects, the HS phase grows from the right corners and edges where the elastic energy is weaker. Then, the domains extend towards the center of the lattice with a clear interface, avoiding the regions with rich concentrations of defects, which convert only at the end of the process. Contrary to the previous case ($R_{00}$ = $R_{LL}$), here the spatiotemporal feature of the HS → LS transition (on cooling branch) is quite similar to that of the heating (LS → HS) branch. Thus, this behavior confirms the neutral role of the defect, which does not favor the nucleation of one state over the other.

<figure><img src="image_11.png" /><figcaption>**Figure 5.** Snapshots of the crystal during the heating and the cooling processes starting from LS state for 6% of defects corresponding to Figure 4a. The red (blue) dots are associated with HS (LS) sites, while the green dots represent the dopant. The nn defect-defect distance is $R_{00} = R_{HL}$.</figcaption></figure>

### 3.3. *Thermal and Mechanical Properties for R<sub>00</sub> = R<sub>HH</sub> and Emergence of HS Residual Fraction*

Now, we analyze the special case where the equilibrium nn dopant-dopant distance is equal to nn bond length of the HS phase ($R_{00} = R_{HH}$), for various defect concentrations, *p*, going from 1 to 30%. Contrary to the case $R_{00} = R_{LL}$, here we start the simulations from the high temperature (i.e., HS) phase where the system is elastically relaxed. The system is then cooled from 200 to 1 K and then heated with a thermal step of 1 K.

The thermoinduced spin transitions for the different concentrations of defects are illustrated in Figure 6 in terms of $n_{HS}$ and average lattice bond lengths, $\langle r \rangle$, respectively.

One can notice from Figure 6a that, once again, the presence of defects reduces the cooperativity of SCO units leading to a gradual transition for defect concentrations bigger than 24%. This reduction in the thermal hysteresis width is now accompanied with the appearance of a residual HS fraction at a low temperature. The presence of this residual HS fraction, leading to incomplete spin transition at low temperatures and its increase with dopant concentration, is in direct line with the thermal variation of the lattice parameters. Indeed, Figure 6b presenting the thermal dependence of $\langle r \rangle$ reveals a significant increase in the average bond length distance at low temperature as *p* increases. The simple explanation is that, within the presently used equilibrium bond length value ($R_{00} = R_{HH}$) of the dopant lattice, the SCO lattice expands locally around the defect positions, particularly in the LS state. This creates a tensile stress which acts on the effective ligand field and favors the local nucleation of HS species around defect regions even at low temperature, when the
local bond distances exceed the $R_{HL}$ value. A careful examination of the p-dependence of the upper and lower switching temperatures, presented in Figure 6c, demonstrates that almost only the HS→LS transition is affected by the dilution.

<figure><img src="image_12.png" /><figcaption>**Figure 6.** (a,b) Thermal dependence of the HS fraction, $n_{HS}$ and average bond lengths $\langle r \rangle$ for different concentration of defects when $R_{00} = R_{HH}$. (c) Evolution of transition temperatures $T_{LS \to HS}$ (resp. $T_{HS \to LS}$) from LS to HS (resp. HS to LS) as a function of defects concentration. (d) Variation of lattice parameter $\langle r \rangle_{LS}$ (resp. $\langle r \rangle_{HS}$) at a low (resp. high) temperature, versus the concentration of defects. (e) Comparison between numerical and analytic equilibrium temperature, $T_{eq}$, as a function of $p$, showing a fair qualitative tendency.</figcaption></figure>

To compare the obtained numerical values of $\langle r \rangle$ at low and high temperatures with the analytical ones, we minimize the total energy in the HS and LS states, using Equations (5) and (6), for $R_{00} = R_{HH}$. The calculations lead to the following expressions of relaxed bond lengths,
$$x_{relax}^{LS} = R_{LL} + \frac{p\delta_R}{2} \cdot \frac{1+p}{p+(1-p)^2} \quad \text{and} \quad x_{relax}^{HS} = R_{HH}. \qquad (12)$$

The expressions of the transition temperature, $T_{eq}$, and corresponding average bond length, $x_{1/2}$, obtained from the relation $E_{tot}^{HS}(x_{1/2}) = E_{tot}^{LS}(x_{1/2})$ (Equations (3) and (4)), are given below:

$$x_{1/2} = R_{HL} + \frac{p(1+p)}{[p+(1-p)^2]} \delta R, \qquad (13)$$

$$T_{eq}(p) = T_{eq}^{(0)} - \frac{(A+4B)}{8k_B \text{Ing}} p \delta R \left[ \frac{(1+p)(2-p)}{[p+(1-p)^2]} \right] \delta R + R_{HH}. \qquad (14)$$

The comparison between the MC results and the analytical predictions show a good qualitative agreement for the general tendencies, as depicted in Figure 6d,e, although some discrepancy exists for the transition temperature (Figure 6e). The latter is mainly due to the limitations of the validity of the homogeneous distribution of dopants assumed in the analytical approach.

## Spatiotemporal Aspects

We have analyzed the presence of defect having the same size as that of HS ion on the spatiotemporal aspects of the thermoinduced spin transition. The obtained results are presented in Figure 7, where we display selected snapshots of the spin state configuration along the thermal hysteresis loop of Figure 6a.

<figure><img src="image_13.png" /><figcaption>**Figure 7.** Snapshots of the crystal showing the nucleation and growth features of a diluted SCO lattice along the heating and cooling processes for $p = 6\%$ of defects. The red (blue) dots are associated with HS (LS). The green sites represent the dopant sites. The nn defect-defect distance is, $R_{00} = R_{HH}$.</figcaption></figure>

One can notice from Figure 7 that, on heating from the LS state, the nucleation of the HS phase mainly starts around the defects thanks to the negative pressure they produce on the LS sites, since, $R_{0L} = \frac{R_{00}+R_{LL}}{2} = \frac{R_{HH}+R_{LL}}{2} = R_{HL} > R_{LL}$. This means that each LS site with a defect neighbor experiences a tensile stress, which helps with the emergence of HS species around the defects. Thus, the nucleation of the HS phase appears more or less as homogeneously appearing everywhere—in the center, the corners and the edges—due
to a random spatial distribution of the defects. In contrast, on cooling, the LS phase tends to appear mostly from corners and edges, far from defects, and coalesce at the end of the process around the defect regions, which are the last to convert into the LS state. In all snapshots, we see that the lattice deforms along the transition, but these deformations remain less severe than those of the pure SCO solid.

### 3.4. Low-Temperature Relaxation of the Photo-induced Metastable HS State

The metastable HS state is produced at low temperature by a photo-excitation through an adequate wavelength by the so-called LIESST effect [67]. In addition, a rapid quenching (or flash cooling) from the high temperature thermodynamic stable HS phase may also lead to the generation of a quenched metastable HS phase at a low temperature. However, the microscopic nature of the latter is different from that of the photo-induced HS state [78–80]. Here, we aim to analyze the effect of dilution on the isothermal HS→LS relaxation process of a SCO system. For that, we investigate the temporal evolution of the HS fraction, $n_{HS}$, the average lattice parameter $\langle r \rangle$ as well as the spatiotemporal properties of the spin organization and lattice deformation along this isothermal relaxation process. Initially, the lattice of size, ($N = 40 \times 40$), is prepared in the HS state from the electronic ($S_i = +1$) and structural ($r_{ij} = R_{HH}$) point of view at $T = 1$ K. At this temperature, the prepared HS state is clearly metastable, and the LS state ($S_i = -1$, $r_{ij} = R_{LL}$) is the thermodynamic stable phase. Then, isothermal MC simulations are performed on both spin and position variables for various dopant concentrations, $p$, and defect sizes ($R_{00}$). The objective here is to enhance the understanding of the combined effects of $p$ and $R_{00}$ on the features of the electro-elastic relaxation mechanism.

#### 3.4.1. Case of $R_{00} = R_{LL}$

In this section, we report the time dependence of the average HS fraction, $n_{HS}$, and lattice parameter, $\langle r \rangle$, during the relaxation from the metastable HS fraction at $T = 1$ K for several defect concentration values ranging from $p = 2\%$ to 10%. The obtained results for the relaxation of the HS fraction (resp. $\langle r \rangle$) are summarized in Figure 8a (resp. Inset of Figure 8a). First, we see that the shape of the relaxation curves changes from being sigmoidal (non-linear) for low $p$ values to almost being exponential for $p = 10\%$. This behavior is clear evidence of the decreasing cooperative character of the system with dilution. On the other hand, the lifetime, $t_{1/2}$, of the metastable HS state, defined as the transit time at $n_{HS} = 1/2$, is expected to grow exponentially as $t_{1/2} \sim t_0 + \alpha e^{\Delta E/kT}$ with the energy barrier $\Delta E$. It is important to notice that $\Delta E$ contains two contributions. The first one is the electronic gap, $\Delta E_{elec} = E_{elec}^{LS} - E_{elec}^{HS} = -2(\Delta - k_B T \ln g)(1 - p)$, which is negative, owing to the stable nature of the LS state at 1 K. The second one is elastic, $\Delta E_{elast} = E_{elast}^{LS} - E_{elast}^{HS}$, which is positive, and constitutes the energy barrier facing the switching between HS and LS states. These two competing contributions are well identified in Equations (3) and (4). Thus, $\Delta E_{elec}$ is the driving force of the relaxation of the HS sites towards the LS state, while $\Delta E_{elast}$ opposes to this spin switching. Overall, the profile of the relaxation curve is due to the competition between these two contributions. In a homogeneous system, with low dopant concentrations (1 to 10%), both energy gaps are expected to linearly depend on the dopant concentration, $p$, as shown in Equations (3) and (4). Interestingly, Figure 8b shows that the $\ln(t_{1/2} - t_0)$ linearly decreases with $p$, which corroborates the assumption of the linear $p$-dependence of $\Delta E_{elast}$. The latter is demonstrated in Figure 8c, which reports the linear decrease in the elastic energy with the dopant concentration, while this is obvious for the electronic energy.
<figure><img src="image_14.png" /><figcaption>**Figure 8.** (a) Time dependence of the metastable HS fraction, *n*<sub>HS</sub>, and average bond lengths ⟨r⟩ (given in inset) throughout the HS→LS relaxation for a square lattice of size 40 × 40. The data have been obtained by MC simulations at T = 1 K for different concentration of defects: *p* = 2% (black), 4% (red), 6% (blue), 8% (pink) and 10% (green). (b) Linear behavior of the logarithm of the lifetime, ln(*t*<sub>1/2</sub> − *t*<sub>0</sub>), of the HS state as a function of *p*. (c) Density of elastic energy as a function of the concentration of defects. (d) Spatiotemporal configurations of the spin and lattice transformations during the HS→LS relaxation. Red (blue) area represents HS (LS) phases, and the green dots stand for the defect. From top to bottom: *p* = 2%, 6% and 10%. Left to right: *n*<sub>HS</sub> = 0.75, 0.5 and 0.25.</figcaption></figure>
Figure 8d illustrates the electro-elastic configuration of the lattice during the relaxation process for three *p* values: 2%, 6% and 10%. In the case of *p* = 10%, the LS phase appears dispersed over the lattice, appearing everywhere, even in the center, which explains the exponential character of the corresponding relaxation curve. This rapid relaxation towards the LS state is induced by the presence of a significant concentration of dopants with size *R*<sub>00</sub> = *R*<sub>LL</sub>, leading to shrink the lattice. This elastic effect acts as an external non-uniform pressure affecting all HS sites. In contrast, for *p* = 6%, we observe that the nucleation of LS domains mainly starts from corners and edges, and then propagates toward the center. Upon decreasing *p*, the number of LS nucleation centers decreases until a single macroscopic domain nucleation mechanism emerges for *p* = 2%. The latter recalls the results obtained for non-diluted SCO systems, leading to sigmoidal relaxation curves of the HS fraction, characteristic of strong cooperative systems.

3.4.2. Case of $R_{00} = R_{HL}$

The relaxation is now performed with defects, assumed to be elastically neutral, since their size, $R_{00} = R_{HL} = \frac{R_{HH}+R_{LL}}{2}$, is between those of HS and LS. The simulations are conducted for different concentrations of defects and we proceed by starting from the metastable HS state, characterized by $\langle r \rangle = R_{HH}$ and $S_i = 1 \forall i$, i.e., $n_{HS} = 1$. The temporal dependence of the HS fraction, $n_{HS}$, and average bond length $\langle r \rangle$, are summarized in Figure 9a,b, respectively. As for the results of Figure 8, the lifetime of the metastable HS state, decreases with the concentration of dopants and the HS fraction reaches the value 0 (that of LS state) at the end of the relaxation. In addition, the shape of the relaxation curves changes from sigmoidal to exponential as *p* increases. On the other hand, the curves of the average bond length, $\langle r \rangle$, displayed in Figure 9b, show a slight difference, compared to those of $n_{HS}$. Indeed, at the end of the relaxation, the stationary value of $\langle r \rangle_{LL}$ remains bigger than the average bond length, $R_{LL}$, of the LS state, although the corresponding value of $n_{HS}$ is zero.

The *p*-dependence of $\langle r \rangle_{LL}$ in the LS state, presented in the inset of Figure 9b shows a quasi-linear dependence, which means that the average lattice bond length in the relaxed LS phase increases with the dilution rate. Consequently, although the spin system converts to the full LS state, a tensile strain does remain in the lattice in the LS phase due to the defect size. To understand this point, let us briefly discuss the analytical expression of $\langle r \rangle_{LL}$. If we consider the SCO lattice as an effective homogeneous medium, then the average lattice parameter in the LS phase is around $\langle r \rangle_{LL} \approx (1-p)R_{LL} + p^2 R_{00} + p(1-p)\frac{R_{00}+R_{LL}}{2}$. For, $R_{00} = R_{HL}$, and neglecting all contributions on $p^2$ (too small), one finds $\langle r \rangle_{LL} \approx R_{LL} + p\frac{\delta R}{4}$, which confirms the linear *p*-dependence of $\langle r \rangle_{LL}$ found in the MC simulations.

<figure><img src="image_15.png" /><figcaption>Figure 9. Cont.</figcaption></figure>
<figure><img src="image_16.png" /><figcaption>**Figure 9.** (a) Time dependence of the metastable HS fraction, $n_{HS}$, along the HS→LS relaxation for a lattice 40 × 40. The data obtained from MC simulations at T = 1 K. The five curves correspond to the following concentrations of defects: p = 2% (black), 4% (red), 6% (blue), 8% (pink) and 10% (green). (b) Corresponding temporal evolution of the average bond length, $\langle r \rangle$. Inset: relaxed average bond length in the LS phase as function of p, showing a quasi-linear behavior. (c) Spatiotemporal configurations of the spin and lattice transformations during the HS→LS relaxation. Red (blue) area represents HS (LS) phases and green dots stand for the defect. From top to bottom: p = 2%, 6% and 10%, respectively. Left to right: $n_{HS}$ = 0.75, 0.5 and 0.25, respectively.</figcaption></figure>

Figure 9c collects some selected snapshots of the electro-elastic configuration of the lattice during the relaxation process for different p values. As for the case of Figure 8d, the relaxation is accompanied with a significant lattice distortion and single-domain nucleation process appears for low dilution rates (2%), while homogenous nucleation is dominant for the higher dilution rate (10%).

### 3.4.3. Case of $R_{00} = R_{HH}$

Now we examine the situation where the size of the dopant is equal to that of the HS ion. In this case, the photo-induced metastable HS state (all SCO spins are set to +1) is free from an elastic strain, since all nn distances are equal to $R_{HH}$. Consequently, the initial diluted HS state has an elastic energy equal to zero. It is then expected that the dopant will
have two roles: first it stabilizes the HS phase, and at some threshold concentration it will destroy the elastic interactions between the SCO atoms, as for the previous cases. The time dependence of the average HS fraction, $n_{HS}$, during the relaxation from the metastable HS is given in Figure 10a and that of $\langle r \rangle$ is presented in the inset of Figure 10a. First, we remark that the shapes of the relaxation curves are very different from those of the two previous studied cases ($R_{00} = R_{LL}$ and $R_{00} = R_{HL}$). Indeed, except the case $p = 2\%$, for all other concentrations, the HS fraction $n_{HS}$ does not reach the value 0 at a long time, which means that a residual HS fraction does subsist in the lattice in the stationary state. Interestingly, this residual fraction, denoted as $n_{HS}^*$, linearly increases with the concentration of dopant, as depicted in Figure 10b. This behavior is also observed on the structural order parameter, $\langle r \rangle$ (see inset of Figure 10a), which does not reach the equilibrium bond length of the LS state, $R_{LL}$, even for weak $p$ values. Consequently, the lifetime of the metastable HS state increases with $p$. This is well explained through the elastic role of the dopant, which expands the lattice, thus preventing its efficient shrinking along the HS → LS relaxation. This negative pressure maintains a global tensile stress that allows the HS state to survive in the immediate vicinity of the dopants. As a result, the SCO atoms close to the dopant stay in the HS state later, compared to those surrounded by SCO sites. Thus, the ligand fields of the SCO sites are renormalized by their nn or nnn dopants, and this renormalization fluctuates from one site to another due to the fluctuations of the local concentration of dopant. At the end, the shape of the relaxation curve of the HS fraction switches from sigmoidal for low $p$ values to stretched exponential, as $p$ increases. The latter is the signature of the decrease in the strength of the interactions combined with disorder arising from the random spatial distribution of the dopants inside the SCO lattice.

Figure 10c illustrates the electro-elastic configuration of the lattice during the relaxation process for different dopant concentration values. Contrary to the previous cases, now the LS phase emerges as a single macroscopic domain mechanism for all concentrations, even for those leading to stretched exponential relaxation shapes (see Figure 10a). In this case, the switching to LS of the HS sites close to defects are delayed to longer relaxation times, due to the negative pressure induced by the dopant, which decreases the effective ligand field of the SCO sites.

<figure><img src="image_17.png" /><figcaption>Figure 10. Cont.</figcaption></figure>
<figure><img src="image_18.png" /><figcaption>**Figure 10.** (a) Time dependence of the metastable HS fraction, *n*<sub>HS</sub>, along the HS→LS relaxation. The data are obtained by MC simulations at T = 1 K, on a square lattice of size 40 × 40, for the following concentrations of defects: *p* = 2% (black), 4% (red), 6% (blue), 8% (pink) and 10% (green). Inset: corresponding time dependence of the average bond lengths ⟨r⟩. (b) Residual HS fraction in the relaxed regime, *n*<sub>HS</sub>*, derived from panel (a), as a function of concentration of defects. (c) Spatiotemporal configurations of the spin and lattice during the HS→LS relaxation. Red (blue) squares represent HS (LS) phases while green square stand for the dopant. From top to bottom: *p* = 2%, 6% and 10%, respectively. Left to right: *n*<sub>HS</sub> = 0.75, 0.5 and 0.25, respectively.</figcaption></figure>

**4. Conclusions**

We investigated by means of MC simulations and an analytical approach the thermal
and the non-equilibrium properties of diluted SCO materials, using an electro-elastic
model combining spins and deformations variables. The present study focuses on the
effects of metal dilution on the bistability of SCO materials. Among them, we examine
the impact of dilution on the thermal hysteresis as well as on the macroscopic nucleation
and growth mechanisms taking place along the thermoinduced phase transition. For
that, we analyze the dependence of the transition temperature and the thermal hysteresis
width along with the concentration of dopant. To prove the important elastic role of
the dopant, we consider three possible sizes for the dopant-dopant distance (R<sub>00</sub>) in the
lattice: R<sub>00</sub> (= R<sub>LL</sub>, R<sub>HL</sub>, R<sub>HH</sub>), corresponding, respectively to those of LS-LS, HS-LS or
HS-HS sites. We find that, when R<sub>00</sub> = R<sub>LL</sub>, the increase in the concentration of defects
changes the thermoinduced spin transition of the crystal from abrupt to gradual, with a
shift of the transition temperature towards lower temperatures [55,61]. It is then found
that the dilution affects both the cooperativity of the system since the hysteresis vanishes
and the effective ligand field, due to the lowering of the transition temperature. When, $R_{00} = R_{HL}$, the thermal hysteresis vanishes on increasing the concentration of defects above the threshold value, *p* = 24%, while the transition temperature remains invariant. In the last case, corresponding to $R_{00} = R_{HH}$, we evidenced the existence of an incomplete SCO transition, with the presence of residual HS fraction at low temperature which linearly depends on the concentration of defects. Here, contrary to the previous cases, the transition temperature gradually decreases with the concentration of dopants, with a vanishing of the thermal hysteresis width due to the collapse of the elastic interactions.

The spatiotemporal aspects of these phase transitions are also investigated as functions of $R_{00}$ and *p*, focusing on the mode of nucleation and growth mechanism along the thermal hysteresis branches. For all studied situations of $R_{00}$, we find that the nucleation starts from the lattice corners for weak *p* values. However, on increasing, *p*, the nature of the nucleation changes significantly according to $R_{00}$ as well as on the considered branch of the thermal hysteresis. Generally, the dopant influences the lattice bond length of its SCO neighbors, which in turn acts on their spin state by creating tensile or compressive stress.

The same reasoning is developed to investigate the non-equilibrium properties of diluted SCO materials by examining the relaxation of the photo-induced metastable HS states at very low temperature. The results of the isothermal relaxation show that the temporal evolution of the HS fraction is reminiscent of that of the thermal dependence of the HS fraction. For low dopant concentrations, all relaxation curves display a sigmoidal behavior, characteristic of cooperative SCO systems. On increasing the rate of dopant, the shape of the relaxation curves changes to regular exponential or stretched exponential according to $R_{00}$ values.

In conclusion, this work demonstrates the ability of the diluted version of the electro-elastic model to mimic the thermal and the temporal behaviors of diluted SCO materials. Several possible extensions and improvements can be realized. For example, the elastic constants between the SCO and the dopant can be taken as dependent on the nature of the connected sites. In addition, dopants with sizes bigger than those of HS spin-crossover atoms are expected to lead (at low concentration) to multi-step or incomplete thermal transitions. Even more importantly, the model can be improved to include structural phase transitions coupled to the spin transition. Within this context, it is important to note that several SCO materials showing structural phase transitions were studied in the literature. One of the well documented case is that of $[Fe(ptz)_6(BF_4)_2]$ [19,81] (abbreviated as Fe(ptz) in the following ptz = 1-propyl-tetrazole). In this compound, the spin transition is coupled to a structural transition between an ordered high temperature structure $R\bar{3}$ and a disordered low-temperature phase. As for the pure material, Fe(ptz), the investigations of the related diluted systems show the presence of a crystallographic phase transition upon slow cooling from the HS state. In a recent study [82] of diluted $Fe_{1-x}Ru_x(ptz)$ single crystals, it has been found that the dilution stabilizes the LS state and induces isothermal multi-stepped relaxation, resulting from the interplay between the spin transition and the structural phase transition. These unexpected and interesting behaviors deserve to be investigated with an adapted elastic model, including lattice symmetry breaking along the spin transition [83] in order to disentangle these two contributions. Furthermore, dinuclear [84] and trinuclear SCO materials [85,86] are well known for exhibiting two- or multi-step transitions due to steric effects, causing elastic frustrations [87] between the SCO units at the transition. The control of the width of the plateau region through metal dilution is of high importance, and the extension of the present model for polynuclear SCO materials, accounting for antagonist elastic effects [33], is also a very promising objective. Finally, the recent progresses in the treatment of 3D versions [88,89] of elastic models opens new opportunities for interesting investigations of the dilution effects in more realistic lattices.
**Author Contributions:** K.A.: software, validation, formal analysis, investigation, data curation, writing the original draft, Y.S.: software, validation, formal analysis, investigation, data curation, writing the original draft. K.B.: designed the model and proposed the subject: formal analysis, investigation, writing, review and editing supervision, project administration and funding acquisition. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by The Agence Nationale de la Recherche, ANR project Mol-CoSM N° ANR-20-CE07-0028-02.

**Data Availability Statement:** Data will be made available on request.

**Acknowledgments:** This work was supported by the ANR project Mol-CoSM N° ANR-20-CE07-0028-02, the Universities of Versailles and Paris-Saclay-UPSAY, and the CNRS (Centre National de la Recherche Scientifique). We thank all of them for their financial support.

**Conflicts of Interest:** The authors declare no conflict of interest.

## References

1. Gütlich, P.; Goodwin, H.A. Spin Crossover—An Overall Perspective. In *Spin Crossover in Transition Metal Compounds I*; Topics in Current Chemistry; Springer: Berlin/Heidelberg, Germany, 2004; Volume 233, pp. 1–47. [CrossRef]

2. Dugay, J.; Giménez-Marqués, M.; Kozlova, T.; Zandbergen, H.W.; Coronado, E.; Van der Zant, H.S.J. Spin Switching in Electronic Devices Based on 2D Assemblies of Spin-Crossover Nanoparticles. *Adv. Mater.* **2015**, *27*, 1288–1293. [CrossRef] [PubMed]

3. Holovchenko, A.; Dugay, J.; Giménez-Marqués, M.; Torres-Cavanillas, R.; Coronado, E.; van der Zant, H.S.J. Near Room-Temperature Memory Devices Based on Hybrid Spin-Crossover@SiO2Nanoparticles Coupled to Single-Layer Graphene Nanoelectrodes. *Adv. Mater.* **2016**, *28*, 7228–7233. [CrossRef] [PubMed]

4. Linares, J.; Codjovi, E.; Garcia, Y. Pressure and Temperature Spin Crossover Sensors with Optical Detection. *Sensors* **2012**, *12*, 4479–4492. [CrossRef] [PubMed]

5. Benaicha, B.; Van Do, K.; Yangui, A.; Pittala, N.; Lusson, A.; Sy, M.; Bouchez, G.; Fourati, H.; Gómez-García, C.J.; Triki, S.; et al. Interplay between spin-crossover and luminescence in a multifunctional single crystal iron(II) complex: Towards a new generation of molecular sensors. *Chem. Sci.* **2019**, *10*, 6791–6798. [CrossRef]

6. Matsuda, M.; Isozaki, H.; Tajima, H. Producible on-off Switching of the Light Emission from the Electroluminescent Device Containing a Spin Crossover Complex. *Thin Solid Film.* **2008**, *517*, 1465–1467. [CrossRef]

7. Abdul-Kader, K.; Lopes, M.; Murgui, C.B.; Kraieva, O.; Hernandez, E.M.; Salmon, L.; Nicolazzi, W.; Carcenac, F.; Thibault, C.; Molnar, G.; et al. Synergistic Switching of Plasmonic Resonances and Molecular Spin States. *Nanoscale* **2013**, *5*, 5288–5293. [CrossRef]

8. Pillet, S. Spin-crossover materials: Getting the most from X-ray crystallography. *J. Appl. Phys.* **2021**, *129*, 181101–181131. [CrossRef]

9. Gallois, B.; Hauw, C.; Real, J.A.; Zarembowitch, J. Structural Changes Associated with the Spin Transition in *Fe(Phen)₂(NCS)₂*: A Single-Crystal X-Ray Investigation. *Inorg. Chem.* **1990**, *29*, 1152–1158. [CrossRef]

10. Tanaka, D.; Aketa, N.; Tanaka, H.; Tamaki, T.; Inose, T.; Akai, T.; Toyama, H.; Sakata, O.; Tajiri, H.; Ogawa, T. Thin Films of Spin-Crossover Coordination Polymers with Large Thermal Hysteresis Loops Prepared by Nanoparticle Spin Coating. *Chem. Commun.* **2014**, *50*, 10074–10077. [CrossRef]

11. Desaix, A.; Roubeau, O.; Jeftic, J.; Haasnoot, J.G.; Boukheddaden, K.; Codjovi, E.; Linares, J.; Nogues, M.; Varret, F. Light-induced bistability in spin transition solids leading to thermal and optical hysteresis. *Eur. Phys. J. B* **1998**, *6*, 183–193. [CrossRef]

12. Müller, E.W.; Ensling, J.; Spiering, H.; Gütlich, P. High-Spin ⇌ Low-Spin Transition in Hexacoordinate Complexes of Iron(II) with Monodentate 1-Alkyltetrazole Ligands: A Variable-Temperature Mossbauer, Magnetic Susceptibility, and Far-Infrared Study. *Inorg. Chem.* **1983**, *22*, 2074–2078. [CrossRef]

13. Vos, G.; Le Febre, R.A.; de Graaff, R.A.G.; Haasnoot, J.G.; Reedijk, J. Unique High-Spin-Low-Spin Transition of the Central Ion in a Linear, Trinuclear Iron (II) Triazole Compound. *J. Am. Chem. Soc.* **1983**, *105*, 1682–1683. [CrossRef]

14. Decurtins, S.; Gütlich, P.; Hasselbach, K.M.; Hauser, A.; Spiering, H. Light-Induced Excited-Spin-State Trapping In Iron(II) Spin-Crossover Systems—Optical Spectroscopic and Magnetic-Susceptibility Study. *Inorg. Chem.* **1985**, *24*, 2174–2178. [CrossRef]

15. Rotaru, A.; Varret, F.; Gindulescu, A.; Linares, J.; Stancu, A.; Létard, J.F.; Forestier, T.; Etrillard, C. Size Effect in Spin-Crossover Systems Investigated by FORC Measurements, for Surfacted [Fe(NH₂ - trz)₃](Br)₂·3H₂O Nanoparticles: Reversible Contributions and Critical Size. *Eur. Phys. J. B Condens. Matter Complex Syst.* **2011**, *84*, 439–449. [CrossRef]

16. Gawali-Salunke, S.; Varret, F.; Maurin, I.; Enachescu, C.; Malarova, M.; Boukheddaden, K.; Codjovi, E.; Tokoro, H.; Ohkoshi, S.; Hashimoto, K. Magnetic and Mössbauer Investigation of the Photomagnetic Prussian Blue Analogue Na₀.₃₂Co[Fe(CN)₆]₀.₇₄·3·4H₂O: Cooperative Relaxation of the Thermally Quenched State. *J. Phys. Chem. B* **2005**, *109*, 8251–8256. [CrossRef]

17. Varret, F.; Chong, C.; Goujon, A.; Boukheddaden, K. Light-Induced Phase Separation (LIPS) in [Fe(Ptz)₆](BF₄)₂ Spin-Crossover Single Crystals: Experimental Data Revisited through Optical Microscope Investigation. *J. Phys. Conf. Ser.* **2009**, *148*, 012036. [CrossRef]
18. Goujon, A.; Varret, F.; Boukheddaden, K.; Chong, C.; Jeftić, J.; Garcia, Y.; Naik, A.D.; Ameline, J.C.; Collet, E. An Optical Microscope Study of Photo-Switching and Relaxation in Single Crystals of the Spin Transition Solid [Fe(Ptz)₆](BF₄)₂, with Image Processing. *Inorg. Chim. Acta* **2008**, *361*, 4055–4064. [CrossRef]

19. Chong, C.; Mishra, H.; Boukheddaden, K.; Denise, S.; Bouchez, G.; Collet, E.; Ameline, J.C.; Naik, A.D.; Garcia, Y.; Varret, F. Electronic and Structural Aspects of Spin Transitions Observed by Optical Microscopy. The Case of [Fe(Ptz)₆](BF₄)₂. *J. Phys. Chem. B* **2010**, *114*, 1975–1984. [CrossRef]

20. Slimani, A.; Varret, F.; Boukheddaden, K.; Chong, C.; Mishra, H.; Haasnoot, J.; Pillet, S. Visualization and Quantitative Analysis of Spatiotemporal Behavior in a First-Order Thermal Spin Transition: A Stress-Driven Multiscale Process. *Phys. Rev. B* **2011**, *84*, 094442. [CrossRef]

21. Varret, F.; Slimani, A.; Boukheddaden, K.; Chong, C.; Mishra, H.; Collet, E.; Haasnoot, J.; Pillet, S. The Propagation of the Thermal Spin Transition of [Fe(Btr)₂(NCS)₂]·H₂O Single Crystals, Observed by Optical Microscopy. *New J. Chem.* **2011**, *35*, 2333–2340. [CrossRef]

22. Chong, C.; Slimani, A.; Varret, F.; Boukheddaden, K.; Collet, E.; Ameline, J.-C.; Bronisz, R.; Hauser, A. The Kinetics Features of a Thermal Spin Transition characterized by Optical Microscopy on the Example of [Fe(bbtr)₃](CLO₄)₂ Single Crystals: Size Effect and Mechanical Instability. *Chem. Phys. Lett.* **2011**, *504*, 29–33. [CrossRef]

23. Sy, M.; Varret, F.; Boukheddaden, K.; Bouchez, G.; Marrot, J.; Kawata, K.; Kaizaki, S. Structure-Driven Orientation of the High-Spin–Low-Spin Interface in a Spin-Crossover Single Crystal. *Angew. Chem. Int. Ed.* **2014**, *53*, 7539–7542. [CrossRef]

24. Pillet, S.; Hubsch, J.; Lecomte, C. Single Crystal Diffraction Analysis of the Thermal Spin Conversion in [Fe(btr)₂(NCS)₂](H₂O) Evidence for Spin-like Domain Formation. *Eur. Phys. J. B-Condens. Matter Complex Syst.* **2004**, *38*, 541–552. [CrossRef]

25. Molnár, G.; Niel, V.; Real, J.A.; Dubrovinsky, L.; Bousseksou, A.; McGarvey, J.J. Raman Spectroscopic Study of Pressure Effects on the Spin-Crossover Coordination Polymers Fe(Pyrazine)[M(CN)₄]·2H₂O (M = Ni, Pd, Pt). First Observation of a Piezo-Hysteresis Loop at Room Temperature. *J. Phys. Chem. B* **2003**, *107*, 3149–3155. [CrossRef]

26. Gütlich, P. Spin Crossover In Iron(II)-Complexes. In *Metal Complexes*; Structure and Bonding; Springer: Berlin/Heidelberg, Germany, 1981; Volume 44, pp. 83–195.

27. Bousseksou, A.; Boukheddaden, K.; Goiran, M.; Consejo, J.; Boillot, M.L.; Tuchagues, J.P. Dynamic Response of the Spin-Crossover Solid Co(H₂(fsa)₂en)(py)₂ to a Pulsed Magnetic Field. *Phys. Rev. B* **2002**, *65*, 172412. [CrossRef]

28. Prins, F.; Monrabal-Capilla, M.; Osorio, E.A.; Coronado, E.; Van der Zant, H.S.J. Room-Temperature Electrical Addressing of a Bistable Spin-Crossover Molecular System. *Adv. Mater.* **2011**, *23*, 1545–1549. [CrossRef]

29. Castro, M.; Roubeau, O.; Piñeiro-López, L.; Real, J.A.; Rodríguez-Velamazán, J.A. Pulsed-Laser Switching in the Bistability Domain of a Cooperative Spin Crossover Compound: A Critical Study through Calorimetry. *J. Phys. Chem. C* **2015**, *119*, 17334–17343. [CrossRef]

30. Mishra, V.; Mukherjee, R.; Linares, J.; Balde, C.; Desplanches, C.; Letard, J.F.; Collet, E.; Toupet, L.; Castro, M.; Varret, F. Temperature-Dependent Interactions and Disorder in the Spin-Transition [Fe<sup>II</sup>(L)₂[CLO₄]₂·C₇H₈ through Structural, Calorimetric, Magnetic, Photomagnetic, and Diffuse Reflectance Investigations. *Inorg. Chem.* **2008**, *47*, 7577–7587. [CrossRef]

31. Raza, Y.; Volatron, F.; Moldovan, S.; Ersen, O.; Huc, V.; Martini, C.; Brisset, F.; Gloter, A.; Stéphan, O.; Bousseksou, A.; et al. Matrix-Dependent Cooperativity in Spin Crossover Fe(Pyrazine)Pt(CN)₄ Nanoparticles. *Chem. Commun.* **2011**, *47*, 11501–11503. [CrossRef]

32. Boukheddaden, K.; Linares, J.; Codjovi, E.; Varret, F.; Niel, V.; Real, J.A. Dynamical Ising-like Model for the Two-Step Spin-Crossover Systems. *J. Appl. Phys.* **2003**, *93*, 7103–7105. [CrossRef]

33. Paez-Espejo, M.; Sy, M.; Boukheddaden, K. Elastic Frustration Causing Two-Step and Multistep Transitions in Spin-Crossover Solids: Emergence of Complex Antiferroelastic Structures. *J. Am. Chem. Soc.* **2016**, *138*, 3202–3210. [CrossRef] [PubMed]

34. Traiche, R.; Sy, M.; Boukheddaden, K. Elastic Frustration in 1D Spin-Crossover Chains: Evidence of Multi-Step Transitions and Self-Organizations of the Spin States. *J. Phys. Chem. C* **2018**, *122*, 4083–4096. [CrossRef]

35. Spiering, H.; Kohlhaas, T.; Romstedt, H.; Hauser, A.; Bruns-Yilmaz, C.; Kusz, J.; Gütlich, P. Correlations of the Distribution of Spin States in Spin Crossover Compounds. *Coord. Chem. Rev.* **1999**, *62*, 190–192. [CrossRef]

36. Kohlhaas, T.; Spiering, H.; Gütlich, P. Monte Carlo Study of the Two-Step Spin Transition in [FeₓZn₁⁻ₓ(2-pic)₃]Cl₂·EtOH. *Z. Phys. B Condens. Matter* **1997**, *102*, 455–459. [CrossRef]

37. Hôo, B.; Boukheddaden, K.; Varret, F. Two-Variable Macroscopic Model for Spin-Crossover Solids: Static and Dynamic Effects of the Correlations. *Eur. Phys. J. B* **2000**, *17*, 449–457. [CrossRef]

38. Gütlich, P.; Hauser, A.; Spiering, H. Thermal and Optical Switching of Iron(II) Complexes. *Angew. Chem. Int. Ed. Engl.* **1994**, *33*, 2024–2054. [CrossRef]

39. Muraoka, A.; Boukheddaden, K.; Linarès, J.; Varret, F.; Linares, J.; Varret, F. Two-Dimensional Ising-like Model with Specific Edge Effects for Spin-Crossover Nanoparticles: A Monte Carlo Study. *Phys. Rev. B* **2011**, *84*, 54119. [CrossRef]

40. Nasser, J.A. First Order High-Spin/Low-Spin Phase Transition Induced by Acoustic-Phonons. *Eur. Phys. J. B* **2001**, *21*, 3–10. [CrossRef]

41. Nasser, J.A.; Boukheddaden, K.; Linares, J. Two-Step Spin Conversion and Other Effects in the Atom-Phonon Coupling Model. *Eur. Phys. J. B* **2004**, *39*, 219–227. [CrossRef]
42. Stoleriu, L.; Chakraborty, P.; Hauser, A.; Stancu, A.; Enachescu, C. Thermal Hysteresis in Spin-Crossover Compounds Studied within the Mechanoelastic Model and Its Potential Application to Nanoparticles. *Phys. Rev. B-Condens. Matter Mater. Phys.* **2011**, *84*, 134102. [CrossRef]

43. Miyashita, S.; Konishi, Y.; Nishino, M.; Tokoro, H.; Rikvold, P.A. Realization of the Mean-Field Universality Class in Spin-Crossover Materials. *Phys. Rev. B* **2008**, *77*, 014105. [CrossRef]

44. Nishino, M.; Boukheddaden, K.; Konishi, Y.; Miyashita, S. Simple Two-Dimensional Model for the Elastic Origin of cooperativity among Spin States of Spin-Crossover Complexes. *Phys. Rev. Lett.* **2007**, *98*, 247203. [CrossRef]

45. Enachescu, C.; Stoleriu, L.; Stancu, A.; Hauser, A. Model for Elastic Relaxation Phenomena in Finite 2D Hexagonal Molecular Lattices. *Phys. Rev. Lett.* **2009**, *102*, 257204. [CrossRef] [PubMed]

46. Slimani, A.; Khemakhem, H.; Boukheddaden, K. Structural Synergy in a Core-Shell Spin Crossover Nanoparticle Investigated by an Electroelastic Model. *Phys. Rev. B* **2017**, *95*, 174104. [CrossRef]

47. Slimani, A.; Boukheddaden, K.; Varret, F.; Oubouchou, H.; Nishino, M.; Miyashita, S. Microscopic Spin-Distortion Model for Switchable Molecular Solids: Spatiotemporal Study of the Deformation Field and Local Stress at the Thermal Spin Transition. *Phys. Rev. B* **2013**, *87*, 14111. [CrossRef]

48. Slimani, A.; Boukheddaden, K.; Yamashita, K. Effect of Intermolecular Interactions on the Nucleation, Growth, and Propagation of like-Spin Domains in Spin-Crossover Materials. *Phys. Rev. B* **2015**, *92*, 014111. [CrossRef]

49. Singh, Y.; Oubouchou, H.; Nishino, M.; Miyashita, S.; Boukheddaden, K. Elastic-Frustration-Driven Unusual Magnetoelastic Properties in a Switchable Core-Shell Spin-Crossover Nanostructure. *Phys. Rev. B* **2020**, *101*, 54105. [CrossRef]

50. Slimani, A.; Boukheddaden, K. An Electro-Elastic Theory for the Mechanically-Assisted Photo-Induced Spin Transition in Core-Shell Spin-Crossover Nanoparticles. *Phys. Chem. Chem. Phys.* **2018**, *20*, 28583–28591. [CrossRef]

51. Enachescu, C.; Nicolazzi, W. Elastic Models, Lattice Dynamics and Finite Size Effects in Molecular Spin Crossover Systems. *Comptes Rendus Chim.* **2018**, *21*, 1179–1195. [CrossRef]

52. Stoleriu, L.; Stancu, A.; Chakraborty, P.; Hauser, A.; Enachescu, C. Analysis of First Order Reversal Curves in the Thermal Hysteresis of Spin-Crossover Nanoparticles within the Mechanoelastic Model. *J. Appl. Phys.* **2015**, *117*, 17B307. [CrossRef]

53. Felix, G.; Nicolazzi, W.; Salmon, L.; Molnar, G.; Perrier, M.; Maurin, G.; Larionava, J.; Long, J.; Guari, Y.; Bousseksou, A. Enhanced Cooperative Interactions at the Nanoscale in Spin-Crossover Materials with a First-Order Phase Transition. *Phys. Rev. Lett.* **2013**, *110*, 235701. [CrossRef] [PubMed]

54. Shepherd, H.J.; Molnar, G.; Nicolazzi, W.; Salmon, L.; Bousseksou, A. Spin Crossover at the Nanometre Scale. *Eur. J. Inorg. Chem.* **2013**, *2013*, 653–661. [CrossRef]

55. Constant-Machado, H.; Linares, J.; Varret, F.; Haasnoot, J.G.; Martin, J.P.; Zarembowitch, J.; Dworkin, A.; Bousseksou, A. Dilution Effects in a Spin Crossover System, Modelled in Terms of Direct and Indirect Intermolecular Interactions. *J. Phys. I* **1996**, *6*, 1203–1216. [CrossRef]

56. Krivokapic, I.; Chakraborty, P.; Bronisz, R.; Enachescu, C.; Hauser, A. Significant Variation of the Singlet–Quintet Intersystem Crossing Rate Constant in an Iron(II) High-Spin Complex as a Function of Temperature. *Angew. Chem. Int. Ed.* **2010**, *49*, 8509–8512. [CrossRef] [PubMed]

57. Jakobi, R.; Spiering, H.; Gütlich, P. Thermodynamics of the Spin Transition in $[Fe_xZn_{1-x}(2-pic)_3]Cl_2 \cdot EtOH$. *J. Phys. Chem. Solids* **1992**, *53*, 267–275. [CrossRef]

58. Baldé, C.; Desplanches, C.; Wattiaux, A.; Guionneau, P.; Gütlich, P.; Létard, J.-F. Effect of Metal Dilution on the Light-Induced Spin Transition in $[Fe_xZn_{1-x}(phen)_2(NCS)_2]$ (phen = 1,10-phenanthroline). *Dalton Trans.* **2008**, *20*, 2702–2707. [CrossRef]

59. Rotaru, A.; Dărtu, M.M.; Enachescu, C.; Tanasa, R.; Linares, J.; Stancu, A.; Garcia, Y. Calorimetric Measurements of Diluted Spin Crossover Complexes $[Fe_xM_{1-x}(btr)_2(NCS)_2] \cdot H_2O$ with MII = Zn and Ni. *Polyhedron* **2009**, *28*, 2531–2536. [CrossRef]

60. Enachescu, C.; Stoleriu, L.; Stancu, A.; Hauser, A. Study of the Relaxation in Diluted Spin Crossover Molecular Magnets in the Framework of the Mechano-Elastic Model. *J. Appl. Phys.* **2011**, *109*, 07B111. [CrossRef]

61. Martin, J.-P.; Zarembowitch, J.; Bousseksou, A.; Dworkin, A.; Haasnoot, J.G.; Varret, F. Solid State Effects on Spin Transitions: Magnetic, Calorimetric, and Moessbauer-Effect Properties of $[Fe_xCo_{1-x}(4,4'-bis-1,2,4-triazole)_2(NCS)_2] \cdot cndot.H_2O$ Mixed-Crystal Compounds. *Inorg. Chem.* **1994**, *33*, 6325. [CrossRef]

62. Sorai, M.; Ensling, J.; Gütlich, P. Mössbauer Effect Study on Low-Spin 1A1 ⇌ High-Spin 5T2 Transition in Tris(2-Picolylamine) Iron Chloride I. Dilution Effect in $[Fe_xZn_{1-x}(2-pic)_3]Ch_2 \cdot C_2H_5OH$. *Chem. Phys.* **1976**, *18*, 199–209. [CrossRef]

63. Zheng, S.; Siegler, M.A.; Costa, J.S.; Fu, W.T.; Bonnet, S. Effect of Metal Dilution on the Thermal Spin Transition of $[Fe_xZn_{1-x}(bapbpy)(NCS)_2]$. *Eur. J. Inorg. Chem.* **2013**, *2013*, 1033–1042. [CrossRef]

64. Ganguli, P.; Gütlich, P.; Mueller, E.W. Effect of Metal Dilution on the Spin-Crossover Behavior in $[Fe_xM_{1-x}(phen)_2(NCS)_2]$ (M = Mn, Co, Ni, Zn). *Inorg. Chem.* **2002**, *21*, 3429–3433. [CrossRef]

65. Martin, J.P.; Zarembowitch, J.; Dworkin, A.; Haasnoot, J.C.; Codjovi, E. Solid-State Effects in Spin Transitions: Influence of Iron(II) Dilution on the Magnetic and Calorimetric Properties of the Series $[Fe_xNi_{1-x}(4,4'-Bis(1,2,4-Triazole))_2(NCS)_2] \cdot Cndot.H_2O$. *Inorg. Chem.* **2002**, *33*, 2617–2623. [CrossRef]

66. Renz, F.; Oshio, H.; Ksenofontov, V.; Waldeck, M.; Spiering, H.; Gütlich, P. Strong Field Iron(II) Complex Converted by Light into a Long-Lived High-Spin State. *Chem. Int. Ed.* **2000**, *39*, 3699–3700. [CrossRef]
67. Decurtins, S.; Gütlich, P.; Köhler, C.P.; Spiering, H.; Hauser, A. Light-Induced Excited Spin State Trapping in a Transition-Metal Complex: The Hexa-1-Propyltetrazole-Iron (II) Tetrafluoroborate Spin-Crossover System. *Chem. Phys. Lett.* **1984**, *105*, 1-4. [CrossRef]

68. Hauser, A. Light-Induced Spin Crossover and the High-Spin→Low-Spin Relaxation. In *Spin Crossover in Transition Metal Compounds*; Springer: Berlin/Heidelberg, Germany, 2004; Volume 234, pp. 155-198. [CrossRef]

69. Hauser, A.; Enachescu, C.; Daku, M.L.; Vargas, A.; Amstutz, N. Low-Temperature Lifetimes of Metastable High-Spin States in Spin-Crossover and in Low-Spin Iron(II) Compounds: The Rule and Exceptions to the Rule. *Coord. Chem. Rev.* **2006**, *250*, 1642-1652. [CrossRef]

70. Hauser, A.; Jeftic, J.; Romstedt, H.; Hinek, R.; Spiering, H. Cooperative Phenomena and Light-Induced Bistability in Iron(II) Spin-Crossover Compounds. *Coord. Chem. Rev.* **1999**, *190-192*, 471-491. [CrossRef]

71. Baldé, C.; Desplanches, C.; Grunert, M.; Wei, Y.; Gütlich, P.; Létard, J.F. Influence of Metal Dilution on the Light-Induced Spin Transition in Two 1D Chain Compounds: [FeₓZn₁₋ₓ(btzp)₃](BF₄)₂ and [FeₓZn₁₋ₓ(endi)₃](BF₄)₂ [btzp = 1,2-Bis(tetrazol-1-yl)Propane and endi = 1,2-Bis(tetrazol-1-yl)ethane]. *Eur. J. Inorg. Chem.* **2008**, *2008*, 5382-5389. [CrossRef]

72. Krivokapic, I.; Chakraborty, P.; Enachescu, C.; Bronisz, R.; Hauser, A. Low-Spin→High-Spin Relaxation Dynamics in the Highly Diluted Spin-Crossover System [FeₓZn₁₋ₓ(bbtr)₃](ClO₄)₂. *Inorg. Chem.* **2011**, *50*, 1856-1861. [CrossRef]

73. Yu, W.C.; Gielisse, P.J. High Pressure Polymorphism in CdS, CdSe and CdTe. *Mater. Res. Bull.* **1971**, *6*, 621-638. [CrossRef]

74. Sorai, M. Heat Capacity Studies of Spin Crossover Systems. In *Spin Crossover in Transition Metal Compounds III*; Springer: Berlin/Heidelberg, Germany, 2006; pp. 153-170. [CrossRef]

75. Bousseksou, A.; Nasser, J.; Linares, J.; Boukheddaden, K.; Varret, F. Ising-like Model for the Two-Step Spin-Crossover. *J. Phys.* **1992**, *2*, 1381-1403. [CrossRef]

76. Jung, J.; Bruchhäuser, F.; Feile, R.; Spiering, H.; Gütlich, P. The Cooperative Spin Transition in [FeₓZn₁₋ₓ(Ptz)₆](BF₄)₂: I. Elastic Properties-An Oriented Sample Rotation Study by Brillouin Spectroscopy. *Z. Phys. B-Condens. Matter* **1996**, *100*, 517-522. [CrossRef]

77. Spiering, H.; Boukheddaden, K.; Linares, J.; Varret, F. Total Free Energy of a Spin-Crossover Molecular System. *Phys. Rev.* **2004**, *70*, 184106. [CrossRef]

78. Felts, A.C.; Slimani, A.; Cain, J.M.; Andrus, M.J.; Ahir, A.R.; Abboud, K.A.; Meisel, M.W.; Boukheddaden, K.; Talham, D.R. Control of the Speed of a Light-Induced Spin Transition through Mesoscale Core-Shell Architecture. *J. Am. Chem. Soc.* **2018**, *140*, 5814-5824. [CrossRef]

79. Ohkoshi, S.I.; Shiraishi, K.; Nakagawa, K.; Ikeda, Y.; Stefanczyk, O.; Tokoro, H.; Namai, A. Reversible photoswitchable ferromagnetic thin film based on a cyanido-bridged RbCuMo complex. *J. Mater. Chem. C* **2021**, *9*, 3081-3087. [CrossRef]

80. Itoi, M.; Maurin, I.; Boukheddaden, K.; Andrus, M.J.; Talham, D.R.; Elkaim, E.; Uwatoko, Y. Sub-micrometer particle size effects on metastable phases for a photoswitchable Co-Fe Prussian blue analog. *J. Appl. Phys.* **2022**, *131*, 085110. [CrossRef]

81. Kusz, J.; Gütlich, P.; Spiering, H. Structural Investigations of Tetrazole Complexes of Iron(II). In *Spin Crossover in Transition Metal Compounds II*; Topics in Current Chemistry; Springer: Berlin/Heidelberg, Germany, 2004; Volume 234, pp. 129-153.

82. Chakraborty, P.; Sy, M.; Fourati, H.; Delgado, T.; Dutta, M.; Das, C.; Besnard, C.; Hauser, A.; Enachescu, C.; Boukheddaden, K. Optical microscopy imaging of the thermally-induced spin transition and isothermal multi-stepped relaxation in a low-spin stabilized spin-crossover material. *Phys. Chem. Chem. Phys.* **2022**, *24*, 982-994. [CrossRef] [PubMed]

83. Slimani, A.; Boukheddaden, K. Microscopic Electroelastic Modeling of a Spin Transition with Symmetry Breaking. *J. Appl. Phys.* **2021**, *129*, 173901. [CrossRef]

84. Dírtu, M.; Schmit, F.; Naik, A.D.; Rusu, I.; Rotaru, A.; Rackwitz, S.; Wolny, J.A.; Schönemann, V.; Spinu, L.; Garcia, Y. Two-Step Spin Transition in a 1D FeII 1,2,4-Triazole Chain Compound. *Eur. J. Chem.* **2015**, *21*, 5843-5855. [CrossRef]

85. Pittala, N.; Thétiot, F.; Charles, C.; Triki, S.; Boukheddaden, K.; Chastanet, G.; Marchivie, M. An unprecedented trinuclear Fe(II) triazole-based complex exhibiting a concerted and complete sharp spin transition above room temperature. *Chem. Commun.* **2017**, *53*, 8356-8359. [CrossRef]

86. Pittala, N.; Cuza, E.; Pinkowicz, D.; Magott, M.; Marchivie, M.; Boukheddaden, K.; Triki, S. Antagonist elastic interactions tuning spin crossover and LIESST behaviours in FeII trinuclear-based one-dimensional chains. *Inorg. Chem. Front.* 2022; in press. [CrossRef]

87. Traiche, R.; Oubouchou, H.; Boukheddaden, K. Elastic Modeling of Thermal Properties in Sterically Frustrated 1D binuclear Spin-Crossover Chains. *Symmetry* **2021**, *13*, 1836. [CrossRef]

88. Di Scala, N.; Belmontouri, N.E.I.; Espejo, M.A.P.; Boukheddaden, K. Three-dimensional electroelastic modeling of the nucleation and propagation of the spin domains in spin-crossover materials. *Phys. Rev. B* **2022**, *106*, 014422. [CrossRef]

89. di Scala, N.; Belmontouri, N.E.I.; Espejo, M.A.P.; Boukheddaden, K. Explaining the origin of the orientation of the front transformation in spin-transition crystals. *Phys. Rev. B* **2022**, *106*, 144107. [CrossRef]