Theory of spin dynamics of magnetic adatoms traced by time-resolved scanning tunneling spectroscopy

The inelastic scanning tunneling microscopy (STM) has been shown recently (Loth et al. Science 329, 1628 (2010)) to be extendable as to access the nanosecond, spin-resolved dynamics of magnetic adatoms and molecules. Here we analyze theoretically this novel tool by considering the time-resolved spin dynamics of a single adsorbed Fe atom excited by a tunneling current pulse from a spin-polarized STM tip. The adatom spin-configuration can be controlled and probed by applying voltage pulses between the substrate and the spin-polarized STM tip. We demonstrate how, in a pump-probe manner, the relaxation dynamics of the sample spin is manifested in the spin-dependent tunneling current. Our model calculations are based on the scattering theory in a wave-packet formulation. The scheme is nonpertubative and hence, is valid for all voltages. The numerical results for the tunneling probability and the conductance are contrasted with the prediction of simple analytical models and compared with experiments.


I. INTRODUCTION
Spin systems are extensively utilized as essential elements for (quantum) information storage or computing devices [1,2]. Thereby, environmental effects play an important role as they generally lead to decoherence and relaxation. The time scales for these processes depend on the underlying coupling mechanisms and exhibit in general a marked temperature and dimensionally dependence.
A detailed insight and a possible control of this relaxation is a key factor for the operation of these devices. Desirably the relaxation time should be larger than the operation (switching) time. In this sense, it is of importance to identify spin systems with suitable relaxation properties and amenable to full fast control of the spins. Prototypical examples are realized in semiconductor nanostructures [3,4] or nanomagnets [5,6]. For a single magnetic atom or molecules adsorbed on a substrate, the surface spin excitations are usually weakly coupled to bulk magnons resulting in relaxation times on the order of nanoseconds, which has been demonstrated to be sufficient for a coherent spin manipulation [7].
The experimental analysis of surface-deposited structures can be performed conveniently by means of the scanning tunneling microscopy, which has an advantage of an atomic spatial res-olution. Using a spin-polarized tip, the spin of a single magnetic adatom can be probed and manipulated [8], since the tunneling probability and the current depend on the relative alignment of the surface and the tip magnetic moments. Very recently the same group demonstrated the potential of the technique for atomic-scale information storage and retrieval [9]. The time resolved STM experiments [10] renders possible the observation of the quantum dynamics. The transient precessional dynamics of the excited spin states is still too fast to be measured. The relaxation process can, however, be mapped onto the current-dependence in a pump-probe experiment. Thus one has a possibility to directly measure the relaxation time of such systems which opens the way for testing different configurations with maximal coherence time.
To our knowledge, the first experiment to measure the relaxation times via STM has been We will consider theoretically the experimental situation (section II) using scattering theory (section III) with the aim to calculate the spin-dependent transport properties. The method is complementary to other formalisms. A widely spread approach is to describe the tunneling by the master equation with transition rates calculated in a standard "golden rule" manner [11,12]. This method is, in essence, a perturbative formalism and hence requires a separate justification for each of the considered experimental configurations. The nonequilibrium Green's function approach has a potential to systematically take into account many-body effects, however, the transition rates are often introduced as adjustable parameters [13]. In the present method we focus on the description of the tunneling starting from a model Hamiltonian including the magnetic anisotropy, the exchange coupling and the coupling to the external magnetic field. Our theory proves that it is indeed possible to directly measure the spin-relaxation times provided they are longer than the precession time. Furthermore, our main findings for the relative number of transmitted electrons with respect to the pump-probe delay (presented in section V) allow for a rigorous determination of the magnetic coupling parameters.

