Rachid Traiche

Hassane Oubouchou

Mourad Zergoug

Kamel Boukheddaden
email: kamel.boukheddaden@uvsq.fr

Monte Mc

; Carlo

Hs

Spatio-temporal aspects of the domain propagation in a spin-crossover lattice with defect

Keywords:

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.

Introduction

Spin-crossover solids (SCO) solids [START_REF] König | Mössbauer effect and magnetism down to 1.2 K in the triplet ground state of an iron(II)-Phenanthroline complex[END_REF][START_REF] Gütlich | Thermal and optical switching of Iron(II) complexes[END_REF][START_REF] Kahn | Molecular bistability: an overview[END_REF], are switchable molecular solids undergoing electronic, magnetic, mechanical and optical changes at the spin transition, and then can be considered as multifunctional materials. They are well known to have potential applications in several areas, like sensing of temperature, pressure [START_REF] Linares | Pressure and Temperature spin-crossover sensors with optical detection[END_REF] or as displays. From the theoretical point view, they have been modelled using regular solutions models [START_REF] Slichter | Pressure-induced electronic changes in compounds of iron[END_REF][START_REF] Wajnflasz | Transitions «Low spin»-«High spin» dans les complexes de Fe 2+[END_REF] which are based on Bragg-Williams approximations, or using Ising-like approaches [START_REF] Doniach | Thermodynamic fluctuations in phospholipid bilayers[END_REF][START_REF] Boukheddaden | Dynamical model for spin-crossover solids. I. Relaxation effects in the mean-field approach[END_REF][9][START_REF] Varret | Light-induced thermal hysteresis and intensity thresholds in molecular switchable solids, by mean-field macroscopic master equation approach: discussion of the experimental data obtained for Co-Fe Prussian Blue Analogues[END_REF], whose mean-field resolutions lead to exactly the same equations as those of regular solutions. Despite of their phenomenological nature, these simple descriptions were able to catch essential aspects of the SCO transition, since they reproduce the conditions of existence of the gradual and first-order transition in these molecular solids. However, since they discard the lattice volume change at the transition [START_REF] Spiering | The elastic interaction of high-spin and low-spin complex molecules in spin-crossover compounds[END_REF], these descriptions remained too qualitative, particularly on the explanation of the true physical origin of the SCO phenomenon. Indeed, in Ising-like models, the first-order thermal transition is obtained through a phenomenological interaction parameter, with electronic nature. In contrast, recent theoretical developments, suggesting a description of the SCO transition using atomistic models based on deformable lattices [START_REF] Miyashita | Realization of the mean-field universality class in spin-crossover materials[END_REF][START_REF] Nishino | Simple two-dimensional model for the elastic origin of cooperativity among spin states of spin-crossover complexes[END_REF][START_REF] Nishino | Intrinsic effects of the boundary condition on switching processes in effective long-range interactions originating from local structural change[END_REF][START_REF] Stoleriu | Thermal hysteresis in spin-crossover compounds studied within the mechanoelastic model and its potential application to nanoparticles[END_REF][START_REF] Enachescu | Model for elastic relaxation phenomena in finite 2D hexagonal molecular lattices[END_REF][START_REF] Nicolazzi | Two-variable anharmonic model for spincrossover solids: a like-spin domains interpretation[END_REF][START_REF] Nishino | Macroscopic nucleation phenomena in continuum media with long-range interactions[END_REF][START_REF] Klinduhov | Vibronic theory of ultrafast intersystem crossing dynamics in a single spin-crossover molecule at finite temperature beyond the Born-Oppenheimer approximation[END_REF][START_REF] Espejo-Parez | Elastic frustration causing two-step and multi-step transitions in spin-crossover solids: emergence of complex antiferroelastic structures[END_REF], were able to retrieve the experimental spatiotemporal behaviors of the SCO transition, observed by optical microscopy [START_REF] Varret | The propagation of the thermal spin transition of [Fe(btr) 2 (NCS) 2 ]•H 2 O single crystals, observed by optical microscopy[END_REF][START_REF] Chong | Electronic and structural aspects of spin transitions observed by optical microscopy. the case of [Fe(ptz) 6 ](BF 4 ) 2[END_REF][START_REF] Slimani | Visualization and quantitative analysis of spatiotemporal behavior in a first-order thermal spin transition: a stress-driven multiscale process[END_REF] on a unique single crystal. These relevant models allowed to reproduce the main features of the spatiotemporal behaviors of SCO solids, which can be summarized as follows: the transition from the high-spin (HS) to low-spin (LS) usually starts from a corner of the crystal, and spreads over the whole material with a well-defined HS/LS interface whose shape and orientation depend on the structural symmetry of the said crystal. In most of the cases, the interface propagates at a very slow velocity (2-10 µm s -1 ), almost nine-orders of magnitude smaller than that of sound (~3000 m s -1 ). This slow velocity of transformation is attributed to the non-elastic character of the spin-transition which is accompanied with a volume change at the macroscopic scale. Our previous optical microscopy investigations demonstrated the presence of stresses which strongly impact the transition temperature and may lead in some cases to irreversible damage of the crystals after several thermal cycling [START_REF] Slimani | Visualization and quantitative analysis of spatiotemporal behavior in a first-order thermal spin transition: a stress-driven multiscale process[END_REF][START_REF] Slimani | Velocity of the high-spin low-spin interface inside the thermal hysteresis loop of a spin-crossover crystal, via photothermal control of the interface motion[END_REF]. The experimental studies clearly underlined that the nucleation and growth process in a spin-crossover single crystal at the thermal first-order transition is a multi-scale (space-time) process driven by the propagation of mechanical stresses ahead of the transformation front-line. This specific character confirms the crucial role of the coupling between the spin state change and the local deformation (molecular volume change) of the lattice, as confirmed by the theoretical investigations performed using the electro-elastic model [START_REF] Slimani | Microscopic spin-distortion model for switchable molecular solids: spatiotemporal study of the deformation field and local stress at the thermal spin transition[END_REF], which has been solved for several crystal shapes [START_REF] Slimani | Microscopic spin-distortion model for switchable molecular solids: spatiotemporal study of the deformation field and local stress at the thermal spin transition[END_REF][START_REF] Slimani | Thermal spin transition of circularly shaped nanoparticles in a core-shell structure investigated with an electroelastic model[END_REF][START_REF] Slimani | Effect of intermolecular interactions on the nucleation, growth, and propagation of like-spin domains in spin-crossover materials[END_REF][START_REF] Boukhedadden | Physical properties of 2D spin-crossover solids from an electroelastic description: effect of shape, size, and spin-distortion interactions[END_REF]. The local volume change is then delocalized over long distances, leading to long-range interactions between the spin-states through the lattice medium. It is also worth mentioning that the role of defects was also pointed out in the field of 1D coordination polymers where the residual high-spin fraction at low temperature is generally attributed to the presence of crystal defects or to the finite size effects of the 1D chains [START_REF] Van Koningsbruggen | Synthesis, crystal structure, EXAFS, and magnetic properties of catena [ < mu > -Tris(1,2-bis(tetrazol-1-yl)propane-N1,N1')iron(II)] Bis(perchlorate). First crystal structure of an Iron(II) spin-crossover chain compound[END_REF].

Thus, optical microscopy observations of the spin transition in single crystals allowed to suspect the presence of mechanical constraints that control the spread of the HS / LS interface. For this purpose, the present work will be devoted to complete the study that we recently published in Ref. [30] on the spatiotemporal behaviors of the HS and LS domains during the relaxation of the photo-induced metastable HS state obtained after photo-excitation through LIESST (Light-Induced Excited Spin State Trapping) effect [START_REF] Decurtins | Light-induced excited spin state trapping in a transition-metal complex-The hexa-1-propyltetrazoleiron(II) tetrafluoroborate spin-crossover system[END_REF] at very low temperature. The investigations are restricted to a 2D rectangular lattice in the presence of a defect, which consists in a hole, injected in the center of the lattice. Here, we are mainly interested in (i) the effect of the hole on the velocity of the HS/LS interface during the relaxation processes, (ii) on the energetic aspects of the HS-LS transformation in the presence of the defect, (iii) on the mechanical properties of this transformation through the analysis the divergence and the rotational of the displacement field, with the goal to clarify the role of the elastic strain in the pinning process of the HS/LS interface around the defect.

The paper is organized as follows: in Section 2,t h em o d e la n d simulation method are described; Section 3 is devoted to the presentation of the obtained results; in Section 4 we summarize the main conclusions and outline some possible developments of the present work.

The model

To describe theoretically the SCO solids, we consider a model which couples electronic and elastic properties of a SCO lattice for which we have already discussed some of its thermodynamic properties in recent works [START_REF] Slimani | Visualization and quantitative analysis of spatiotemporal behavior in a first-order thermal spin transition: a stress-driven multiscale process[END_REF][START_REF] Slimani | Microscopic spin-distortion model for switchable molecular solids: spatiotemporal study of the deformation field and local stress at the thermal spin transition[END_REF]. This model is based on the description of the HS and LS states of the SCO molecule by a two-states fictitious spin, S , . The molecules (here, the sites) interact via elastic springs, as depicted in Fig. 1.

The Hamiltonian of the elastic model including both electronic and elastic degrees of freedom of SCO solids writes as follows

∑∑ ∑ H kT g S A r R S S Br R S S = 1 2 [∆ - ln ( )] + [ -( , )] +[ -′ (, ) ] i Bi ij ij i j ik ik i k (,) 0 2 (, ) 0 2 (1)
The first contribution in (1) is the effective ligand field energy which splits the degeneracy between the LS and HS states of an isolated molecule, where, Δ, comes the difference of ligand fields in the HS and LS states, and kT g -l n

B is an entropic contribution (electronic and vibrational) arising for the difference of degeneracies between the LS and the HS states, whose ratio is g. The second and third contributions account for elastic interactions, between nearest-neighbors (nn) and next-nearest-neighbors (nnn) of spin-crossover units, respectively, which can be easily identified in Fig. 1. There, RSS (,

) ij 0 (resp. RSS ′ (, ) ik 0
) is the equilibrium distance between the nn sites i and j (respectively, i and k). A (resp. B) are the elastic strengths of springs connecting nn and nnn, respectively. The quantity r ij (resp. r ik ) is the instantaneous distance between the nn (resp. nnn) nodes i and j (i and k).

Since two connected SCO particles may have only four configurations, we denote by, R HH 0 , R HL 0 ( R = LH 0 ) and R LL 0 , the respective nn distances between two HS-HS, HS-LS (LS-HS), and LS-LS interacting neighboring SCO sites. In addition, we assume for simplicity that, RR R =( + )/2

Here, δ R is the lattice parameter misfit between the LS and HS states. According to these expressions, we then have:

RR =( + 1 , + 1 ) HH 0 0 , RR =( -1 , -1 ) LL 0 0
, and RR =( + 1 , -1 ) .

HL 0 0

In Ref. [START_REF] Nicolazzi | Two-variable anharmonic model for spincrossover solids: a like-spin domains interpretation[END_REF], we have demonstrated that the present elastic model can be mapped as a generalized Ising-like model with local interactions and local fields which depend on the positions of the spin states. The main results of this analysis can be quoted in the following. This feature allowed to demonstrate the following points:

(i) the equilibrium transition temperature, denoted T eq , is found to be equal to that of the usual Ising-like model, i.e.,

T kl n g = ∆ ;

eq B

(3)

(ii) the condition of the existence of a first-order transition is

k TA >( ) Be q δ 4 2 R , where δR R =- R HH LL 00
is the lattice misfit between the HS and LS states.

The simulations performed in the present work have been done using the following model parameters:

∆ = 450 K, g = 150,, R =1.2 HH 0 nm, R =1 LL 0 nm, R =1.1 HL 0 nm (nnn equilibrium distances=nn equilibrium distances × 2 ) and AB = = 23000 Knm -2 .
From these parameter values, the derived bulk modulus value, is ~6 GPa, and the transition temperature is evaluated as T ≈90 eq K.

Results and discussion

We investigate the relaxation of the photo-induced HS metastable state at low temperatures of a planar rectangular lattice with NN N =× = 6 0 × 2 0 xy particles with a defect (a hole) in the center. We consider free boundary conditions in the simulations, and the calculations are performed at T=10 K. For so doing, we prepare all sites in the HS state, by initializing all spin values at S =+1 i and all nn distances r R == 1 . 2 ij HH 0 nm. We then perform usual MC simulations (using Metropolis algorithm) on the spin states and also on lattice positions to relax the total elastic energy. After each MC step, we measure the average value of the HS fraction as well as the average nn distances.

In practice, the MC simulations are conducted on both electronic and structural degrees of freedom, which corresponds to change the fictitious spins variables (spin flip) and to move the atomic positions. The technical procedure is as follows: a site of coordinates ij (, ) is randomly chosen, its associated spin state is flipped using the conventional MC Metropolis algorithm. After this first step, we visit each lattice site and we change its coordinates in a random direction with a distance δu =0 . 0 1 5 nm, which is very small compared to the lattice parameter values in the HS (1.2 nm) and LS (1 nm) states. The new coordinates are accepted or rejected following the usual Monte Carlo Metropolis algorithm, after calculating the elastic energy change associated to the displacement of the atom. If the process is rejected, the atom is left in its original position. This mechanical relaxation continues for all the sites of the lattice. This step corresponds to the mechanical relaxation of the lattice positions and it is repeated 10 times, for each spin update which ensures reaching the stable mechanical state. The previous procedure is repeated for each of the NN × x spins, which defines one MC step (MCS). At the end, in 1 MCS, each atomic positions has been updated 12,000 (60 × 20 × 10) times. To reach the thermodynamic equilibrium state, 10 3 MCS are performed at each temperature, followed by a second 10 3 MCS which are utilized to evaluate the averages of physical quantities of interest, like the highspin fraction, the average lattice parameter, elastic energy, etc.

Fig. 2a displays some selected snapshots of the electro-elastic configuration of the lattice, summarizing the spatiotemporal features of the spin and structural transformations during the thermal relaxation of the lattice with defect published in [START_REF] Van Koningsbruggen | Synthesis, crystal structure, EXAFS, and magnetic properties of catena [ < mu > -Tris(1,2-bis(tetrazol-1-yl)propane-N1,N1')iron(II)] Bis(perchlorate). First crystal structure of an Iron(II) spin-crossover chain compound[END_REF]. One can easily see that the HS (red) to LS (blue) transformation starts from bottom left corner of the lattice (snapshot B) with a clearly curved interface, obtained after 1944 MCS. This interface shape is the ones which minimizes the elastic energy. Concomitantly, we may remark the presence of some nucleus of the stable LS phase (red points) on the surface of the lattice as well as around the defect. Their existence can be explained as due to the lack of neighbors which then facilitates the stochastic nucleation of the LS state. In the course of the process, the transformation spreads out over the whole lattice (snapshots C-F). First, the relaxation is dominated by a transformation along the transverse direction following the shortest distance to the surface. Then a longitudinal propagation proceeds in the lattice with a characteristic arrow-shaped HS/LS interface, displayed in snapshot C, which is obtained after 2986MCS. The particular interface's shape is due to the presence of the defect, which accelerates the relaxation of the HS sites located on both sides of y = N 2 y , because they are situated at the shortest distances from the hole, which then plays the role of an additional surface, where the excess of elastic energy is relaxed. Snapshot D is reached after 6709 and from D to E, the HS-LS interface rides very slowly between the hole and the lattice border, since snapshot E is obtained only after 11529 MCS. In addition, and probably due to the dynamic slowing down of the interface velocity, the nucleation of the stable LS phase has time to start from ahead the front as shown in snapshots D-F. Finally a second new nucleation site appears in the top right corner of the lattice (snapshot G) which complicates the final steps of the relaxation, as well as the shape of the lattice. In contrast, simulations performed on a perfect lattice, using the same model parameter, shown in Fig. 2b, display a regular and straight interfaces in the region of the cruise climb.

To estimate the interface velocity, we determined the elastic profiles of the interface along a reference line, y = N 2 y , in the direction of propagation. In practice, the elastic, interface profile is obtained by plotting the spatial dependence, along the interface propagation direction (here, x axis), of the distance, d between two successive sites, located along the reference line, y = N 2 y .

We illustrate in Fig. 3a some spatial profiles corresponding to the interface propagation of Fig. 2a at several selected times. The three dimensional profile of the elastic HS-LS interface is presented in Fig. 3b at time t = 3086 MCS. Far from the hole, the interface is sharp and its velocity (see Fig. 4) is constant. When approaching the defect region, which can be considered here as an additional surface at which the system can relax the excess of elastic energy, we do not see any significant effect on the curves of dx ( ), which globally keep the same shape, except at t = 3241 MCS, where the HS-LS interface crosses the hole's region.

On the other hand, one can notice in Fig. 3, the presence in the left (resp. right) side of the elastic interface, of a region with dR < LS 0 (resp. dR > HS 0 ) which clearly indicates the specific nature of the elastic interface which stores the elastic strain due to the misfito f lattice parameter between the LS and the HS states. Far from the interface, this misfit becomes zero and the system recovers the equilibrium distances.

The interface positions of Fig. 2 as function of the Monte Carlo time are reported in Fig. 4. We note the presence of five regimes of propagation. At the beginning, the growth starts from one corner of the rectangle along the width; this regime, which is quite short and was not depicted in Fig. 2, is the regime of the interface creation. So, we denote as the regime (1), the flow regime of longitudinal propagation at almost constant velocity V =0.01 1 nm/MCS. Then a second regime [START_REF] Gütlich | Thermal and optical switching of Iron(II) complexes[END_REF] where the propagation significantly slows down, takes place when the HS/LS interface's gets the hole. There, the velocity is nonzero but very small, V = 0.00025 2 , nm/MCS, and one can say that the interface can be considered as stapled around the hole. This phenomenon takes place at t~3000 MCS and ends at t = 7000 MCS. The third regime (3) corresponds to the time region where the interface comes off from the hole, at relatively low velocity V =0.001 3 nm/MCS. Finally, the last regime (4) corresponds to a rapid motion of the interface, V =0.018 4 nm/MCS, when the elastic interface "feels" the lattice surface border. Now, we discuss the energetic aspects of this transformation. In Fig. 5, we present the spatial distribution of the local elastic energy at several time values along the relaxation process of the HS fraction of Fig. 2. The local elastic energy, at a chosen site i, is given by:

∑∑ E Ar R S S BrR S S =[ -( , ) ] +[ -′ (, ) ]. i j ij i j k ik i k () 0 2 () 0 2 (4)
The inspection of the maps of Fig. 5 indicates that the elastic energy is stored around the HS/LS interface regions (maps B and C), where the lattice parameter misfit reaches it maximum value. Moreover, around the defect, the elastic energy is easily relaxed (picture D), since the hole introduces an additional surface (or border) in the lattice where the energy excess can be released. When the interface overpasses the defect, the elastic energy increases again (map F) and enhances in snapshot G where the two interfaces coming from the left and the right tips of the lattice meet each other. At this point again, the HS domain, surrounded from both sides by LS phases, experiences a strong elastic strain. To illustrate quantitatively these facts, we have calculated and plotted in Fig. 6 the time-dependence of the total elastic energy during the HS to LS relaxation process, for (i) a the lattice with a defect (red curve) and also for (ii) the perfect lattice (blue curve). It is easily remarked that the total elastic energy of the red curve is smaller than that of the blue curve around the defect region, and waves in the course of the process, due to the effect of the hole. In addition, the elastic energy of the lattice with defect shows two peaks, clearly related to the elastic energy excess, first in the beginning of the process, and last when the interface gets and leaves the hole, which slows down the interface motion. Let us briefly discuss the elastic energy of the perfect lattice in the metastable HS state. When producing the metastable HS state at low-temperature, we prepare all sites with a spin S =+ 1 i and in addition all distances between neighboring sites are set equal to the HS equilibrium distance, that is, R HH 0 . Within these spin and lattice configurations, the total elastic energy is zero because there is no elastic misfit, however the electronic energy is,

kT g S kT g ∑ [∆ - ln ( )] = [∆ - ln ( )] i N Bi N B =1 1 22
. Since the electronic energy of the uniform LS state (S =+ 1

i ) with R LL 0 distances is kT g -[ ∆ - l n () ] N B 2
and is lower, the HS lattice will then irremediably relax to the LS state. During the electronic relaxation, any spin flip from S =+1 i to S =-1 i at fixed lattice coordinates (those of the HS state) will cause an elastic misfit since the equilibrium distances between the site i and its closest neighbors will be will be R HL 0 and not R LL 0 . As a result the lattice positions will move in order to satisfy this new condition which minimizes the local elastic energy. After the total electronic and mechanical relaxation, the system reaches the LS state, where S =-1 i i ∀ and all nearest-neighbors distances are equal to, R LL 0 . Since all positions are at equilibrium, the elastic energy of the LS state is also zero.

To go further ahead in the analysis of the dynamics of this process, we have also studied the distribution of the local pressure in the lattice, calculated through the following expression,

∑∑ PA r R S SB r R S S =- [ -( , )]- [ -( , )], i j ij ij i j k ik ik i k 00 (5)
In Eq. ( 5) the local pressure can be viewed as the gradient of the local elastic energy at site i. We see in Eq. ( 5), that the local pressure can be positive or negative, according to the stress exerted on the site. If the distance to the neighbors is smaller than the expected equilibrium distance, the system will experience a positive strain and the site will be under pressure. In contrast, if r RSS >( , )

ij i j 0
, the site i will be under tension and will feel a negative local pressure. Fig. 7 displays the Fig. 5. Maps of the local elastic energy at several positions along the relaxation process of Fig. 2, for the lattice with defect (a), and for a perfect lattice (b). In both cases, the elastic energy is mostly stored around the HS/LS interface, except for the case (a) where it also appears in the neighborhood of the defect region where the latter can be released. Fig. 6. After [START_REF] Slimani | Velocity of the high-spin low-spin interface inside the thermal hysteresis loop of a spin-crossover crystal, via photothermal control of the interface motion[END_REF]. Time dependence of the total elastic energy, stored in the lattice during the HS to LS relaxation of the HS fractions of Fig. 1 (red curve) and that of a perfect lattice (blue curve). In the middle region, 4000-10,000 MCS, corresponding to the longitudinal propagation of the interface, the presence of the hole lowers the total strain in the lattice. Parameter values are the same as in Fig. 2. See text for more explanations. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article). Fig. 7. Maps of the local pressure of the lattice with a defect at several times along the relaxation process of Fig. 2a.

snapshots of spatial distribution of local pressures during the propagation process of the LS/HS interface of Fig. 2. We see that the elastic strain is delocalized over long distances from the HS-LS interfaces, which demonstrates the crucial role of the long-range interactions in this problem.

Study of the displacement field of the lattice with defect

To reveal the lattice strain in the relaxation process, we focus on the study of the displacement field of lattice sites A→H in the configurations of Fig. 2. We denote by uij ⎯ → (, ) the displacement field associated with the site coordinates ij (, ) and defined by:

uij rij r ij ⃗ ⃗ ⎯ → (, )= (, )-(, ) 0 (6)
Where, r ij ⃗ (, ) 0 and r ij ⃗ (, ) are the initial and final positions of the site vectors ij (, ) in the lattice. In this analysis, we take the positions of the HS state as reference positions. Therefore, we have,

r i j iR jR ⃗ (, )=( , )
.

HH HH 0 00 [START_REF] Doniach | Thermodynamic fluctuations in phospholipid bilayers[END_REF] To reveal the lattice strain, we calculate the derivatives of the various components of the displacement field. For example, the derivatives of x and y components of the displacement field along the x direction are written as follows in the approximation of continuous medium:

u x ui j uij R u x ui j uij R ∂ ∂ = ( +1
⎛ ⎝ ⎜ ⎞ ⎠ ⎟ ε u β u α = 1 2 ∂ ∂ + ∂ ∂ αβ α β (9)
where, α and β correspond to x and y respectively. The composition of the various elements of the strain tensor elucidates the mechanical effects induced by the growth areas. The divergence and rotational of the displacement field are interesting quantities to study and their expressions are given by: Figs. 8 and9, display the spatial distribution of the divergence and the rotational of the displacement field during the relaxation of the HS fraction, given in Fig. 2. In the picture (a) of Fig. 8, the relative volume expansion is zero in the HS state, because the latter is used as the reference state. One can see that div u r ⃗ [ ⎯ → () ]becomes negative when the LS domains develops on time and the most important changes are observed in the interface regions. In particular, this effect is enhanced around the hole (snapshot C and D), due to the significant contraction of lattice, when the HS/LS interface crosses the defect region. In contrast, the spatial distribution of rot u r ⃗ ⎯→ ⎯ [ ⎯ → () ], given in Fig. 9, reaches its maximum values at the most fragile regions of the lattice, which are mainly located at the intersection between the HS/LS interface and the edges of the lattice, as well as around the defect (snapshot C and D). This behavior clearly allows to postulate that this region may play the role of a nucleation center to initiate dislocations and/or cracks during the transition, as often observed by optical microscopy.

div u r ε r ε r rot u r ε r ε r ⃗ ⃗ ⃗ ⃗ ⃗ ⃗ [ ⎯ →

Conclusion

The effect of the structural defects on the process of propagation of HS and LS domains is still an open problem in spin-crossover solids; and was not yet considered from the theoretical point of view, although it is largely admitted that the elastic interactions are at the heart of the spin-crossover phenomenon. Knowing that elastic interactions generates naturally both short-and long-range interactions, the presence of a punctual defect (structural defect, topological defect, etc.) may significantly alter the propagation process of the HS fraction during the thermo-induced spin transition or during the relaxation of the photoinduced metastable HS state. At this end, we investigated in this work the spatiotemporal properties of a 2D lattice made of SCO sites interacting elastically. To consider the problem of the defect, we introduced in the center of the SCO lattice a hole with a fixed size and shape. Thus all atoms surrounding the hole have lost at least one neighbor for those situated in the edge and two neighbors for those in the corner. As a consequence, the presence of this particular defect introduces an additional surface to the lattice. The elastic and electronic properties of this lattice have been studied using the elastic spindistortion model, which is solved here by Monte Carlo technique. The model accounts for the spin state and volume changes accompanying the spin transition. A major issue of this work is the fact that the defect plays the role of the pinning site for the front interface. Indeed, results have shown that the interface propagation is significantly altered around the defect: the interface velocity slows down by one order of magnitude, which then can be used as a way to stabilize or to control [START_REF] Bousseksou | Ising-like model for the 2-step spin-crossover[END_REF][START_REF] Bedoui | Triggering a phase transition by a spatially localized laser pulse: role of the strain[END_REF][START_REF] Sy | Reversible control by light of the high-spin low-spin elastic interface inside the bistable region of a robust spin transition single crystal[END_REF] the dynamic HS/LS elastic interface. Furthermore, we observed that the interface shape changes well in advance to get the defect. Indeed, whereas for a perfect lattice the interface shape is straight and perpendicular to the lattice border, here a specific pointed shape is observed far from the hole, around which it is accentuated. This effect is really interesting because it can be used as sensor allowing to reveal the existence of micro-structural defects. Its understanding constitutes an interesting objective for a future work in which the position and the size of the defect on the shape of the HS/LS interface will be analyzed. The idea here is get some universal feature correlating the deformation of the HS/LS interface as function of its distance to the defect and the strength of elastic interactions between SCO sites.

Fig. 1 .

 1 Fig. 1. Adapted from Ref. [25]. The electro-elastic configuration of a SCO molecule (red ball in the center) in a 2D lattice, elastically-coupled to surrounding SCO molecules by four nearest-neighbors and four next-nearest-neighbors.

Fig. 2 .

 2 Fig. 2. Spatiotemporal configurations of the spin and lattice transformations during the HS to LS relaxation. Red points denote LS sites and blue dots correspond to HS sites, (a) fora lattice with defect and (b) for a perfect lattice. Remark in (a) the particular shape of the HS/LS interface, which is affected by the presence of the defect. From A to H, Monte Carlo times were: 0, 1944, 2986, 6709, 11529, 12949, 13580, 16451 MCS. For a perfect lattice (calculations made with the model parameters as for (b) b) the interface has a straight character and is perpendicular to the lattice borders. From A ' to F ' Monte Carlo times are: 0, 3785, 4476, 9255, 10901 and 13610 MCS. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Fig. 3 .

 3 Fig. 3. a) Elastic interface profiles along the propagation direction (x axis) for the atomic line, N =/ 2 y , plotted for various MC times during the relaxation process of Fig. 2. See text for more explanations. b) Spatial distribution of the nearest neighboring distances showing the spatial profile of the elastic interface around the defect region.

Fig. 4 . 2 y

 42 Fig. 4. time dependence of the HS:LS interface position along the propagation direction (x axis) for the reference atomic line y N =/ 2 y showing four propagation regimes. a) For the lattice with defect: the first regime (blue line) corresponds to the interface propagation far from the defect. Around the defect, the interface is fixed to the hole, and its velocity is extremely small. It comes off in the third regime. The last regime is driven by the presence of the lattice border. b) For a perfect lattice, showing a linear plot indicating that the interface speed is constant overall the process, except at the border. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

  in the same way. These terms provide access to the strain tensor ε αβ defined by:

  ] is the trace of the strain tensor. It describes the relative volume expansion, while rot u r ⃗ ⎯→ ⎯ [ ⎯ → () ] gives information about the shear stresses which causes the distortion of the lattice during the spin transition.

