Coupled spin-lattice dynamics from the tight-binding electronic structure

We developed a method which performs the coupled adiabatic spin and lattice dynamics based on the tight-binding electronic structure model, where the intrinsic magnetic field and ionic forces are calculated from the converged self-consistent electronic structure at every time step. By doing so, this method allows us to explore limits where the physics described by a parameterized spin-lattice Hamiltonian is no longer accurate. We demonstrate how the lattice dynamics is strongly influenced by the underlying magnetic configuration, where disorder is able to induce significant lattice distortions. The presented method requires significantly less computational resources than ab initio methods, such as time-dependent density functional theory (TD-DFT). Compared to parameterized Hamiltonian-based methods, it also describes more accurately the dynamics of the coupled spin and lattice degrees of freedom, which becomes important outside of the regime of small lattice and spin fluctuations.


I. INTRODUCTION
The interplay between spin and lattice degrees of freedom holds paramount significance in various fields of condensed matter physics, notably in ultrafast dynamics [1][2][3][4][5].At extremely low temperatures, ionic motion has been demonstrated to be negligible compared to spin motion [6][7][8][9].Hence, in this regime, the spin and lattice degrees of freedom are treated independently, which effectively describes a wide range of applications [10].However, a growing number of phenomena of interest to the scientific community necessitate consideration of magnon-phonon coupling when magnon and phonon frequencies are of similar magnitude: i) In multiferroics, the polarization induced by the combination of the charge and/or magnetic non-collinearity distorts the lattice, making the spin-lattice coupling the pivotal mechanism behind the magnetoelectric effect [4,[11][12][13]; ii) In spintronics applications and for terahertz applications, angular momentum currents can be converted between lattice and spin degrees of freedom via the spin-lattice coupling [14][15][16]; iii) Moreover, in femto-second ultrafast demagnetization processes it is debated that the angular momentum transfer is mediated by a coupling between lattice, spin, and electronic degrees of freedoms [17,18].
In the adiabatic approximation, both forms of dynamics are well understood [10,19].The magnetization dynamics usually takes the form of the Landau-Lifshitz-Gilbert (LLG) equation [20,21].In such simulations, one often introduces a coupling to a thermal reservoir, through the Langevin dynamics (LD) approach, which ensures that the distribution of the energy of atomic spins follows a Boltzman distribution [10,22].Regarding the lattice dynamics, it follows the principles of Newtonian dynamics [19].In this case, the time evolution of the ionic positions is driven by the atomic forces at every time step.Different techniques to minimize the atomic forces can be found in the literature [23], in order to lead the system to its respective relaxed structure.Similar techniques are also applied to minimize torques in the magnetic system [24].
In recent works, progress has been made towards coupling of these two dynamical processes.One way is to consider that the spin-spin exchange parameters, e.g. the isotropic (J ij ) and Dzyaloshinskii-Moriya ( ⃗ D ij ) interaction, depend on the atomic displacements, where i, j are the sites indices.This approach has been considered and applied in Ref. [25], where a Taylor expansion of an effective spin Hamiltonian was considered in order to collect these contributions, always assuming the adiabatic limit.Such expansion breaks the rotational invariance and conservation of angular momentum [26,27], which however may be restored by expressing the expansion using instantaneous differences between atomic positions instead of displacements [27,28].However, the use of displacements is a practical way to compute the parameters of various terms in the spin-lattice Hamiltonian, e.g. the force constants, and, therefore, this approach will be used in the present work.
Moreover, in Ref. [2] classical spin dynamics and arXiv:2311.00765v1[cond-mat.mtrl-sci] 1 Nov 2023 ab initio molecular dynamics were coupled in order to study CrN, where the authors were able to identify nonadiabatic effects, beyond an effective Hamiltonian description, from both the dynamics of the lattice as well as from the spin.When using such an effective Hamiltonian to describe the time evolution of the magnetic moments, it is assumed that their motion do not affect the parameters calculated from the initial state.Recent literature, however, has shown that this assumption only holds for small fluctuations around the ground state [29][30][31][32][33][34][35][36][37][38].It can be understood that in strongly non-collinear magnetic states, the electrons cause a spin current by hop between the atomic sites of the sample, giving rise to high-order magnetic interactions (or multi-spin parameters) [39][40][41].Another way to look at this phenomenon is that the interactions between two sites only depend on their current magnetic configuration [31,32,34].Therefore, it is expected that for systems where it becomes highly relevant to treat strong non-collinear magnetic configurations, the spin-Hamiltonian assumption might break or miss relevant features in magnon and phonon spectra compared to experiment [34,42].Therefore, as far as computational tools go to describe spin-lattice dynamics of a given system, Hamiltonian based models rely on the correct description of the system electronic structure throughout the parameters, which can be a delicate job for materials with non-collinear magnetism, while ab initio methods can only treat a limited number of atoms, due to the heavy computational cost.This paper presents a solution that involves computing the key components of the interconnected spin-lattice dynamics, including molecular forces and effective magnetic fields, directly from the electronic structure.The approach utilized here is based on a Slater-Koster parameterized tight-binding model using the National Research Lab (NRL) formalism [43][44][45], to calculate the ground state of a given system at every time step, from which we evaluate the molecular forces via the Hellmann-Feynman theorem [46][47][48][49] and the effective magnetic field as the field opposing the constraining magnetic field [50,51].Based on this, we performed coupled spin-lattice dynamics of Cr and Fe clusters (triangular and nanochains with 10 atoms) to find the ground state, when varying both but also only one degree of freedom.From this ground state, we explicitly calculate the parameter for an effective Hamiltonian and analyze how these parametersforce constants and exchange coupling parameters -vary in different situations and magnetic configurations.Our results show that for small variations near the magnetic ground state, these parameters are fairly constant.Far away from the ground state, particularly for the Cr triangle case where the magnetic ground state is non-collinear, there is a significant difference.As to the magnetic properties, Cr tends to stabilize in a AFM fashion, which results in the Néel AFM magnetic structure for the triangular trimer and a collinear AFM for the nanochains.The Fe based systems present a FM ordering for the triangular trimer and a spin-spiral like configuration for the nanochains, resulted from the interplay between short and long range interactions.During the relaxation process of the atomic positions, we show how the magnetic moment length varies and how the lattice dynamics is affected by the motion of the magnetic moments.In fact, for the Cr trimer, our results reveal that a distorted magnetic structure can lead to a distorted lattice structure.Lastly, we compared our tight-binding spin-lattice dynamics with semi-classical spin-lattice dynamics from the effective Hamiltonian, which revealed missing contributions to the effective Hamiltonian in the highly disordered case.
This paper is structured as the following: in Sec.II, we discuss about the adiabatic spin-lattice dynamics and introduce the equations regarding the time evolution of the magnetic moments and ionic positions cause by constraining fields and forces, directly obtained from the underlying electronic structure.In Sec.III, we introduce the semi-classical effective spin-lattice Hamiltonian which is used to compare the tight-binding with the semi-classical spin-lattice dynamics and for which we calculated respective parameter.These parameters are discussed and analyzed for the cluster of Fe and Cr in Sec.IV.Finally, in Sec.V, the spin-lattice dynamics of these nanoclusters obtained from tight binding and from an effective Hamiltonian are presented.We summarize our findings in Sec.VI.

