Impact of electron correlation on the light-induced demagnetization of elemental ferromagnetic metals

The local spin-density approximation (LSDA) is known to describe poorly the electronic structure of 3d transition metals, yet most density-functional-based ab-initio studies of ultra-fast demagnetization rely on it. One way to account for Coulomb correlations among the localized d electrons and go beyond LSDA is to include the effective correlation energy (or Hubbard) U. By doing so, we show here that electronic correlations lead to sizable changes of the laser-induced demagnetization of iron, cobalt, and nickel. We study how the various laser parameters, such as pulse duration or intensity, change the magnetization dynamics. It turns out that the total laser fluence is not suitable to quantify how much a laser pulse demagnetizes a material, as changes in pulse duration and shape influence significantly the outcome. The findings are traced back to the electronic structure of the material, and explained based on phase space for optical transitions.


Introduction
In 1996, Beaurepaire et al. [1] showed that a 60 fs optical laser demagnetizes a film of nickel on a sub-picosecond timescale triggering a vast number of experimental and theoretical studies [2][3][4][5][6][7][8][9] aiming at the ultra-fast control of magnetism without the need for magnetic fields. Despite being the subject of intense research for over two decades, the underlying mechanisms that govern the initial state of laser-induced demagnetization and the subsequent paths to equilibrium are not fully clarified yet. Clearly, different (electronic, phononic, spinrelated) degrees of freedom are dominantly active at particular timescales [1,3,[5][6][7][8][9][10]. Generally, one may distinguish between thermal effects upon laser heating and non-thermal effects, the latter corresponding to faster all-optical demagnetization and switching. The dynamics in this case (typically atto-to femtoseconds in metals) is dominated mainly by coherence, charge transfer and dephasing effects among the laser-excited electronic states. These aspects can be captured quantum mechanically. However, a suitable framework is needed to describe highly excited spin-dependent manyelectron systems. Clearly, theory is indispensable for clarifying the non-equilibrium spin-dependent charge quantum dynamics by providing simulations for key a e-mail: tomas.perie-de-barros@physik.uni-halle.de b e-mail: miguel.marques@physik.uni-halle.de (corresponding author) quantities such as the electronic charge and the charge and spin current densities. Ab initio simulations have addressed the role of spin-orbit interaction for Ni and Co films [3]. Other ab initio studies demonstrated the influence of finite size effects [11], or the possibility to demagnetized correlated magnetic insulators [12].
To date, most ab initio studies are performed using time-dependent spin-density functional theory (TDDFT) considering only the electronic degrees of freedom and neglecting the phonons. These calculations are usually done within the adiabatic local spin-density approximation (LSDA), which is unfortunately known to suffer from several weaknesses, including its incapability to describe correctly the electronic structure of 3d transition metals such as bulk iron or nickel [13]. The question of how well the electronic properties of the initial magnetic state need to be described by ab initio simulations, and which implications the quality of the initial state means for the light-induced demagnetization remain so far elusive. Beyond the problem of properly describing the electronic property of materials, the LSDA is also known to produce an exchange-correlation magnetic field parallel to direction of the local magnetization, which leads to the absence of the so-called exchange torque, i.e., a torque exerted by the local magnetic field on the magnetization [14]. This absence might hamper the description of the optical-spintransfer torque [15]. Another problem of the description of the demagnetization using the LSDA is that the predicted values of demagnetization are much smaller than the experimental ones, when considering similar laser parameters. Indeed, as reported by several prior studies [16][17][18], TDDFT within the adiabatic LSDA strongly underestimates the amount of demagnetization, and more intense laser pulses, having an intensity of up to two orders of magnitude larger than the experimental ones, need to be employed to reach similar demagnetizations. A recent study based on a combination of TDDFT and dynamical mean-field theory [18] demonstrated that the usually missing memory effects in exchange-correlation functionals play a key role in the proper description of the amount of demagnetization induced by light, and that including them leads to a larger demagnetization. The role of phonons is also starting to be investigated [19]. These studies partly fill the gap between theoretical predictions and experiments but the picture is far from complete.
In view of the wide classes of magnetic compounds, to develop a picture of laser-induced magnetism at the initial excitation stage one has to establish the link between electronic structure, exchange-correlation functionals, and their implication for the magnetic dynamics. Our aim here is to investigate, with the established tools of TDDFT, the role of electronelectron interaction on the early stage of the demagnetization process, which occur while the laser field is still driving the material, leaving aside the problem of missing memory effects in the exchange-correlation functional. For this, we employ the DFT+U method [20][21][22][23], which is known to improve substantially the electronic properties of 3d transition metals compared to the LSDA [24], while retaining a tractable numerical cost. How the localization of electrons in the d orbitals of the transition metals, and the resulting change in the electronic structure of simple ferromagnets affects their light-induced demagnetization is a central question of this work. Moreover, as discussed later in the paper, the use of Hubbard U helps us unravel the central role of the joint density of state in the demagnetization process.
Apart from the choice of the exchange-correlation functional, TDDFT simulations have the potential to guide experiments. To be predictive, the effect of laser field parameters, such as pulse duration, intensity, or wavelength, need to be carefully understood. Given the numerical burden associated with the simulations of all-optical demagnetization, the effect of these experimental knobs on the magnetization dynamics is still not fully elucidated, and we also aim here at a better understanding of how these parameters affect the demagnetization of elemental transition metals.
The outline of this paper is as follows. First, in Sect. 2, we give details of the ab initio method used to model the light-induced electron dynamics in the elementary ferromagnetic transition metals Fe, Ni, and Co considered in this work. Then, in Sect. 3, we study effect of the electron-electron interaction, as well as of the laser parameters on the ultrafast demagnetization of these materials. We then present our conclusions in Sect. 4.

Method
All simulations presented in this work rely on a realtime TDDFT framework. The time-dependent wavefunctions are described by Pauli spinors and are computed by propagating the Kohn-Sham equations, as provided by the Octopus package. [25] In atomic units, used throughout the text, the time-dependent Kohn-Sham equations within the adiabatic approximation reads where |ψ n,k is a Pauli spinor representing the Bloch state with the band index n at the point k in the Brillouin zone, A(t) is the external vector potential describing the laser field,v ext is the electron-ion potential containing the non-local pseudopotential, also describing spin-orbit coupling,v H is the Hartree potential, and v xc is the exchange-correlation potential given by the adiabatic LSDA. Finally,v U is the potential coming from the DFT+U term [26], and that is added to the d orbitals to improve the electronic localization around the transition metal atoms. The potentialv U is given in the case of noncollinear spin by [26] v U |ψ n,k = where I refers to atomic site, m and m to the localized d orbitals |φ I m and |φ I m of this site, and σ, σ are spin components. We defined where U eff I is a parameter and n I,σσ is the so-called occupation matrix of the localized orbitals attached to the atomic site I. It reads where f nk is the occupation of the spinor state |ψ nk , and w k is the weight associated to the point k. In order to better understand the early stages of the demagnetization process, the ions are clamped in all our simulations.
Calculations for body-centered cubic (BCC) iron were performed using a lattice parameter of 2.856Å, a uniform spacing for the real-space grid used in Octopus of 0.111Å, and a 12 × 12 × 12 k-point grid to sample the Brillouin zone. For the case of face-centered cubic (FCC) cobalt, a lattice parameter of 3.544Å, a spacing of 0.132Å, and a 16 × 16 × 16 k-point grid were used. Finally, for FCC nickel, calculations were performed using a lattice parameter of 3.520Å, a spacing of 0.132Å, and a k-point grid of 14 × 14 × 14. We employed norm-conserving fully relativistic Hartwigsen-Goedecker-Hutter pseudo-potentials [27]. Time-dependent calculations were performed using a time-step of 0.73 as for Fe, 1.21 as for Co and 0.98 as for Ni.
In all cases the laser pulse used to induce the demagnetization dynamics was linearly polarized along the xaxis and perpendicular to the direction of the magnetic moment. Unless stated differently, we employed a laser pulse with a cosinoidal envelope of 6 fs duration, a peak intensity of 10 15 W/cm 2 and a central wavelength of 457.1 nm which corresponds to a carrier photon energy of 2.712 eV.

LSDA calculations
We start our analysis by looking at the all-optical demagnetization, as described at the level of the standard LSDA. The dynamics of the magnetic moments for bulk Fe, Co and Ni induced by this short and intense laser pulse is presented in Fig. 1. In all three metals demagnetization is observed, with Fe demagnetizing around 9%, Co around 6% and Ni around 8%. We note that the values for the demagnetization, though smaller than what is observed experimentally [28,29] where the demagnetization measured is much more pronounced, are in accordance with other similar theoretical studies [18].
Before commenting on the effect of adding correlations in the form of an Hubbard U , let us comment on the part of the demagnetization after the end of the laser pulse. In sharp contrast with experimental results, [30,31] we find that the demagnetization stops after the end of our laser pulse. This result is in agreement with prior theoretical results where only electronic effects have been considered, see for instance [12,18,32]. It is in fact possible to understand this for the equation of motion of the magnetization, when ionic motion is neglected. One can derive the equation of motion of the magnetizationm(x, t), starting from a locally gauge-invariant U(1)xSU(2) Hamiltonian, including spin-orbit coupling (see Ref. [33] for the derivation). One obtains that  where B is the (Kohn-Sham) magnetic field, E is the total electric field acting on the electrons. Here e is the elemental charge, m is the electron mass, and c is the speed of light. We used the physical spin current J phys,pk , as defined in Ref. [33]. Here ijk denotes the Levi-Civita symbol.
From Eq. (5), we can analyze the equation of motion of the macroscopic magnetization of the Kohn-Sham system. The divergence of the Kohn-Sham spin current vanishes as we integrate over space, and we get the equation of motion of the total magnetization Without SOC, the second term is absent and we find the usual result that, in absence of SOC, there is no change in the magnetization, unless the functional violates the zero-torque theorem (which is not the case for LSDA, as the LSDA exchange-correlation torque is uniformly zero).
Let us now analyze the term that originates from the SOC, which is the second term on the right-hand side of Eq. (6). It is made of two ingredients, the electric field acting on the electrons E, which is given by and the physical spin current J phys . Now, when the laser ends, i.e., when A laser (t) is zero, we are left with the spin current, and the gradient of the Kohn-Sham potential v s . Therefore, if no spin current is flowing in the system (on average) at the end of the laser pulse, there is therefore no demagnetization occurring. This is clearly the case in the linear response regime and when ions are not moving, as the average value of the (spin) currents are simply proportional to the perturbing electric field. As we shown later, changing the intensity of the perturbation, which should enhance the induced spin current, does not show signatures of demagnetization at the end of the laser pulse. The oscillations of the magnetization after the end of the laser pulse are attributing to the fluctuations in time of the spin-current in time around a vanishing averaged value.
The time-average of the spin current could still be non-zero at the end of the laser pulse, due to nonlinear effects like the shift current, but this is usually forbidden by symmetries in cubic materials and for a linearlypolarized driving field. This could only occur if a ultrashort pulse is used, for which the time integral of the vector potential is non-zero, leading to a population imbalance in momentum space and to a macroscopic finite current, but this is not the case here.
Moreover, we note that if the spin current is only macroscopic (uniform over the unit cell), what remains as a contribution to the macroscopic magnetization is the integral over space of the gradient of the Kohn-Sham potential, which leads to zero contribution. Therefore, in any system, when the laser pulse is over, only the time-average of the microscopic part (i.e. the local fluctuation of the current within the unit cell) of the spin current is playing a role in the demagnetization process, and this is expected to be a small contribution. The microscopic spin current also most likely becomes zero at the end of the pulse, as one would expect from linear response, leading to no change in the magnetization. As we will see later in Fig. 6, there is no indication that for higher intensities demagnetization occurs after the end of the laser pulse, so it seems to be valid to assume no remaining spin current at the end of the laser pulse.

LSDA + U calculations
To better understand the microscopic origin of the process of all-optical demagnetization, we now address at first the role of the equilibrium bandstructure in the demagnetization process. Conceptually, at the early stage of the coherent optical excitation, changing the equilibrium bandstructure must affect the process of laser excitation of carriers to the conduction bands, Fig. 2 Excess energy for each material, using the same laser pulse as in Fig. 1, for both the LSDA and LSDA + U cases. This is computed as the energy difference between the beginning and the end of the simulation as the optical transitions available at a given photon energy are affected by a change in the bandstructure. To investigate how the description of electron localization around the d orbitals in these metals affects the demagnetization process, we repeated the same calculations but applying an effective Hubbard correction U eff = U − J to the d orbitals of 2 and 4 eV for Fe, 1.6 and 3.2 eV for Co and 1.2 and 2.4 eV to Ni. While the inclusion of an Hubbard U is not necessarily sufficient to really improve the description of the density of states of these metals, it offers us a nice tool to study the effect of the bandstructure of the material on the demagnetization process.
The results are also shown in Fig. 1. Overall, we find that the dynamics of the magnetic moments for all the elemental ferromagnets is similar to the one obtained using the LSDA. However, we note that in the case of Fe and Ni the amount of demagnetization is reduced in the LSDA + U case compared with the standard LSDA case while in Co the Hubbard U increases the amount of demagnetization.
To understand the reason for this change in demagnetization, we plot in Fig. 2 the energy absorbed by the materials (the difference in the system's total energy before and after the laser pulse) depending on the energy functional used to describe the electron-electron interaction. We find that these elemental ferromagnets demagnetize less when they absorb less energy. This directly points towards the direction of less electrons being excited by the laser field, which we should be able to relate to the bandstructure of the material at equilibrium. The Hubbard U correction introduces a shift in the d orbitals of the metals which affects their bandstructure, and hence the number of optical transitions available at the laser frequencies. This, in turn, is directly reflected in the amount of energy absorbed and in the demagnetization process.
To corroborate our interpretation, we computed the joint density of states (JDOS) for each of the elemental solids, at both the LSDA and LSDA+U level, and we compared it with the energy spectrum of the considered laser pulse. Our results, shown in Fig. 3, confirm our interpretation. Indeed, we observe that the presence of the Hubbard U correction can significantly change the JDOS. For the Fe case the JDOS is found to be smaller for the LSDA + U case in the region of the laser excitation, which is in accordance with our interpretation that with the Hubbard U correction there are less transitions available to the system at the laser frequencies. In the Co and Ni cases things are more complex as the relative strength of the LSDA and LSDA + U JDOS changes over the region of the laser frequencies. We can make things clearer by choosing a laser with a central frequency in a region where the JDOS is markedly larger in one case over the other. For that we choose a laser with a carrier photon energy of 1.2 eV, the spectrum of which can be seen plotted in Fig. 4 over the Co JDOS in the LSDA and U = 3.2 eV cases. In this frequency range the JDOS is clearly larger at the LSDA level. The corresponding demagnetization dynamics for this laser field are shown in Fig. 5, where we note the system demagnetizes more in the LSDA case, unlike with the original laser pulse in Fig. 1 but as expected based on the JDOS. This again indicates that the magnitude of demagnetization observed here is directly related to the amount of optical transitions available. This result is important and useful in the search for materials suitable for alloptical demagnetization.

Effect of the laser parameters
We now turn our attention to the effect of the laser parameters on the demagnetization process. We start by investigating the effect of the peak field intensity, keeping the other laser parameters fixed. If the results are merely due to electronic transitions between valence and conduction bands, then an higher peak intensity should imply more transitions and so a larger demagnetization. In particular, in a picture of demagnetization being initiated by excitation of carriers to the conduction bands, the amount of demagnetization should be dictated by the number of electrons being excited to the conduction bands. Assuming a single-photon absorption, this should therefore scale linearly with the intensity of the laser field. To check this, we performed simulations for Ni, using laser fields of different intensities, while keeping the pulse duration and the frequency fixed (see Fig. 6).  We can see that, as expected, the lower the peak intensity, the less the system demagnetizes. The amount of demagnetization for each given intensity is plotted in Fig. 7. It is clear that demagnetization scales roughly linearly with the intensity of the laser field, as we expected. We can also see that as the intensity increases the demagnetization starts to stray from this linear scaling, possibly as a result of non-linear effects given the high intensity of the laser pulses used.
This again suggests that the demagnetization, as described at the LSDA level, is tightly linked to the optical transitions available at the laser frequency.
One question naturally rises from this result, namely whether the laser intensity or the laser fluence is the suitable quantity to understand demagnetization? Therefore, we considered different laser fields with the same central frequency, but different pulse duration T and peak intensities I 0 , such that the fluence for all the considered pulse remained equal. Clearly, our results, shown in Fig. 8, indicate that the laser fluence is not a good parameter to quantify how much a given laser pulse can demagnetize the material, as we can tune the  Excess energy for Co under the influence of different laser fields with the same fluence. The color coding for the different laser pulses is the same as in Fig. 8 amount of demagnetization, as well as the timescale at which it occurs, by changing the intensity and pulse duration, while keeping the fluence fixed. This is true, even though the total amount of energy carried by the different laser fields is the same. The energy absorbed by the system also changes under the different laser fields, as depicted in Fig. 9, and it is smaller in the cases where the system demagnetizes less.
We note, as in previous cases, that the demagnetization is only present while the laser field is on. These new longer simulations reveal an important facet of laser induced demagnetization, which is not clearly visible for shorter laser pulses: the magnetization seems to occur mostly within the first half of the laser pulse and then saturates slightly after the middle of the laser pulse. Indeed, in the case of the laser field with an intensity I 0 = 2 × 10 14 W/cm 2 and a total duration of T = 30 fs, the relative magnetization oscillates around a fixed value from approximately 20 fs onward. The last part of the driving field has no meaningful effect on the magnetization. This effect is related to the excitation to the conduction bands. Indeed, in the view of singlephoton absorption, the number of electrons excited by a Fig. 10 Detail of the joint density of states for Co with the frequency spectrum of the different laser pulses used superimposed. The laser spectrums are color coded the same way as in Fig. 8 laser pulse resemble the time integral of the pulse. Said differently, most of the ionization takes place around the field maximum, and not towards the end of the pulse. The fact that the demagnetization process resembles the time integral of the laser pulse is yet another indication of the tight connection to the number of electrons excited in the conduction bands.
Another important aspect of the numerical results shown in Fig. 8 is that shorter pulses lead to more significant demagnetization in Co, at least for the laser parameters considered here. Clearly, as we employed a more intense laser field, we are exciting more electrons within each half cycle of the laser field. However, if we assume this process to be linear in intensity, we should get for the same fluence the same demagnetization, which is clearly not the case from the simulations, especially for ultrashort pulses. There is therefore another effect coming into play that affects most the shortest laser pulses. Looking at the frequency spectrum of these different laser fields compared to the JDOS of cobalt (Fig. 10), it is clear that the shortest pulse we considered is so short that it covers a wide range of frequencies, with a full width at half maximum of roughly 1.5 eV. This implies that this laser field contains a large set of frequency components, and therefore has access to a larger range of optical transitions, given its broader energy spectrum. When the number of optical transitions varies significantly over the energy range covered by the laser spectrum, the fluence stops being a valid quantity to describe demagnetization, and one should consider instead pulse duration and intensity as two independent and equally relevant parameters having different effects on the light-induced demagnetization. On the other end, the two longer laser fields only differ mildly on the bandwidth, and the total JDOS does not vary much over this energy range. We therefore expect that for even longer pulses, the fluence is a good quantity to quantify how much a laser demagnetizes a given material, as such pulses would probe an energy region where the JDOS is essentially constant.

Conclusions
To conclude, in this work, we investigated how the description of electronic correlations, and the resulting bandstructure, influences the light-induced demagnetization of elementary ferromagnetic transition metals Ni, Co, and Fe. For this, we employed the LSDA + U method, which is known to improve the description of the electronic properties of these materials with respect to the simple LSDA. We showed that these materials generally demagnetize less at the LSDA + U level than when described using the LSDA, and we related this to the fact that at the LSDA + U level, the material absorbs less energy, for the optical frequency considered here, because of the change in the electronic structure of the material.
We also performed a detailed analysis on how the laser parameters such as pulse duration and intensity influence the demagnetization process. Importantly, we showed that the amount of demagnetization for the various laser conditions can be understood on the basis of the possible optical transitions, given by the jointdensity-of-states of the material and the energy position and bandwidth of the considered laser field. In particular, our study reveals that shorter laser pulses have the potential to demagnetize more than longer pulses, the latter having an energy bandwidth too narrow to cover enough optical transitions to induce a significant demagnetization. It also appears, from our analysis, that the laser fluence is not a valid quantity to determine whether a laser pulse can demagnetize a material or not when considered ultrashort pulses.
It should be mentioned that many aspects of the demagnetization still need to be further investigated. In the context of the DFT + U calculations, the effect of dynamical correlations, in the form of a time-dependent Hubbard U [34] still needs to be explored. More generally, the effect of laser-induced band renormalization, Pauli blocking, and electron-electron scattering needs to be investigated. Exchange-correlation functionals with non-zero exchange torque need to be developed for non-colinear magnetism and tested in the context of all-optical demagnetization, in order to better understand the potential role of optical torques in this process. Other important aspects are the choice of exchange-correlation functionals and how they affect the material-specific magnetic disordering in antiferromagnets, spin glasses, or rare-earth ferromagnets.

Author contributions
TB performed the calculations. All authors contributed equally to the interpretation of the results and the writing of the article. Data availability The data that support the findings of this study are available from the corresponding author upon reasonable request

Declarations
Conflict of interest The authors declare no competing financial or non-financial interests.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder.