We study the influence of transport effects on time- and space-resolved magnetization dynamics in a laser-excited thick nickel film. We explicitly include diffusive heat transport and spin-resolved charge transport as well as Seebeck and Peltier effects and calculate the dynamics of spin-dependent electronic temperatures, chemical potentials, lattice temperatures, and magnetization. We find that transport has an influence on the magnetization dynamics closer to the excited surface as well as in regions deeper than the penetration depth of the laser. We reveal that, for higher absorbed fluences and in the presence of transport, thick magnetic films show a quenching time nearly independent of depth, though the magnitude of quenching is depth-dependent.

Ultrafast demagnetization dynamics combines the utility of conventional electronics and the energy efficiency of spintronics with greater operational speeds.^{1} Since it was discovered by Beaurepaire *et al.*,^{2} enormous progress has occurred in the study and applications of light-driven,^{3–5} spin-driven,^{6–9} and current-driven^{10,11} magnetization dynamics. Recent developments also link the phenomenon with ultrafast magnonics,^{12–15} spincaloritronics,^{16,17} and terahertz physics.^{18}

With the ubiquity of transport mechanisms relevant for ultrafast spintronic devices, there is a need to better understand the role of transport in ultrafast demagnetization dynamics. To this order, the influence of several mechanisms, such as diffusive, superdiffusive, and ballistic transport, have been studied theoretically as well as experimentally.^{19–21} In the case of magnetic heterostructures as well as thin magnetic films, atomistic simulations based on Landau -Lifshitz -Gilbert models,^{13,22} time-dependent Density Functional Theory (DFT) calculations,^{23,24} and Monte Carlo simulations^{25,26} have been proven successful in predicting, modeling, and interpreting experimental observations. In contrast to a thin film, exciting a thick magnetic film with an ultrafast laser pulse leads to nonuniform heating and results in gradients of temperature as well as of particle densities or chemical potentials. Understanding the mutual dependence of different transport effects on magnetization dynamics is an urgent challenge for the development of spintronic devices.

In this work, we study the influence of transport on the spatially resolved magnetization dynamics. To that end, we apply the thermodynamic *μ*T-model^{27} and introduce transport effects caused by the laser-induced transient gradients of temperatures and chemical potentials. Namely, we explicitly include diffusive heat transport, spin-resolved charge transport as well as Seebeck and Peltier effects and determine the temporal and spatial changes of spin-resolved electronic temperatures, chemical potentials, lattice temperature as well as magnetization. In order to disentangle the effect of the depth-dependent laser absorption profile and the effect of transport on spatiotemporal magnetization, we consider the extreme cases when (i) all transport effects are set to zero (without-transport) and (ii) all transport effects are included simultaneously (with-transport). The scenarios of without- and with-transport are also comparably applicable to magnetic metals with weak and high conductive properties. Our results show a strong depth-dependence of the measurable characteristics like maximum magnitude of quenching as well as quenching time, which we define as the time to reach the minimum of magnetization.

The dynamics of spin-resolved electronic $ue\sigma $ and phononic energy densities *u _{p}* as well as the spin-resolved particle densities $n\sigma $ are described by the continuity equations given by

The electronic energy and particle densities are determined by the space (*z*)- and time (*t*)-dependent transient variables of electronic temperatures $Te\sigma (z,t)$, chemical potentials $\mu \sigma (z,t)$, and magnetization *m*(*z*, *t*). They are connected to each other through the zeroth and first moment of the Fermi distribution of the spin-resolved electrons. The energy density of phonons depends on phonon temperature $Tp(z,t)$. In order to simplify the heat-flux and particle-flux terms, they are Taylor-expanded with respect to these transient variables. We consider the diffusive transport to be predominant in metallic ferromagnets;^{21} therefore, only the expansion up to first order derivatives are considered. Therefore, the *Ansätze*,

*Ansatz*(2). Figure 1 illustrates the electronic transport channels included in our approach. Additionally, the energy transport due to gradients in the phononic temperature is considered via thermal conductivity of the lattice

*κ*, see Eq. (2b).

_{p}In the following, we determine the transport parameters mentioned above. The particle flux density of electrons, $j\u2192e$, can be identified with the electric current per elementary charge *e*. With the electrical conductivity $\sigma cond$, an internal electric field $E\u2192int$, and Ohm's law, we obtain