II. ADIABATIC SPIN-LATTICE DYNAMICS
Following our assumption in Ref. [52], we are interested to describe physics of magnetic moments e and lattice displacements u on a timescale above 1 fs.Here, adiabatic dynamics operates under the premise that electrons evolve significantly faster than the dynamics of the respective degrees of freedom [9,53].This premise is solidly grounded, particularly for quasi-particle excitations possessing energies considerably below intrinsic electron energies like Stoner spin splitting and variations in electron hopping.On comparable energies, the dynamics of the electron need to be addressed from first principle methods, e.g., time-dependent density functional theory (TD-DFT) [8,54].
In the adiabatic approximation, the total energy of the system is a function of only the magnetic moment directions {e} and lattice displacements {u}, where the electronic degrees of freedom can be regarded as being in a quasi-equilibrium state, E = E({e}, {u}). (2.1) For precise computation of electronic states and, consequently, the total energy within first-principles methods, it is imperative to rigorously constrain the calculated magnetic moments with specified 'input' moment directions denoted as e.Without this constraint, the system would naturally experience a finite torque and revert to the absolute ground state due to relaxation processes [50,55,56].Adding an extra term representing the constraining field B con i at a lattice site i to the electron Hamiltonian [51,52], we conserve calculated magnetic moments in the direction e.Here, B con i is obtained according to B con i • e i = 0 in an iterative, self-consistent algorithm [50,55].This criterion for the constraining field shows that its acts perpendicularly to the magnetic moment direction e i and, thus, do not change the magnetic moment length.
In the constrained electronic solution, the magnetic moment experiences an intrinsic effective field B eff i = −B con i [52].This field drives the time evolution of the magnetic moment m i = m i e i following the Landau Lifshitz Gilbert equation [10,53]: where γ is the gyromagnetic ratio, α is the phenomenological Gilbert damping.We keep the Gilbert damping finite due to dissipation effects by spin disorder or electron correlation effects [57].For instance, the origin of the Gilbert damping has been discussed to originate from spin-orbit coupling [58].
Regarding the dynamics of the lattice degree of freedom, the Born-Oppenheimer approximation can be considered, where the electronic and ionic degrees of freedoms can be treated separately.Here, forces F k acting on the atom k are well defined by the Hellmann-Feynman forces [46][47][48][49] and the adiabatic lattice dynamics simply consists in solving numerically Newton's equation of motion to produce the trajectory u k of the atom with mass M k and momentum p k .Also here, a friction force proportional to the momentum and scaled by the parameter ν [23] is added, analogous to the damping constant α in the case of spin-dynamics, where ν is the friction force.It should be noted that we ignore the influence of temperature in both equations of motion, that e.g. in Langevin dynamics is introduced via a stochastic field.
In our work we represent the electronic structure by a non-orthogonal Slater-Koster parametrized, realspace tight-binding model implemented in the software package Cahmd [59].The Slater-Koster parameters are expanded in terms of the hopping distance of the electrons between the atoms according to the NRL tight-binding method.Magnetism is included by a Stoner term proportional to the magnetic moments e.A local charge neutrality term conserves the total charge of the system.Since the orbital overlap matrix enters via Mulliken transformation [60] also into the magnetic as well as charge neutrality terms of the Hamilonian, not only the hopping but also the exchange splitting and charge distribution in the system become dependent on the lattice displacement.More details about the model are given in the supplemental material A and in Ref. [44,45].The electronic structure is self-consistently iterated with respect to the charge, the magnetic moment length, and the constraining field up to the accuracy of 1 • 10 −10 , i.e. the relative error of each convergence parameter (charge, magnetic moment length and constraining field) is of that order.From the solution we extract the effective field [51] and calculate from the Hellman-Feynman theorem the lattice forces [47,61].More details in how the effective magnetic field and ionic forces, used in Eqs.(2.2) and (2.4), are calculated can be found in Ref. [51] and in Appendix.C, respectively.It should be noted that we neglect terms in the effective fields and forces that are proportional to the variation of the self-consistent parameters, say changes of the charges with the respective degree of freedom ∂n /∂⃗ e, ∂n /∂⃗ u as well as changes of the magnetic moment length with the respective degree of freedom ∂m /∂⃗ e, ∂m /∂⃗ u.These gradients are typically small and negligible.
Having the effective fields and forces, we integrate Eq. (2.2) and Eq.(2.4) using the implicit mid-point method as described in Ref. [25] using parameters α = 0.1 and ν = 0.1 fs −1 .Note that a common way to solve Eqs.(2.3-2.4) is to use the so-called Verlet-type based algorithms [23].However, these methods, for the purpose of spin-lattice dynamics, produce numerical instabilities which lead to the non-conservation of the total energy of the system.Our numerical integration step is ∆t = 1 as

