text:
# Magnetic field effects in energy relaxation mediated by Kondo impurities

## I Motivation and Overview

Energy relaxation of “hot” electrons in disordered conductors at low enough energies was for a long time believed to be determined by direct interaction between electrons. Recent experiments on mesoscopic wires (see Ref.  and references therein) have shown that Kondo impurities (localized spins) lead to much higher energy relaxation rates than predicted by the standard theory.

On the other hand the inelastic collisions of electrons and their spin–flips are directly related to the phase coherence time probed by weak localization effects, such as low field magnetoresistance. The problem of decoherence in weak localization was recently revisited and intensively discussed.

Theoretical studies by various groups  lead to a satisfactory and consistent explanation of energy relaxation experiments by Kondo impurity mediated electron–electron interaction in the gold wires by Pierre et al.  Those gold wires were contaminated by iron impurities, and the concentration of the impurities could be independently estimated from the magnetoresistance as well as from the temperature dependence of the resistivity. As a result, it was possible to carry out a parameter free comparison of theory an experiment.

At the same time the copper samples in the first energy relaxation experiment by Pothier et al.  were fitted using the concentration and Kondo temperature of the paramagnetic impurities as free parameters. The parameters obtained from energy relaxation disagree with those obtained from the magnetoresistance experiments of the same sample.  Therefore, if magnetic impurities are responsible for these effects is dubious.

Since the behavior of magnetic impurities is sensitive to the applied magnetic field, studies of energy relaxation in presence of the magnetic field could either rule out or validate magnetic impurities as relevant scattering process in energy relaxation. Recently, Anthore et al.  reported results of such experiments in Cu wires that indicate a strong dependence of the energy relaxation on the magnetic field suggesting that magnetic impurities indeed play a role for the copper wires as well.

In this article we perform a theoretical study of transport and energy relaxation in a mesoscopic wire in dependence on an applied magnetic field. We use a diffusive Boltzmann equation to account for the static scatterers and focus in the inelastic collision integral on magnetic impurities. The findings are in qualitative agreement with the experimental data in Ref.  supporting the presumption that scattering with magnetic impurities is the essential mechanism of energy transfer at low temperatures. However, the apparent inconsistency between the values of the experimentally observed energy relaxation rate and the dephasing rate extracted from the magnetoresistance in Cu wires remains puzzling.

Starting with a brief discussion of the experiment we propose in Sec. II a simple physical picture to explain the anomalous dependence of the energy relaxation on the magnetic field. In Sec. III follows a theoretical description in terms of a renormalized Hamiltonian restricting the interaction processes to the coupling to electron–hole pairs only. We then present in Sec. IV the numerical procedure and the comparison with experimental results. Sec. V is devoted to a discussion of the interpretation and validity of the approach as well as its possible extensions.

## II Experimental situation and physical picture

Here, we briefly describe the experimental situation in Ref.  and a possible qualitative explanation of their findings.

The experimental setup in Ref.  consists of a thin copper wire of about 45 nm thickness, 105 nm width, and 5 $`\mu `$m length connected with two metallic leads. The leads are biased by an external voltage source $`U=0.1`$ mV and $`U=0.3`$ mV imposing a steady state current through the wire. The setup is placed in a dilution refrigerator with a temperature of $`25`$ mK and a magnetic field up to $`2.2`$ T is applied. The elastic mean free path can be estimated to be much smaller than the length of the wire, $`L`$, so that the transport of electrons between the contacts is diffusive. The diffusion constant $`D`$ = 90 cm2/s, estimated from the low temperature resistance, leads to a diffusion time of $`\tau _D=L^2/D=2.8`$ ns.

The aim of the experiment was to study energy distribution of electrons. The distribution function was determined by tunneling to an underlying aluminum probe electrode. The differential tunneling conductance is given by a convolution over the electron distribution functions in the wire and in the probe electrode, both for $`B=0`$ and $`B0`$, cf. with Refs.  and , respectively. In the case of vanishing magnetic field the aluminum probe electrode is in the superconducting state, and the peaked density of states of the probe electrode allows one to straightforwardly extract the energy distribution in the wire from the I-V characteristic of the junction.

For finite magnetic fields, however, the aluminum probe electrode is in the normal state. In this case the deconvolution has been made using a zero bias anomaly. Unfortunately, the latter procedure is less accurate. Considering the numerical transformation depicted in Fig. 2 of Ref. , we expect that the strongest variations due to uncertainties in the shape and the depth of the zero bias anomaly arise at the “Fermi points” $`ϵ=\pm eU/2`$. Furthermore, the tunnel probe experiment in equilibrium, $`U0`$, gave the temperature of $`65`$ mK which differs from the actual refrigerator temperature. This might be due to an oversimplification of the environmental impedance responsible for the zero bias anomaly.  Consequently, the resulting distribution function has to be taken with some care, in particular near the Fermi points.

The electron distribution function in the absence of inelastic scattering is just a linear combination of the distributions in the left and right electrode

$$f_0(ϵ,x)=(1x)f_F(ϵeU/2)+xf_F(ϵ+eU/2).$$
(1)
Here, $`f_F(ϵ)`$ is the Fermi function while $`x`$ is the longitudinal coordinate of the observation point in units of the total length, $`L`$. At low temperatures the distribution $`(\text{1})`$ has steps in both $`ϵ`$ and $`x`$ dependencies. These steps are smeared by inelastic processes, such as electron–phonon and electron–electron interaction. As a result, the distribution in the middle of the wire turns out to be almost insensitive to the bath temperature.  The smearing depends on the effective inelastic relaxation time, and the latter can be estimated from the experimentally observed distribution function.

The first experiments  for $`B=0`$ have clearly shown that the smearing is too strong to be attributed to the electron–electron or electron–phonon interaction. Therefore, in the following we do not take these interactions into account.

The qualitative outcome of the experiment  is that the behavior of the inelastic relaxation rate is a non–monotonous function of the magnetic field $`B`$. For $`BB_1=eU/4.3\mu _B`$, the relaxation rate *increases* with magnetic field and reaches a maximum at $`BB_1`$. At stronger fields, $`B>B_1`$, it decreases with further increasing magnetic field, and at $`BB_2=eU/2\mu _B`$ it reaches almost the same value as it had at $`B=0`$ and then decreases further.

The explanation of such a complex behavior given in Ref.  is based on electron scattering by magnetic impurities. For vanishing magnetic field the spin system is degenerate and only second or higher order scattering processes contribute to energy relaxation. For finite magnetic fields, there exists also a first order contribution with energy $`E_H=g\mu _BB`$ (equal to the Zeeman splitting) transferred to or from the spin system. Therefore, the energy relaxation rate increases. However, for $`E_H>eU`$ the spins are completely polarized and can no longer contribute to energy relaxation. Consequently, again only higher order processes are effective. Comparing this explanation with the experiment, one estimates the gyromagnetic factor for the impurity spins as $`g2`$.

Thus, electron–spin interaction taken into account in the lowest order of the perturbation theory explains, in principle, the main experimental features. However, from the theoretical point of view, there appears a subtlety. The problem is that the higher order terms, estimated within the framework of the $`t`$–matrix approach in Refs. , lead to a divergent contribution $`J^4/(ϵE_H)^2`$ to be integrated over. Here $`J`$ is the renormalized coupling constant which defines the strength of the electron–impurity interaction. The suggestion  to introduce a cutoff at the Korringa width , $`KJ^2`$ gives a result that is comparable with the first order contribution. Consequently, there is no systematic expansion in powers of the coupling constant, and one needs a generalized approach which is not based on expansion in terms of the interaction strength.

The aim of the present paper is to develop an approach able to treat lowest and higher order contributions within a unified scheme. We show that the problem can be described as electrons coupled to a dissipative spin system. The dissipation of impurity spins is, in turn, caused by creation/annihilation of electron–hole pairs due to electron–spin coupling. We show that this scheme includes the divergent higher–order contributions appearing in the $`t`$–matrix approach. Using our approach, we derive an electron–spin collision integral expressed through spin–spin correlation functions. An important feature of these correlation functions is that their dependences on the electron energy is automatically broadened by Korringa–type processes. These processes, however, depend on the actual electron distribution rather than on the thermal equilibrium distribution as in Ref. . As a result, the applied voltage plays the role of an effective temperature.

The crucial difference between the cases $`B=0`$ and $`B0`$ is the following. At zero magnetic field all spin–spin correlation functions are centered at zero energy and behave as $`K/(ϵ^2+K^2)`$. It is important that in the collision integral these correlators are multiplied by a combination of the electron distribution functions which at small energies is proportional to $`ϵ`$. Consequently, the broadening turns out to be unimportant , and one can omit $`K^2`$ in the denominator. This way the resulting collision integral becomes proportional to $`J^4`$ and one recovers the results  of the $`t`$–matrix approach.

At finite magnetic fields, however, the spin–spin correlation function decomposes into three contributions. The non spin–flip part is still centered at zero energy, while the two spin–flip contributions are peaked at $`ϵ=\pm E_H`$. Let us for a moment accept the above simplified form of the correlation function and assume $`K0`$. Since in this case $`K/[(ϵE_H)^2+K^2]\pi \delta (ϵE_H)`$, the two spin–flip correlation functions with finite energy transfer lead indeed to the desired first order in $`J^2`$ contribution. In this simplified case only the non spin–flip correlation function contributes to the order $`J^4`$ in the collision integral.

Using the $`t`$–matrix approach with a finite cutoff $`K`$, introduced by hand, would lead to a double counting of the first order in $`J^2`$ contribution. Here, on the other hand, the cutoff is included automatically and the first order is accounted for correctly.

In this simplified case the reasoning to explain the experimental data follows the lines of the experimentalists. The only difference is that for $`E_H>eU`$ where the spin–flip contribution is already frozen out there remains just $`1/3`$ of the $`B=0`$ energy relaxation because only the non spin–flip component contributes to energy relaxation in order $`J^4`$. In practice, $`K`$ is not constant but depends on frequency that the correlation functions are not Lorentzian shaped and therefore the spin–flip terms will also contribute to the $`J^4`$ term. Further, the width decreases with increasing magnetic field for the non spin–flip component and therefore energy relaxation mediated by magnetic impurities dies out for $`E_HeU`$.

We believe that our approach provides a consistent explanation of the magnetic field dependence of the non–equilibrium electron distribution in diffusive wires with magnetic impurities.

## III Theoretical description

Here we present a simplified version of the theory with isotropically renormalized coupling constant $`J`$ independent of energy. In general, the renormalization can be anisotropic and energy dependent. These generalizations which do not alter the underlying physics are discussed in the Appendix.

We assume that the metallic wire is in the diffusive limit, i.e. the elastic relaxation time is much smaller than other time scales. We also assume that the distribution function of electrons does not depend on the spin. The energy distribution of the electrons is governed by the Boltzmann equation

$`{\displaystyle \frac{f(ϵ,x)}{t}}{\displaystyle \frac{1}{\tau _D}}{\displaystyle \frac{^2f(ϵ,x)}{x^2}}+I\{f\}=0,`$ (2)
$`I\{f\}={\displaystyle }d\omega \{f(ϵ)[1f(ϵ\omega )]W(\omega )`$ (3)
$`[1f(ϵ)]f(ϵ\omega )W(\omega )\}.`$ (4)
Here, we include the density of states $`\rho `$ in the scattering rate $`W(\omega )`$ and omit for convenience the explicit spatial dependence of the distribution function. The rate $`W`$ describes the transitions between two electron states with energies $`ϵ`$ and $`ϵ\omega `$ mediated by coupling to the dissipative spin system. Its explicit form is given by

$$W(\omega )=(c_{\text{imp}}/\rho \mathrm{})(\rho J/2)^2C(\omega )$$
(5)
where $`c_{\text{imp}}`$ is the impurity density. Further, $`C(\omega )`$ is the Fourier transform of a spin–spin correlation function. The latter can be split as