Fig. 8 .

 8 Fig.8. Maps of the divergence of the displacement field at corresponding to snapshots of Fig.2. The reference state is the HS phase, thus leading to a negative value of the divergence in the LS state. The significant changes are observed around the HS/LS interfaces.

Acknowledgments

The authors acknowledge by the Université Paris-Saclay, Université de Versailles, CNRS and the "Agence Nationale de la Recherche" (ANR project BISTA-MAT: ANR-12-BS07-0030-01) for their financial support. The authors are indebted to M. Sy for reading the manuscript and for the technical help. R.T. is indebted to the government of Algeria for the fellowship during his stay in GEMaC laboratory.


A Pron

L Brun

C Deruelle

O Coulon

DENSE AND STRUCTURED REPRESENTATIONS OF U-SHAPE FIBER CONNECTIVITY IN THE CENTRAL SULCUS

Keywords:

U-shape fibers are superficial white matter fibers connecting adjacent gyri. In this paper, we present a method to characterize the connectivity of U-shape fibers coursing around the central sulcus. Pre-and post-central gyral crests are semi-automatically drawn and used to build a connectivity space that is identical between subjects. This group space provides a dense representation of the shortrange connectivity between the two gyri, as well as a structured representation after clustering.

INTRODUCTION

Although main fiber bundles of the brain have been intensively studied during the last ten years, less is known about superficial white matter and more precisely about short associations fibers. Recent technical progress in diffusion MRI, the modality commonly used to probe white matter, have made feasible to access the structure of these fibers. Few studies have investigated such fibers and most of them concentrate on the central sulcus area, a sulcus with little variability across subjects and a well-known functional organization [START_REF] Coulon | Two new stable anatomical landmarks on the Central Sulcus: definition, automatic detection, and their relationship with primary motor functions of the hand[END_REF], [START_REF] Yousry | Localization of the motor hand area to a knob on the precentral gyrus. A new landmark[END_REF]. Such studies have tried to characterize the U-shape fibers connecting the pre-and post-central gyri around the central sulcus, either at the entire sulcus level, or locally based on predefined regions of interests (ROIs). In [START_REF] Catani | Short frontal lobe connections of the human brain[END_REF], tracks were selected from a whole-brain deterministic tractogram using nine ROIs in the pre-and post-central gyri. A similar approach was used in [START_REF] Magro | Characterization of short white matter fiber bundles in the central area from diffusion tensor MRI[END_REF] with six ROIs. The connectivity between ROIs was assessed through fiber count, normalized by the white matter volume within each ROI. In [START_REF] Zhang | Characterization of U-shape streamline fibers: Methods and applications[END_REF], associations fibers are clustered into four categories at the whole brain level. Connectivity profile of these fibers was assessed with respect to drawn sulcus fundi and gyral crest and the connectivity of sharp U-shape fibers was assessed globally between gyri as the number of fibers connecting pairs of gyri. Recently, a clustering technique dedicated to short association fibers was introduced in [START_REF] Román | Clustering of Whole-Brain White Matter Short Association Bundles Using HARDI Data[END_REF] and used to build a specific atlas from seventy subjects.