III. THE SPIN-LATTICE HAMILTONIAN
A well defined procedure is to project the energy in Eq. (2.1) to a parametrized Hamiltonian.Following Ref. [25], such a Hamiltonian can be defined as where J αβ ij is the exchange tensor, Γ αβµ ijk is the derivative of the exchange tensor with respect to the lattice dis- , Φ µν kl is the force constant and the last term is the kinetic energy.In our study, no spin-orbit coupling is considered such as that the exchange tensor J αβ ij has only diagonal, identical terms.It is important to highlight that in the exchange striction constants Γ ijk for nanostructure, the term i = j is finite and is required to fulfill Newton's third law for the exchange striction term in the ferromagnetic collinear state: ij Γ µ ijk = 0 for all k.From Eq. (3.1) and for performing semi-classical spin-lattice dynamics, one can calculate the effective magnetic field B eff i = − 1 /mi ∂H SLD /∂ei and the ionic force F k = − ∂H SLD /∂u k by taking the derivative of the spin-lattice Hamiltonian with respect to the magnetic moment of site i and displacement of site k, respectively [2,10,30,52].
In order to compare the dynamics with effective fields and forces from tight binding electronic structure and from semi-classical parametrization, we calculate the parameters in Eq. (3.1) also from the tight binding method.In Ref. [33], we solved initial challenges to calculate the spin exchange for an arbitrary non-collinear magnetic configuration.This study is based on the original work by Lichtenstein et al [62,63].The exchange striction term is calculated via the numerical procedure described in Ref. [25].On the other hand, the exchange striction parameter can also be calculated explicitly via Green's function methods as demonstrated in Ref. [64].The numerical procedure, though, allows to calculate also the exchange striction for any arbitrary magnetic state due to corrections coming from the constraint methodology.
The force constants are obtained from the numerical derivation of the forces where {µ, ν} are the cartesian coordinates and {k, l} are the atomic sites.The Φ's are symmetrized according to point group symmetrical relations.For instance, as a consequence of the Newton's third law, the force constants have the following sum rule Also here, the force constants can be explicitly calculated from the Green's function of the underlying electronic structure as proposed by Fransson et al [65].In Eq.(3.1), we explicitly neglect the expansion of the force constants with respect to the magnetic moment configuration, following Ref.[66].These are terms of the form ) To nevertheless prove the existence of such high-order terms, we calculate the changes of the force constants when rotating the magnetic moment.Due to the explicit form of the constrained isotropic spin-spin exchange [33], we are able to also time-resolve J αβ ij .But the numerical treatment of the force constants as well as exchange striction terms, do not allow for a time resolved output during the spin-lattice tight binding dynamics.Calculating the parameters of the spin-lattice Hamiltonian (3.1) in the ground state of the system, we are also able to perform semi-classical dynamics based on the H SLD .It should be noticed that calculations of the parametrization as well as the semi-classical dynamics are also done via the package Cahmd [59].Thus, the entire methodology presented is on equal footing.