II. MODEL HAMILTONIAN
Our model Hamiltonian matching the experimental situation in figure 1 includes which describes the tunneling electrons (with a mass m and momentum operator p) emanating from the tip and are subject to the magnetic field B. They are coupled to the surface spin subsystem via a spin-dependent term V B . The spin of the tip electrons is described by the vector operatorŝ τ = (τ x , τ y , τ z ) containing the Pauli matrices in their standard representation. In a similar way, S = (S x , S y , S z ) are the operators for the surface spin. The tunneling barrier V B (x) is measured relative to the Fermi levels of the tip or the sample (for simplicity we assume both to be equal).
which indicates that the easy axis is parallel to the (strong) magnetic field. We have chosen typical values for a rather strong magnetic anisotropy. This does not directly resemble the experimental situation, but serves to illustrate the effects arising from the anisotropic environment more clearly.
In particular, we take δ ≈ −5 meV and ε ≈ 1.5 meV for the anisotropy constants. The value of the magnetic field is kept with 7 T. The Landé factor g is experimentally [14] found to be approximately 2, which corresponds to the total spin of the free Fe atom (S = 2). The Hamiltonian (2) has 5 eigenstates |χ ν with eigenenergies E S ν (table I lists their properties). The ground state in orbital representation is depicted in figure 2.
The third term on the left side in equation (2) results in a mixing of the different eigenstates |m of the free spin in the magnetic field. Thus, the expectation value S z in the ground state deviates from the noninteracting value -2. The ground state contains a large amount of free spin configurations (see table I). The first excited state corresponds to the spin pointing in the opposite direction. Hence, the relaxation |χ 2 → |χ 1 is dominated by the transitions |m = +2 → |m = −2 .
For the spin-dependent interaction of the tunneling electrons with the adsorbate we adopt the picture of extended states (tunneling electrons) coupled to a localized magnetic moment. In brief, such a process is known to be determined by the Hunds rules and the local electronic correlation.
E.g. for s states (with energy ǫ k ) coupled to a pinned magnetic impurity (with energy ǫ d ) one may utilize the Anderson model [15]: H And = k,σ ǫ k n kσ + σ ǫ d n dσ +Un d↑ n d↓ +H ds , where σ labels the spins and n kσ (n kσ ) is the corresponding particle density operator. U stands for the Coulomb repulsion between electrons with antiparallel spins. The term H ds = k,σ V kd c † kσ c dσ + h.c. accounts for the s-d hybridization with the mixing matrix element V kd . The Anderson Hamiltonian Typical values [17] are in the order of |V k F d | 2 ≈ 2 eV, U ≈ 8 eV and ǫ d ≈ 5 eV, meaning that J ≈ 1 eV.
In STM experiments the tip-sample distance (we take 5Å) is usually not negligible to the typical extent of d-orbitals in Fe adatoms (which is about 1Å). In the calculations, we use therefore an interaction H int of the Kondo type with a spatial-dependent Kondo coupling The u 0 > 0 is the coupling strength and the function w(x) is a dimensionless form factor defining the width of the interaction region. We choose a strongly localized function with a width about

1Å.
The total Hamiltonain of the system is thus The total effective potential experienced by the tunneling electrons V B (x) + χ ν , ↑ |H int |χ ν , ↑ is thus different for the ground state ν = 1 and the first excited state ν = 2, as sketched in figure 1 (b) (the tip spin is ↑ in both cases). The anti-ferromagnetic coupling of the tip spin and the surface spin effectively raises the tunneling barrier for the tip electrons, for they have to overcome the additional spin-interaction energy. On the other hand, any excitation of the surface spin lowers the effective potential barrier (c. f. solid and dashed lines in figure 1 (b)). Therefore the excitation |χ 1 → |χ 2 is favored as long as the applied voltage is higher than the first excitation energy, i. e.
For a smaller voltage this excitation is not reached energetically. All higher excited states also lower the effective potential barrier, but the effect is maximal for |χ 2 . This is reflected in the enhanced excitation probability for the first excited state.
If we on the other hand consider the case of | ↓ -electrons beeing emitted from the tip, the effective potential V B (x) + χ ν , ↓ |H int |χ ν , ↓ undergoes a change of the sign in the last term.
Therefore, χ 1 , ↓ |H int |χ 1 , ↓ becomes a negative contribution and thus enhances the tunneling current, whereas χ 2 , ↓ |H int |χ 2 , ↓ acts as an additional barrier. The other consequence is the smaller excitation probability because the transition |χ 1 → |χ 2 does not reduce the tunneling barrier. Physically, this picture reflects the anti-ferromagnetic coupling again.
We remark that the existence of a direct transition |χ 1 → |χ 2 is governed by the transition matrix element χ µ , ↑ |H int |χ ν , ↑ in the surface spin space, which is connected to the anisotropy constant ε. For the rather large value of ε ≈ 5 meV, the direct excitation |χ 1 → |χ 2 is the dominant process [18].