1 https://surfer.nmr.mgh.harvard.edu/fswki/ In this paper, we introduce a method to characterize the structure of U-shape fibers around the central sulcus in a local and continuous manner, providing a dense description of the U-shape fiber distribution, independent of any ROI-based seed. This distribution is characterized at the group level after extracting U-shape bundles and representing each tract as a point in a local connectivity space based on its extremities. This group-level connectivity space is defined from a parameterization of the pre-and post-central gyral crest lines.

This new representation not only provides a dense description of the short association fibers connecting both gyri but can also be clustered into bundles in order to provide a more structured description. The method was applied on the left central sulcus of ten subjects from the Human Connectome Project dataset.

MATERIAL AND METHODS

Dataset

Ten subjects from the publicly available Human Connectome Project dataset were used in this study. All individual were right-handed males aged of 25-30 years. Scans were acquired using a modified version of Siemens Skyra 3T scanner (Siemens, Erlangen, Germany) with a maximum gradient strength of 100mT/m, slew rate of 200T/m/s and a 32-channel head coil.

T1-weighted scans were acquired using 3D MPRAGE sequence (TR/TE=2400/2.14ms, flip angle=8°, FOV=224x224 mm², resolution=0.7 mm isotropic). Diffusion-weighted images were acquired using a spin-echo EPI sequence consisting of 3 shells (b=1000, 2000 and 3000 s/mm²) of 90 diffusion weighted volumes each and 6 interleaved non-weighted b0 volumes (TR/TE = 5520/89.5ms, resolution=1,25mm isotropic, FOV= 210x180mm², 111 axial slices, multiband factor=3, partial Fourier = 6/8, echo spacing=0.78ms). Each direction was acquired twice with reversed-phase encoding.