IV. PARAMETRIZATION OF THE SPIN-LATTICE HAMILTONIAN
In the following we will first discuss the parameterization of the spin-lattice Hamiltonian (3.1) for trimers and nanochains of Fe and Cr.These parameters are required to calculate the effective fields and forces from the Hamiltonian (3.1) and to perform semi-classical dynamics.Besides analyzing the parameters, we will also discuss the magnetic moment dependent force constants just for the trimer case, to confirm whether the terms (3.4) and (3.5) can be neglected.

A. Triangular trimer
Here, we calculate the spin-lattice Hamiltonian parametrization using the magnetic and lattice ground state as reference state.The ground state was obtained by performing the tight binding spin-lattice dynamics starting from an arbitrary state (see next Section).Due to the large energy dissipation, the system converges rapidly (≈ 500 fs) to its zero-temperature ground state.In general, the dissipation time of each channel (lattice or spin) will depend on how far each channel is from its respective ground state.For instance, for the Cr trimer case, since the ground state is non-collinear, the magnetic moments takes longer to reach equilibrium compared to the lattice counterpart.Here the magnetic ground state found is ferromagnetic (FM) and Néel antiferromagnetic (N-AFM) state for Fe and Cr, respectively.In the case of Fe, the isotropic spin-spin exchange is J 12 = 362.8meV, force constants are Φ xx 12 = −5.92eV/a.u. 2 ,Φ xy 12 = −0.40eV/a.u. 2 , and Φ yy 12 = 0.23 eV/a.u. 2 in an irreducible representation.For the equilateral triangle sitting in the xy-plane, the elements Φ xz , Φ zx and Φ yz , Φ zy are zero by symmetry.The exchange striction terms are calculated (data not shown, but used in Section IV).The magnitude of Γ is comparable to the values of the spin-spin exchange or larger as also found in Ref. [25].Other symmetries mentioned in Ref. [25] are also confirmed.
For Cr, we have J 12 = −143.9meV, force constants are Φ xx 12 = −0.81eV/a.u. 2 ,Φ xy 12 = 0.13 eV/a.u. 2 , and Φ yy 12 = −0.07eV/a.u. 2 .It is important to point out that representing these interaction in an irreducible form is only possible due to symmetries available in the ground state.It can be shown in this case that the pair interactions can be transformed into another.Breaking the symmetries in the magnetic and lattice ground state breaks also the symmetry relations of the pair couplings.The same holds also for the exchange striction terms, which are shown for two high-symmetry cases in Fig. 2. Important to notice are the points where i = j, which are finite for these nanostructures.We checked also nanostructures with periodic boundary conditions and obtained zero for parameters where i = j (data not shown here).The magnetic configuration cases are the FM state and the N-AFM state.Remarkably, both magnitude and relation between different sites are significantly different.As already shown in Ref. [52], magnetic states far away from the actual ground state (here the Néel state) show stronger fluctuations when the internal degrees of freedom are varied.This causes larger values of Γ in the ferromagnetic Cr trimer.Furthermore, the distance between the atoms turned out to be different for the FM state compared to the Néel state, after the lattice relaxation process, which in turn impact electron orbital overlaps that impacts the strength of the spin-spin coupling and variations of the spin-spin coupling.These strong changes of Γ have an impact also on the dynamics, however, in semi-classical dynamics often only the values of the ground state are used.This aspect will be discussed further in more details in the next section where these parameters are used to perform spin-lattice dynamics and compared to our tight binding spin-lattice dynamics implementation.
It has been heavily discussed in the literature how magnetic interactions are dependent on the underlying magnetic configuration of the system [31,33,34,67].The same question arises with respect to the force constants: do the force constants depend on the magnetic state?For instance, Ref. [66] concluded that this dependence is neglegible.In order to investigate such dependence, we calculated the onsite force constants as a function of the magnetic configuration.The magnetic moment of a single atom, sitting along the y-axis, is rotated around the z-axis while the onsite force constant matrix is calculated (see Fig. 1).Given the nature of the onsite term, Eq. (3.3), the results can be interpreted as how the pairwise force constants of the other atoms react to the mag- netic moment rotation of the atom considered.The rotation is done taking the ground state as reference, i.e.FM and N-AFM states for Fe and Cr, respectively.That magnetic dependence can be seen as the emergence of high-order terms in the magnetic moment, as the ones described in Eqs.(3.4)-(3.5).In the cases studied here, the only components that are not zero are Φ xx , Φ yy and Φ zz .Particularly for Cr, given the noncollinear magnetic ground state, the Φ xy = Φ yx is also not zero and present a linear dependence as opposed to the quadratic one seen in the other terms (Fig. 3).It is worth mentioning that components Φ zz and Φ xy = Φ yx (for the Cr case) are zero if the system is in a magnetic ground state.The changes are approximately 2% for the Φ xx and 1% for the Φ yy components for Fe.For Cr, the values are 8% and 5% for Φ xx and Φ yy , respectively.These percentages were obtained comparing the Φ's between θ = 0 rad and θ = 0.87 rad.The element dependence of these interactions is unclear since it depends on a variety of factors.Nevertheless, our examples show clearly that terms like (3.4) and (3.5) should be considered in semi-classical dynamics.