III. SCATTERING FORMULATION
In a wave-packet scattering picture of the process depicted in figure 1 one identifies an undisturbed part H free = p 2 /(2m) + g e µ B Bτ z + H S and a scattering part H sc = V B (x) + H int . The unperturbed (asymptotic) scattering states are cast as (± stand for incoming/outgoing waves) is the spin ground state of the surface (tip). | ± k is the orbital part of the wave with a wave vector k which is determined by the direction and by the amplitude of the applied bias voltage V , i.e. k = 2meV / 2 (and e is the electron charge). By propagating |Φ + ν,τ the state of the eigenstates |ψ ± ν,τ of the full Hamiltonian can be obtained (up to a phase) [19][20][21], i.e. |ψ + ν,τ = U(0, −∞)|Φ + ν,τ and |ψ − ν,τ = U(+∞, 0)|Φ − ν,τ , where U(t, t 0 ) is the propagator in the interaction representation. Tunneling is governed by the probability amplitude to scatter from the incoming eigenstate |ψ + µ,σ into the outgoing state |ψ − ν,τ is given by we find the tunneling probability to scatter from the tip electron with a spin state |τ into |σ and from the surface spin state |χ ν into |χ µ by taking the absolute value of the S-matrix (equation (5)). For a given initial spin-spin density matrix ρ 0 (of the dimension 10): we obtain the total tunneling probability T (V ) and, thus, the normalized density matrix elements by The probability to scatter into a specific channel (ν, σ) is given as Whereas the total tunneling probability is T (V ) = ν,τ T ντ (V ) = Trρ.
In the setup figure 1 we are interested in the surface spin dynamics only. Hence, we trace out the spin |σ of the tunneling tip electron. This yields the reduced density matrix of the surface spin. Furthermore, we fix the initial spin state of the tip electron as |τ = | ↑ as determined by the polarization of the tip and |χ ν = |χ 1 , i. e. the surface spin is in the ground state when we switch on the voltage. We simplify the notations and write the resulting density matrix as ρ S µν (V). The population of the final surface spin state |χ ν as a function of the applied voltage is then given by , whereas the phase information can be extracted from the off-diagonal elements. The assumption |τ = | ↑ corresponds to the case of a full polarization of the tip. A partial polarization β can be included in the formalism via the initial spin density matrix of the tip, i.e.
For clarity however, we employ β = 1 throughout this work and comment on the change of the results for β < 1.
To calculate the individual S-matrix elements we utilize the wave packet approach [22,23].
This method exploits the fact that in the Schrödinger picture the full eigenstates can be expressed as a superposition of propagated wave packets, i.e.
where the normalization factor η ± ντ is defined by The wave vector k for η + ντ , i. e. for the incoming waves, is given by k = 2meV / 2 , whereas for the outgoing state, the wave vector may be reduced due to inelastic tunneling. In this case, we set is the energy of the tip (surface) with a spin state γ. The wave packets |φ ± ν,τ are superpositions of full eigenstates. However, by placing the wave packets far away from the interaction region U(t, t 0 ) becomes the unit operator in which case the packets by composed of the eigenstates of H free , i. e. we can write |φ ± ντ as the product state |g ± |χ ν |τ with a Gaussian function g ± . The subscripting ± now corresponds to centering the Gaussian to the left and to the right of the interaction region.
Using the spectral representation (equation 9) we find the S-matrix elements as The formula (11) is then evaluated by (i) propagating the initial product state, (ii) calculating the overlap with the final state on the right side of the barrier, and (iii) by performing a Fourier transformation.