A gradient in chemical potential introduces a spin-voltage, which evokes an electrical field of $E\u2192int\sigma =\u2207\mu \sigma e$. Thus, we obtain the particle flux density of electrons as

Comparing with Eq. (2c) for the case of vanishing temperature gradient, we identify

A gradient in the chemical potential also induces a heat flux density through the Peltier effect, see Fig. 1. It is given as $q\u2192e=\u2212e\Pi j\u2192e,$ where Π is the Peltier coefficient. Consequently, we find through Eqs. (2)

The inverse effect, where a temperature gradient leads to the electric current, is called the Seebeck effect. The induced internal Seebeck field is connected to the temperature gradient as $ESeebeck=\Sigma \u2207Te$, where Σ is the Seebeck coefficient. With (3) and (2c) for the case of vanishing gradient of chemical potential, we directly obtain

Finally, we determine the heat conduction parameter *κ _{uT}*. In the absence of particle currents, the heat flux density is proportional to the temperature gradient via Fourier's law

^{28}

where *κ* is the common electronic thermal conductivity. Setting $j\u2192e=0$ in the equation system (2), we obtain

and thus,

which introduces a small correction to the common thermal conductivity *κ*, caused by all other considered transport effects. Note that Peltier- and Seebeck coefficients are connected through the Onsager relation, $\Pi =\Sigma Te$.^{28} The calculated values of transport coefficients at $Te=300\u2009K$ are given in the supplementary material.

Following Ref. 27, the partial derivatives of the spin-resolved electronic energy densities $ue\sigma $ and particle densities $n\sigma $ as well as phononic energy density *u _{p}* are expressed in terms of the temperatures, chemical potentials, magnetization, electron-lattice coupling parameter $g\sigma $, coupling parameter between the chemical potentials,

*ν*, and the inner-electronic temperature coupling,

*γ*,

^{27}as

where $S\sigma (r\u2192,t)$ is the energy input to the electrons^{29} of the majority- and minority electronic bands. The dynamics of the magnetization $m=(n\u2191\u2212n\u2193)$ is determined by

The temporal derivatives of the energy- and particle density are expressed in terms of capacity-equivalents for energy densities $cx=\u2202ue\sigma \u2202x$ and analogously for particle densities $px=\u2202n\sigma \u2202x$ with $x={Te\sigma ,\mu \sigma ,m}$, permitting the left-hand sides of Eq. (1) to be written in terms of transient variables ${Te\sigma ,\mu \sigma ,m}$. The calculations of these capacity-equivalents are described in the supplementary material.

Inserting the total differentials

and Eqs. (11) into Eq. (1) and rearranging (see the supplementary material), we obtain a coupled nonlinear partial differential equation of the form

with $X\u2192T(z,t)=(Te\u2191,Te\u2193,Tp,\mu \u2191,\mu \u2193,m)$ being the time- and space-dependent vector of transient variables; $C,\u2009G$, and $K$ are the matrices of capacity-equivalents, coupling parameters, and transport coefficients, respectively (see the supplementary material), and $S\u2192(z,t)=(S\u2191(z,t),S\u2193(z,t),0,0,0,0)$ is the source term corresponding to excitation of the electronic subsystems.

We solve the thermodynamic *μ*T-model, Eq. (14) for a 50 nm thick nickel film excited with a Gaussian light pulse with a wavelength of 800 nm and an absorbed fluence of $F0=5\u2009mJ/cm2$. The full width at half maximum $tFWHM=$ 50 fs, centered at $t0=$ +50 fs was chosen. The decay of absorbed energy within the material is modeled following Beer–Lambert's law with a laser penetration-depth of $\lambda pen=$ 15 nm, which is typical for laser pulses in visible wavelengths.^{30} Assuming equal absorption coefficients for majority and minority electrons, the energy of the laser pulse is divided equally between the spin-resolved electronic systems $S\u2191=S\u2193=S/2$ with

being the total laser source. The capacity matrix $C$ is calculated dynamically as a function of the transient spin-resolved electronic temperatures and chemical potentials as described in the supplementary material by applying the DFT-calculated density of states for nickel.^{31} The elements of the coupling matrix $G$ are taken from Ref. 27, and the elements of the transport coefficient matrix are given in the supplementary material. The initial Stoner exchange at $T=300\u2009K$ is taken to be $\Delta (z,t=0)=0.5\u2009eV$.^{32}