B. Nanochain
Here, we considered a nanochain with 10 atoms along the x-axis, as shown in Fig. 1.Although being a simple one-dimensional example, there are a set of experimental data of nanochains supported on a class of metallic surfaces whose ground state is noncollinear [68][69][70].
On Fig. 4 we show the calculated exchange striction Γ µ ijk and the isotropic spin-spin exchange parameter J ij for the nanochain.In the nanochain case, only the xcomponent of the spin-lattice parameter, Γ x ijk , is not zero by symmetry.On the left side of the figure, we keep i and j fixed and vary k, which is the atom being displaced, from the i−th site (Γ x 565 ) itself until the end of the chain, k = 10.For that case, Fe presents a bigger sensitivity to lattice distortions far away from the local pair {i, j}, as oppose to the Cr nanochain case which the interaction shows itself more localized.On the right side, we keep i and k fixed to 5 and vary j until 10.For such case, the spin-lattice interaction for both Fe and Cr systems are localized with only relevant terms being the next-nearest and next next-nearest neighbors.For the isotropic spin-spin interaction, the relevant terms are also for the next-nearest and the next next-nearest neighbors.Regarding the force-constants Φ µν lk , only the components xx, yy and zz are not zero in the nanochain case, with yy and zz having the same magnitude.The values are Φ xx 55 = 4.72 eV/a.u. 2 for the Fe and Φ xx 55 = 3.39 eV/a.u. 2 for Cr, while the other components are Φ yy 55 = Φ zz 55 = −0.12eV/a.u. 2 and −0.34 eV/a.u. 2 for Fe and Cr, respectively.It should be noted that the force constants are calculated for the respective systems magnetic ground state.The magnetic dependence of the force constants having the ground state as reference is also verified for this case (data not shown).
The above analysis of the interactions shows clearly that a precise description of the parametrized total energy (2.1) by a defined Hamiltonian, e.g.Eq. (3.1), is complex and can vary from material to material.Based solely on the values of J's, ϕ's, and Γ's, it is challenging to gauge the extent of the real-world impact that, e.g. the orientation of magnetic moments has on lattice dynamics or the lattice displacements have on the magnetic ground state.To delve deeper into assessing the consequences of these alterations, we conducted semi-classical spin-lattice dynamics analysis in Sec.V B. This involved employing a classical spin-lattice Hamiltonian with fixed parameters and contrasting these outcomes with those obtained from a tight-binding spin-lattice approach.In the following section, we present the results of this comparative analysis.

V. SPIN-LATTICE DYNAMICS OF NANOCLUSTERS
We conducted a comparative study by performing spin-lattice dynamics using two different approaches: the tight-binding (TB) model described above, and a spinlattice Hamiltonian.In the former, we directly computed forces and the effective magnetic field from the electronic structure.In the latter, we computed the parameters of the spin-lattice Hamiltonian, assuming they remained constant during the dynamical process.
As discussed above, the semi-classical dynamics requests the energy parametrization from the ground state and, thus, misses effects that arbitrary configuration {e} and {u} are doing to the parameters, when not considering even larger expansions of the total energy in a parametrized Hamiltonian.The tight binding based dynamics, on the other hand, is able to calculate precisely fields and forces for these arbitrary configuration.Since the magnetic moment lengths are self-consistently determined at every time-step, the TB spin-lattice dynamics includes also longitudinal moment fluctuations, which are not addressed at all in typical semi-classical  ijk (exchange striction) and isotropic spin-spin exchange parameter Jij for the nanochains of Fe (top) and Cr (bottom).On the left panel, the atom 5 is taken as the reference atom while varying the k-th atom.On the right panel, the atom 5 is also taken as reference while the j-th atom is the one varying.Still on the right panel, for the spin-lattice parameter, the k-th atom is also being fixed.
dynamics [10].These fluctuations, however, will mainly influence the precession frequency of the moments, rather then the actual ground state.