Preprocessing

T1-images were segmented using Freesurfer 1 and imported into Brainvisa 2 to obtain grey/white matter interface meshes. Diffusion-weighted scans were preprocessed using the Diffuse toolbox 3 , an in-house Brainvisa toolbox [START_REF] Pron | Diffusion MRI: assessment of the impact of acquisition and preprocessing methods using Brainvisa-Diffuse toolbox[END_REF] wrapping FSL 4 and Dipy [START_REF] Garyfallidis | Dipy, a library for the analysis of diffusion MRI data[END_REF] functions. Scans were corrected for subject motion and susceptibility artefacts. Corrected diffusion images and T1-images were non-linearly registered, then the white matter meshes could be brought into the diffusion space. Local Multi Shell Multi Tissue (MSMT) model parameters were estimated in Mrtrix 5 . MSMT impulsionnal response was first estimated for each subject using dhollander method [START_REF] Dhollander | Unsupervised 3-tissue response function estimation from single-shell or multi-shell diffusion MR data without a co-registered T1 image[END_REF] with default parameters. Whole brain deterministic tractography based on the MSTM model was performed using the Mrtrix SD_Stream algorithm [START_REF] Tournier | MRtrix: Diffusion tractography in crossing fiber regions[END_REF] with the following parameters: 50 million seeds, random placement in grey and white matter volume, Runge-Kutta 4 th order integration, 45° angle, constrained by anatomical tissues and cropping at grey-white matter interface [START_REF] Smith | Anatomically-constrained tractography: Improved diffusion MRI streamlines tractography through effective use of anatomical information[END_REF]. Resulting tractograms were filtered with the SIFT algorithm [START_REF] Smith | SIFT: Spherical-deconvolution informed filtering of tractograms[END_REF] .