The coupled nonlinear partial differential equation (14) is converted to a set of coupled nonlinear difference equations using the “forward in space—centered in time” discretization scheme.^{33,34} The initial electronic and lattice temperatures are 300 K, and the initial spin-resolved chemical potentials are determined to be 0.0047 eV below the Fermi level (see the supplementary material). The time stepping is done through the fully implicit discretization scheme and the Thomas algorithm^{34} generalized to coupled difference equations.^{33,35}

We study the influence of spin-resolved charge and heat transport on spatially resolved dynamics of temperatures, chemical potentials, and magnetization. In the following, we compare the dynamics at the surface, where the highest energy absorption from the laser pulse takes place with the dynamics in the center of our simulated film, where the laser excites the material only weakly.

Figures 2(a) and 2(b) show the dynamics of spin-resolved electronic and lattice temperatures at the surface and at a depth of 25 nm in the absence (dashed lines) and presence (solid lines) of heat and particle transport. The electronic temperatures sharply rise in a few femtoseconds following laser irradiation until maxima are reached. The maximum temperatures as well as the time at which maxima occur depends on depth, spin and on transport. Though the energy absorbed from the laser pulse is equally distributed to the majority and minority electrons, the maximum temperature of minority electrons is higher than that of majority electrons. This is the consequence of the lower heat capacity of the minority electrons as compared to the majority electrons. The maximum electronic temperatures at the surface are higher than the ones at 25 nm due to the decreasing laser energy absorption profile.

In the absence of transport, the absorbed energy at the surface remains trapped close to the surface, whereas in the presence of transport, the colder regions in the material act as an energy and particle sink. Therefore, the maxima of temperatures at the surface are higher when transport is absent in comparison to the case when it is permitted. At 25 nm, the heat- and particle flow from the hotter surface into the colder bulk acts as an additional driver for temperature dynamics. Therefore, the maxima of electronic temperatures in the bulk are higher when transport is present.

Energy and particle transport, implicitly and explicitly, influence the dynamics of spin-resolved chemical potentials. They are shown in Figs. 2(c) and 2(d) at the surface and at a depth of 25 nm, respectively, in dependence on time. The curves are depicted for the case without transport with dashed lines and for the case of included transport effects with solid lines. At both depths, the chemical potential of the majority electrons (blue) initially steeply increases following the absorption of the laser pulse whereas the chemical potential of minority electrons falls. This is a direct consequence of the spin-resolved density of states. For the majority-electrons, the d-band is fully occupied, and the density of states decreases at the initial chemical potential. Therefore, for increasing temperatures, the chemical potential shifts to higher energies, ensuring particle conservation. For minority electrons, the case is vice versa: The d-band is only partially occupied, and the density of states is increasing at the initial chemical potential. Thus, upon excitation, particle conservation leads to a shift of the chemical potential to lower energies. Two coupled reservoirs with different chemical potentials exchange particles until the chemical potentials equilibrate. Since the majority electrons have a larger chemical potential, spin-flip scattering mainly transfers electrons from the majority to minority subsystem. This process was identified in Ref. 36 as one of the driving forces behind ultrafast demagnetization. Similar to temperature-dynamics, the dynamics of chemical potentials are weakened at the surface and strengthened at the bulk when transport is included. Details of the dynamics of chemical potentials are explained in Ref. 27. Note the transient equilibrium of the chemical potentials at the crossover point of $\mu \u2191$ and $\mu \u2193$, which marks the minimum of the transient magnetization.^{27} At the surface, this equilibration is accelerated, when electronic transport is active.

The resulting magnetization dynamics is shown in Fig. 3 as a function of time at the surface and at a depth of 25 nm, respectively, in the absence (dashed lines) as well as the presence (solid lines) of transport. In all four curves, the magnetization dynamics shows rapid drop following laser excitation and slower recovery. Transport influences the magnitude of quenching as well as the time when maximal quenching is reached (quenching time). In the presence of transport, quenching at the surface (red solid line) is weaker and faster compared to the case where transport was absent (red dashed line). This is due to loss of energy into colder bulk and lowering of electronic temperatures at the surface, compare Fig. 2(a). As in the case of the transient temperatures and chemical potentials, see Fig. 2, transport enhances the effect in the bulk and, thus, leads to stronger quenching at a depth of 25 nm as compared to the case without transport.