A. Dynamics via tight-binding model
When focus extends beyond dynamics near an equilibrium state, as is often the case, utilizing tight-binding spin-lattice dynamics becomes important, at least for the systems considered here.This is especially pronounced when commencing from a randomly chosen initial state, a scenario frequently encountered when seeking to ascertain the precise ground state of a nanostructure.The TB spin-lattice dynamics calculations began with initially disordered magnetic configurations for both the triangular trimer and nanochain systems.Specifically, for the trimer, we investigated spin-lattice dynamics in various scenarios.These scenarios included comparing dynamics with and without fixing the magnetic moment orientation, which was initiated from a disordered magnetic configuration.
When the magnetic moments are fixed (lattice dynamics only performed), the results show that the magnetic disorder induces a distortion in the systems as the atomic sites relax to their lowest energy positions, as shown in Fig. 5, lower panel.This is due to the force induced by the non-collinearity.Such forces are, e.g., due to the exchange striction term in Eq.(3.1).However, forces calculated from the tight binding method do not fully agree to them (data not shown here) since they contain also higher order terms, such like forces coming from terms of the form (3.4) and (3.5).This distortion effect is more dominant for Cr than for Fe (data not shown here) with a change in the distance up to 1 a.u. and has also impact on the magnetic exchange coupling, making it spatially anisotropic.Now, if the magnetic moments are allowed to relax together with the atomic sites, i.e. the complete spin-lattice dynamics is performed, the systems (triangular trimers) find their true magnetic and molecular ground states, which in all cases correspond to the atomic sites sitting in the vertices of an equilateral triangle and the magnetic moments having a ferromagnetic ordering for Fe, while having an antiferromagnetic Néel state for the case of Cr (see the spin configurations in Fig. 6).The trajectory of the spins towards this ground state is shown in Fig. 5, upper and center panel.Although we used the same energy dissipation parameter for both materials, the relaxation time of the magnetic moments is much larger for Cr compared to Fe and in the order of hundreds of femtoseconds.Similar isotropic spin-spin exchange constants J ij argue also for comparable timescales.It can be debated that during the relaxation process, the Cr magnetic moments cover a wider range of the configurational space, making higher order coupling mechanism to most likely play a more crucial role during the process for the Cr case as opposed to the Fe case.However, these couplings are hard to resolve for the constrained dynamics.
Additional insight into the distinct timescales observed can be obtained by examining both the magnetic moment lengths and lattice displacements, as depicted in Fig. 6.For Cr, the magnetic moment length exhibits more significant variations, reaching up to 0.4µ B , in contrast to Fe, where the variations are limited to 0.04 µ B .Moreover, Cr inherently possesses larger magnetic moments, which in turn results in a lower precession magnetic field and, consequently, longer relaxation times.Analyzing lattice displacements reveals another disparity between Fe and Cr.In the case of Cr, the relaxation times of the lattice closely match those of the spins.This suggests a stronger spin-lattice coupling in Cr compared to Fe.Furthermore, around the 100 fs mark in the case of Cr, there's a notable divergence in the amplitudes of atomic distances, despite an initially symmetric atomic position.This discrepancy hints at an energy transfer between the magnetic moments and the lattice in Cr, likely due to high-order terms naturally considered in the TB approach and not taken into account via spin-lattice Hamiltonian, a phenomenon much weaker in the case of Fe.
The nanochain case is a slightly more complicated example compared to the trimer.This is due to the fact that now there are long range interactions that can compete with each other and may lead to a set of complex magnetic configurations, in case of spin-dynamics, Cr trimer with fixed moments Figure 5.The angle between the magnetic moments of the Fe and Cr triangular trimer as a function of time in femtoseconds, as well as their respective exchange coupling parameters calculated at every time step.In the bottom panel, we show the results when the magnetic moments are considered fixed (lattice dynamics only), for the Cr triangular trimer case.The disordered magnetic moments induce a distortion in the lattice.
and complex structural composition, in case of lattice-We started the nanochains with a random magnetic configuration and the atoms sitting along the x-axis with a spacing of 4.6 a.u between them, as can be seen in Fig. 7.In the case of the nanochain with 10 atoms, the first phenomenon which can be observed is the tendency of nearest neighbors atoms to occupy positions closer to each other.This is the so-called Peierls transition or Peierls distortion [71].Because electrons are free to move, it creates charge density waves which induce distortion to atomic positions.This has been observed before in weakly coupled molecular chains [72].As these atoms get closer to each other, their magnetic moments start to behave similarly, i.e. they act as a magnetic sub-lattice.For instance, for both Fe and Cr, this distortion leads to a dimerization of the atomic sites and their magnetic moment directions are then collinear between themselves, assuming a ferromagnetic alignment for the Fe and an antiferromagnetic alignment for the Cr case.Moreover, in the Fe case, these dimerized magnetic moments form a spin-spiral (with a long wave-length, as seen in Fig. 7) as the magnetic ground state.The phenomenon of dimerization was also discussed in Ref. [65], were the authors highlight the presence of a dominant bilinear spin-lattice coupling term as the primary mechanism linking the spin and lattice behavior.Notably, this coupling term, which plays a crucial role in understanding the interplay between spin and lattice dynamics, is not considered in the semi-classical Hamiltonian (3.1).The omission of this term is due to the global frame defined in Ref. [33], where its inclusion would disrupt the time reversal symmetry of the system.Consequently, the semi-classical description provided by (3.1) would fail in capturing the dimerization effect.This underscores the distinct advantage of the here proposed tight-binding method, which accommodates the bilinear and other possible spin-lattice coupling term, enabling a more comprehensive analysis of the dimerization phenomenon.An infinitely long one-dimensional chain must transform under a Peierls distortion, but it is unclear if one-dimensional chains of limited size must distort as well.Most likely this is system dependent, but we note that the results presented here are consistent with what one would expect for infinite sized chains.Since the nanochain system has an even number of atoms, we can compare the left half with the right half.In Fig. 7, we compare the different sides with the same color, but different linestyles.While the distance between pairs of atoms and their respective magnetization is shown with a full line for the left side, the right side is represented with a dashed line.Concerning the time evolution of the distances between different pairs of atoms, it is expected that both full and dashed lines, representing equivalent pairs of atoms from different sides, to behave equally.Instead, one can verify that full and dashed lines can be distinguished, specially in the beginning of the simulation.This is due to the fact that although both sides have equivalent spacing between the atoms, their magnetic moments are not equivalent.It suggests a transfer of angular momentum between the magnetic moments and the lattice which leads to a slightly different lattice dynamics on each side.That effect is, as it was for the trimer system, stronger on Cr than on Fe.Cr trimer Figure 6.Spin-lattice dynamics for the triangular trimers of Fe and Cr.For these calculations, both the spin and lattice damping parameters are set to 0.05.On the left, we show the Cr triangular trimer relaxation process for the magnetic moment length (top) and the atomic positions (down).The figures in the middle denote the start configuration (left) and the final configuration (right).Note that the relaxation processes have different time scales, where the atomic positions relaxes at a faster rate than the magnetic moment lengths.The Fe triangular trimer is shown on the left, where magnetic moment lengths and atomic positions relaxes at a similar rate.