Pre-and post-central gyrus crest lines delineation

The left pre-central and post-central gyral crest were drawn semi automatically onto each subject's left white matter mesh using the SurfPaint module [START_REF] Le Troter | An Interactive Sulcal Fundi Editor in BrainVisa[END_REF] in BrainVisa after manually choosing their end points (see Fig. 2). Each of the two lines was isometrically parametrized from the ventral to the dorsal extremity, such as to obtain a reference 1-D coordinate system [0-100], as described in [START_REF] Bodin | Anatomo-functional correspondence in the superior temporal sulcus[END_REF].

Hand-knob alignement reparameterization

In order to provide a better correspondence of the pre-and post-central gyri across subjects than the simple [0-100] parameterization, we reparameterized both crest lines to enforce the matching of a specific structure, the hand-knob, following the same procedure as in [START_REF] Bodin | Anatomo-functional correspondence in the superior temporal sulcus[END_REF]. The hand-knob is a part of the central sulcus with a characteristic inverted omega shape, known to locate the primary hand motor area reliably across subjects [START_REF] Yousry | Localization of the motor hand area to a knob on the precentral gyrus. A new landmark[END_REF].

The location of the most posterior extension of the handknob was located manually on the white matter mesh of each individual. The plane orthogonal to the inertia plane of the central sulcus and passing through the hand-knob extension was computed. The pre-central and post-central hand-knob landmarks were then chosen as the closest vertex to this plane on each of the two crest lines.