The spatial profiles of the spin-resolved temperatures, chemical potentials, and magnetization are shown in Fig. 4 for the time-instants of 100, 500, and 3000 fs. Figure 4(a) depicts the spin-resolved electronic temperatures in the presence of transport. The temperatures of minority electrons are higher than that of majority electrons at 100 fs, independent of depth. At 500 fs, this trend has reversed. At 3000 fs, the spin-resolved electronic temperatures appear to be equal. Also the spatial gradients have decreased considerably.

Figure 4(b) shows spatial profiles of the spin-resolved chemical potentials, when all considered transport processes are included. The horizontal black line represents the initial chemical potentials at 300 K, located at −0.004 73 eV below the Fermi energy. Following laser irradiation, the chemical potentials show a strong separation at 100 fs. The stronger gradient of the majority profile as compared to the gradient of the minority profile suggests that the transport of majority electrons plays a stronger role in depth-dependent magnetization dynamics. This is in agreement with Ref. 37. At 500 fs, both chemical potentials are below the Fermi level. At 3000 fs, the chemical potentials show nearly no gradients anymore but are not yet equilibrated.

Figure 4(c) shows the magnetization, normalized to its initial value, as a function of depth at various time instances after laser excitation in the presence (solid lines) and the absence (dashed lines) of transport. The quenching is given as the difference to the initial value. In the absence of transport, the magnetization profiles are monotonic, and the quenching decreases exponentially with depth due to Beer–Lambert's absorption profile. The maximal quenching occurs at the surface while it is smallest at the depth of 50 nm. When transport is present, this exponential behavior is lost. The weakening of demagnetization (in comparison to the no-transport cases) is strongest at the surface and gradually decreases proceeding inwards until time dependent inflection points are reached. Beyond a certain time-dependent depth, transport leads to strengthening of quenching at all depicted time-instances.

Figure 5 shows the relation between the magnitude of quenching and the quenching time (time to reach the maximum quenching) in the presence and the absence of transport for several absorbed fluences at various depths of the material. The depths are color-coded and range from the surface to the total thickness of the simulated film. In the absence of transport, solely the local energy input determines the magnetization dynamics. Thus, different depths are equivalent to different irradiation strengths, and we obtain a joint curve for all fluences. In contrast, the relation between maximum quenching and quenching time becomes strongly dependent on the absorbed fluence, when transport processes influence the dynamics. Figure 5 shows for each applied fluence, a narrower spread of the quenching times in the presence of transport. For a given fluence, the quenching time appears to be rather independent of depth, while the maximum magnitude of quenching is depth-dependent.

As compared to the case without transport, quenching is weakened but faster in the vicinity of the surface, whereas it is strengthened and delayed within the bulk. This is due to the heat and particle currents affecting the transient variables in opposing ways. The microscopic mutual interplay of local and nonlocal processes needs further investigation. Quenching times and magnitudes of quenching are experimentally accessible quantities. The relations shown in Fig. 5 could be, in principle, experimentally verified in magnetic materials using depth- and time-resolved magneto-optical Kerr Effect measurements.^{38}

In summary, we explicitly included spin-resolved charge- and heat-transport as well as Seebeck and Peltier effects into the thermodynamic *μ*T-model. We studied the spatially resolved dynamics of temperatures, chemical potentials, and magnetization. We then studied depth-resolved quenching and demagnetization times for several absorbed fluences. Additionally, it was revealed that, for higher absorbed fluences and in the presence of transport, thick magnetic films show depth-independent demagnetization time though the quenching is depth-dependent. Our theoretical model is sufficiently complex to include several underlying transport mechanisms and at the same time computationally simple enough to be extended to heterostructures and magnetic multilayers.

See the supplementary material for calculation details of the partial heat capacities, numerical values of individual interaction and transport coefficients, and additional computational details.

The research leading to this study was funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation)—TRR 173—268565370 Spin+X (Project No. B03). We thank Jonas Hoefer and Martin Stiehl for helpful discussions.

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

## DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding author upon reasonable request.