$$C(t)=[C_+(t)+C_{}(t)]/2+C_z(t)$$
(6)
where

$$C_\pm (t)=S^\pm (t)S^{}(0),C_z(t)=S^z(t)S^z(0).$$
(7)
The averages here mean the spin and electron trace weighted with the unknown non–equilibrium density. The time evolution is governed by the Hamiltonian $`H=H_0+H_I`$ where

$$H_0=\underset{k\sigma }{}ϵ_{k\sigma }C_{k\sigma }^{}C_{k\sigma }E_HS^z$$
(8)
describes free electrons. Here, operators $`C_{k\sigma }^{}`$ and $`C_{k\sigma }`$ create and annihilate an electron in a given orbital, $`k`$, and spin, $`\sigma `$, state. $`ϵ_{k\sigma }`$ is the energy of this state. The second term in Eq. $`(\text{8})`$ describes a spin $`1/2`$ impurity with Zeeman splitting $`E_H=g\mu _BB`$. The interaction Hamiltonian

$$H_I=J\underset{kk^{}\sigma \sigma ^{}}{}𝐒𝐬_{\sigma ^{}\sigma }C_{k^{}\sigma ^{}}^{}C_{k\sigma }$$
(9)
couples electrons to the impurity spin system via the renormalized coupling strength, $`J`$, rather than the bare one, $`J_0`$. Further, the electron is coupled only to one impurity spin since we assume that the impurity density $`c_{\text{imp}}`$ is small enough to neglect higher order terms.

Using this renormalized Hamiltonian we will restrict our calculation for the time evolution of $`C(t)`$ to coupling to simple electron–hole pair excitations only, see Appendix A for details. A similar procedure was already used in Ref.  to discuss the impurity spin resonance linewidth in equilibrium using Baym and Kadanoff’s kinetic equations. Here, the spin–spin correlation functions are calculated using the conventional projection operator technique. The results read

$`C_z(\omega )`$ $`=`$ $`{\displaystyle \frac{1}{2}}{\displaystyle \frac{\nu _z(\omega )}{\omega ^2+\nu _z(\omega )^2}}`$ (10)
$`C_\pm (\omega )`$ $`=`$ $`{\displaystyle \frac{2𝒫_\pm \nu _\pm (\omega )}{[\omega E_H]^2+\nu _\pm (\omega )^2}}.`$ (11)
The functions $`\nu _z,\nu _\pm `$ describing damping of spin fluctuations can be expressed through an auxiliary function

$$\zeta (\omega )=𝑑ϵ^{}f(ϵ^{})[1f(\omega +ϵ^{})]$$
(12)
stemming from the coupling to electron–hole pairs. They read:

$`\nu _z(\omega )`$ $`=`$ $`\pi (\rho J)^2[𝒫_+\zeta (\omega E_H)+𝒫_{}\zeta (\omega +E_H)]`$ (13)
$`\nu _\pm (\omega )`$ $`=`$ $`(\pi /4)(\rho J)^2\left[2\zeta (\omega E_H)+\zeta (\omega )/𝒫_\pm \right].`$ (14)
Further, $`𝒫_\pm `$ is the occupation probability for impurity spin up or down, respectively. These probabilities are determined by a master equation

$$d𝒫_\pm /dt=\mathrm{\Gamma }_\pm 𝒫_\pm +\mathrm{\Gamma }_{}𝒫_{},𝒫_++𝒫_{}=1$$
(15)
which can be solved in the steady state leading to

$$𝒫_\pm =\mathrm{\Gamma }_{}/(\mathrm{\Gamma }_++\mathrm{\Gamma }_{}).$$
(16)
Here, the inverse lifetime for the spin up/down state, $`\mathrm{\Gamma }_\pm `$, is determined by the expression

$$\mathrm{\Gamma }_\pm =\frac{(\rho J)^2}{4\mathrm{}𝒫_\pm }𝑑\omega \zeta (\omega )C_\pm (\omega ).$$
(17)
For a given electron distribution, the set of Eqs. $`(\text{16})`$ and $`(\text{17})`$ determines the occupation probabilities of the spin system. From Eqs. $`(\text{6})`$ and $`(\text{7})`$ one can prove the sum rule for the correlation function,

$$C(t=0)=(d\omega /2\pi )C(\omega )=3/4=S(S+1),$$
(18)
which is independent of magnetic field. In the weak coupling limit

$$C_\pm (\omega )=2\pi 𝒫_\pm \delta (\omega E_H),C_z(\omega )=\pi \delta (\omega )/2.$$
(19)
Inserting this results into the rates $`(\text{5})`$ and $`(\text{17})`$, we recover the Fermi’s Golden Rule expressions for the electron collision operator for the case of interaction with a single spin. The $`\delta `$–functions here signalize energy conservation. The dissipation leads to an energy uncertainty that, in turn, results in a broadening of the $`\delta `$–functions. The final result differs from that obtained by the $`t`$–matrix approach  by allowing for a finite energy uncertainty.

For vanishing magnetic field, $`E_H=0`$, the occupation numbers $`𝒫_\pm =1/2`$, and the correlation function simply reads $`C(\omega )=3C_z(\omega )`$. Inserted in the collision integral, our expression with the renormalized coupling constant $`J`$ from Eq. (B9) is consistent with the results of Refs. . The advantage here is that the cutoff, $`\nu _z(0)=\pi (\rho J)^2\zeta (0)`$, suggested in Refs.  is naturally included. This cutoff equals to the Korringa width $`K`$. For weak coupling only elastic scattering survives and we get $`W(ϵ)=\tau _{\text{sf}}^1\delta (ϵ)`$ where the time $`\tau _{\text{sf}}`$ is usually referred to as spin–flip time . In this case the collision integral vanishes identically and one has to go beyond the lowest order.

Further, our expression for the correlation function in equilibrium and vanishing magnetic field coincides with the results of Wöger and Zittartz  based on a Nagaoka–like decoupling scheme. Using Suhl’s $`t`$–matrix for the renormalized coupling constants, the spin susceptibilities follow from proper kinetic equations. The procedure is actually similar to the one used in Ref.  for the analysis of the impurity spin resonance linewidth.

## IV Numerical procedure and comparison with Experiment

For the numerical procedure we use anisotropic and energy dependent coupling constants $`J_\pm ^z(ϵ)`$ and $`J^\pm (ϵ)`$ given in the Appendix B, Eqs. $`(\text{B9})`$ and $`(\text{B11})`$, respectively. They have to be calculated self–consistently as functionals of the final non–equilibrium distribution function. This procedure leads to a slight complication of the formulas in the previous section but does not alter their structure. The detailed changes to be made are listed in the Appendix D. Further, in Eqs. $`(\text{10})`$ and $`(\text{11})`$ we fix the lifetimes in the denominator at resonance, $`\nu _z(0)`$ and $`\nu _\pm (\pm E_H)`$, respectively, and neglect $`\omega `$–dependence of $`\nu _z`$ and $`\nu _\pm `$. This approximation formally violates the sum rule $`(\text{18})`$. However, this violation appears in higher orders in $`J`$ and thus is not worrisome. On the other hand in our approximation the Boltzmann equation $`(\text{2})`$ leads to the effective Fermi distribution in the middle of a sufficiently long wire, i.e., gives the correct “hot electron” limit.

We have started with the solution (1) of the impurity–free problem inserted in Eqs. (B9) and (B11) for the coupling constants, Eqs. (D9) and (D10) for the auxiliary functions $`\zeta _z(\omega )`$ and $`\zeta _\pm (\omega )`$, and Eq. (4) for the collision integral. Then the distribution function was evolved iteratively using the Boltzmann equation (2) with collision operator (4). At each iteration both the coupling constants and the correlation functions were updated. After about 120 iterations a stationary solution has been reached.

To make an independent comparison of the finite field data, we have fitted the impurity density with data at $`B=0`$. The resulting impurity concentration $`c_{\text{imp}}=8`$ ppm is in accordance with the experimental purity of copper of $`99.999\%`$ The density of states is chosen to be $`\rho =0.21/`$(site$``$eV).  The Kondo temperature, not known so far, has been assumed as $`T_K=0.4`$ K. Further, the gyromagnetic factor was chosen as $`g2`$, see review of the experimental results in Sec. II.

The comparison of the distribution functions for $`U=0.1`$ mV is shown in Fig. 1 and for $`U=0.3`$mV in Fig. 2. The symbols are the distribution functions for several magnetic fields obtained from experimental data by a deconvolution procedure. Note, that a magnetic field of $`B=1`$ T leads to a Zeeman splitting $`E_H0.12`$ meV. Thick lines are the outcome of our numerical procedure allowing for the non–equilibrium occupation numbers (16) for the impurity spins. For comparison, thin broken lines are the results obtained from the Boltzmann equation for equilibrium impurity spins. It is clear that the results for relatively weak magnetic fields, $`E_H<eU`$, agree with theoretical predictions only if the non–equilibrium spin population is taken into account. Consequently, we conclude that the impurity states are indeed out of equilibrium. Since the only spin relaxation mechanism taken onto account is an interaction with electrons it follows that “hot” impurity centers serve as *mediators* for the electron–electron interaction.

In large magnetic fields, $`E_H>eU`$, according to the theory, the impurity–induced energy relaxation is frozen out. Consequently, to obtain a quantitative agreement with the experimental results shown in Fig. 1 one should take into account other, though weak, scattering mechanisms.

Further, we observe that the numerical data for $`U=0.3`$mV and $`B=1.2`$T show an additional step in the middle of the energy distribution function. This can be explained in terms of lowest order, $`J^2`$, scattering with impurity spins where the transferred energy between electron and impurity spin is always $`\pm E_H`$. The origin is that the Zeeman splitting, $`E_H`$, at this magnetic field is almost $`eU/2`$. Electrons in the energy region slightly below $`eU/2`$ are scattered slightly below $`ϵ=0`$ because the distribution function differs strongly in these two regions. The same statement is correct for electrons which scatter from slightly above $`ϵ=0`$ to slightly above $`eU/2`$. Scattering between other regions is suppressed by the small difference in the distribution function at the contributing regions. Using only the first order in $`J^2`$ description for energy relaxation this feature would be much more pronounced as the one shown in Fig. 2. The finite width in the correlation functions $`(\text{11})`$ responsible for these inelastic processes smear already out most of the sharp features. We assume that additional energy relaxation processes, such as direct electron–electron interaction are responsible for the missing step in the corresponding experimental data. This feature is therefore assumed to be an artifact of our numerical calculation caused by neglecting other scattering mechanisms.

Since in our theoretical approach additional scattering mechanisms which lead to small energy transfer are missing and further the experimental data uncertainties are most pronounced at the “Fermi points” $`ϵ=\pm eU/2`$ it is not reasonable to consider these energy regions in more detail. To characterize the quality of our fits we focus the comparison between our numerical data and the experiment to the “plateau” region around $`ϵ0`$ in the middle of the two–step distribution function. The negative slope of the distribution function is a good indicator of the energy relaxation strength as long as it is small. For vanishing energy relaxation it is zero and it increases with increasing energy relaxation. In Fig. 3 we show the averaged negative slope of the energy distribution functions at the plateau near $`ϵ0`$ in dependence of the magnetic field $`B`$. We find that the theoretical data (solid lines) meet the qualitative outcome in Sec. II of increasing energy relaxation up to the maximum value at about $`E_HeU/2`$ and then decreasing again showing at $`E_HeU`$ about the same value as at $`E_H=0`$. The data, however, lack to explain the experimental outcome (dashed lines) in quantitative detail, which is not further surprising taking into account the simplifications we made.

## V Discussion

Let us first discuss the outcome of our findings, then focus on the justification of various assumptions not addressed in previous sections, and conclude with an outlook to further work.