IV. SPIN STATE POPULATION
We now proceed with computing the population P S ν of the surface spin eigenstates |χ ν as a function of the applied voltage V for different coupling constants u 0 . We assume the surface spin to be initially in the ground state and the spin of the tip's electron is being aligned with the magnetic field. Figure 3 figure 3 by vertical dashed lines.
The basic trend of the voltage-dependence of the population is to reduce the coupling energy by exciting the surface spin. As long as the energy of the tunneling electrons is high enough, the population of the ground state decreases. For voltages in the interval between the first two dashed lines, the first excited state is the only accessible tunneling channel, so that the population is solely transferred to |χ 2 . The efficiency of this transition is determined by the coupling constant u 0 .
If the voltage increases further, the next tunneling channels open and the spin interaction energy can again be lowered by exciting the surface spin. Therefore, P S 2 decreases and P S 3 rises by the same amount. The situation is similar between the last two dashed lines, where P S 3 is decreased in order to increase P S 4 . The highest state |χ 5 corresponds to a very small projection parallel to the magnetic field as the major contribution is |m = 0 (see table I). The transition probability is, thus, much smaller than a transition into the other states. For the surface spin it is therefore more convenient to maintain a higher population of |χ 4 instead of decreasing P S 4 when the voltage exceeds the last threshold.
For larger coupling constants u 0 the probability to scatter to the excited states of the surface spin increases. However, this holds true as long as u 0 < u c 0 , where u c 0 is about 7 eV. For larger u 0 , the transition efficiency decreases again. In order to reveal the dependence on u 0 , we calculated the population of the ground state as a function of V and u 0 (figure 4). Since the effects are especially pronounced for the transition |χ 1 → |χ 2 , we will only consider the corresponding voltage interval. As it turns out, the excitation of the surface spin is only probable for certain values of the coupling constant. For large u 0 , the excitation probability even approaches zero.
In order to explain this behavior, we have to apply a further simplification that renders possible an analytical treatment (see appendix). Analyzing the simplified model we infer that with C 22 = χ 2 |S z |χ 2 and q = 2mV B / 2 for the optimal value of the coupling strength. Physically, the occurrence of this optimal value can simply be interpreted as a resonance phenomenon of scattering on the narrow quantum well of the spin-spin-interaction potential. The value of u c 0 is represented by the horizontal white dashed line in figure 4.
Generally, for values for the tip polarization β < 1, the excitation probability decreases. We have already noted that the transitions induced by the tunneling electrons with the opposite spin direction (|τ = | ↓ ) are weaker in comparison to | ↑ -electrons. For β decreasing from β = 1, the strong (| ↑ ) and the weak (| ↓ ) contribution are mixed together [24]. Despite from beeing reduced to a smaller intervall of the probability, the shape of the curves in figure 3 remains the same for realistic values of β ≈ 30 %.
The resonance phenomenon arising from the quantum well of the effective potential becomes also less visible for a smaller tip polarization, as the quantum well is inverted to a barrier for electron spins of the opposite direction. The value of u c 0 is however hardly influenced.