B. Dynamics via spin-lattice Hamiltonian
It has been widely discussed in the literature the concept between local and global Hamiltonian [33,67].While the former relies on describing locally the energy of a given system, the latter uses a Hamiltonian with a complete set of parameters capable of describing the total energy of the system at any point in the phase space.Although convenient to use the global approach, the Hamiltonian cannot be known a priori.A way of avoiding this issue is to extract both ionic forces and effective magnetic field directly from the electronic structure as shown in the previous Sections.
In this Section, instead of calculating the ionic forces and effective magnetic field directly from the electronic structure, we use the parameters described in Sec.IV as an input for a semi-classical spin-lattice Hamiltonian.From this Hamiltonian, the forces and effective field are calculated and used in the coupled spin-lattice dynamics.The goal is to determine validity of the semi-classical approach for a magnetically disordered system when comparing it to the tight binding spin-lattice dynamics.We started the simulation from the magnetic and structural ground state obtained from above described tight binding spin-lattice dynamics and by inducing an initial magnetic disorder while atoms sit in their relaxed atomic positions.If the system is in a magnetic ordering, since the ionic forces coming from the the exchange striction will be zero, the atoms will not move during the spin-lattice Cr nanowire Figure 7. Spin-lattice dynamics for the nanochains with 10 atoms of Fe and Cr.For these calculations, both the spin and lattice damping parameters are set to 0.5.The figures in the middle denote the start configuration (left) and the final configuration (right).The full lines concern the atoms in the left half of the chain, while the dashed lines stand for the right half.One can note that during the first femtoseconds of simulation, they can be easily distinguished due to the disordered magnetic configuration.
dynamics.However, the initial magnetic disorder gives rise to a net force and the atoms evolve.That can be seen in Fig. 8.That movement can be understood as the angular momentum of the magnetic moments is transfer to the lattice, causing the atoms to move initially and going back to their ground state after some relaxation time.Surprisingly, we also observed such large discrepancies in the relaxation time between the Fe and Cr trimer as in the above discussed tight-binding dynamics.Since the semiclassical Hamiltonian includes the parametrization from the magnetic and atomic ground state, both the magnetization and lattice displacements relax toward this state.
When comparing the tight-binding (dashed lines) with the semi-classical (solid lines) in Fig. 8, difference between both can be seen already in the early femtosec-onds of the simulation and mostly in the displacements of the atoms.This might be due to the emergence of highorder terms of the form (3.4) and (3.5) which we discussed in Sec.IV but not explicitly included in the spin-lattice Hamiltonian (3.1).Surprisingly, the magnetic moments in the Fe trimer do not show much difference between the different methods.The discrepancy for the Cr trimer case is much more apparent.This is both in the relaxation time of the respective degree of freedom as well as in the magnitude of the perturbation.The lattice displacements are more significant in the tight binding case compared to the semi-classical model.This behavior can be comprehended by examining the exchange striction constants Γ, as demonstrated in Fig. 2. It is evident that Γ exhibits an upward trend when transitioning from the Néel state to a ferromagnetic state.This increase in Γ . Comparison between the semi-classical (full line) and the tight-binding (dashed line) spin-lattice dynamics.The damping for both lattice and spin were set to 0.1.The parameters to the semi-classical spin-lattice Hamiltonian were calculated from the systems respective ground state, which are FM to the Fe trimer and Néel AFM for the Cr.Each system starts with disordered magnetic moments and in their respective ionic ground state.
facilitates from the movement of spins that the atoms more readily displaced in comparison to the Néel state.
An increase in the parameters also causes a growth in the precession frequency and, thus, relaxation times become smaller compare to the semi-classical case where the Γ's of the Néel state are used.Furthermore, besides the pronounced variations in exchange striction, it is worth noting that other higher-order terms may also contribute to this disparity observed between tight-binding and semiclassical dynamics.