We have succeeded to explain the main features of the inelastic relaxation due to impurity–mediated electron–electron interaction – effective relaxation in weak fields and its decrease in large magnetic fields, $`E_H>eU`$. Further, we found that the impurity spins are *out of equilibrium*. Namely, for $`E_H<eU`$ both spin states are occupied even at zero temperature, and the role of the effective temperature is played by the applied voltage. In this sense the impurity centers turn out to be “hot”. Contrarily, when the magnetic field is strong, $`E_H>eU`$, only the lowest spin level should be occupied, the occupation of the “unfavorite” spin state is determined by high–order processes. Quantitative agreement in this region, however, lacks for small energy transfer, responsible for the smearing at the “Fermi points” $`ϵ=\pm eU/2`$, where the experimentally determined distribution functions show much stronger energy relaxation than the one expected from the theory. It is so far not clear if the small energy transfer discrepancy comes from a failure of our approach in this regime, from other scattering mechanisms, or is an artifact of the deconvolution procedure which was not ruled out in Ref. . In the present work we interpret the results using a spin $`1/2`$ model which, as we believe, is valid at least when the renormalized coupling constant is small, $`\rho J1`$. This regime is usually referred to as “above the Kondo temperature” and its meaning in non–equilibrium situations will be explained in the discussion below.

For vanishing magnetic field the experiments by Pothier et al.  and Pierre et al.  have shown that the distribution function depends only on the ratio $`ϵ/eU`$ rather than on $`ϵ/k_BT`$. This observation is compatible with the fact that the Korringa inverse time, $`K`$, is proportional to $`eU`$ at $`eUk_BT`$. Moreover, using the results in Sec. III we would even get a distribution function which depends on the pair of dimensionless energies, $`ϵ/eU`$, and, $`E_H/eU`$, only which qualitatively meets the experimental findings. This relation, however, relies on the assumption that the coupling constant is isotropic and energy independent.

### A Inelastic relaxation rate

Since the distribution functions both of electrons and spins are out of equilibrium, the inelastic relaxation cannot, in general, be discussed in terms of a single relaxation rate. In the following we use two quantities to describe the inelastic relaxation.

First, we consider the collision integral $`(\text{4})`$ in relaxation time approximation. We find the corresponding rate to be

$`{\displaystyle \frac{1}{\tau _{\text{rt}}}}`$ $``$ $`{\displaystyle 𝑑ϵ^{}W(ϵ^{})[1f(ϵϵ^{})+f(ϵ+ϵ^{})]}`$ (20)
$`=`$ $`{\displaystyle \frac{1}{\tau _{\text{sf}}}}{\displaystyle 𝑑ϵ^{}f(ϵϵ^{})[W(ϵ^{})W(ϵ^{})]},`$ (21)
with the spin–flip time defined by the exact relation following from the sum rule $`(\text{18})`$ for the correlation function, $`C(\omega )`$,

$$𝑑ϵW(ϵ)=\frac{\pi }{2\mathrm{}}\frac{c_{\text{imp}}}{\rho }S(S+1)\left(\rho J\right)^2\frac{1}{\tau _{\text{sf}}}.$$
(22)
This relaxation time equals the imaginary part of the electron self energy considering only interaction with localized impurity spins. To discuss this expression at finite magnetic fields we approximate $`W(ϵ)`$ by its weak coupling form to get

$$\frac{1}{\tau _{\text{rt}}}=\frac{1}{\tau _{\text{sf}}}\left\{1\frac{2}{3}\left[f(ϵE_H)f(ϵ+E_H)\right](𝒫_+𝒫_{})\right\}$$
(23)
where $`(𝒫_+𝒫_{})`$ is the mean polarization of the localized spin, see Eqs. $`(\text{15})`$ and $`(\text{16})`$. For monotonous distribution functions the additional term lead to a non–positive contribution and we have always $`1/\tau _{\text{rt}}1/\tau _{\text{sf}}`$. At small magnetic fields the combination of occupation numbers and distribution functions in $`(\text{23})`$ decreases and becomes of order $`(E_H/eU)^2`$. In equilibrium we find this contribution to be $`\mathrm{tanh}^2(\beta E_H/2)`$ that leads to an exact cancellation of the spin–flip term, $`1/\tau _{\text{rt}}=1/3\tau _{\text{sf}}`$, for large magnetic fields, $`E_Hk_BT`$. This is in accordance to the fact that for large magnetic fields the spin–flip contribution is completely frozen out. Out of equilibrium, the contribution stemming from the distribution functions behaves similarly and cannot explain the increase of energy relaxation for small magnetic fields.

Whereas in the collision operator $`(\text{4})`$ small energy transfer is canceled by the distribution functions the relaxation rate in Eq. $`(\text{21})`$ includes terms stemming from elastic scattering for $`ϵ^{}=0`$. To consider the pure inelastic part we define an inelastic scattering rate by

$$\frac{1}{\tau _{\text{in}}}\frac{1}{\tau _{\text{rt}}}\frac{1}{\tau _{\text{el}}}=\frac{1}{\tau _{\text{rt}}}_K^{}^K^{}𝑑ϵ^{}W(ϵ^{})$$
(24)
with some given cutoff, $`K^{}`$, not further specified.

To understand the behavior of the inelastic relaxation time at small magnetic field one has to recall that the scattering processes at $`B=0`$ are almost elastic, a typical energy transfer being smaller or of the order of the Korringa rate, $`K`$. Thus, subtraction of the elastic processes leads to a drastic decrease of the inelastic relaxation rate, or an increase of the inelastic relaxation time. For finite, small magnetic fields, however, the inelastic relaxation is of order $`K+E_H`$ and therefore will lead to an increase of order $`(E_H/K)^2`$. This term has to be compared with the monotonous decrease of order $`(E_H/eU)^2`$ and dominates in the regime $`\rho J1`$ explaining the increase of the inelastic relaxation rate for small magnetic fields.

A second way to discuss inelastic relaxation is to introduce the so–called energy relaxation time defined as an average energy transfer rate. It is defined as

$$\frac{1}{\tau _ϵ}\frac{1}{\overline{ϵ}}𝑑ϵ\rho (ϵ)ϵI\{f\}.$$
(25)
Here $`\rho (ϵ)`$ is the electron density of states, while $`\overline{ϵ}`$ is a constant reference energy to normalize the particle energy. A natural scale for $`\overline{ϵ}`$ at $`eUk_BT`$ is $`eU`$. Since $`eUD`$, where $`D`$ is the electronic bandwidth, the electron density of states can be regarded as energy–independent, and after some algebra we get

$$\frac{1}{\tau _ϵ}=\frac{\rho }{\overline{ϵ}}𝑑ϵ𝑑\omega W(\omega )\omega f(ϵ)[1f(ϵ\omega )].$$
(26)
One notices that this expression does not contain elastic scattering.

In Fig. 4 we show the magnetic field dependence of the inverse energy relaxation time $`1/\tau _ϵ`$ in arbitrary units for $`U=0.3`$mV (solid line). We observe a monotonous decrease of the total energy transferred with increasing magnetic field. Considering Eq. $`(\text{26})`$ as the average $`(\rho /\overline{ϵ})𝑑ϵ𝑑\omega \mathrm{}`$ of the energy transfer, one finds in the integrand, depending on the energy, both positive and negative contributions. Although the average $`(\text{26})`$ monotonously decreases when $`B`$ increases, the positive and negative contributions in fact increase for small magnetic fields $`E_HeU`$. For comparison we show by the dashed line the average $`(\text{26})`$ where all transferred energies $`\omega `$ are added positively. This quantity, in contrast, show the non–monotonous behavior in dependence of magnetic field. The physical relevance of this quantity, however, is not clear. This observation shows that the notion of a single energy relaxation time to describe energy relaxation is quite misleading. What happens is that non–equilibrium spins just *redistribute* the electron energy between different energy regions, the extra energy being transferred due to diffusion of “hot” electrons in real space.

In the experimental papers Refs. , as well as in the theoretical ones based upon the $`t`$–matrix approach, Refs. , a contribution of two–particle scattering processes was also discussed. The results are expressed in terms of an effective interaction strength $`\gamma `$ which in our notation reads $`\gamma =(\rho J)^2/\tau _{\text{sf}}`$. In the experimental papers the results for $`B=0`$ and for finite magnetic fields were fitted assuming the same shape of the electron–electron interaction kernel $`1/\omega ^2`$. According to our findings, the behavior of the interaction kernel differs from $`\omega ^2`$. Consequently, it is difficult to compare their and our results.

### B Parameters obtained from experiment

Now let us briefly discuss the parameters obtained from the fit of the experimental data. Similarly to Ref.  where the first experiments with Cu were fitted, we determine the Kondo temperature and the impurity density (see Sec. IV) from the distribution function at $`B=0`$. These values yield $`T_K0.4`$ K and $`c_{\text{imp}}=8`$ ppm, respectively. The procedure proposed in Ref. , based on the analysis of the saturated phase breaking time, would suggest a Kondo temperature below $`T0.1`$ K, see also Ref. . However, assuming such a small $`T_K`$ in our formalism, we would disagree with the experimentally observed distribution function at $`B=0`$ in Ref. . Here, we want to mention that we used the Kondo temperature in its simplest form based on the leading logarithmic approximation. Corrections may still be large for this quantity and this disagreement could possibly be resolved using a refined version.

The procedure based on the estimate of the decoherence time, $`\tau _\varphi `$, leads to some other inconsistencies. Namely, an estimation of the phase coherence time with our parameters would lead to $`\tau _\varphi \tau _{\text{sf}}/20.06`$ ns while, according to Ref. , $`\tau _\varphi 1`$ ns. In return, a fit of the decoherence time in Ref.  using the same $`\tau _{\text{sf}}`$ has lead to an impurity density of $`c_{\text{imp}}0.15`$ ppm much lower than our estimated impurity density. This means that using the standard theory of weak localization  the impurity spins should produce much higher decoherence rates than observed experimentally. The main problem in comparing these two quantities may be the following – the decoherence time is measured in equilibrium and it saturates in the case of Cu below the Kondo temperature if we assume that both $`\tau _\varphi `$ and the energy relaxation are determined by the Kondo impurities. In that region, however, our treatment, as well as the standard theory of weak localization using the spins as random classical objects, become inaccurate. It is the aim of future work to develop a theory valid below $`T_K`$.

In our treatment we used some assumptions which were not addressed so far. First, we used the same distribution function for the electrons with different spin projections and the standard diffusive form of the Boltzmann equation. The latter assumption can be justified because the elastic relaxation time $`\tau _{\text{imp}}0.01`$ ps is much smaller than that for inelastic scattering. The usage of a single distribution function is certainly fulfilled when the spin–orbit relaxation rate, $`\tau _{\text{so}}^1`$, is larger than the rate $`\tau _{\text{sf}}^1`$ of processes which tend to violate the electron spin symmetry. This assumption is supported by the results of the experiment  yielding $`\tau _{\text{so}}39`$ ps whereas $`\tau _{\text{sf}}`$ is of the order of few nanoseconds.

We have disregarded some scattering mechanisms to focus on the magnetic dependence and to discuss the influence of this sort of impurities only. Direct electron–electron interaction leads to a pronounced smearing of the distribution function at low energies. That can be important at rather large magnetic fields when the energy relaxation mediated by magnetic impurities vanishes. Since this additional effect at $`B=0`$ is much weaker than the contribution by magnetic impurities it cannot explain the deviations from experiment in the region of small energy transfer.