These two landmarks were then forced the same coordinate across subjects following the procedure in [START_REF] Coulon | Two new stable anatomical landmarks on the Central Sulcus: definition, automatic detection, and their relationship with primary motor functions of the hand[END_REF]. Specifically, each crest line was reparameterized using the following affine transformation: where 𝑦 𝑖 is the hand-knob coordinate for subject 𝑖 and 𝑦 𝑚 is the corresponding mean coordinate across subjects. After reparametrization, the hand-knob landmark has coordinate 𝑦 𝑚 across all subjects (Fig. 1). Both pre-central and postcentral crest lines are processed independently.

𝑖𝑓 0 ≤ 𝑦 ≤

𝑈 𝑠ℎ𝑎𝑝𝑒 𝑓𝑖𝑏𝑒𝑟𝑠 = { 𝑓 ∈ 𝐹 / (𝑑 𝑠 𝑝𝑟𝑒 < 𝑑 𝑠 𝑝𝑜𝑠𝑡 ) . (𝑑 𝑒 𝑝𝑟𝑒 > 𝑑 𝑒 𝑝𝑜𝑠𝑡 ) + (𝑑 𝑠 𝑝𝑜𝑠𝑡 ≤ 𝑑 𝑠 𝑝𝑟𝑒 ) . (𝑑 𝑒 𝑝𝑟𝑒 ≤ 𝑑 𝑒 𝑝𝑜𝑠𝑡 ) }