V. RELAXATION DYNAMICS
STM experiments are not only suitable for determining the transport properties, but also offers an insight in the relaxation dynamics of the excited surface spin by using a pump-probe scheme.
The time resolution is not sufficient to resolve the precession or the transitions between neighboring |m states, but it is still high enough to measure the relaxation from |m = +2 to |m = −2 .
This corresponds in a good approximation to the transition |χ 2 → |χ 1 .
The idea of determining the relaxation dynamics is to apply a pump voltage pulse in order to excite the surface spin, i. e. to increase the population P S 2 and to detect the time evolution by a second probe pulse. The voltage of the latter has to be chosen according to the excitation threshold.
We set the maximal values V ing, there is more time for the surface spin to relax back to the ground state, so N approaches N 0 again for a large ∆t.
For a quantitative analysis, we note that all quantities are only parametrically time-depended, since the time-dependence of V pump (t) and V probe (t) is very slow on the typical time scales of the system. This means we can use the results of the steady-state considerations from the section III and IV. From the tunneling probability for the different scattering channels T νσ (V ) we compute the corresponding currents as j 1,2 = (e k 1,2 /m)T (1,2)↑ (V ) with k 1 = 2meV / 2 and k 2 = 2m(eV + E S 1 − E S 2 )/ 2 for eV > E S 2 −E S 1 and k 2 = 0 otherwise. The total current is expressed as the sum j = j 1 + j 2 proportional to the total number of tunneling electrons dN in the given time interval dt. As the calculations are based on the quasi-static regime, it is still possible to define the conductance G(V ) ∝ dj/dV , which then is parametrically time-dependent, as well. As we have pointed out in section IV, the excitation of the state |χ 2 is favored as long as the voltage is sufficiently high. For this reason, the average slope of j 2 (V ) is larger than the average slope of j 1 (V ), i. e. the conductance jumps to a higher value if eV exceeds E S 2 − E S 1 . The nonlinear conductance profile is depicted in figure 6.
We note the population of the ground state and the first excited state evolve only adiabatically with resultant voltage V pump (t) + V probe (t), because the excitation |χ 1 → |χ 2 occurs instantly on the time scale of the voltage pulses. With the surface spin driven into the excited state, it undergoes i) a precession (the characteristic time is about 2 ps) and ii) relaxes. Since the computation of N involves an integration over the time interval of several hundred nanoseconds, the oscillatory part due to the precession does not play a role. We can, therefore, describe the dynamics of the surface spin in the case of a falling voltage by using the Bloch equations including damping, but at a zero frequency. The time scale for this relaxation process, chosen according to the experiment [8], is again very slow for our system.
For a pump-probe delay ∆t < 0 such that the probe pulse still advances the pump pulse, we therefore obtain N probe = 1 eL dt P S 1 (t)j 1 (V probe (t)) is, thus, lowered in comparison to the probe pulse advancing the pump pulse. Figure 7 shows the total dependence of N on the time delay ∆t. Again the curves in figure 7 evidence a strong dependence on the coupling of the tunneling electron to the local moment u 0 .
The peak at ∆t = −200 ns (corresponding to a complete overlap of the pump and the probe pulses) is due to the nonlinearity of current-voltage characteristic, as j(V probe ). The magnitude of this effect can be explained by means of the the conductance profile ( figure 6). For small values of u 0 , the relative jump of the conductance is small, so j(V ) is nearly probe ). Hence, N only slightly exceeds N 0 when both pulses coincide. For this reason, the first peak in figure 7 is small for smaller u 0 . The maximum becomes more pronounced when approaching the optimal value of u 0 for the depopulation of the ground state. From the conductance profile (figure 6) we conclude that the relative jump of G when exceeding the threshold voltage attains the highest value for u 0 = u c 0 .
The nonlinearity then plays a dominant role, resulting in an especially strong increase of N with respect to N 0 . For larger coupling strength, the current-voltage dependence resembles again the linear case. Note that the averaged conductance is much smaller for u 0 > u c 0 than in the case u 0 < u c 0 . The minimum at ∆t = 0 is caused by the reduction of the ground state population P S 1 due to the advancing pump pulse. We have already discussed the efficiency of this process ( figure 4). For values u 0 ≪ u c 0 or u 0 ≫ u c 0 , P S 1 1 holds. Hence, N probe and for this matter N decrease only by a small amount. The depth of the minimum is maximal for u 0 = u c 0 . We remark that the relative jump in the conductance profile also decreases for β < 1. The tunneling channel connected to the excitation |χ 1 → |χ 2 is much less effective for | ↓ -electrons emitted from the tip. Hence, the major contribution for this electrons arises from the adsorbate spin remaining in the ground state. Hence, the conductance for voltages eV > E S 2 − E S 1 differs only slightly from the value left to this threshold. This effect again shows a linear dependence on the tip polarization. The total value of the conductance on the other hand increases with smaller β, which is due to the reduction of the effective potential for | ↓ -electrons with the adsorbate spin in the ground state.
Similarly to the discussion in section IV, the special features with respect to u c 0 become less visible if we reduce β. Hence, the height of maximum in figure 7 decreases. This can be explained in terms of the conduction profile which more and more resembles the linear case. On the other hand, the transition of the surface spin to the first excited state effectively blocks the | ↓ -electrons from the tips and therefore reduces the current more than in the case β = 1. We thus conclude that the minimum becomes deeper for a realistic tip polarization.