In the Appendix B we determine the renormalized coupling constant $`J`$ which replaces the bare coupling constant $`J_0`$. The calculations are performed within a RPA–like approximation. All higher order corrections are supposed to be taken into account by the renormalization of the coupling constant $`J`$. This procedure, as known, leads to the so–called overcounting problem which has been considered both for the electron  and pseudo Fermion  self energy. In these papers an attempt to cure this problem by setting the energy variable of the renormalized quantity to some special value has been made. No general method that would allow to avoid the overcounting has been suggested so far. The only message is that as long as the renormalization is weak and the renormalized vertices are changing slowly with energy one can use the “double dressed” vertices instead of the “single dressed” ones, using the language of Ref. . Here, we simply adopt this “rule”, however, having in mind that if the renormalized quantities get too much structure, i.e. in the Kondo regime, this approach may become inaccurate. Our renormalized vertices are similar to those obtained using the Hamann approximation  in the $`t`$–matrix approach by Keiter , which is a generalization of Suhl’s dispersion approach  to finite magnetic fields. We further neglected the non spin–flip contribution which may be put into the free Hamiltonian and does not influence energy relaxation.

In the calculation of the correlation functions (10) and (11) we neglected the Knight shift as well as the frequency and Zeeman energy renormalizations. The Knight shift just leads to a simple addition to the “not really known” Zeeman splitting of the impurity spin and can be disregarded. The other terms can be determined by the Kramers–Kronig relation from the Korringa widths (13) and (14). They vanish for $`ϵ,E_H=0`$ and give just a renormalization of higher order in $`(\rho J)^2`$ which is usually disregarded well above the Kondo temperature. However, below the Kondo temperature these corrections would lead to an anomalous behavior and one has to go beyond the lowest order in the memory function even if one uses already renormalized coupling constants, see Ref.  for a discussion.

### C Kondo Temperature

The regime where all our calculations are actually valid is defined here as the regime where the lowest order in the logarithmic corrections, such as in Eqs. $`(\text{B6})`$ and $`(\text{B7})`$, are still small compared to $`1`$. This is what we actually mean when speaking about a “non–equilibrium” regime or about a situation “above the Kondo temperature”. The Kondo temperature is usually defined by the electron temperature where perturbation theory starts to fail. The terms of the perturbative expansion explicitly depend on the electron distribution function. As a result even when the bath temperature is reduced below $`T_K`$ the perturbative approach can still be applied provided that the deviation from equilibrium is strong enough.

In Refs. , , and  the use of an effective spatially dependent Kondo temperature $`T_K^{}(x)=T_K^{1/x}/(eU)^{(1x)/x}`$ was proposed for electrons energetically near the Fermi point $`ϵ=eU/2`$ using the free solution $`(\text{1})`$ as a basis of the renormalization. For electrons near $`ϵ=eU/2`$ one has simply to replace $`x1x`$. In these terms the effective Kondo temperature in the non–equilibrium situation is always smaller than the equilibrium Kondo temperature.

For an arbitrary distribution function it is not as straightforward to calculate the renormalized quantities as it is for the free solution $`(\text{1})`$. Equivalently to the notion of a renormalized Kondo temperature we propose to use a local temperature $`T^{}(x)`$ describing the non–equilibrium situation. The Kondo temperature then equals its bulk equilibrium value and the local temperature is a spatially dependent functional of the distribution function. For strong electron–electron interaction or equivalently in the middle of a very long wire $`T^{}(x)`$ equals the analytically determined “hot electron” temperature  independent of the scattering mechanism. Further, this temperature fulfills the boundary conditions $`T^{}(x=0,1)=T_{\mathrm{bath}}`$. Therefore, at the ends of the wire and for temperatures explored experimentally we are below Kondo temperature and our approach is inadequate. Here, we assume that the boundary regions where our approach fails are narrow and do not influence the distribution function in the middle of the wire. We believe that this assumption is fulfilled for long wires that were studied in the experiments.  Further, it is obvious that our approach becomes more accurate for larger applied voltages when the distribution function at $`B=0`$ is a function of $`ϵ/eU`$ only. This is the reason why we focus more on $`U=0.3`$ mV. It turns out (see below for details) that for $`U=0.1`$ mV and smaller voltages our theoretical approach reaches its limit of validity within the experimentally explored temperature range.

To be more specific we give an analytic expression for the effective temperature in the case of weak smearing. The distribution function then equals the free solution $`(\text{1})`$ and the smearing is just included in an effective temperature $`T_{\mathrm{eff}}`$ characterizing the energetic width of the smeared steps. (Here, the Boltzmann constant $`k_B`$ is set to one.) For simplicity we restrict ourselves to $`B=0`$ but the generalization is straightforward. As a basis of the renormalization we use Eqs. $`(\text{B5})`$ and $`(\text{B7})`$. Following the usual poor man’s scaling procedure we get

$$\rho J(ϵ)=\frac{\rho J_0}{1\rho J_0g(ϵ)}$$
(27)
whereby we adiabatically integrated out all high energy contributions from $`D`$ to $`ϵ\stackrel{~}{D}`$ and $`D`$ to $`ϵ+\stackrel{~}{D}`$ letting $`\stackrel{~}{D}`$ tend to zero. This procedure at finite temperature or finite effective temperature leads to an effective lower energy cutoff which is either $`ϵ`$, $`eU/2`$, or $`T_{\mathrm{eff}}`$ and follows immediately from

$`g(ϵ)`$ $``$ $`(1x)\mathrm{ln}\left[{\displaystyle \frac{D}{|ϵeU/2|+T_{\mathrm{eff}}}}\right]`$ (28)
$`+`$ $`x\mathrm{ln}\left[{\displaystyle \frac{D}{|ϵ+eU/2|+T_{\mathrm{eff}}}}\right].`$ (29)
Inserted into $`(\text{27})`$ we get for the renormalized coupling constant

$$\rho J(ϵ)=\mathrm{ln}\left[\frac{(|ϵeU/2|+T_{\mathrm{eff}})^{1x}(|ϵ+eU/2|+T_{\mathrm{eff}})^x}{T_K}\right]^1$$
(30)
with

$$T_K=De^{1/\rho J_0}.$$
(31)
For sufficiently weak smearing $`T_{\mathrm{eff}}eU`$ mostly electrons near the Fermi points $`ϵ\pm eU/2`$ contribute. We may define a local temperature for electrons depending on the Fermi point. For electrons at $`ϵeU/2`$ we get $`T_+^{}(x)=T_{\mathrm{eff}}^{1x}(eU)^x`$. The other local temperature for electrons at $`ϵeU/2`$ reads $`T_{}^{}(x)=T_{\mathrm{eff}}^x(eU)^{1x}`$. Here, of course, the effective smearing $`T_{\mathrm{eff}}`$ depends on $`x`$ and coincides at the ends of the wire $`(x=0,1)`$ with the bath temperature $`T_{\mathrm{bath}}`$. Since electrons of both Fermi points contribute to inelastic processes one has to average the effective temperatures. The simplest average which yields the correct expression in the middle and the ends of the wire reads $`T^{}(x)=(1x)T_{\mathrm{eff}}^{1x}(eU)^x+xT_{\mathrm{eff}}^x(eU)^{1x}`$, because at $`x=0,1`$ only electrons at one Fermi point contribute.

In the middle of the wire we get simply the requirement $`T_{\mathrm{eff}}eUT_K^2`$ for the calculations to be valid. If the smearing becomes stronger and $`T_{\mathrm{eff}}`$ is of the order of $`eU`$ or even larger we get from Eq. $`(\text{30})`$ the requirement $`T^{}=T_{\mathrm{eff}}T_K`$.

On the other hand, we may define an effective Kondo temperature by the failure of the perturbation theory. For electrons near the lower Fermi point $`ϵeU/2`$ we write

$$\rho J(ϵ)=1/x\mathrm{ln}\left[T_{\mathrm{eff}}\frac{|eU|^{(1x)/x}}{T_K^{1/x}}\right]$$
(32)
and find the effective Kondo temperature to be $`T_K^{}(x)=T_K^{1/x}/(eU)^{(1x)/x}`$ as proposed in Ref. . The Kondo temperature for electrons near the upper Fermi point follows analogously. To our opinion the introduction of a local temperature is much more intuitive than a modification of the Kondo temperature. It seems to us that the latter concept can sometimes be misleading.

### D Conclusion

In this article we studied energy relaxation mediated by magnetic impurities above Kondo temperature. We have shown that the effective electron–electron interaction is already included in the dissipative nature of the spin system. The finite line width in non–equilibrium proportional to the applied voltage lead to strong energy relaxation even for vanishing magnetic field where the spin states are degenerate. We succeeded to derive a non–perturbative description valid for arbitrary magnetic fields where the perturbative t–matrix approach fails due to divergences in the two electron scattering rate.

We characterized the efficiency of the energy relaxation by the slope of the distribution function at $`ϵ0`$. The magnetic field dependence of this quantity shown in Fig. 3 shows a non–monotonous behavior. For small magnetic field the energy relaxation strength is increased and for larger magnetic fields it decreases again. These features meet qualitatively recent experimental data on energy relaxation in mesoscopic copper wires. Although the detailed comparison of the energy distribution functions lacks to fit quantitatively all energy regions, this gives a strong indication that magnetic impurities are indeed responsible for energy relaxation in copper wires.

Several things remain to be done in order to achieve a quantitative description for the energy relaxation experiment in the broad regions of temperature and magnetic field. To go below the Kondo temperature the overcounting problem has to be addressed. High–order corrections in impurity density might require a more complicated kinetic equation including quantum effects. Finally, for larger magnetic fields and/or low concentration of the localized spins other relaxation mechanisms have to be accounted for.

To show that magnetic impurities are indeed present in copper wires one has to describe all the available experiments using the same set of parameters. Since there are strong deviations to the magnetoresistance experiments this is still an open question.

## Acknowledgments

The authors would like to thank I. L. Aleiner, L. I. Glazman, F. Pierre, and H. Pothier for valuable discussions. Financial support was provided by the Deutsche Forschungsgemeinschaft (DFG), as well as by Princeton University, NEC Research Institute and Argonne National Laboratory. Further, this research is partly supported by the US DOE Office of Science under contract No. W-31-109-ENG-38.

## A Derivation of the collision integral

We start with the ‘bare’ Hamiltonian $`H=H_0+H_I`$ where the free Hamiltonian is given in Eq. (8). In terms of pseudo Fermions it reads:

$$H_0=\underset{𝐤\sigma }{}ϵ_{𝐤\sigma }C_{𝐤\sigma }^{}C_{𝐤\sigma }+\underset{\beta }{}E_\beta a_\beta ^{}a_\beta .$$
(A1)
Here the creation/annihilation operators $`C_{k\sigma }^{},C_{k\sigma }`$ describe electrons while the operators $`a_\beta ^{},a_\beta `$ describe pseudo Fermions. The pseudo Fermion states are specified by $`\beta =\pm `$ having the energies $`E_\pm =g\mu _BB/2`$. The interaction Hamiltonian

$$H_I=J_0\underset{\mathrm{𝐤𝐤}^{}\sigma \sigma ^{}}{}𝐒𝐬_{\sigma ^{}\sigma }C_{𝐤^{}\sigma ^{}}^{}C_{𝐤\sigma }$$

is the $`sd`$ exchange Hamiltonian with the bare coupling $`J_0`$. Expressing spin operators through pseudo Fermions one can rewrite it as

$$H_I=\underset{\mathrm{𝐤𝐤}^{}}{}J_{\sigma ^{}\sigma }^{\beta ^{}\beta }C_{𝐤^{}\sigma ^{}}^{}C_{𝐤\sigma }a_\beta ^{}^{}a_\beta $$
(A2)
(summation over repeated indices is expected). In addition, we have to take into account the operator constraint $`a_+^{}a_++a_{}^{}a_{}=1`$ which can be formally done by a projection with the help of a complex chemical potential .

Linear in the impurity density $`c_{\text{imp}}`$ the angular averaged collision integral for the classical Boltzmann equation can be expressed as