where 𝐹 stands for whole-brain tractogram, + and . stands for the OR and AND boolean operators. This first criterion insures that tracts are indeed going from one gyrus to another. In addition, tracks extremities (s and e) were associated with their respective closest mesh vertices. A geodesic distance map from the two gyral crest was then computed on the mesh. To insure only short U-shape streamlines connecting adjacent gyri are kept, a threshold of t =16 mm was set on the geodesic distance map. Streamlines whose closest vertices were further than t from the gyral crests were discarded, insuring that remaining streamlines go up in the gyri close enough to their crest line. This threshold was empirically found to include most fibers connecting the two gyri while removing spurious streamlines in the tractogram (Fig. 2.). Fig. 3. Group connectivity profile into the dense group connectivity space. Each marker aggregates the streamlines over unitary bins to ease visualization. Size of the marker is proportional to the number of streamlines. Group hand-knob position is plotted by green cross.

Group connectivity space

Each U-shape fiber selected during the previous step is associated with the pre-and post-central coordinates corresponding to its extremities (closest vertices on the gyral crest lines). Therefore each streamline corresponds to a point in the [0;100]x[0;100] space defined by the parametrization (section 2.2), and hereafter called 'group connectivity space'. Because the parameterization is common to all subjects, so is the group connectivity space, and all U-shape streamlines of all subjects can be plotted in this space, with a total of 361241 streamlines for the ten subjects.