VI. CONCLUSIONS
We analyzed theoretically the dynamics of a surface-deposited localized magnetic moment, triggered by a spin-polarized current pulse from an STM tip in the presence of a magnetic field.
We showed and explained how and why the spin population of the surface-deposited structure can be manipulated by applying appropriate voltages and investigated the dependence on the coupling strength u 0 of the tunneling electrons to the localized moments. A simplified analytical model allowed us to predict a value of u 0 for which the excitation probability of the surface spin is particularly high. The obtained value is in excellent agreement with the numerical calculations.
The relaxation was incorporated by using the Bloch equations, as our model does not provide a relaxation mechanism. For future work, this point deserves a special consideration. In this respect, the main advantage of the Cu 2 N-overcoat is the magnetic decoupling of the adatom and the substrate, the excitation of magnons is strongly suppressed leading to prolonged relaxation times. In general, we suggest the spin-orbit interaction, that transfers the energy to the substrate, as the major dissipation channel.
The above treatment is based on an adiabatic approximation in the sense, that the voltage pulses durations are very long on the time scale of the electron and spin dynamics, which is a very well justified assumption for the experiments under consideration (figure 1). In future experiments one may envisage however to go beyond the relaxation dynamics with the aim to image the precession of the adsorbate spin as well. In this case picosecond pulses are required. Such pulses have already been utilized for tracing the ultrafast spin dynamics in a proof of principle experiments [25,26].
Alternatively, one may consider laser-pulse induced tunneling currents, as in the recent experiments [27,28]. Generally, the characteristic times for the tunneling process (some femtoseconds) and the spin dynamics (picoseconds) are still very different. Hence, the dynamics driven by the voltage pulses with picosecond durations can be adiabatically decomposed into (i) the time evolution of the electron states, (ii) the electron pulse driving the surface spin, and (iii) an adiabatic dependence of the current on the surface-deposited spin dynamics (which images of the spin dynamics). This time scale separation is the key for future experiments and calculations on road to full insight into the spatiotemporal evolution of the quantum spin dynamics.

VII. APPENDIX: ANALYTICAL MODEL
In section IV that for the parameters of interest the excitation probability exhibits a maximum for u 0 ≈ 7 eV. To explain the specific value of the maximum and to gain further insight in the underlying physics, we introduce here some reasonable simplifications that allow an analytically solvable model. Since the scattering channels are limited to the first two states of the surfacedeposited magnetic moment and the tip electron spin can not flip within this transition, we can reduce our considerations to a two-channel model. The operatorτ ·Ŝ can be replaced by S z as the electron spin is fixed in the z-direction. For small voltages, the gradient of the potential barrier can be neglected. We can, thus, assume a rectangular shape of the barrier. The function w(x) describing the spin-interaction area is localized around the Fe atom. Hence we approximate the coupling Hamiltonian by where x 0 is a typical length associated with w(x), L is the width of the vacuum barrier and C is the remaining coupling matrix that contains the transition amplitudes, i. e. C αβ = χ α |S z |χ β .
The wave functions ψ + ν (x) are then linear combinations of plane waves e ±ik(x−L) for x < 0 or x > L and e ±q(x−L) for x ∈ (0, L). The respective wave vector follows as k ν = 2m(eV + E S 1 − E S ν )/ 2 , whereas q ν = 2m(V B − eV − E S 1 + E S ν )/ 2 . Since the height of the vacuum barrier V B is much greater than the electron energy, we can approximately set q ν = q = 2mV B / 2 . For 0 < x < L we write ψ + ν (x) = α < ν e qx + β < ν e −qx and for x > L we define ψ + ν (x) = α > ν e ikν x + β > ν e −ikν x . The coefficients are linked by the corresponding transfer matrices. Imposing the standard matching conditions for the wave-function and its derivative we obtain for the reflection coefficient of the first excited spin state in the interval (0, L): where > and < labels the coefficients to the left and to the right of x = L. We can furthermore neglect the term k 2 /q because of the small voltages. The remaining terms are then real, so that the absolute value of β < 2 is minimal if the right side in equation (15) is minimal. Hence, α > 1 becomes small near the optimal value u c 0 . We conclude that β < 2 is small for q = 2mx 0 u 0 C 22 / 2 , yielding for the optimal value u c 0 = 2 q/(2mx 0 C 22 ).