$`I\{f\}={\displaystyle \frac{i}{\mathrm{}}}\left\{f(ϵ)\mathrm{\Sigma }^>(ϵ)+[1f(ϵ)]\mathrm{\Sigma }^<(ϵ)\right\}`$ (A3)
with $`\mathrm{\Sigma }^{>/<}(ϵ)=_\sigma \mathrm{\Sigma }^{>/<}(𝐤\sigma ,ϵ)/2`$ where $`ϵ=ϵ_{𝐤\sigma }`$ is the spin averaged self energy assumed to be independent of the angular momentum. $`f(ϵ)`$ is the angular averaged distribution function for electrons of energy $`ϵ`$ where we suppressed the spatial dependence for convenience. Since the self energy is proportional to the impurity density we already replaced the electron Green’s functions by their unperturbed form and integrated over frequency to get the classical form of the Boltzmann equation.

The electron self energy in lowest nonvanishing order, depicted in Fig. 5, is given by a pseudo Fermion bubble and an electron or hole line in between

$`\mathrm{\Sigma }_2^>(𝐤\sigma ,ϵ)`$ $`=`$ $`c_{\text{imp}}J_{\sigma ^{}\sigma }^{\beta ^{}\beta }J_{\sigma \sigma ^{}}^{\beta \beta ^{}}{\displaystyle \underset{𝐤^{}}{}}{\displaystyle \frac{d\omega }{2\pi }\frac{d\omega ^{}}{2\pi }G_0^>(𝐤^{}\sigma ^{},ϵ\omega ^{})}`$ (A5)
$`𝒢_0^>(\beta ^{},\omega +\omega ^{})𝒢_0^<(\beta ,\omega )`$
where summation over internal spins is implied. The self energy $`\mathrm{\Sigma }_2^<`$ is given by the obvious change of the $`<`$ and $`>`$ signs in (A5). Since all corrections in the electron Green’s functions are of higher orders in $`c_{\text{imp}}`$ they can always be used in the free form $`G_0^<(𝐤\sigma ,\omega )=2\pi if(\omega )\delta (\omega ϵ_{𝐤\sigma })`$ and $`G_0^>(𝐤\sigma ,\omega )=2\pi i[1f(\omega )]\delta (\omega ϵ_{𝐤\sigma })`$, respectively. The pseudo Fermion Green’s function in lowest order reads $`𝒢_0^<(\beta ,\omega )=2\pi i𝒫_\beta \delta (\omega E_\beta )`$ and has to be renormalized in the following. Note that when using the Abrikosov technique , the occupation numbers $`𝒫_\pm `$ acquire an additional chemical potential $`\mathrm{exp}(i\lambda )`$ and therefore $`𝒢_0^>(\beta ,\omega )=2\pi i\delta (\omega E_\beta )`$, see Ref. . In contrast to Ref.  we assume that the occupation probabilities $`𝒫_\pm `$ are determined by the non–equilibrium electron distribution. Therefore they have to be found selfconsistently.

For further proceeding we have to classify the appearing graphs. Graphs including two pseudo fermion bubbles like shown in Fig. 6 are of second order in the impurity density and are therefore neglected. It is obvious that we may include an arbitrary number of additional electron–hole bubbles where each one is at least of order $`J_0^2`$. The two lowest order graphs with one additional electron–hole bubble are depicted in Fig. 7. Considering only these two topologically different graphs in the self energy we obtain the $`t`$–matrix expression for the collision integral derived in Refs. , however, with bare coupling constants. The electron–hole pairs do not only renormalize the pseudo Fermion propagator, which has been considered in Ref. , but also include bubbles connecting the upper and lower pseudo Fermion line. The latter graphs are so called crossed diagrams and lead to vertex corrections  not included in a simple non–crossing approach.

All other corrections lead to higher orders in $`J_0`$ for a single electron–hole bubble or outer electron line. Above the Kondo temperature it is conventionally accepted that one may write these corrections as a renormalization of the corresponding vertices and that in the leading logarithmic order the vertices are renormalized independently . As a result, two dressed vertices instead of one are used. This procedure has been explicitly proven in the leading logarithmic approximation both for the electron and pseudo Fermion  self energy. An alternative derivation of the independent vertex renormalization has been given in Ref. .

Assuming this independent renormalization of vertices to be correct also for more involved graphs and that this renormalization depends only on the electron energies, one can equivalently start with a renormalized Hamiltonian restricting the appearing graphs to simple electron–hole bubbles of the order of $`J^2`$ in the *renormalized* matrix elements $`\stackrel{~}{J}_{\sigma ^{}\sigma }^{\beta ^{}\beta }`$.

Neglecting the energy dependence of $`\stackrel{~}{J}_{\sigma ^{}\sigma }^{\beta ^{}\beta }`$ one obtains

$`\mathrm{\Sigma }^>(𝐤\sigma ,\omega )`$ $`=`$ $`c_{\text{imp}}\stackrel{~}{J}_{\sigma ^{}\sigma }^{\beta ^{}\beta }\stackrel{~}{J}_{\sigma \sigma ^{}}^{\beta \beta ^{}}{\displaystyle \underset{k^{}}{}}{\displaystyle \frac{d\omega ^{}}{2\pi }\frac{d\widehat{\omega }}{2\pi }G_0^>(𝐤^{}\sigma ^{},\omega \omega ^{})}`$ (A7)
$`𝒢_0^>(\beta ^{},\widehat{\omega }+\omega ^{})𝒢_0^<(\beta ,\widehat{\omega })_{\mathrm{eh}}`$
where the average $`\mathrm{}_{\text{eh}}`$ means the dressing of the pseudo Fermion bubble with all possible electron–hole pairs using renormalized coupling constants. Rephrasing the pseudo Fermions in terms of spin operators and assuming isotropic coupling $`J`$ we may write

$`\mathrm{\Sigma }^>(\omega )=i{\displaystyle \frac{c_{\text{imp}}}{4\rho }}(\rho J)^2{\displaystyle 𝑑\omega ^{}C(\omega ^{})[1f(\omega \omega ^{})]}`$ (A8)
with the time–dependent correlation function $`C(t)`$ given by Eqs. $`(\text{6})`$ and $`(\text{7})`$. The self energy $`\mathrm{\Sigma }^<`$ is given by the replacement $`f1f`$ and $`C(\omega )C(\omega )`$ which is obvious since $`C(t)`$ is an autocorrelation function. Inserting this expression in (A3) we find the collision integral in the desired form (4). The spin–spin correlation functions are calculated then using a projection operator technique . The actual procedure follows the lines outlined in Ref.  for the case of two–level systems.

## B Vertex renormalization

We proceed to determine the renormalized coupling constants. They can be derived in different ways: by a poor man’s scaling procedure , Suhl’s $`t`$–matrix approach , following Abrikosov and solving a generalized vertex equation , or considering lower order corrections and assuming similar resummation procedure as in equilibrium. According to Refs.  and  the renormalized vertex reads

$$\mathrm{\Gamma }(ϵ^{}\sigma ^{},\omega ^{}\beta ^{}|ϵ\sigma ,\omega \beta )=\mathrm{\Gamma }_0+\mathrm{\Gamma }_e(ϵ+\omega )+\mathrm{\Gamma }_h(ϵ\omega ^{})$$
(B1)
where

$$\mathrm{\Gamma }_0=J_{\sigma ^{}\sigma }^{\beta ^{}\beta }=J_0𝐒_{\beta ^{}\beta }𝐬_{\sigma ^{}\sigma }$$
(B2)
is the bare vertex of the $`sd`$ exchange Hamiltonian. $`\sigma ,\beta `$ are the incoming electron and pseudo Fermion spins and $`ϵ,\omega `$ are the incoming electron and pseudo Fermion energies. The primed quantities are the outgoing spins and energies, respectively. Energy conservation is fulfilled at each vertex meaning $`ϵ^{}=ϵ+\omega \omega ^{}`$. The electron and hole vertex parts can be assumed to depend on a single energy variable.  The electron vertex $`\mathrm{\Gamma }_e`$ depends on the sum of incoming electron and pseudo Fermion energies $`ϵ+\omega `$, and the hole vertex part $`\mathrm{\Gamma }_h`$ on the difference of the incoming electron and outgoing pseudo Fermion energies $`ϵ\omega ^{}`$. In the lowest order in $`J_0`$ the retarded quantities read

$`\mathrm{\Gamma }_{e/h}^r(ϵ)={\displaystyle \frac{i\rho J_{\sigma ^{}\kappa }^{\beta ^{}\gamma }J_{\kappa \sigma }^{\gamma \beta }}{2}}{\displaystyle }{\displaystyle \frac{d\xi }{2\pi }}{\displaystyle }dϵ_{q\kappa }[G_0^r(q\kappa ,ϵ\xi )𝒢_0^K(\gamma ,\xi )`$ (B3)
$`+G_0^K(q\kappa ,ϵ\xi )𝒢_0^r(\gamma ,\xi )]J^{\beta ^{}\gamma }_{\sigma ^{}\kappa }J^{\gamma \beta }_{\kappa \sigma }\rho g(ϵE_\gamma ).`$ (B4)
Here the upper sign stands for the electron vertex and the lower sign for the hole vertex, respectively. The auxiliary function

$$g(ϵ)=_D^D𝑑ϵ^{}\frac{f(ϵ^{})1/2}{ϵϵ^{}+i\delta }.$$
(B5)
leads in equilibrium to the usual logarithmic corrections.

Considering the lowest order vertex corrections to the self energy in Fig. 5 and assuming that the pseudo Fermion energies are fixed to resonance, we get

$$\mathrm{\Gamma }^\pm =\frac{J_0}{2}\left\{1+\frac{\rho J_0}{2}[g(ϵ\pm E_H)+g(ϵ)]\right\}$$
(B6)
for the spin–flip component proportional to $`S^\pm `$ and

$$\mathrm{\Gamma }_\pm ^z=(J_0/4)\left[1+\rho J_0g(ϵ\pm E_H)\right]$$
(B7)
for the component proportional to $`S^z`$ applying to a spin up/down electron. This means that the leading renormalization of the ‘up’ electron is given by the ‘down’ electrons and vice versa. This anisotropic decomposition of the renormalized coupling constants has also been found by Keiter . Setting the pseudo Fermion energy to resonance is correct only for the lowest order correction in the electron self energy in Fig. 5. However, considering corrections to the graphs in Fig. 7 we find that there the pseudo Fermion energies are weighted with a divergent term $`(\omega E_\beta )^2`$ justifying the use of the resonance energy $`\omega E_\beta `$ in the renormalization. Assuming that this is also justified for more involved graphs we use these corrections as a basis of our renormalization. To write the renormalized Hamiltonian in terms of effective anisotropic and energy dependent coupling constants we reintroduce the spin operators and define $`J_\pm ^z(ϵ)=4\mathrm{\Gamma }_\pm ^z(ϵ)`$ and $`J^\pm (ϵ)=2\mathrm{\Gamma }^\pm (ϵ)`$.

To define a ‘Kondo scale’ we consider the renormalized coupling constants as functions of the incoming electron energy and use the Hamann approximation  for the renormalized quantities. This approximation yields the same leading–logarithm expansion as the poor man’s scaling procedure based on the expressions (B6) and (B7), however, taking care of the unitarity condition in each order. In this approximation the $`S^z`$ coupling constant reads

$`J_\pm ^z(ϵ)/J_0`$ $`=`$ $`\left\{\right|1(\pi \rho J_0)^2S(S+1)/4\rho J_0g(ϵE_H)|^2`$ (B9)
$`+(\pi \rho J_0)^2S(S+1)\}^{1/2}`$
in accordance with a high temperature expansion of Keiter . Analogously the spin–flip process renormalization reads

$`J^\pm (ϵ)/J_0`$ $`=`$ $`\left\{\right|1(\pi \rho J_0)^2S(S+1)/4\rho J_0[g(ϵ)`$ (B11)
$`+g(ϵ\pm E_H)]/2|^2+(\pi \rho J_0)^2S(S+1)\}^{1/2}.`$
These quantities coincide for $`B=0`$ with those used in Refs.  and cover the results obtained by the poor man’s scaling procedure, cf. discussion about Kondo temperature in Sec. V. Note that the renormalized coupling constants remain finite also below the Kondo temperature where their meaning is questionable.