Fiber clustering

The group connectivity space provides a multi-subject, continuous representation of U-shape fibers. However, in order to get a structured and simplified representation of this space we then cluster the group connectivity space using the k-medoids [START_REF] Kaufman | Finding Groups in Data[END_REF] algorithm with 1000 initializations and 300 iterations each, to ensure clusters stability. The number of clusters was empirically set to k=5. In this way, fibers that have similar starting and end points are clustered together in bundles defined by the cortical areas they connect.

RESULTS AND DISCUSSION

We presented a method to study connectivity profile of Ushape association fibers underneath the central sulcus. We introduced a continuous square domain common to all subjects by parameterization of the pre-and post-central gyral crests. In this continuous domain, a U-shape fiber connecting the two gyri is represented by a point. This group connectivity space is presented in Fig. 3. We observed a group connectivity variable density along the central sulcus from the ventral to the dorsal extremities. The denser and most numerous connections represented by larger disks are essentially located near to the hand-knob position, i.e. in the [60;80]x[60;80] domain.

A dense but smaller connectivity is also observed in the dorsal area i.e. in the [90;100]x[90;100] domain, while we observed sparse connectivity in the ventral area i.e.

[0;40]x[0;40] domain. These observations are similar to the functionally coherent (hand, foot and face area) tracts obtained and confirmed by dissection in [START_REF] Catani | Short frontal lobe connections of the human brain[END_REF]. In addition, we observe connections between the ventral pre-central gyri and the middle part of the post-central gyri. It may correspond to an extension of face superior U-tract described in [START_REF] Catani | Short frontal lobe connections of the human brain[END_REF] but also to previously undetected fibers or tractography artefacts. This connectivity organization is confirmed by the clustering shown in Fig. 4. Obtained clusters, three around the handknob area, one in the extreme ventral region and one in the extreme dorsal region define bundles that are similar to those described in [START_REF] Catani | Short frontal lobe connections of the human brain[END_REF].

By using a continuous parametric domain based on adjacent gyral crest lines, we provide a continuous description of structural connectivity along central sulcus without introducing priors by choosing region of interest location along sulcus. This method, applied on high diffusion data quality, suggests a discrete organization of the connectivity along the central sulcus which is coherent with the literature. However, we found some additional short-range associations connections that need to be further investigated.

Further work will concentrate on several points: the use of a micro-structure informed filtering method such as COMMIT [START_REF] Daducci | COMMIT: Convex Optimization Modeling for Microstructure Informed Tractography[END_REF] in order to reduce the number of biologically erroneous tracts and improve the quality of the tractogram; a study on a much larger group of subjects, such as the entire HCP900 database; and finally a groupwise clustering technique based on individual connectivity space representations in order to takes into account inter-subject variability. Finally, other sulci such as the temporal superior sulcus or the Calcarine sulcus.

Fig. 1 .

 1 Fig. 1. Schematic 2D representation of affine transform on a gyral crest line: subjects hand-knob coordinates 𝒚 𝒊 and 𝒚 𝒋 are transported to group mean hand-knob coordinate 𝒚 𝒎 .

Fig. 2 . 2 . 5 .

 225 Fig. 2. Dorso-ventral view of one subject: pre-(red), post-(blue) central gyri crests and resulting U-shape streamlines (green) behind the left white matter mesh. Hand-knob position is marked with the dotted line.

Fig. 4 .

 4 Fig. 4. Discrete and structured representation of group connectivity space achieved with k-medoids, k=5. Centroids are represented as orange dots, group handknob position is plotted by grey cross.

http://brainvisa.info

https://github.com/MecaLab/Brainvisa-Diffuse

https://fsl.fmrib.ox.ac.uk/fsl/fslwiki

www.mrtrix.org

4. ACKNOWLEDGEMENTS Data were provided by the Human Connectome Project, WU-Minn Consortium (Principal Investigators: David Van Essen and Kamil Ugurbil; 1U54MH091657) funded by the 16 NIH Institutes and Centers that support the NIH Blueprint for Neuroscience Research; and by the McDonnell Center for Systems Neuroscience at Washington University. A.P is supported by a grant from Aix-Marseille University.