VI. SUMMARY AND DISCUSSIONS
We have developed a fully integrated spin-lattice dynamics approach rooted directly in the tight-binding electronic structure.Unlike traditional methods that relies on a parametrized classical spin-lattice Hamiltonian, ours extracts molecular forces and the effective magnetic fields straight from the tight-binding structure via the Hellman-Feynman theorem and a self-consistent constraining field, respectively.
Our investigation into semi-classical Hamiltonians for Fe and Cr trimers and nanochains revealed significant variations in the Hamiltonian parameters based on the magnetic reference state.This aligns with literature findings about the magnetic texture influence on the magnetic interactions [29,31,33,36,38].More specifically, our tight-binding spin-lattice dynamics suggests the profound influence of magnetic texture on the Cr trimer case, whose ground state is non-collinear.In this case, parameters calculated from the ferromagnetic configuration as reference cannot correctly describe the system's dynamics.Still, the parameters calculated from the ground state describe the spin-lattice dynamics correctly only locally.i.e. for small variations of the magnetic moment orientation around the ground state.We argue that this scenario can be strongly present in sytsems whose magnetic ground state is non-collinear, such as 2D Kagomé magnets [73].Moreover, small differences can also be seen for the other systems.
The advantages of our present method is evident.It treats spin and lattice dynamics simultaneously, directly from the electronic structure, differing from previous methods like those in Ref. [2].The self-consistent electronic structure calculation at each step enables consideration of magnetic moment length relaxation, although certain abrupt changes can only be approximated at this stage.
Looking ahead, we are focusing on incorporating spinorbit coupling effects, anticipating even more pronounced variances between tight-binding and semi-classical dynamics, especially with non-collinear effects at finite temperatures.In essence, our new approach offers a complete tool for probing spin-lattice dynamics from first principles across diverse materials.According to the Eqs.(C3) and (C7), we can finally write where ℓ, ℓ ′ are the orbital indexes and {k, j} are the site indexes.Given the fact that both the Hamiltonian H kℓ,jℓ ′ and the overlap matrix S kℓ,jℓ ′ are a function of the Slater-Koster table in the NRL tight-binding approach, the derivative with respect to the atomic positions ⃗ R k have an analytical expression which is essentially the derivative of the Slater-Koster table with respect to the direction cosines.This was done in Ref. [61].More details on the tight-binding model can be seen fully in Refs.[30,45,52,60].

Figure 1 .
Figure 1.Schematic representation of (a) the triangle trimer and (b) the nanowire with 10 atoms, with their respective coordinate system represented by the tripod.

Figure 2 .
Figure2.The spin-lattice parameter Γ µ ijk for the triangular trimer of Fe (top panel) and Cr (middle and bottom panel).For the Cr, we present the spin-lattice parameter calculated from two different magnetic configuration: ferromagnetic (middle panel) and Néel state (bottom panel).

2 )Figure 3 .
Figure 3. Onsite force constant matrix for the rotated atom in the Fe and Cr triangular trimer, top and bottom, respectively.

Figure 4 .
Figure 4. Spin-lattice parameters Γ µijk (exchange striction) and isotropic spin-spin exchange parameter Jij for the nanochains of Fe (top) and Cr (bottom).On the left panel, the atom 5 is taken as the reference atom while varying the k-th atom.On the right panel, the atom 5 is also taken as reference while the j-th atom is the one varying.Still on the right panel, for the spin-lattice parameter, the k-th atom is also being fixed.