## C Calculation of the correlation function

We briefly present the method used to derive the results in Sec. III. Our main interest is the spin–spin correlation functions (7). Since $`C(t)`$ is an autocorrelator its Fourier transform is given by the expression $`C(\omega )=\stackrel{~}{C}(i\omega +\delta )+\stackrel{~}{C}(i\omega \delta )`$ where $`\stackrel{~}{C}(z)`$ the Laplace transform of the correlator. Using projection operators $`P^zX=S^zXS^z/S^zS^z`$ for $`C_z`$ and $`P^\pm X=S^\pm XS^{}/S^\pm S^{}`$ for the $`C_\pm `$ component, one can derive a formally exact integro–differential equation

$$\dot{C}_a(t)=\mathrm{\Phi }_aC_a(t)_0^t𝑑u\varphi _a(tu)C_a(u)$$
(C1)
with the solution in terms of the Laplace transform

$$\stackrel{~}{C}_a(z)=\frac{C_a(t=0)}{z\mathrm{\Phi }_a+\stackrel{~}{\varphi }_a(z)}$$
(C2)
where $`a=z,\pm `$. Here, $`\mathrm{\Phi }_z=\dot{S}^zS^z/S^zS^z=0`$ and $`\mathrm{\Phi }_\pm =\dot{S}^\pm S^{}/S^\pm S^{}=i\stackrel{~}{E}_H`$ leads to the free propagation of the correlators, where $`\stackrel{~}{E}_H`$ includes the Knight shift neglected throughout the article. The averages are calculated with the proper steady state density defined by the stationary solution of the Boltzmann equation. The memory kernel $`\varphi _a(t)`$ for the correlation function plays a similar role as the self energy for the Green’s function. For the $`C_\pm `$ term we find

$$\varphi _\pm (t)=\frac{\dot{S}_r^\pm (t)\dot{S}^{}}{S^\pm S^{}}+\mathrm{\Phi }_\pm \frac{\dot{S}_r^\pm (t)S^{}}{S^\pm S^{}}.$$
(C3)
Here, the index $`r`$ in $`S_r^\pm (t)`$ indicates that the dynamics of the spin operator is reduced by the projection. It is determined by the expression $`\dot{S}_r^\pm (t)=\mathrm{exp}[i\widehat{L}(1P^\pm )t]\dot{S}^\pm `$ with the Liouville operator $`\widehat{L}`$ acting as $`\widehat{L}\widehat{X}=[H,\widehat{X}]/\mathrm{}`$. The memory kernel for the $`C_z`$ correlation function is simply given by (C3) with the replacement $`\pm ,z`$. Though formally exact, the above equations do not allow a simple calculation of the correlators. Here, we expand the kernel up to second order in the renormalized coupling $`J`$. Since the dynamics of the expanded kernel function is oscillatory the Fourier transformed correlation function has always the simple form

$$C_a(\omega )=\frac{2C_a(t=0)\mathrm{Re}\varphi _a(\omega )}{[\omega i\mathrm{\Phi }_a+\mathrm{Im}\varphi _a(\omega )]^2+[\mathrm{Re}\varphi _a(\omega )]^2}$$
(C4)
with $`a=z,\pm `$. Similar to Green’s functions, this reflects the fact that in the steady state the retarded and advanced self energies are complex conjugate, leading to a similar structure of the “$`>`$” or “$`<`$” Green’s functions. Note that the roles of the imaginary and real parts of the memory kernel are opposite to those of the Green’s functions because of different set of definitions. Further, we define $`\mathrm{Re}\varphi _a(\omega )\mathrm{Re}\{\stackrel{~}{\varphi }_a(i\omega +\delta )\}`$ and the imaginary part $`\mathrm{Im}\varphi _a(\omega )`$ follows from a Kramers–Kronig relation. With $`C_z(t=0)=S^zS^z=1/4`$ and $`C_\pm (t=0)=S^\pm S^{}=𝒫_\pm `$ and further neglecting the imaginary parts in the denominators which lead to a frequency and Zeeman energy renormalization, we find the expressions $`(\text{10})`$ and $`(\text{11})`$ for the correlation functions.

## D Detailed theoretical description

Starting with anisotropic and energy dependent coupling constants derived in Appendix B the interaction Hamiltonian in terms of impurity spin operators reads

$`H_I`$ $`=`$ $`{\displaystyle \frac{1}{2}}{\displaystyle \underset{kk^{}}{}}\{S^+J^+(ϵ_k)C_k^{}^{}C_k+S^{}J^{}(ϵ_k)C_k^{}^{}C_k`$ (D2)
$`+S^z[J_+^z(ϵ_k)C_k^{}^{}C_kJ_{}^z(ϵ_k)C_k^{}^{}C_k]\}.`$
Considering the electron self energy in the Boltzmann equation the scattering rate $`W`$ in the collision integral $`(\text{4})`$ depends on both, incoming electron energy $`ϵ`$ and transferred energy $`\omega `$

$`W(ϵ,\omega )`$ $`=`$ $`(c_{\text{imp}}\rho /8\mathrm{})\{J^{}(ϵ)J^+(ϵ^{})C_+(\omega )`$ (D5)
$`+J^+(ϵ)J^{}(ϵ^{})C_{}(\omega )`$
$`+[J_+^z(ϵ)J_+^z(ϵ^{})+J_{}^z(ϵ)J_{}^z(ϵ^{})]C_z(\omega )\}`$
with the outgoing electron energy $`ϵ^{}=ϵ\omega `$, cf. Appendix A and Eqs. $`(\text{6})`$ and $`(\text{7})`$. The calculation of the correlation functions follows the lines in Appendix C and the general form $`(\text{C4})`$ leading the Eqs. $`(\text{10})`$ and $`(\text{11})`$ remain the same. The changes are only in the damping

$`\nu _z(\omega )`$ $`=`$ $`\pi \rho ^2[𝒫_+\zeta _+(\omega E_H)+𝒫_{}\zeta _{}(\omega +E_H)]`$ (D6)
$`\nu _\pm (\omega )`$ $`=`$ $`(\pi /4)\rho ^2\left[\zeta _z(\omega E_H)+\zeta _{}(\omega )/𝒫_\pm \right]`$ (D7)
with the auxiliary functions

$`\zeta _z(\omega )`$ $`=`$ $`{\displaystyle 𝑑ϵ^{}\left[J_+^z(ϵ^{})J_+^z(\omega +ϵ^{})+J_{}^z(ϵ^{})J_{}^z(\omega +ϵ^{})\right]}`$ (D9)
$`\times f(ϵ^{})[1f(\omega +ϵ^{})]`$
$`\zeta _\pm (\omega )`$ $`=`$ $`{\displaystyle 𝑑ϵ^{}J^{}(ϵ^{})J^\pm (\omega +ϵ^{})f(ϵ^{})[1f(\omega +ϵ^{})]}.`$ (D10)
The corresponding master equation $`(\text{15})`$ for the occupation probabilities $`𝒫_\pm `$ remains the same but the rate changes to

$$\mathrm{\Gamma }_\pm =\frac{\rho ^2}{4\mathrm{}𝒫_\pm }𝑑\omega \zeta _\pm (\omega )C_\pm (\omega ).$$
(D11)
We used these formulae for the numerical determination of the distribution function in Sec. IV, however, the main features do not alter even if one uses the simplified description in Sec. III.


format: text

text:
# Continuous Unitary Transformations

## I Introduction

Unitary transformations play an extremely important role in physics. For example, with the help of unitary transformation sometimes it is possible to simplify a problem of interacting quasiparticles by eliminating “inconvenient” interaction terms from initial Hamiltonian and construct corresponding effective interaction terms.

In the present time we observe a growing number of publications (see, ,,,,,, ,,,,,, ,), where the, so-called, flow equations (FE) method is successfully used to diagonalize Hamiltonians by means of an appropriate unitary transformation. This method of continuous unitary transformations looks very convenient and universal. The transformed Hamiltonian appears as a result of solving differential equations on some formal parameter.

It should be noted that quite the same idea and an analogous method (Safonov, 1982 , 1983 ) to construct nonlinear unitary transformation was published far before the original paper on FE method (Wegner, 1994 ) and just did not attract enough attention of researches. For convenience, we shall call the approach of Refs. , as a one-step continuous unitary transformation (OS CUT) method.

The goal of the present paper is to discuss and compare FE and OS CUT methods. As it will be shown later, the FE method is a generalization of OS CUT approach. The only difference between the methods is that the unitary transformation generator is assumed to be dependent on a formal parameter in the FE method and the generator does not depend on the parameter in the OS CUT method. In the CUT method we obtain linear differential equations. In the FE method we should solve a nonlinear system of differential equations. Below we shall demonstrate both methods by several examples.

It is interesting that the OS CUT method idea can be applied to the case of quantum equations of motion. Recently Mišta and Filip  have developed a method of non-perturbative solutions of nonlinear Heisenberg equations. We shall briefly discuss this approach at the end of this paper.

## II Continuous Unitary Transformation method

### A One step continuous unitary transformation

Here we shall formulate the key idea of the OS CUT method , . A general form of unitary transformation can be written as

$$H(\theta )=e^{\theta R}He^{\theta R}$$
(1)
where$`\theta `$ is a formal parameter, $`R`$is an anti-Hermitian operator ($`R^{}=R`$), $`H`$ and $`H(\theta )`$ are the initial and transformed Hamiltonians, respectively. The expression (1) is the solution of the equation

$$\frac{d}{d\theta }H(\theta )=[R,H(\theta )]$$
(2)
with the initial condition$`H(0)=H`$. In this form the Hamiltonian is “rotated in an operator space around $`R`$.”

In order to solve this equation we should write the most general form of $`H(\theta )`$ as the expansion in terms of operator combinations with unknown $`\theta `$-dependent coefficients. The most general form of $`R`$ as anti-Hermitian operator combinations also should be used.

After substituting these general (for $`H(\theta )`$ and $`R`$) expressions into (2), one can obtain a set of linear differential equations by comparing coefficients in analogous operator compositions in both sides. Solving these equations with the initial conditions, we obtain the transformed Hamiltonian $`H(\theta )`$. In order to eliminate ‘inconvenient’ terms one needs to put their coefficients (for example, for $`\theta =1`$) equal to zero. This condition defines the choice of $`R`$.

Usually it is simple to find the anti-Hermitian operator, which eliminates given non-diagonal term. The OS CUT method first has been applied for the spin Hamiltonian diagonalization ,. Then it was successfully used in physics of nuclear spin waves , magnon-impurity interactions , theory of superconductivity  and for eliminating three-boson interactions .

### B Flow equations for Hamiltonians

The flow equations method  begins with a general form of a unitary transformation:

$$H(\mathrm{})=U(\mathrm{})HU^{}(\mathrm{}),$$
(3)
where $`\mathrm{}`$ is a continuous flow parameter. The Hamiltonian $`H(\mathrm{})`$ is transformed from the initial Hamiltonian $`H(0)=H`$. Differentiation (3) yields

$$\frac{d}{d\mathrm{}}H(\mathrm{})=[\eta (\mathrm{}),H(\mathrm{})]$$
(4)
with the generator

$$\eta (\mathrm{})=\frac{dU(\mathrm{})}{d\mathrm{}}U^{}(\mathrm{})=\eta ^{}(\mathrm{}).$$
(5)
The generator $`\eta (\mathrm{})`$ should be chosen in such a way, that the off-diagonal matrix elements decay. A good choice is defined by

$$\eta (\mathrm{})=[H_d(\mathrm{}),H(\mathrm{})],$$
(6)
where $`H_d(\mathrm{})`$ is the diagonal part of $`H(\mathrm{})`$ (comments ***A rigorous criterion that this choice is most optimal seems has not been obtained.,The condition of diagonality and the condition that $`\eta (\mathrm{})=[H_d(\mathrm{}),H(\mathrm{})]=0`$ are not equivalent to each other.).

One can see that the FE method utilizes the same key tool (first proposed in the OS CUT method): the unitary transformation is represented as a solution of differential equations on a formal parameter. With an accuracy of notations the only difference between FE procedure and the OS CUT method is the parameter dependence of the generator $`\eta (\mathrm{})`$. If $`\eta `$ does not depend on $`\mathrm{}`$ we have $`U(\mathrm{})=\mathrm{exp}(\eta \mathrm{})`$. This case, obviously, exactly corresponds to the OS CUT method. If $`U(\mathrm{})=\mathrm{exp}(RF(\mathrm{}))`$, where $`R`$ is the parameter independent anti-Hermitian operator and $`F(\mathrm{})`$ some function on $`\mathrm{}`$, then $`\eta (\mathrm{})=(dF(\mathrm{})/d\mathrm{})R`$. The Eq.(4) in this case can be written as

$$\frac{d}{(dF(\mathrm{})/d\mathrm{})d\mathrm{}}H(\mathrm{})=[R,H(\mathrm{})].$$
(7)
This equation is reduced to Eq.(2) using a simple change of variable $`d\theta =(dF(\mathrm{})/d\mathrm{})d\mathrm{}`$. For $`U(\mathrm{})=\mathrm{exp}\left(R_1F_1(\mathrm{})+R_2F_2(\mathrm{})\right)`$ with commuting $`R_1`$ and $`R_2`$ anti-Hermitian operators one has $`\eta (\mathrm{})=(dF_1(\mathrm{})/d\mathrm{})R_1+(dF_2(\mathrm{})/d\mathrm{})R_2`$. The diagonalization procedure in this case can be considered as a two transformations with $`U_1(\mathrm{})=\mathrm{exp}\left(R_1F_1(\mathrm{})\right)`$ and $`U_2(\mathrm{})=\mathrm{exp}\left(R_2F_2(\mathrm{})\right)`$ unitary operators. Each transformation is described by (7) with $`R_1`$, $`F_1(\mathrm{})`$, or $`R_2`$, $`F_2(\mathrm{})`$.

In the most general case the $`\mathrm{}`$-dependent generator transformation can not be reduced to the parameter independent case. In such a situation the FE method can give principally different results than the OS CUT method. Let us consider several examples.

## III Quadratic form diagonalization

Stein  has considered the following quadratic Hamiltonian:

$$H_1(\mathrm{})=\underset{𝐪}{}\left(f_𝐪(\mathrm{})(a_𝐪^{}a_𝐪+b_𝐪^{}b_𝐪)+g_𝐪(\mathrm{})(a_𝐪^{}b_𝐪^{}+a_𝐪b_𝐪)\right),$$
(8)
where $`a_𝐪^{}`$, $`a_𝐪`$ and $`b_𝐪^{}`$, $`b_𝐪`$ are the creation and annihilation Bose operators. Applying flow equation

$$\frac{d}{d\mathrm{}}H_1(\mathrm{})=[\eta _0(\mathrm{}),H_1(\mathrm{})]$$
(9)
with

$$\eta _0(\mathrm{})=\frac{1}{2}\underset{𝐪}{}g_𝐪(\mathrm{})(a_𝐪^{}b_𝐪^{}a_𝐪b_𝐪),$$
(10)
the following nonlinear differential equations have been obtained:

$$\frac{d}{d\mathrm{}}f_𝐪(\mathrm{})=g_𝐪^2(\mathrm{}),\frac{d}{d\mathrm{}}g_𝐪(\mathrm{})=f_𝐪(\mathrm{})g_𝐪(\mathrm{}).$$
(11)
The solutions of (11) are

$$f_𝐪(\mathrm{})=ϵ_𝐪\mathrm{coth}(\mathrm{}ϵ_𝐪+\mathrm{}_0(𝐪)),g_𝐪(\mathrm{})=\frac{ϵ_𝐪\mathrm{sgn}(\gamma _𝐪)}{\mathrm{sinh}(\mathrm{}ϵ_𝐪+\mathrm{}_0(𝐪))},$$
(12)
where $`ϵ_𝐪=(1\gamma _𝐪^2)^{1/2}`$, $`\mathrm{}_0(𝐪)=\frac{1}{2}\mathrm{ln}\left(\frac{1+ϵ_𝐪}{1ϵ_𝐪}\right)`$ and

$$ϵ_𝐪=\sqrt{f_𝐪^2(\mathrm{})g_𝐪^2(\mathrm{})}=\sqrt{f_𝐪^2(0)g_𝐪^2(0)}$$
(13)
is the spectrum of the diagonal (at $`\mathrm{}\mathrm{}`$) Hamiltonian (8).

Let us now consider the analogous diagonalization procedure in the framework of the OS CUT method. For convenience, the same notations will be used. Simple analysis shows that the only anti-Hermitian operator that can eliminate $`a_𝐪^{}b_𝐪^{}+a_𝐪b_𝐪`$ in (8) is $`a_𝐪^{}b_𝐪^{}a_𝐪b_𝐪`$. Thus, the generator (independent on $`\mathrm{}`$) should be taken in the form:

$$R=\frac{1}{2}\underset{𝐪}{}G_𝐪(a_𝐪^{}b_𝐪^{}a_𝐪b_𝐪),$$
(14)
where $`G_𝐪`$ is unknown function independent on $`\mathrm{}`$. From $`dH_1(\mathrm{})/d\mathrm{}=[R,H_1(\mathrm{})]`$ it is simple to obtain linear differential equations

$$\frac{d}{d\mathrm{}}f_𝐪(\mathrm{})=G_𝐪g_𝐪(\mathrm{}),\frac{d}{d\mathrm{}}g_𝐪(\mathrm{})=G_𝐪f_𝐪(\mathrm{}).$$
(15)
The solution of (15) is simple:

$`f_𝐪(\mathrm{})`$ $`=`$ $`{\displaystyle \frac{f_𝐪(0)+g_𝐪(0)}{2}}e^{G_𝐪\mathrm{}}+{\displaystyle \frac{f_𝐪(0)g_𝐪(0)}{2}}e^{G_𝐪\mathrm{}},`$ (16)
$`g_𝐪(\mathrm{})`$ $`=`$ $`{\displaystyle \frac{f_𝐪(0)+g_𝐪(0)}{2}}e^{G_𝐪\mathrm{}}{\displaystyle \frac{f_𝐪(0)g_𝐪(0)}{2}}e^{G_𝐪\mathrm{}}.`$ (17)
Let the Hamiltonian (8) be diagonal at $`\mathrm{}=1`$: $`g_𝐪(1)=0`$. From this condition we immediately find the unknown function $`G_𝐪`$:

$$G_𝐪=\frac{1}{2}\mathrm{ln}\left(\frac{f_𝐪(0)+g_𝐪(0)}{f_𝐪(0)g_𝐪(0)}\right)$$
(18)
and the spectrum

$$ϵ_𝐪=f_𝐪(1)=\sqrt{f_𝐪^2(0)g_𝐪^2(0)}.$$
(19)
Thus, we see that in the case of exact diagonalization of quadratic Hamiltonian the results of FE and OS CUT methods coincide with each other.

## IV Electron-phonon coupling

Consider now the elimination of electron-phonon interaction to obtain the effective electron-electron scattering in the framework of OS CUT and FE methods. The Hamiltonian is (see, ):

$$H=H_0+H_{eph}+H_{ee},$$
(20)
where

$$H_0=\underset{𝐪}{}\omega _𝐪a_𝐪^{}a_𝐪+\underset{𝐤}{}\epsilon _𝐤c_𝐤^{}c_𝐤$$
(21)
describes free phonons ($`a_𝐪^{}`$, $`a_𝐪`$) and electrons ($`c_𝐤^{}`$, $`c_𝐤`$),

$$H_{eph}(\mathrm{})=\underset{𝐤,𝐪}{}M_{𝐤,𝐪}(\mathrm{})a_𝐪^{}c_{𝐤+𝐪}^{}c_𝐤+\mathrm{h}.\mathrm{c}.$$
(22)
describes the electron-phonon interaction and

$$H_{ee}(\mathrm{})=\underset{𝐤,𝐤^{},𝐪}{}V_{𝐤,𝐤^{},𝐪}(\mathrm{})c_{𝐤+𝐪}^{}c_{𝐤^{}𝐪}^{}c_𝐤^{}c_𝐤$$
(23)
is the Hamiltonian of electron-electron interaction. The initial conditions are:

$$M_{𝐤,𝐪}(\mathrm{}=0)=M_𝐪,V_{𝐤,𝐤^{},𝐪}(\mathrm{}=0)=V_{𝐤,𝐤^{},𝐪}(0).$$
(24)
Let us first consider the diagonalization (elimination of (22)) in the framework of OS CUT method. From a simple analysis it follows that the parameter independent generator should be taken in the form:

$$R=\underset{𝐤,𝐪}{}R_{𝐤,𝐪}a_𝐪^{}c_{𝐤+𝐪}^{}c_𝐤\mathrm{h}.\mathrm{c}.$$
(25)
The diagonalization procedure gives the following linear differential equations:

$$\frac{M_{𝐤,𝐪}(\mathrm{})}{\mathrm{}}=\alpha _{𝐤,𝐪}R_{𝐤,𝐪}$$
(26)
and for the electron pairs (with zero momentum)

$$\frac{V_{𝐤,𝐤,𝐪}(\mathrm{})}{\mathrm{}}=R_{𝐤,𝐪}M_{𝐤𝐪,𝐪}(\mathrm{})R_{𝐤𝐪,𝐪}M_{𝐤,𝐪}(\mathrm{}).$$
(27)
From (26) we obtain

$$M_{𝐤,𝐪}(\mathrm{})=M_{𝐤,𝐪}(0)\alpha _{𝐤,𝐪}R_{𝐤,𝐪}\mathrm{}.$$
(28)
Let the electron-phonon Hamiltonian (22) vanishes at $`\mathrm{}=1`$$`M_{𝐤,𝐪}(1)=0`$. From this condition we immediately find

$$R_{𝐤,𝐪}=M_{𝐤,𝐪}(0)/\alpha _{𝐤,𝐪}.$$
(29)
Simple solution of (27) with (28), (29) and initial conditions (24) gives

$$V_{𝐤,𝐤,𝐪}(1)=V_{𝐤,𝐤,𝐪}(0)\frac{M_𝐪^2\omega _𝐪}{\omega _𝐪^2(\epsilon _{𝐤+𝐪}\epsilon _𝐤)^2}.$$
(30)
This result corresponds to the Fröhlich’s effective interaction obtained in the second order of the perturbation theory (see, e.g., ).

Let us now consider the FE method. The generator is taken in the form :

$`\eta (\mathrm{})`$ $`=`$ $`[H_0,H_{eph}(\mathrm{})]`$ (31)
$`=`$ $`{\displaystyle \underset{𝐤,𝐪}{}}M_{𝐤,𝐪}(\mathrm{})\alpha _{𝐤,𝐪}a_𝐪^{}c_{𝐤+𝐪}^{}c_𝐤\mathrm{h}.\mathrm{c}.`$ (32)
with the energy difference

$$\alpha _{𝐤,𝐪}=\epsilon _{𝐤+𝐪}\epsilon _𝐤+\omega _𝐪.$$
(33)
This generator yields several contributions to $`dH(\mathrm{})/d\mathrm{}=[\eta (\mathrm{}),H(\mathrm{})]`$. The contribution to the change of $`M_{𝐤,𝐪}(\mathrm{})`$ results from $`[\eta (\mathrm{}),H_0]`$

$$\frac{M_{𝐤,𝐪}(\mathrm{})}{\mathrm{}}=\alpha _{𝐤,𝐪}^2M_{𝐤,𝐪}(\mathrm{})$$
(34)
with the solution

$$M_{𝐤,𝐪}(\mathrm{})=M_𝐪\mathrm{exp}(\alpha _{𝐤,𝐪}^2\mathrm{}),$$
(35)
where $`M_𝐪`$ is the initial electron-phonon coupling. For $`\alpha _{𝐤,𝐪}0`$ the solution (35) eliminates the electron-phonon interaction from the Hamiltonian at $`\mathrm{}\mathrm{}`$. The contribution to the effective electron-electron interaction is obtained from $`[\eta (\mathrm{}),H_{eph}(\mathrm{})]`$. For the electron pairs (with zero momentum) we have the following nonlinear differential equation:

$$\frac{V_{𝐤,𝐤,𝐪}(\mathrm{})}{\mathrm{}}=(\alpha _{𝐤,𝐪}+\alpha _{𝐤𝐪,𝐪})M_{𝐤,𝐪}(\mathrm{})M_{𝐤𝐪,𝐪}(\mathrm{})$$
(36)
with the solution

$$V_{𝐤,𝐤,𝐪}(\mathrm{})=V_{𝐤,𝐤,𝐪}(0)\frac{M_𝐪^2\omega _𝐪}{\omega _𝐪^2+(\epsilon _{𝐤+𝐪}\epsilon _𝐤)^2}.$$
(37)
This result differs from (30) (plus sign between the two squares in the denominator). The difference is explained as: “perturbation theory for Hamiltonians is not uniquely defined” .

In the FE method a choice of $`\eta (\mathrm{})`$ is motivated to be defined by $`[H_0,H_{eph}(\mathrm{})]`$ (31). However, a rigorous criterion that this choice (or, another one) is most optimal for a given problem seems has not yet been formulated.

## V Three-boson Hamiltonian

In this example we consider the elimination of forbidden three-boson terms from the Hamiltonian

$$H(\mathrm{})=\underset{𝐤}{}\epsilon _𝐤b_𝐤^{}b_𝐤+H_3(\mathrm{})+H_4(\mathrm{}),$$
(38)
where

$$H_3(\mathrm{})=\frac{1}{3}\underset{𝐤,𝐪}{}\{\mathrm{\Psi }_{𝐤,𝐪}(\mathrm{})b_𝐤^{}b_𝐪^{}b_{𝐤𝐪}^{}+\mathrm{h}.\mathrm{c}.\}$$
(39)
and

$`H_4(\mathrm{})`$ $`=`$ $`{\displaystyle \frac{1}{2}}{\displaystyle \underset{𝐤_1,𝐤_2;𝐤_3,𝐤_4}{}}\mathrm{\Phi }_{𝐤_1,𝐤_2;𝐤_3,𝐤_4}(\mathrm{})b_{𝐤_1}^{}b_{𝐤_2}^{}b_{𝐤_3}b_{𝐤_4}`$ (41)
$`\times \mathrm{\Delta }(𝐤_1+𝐤_2𝐤_3𝐤_4).`$
Here $`\mathrm{}`$ is a formal parameter and $`\mathrm{\Delta }(𝐤)`$ is the Kronecker delta function. For simplicity, we shall neglect the energy $`\epsilon _𝐤`$ dependence on $`\mathrm{}`$.

Using OS CUT method, we write the parameter-independent operator in the form :

$$R=\underset{𝐤,𝐪}{}(R_{𝐤,𝐪}b_𝐤^{}b_𝐪^{}b_{𝐤𝐪}^{}\mathrm{h}.\mathrm{c}.).$$
(42)
Calculating the commutators $`[R,_𝐤\epsilon _𝐤b_𝐤^{}b_𝐤]`$ and $`[R,H_3(\mathrm{})]`$, we obtain the following linear differential equations:

$$\frac{d}{d\mathrm{}}\mathrm{\Psi }_{𝐤,𝐪}(\mathrm{})=3(\epsilon _𝐤+\epsilon _𝐪+\epsilon _{𝐤𝐪})R_{𝐤,𝐪},$$
(43)
$`{\displaystyle \frac{d}{d\mathrm{}}}\mathrm{\Phi }_{𝐤_1,𝐤_2;𝐤_3,𝐤_4}(\mathrm{})`$ $`=`$ $`6R_{𝐤_1,𝐤_2}\mathrm{\Psi }_{𝐤_3,𝐤_4}^{}(\mathrm{})`$ (45)
$`6R_{𝐤_3,𝐤_4}^{}\mathrm{\Psi }_{𝐤_1,𝐤_2}(\mathrm{}).`$
Solving these equations with the initial conditions $`\mathrm{\Psi }_{𝐤_1,𝐤_2}(\mathrm{}=0)=\mathrm{\Psi }_{𝐤_1,𝐤_2}`$, $`\mathrm{\Phi }_{𝐤_1,𝐤_2;𝐤_3,𝐤_4}(\mathrm{}=0)=0`$, and demanding $`\mathrm{\Psi }_{𝐤_1,𝐤_2}(\mathrm{}=1)=0`$, we obtain

$$R_{𝐤,𝐪}=\frac{\mathrm{\Psi }_{𝐤,𝐪}}{3(\epsilon _𝐤+\epsilon _𝐪+\epsilon _{𝐤𝐪})}$$
(46)
and

$`\mathrm{\Phi }_{𝐤_1,𝐤_2;𝐤_3,𝐤_4}(1)`$ $`=`$ $`\mathrm{\Psi }_{𝐤_1,𝐤_2}\mathrm{\Psi }_{𝐤_3,𝐤_4}^{}({\displaystyle \frac{1}{\epsilon _{𝐤_1}+\epsilon _{𝐤_2}+\epsilon _{𝐤_1𝐤_2}}}`$ (48)
$`+{\displaystyle \frac{1}{\epsilon _{𝐤_3}+\epsilon _{𝐤_4}+\epsilon _{𝐤_3𝐤_4}}}).`$
It should be noted that this result for the effective four-boson interaction coincides (with the additional condition $`\mathrm{\Psi }_1=\mathrm{\Psi }_1^{}`$) with the result of classical diagonalization of (39) (when $`b_𝐤^{}`$ and $`b_𝐤`$ are $`c`$-numbers) by canonical quasi-linear transformation .

Let us now consider the same problem in the framework of FE method. The (parameter-dependent) generator in this case is:

$$\eta (\mathrm{})=\frac{1}{3}\underset{𝐤,𝐪}{}\{\beta _{𝐤,𝐪}\mathrm{\Psi }_{𝐤,𝐪}(\mathrm{})b_𝐤^{}b_𝐪^{}b_{𝐤𝐪}^{}\mathrm{h}.\mathrm{c}.\},$$
(49)
where

$$\beta _{𝐤,𝐪}=\epsilon _𝐤+\epsilon _𝐪+\epsilon _{𝐤𝐪}.$$
(50)
Corresponding differential equation for the three-boson amplitude has the form:

$$\frac{d}{d\mathrm{}}\mathrm{\Psi }_{𝐤,𝐪}(\mathrm{})=\beta _{𝐤,𝐪}^2\mathrm{\Psi }_{𝐤,𝐪}(\mathrm{})$$
(51)
with the solution (compare with (35)):

$$\mathrm{\Psi }_{𝐤,𝐪}(\mathrm{})=\mathrm{\Psi }_{𝐤,𝐪}\mathrm{exp}(\beta _{𝐤,𝐪}^2\mathrm{}).$$
(52)
The differential equation for the four-boson effective amplitude in this case is

$$\frac{d\mathrm{\Phi }_{𝐤_1,𝐤_2;𝐤_3,𝐤_4}(\mathrm{})}{d\mathrm{}}=2(\beta _{𝐤_1,𝐤_2}+\beta _{𝐤_3,𝐤_4})\mathrm{\Psi }_{𝐤_1,𝐤_2}(\mathrm{})\mathrm{\Psi }_{𝐤_3,𝐤_4}^{}(\mathrm{}).$$
(53)
The solution at $`\mathrm{}\mathrm{}`$ (when $`\mathrm{\Psi }_{𝐤,𝐪}(\mathrm{})`$ vanishes) has the form

$$\mathrm{\Phi }_{𝐤_1,𝐤_2;𝐤_3,𝐤_4}(1)=2\mathrm{\Psi }_{𝐤_1,𝐤_2}\mathrm{\Psi }_{𝐤_3,𝐤_4}^{}\frac{\beta _{𝐤_1,𝐤_2}+\beta _{𝐤_3,𝐤_4}}{\beta _{𝐤_1,𝐤_2}^2+\beta _{𝐤_3,𝐤_4}^2}.$$
(54)
We see that this formula differs from the OS CUT result (48) and therefore (54) does not correspond to the classical case .

## VI Non-perturbative solution of quantum equations

### A Heisenberg picture

The Heisenberg equation for an operator $`A(t)`$

$$i\mathrm{}\frac{d}{dt}A(t)=[A(t),H(t)]$$
(55)
can be rewritten in the form:

$$\frac{d}{dt}A(t)=[iH(t)/\mathrm{},A(t)].$$
(56)
Formally this equation is similar to (2) (or, (4) for time-dependent Hamiltonian): time $`t`$ plays a role of a formal parameter and $`iH(t)/\mathrm{}`$ is an anti-Hermitian operator. As in the CUT method, we can write the most general form of $`A(t)`$ as the expansion in terms of operator combinations with unknown time-dependent coefficients. After substituting these form into (56) with a given $`iH(t)/\mathrm{}`$, one can obtain a set of differential equations by comparing coefficients in analogous operator compositions in both sides. These equations will be linear for time-independent $`H`$ and nonlinear for time-dependent $`H(t)`$. Then we solve these equations with the initial conditions defined by $`A(t)=A(0)`$. Recently Mišta and Filip  have developed this idea as a method of non-perturbative solutions of nonlinear Heisenberg equations with an illustrative example of two coupled harmonic oscillators.

### B Density matrix

It is obvious that a quite analogous idea can be applied for the density matrix. In this case we can write the density matrix equation as

$$\frac{d}{dt}\rho (t)=[iH(t)/\mathrm{},\rho (t)],$$
(57)
where $`iH(t)/\mathrm{}`$ plays a role of “transformation generator”.

As an illustration consider a system of nuclear spins $`𝐈_j`$ with a Hamiltonian $`H=H_0+H_{int}`$, where $`H_0=\mu B_jI_j^z`$ is the Zeeman energy and $`H_{int}`$ describes spin-spin interactions. As an initial state assume the transverse magnetization. In the high temperature approximation one can write $`\rho (0)=1\alpha _jI_j^x`$. Then the density matrix can be represented as $`\rho (t)=1_{n=1}^{\mathrm{}}\mathrm{\Delta }\rho _n(t)`$, where $`\mathrm{\Delta }\rho _n(t)`$ describes compositions of $`n`$ spins with time dependent coefficients, e.g., $`\mathrm{\Delta }\rho _1(t)=_j\left(X_j(t)I_j^x+Y_j(t)I_j^y\right)`$. From $`[H_0,\mathrm{\Delta }\rho _1(t)]`$ we obtain differential equations $`dX_j(t)/dt=\omega _0Y_j(t)`$, $`dY_j(t)/dt=\omega _0X_j(t)`$ and solutions $`X_j(t)=\alpha \mathrm{cos}\omega _0t`$, $`Y_j(t)=\alpha \mathrm{sin}\omega _0t`$, where $`\omega _0=\mu B/\mathrm{}`$ is the nuclear magnetic resonance frequency Contributions to two spin motions can be found from $`[H_0,\mathrm{\Delta }\rho _2(t)]+[H_{int},\mathrm{\Delta }\rho _1(t)]`$ and so on. This example can demonstrate how the order stored in magnetization is transformed to higher-order spin correlations.

The author wish to thank G. S. Uhrig and A. K. Khitrin for valuable comments.