text:
# Density and current response functions in strongly disordered electron systems: Diffusion, electrical conductivity and Einstein relation

## I Introduction

Low-energy physics of equilibrium systems with weakly interacting electrons is well understood both qualitatively and quantitatively. The relevant information about the equilibrium system is contained in one- and two-particle Green functions, in particular in their behavior around the Fermi energy. A number of reliable approximate methods have been developed for the calculation of these functions. Among them, systematic renormalized perturbation expansions based on the many-body Feynman diagrammatic technique have proved most effective. When impurities or quenched configurational randomness are added we still are able to describe equilibrium properties of such a system quite reliably within diagrammatic schemes and a mean-field-type coherent-potential approximation (CPA).Elliot74

The situation gets less straightforward if we come out of equilibrium. This is the case when the system is disturbed by a time-dependent external force but does not manage to reach a new equilibrium within the relaxation times of the experimental setup. We then have to determine the response of the system to the external perturbation to obtain experimentally relevant data. However, we do not have at our disposal many established methods to calculate nonequilibrium response functions. When the perturbation is weak, which is the most common situation in practice, we can use linear response theory to calculate effectively the response functions.Plischke94  A fundamental tool for the calculation of response functions within the linear-response theory is the Kubo formalism.Kubo57  It determines how the response functions of the disturbed nonequilibrium system can be expressed in terms of two-particle Green (correlation) functions that are characteristics of an equilibrium state. Kubo formulas are means for the description of nonequilibrium systems with equilibrium functions. Unlike the response functions, the equilibrium Green functions obey equations of motion and are suitable for developing systematic approximations within many-body diagrammatic techniques.

Kubo formalism, however, provides different representations for different response functions that we have to approximate separately. The problem that emerges with independent approximations for different response functions via different Kubo formulas is that we are usually unable to keep exact relations between response functions that may hold due to special symmetries of the system under investigation. Most pronounced case with “hidden” symmetries is the response to a weak electromagnetic perturbation where we have to satisfy gauge invariance and charge and current conservations. Instantaneous field-dependent deviations of charge and current densities from their equilibrium values are described by density and current response functions. Within linear-response theory these functions are obtained from independent Kubo formulas with density-density and current-current correlation functions, respectively. Although both correlation functions are derived from the same two-particle Green function, they are normally approximated independently according to the purpose they serve to.

When quantum coherence is negligible we can calculate the transport properties from the coherent-potential approximation.Elliot74  There is no contribution to the homogeneous current-current correlation function beyond the single electron-hole bubble in the single-band coherent-potential approximation.Velicky69  Hence most theories beyond the CPA, either on model or realistic level, use the Kubo formula for the electrical conductivity with the current-current correlation function to determine transport properties of disordered solids.Gonis92 ; Turek97

On the other hand, when quantum coherence effects are substantial and we expect the Anderson metal-insulator transition, electrical conductivity is usually calculated from the electron-hole correlation function with the aid of the diffusion constant controlling its low-energy behavior near the diffusion pole.Vollhardt92 ; Belitz94  The diffusion pole is crucial and enables application of scaling arguments and the renormalization group approach to Anderson localization.Lee85

There is a number of more or less heuristic arguments in the literature that relate the density response with the conductivity.Rammer98  They are essentially based on gauge invariance and charge conservation of a semiclassical equilibrium description of an electron gas exposed to an electromagnetic field. In weakly disordered quantum systems (described by continuum models) a relation between the density response and conductivity can formally be derived as follows.Rammer98  Gauge invariance is used to relate the external scalar potential with the electric field $`𝐄=\mathbf{}\phi `$. The current density generated by the external field then is

$$𝐣(𝐪,\omega )=𝝈(𝐪,\omega )𝐄(𝐪,\omega )=i𝝈(𝐪,\omega )𝐪\phi (𝐪,\omega )$$
(1)
where $`𝝈(𝐪,\omega )`$ denotes tensor of the electrical conductivity. Charge conservation is expressed by a continuity equation. In equilibrium we can use its operator form that follows from the Heisenberg equations of motion for the current and density operators. For Hamiltonians with quadratic dispersion relation we have

$$e_t\widehat{n}(𝐱,t)+\mathbf{}\widehat{𝐣}(𝐱,t)=0.$$
(2)
Energy-momentum representation of the continuity equation in the ground-state solution is

$$i\omega e\delta n(𝐪,\omega )+i𝐪𝐣(𝐪,\omega )=0.$$
(3)
We have to use a density variation of the equilibrium density, i. e., the externally induced density $`\delta n(𝐪,\omega )=n(𝐪,\omega )n_0`$ in the continuity equation with averaged values of operators. From the above equations and for linear response $`\delta n(𝐪,\omega )=e\chi (𝐪,\omega )\phi (𝐪,\omega )`$ we obtain in the isotropic case

$$\sigma (𝐪,\omega )=\frac{ie^2\omega }{q^2}\chi (𝐪,\omega ).$$
(4)
The derived equality formally holds for complex functions without restrictions on momenta or frequencies. Frequencies can, in principle, be even complex. Relation (4) is often taken as granted for the whole range of the disorder strength and used for the definition of the zero-temperature conductivity when describing Anderson localization transition.Vollhardt92 ; Belitz94

Although the above derivation may seem very general it suffers from a few flaws. First, the operator continuity equation (2) cannot be directly used out of equilibrium. A new term due to causality of the response functions enters the continuity equation when the equilibrium is disturbed. Moreover the nonequilibrium density and current operators in linear-response theory no longer obey Heisenberg equations of motion with the perturbed Hamiltonian. The perturbation is decoupled from the equilibrium Hamiltonian and is treated only to linear order. Second, strongly disordered electron systems can be desribed only by lattice models with a nonquadratic dispersion relation. Continuity equation (2) is then to be modified even in equilibrium. Hence the above, in the literature broadly disseminated reasoning leading to Eq. (4) cannot be fully trusted in quantitative studies of strongly disordered systems out of equilibrium.

Presently there is no reliable theory of strongly disordered electrons beyond the mean-field CPA. Since mean-field approximations do not include vertex corrections to the electrical conductivity they are not suitable for the investigation of localization effects in three spatial dimensions. Unlike low dimensions ($`d2`$), Anderson localization in $`d=3`$ may occur only in strongly disordered systems. Anderson localization in bulk systems has not yet been understood or quantitatively described in a satisfactory manner. Its quantitative description demands bridging the gap between the mean-field (CPA) transport theory and theories for weakly and strongly localized electron states. It is clear that such an interpolating scheme should be based on advanced approximations for two-particle functions.

Recently a diagrammatic method for summations of classes of two-particle diagrams has been proposed.Janis01b  This theory has a potential to interpolate between the mean-field and localization theories provided exact relations between the density and current response functions have been proved for this scheme in the whole range of the disorder strength. Criteria for validity of relations between density-density and current-current correlation functions in approximate treatments of strongly disordered systems have not yet been established. When we want to go beyond the CPA we still have to establish such relations in the metallic regime of tight-binding models with extended electron states. Simultaneously we have to formulate conditions under which the derived relations hold or may be broken.

The aim of this paper is threefold. First, we derive various relations between current-current, current-density, and density-density correlation functions for strongly disordered electrons out of equilibrium, i. e., beyond the reach of the existing derivations. Second, we articulate conditions under which these relations hold or to which extent they may be broken in quantitative treatments and show when the conductivity can be calculated from the density response function and vice versa. Third, we introduce a generalization of diffusion via a quantum response function, test validity of the Einstein relation for this quantum diffusion, and relate quantum diffusion with its classical counterpart and with the diffusion pole in the electron-hole correlation function.

To reach this goal we use Kubo formalism and averaged many-body Green functions. In this approach a weaker form of the continuity equation expressed in terms of one- and two-particle Green functions replaces the operator identity. It is derived from equations of motion for Green functions and Ward identities between them. We use two Ward identities that are not fully equivalent. One expresses conservation of probability and the other reflects charge conservation. Each Ward identity is used in a different manner to relate the conductivity with the density response.

The layout of the paper is as follows. In Sec. II we summarize the definitions and useful representations of the density and current response functions. In Sec. III we show how the current-current correlation function emerges from a momentum expansion of the density response function. Then in Sec. IV we derive continuity equations expressed in Green functions and use them to show how the density response function can be revealed from the conductivity. Based on the exactly derived relation between the density and current response function we introduce in Sec. V quantum diffusion, derive the diffusion pole in the hydrodynamic limit of the density response function and relate the diffusion constant to the conductivity. In Sec. VI we illustrate general results on an exactly solvable limit of infinite spatial dimensions. In Appendix we discuss assumptions for and the range of validity of the Vollhardt-Wölfle-Ward identity used in the derivation of the continuity equation for Green functions.

## II Density and current response functions in disordered systems

A simplest description of the electron motion in impure, weakly correlated metals is provided by a tight-binding Anderson model. It assumes noninteracting spinless electrons moving in a random, site-diagonal potential $`V_i`$. Its Hamiltonian reads

$`\widehat{H}_{AD}`$ $`=`$ $`{\displaystyle \underset{<ij>}{}}t_{ij}c_i^{}c_j+{\displaystyle \underset{i}{}}V_ic_i^{}c_i.`$ (5)
The values of the random potential $`V_i`$ are site-independent and obey a disorder distribution $`\rho (V)`$. I. e., a function depending on the random potential $`V_i`$ is averaged as

$`X(V_i)_{av}`$ $`=`$ $`{\displaystyle _{\mathrm{}}^{\mathrm{}}}𝑑V_i\rho (V_i)X(V_i).`$ (6)
The averaged two-particle propagator (resolvent) is defined as the averaged product of one-particle propagators

$$\begin{array}{c}G_{ij,kl}^{(2)}(z_1,z_2)=\hfill \\ \hfill \left[z_1\widehat{1}\widehat{t}\widehat{V}\right]_{ij}^1\left[z_2\widehat{1}\widehat{t}\widehat{V}\right]_{kl}^1_{av}.\end{array}$$
(7)
Fourier transform to momentum space is not uniquely defined, since due to momentum conservation we have only three independent momenta. For our purposes we choose the following notation and definition of the Fourier transform

$$\begin{array}{c}G_{\mathrm{𝐤𝐤}^{}}^{(2)}(z_1,z_2;𝐪)=\frac{1}{N}\underset{ijkl}{}e^{i(𝐤+𝐪/2)𝐑_i}e^{i(𝐤^{}+𝐪/2)𝐑_j}\hfill \\ \hfill \times e^{i(𝐤^{}𝐪/2)𝐑_k}e^{i(𝐤𝐪/2)𝐑_l}G_{ij,kl}^{(2)}(z_1,z_2).\end{array}$$
(8)
The density response function is determined from the following thermodynamic representation of the Kubo formulaEnz92

$$\begin{array}{c}\chi (𝐪,i\nu _m)=\hfill \\ \hfill \frac{1}{N^2}\underset{\mathrm{𝐤𝐤}^{}}{}k_BT\underset{n=\mathrm{}}{\overset{\mathrm{}}{}}G_{\mathrm{𝐤𝐤}^{}}^{(2)}(i\omega _n,i\omega _n+i\nu _m;𝐪).\end{array}$$
(9)
Here $`\omega _n=(2n+1)\pi T`$ are fermionic and $`\nu _m=2m\pi T`$ bosonic Matsubara frequencies at temperature $`T`$. Due to the analytic properties of the two-particle Green function we can perform analytic continuation from imaginary Matsubara to real frequencies using contour integrals. Only the contribution along the real axis remains and we obtain

$$\begin{array}{c}\chi (𝐪,\omega +i0^+)=\frac{1}{N^2}\underset{\mathrm{𝐤𝐤}^{}}{}_{\mathrm{}}^{\mathrm{}}\frac{dE}{2\pi i}\{[f(E+\omega )f(E)]G_{\mathrm{𝐤𝐤}^{}}^{AR}(E,E+\omega ;𝐪)\hfill \\ \hfill +f(E)G_{\mathrm{𝐤𝐤}^{}}^{RR}(E,E+\omega ;𝐪)f(E+\omega )G_{\mathrm{𝐤𝐤}^{}}^{AA}(E,E+\omega ;𝐪)\}.\end{array}$$
(10)
We used the following definition

$$G_{\mathrm{𝐤𝐤}^{}}^{AR}(E,E+\omega ;𝐪)=G_{\mathrm{𝐤𝐤}^{}}^{(2)}(Ei0^+,E+\omega +i0^+;𝐪)$$

and analogously for the functions $`G^{RR}`$ and $`G^{AA}`$, where both energies are from the upper and lower complex energy half-plane, respectively. We denoted the Fermi function$`f(E)=[1+\mathrm{exp}\{\beta (E\mu )\}]^1`$ with the chemical potential $`\mu `$.

The two-particle Green function is generally determined from a Bethe-Salpeter equation in which the input is a two-particle irreducible vertex. The latter is known only approximately, except for special limits. However, a few specific elements of the full two-particle Green function can be evaluated without knowing the two-particle irreducible vertex. It is sufficient to know only the one-particle Green function to determine them. Determination of two-particle functions from one-particle ones is enabled by Ward identities. The Ward identity relating the averaged one- and two-particle Green functions reads

$$\begin{array}{c}\frac{1}{N}\underset{𝐤^{}}{}G_{\mathrm{𝐤𝐤}^{}}^{(2)}(z_1,z_2;\mathrm{𝟎})=\frac{1}{z_2z_1}[G(𝐤,z_1)\hfill \\ \hfill G(𝐤,z_2)].\end{array}$$
(11)
and was proved for the first time within the coherent-potential approximation by Velický.Velicky69  It is a nonperturbative identity valid quite generally beyond the CPA.Janis01b  As a consequence of this identity we obtain vanishing of the density response function for a homogeneous global perturbation, i. e., for $`𝐪=0`$. Using (11) in (10) we easily find

$$\begin{array}{c}\chi (\mathrm{𝟎},\omega +i0^+)\hfill \\ \hfill =\frac{1}{\omega }_{\mathrm{}}^{\mathrm{}}\frac{dE}{2\pi i}\{f(E+\omega )[G^A(E+\omega )\\ \hfill G^R(E+\omega )]f(E)[G^A(E)G^R(E)]\}=0.\end{array}$$
(12)
Another situation where we do not need to know the two-particle irreducible vertex is the static limit, $`\omega =0`$. First, the static density response function is real and reads

$$\chi (𝐪,0)=\frac{1}{N^2}\underset{\mathrm{𝐤𝐤}^{}}{}_{\mathrm{}}^{\mathrm{}}\frac{dE}{\pi }f(E)\mathrm{}G_{\mathrm{𝐤𝐤}^{}}^{RR}(E,E;𝐪)$$
(13a)
that in the limit $`𝐪0`$ goes over to
$$\chi (\mathrm{𝟎},0)=_{\mathrm{}}^{\mathrm{}}\frac{d\omega }{\pi }\frac{f(\omega )}{\omega }\mathrm{}G^R(\omega )\underset{T0}{\overset{}{}}n_F$$
(13b)
with $`n_F`$ as the density of states at the Fermi level and $`G^R(\omega )=N^1_𝐪G^R(𝐪,\omega )`$. Note that Eqs. (13) hold only in the limit $`\omega /q0`$. In the inverse case, $`q/\omega 0`$, Eq. (12) applies. Noncommutativity of the limits $`\omega 0`$ and $`q0`$ indicates that the point $`q=0,\omega =0`$ is not analytic in the same manner as in the Fermi liquid theory.Rickayzen80

The current-current correlation function $`\mathrm{\Pi }_{\alpha \beta }(𝐪,t)`$ relates the current of the system $`𝐣(𝐪,t)`$ perturbed by an external vector potential $`𝐀(𝐪,t)`$. Tensor of the electrical conductivity $`\sigma _{\alpha \beta }(𝐪,t)`$ determines the current response to an external electric field. It can be obtained from the current response function by using a relation $`𝐄(𝐱,t)=\dot{𝐀}(𝐱,t)`$. In the Fourier representation we obtainMahan00

$$\sigma _{\alpha \beta }(𝐪,\omega _+)=\frac{i}{\omega }\left[\mathrm{\Pi }_{\alpha \beta }(𝐪,\omega _+)\mathrm{\Pi }_{\alpha \beta }(𝐪,0)\right]$$
(14)
with $`\omega _+=\omega +i0^+`$. The compensation term $`\mathrm{\Pi }_{\alpha \beta }(𝐪,0)`$ on the right-hand side of Eq. (14) is real and was introduced to warrant finiteness of the complex conductivity. This additive term is not generated from the Kubo formula with a commutator of current operators but originates from gauge invariance of the electron gas in an external electromagnetic field.Note1  According to Eq. (14) there is no current in the system with a static vector potential and the static (complex) conductivity can be defined only via a dynamical one from the limit $`\omega 0`$.

Unlike the density response function the current response function cannot be simplified in any limit and to determine the conductivity from Eq. (14) we have to know the two-particle irreducible vertex. The current response function alike the density response function is a complex quantity with the same analytic properties. Of interest for us are the parallel components that for real frequencies read

$$\begin{array}{c}\sigma _{\alpha \alpha }(𝐪,\omega )=\frac{e^2}{4}\frac{1}{N^2}\underset{𝐤,𝐤^{}}{}\left[v_\alpha (𝐤^{}+𝐪/2)v_\alpha (𝐤+𝐪/2)\right]\left[v_\alpha (𝐤^{}+𝐪/2)v_\alpha (𝐤^{}+𝐪/2)\right]\hfill \\ \hfill \times _{\mathrm{}}^{\mathrm{}}\frac{dE}{2\pi \omega }\{[f(E+\omega )f(E)]G_{\mathrm{𝐤𝐤}^{}}^{AR}(E,E+\omega ;𝐪)+f(E)G_{\mathrm{𝐤𝐤}^{}}^{RR}(E,E+\omega ;𝐪))\\ \hfill f(E+\omega )G_{\mathrm{𝐤𝐤}^{}}^{AA}(E,E+\omega ;𝐪)f(E)[G_{\mathrm{𝐤𝐤}^{}}^{RR}(E,E;𝐪)G_{\mathrm{𝐤𝐤}^{}}^{AA}(E,E;𝐪)]\}.\end{array}$$
(15)
We denoted the group velocity $`v_\alpha (𝐤)=m^1/k_\alpha ϵ(𝐤)`$, where $`ϵ(𝐤)`$ is the dispersion relation of the underlying lattice.

## III Dynamical conductivity calculated from the density response

As first we show how the dynamical conductivity can be obtained from the long-range behavior of the density response function (hydrodynamic limit of small transfer momentum). To this purpose we use the Velický-Ward identity, Eq. (11). It enables us to evaluate averaged matrix elements of the two-particle Green function with zero transfer momentum in terms of the one-particle propagator. Zero transfer momentum is a severe restriction on the applicability and utilization of this Ward identity. Physically we are interested in the hydrodynamic limit, i. e., the asymptotics $`q0`$. We are unable to extend identity (11) beyond the homogeneous case. However, if we assume analyticity of the asymptotics $`q0`$, we can use momentum $`q`$ as an expansion parameter and investigate the hydrodynamic limit perturbatively. The hydrodynamic limit is analytic if frequency $`\omega 0`$, i. e., we are in the regime $`q/\omega 1`$. We then can expand the density response function in powers of momentum $`q`$.

We define two configuration-dependent resolvent operators

$$\widehat{G}_\pm (z)=\left[z_\pm \widehat{1}\widehat{t}\mathrm{\Delta }_1\widehat{t}\mathrm{\Delta }_2\widehat{t}\widehat{V}\right]^1$$
(16)
where $`\mathrm{\Delta }_1\widehat{t}\pm \mathrm{\Delta }_2\widehat{t}`$ is a difference in the dispersion relation of the two resolvents and is defined via its matrix elements $`𝐤|\mathrm{\Delta }_1\widehat{t}|𝐤^{}=\delta (𝐤^{}𝐤)𝐯(𝐤)𝐪/2`$, $`𝐤|\mathrm{\Delta }_2\widehat{t}|𝐤^{}=\delta (𝐤^{}𝐤)(𝐪_k)^2ϵ(𝐤)/8`$. We will be interested in the following function of two energies

$$I(z_+,z_{};𝐪)=\frac{1}{N}\text{Tr}\left[\widehat{G}_+(z_+)\widehat{G}_{}(z_{})\right]$$
(17)
that can be, after averaging over the configurations of the random potential, expressed via the two-particle Green function

$$I(z_+,z_{};𝐪)_{av}=\frac{1}{N^2}\underset{\mathrm{𝐤𝐤}^{}}{}G_{\mathrm{𝐤𝐤}^{}}^{(2)}(z_+,z_{};𝐪).$$
(18)
We assume that the expansion in the operators $`\mathrm{\Delta }_{1,2}\widehat{t}`$ commutes with the configurational averaging. We first expand the quantity $`I`$ and then average the series term by term. It is sufficient for our purposes to expand only to second order in momentum $`𝐪`$. This precision determines the leading small-momentum behavior.

Expanding quantity $`I`$ in $`\mathrm{\Delta }_{1,2}\widehat{t}`$ we have to keep the order of operators in the products, since $`\mathrm{\Delta }_{1,2}\widehat{t}`$ and the resolvent $`\widehat{G}(z)`$ do not commute. The expansion to second order reads

$$\begin{array}{c}I_2(z_+,z_{}𝐪)=\text{Tr}\{\widehat{G}(z_+)\widehat{G}(z_{})+\widehat{G}(z_+)[\mathrm{\Delta }_1\widehat{t}+\mathrm{\Delta }_2\widehat{t}]\widehat{G}(z_+)\widehat{G}(z_{})\hfill \\ \hfill \widehat{G}(z_+)\widehat{G}(z_{})\left[\mathrm{\Delta }_1\widehat{t}\mathrm{\Delta }_2\widehat{t}\right]\widehat{G}(z_{})\widehat{G}(z_+)\mathrm{\Delta }_1\widehat{t}\widehat{G}(z_+)\widehat{G}(z_{})\mathrm{\Delta }_1\widehat{t}\widehat{G}(z_{})\\ \hfill +\widehat{G}(z_+)[\mathrm{\Delta }_1\widehat{t}\widehat{G}(z_+)\mathrm{\Delta }_1\widehat{t}\widehat{G}(z_+)+\widehat{G}(z_{})\mathrm{\Delta }_1\widehat{t}\widehat{G}(z_{})\widehat{\mathrm{\Delta }_1t}]\widehat{G}(z_{})\}.\end{array}$$
(19)
Each direct product of the resolvent operators $`\widehat{G}(z_\pm )\widehat{G}(z_{})`$ can be simplified using identity (11). Doing this consequently we end up with a sum of products of two resolvents. We have three different terms to analyze:

$$\begin{array}{c}\text{Tr}\left\{\widehat{G}(z_+)\left[\left(\mathrm{\Delta }_1\widehat{t}+\mathrm{\Delta }_2\widehat{t}\right)\widehat{G}(z_{})\widehat{G}(z_+)\left(\mathrm{\Delta }_1\widehat{t}\mathrm{\Delta }_2\widehat{t}\right)\right]\widehat{G}(z_{})\right\}=\hfill \\ \hfill \frac{1}{z_+z_{}}\{\text{Tr}\left[(\widehat{G}^{}(z_+)+\widehat{G}^{}(z_{}))\mathrm{\Delta }_1\widehat{t}\right]+\text{Tr}\left[(\widehat{G}^{}(z_+)\widehat{G}^{}(z_{}))\mathrm{\Delta }_2\widehat{t}\right]\\ \hfill \frac{2}{z_+z_{}}\text{Tr}\left[\mathrm{\Delta }_1\widehat{t}(\widehat{G}(z_+)\widehat{G}(z_{}))\right]\},\end{array}$$
(20a)
$$\begin{array}{c}\text{Tr}\left\{\widehat{G}(z_+)\mathrm{\Delta }_1\widehat{t}\widehat{G}(z_+)\widehat{G}(z_{})\mathrm{\Delta }_1\widehat{t}\widehat{G}(z_{})\right\}=\frac{1}{(z_+z_{})^2}\{\text{Tr}\left[\widehat{G}(z_+)\mathrm{\Delta }_1\widehat{t}\widehat{G}(z_+)\mathrm{\Delta }_1\widehat{t}\right]\hfill \\ \hfill +\text{Tr}\left[\widehat{G}(z_{})\widehat{\mathrm{\Delta }_1t}\widehat{G}(z_{})\mathrm{\Delta }_1\widehat{t}\right]2\text{Tr}\left[\widehat{G}(z_+)\mathrm{\Delta }_1\widehat{t}\widehat{G}(z_{})\mathrm{\Delta }_1\widehat{t}\right]\},\end{array}$$
(20b)
$$\begin{array}{c}\text{Tr}\left\{\widehat{G}(z_+)\left[\mathrm{\Delta }_1\widehat{t}\widehat{G}(z_+)\mathrm{\Delta }_1\widehat{t}\widehat{G}(z_+)+\widehat{G}(z_{})\mathrm{\Delta }_1\widehat{t}\widehat{G}(z_{})\mathrm{\Delta }_1\widehat{t}\right]\widehat{G}(z_{})\right\}=\hfill \\ \hfill \frac{1}{z_+z_{}}\{\text{Tr}\left[\widehat{G}^{}(z_+)\mathrm{\Delta }_1\widehat{t}\widehat{G}(z_+)\mathrm{\Delta }_1\widehat{t}\right]\text{Tr}\left[\widehat{G}^{}(z_{})\mathrm{\Delta }_1\widehat{t}\widehat{G}(z_{})\mathrm{\Delta }_1\widehat{t}\right]\\ \hfill \frac{1}{z_+z_{}}\text{Tr}\left[(\widehat{G}(z_+)\widehat{G}(z_{}))\mathrm{\Delta }_1\widehat{t}(\widehat{G}(z_+)\widehat{G}(z_{}))\mathrm{\Delta }_1\widehat{t}\right]\}.\end{array}$$
(20c)
We insert Eqs. (20) in expansion (19) and average over the configurations of the random potential $`V_i`$. To recover the density response function we have to limit the complex energies to the real axis from above or below with which we distinguish the causality. We define correlation functions that can be represented as a trace of the averaged two-particle Green function

$$\mathrm{\Phi }_E^{ab}(𝐪,\omega )=\frac{1}{N^2}\underset{\mathrm{𝐤𝐤}^{}}{}G_{\mathrm{𝐤𝐤}^{}}^{ab}(E,E+\omega ;𝐪).$$
(21)
where $`a,b`$ stand for $`A,R`$ depending on whether the imaginary part of the corresponding frequency argument is positive or negative. We then can represent the leading asymptotics in small momenta $`𝐪`$ of the averaged electron-hole correlation function

$$\begin{array}{c}\mathrm{\Phi }_E^{AR}(𝐪,\omega )\mathrm{\Phi }_E^{AR}(\mathrm{𝟎},\omega )=\frac{1}{2\omega }\left[_E\left(G^R(E+\omega )+G^A(E)\right)𝐪𝐯\right]\hfill \\ \hfill +\frac{1}{8\omega }\left[_E\left(G^R(E+\omega )G^A(E)\right)(𝐪)^2ϵ\right]+\frac{1}{\omega ^2}\left[\left(G^R(E+\omega )G^A(E)\right)𝐪𝐯\right]\\ \hfill \frac{1}{4\omega N^2}\underset{\mathrm{𝐤𝐤}^{}}{}𝐪𝐯(𝐤)𝐪𝐯(𝐤^{})\{_\epsilon [G_{\mathrm{𝐤𝐤}^{}}^{AA}(E+\epsilon ,E;\mathrm{𝟎})G_{\mathrm{𝐤𝐤}^{}}^{RR}(E+\omega ,E+\omega +\epsilon ;\mathrm{𝟎})]_{\epsilon =0}\\ \hfill +\frac{2}{\omega }[G_{\mathrm{𝐤𝐤}^{}}^{AA}(E,E;\mathrm{𝟎})+G_{\mathrm{𝐤𝐤}^{}}^{RR}(E+\omega ,E+\omega ;\mathrm{𝟎})2G_{\mathrm{𝐤𝐤}^{}}^{AR}(E,E+\omega ;\mathrm{𝟎})]\}\end{array}$$
(22a)
and analogously of the electron-electron ($`\mathrm{\Phi }^{RR}`$) and hole-hole ($`\mathrm{\Phi }^{AA}`$) ones that can be written generically as
$$\begin{array}{c}\mathrm{\Phi }_E^{aa}(𝐪,\omega )\mathrm{\Phi }_E^{aa}(\mathrm{𝟎},\omega )=\frac{1}{2\omega }\left[_E\left(G^a(E+\omega )+G^a(E)\right)𝐪𝐯\right]\hfill \\ \hfill +\frac{1}{8\omega }\left[_E\left(G^a(E+\omega )G^a(E)\right)(𝐪)^2ϵ\right]+\frac{1}{\omega ^2}\left[\left(G^a(E+\omega )G^a(E)\right)𝐪𝐯\right]\\ \hfill \frac{1}{4\omega N^2}\underset{\mathrm{𝐤𝐤}^{}}{}𝐪𝐯(𝐤)𝐪𝐯(𝐤^{})\{_\epsilon [G_{\mathrm{𝐤𝐤}^{}}^{aa}(E+\epsilon ,E;\mathrm{𝟎})G_{\mathrm{𝐤𝐤}^{}}^{aa}(E+\omega ,E+\omega +\epsilon ;\mathrm{𝟎})]_{\epsilon =0}\\ \hfill +\frac{2}{\omega }[G_{\mathrm{𝐤𝐤}^{}}^{aa}(E,E;\mathrm{𝟎})+G_{\mathrm{𝐤𝐤}^{}}^{aa}(E+\omega ,E+\omega ;\mathrm{𝟎})2G_{\mathrm{𝐤𝐤}^{}}^{aa}(E,E+\omega ;\mathrm{𝟎})]\}.\end{array}$$
(22b)
In the above equations we used angular brackets to denote summation over fermionic momenta in the first Brillouin zone of one-electron functions, that is

$$G(\omega )f_𝐪=\frac{1}{N}\underset{𝐤}{}G(𝐤,\omega )f_𝐪(𝐤).$$
(23)
Note that only the term with $`\omega ^2`$ from the electron-hole correlation function diverges in the limit $`\omega 0`$ while in the electron-electron (hole-hole) functions it remains finite even in the static limit. This is a manifestation of the diffusion pole missing the latter two correlation functions.

In Eq. (22) different powers of momentum $`𝐪`$ and frequency $`\omega `$ appear. There is no restriction on validity of Eqs.(22) in frequency but it holds only for small momenta, more precisely only perturbatively up to second order. Equations (22) establish relations between the density-density and the current-current correlation functions in the asymptotic limit $`q=0`$. Generally, however, these two functions are not directly proportional, since the right-hand sides of Eqs.(22) contain one-particle contributions. They cancel each other if we combine the electron-hole with electron-electron and hole-hole correlation functions appropriately to build up the density response function.

For isotropic situations we define in the hydrodynamic limit

$$g(\omega )=i\underset{q0}{lim}\frac{\omega }{q^2}\chi (𝐪,\omega ).$$
(24a)
This quantity can be generalized to anisotropic cases as
$$g_{\alpha \beta }(\omega )=i\omega \frac{^2}{q_\alpha q_\beta }\chi (𝐪,\omega )|_{q=0}.$$
(24b)
Function $`g(\omega )`$ measures leading long-range correlations of the density response function. Its real part can be identified with the diffusive conductivity or mobility of the system. It is now easy to find an explicit representation of this function using Eqs. (22). All contributions except for the last terms on the right-hand sides of Eqs. (22) cancel and we obtain equality $`\sigma _{\alpha \alpha }(\mathrm{𝟎},\omega )=e^2g_{\alpha \alpha }(\omega )`$ being just Eq. (4) in the limit $`𝐪=\mathrm{𝟎}`$. Note that in the course of the derivation we had to sum contributions from the electron-hole, electron-electron, and hole-hole correlation functions and integrate over energies.

If we resort to the low-frequency limit $`\omega 0`$, the diffusion pole in the electron-hole correlation function dominates in the density response and we recover the conductivity directly from $`\mathrm{\Phi }^{AR}`$. It is interesting to note that in this static limit we derive the mobility ($`dc`$-conductivity) with contributions from $`G^{AR}`$ and $`G^{RR}`$ solely from the electron-hole correlation function $`\mathrm{\Phi }^{AR}`$.

## IV Density response function calculated from the conductivity

We showed in the preceding section how the conductivity can be revealed from the small-momentum behavior of the density response function. In this section we express the current-current correlation function in terms of the density-density correlation function. The aid we use for this task is a continuity equation and another Ward identity

$$\begin{array}{c}\mathrm{\Sigma }(𝐤_+,z_+)\mathrm{\Sigma }(𝐤_{},z_{})=\frac{1}{N}\underset{𝐤^{}}{}\mathrm{\Lambda }_{\mathrm{𝐤𝐤}^{}}(z_+,z_{};𝐪)\hfill \\ \hfill \times \left[G(𝐤_+^{},z_+)G(𝐤_{}^{},z_{})\right].\end{array}$$
(25)
proved for retarded and advanced functions by Vollhardt and Wölfle in Ref. Vollhardt80b, . We denoted $`𝐤_\pm =𝐤\pm 𝐪/2`$. Note that identity (25) for $`𝐪=\mathrm{𝟎},z_+=E_F+i0^+,z_{}=E_Fi0^+`$ was already used in Ref. Maleev75, . Identity (25) is related to Eq. (11), however, the two Ward identities are identical neither in the derivation nor in the applicability and validity domains. The latter holds for nonzero transfer momentum $`𝐪`$, i. e., for an inhomogeneous perturbation while the former only for $`q=0`$. On the other hand, we show in the appendix that Eq. (25), unlike Eq. (11), holds only perturbatively within an expansion in powers of the random potential. I. e., we can prove this identity for nonzero transfer momentum $`𝐪`$ only by assuming that perturbation expansions for the self-energy $`\mathrm{\Sigma }`$ and simultaneously for the two-particle irreducible vertex $`\mathrm{\Lambda }`$ converge. It was shown earlier that in the homogeneous case, $`q=0`$, the Vollhardt-Wölfle-Ward identity follows from the Velický-Ward one.Janis01b

To derive a continuity equation for Green functions we start with an equation of motion for the two-particle Green function. It is a Bethe-Salpeter equation where the input is a two-particle irreducible vertex. We have three possibilities (topologically distinct scattering channels) how to construct a Bethe-Salpeter equation.Janis01b  For our purposes the electron-hole channel is the relevant or most suitable one. There the Bethe-Salpeter equation in momentum representation reads

$$G_{\mathrm{𝐤𝐤}^{}}^{(2)}(z_+,z_{};𝐪)=G(𝐤_+,z_+)G(𝐤_{},z_{})\left[\delta (𝐤𝐤^{})+\frac{1}{N}\underset{𝐤^{\prime \prime }}{}\mathrm{\Lambda }_{\mathrm{𝐤𝐤}^{\prime \prime }}(z_+,z_{};𝐪)G_{𝐤^{\prime \prime }𝐤^{}}^{(2)}(z_+,z_{};𝐪)\right]$$
(26)
where $`\mathrm{\Lambda }_{\mathrm{𝐤𝐤}^{}}(z_+,z_{};𝐪)`$ is the two-particle irreducible vertex from the electron-hole channel. One-electron propagators $`G(𝐤_\pm ,z_\pm )=[z_\pm \mathrm{\Sigma }(𝐤\pm 𝐪/2,z_\pm )ϵ(𝐤\pm 𝐪/2)]^1`$ are the averaged resolvents.

The product of the one-electron propagators from the right-hand side of Eq. (26) can be decomposed into

$$\begin{array}{c}G(𝐤_+,z_+)G(𝐤_{},z_{})=\hfill \\ \hfill \frac{\mathrm{\Delta }_𝐪G(𝐤;z_+,z_{})}{\mathrm{\Delta }z\mathrm{\Delta }_𝐪\mathrm{\Sigma }(𝐤;z_+,z_{})\mathrm{\Delta }_𝐪ϵ(𝐤)}\end{array}$$
(27)
where we denoted

$$\mathrm{\Delta }_𝐪ϵ(𝐤)=ϵ(𝐤_+)ϵ(𝐤_{}),$$
(28a)
$$\mathrm{\Delta }_𝐪G(𝐤;z_+,z_{})=G(𝐤_+,z_+)G(𝐤_{},z_{}).$$
(28b)
and analogously the difference $`\mathrm{\Delta }_𝐪\mathrm{\Sigma }(𝐤;z_+,z_{})`$.
We multiply both sides of the Bethe-Salpeter equation (26) by the denominator from the right-hand side of Eq. (27) and assume validity of Eq. (25) relating the two-particle irreducible vertex $`\mathrm{\Lambda }`$ with the one-particle irreducible self-energy $`\mathrm{\Sigma }`$. We then obtain a ”difference” equation of motion

$$\begin{array}{c}\left[\mathrm{\Delta }_𝐪ϵ(𝐤)\mathrm{\Delta }z\right]G_{\mathrm{𝐤𝐤}^{}}^{(2)}(z_+,z_{};𝐪)=\mathrm{\Delta }_𝐪G(𝐤;z_+,z_{})\delta (𝐤𝐤^{})+\frac{1}{N}\underset{𝐤^{\prime \prime }}{}\mathrm{\Lambda }_{\mathrm{𝐤𝐤}^{\prime \prime }}(z_+,z_{};𝐪)\hfill \\ \hfill \times \left[\mathrm{\Delta }_𝐪G(𝐤;z_+,z_{})G_{𝐤^{\prime \prime }𝐤^{}}^{(2)}(z_+,z_{};𝐪)G_{\mathrm{𝐤𝐤}^{}}^{(2)}(z_+,z_{};𝐪)\mathrm{\Delta }_𝐪G(𝐤^{\prime \prime };z_+,z_{})\right].\end{array}$$
(29)
We now define correlation functions generalized to complex frequencies describing density-density and density-current correlations by summing over the fermionic momenta in Eq. (29)

$`\mathrm{\Phi }(z_1,z_2;𝐪)`$ $`={\displaystyle \frac{1}{N^2}}{\displaystyle \underset{\mathrm{𝐤𝐤}^{}}{}}G_{\mathrm{𝐤𝐤}^{}}^{(2)}(z_1,z_2;𝐪),`$ (30a)
$`\mathrm{\Phi }_ϵ(z_1,z_2;𝐪)`$ $`={\displaystyle \frac{1}{N^2}}{\displaystyle \underset{\mathrm{𝐤𝐤}^{}}{}}\mathrm{\Delta }_𝐪ϵ(𝐤)G_{\mathrm{𝐤𝐤}^{}}^{(2)}(z_1,z_2;𝐪),`$ (30b)
$`\overline{\mathrm{\Phi }}_ϵ(z_1,z_2;𝐪)`$ $`={\displaystyle \frac{1}{N^2}}{\displaystyle \underset{\mathrm{𝐤𝐤}^{}}{}}G_{\mathrm{𝐤𝐤}^{}}^{(2)}(z_1,z_2;𝐪)\mathrm{\Delta }_𝐪ϵ(𝐤^{}).`$ (30c)
The contributions from the two-particle irreducible vertex $`\mathrm{\Lambda }`$ in Eq. (29) cancel each other provided the two-particle irreducible vertex is symmetric, i. e., $`\mathrm{\Lambda }_{\mathrm{𝐤𝐤}^{}}=\mathrm{\Lambda }_{𝐤^{}𝐤}`$. If so, we end up with a continuity equation relating the generalized density-density and density-current correlation functions

$$\begin{array}{c}\mathrm{\Phi }_ϵ(z_+,z_{};𝐪)\mathrm{\Delta }z\mathrm{\Phi }(z_+,z_{};𝐪)\hfill \\ \hfill =\frac{1}{N}\underset{𝐤}{}\mathrm{\Delta }_𝐪G(𝐤;z_+,z_{}).\end{array}$$
(31a)
Another continuity equation can be derived by multiplying Eq. (29) with the energy difference and summing over the fermionic momenta. We obtain an equation relating the generalized current-current and density-current correlation functions
$$\begin{array}{c}\mathrm{\Phi }_{ϵϵ}(z_+,z_{};𝐪)\mathrm{\Delta }z\overline{\mathrm{\Phi }}_ϵ(z_+,z_{};𝐪)\hfill \\ \hfill =\frac{1}{N}\underset{𝐤}{}\mathrm{\Delta }_𝐪G(𝐤;z_+,z_{})\mathrm{\Delta }_𝐪ϵ(𝐤).\end{array}$$
(31b)
Combining the above two continuity equations we obtain the resulting relation between the generalized current-current and density-density correlation functions replacing the operator continuity equation (3)
$$\begin{array}{c}\mathrm{\Phi }_{ϵϵ}(z_+,z_{};𝐪)(\mathrm{\Delta }z)^2\mathrm{\Phi }(z_+,z_{};𝐪)\hfill \\ \hfill =\frac{1}{N}\underset{𝐤}{}\mathrm{\Delta }_𝐪G(𝐤;z_+,z_{})\left[\mathrm{\Delta }z+\mathrm{\Delta }_𝐪ϵ(𝐤)\right].\end{array}$$
(31c)
Actually, correlation function $`\mathrm{\Phi }_{ϵϵ}`$ is strictly speaking not the current-current correlation function, since it is not the momentum or velocity that is summed with the two-particle Green function. We used an energy difference $`\mathrm{\Delta }_𝐪ϵ(𝐤)`$ from Eq. (28a). Only in the case of quadratic dispersion relation we have $`\mathrm{\Delta }_𝐪ϵ(𝐤)=𝐪𝐤/m`$ and the energy difference is proportional to the group velocity. Otherwise $`\mathrm{\Phi }_{ϵϵ}`$ equals the current-current correlation function only in the small-momentum limit, $`q0`$. We use

$$\mathrm{\Delta }_𝐪ϵ(𝐤)𝐪𝐯(𝐤)$$
(32)
to convert the energy difference in the hydrodynamic limit to a multiple of the group velocity $`𝐯(𝐤)`$ and denote

$$\begin{array}{c}\mathrm{\Phi }_{\alpha \beta }(z_+,z_{};𝐪)\hfill \\ \hfill =\frac{1}{N^2}\underset{\mathrm{𝐤𝐤}^{}}{}v^\alpha (𝐤)v^\beta (𝐤^{})G_{\mathrm{𝐤𝐤}^{}}^{(2)}(z_+,z_{};𝐪).\end{array}$$
(33)
With the aid of Eq. (33) we rewrite continuity equation (31c) to

$$\begin{array}{c}\underset{\alpha \beta }{}q^\alpha q^\beta \mathrm{\Phi }_{\alpha \beta }(z_+,z_{};𝐪)\underset{\alpha }{}q^\alpha v^\alpha \mathrm{\Delta }_𝐪G(z_+,z_{})\hfill \\ \hfill =\mathrm{\Delta }z\left[\mathrm{\Delta }z\mathrm{\Phi }(z_+,z_{};𝐪)+\mathrm{\Delta }_𝐪G(z_+,z_{})\right]\end{array}$$
(34)
where the angular brackets stand for the summation over the fermionic momenta from the first Brillouin zone as defined in Eq. (23). In the isotropic case we obtain

$$\begin{array}{c}\mathrm{\Phi }_{\alpha \alpha }(z_+,z_{};𝐪)=\frac{1}{q^2}𝐪𝐯\mathrm{\Delta }_𝐪G(z_+,z_{})\hfill \\ \hfill +\frac{\mathrm{\Delta }z}{q^2}\left[\mathrm{\Delta }z\mathrm{\Phi }(z_+,z_{};𝐪)+\mathrm{\Delta }_𝐪G(z_+,z_{})\right].\end{array}$$
(35)
It is the most general relation between the current-current and density-density correlation functions for arbitrary complex frequencies and non-zero momenta. It holds for finite momenta as far as the dispersion law remains quadratic. For general dispersion relations the validity of Eq. (35) remains in the asymptotic regime of small momenta.

To return to measurable quantities we limit the complex frequencies to the real axis. To derive exact relations we first define a generalization of the conductivity with the correlation functions $`\mathrm{\Phi }_{ϵϵ}^{ab}`$:

$$\begin{array}{c}\sigma _{ϵϵ}(𝐪,\omega )=e^2_{\mathrm{}}^{\mathrm{}}\frac{dE}{2\pi \omega }\{[f(E+\omega )f(E)]\mathrm{\Phi }_{ϵϵ}^{AR}(E,E+\omega ;𝐪)+f(E)\mathrm{\Phi }_{ϵϵ}^{RR}(E,E+\omega ;𝐪)\hfill \\ \hfill f(E+\omega )\mathrm{\Phi }_{ϵϵ}^{AA}(E,E+\omega ;𝐪)f(E)[\mathrm{\Phi }_{ϵϵ}^{RR}(E,E;𝐪)\mathrm{\Phi }_{ϵϵ}^{AA}(E,E;𝐪)]\}.\end{array}$$
(36)
We use continuity equation (31c) to relate $`\sigma _{ϵϵ}`$ and $`\chi `$ and find

$$\begin{array}{c}\sigma _{ϵϵ}(𝐪,\omega )+ie^2\omega \chi (𝐪,\omega )=\frac{e^2}{\omega }_{\mathrm{}}^{\mathrm{}}\frac{dE}{2\pi }f(E)\left[\mathrm{\Phi }_{ϵϵ}^{RR}(E,E;𝐪)\mathrm{\Phi }_{ϵϵ}^{AA}(E,E;𝐪)\right]\hfill \\ \hfill \frac{ie^2}{\omega }_{\mathrm{}}^{\mathrm{}}\frac{dE}{\pi }f(E)\frac{1}{N}\underset{𝐤}{}\mathrm{\Delta }_𝐪ϵ(𝐤)\left[\mathrm{}G^R(𝐤_+,E)\mathrm{}G^R(𝐤_{},E)\right]=0.\end{array}$$
(37)
Vanishing of the right-hand side in Eq. (37) can be explicitly manifested when continuity equation (31c) is applied to the two-particle correlation functions $`\mathrm{\Phi }_{ϵϵ}^{RR}`$ and $`\mathrm{\Phi }_{ϵϵ}^{AA}`$. The two terms on the right-hand side of Eq. (37) add to zero.

We can extract the electrical conductivity from $`\sigma _{ϵϵ}(𝐪,\omega )`$ in the limit $`q0`$ or for quadratic dispersion law by replacing $`\mathrm{\Phi }_{ϵϵ}=_\alpha q_\alpha ^2\mathrm{\Phi }_{\alpha \alpha }`$ leading to $`\sigma _{ϵϵ}(𝐪,\omega )=_\alpha q_\alpha ^2\sigma _{\alpha \alpha }(𝐪,\omega )`$. In the isotropic case we then reveal Eq. (4) from Eq. (37). It is important to state that this proof of validity of Eq. (4) is strongly based on the Ward identity (25) for nonzero momenta. The latter can be proved only perturbatively, see Appendix, and hence it is not clear whether Eq.(4) holds beyond the perturbative regime near and beyond the Anderson localization transition.

## V Diffusion, diffusion pole, and Einstein relation in quantum systems

In the preceding sections we proved validity of Eq. (4) for arbitrary frequencies and small momenta from the region where quadratic dispersion relation is accurate in strongly disordered lattice systems. It is a useful relation in particular in approximate theories. It enables one to approximate only either the density response function or the conductivity and to determine the other one from Eq. (4). However, both quantities in Eq. (4) are response functions that are not well suited for systematic approximations. We would prefer exact relations or identities for equilibrium correlation or Green functions. This can be achieved when we introduce diffusion. Diffusion in classical systems is defined by Fick’s law relating the current with the negative gradient of the charge density.Kubo57 ; Enz92  The diffusion constant introduced in this way, i. e. via Fick’s law, has in classical physics an intuitive phenomenological character without a proper microscopic backing that should come from quantum physics. The study of diffusion in quantum disordered systems was launched by the seminal work of Anderson analyzing destructive effects of quantum coherence on diffusion of disordered electrons.Anderson58  Since then diffusion has become part of transport studies of disordered and interacting quantum itinerant systems.Altshuler85 ; Kramer93 ; Belitz94  However even in these quantum approaches diffusion is introduced either via (self-consistent) perturbation expansions leading to the diffusion pole in the electron-hole correlation function or via the semiclassical limit where the electron-hole correlation function becomes a Green function of the classical diffusion equation.Rammer98

Here we propose a nonperturbative way to define diffusion from first principles. We introduce a quantum generalization of Fick’s law with a diffusion response function relating the averaged current density with the negative gradient of the charge density of a perturbed system. We hence use

$$\begin{array}{c}𝐣(𝐱,t)_{av}=e_{\mathrm{}}^{\mathrm{}}𝑑t^{}d^dx^{}D(𝐱𝐱^{},tt^{})\hfill \\ \hfill \times \mathbf{}\delta \stackrel{~}{n}(𝐱^{},t^{})_{av}\end{array}$$
(38a)
as a definition of the diffusion response function $`D(𝐱,t)`$ being an extension of the diffusion constant to which it should reduce in the semiclassical limit. In the above definition we used a perturbed charge density $`\delta \stackrel{~}{n}(𝐱,t)`$ that is not the total change of the charge density $`\delta n(𝐱,t)=n(𝐱,t)n_0`$. It is a “dynamical” density change obeying a condition
$$_{\mathrm{}}^{\mathrm{}}𝑑t\delta \stackrel{~}{n}(𝐱,t)_{av}=0$$
(38b)
saying that the total change at any place during the entire time evolution is zero. This restriction is dictated by the same condition obeyed by the averaged current density defined from the conductivity (14).

Definition of diffusion (38) is nonperturbative, has quantum character and does not demand any assumption about the behavior in the semiclassical limit or the existence of the diffusion pole in the electron-hole correlation function. Moreover, it offers a microscopic determination of the diffusion function from first principles, i. e., quantum response functions of the perturbed system and via the Kubo formalism from correlation and Green functions of an equilibrium solution. It is hence directly related to measurable quantities.

In Fourier representation we obtain for the averaged current

$$𝐣(𝐪,\omega _{av}=ie𝐪D(𝐪,\omega )\delta \stackrel{~}{n}(𝐪,\omega )_{av}.$$
(39)
The dynamical charge density $`\delta \stackrel{~}{n}`$ is induced by an external dynamical scalar potential, which can be expressed in Fourier representation as

$$\delta \stackrel{~}{n}(𝐪,\omega )_{av}=e\left[\chi (𝐪,\omega )\chi (𝐪,0)\right]\phi (𝐪,\omega ).$$
(40)
It is easy to verify that with this definition Eq. (38b) is fulfilled. We use Eq. (1) together with Eqs. (39) and (40) to find an identity between the conductivity, diffusion and density response

$$\sigma (𝐪,\omega )=e^2D(𝐪,\omega )\left[\chi (𝐪,\omega )\chi (𝐪,0)\right].$$
(41)
We have derived a general relation between quantum diffusion and the electrical conductivity. It is now to show that quantum diffusion response function $`D(𝐪,\omega )`$ behaves in the semiclassical limit as expected from the diffusion constant in the diffusion equation. Namely, we have to show that the homogeneous static diffusion $`D(\mathrm{𝟎},0)`$ obeys Einstein’s relation and enters the diffusion equation in the semiclassical limit. To come to it we use Eq. (4) to exclude the density response $`\chi (𝐪,\omega )`$ from Eq. (41) and Eq. (13a) to determine $`\chi (𝐪,0)`$. Doing so we obtain in the isotropic case

$$\begin{array}{c}\left[1+\frac{iq^2}{\omega }D(𝐪,\omega )\right]\sigma (𝐪,\omega )\hfill \\ \hfill =e^2D(𝐪,\omega )_{\mathrm{}}^{\mathrm{}}\frac{dE}{\pi }f(E)\mathrm{}\mathrm{\Phi }_E^{RR}(𝐪,0).\end{array}$$
(42)
If we define the homogeneous dynamical conductivity $`\sigma (\omega )=\sigma (\mathrm{𝟎},\omega )`$ (in the same way also the dynamical diffusion) and use the Velický-Ward identity (11) in the limiting case $`z_1z_20`$ we obtain a dynamical generalization of the well known Einstein relation

$$\begin{array}{c}\sigma (\omega )=e^2D(\omega )_{\mathrm{}}^{\mathrm{}}\frac{dE}{\pi }f^{}(E)\mathrm{}G^R(E)\hfill \\ \hfill =e^2D(\omega )\left(\frac{n}{\mu }\right)_T\end{array}$$
(43a)
that at zero temperature reduces to
$$\sigma (\omega )=e^2n_FD(\omega ).$$
(43b)
A classical version of this formula was proved by Einstein for the Brownian motion of a particle in a random mediumEinstein56  and later on it was re-derived in the framework of nonequilibrium statistical mechanics by Kubo.Kubo57  Here we derived Einstein’s relation for quantum response functions and showed that it is a consequence of gauge invariance of the system, namely Eq. (4).

It is not yet clear from the above reasoning how the diffusion function $`D(𝐪,\omega )`$ is related to the classical diffusion equation and to the diffusion pole. The existence of the diffusion pole in quantum systems is usually deduced from the semiclassical limit leading to the diffusion equation. We can prove the existence of the diffusion pole entirely from first quantum principles and from the quantum generalization of Fick’s law, Eq. (38). We use Eq. (4) to exclude now conductivity $`\sigma (𝐪,\omega )`$ from Eq. (41). It is easy to come to a representation of the density response function

$$\chi (𝐪,\omega )=\frac{D(𝐪,\omega )\chi (𝐪,0)q^2}{i\omega +D(𝐪,\omega )q^2}$$
(44)
known from other treatments of diffusion.Forster90  Eq. (44) holds for arbitrary frequencies and momenta within the range of quadratic dispersion relation. This representation directly indicates a singularity in the density response function in the limit $`\omega 0,q0`$. Since the order of the limits is relevant for the result, we have to specify it explicitly. To single out the singular contribution in the density response function we have to choose $`\omega /q1`$, that is, the inverse situation to Sec. III and to Eqs. (43).

The singularity in the density response function leads to a pole in the two-particle Green function. It is demonstrated in the low-energy and small-momentum asymptotics of the density response function at zero temperature. We find from Eq. (10)

$$\begin{array}{c}\chi (𝐪,\omega )=\chi (𝐪,0)+\frac{i\omega }{2\pi }\left(\mathrm{\Phi }_{E_F}^{AR}(𝐪,0)+O(q^0)\right)\hfill \\ \hfill +O(\omega ^2).\end{array}$$
(45)
Using this asymptotics in Eq. (44) we obtain an explicit manifestation of the diffusion pole in the zero-temperature electron-hole correlation function

$$\mathrm{\Phi }_{E_F}^{AR}(𝐪,\omega )\frac{2\pi n_F}{i\omega +Dq^2}$$
(46)
representing the singular part of the electron-hole correlation function. We denoted a diffusion constant $`D=lim_{q0}D(𝐪,0)`$. This representation was derived for zero temperature in the asymptotic limit $`q0,\omega 0`$ with the restriction $`\omega /q1`$. To extend it to the opposite limit, $`q/\omega 1`$, we have to show that the diffusion function $`D(𝐪,\omega )`$ is analytic in the limit $`q0,\omega 0`$. It is not a priori clear that the diffusion constant in Eq. (46) equals the diffusion constant obtained from the inverse order of limits $`D=lim_{\omega 0}lim_{q0}D(𝐪,\omega )`$ used in the Einstein relation, Eq. (43b). However, the diffusion response function $`D(𝐪,\omega )`$ was introduced in such a way that the limits $`\omega 0`$ and $`q0`$ commute and we have only one definition for the static diffusion constant from the density response function

$$\begin{array}{c}D=\frac{i}{2}\frac{lim_{\omega 0}_q^2\chi (𝐪,\omega )|_{q=0}}{lim_{q0}_\omega \chi (𝐪,\omega )|_{\omega =0}}\hfill \\ \hfill =\frac{2\pi n_F}{lim_{q0}q^2\mathrm{\Phi }^{AR}(𝐪,0)}.\end{array}$$
(47)
The left equality in Eq. (47) holds quite generally while the right one only at zero temperature. Note that representation (46) holds only at zero temperature for small momenta and in the low-frequency limit.

We can see from Eq. (46) that the diffusion function introduced in the quantum version of Fick’s law reduces to the diffusion constant determining the long-range fluctuations of the electron-hole correlation function and hence in the semiclassical limit also the diffusive behavior in the diffusion equation. However, only the static diffusion constant $`D`$ enters the denominator of the electron-hole correlation function. A generalization of the diffusion constant from the electron-hole correlation function to a frequency-dependent quantity via

$$n_F\stackrel{~}{D}(\omega )=\frac{\omega ^2}{4\pi }_q^2\mathrm{\Phi }_{E_F}^{AR}(𝐪,\omega )|_{q=0},$$
(48)
does not have a rigorous justification. It is nevertheless often used in the literature for a quantitative treatment of Anderson localization.Vollhardt92  Beware that the dynamical diffusion defined in this way does not fulfill the Einstein relation and hence deviates from the dynamical conductivity. We show in the next section that in an exactly solvable mean-field limit the left-hand side of Eq. (48) related to the conductivity via the Einstein relation differs at finite frequencies from the conductivity obtained from the Kubo formula with the electrical current.

## VI Infinite-dimensional model: Explicit exact solution

We derived exact relations between density and current correlation functions for a lattice electron gas in a random potential. We now demonstrate the generally derived formulas explicitly on an exactly solvable limit of infinite spatial dimensions. This limit serves as a mathematical tool for the definition of a mean-field theory not only for classical spin systems but also for itinerant disordered and interacting models.Janis89 ; Janis92 ; Georges96  In case of the Anderson model of disordered electrons the mean-field theory, i. e., the limit of infinite spatial dimensions, equals the coherent potential approximation. We hence resort in this section to this solution.

Equation determining the local self-energy in the coherent-potential approximation is Soven’s equation that can be written asJanis01b

$$G(z)=\left[G^1(z)+\mathrm{\Sigma }(z)V_i\right]^1_{av}$$
(49a)
where $`G(z)=N^1_𝐤G(𝐤,z)`$. The two-particle irreducible vertex then is
$$\begin{array}{c}\lambda (z_+,z_{})=\frac{1}{G(z_+)G(z_{})}\left[1\frac{1}{1+\left(\mathrm{\Sigma }(z_+)V_i\right)G(z_+)}\frac{1}{1+\left(\mathrm{\Sigma }(z_{})V_i\right)G(z_{})}_{av}^{}{}_{}{}^{1}\right]\hfill \\ \hfill =\frac{\mathrm{\Sigma }(z_+)\mathrm{\Sigma }(z_{})}{G(z_+)G(z_{})}=\frac{\mathrm{\Delta }\mathrm{\Sigma }}{\mathrm{\Delta }G}.\end{array}$$
(49b)
The angular brackets were defined in Eq. (23). The second equality in Eq. (49b) is the Velický-Ward identity.

It is straightforward to find an explicit form of the two-particle Green function in the coherent-potential approximation

$$\begin{array}{c}G_{\mathrm{𝐤𝐤}^{}}^{(2)}(z_+,z_{};𝐪)=G_+(𝐤)G_{}(𝐤)[\delta (𝐤𝐤^{})\hfill \\ \hfill +\frac{\lambda (z_+,z_{})G_+(𝐤^{})G_{}(𝐤^{})}{1\lambda (z_+,z_{})G_+G_{}}]\end{array}$$
(50)
where $`G_\pm (𝐤)=G(𝐤\pm 𝐪/2,z_\pm )`$ and $`G_+G_{}=N^1_𝐤G_+(𝐤)G_{}(𝐤)`$.

Using the two-particle Green function, Eq. (50), we easily obtain explicit representations for the density-density and density-current correlation functions

$`\mathrm{\Phi }(z_+,z_{};𝐪)={\displaystyle \frac{G_+G_{}}{1\lambda (z_+,z_{})G_+G_{}}},`$ (51a)
$`\mathrm{\Phi }_ϵ(z_+,z_{};𝐪)={\displaystyle \frac{\mathrm{\Delta }_𝐪ϵG_+G_{}}{1\lambda (z_+,z_{})G_+G_{}}}.`$ (51b)
We now insert Eqs. (51) together with the explicit representation for the two-particle Green function, Eq. (50), to the right-hand side of Eq. (31a) and find

$$\begin{array}{c}\mathrm{\Phi }_ϵ(z_+,z_{};𝐪)\mathrm{\Delta }z\mathrm{\Phi }(z_+,z_{};𝐪)\hfill \\ \hfill =\frac{\left(\mathrm{\Delta }_𝐪ϵ\mathrm{\Delta }z\right)\frac{\mathrm{\Delta }_𝐪G}{\mathrm{\Delta }_𝐪ϵ\mathrm{\Delta }z\mathrm{\Delta }\mathrm{\Sigma }}}{1\lambda (z_+,z_{})G_+G_{}}\\ \hfill =\frac{\mathrm{\Delta }G\mathrm{\Delta }\mathrm{\Sigma }G_+G_{}}{1\lambda (z_+,z_{})G_+G_{}}=\mathrm{\Delta }G,\end{array}$$
(52)
that is, continuity equation (31a) is fulfilled. Analogously we can explicitly verify Eqs. (31b) and (31c). Actually, validity of continuity equations (31) follows from the Velický-Ward identity, since the two-particle irreducible vertex $`\lambda (z_+,z_{})`$ is momentum independent. In such a case the Vollhardt-Wölfle and the Velický identities are equivalent.

As next we show that the $`dc`$-conductivity in the coherent-potential approximation can be calculated using the diffusion constant $`D`$ from Eq. (46). To this purpose we need to evaluate $`_q^2G_+G_{}`$ at zero transfer frequency $`\omega =0`$. We have

$$\begin{array}{c}_q^2G_+G_{}=\frac{1}{2}_q\left(G_{}_kG_+G_+_kG_{}\right)\hfill \\ \hfill =\frac{1}{4}G_{}_k^2G_++G_+_k^2G_{}2_kG_+_kG_{}\\ \hfill =_kG_+_kG_{}\end{array}$$
(53)
where in the last equality we used integration per parts in momentum space. We utilize this result and Eq. (11) to obtain the electron-hole correlation function

$$\mathrm{\Phi }(𝐪,\mathrm{\Delta }z)=\frac{\mathrm{\Delta }G}{\mathrm{\Delta }z{\displaystyle \frac{1}{2}}{\displaystyle \frac{\mathrm{\Delta }\mathrm{\Sigma }}{\mathrm{\Delta }\mathrm{\Sigma }\mathrm{\Delta }z}}{\displaystyle \frac{v_k^2(\mathrm{\Delta }G)^2}{\mathrm{\Delta }G}}q^2},$$
(54)
that in the low-energy limit reduces for real frequencies to

$$\mathrm{\Phi }^{AR}(𝐪,\omega )=\frac{2\pi n_F}{i\omega +\sigma ^{CPA}q^2/e^2n_F}$$
(55)
proving the Einstein relation for the $`dc`$-conductivity $`\sigma ^{CPA}`$ of the coherent-potential approximation. It is evident from Eq. (54) that the frequency-dependent coefficient in the denominator of the electron-hole correlation function does not equal the dynamical diffusion from Eq. (43b). To show this explicitly we use definition (48) and compare the result with the real-part of the conductivity at zero temperature given by

$$\begin{array}{c}\mathrm{}\sigma _{\alpha \alpha }(\omega )=e^2_{E_F\omega }^{E_F}\frac{dE}{2\pi \omega }[\mathrm{}\mathrm{\Phi }_{\alpha \alpha }^{AR}(E,E+\omega ;\mathrm{𝟎})\hfill \\ \hfill \mathrm{}\mathrm{\Phi }_{\alpha \alpha }^{RR}(E,E+\omega ;\mathrm{𝟎})]=\sigma ^{AR}_{\alpha \alpha }(\omega )+\sigma ^{RR}_{\alpha \alpha }(\omega ).\end{array}$$
(56a)
We separated contributions to the conductivity from the electron-hole and electron-electron current-current correlation functions, $`\sigma ^{AR}`$, $`\sigma ^{RR}`$. It is useful to introduce also a simplified asymptotic form of the low-frequency conductivity
$$\begin{array}{c}\mathrm{}\overline{\sigma }_{\alpha \alpha }(\omega )=e^2\frac{1}{2\pi }[\mathrm{}\mathrm{\Phi }_{\alpha \alpha }^{AR}(E_F\omega ,E_F;\mathrm{𝟎})\hfill \\ \hfill \mathrm{}\mathrm{\Phi }_{\alpha \alpha }^{RR}(E_F\omega ,E_F;\mathrm{𝟎})]\end{array}$$
(56b)
where the integrand in the energy integral is replaced by its initial value. This representation is asymptotically exact in the limit $`\omega 0`$ for smooth current-current correlation functions $`\mathrm{\Phi }_{\alpha \alpha }`$ near the Fermi energy $`E_F`$ and carries the same frequency dependence as the diffusion constant $`\stackrel{~}{D}(\omega )`$ defined in Eq. (48).

For the numerical calculations we use a binary alloy with two values of the random potential $`V_i=\pm \mathrm{\Delta }`$ weighted with probability $`x`$ and $`1x`$. Fig. 1 shows the two dynamical conductivities and the diffusion constant from Eq. (48) on a simple cubic lattice with parameters $`x=0.1`$, $`\mathrm{\Delta }=0.6w`$ for frequencies $`\omega =0.2w`$ and $`\omega =0.4w`$ with $`w`$ being the half-bandwidth. Note that the split-band value of the disorder strength is $`\mathrm{\Delta }_c0.4w`$.

Although the three quantities differ for low frequencies only slightly inside the band, they behave differently near the band edges and for higher frequencies. In particular the diffusion constant shows anomalous behavior when the Fermi energy approaches a band edge. The compensating terms from the electron-electron (density-density) correlation function get relevant there. Anomalous behavior of the diffusion constant is more transparent for higher frequencies (lower pane). We cannot evidently rely on $`n_F\stackrel{~}{D}(\omega )`$ as a good approximation to the conductivity for finite frequencies except for Fermi energies deep inside the energy band. The smoothing impact of the integral over frequencies on the behavior of the dynamical conductivity gets clear from our numerical results.

Recently a discussion was renewed about the proportion of contributions to the Kubo formula for the electrical conductivity from the electron-hole and electron-electron current-current correlation functions, $`\sigma ^{AR}`$ and $`\sigma ^{RR}`$ in Eq. (56a), respectively.Nikolic01  In Fig. 2 these contributions are compared for the same setting of the binary alloy on a simple cubic lattice for $`\omega =0`$. The contribution from $`\sigma ^{AR}`$ dominates inside the band far from the band edges. Outside the central band and in the satellite impurity band the importance of the compensating effects of $`\sigma ^{RR}`$ is evident.

The situation worsens when we go over to the dynamical conductivity, Fig. 3. We can see that there is no region where the electron-hole contribution would dominate or approximate the full conductivity reliably.

It is important to note that conductivity $`\sigma ^{AR}`$ differs from the (dynamical) diffusion constant defined from the electron-hole density-density correlation function, Eq. (48). The latter contains even in the static limit both terms, $`\sigma ^{AR}`$ and $`\sigma ^{RR}`$. The difference is evident in Fig. 3. The conductivity calculated from $`\stackrel{~}{D}(\omega )`$ is much closer to the full electrical conductivity than $`\sigma ^{AR}`$. With increasing frequency the role of the term $`\sigma ^{RR}`$ increases. We hence cannot interchange contributions to the electrical conductivity from the electron-hole density-density correlation function $`\mathrm{\Phi }^{AR}(𝐪,\omega )`$ (its leading $`q`$-dependent term) and the conductivity $`\sigma ^{AR}(𝐪,\omega )`$. The former is generally better approximation than the latter except for Fermi energies near band edges and in band tails where both the approximations equally fail.

## VII Discussion and conclusions

In this paper we derived a number general relations between density-density, density-current and current-current correlation functions. With the aid of these relations we proved validity of an equation relating the electrical conductivity to the density response function, Eq. (4), in the hydrodynamic limit of strongly disordered electron lattice systems and derived a dynamical generalization of the Einstein relation between the conductivity and diffusion, Eq. (43). We found to what extent Eq. (4) can be used in quantitative calculations and when it may be broken. We showed that the dynamical diffusion from the Einstein relation cannot be related to the diffusion pole beyond the static limit. We used only exact reasoning and equalities derived directly from equations of motion for Green functions. We avoided extrapolations of semiclassical results or conclusions derived within equilibrium theory.

The electrical conductivity can be related to the density response function only if Ward identities are fulfilled. The same holds for the existence of a diffusion pole in the electron-hole correlation function with the static optical conductivity governing its long-range spatial fluctuations. There are two, not fully equivalent, Ward identities for noninteracting electrons in a random potential. Although they reflect conservation laws they can be proved only if perturbation expansion in the strength of the random potential for the irreducible two-particle vertex converges and results in an analytic function.

In Sec. III we used Velický-Ward identity (11) to relate mobility with the density response. This Ward identity expresses conservation of probability in the state space and holds only for the homogeneous case, i. e., zero transfer momentum in the two-particle Green function. To be of use in the calculation of the conductivity we had to assume analyticity of the hydrodynamic limit for all frequencies and to use the transfer momentum as an expansion parameter. Vollhardt-Wölfle-Ward identity (25) used in Sec. IV to reveal the density response function from the conductivity expresses charge conservation. It holds, unlike the Velický-Ward identity, for arbitrary transfer momenta in the two-particle function, but it can be proved only via a perturbation (diagrammatic) expansion for the two-particle vertex function. Moreover, to warrant this identity we are not allowed to sum selected classes of relevant diagrams. A single relevant diagram for the two-particle irreducible vertex of order $`n`$ produces $`n1`$ irrelevant diagrams in the proof of the Ward identity. Selection of relevant diagrams may be, however, dictated by the existence of a singularity in the two-particle function or by analyticity (causality) of the theory. There is hence no general guarantee that Ward identities and consequently the relation between the conductivity and the density response hold beyond the perturbative regime near the Anderson localization transition.

When the conductivity can be related to the density response via Eq. (4) we can derive a quantum version of the Einstein relation. To that purpose we introduced a diffusion response function relating the averaged current with the negative gradient of the charge density as a quantum generalization of Fick’s law. Diffusion defined in this way enters the dynamical Einstein relation, Eq. (43), but it is not à priori clear whether this quantum diffusion reduces to the classical one in the semiclassical limit. The solution of the classical diffusion equation pops up from a quantum theory via the diffusion pole and the long-range asymptotics of the electron-hole correlation function. We showed that the “classical” diffusion from the diffusion pole equals the quantum one from the Einstein relation and hence the conductivity only in the static limit. A dynamical diffusion constant defined from the electron-hole correlation function is no longer proportional to the conductivity beyond the static limit. This was demonstrated generally as well as explicitly on an exact solution in infinite spatial dimensions.

To conclude, we showed that gauge invariance of an electron gas in a random potential responding to an electromagnetic perturbation is guaranteed within linear-response theory if Ward identities (11) and (25) are fulfilled. In this case an (approximate) averaged two-particle Green function generates density and current response functions consistent in the hydrodynamic limit with relation (4). Or, what is more important, we can safely apply the approximate two-particle Green function only in either Kubo formula and use Eq. (4) to determine the other response function. We explicitly demonstrated that Ward identities are fulfilled only in the “diffusive” regime characterized by convergence of the perturbation expansion for the two-particle irreducible vertex and by analyticity of the hydrodynamic limit for all frequencies. We, however, know that in solutions with weakly or strongly localized electrons a pole in the two-particle irreducible vertex emerges and validity of the Ward identities cannot be generally proved. One has to bear this aspect in mind when interpolating between the diffusive and localized regimes.

###### Acknowledgements.

This work was supported in part by Grant No. 202/01/0764 of the Grant Agency of the Czech Republic.

*

## Appendix A Proof and validity of the Vollhardt-Wölfle-Ward identity

Vollhardt-Wölfle-Ward identity (25) is a primary tool for proving the continuity equations relating the current and density correlation functions. Understanding its proof and the domain of its validity is hence very important. The original proof of Vollhardt and WölfleVollhardt80b  is restricted to real frequencies with different small imaginary parts only. However, it is valid for arbitrary complex frequencies. To demonstrate this we summarize the principal assumptions and steps of its proof from which we conclude on the domain of validity of Eq. (25).

First, we have to assume that the self-energy $`\mathrm{\Sigma }`$ is a functional of the full one-electron propagator $`G`$

$$\mathrm{\Sigma }(𝐩,z)=\mathrm{\Sigma }\left[G\right](𝐩,z)$$
(57)
so that the perturbation expansion in the random potential consists of irreducible diagrams only. We consider here only noninteracting electrons in a random potential although the reasoning can be extended to correlated electrons as well. Perturbation expansion of the self-energy $`\mathrm{\Sigma }(𝐩,z)`$ in the random potential reads

$$\begin{array}{c}\mathrm{\Sigma }(𝐩,z)=V_{av}+\underset{n=1}{\overset{\mathrm{}}{}}\underset{k=1}{\overset{[(n+1)/2]}{}}\underset{j_1,\mathrm{},j_k=1}{\overset{n+1}{}}\underset{l_1,\mathrm{},l_k}{}\frac{1}{l_1!\mathrm{}l_k!}\delta (l_1+\mathrm{}l_kn1)V_{j_1}^{l_1}_c\mathrm{}V_{j_k}^{l_k}_c\hfill \\ \hfill \underset{𝐑_1,\mathrm{},𝐑_{n+1}}{}\underset{𝒫\{i_1,\mathrm{},i_{n+1}\}}{}\delta _{j_1,i_1}\mathrm{}\delta _{j_1,i_{l_1}}\delta _{j_2,i_{l_1+1}}\mathrm{}\delta _{j_k,i_{n+1}}\\ \hfill \frac{1}{N^n}\underset{𝐩_1,\mathrm{},𝐩_n}{}e^{i𝐑_{i_1}(𝐩𝐩_1)}e^{i𝐑_{i_2}(𝐩_1𝐩_2)}e^{i𝐑_{i_{n+1}}(𝐩_n𝐩)}G(𝐩_1,z)\mathrm{}G(𝐩_n,z)\end{array}$$
(58)
where $`[x]`$ is the integer part of $`x`$, $`𝒫\{i_1,\mathrm{},i_{n+1}\}`$ denotes a permutation of the indices $`\{1,2,\mathrm{},n+1\}`$, and the angular brackets stand for cumulant averages defined from

$$\mathrm{exp}V_{av}=\mathrm{exp}\left\{\underset{n=1}{\overset{\mathrm{}}{}}\frac{1}{n!}V^n_c\right\}.$$
(59)
Only the momentum variables and momentum-dependent functions are of importance. We denote $`𝒟_n`$ the sum of all diagrams with $`n`$ internal fermionic lines with internal momenta $`𝐩_1,\mathrm{},𝐩_n`$. Then representation (58) can be simplified to

$$\mathrm{\Sigma }(𝐩,z)=V_{av}+\underset{n=1}{\overset{\mathrm{}}{}}\frac{1}{N^n}\underset{𝐩_1,\mathrm{},𝐩_n}{}\underset{𝒟_n}{}X_{𝒟_n}(V;𝐩,𝐩_1,\mathrm{},𝐩_n)G(𝐩_1,z)\mathrm{}G(𝐩_n,z).$$
(60)
We use notation $`G_\pm (𝐩)=G(𝐩\pm 𝐪/2,z_\pm )`$ and $`\mathrm{\Sigma }_\pm (𝐩)=\mathrm{\Sigma }(𝐩\pm 𝐪/2,z_\pm )`$ and apply expansion (60) for the self-energy difference. We obtain

$$\begin{array}{c}\mathrm{\Delta }\mathrm{\Sigma }(𝐩)=\mathrm{\Sigma }_+(𝐩)\mathrm{\Sigma }_{}(𝐩)=\underset{n=1}{\overset{\mathrm{}}{}}\frac{1}{N^n}\underset{𝐩_1,\mathrm{},𝐩_n}{}\underset{𝒟_n}{}X_{𝒟_n}(V;𝐩,𝐩_1,\mathrm{},𝐩_n)\hfill \\ \hfill \left[G_+(𝐩_1)\mathrm{}G_+(𝐩_n)G_{}(𝐩_1)\mathrm{}G_{}(𝐩_n)\right].\end{array}$$
(61)
The difference of the products of one-electron propagators can further be rewritten to a sum

$$\begin{array}{c}G_+(𝐩_1)\mathrm{}G_+(𝐩_n)G_{}(𝐩_1)\mathrm{}G_{}(𝐩_n)=\underset{i=1}{\overset{n}{}}G_+(𝐩_1)\mathrm{}G_+(𝐩_{i1})\hfill \\ \hfill \times \mathrm{\Delta }G(𝐩_i)G_{}(𝐩_{i+1})\mathrm{}G_{}(𝐩_n)\end{array}$$
(62)
where left (right) to the difference $`\mathrm{\Delta }G(𝐩_i)`$ only $`G_+(𝐩_j)`$ ($`G_{}(𝐩_j)`$) appear. We sum all diagrams for the fixed difference of the one-electron propagators, being now two-particle irreducible diagrams from the electron-hole channel. We then come to a new representation

$$\mathrm{\Delta }\mathrm{\Sigma }(𝐩)=\underset{n=1}{\overset{\mathrm{}}{}}\underset{i=1}{\overset{n}{}}\frac{1}{N}\underset{𝐩_1}{}\mathrm{\Lambda }_{\mathrm{𝐩𝐩}_i}^{(i,ni)}(z_+,z_{};𝐪)\mathrm{\Delta }G(𝐩_i)$$
(63)
where $`\mathrm{\Lambda }_{\mathrm{𝐩𝐩}_i}^{(n,i)}(z_+,z_{};𝐪)`$ is a sum of two-particle irreducible diagrams with $`n`$ internal one-electron lines of which $`i`$ lines carry energy $`z_+`$ and $`ni`$ lines energy $`z_{}`$.

Last but very important step in the proof of identity (25) is an assumption that removing the difference $`\mathrm{\Delta }G`$ from Eq. (63) and fixing the internal variable $`𝐩_i=𝐩^{}`$ does not change summability of the perturbation expansion. If so, we can write

$$\mathrm{\Lambda }_{\mathrm{𝐩𝐩}^{}}(z_+,z_{};𝐪)=\underset{n=1}{\overset{\mathrm{}}{}}\underset{i=1}{\overset{n}{}}\mathrm{\Lambda }_{\mathrm{𝐩𝐩}^{}}^{(i,ni)}(z_+,z_{};𝐪).$$
(64)
The equality holds if and only if the perturbation expansion in the random potential for the two-particle irreducible vertex $`\mathrm{\Lambda }_{\mathrm{𝐩𝐩}^{}}(z_+,z_{};𝐪)`$ converges point-wise for the chosen values of independent variables $`z_+,z_{},𝐩,𝐩^{},𝐪`$. This is a rather severe restriction on applicability of the Ward identity (25). It says that the Ward identity is valid only in cases where no selective rules for sums of two-particle irreducible diagrams apply, i. e., all classes of diagrams are equally important. This happens when the irreducible two-particle vertex is regular (bounded). However, physically the most interesting situation occurs when due to backscatterings a Cooper pole appears in the irreducible vertex $`\mathrm{\Lambda }`$. Thanks to this pole only a class of relevant (crossed electron-hole) diagrams determines the low-energy behavior of the two-particle vertex. In the asymptotic region of the Cooper pole in the electron-hole irreducible vertex, the Ward identity (25) cannot be proved. In fact it gets violated whenever we take into account a selected series of diagrams dictated by causality of the approximation or leading to a divergence and we have to consider a nonperturbative solution for the two-particle irreducible vertex. Each selected (relevant) diagram of order $`n`$ from the formal expansion of the singular two-particle vertex generates in the Ward identity $`n1`$ irrelevant diagrams to complete the sum from the right-hand side of Eq. (62) so that the self-energy difference can be represented via the vertex $`\mathrm{\Lambda }`$. The irrelevant diagrams are neglected in approximations with selected dominant classes of diagrams violating then identity (25). Hence summability of the perturbation expansion for the two-particle irreducible vertex is an additive assumption in the proof of the Ward identity (25). There is no nonperturbative proof of this identity except for the homogeneous case, $`q=0`$, where it is a consequence of the Velický-Ward identity (11).


format: text

text:
# Hadron formation in high energy photonuclear reactions

## I Introduction

In a high energy collision between two hadrons or a photon and a hadron it takes a finite amount of time for the reaction products to evolve to physical particles. During the collision process some momentum transfer between the hadrons or some hard scattering between two of the hadrons’ constituents leads to the excitation of hadronic strings. The time that is needed for the creation and fission of these strings as well as for the hadronization of the string fragments cannot be calculated within perturbative QCD because the hadronization process involves small momentum transfers of typically only a few hundred MeV. One can perform an estimate of the formation time $`\tau _f`$ in the rest frame of the hadron. It should be of the order of the time that the quark-antiquark (quark-diquark) pair needs to reach a separation, that is of the size of the produced hadron ($`r_h0.60.8`$ fm):

$$\tau _f\frac{r_h}{c}.$$
(1)
During their evolution to physical hadrons the reaction products will react with reduced cross sections. This is motivated by means of color transparency: the strings and the substrings created during the fragmentation are in a color singlet state and therefore react with a cross section that increases with their transverse size. As a consequence the produced hadrons travel inside the nuclear medium with a reduced scattering probability during their formation time. Hence the formation time plays an important role in the dynamics of nuclear reactions, e.g. heavy ion collisions, proton and pion induced reactions as well as photon and electron induced reactions on nuclei.

The latter two are of special interest because they are less complex than heavy ion collisions and, in contrast to hadron induced reactions, the primary reaction does in general not only take place at the surface of the nucleus but also at larger densities. Experiments at TJNAF and DESY, for example, deal with exclusive and semi-inclusive meson photo- and electroproduction at high energies. Large photon energies $`E_\gamma `$ are of special interest because of time dilatation the formation length $`l_f`$ in the rest frame of the nucleus can exceed nuclear dimensions:

$$l_f=v_h\gamma \tau _f=\frac{p_h}{m_h}\tau _f.$$
(2)
If one chooses the formation time to be $`\tau _f=0.8`$ fm/c, the formation length in the rest frame of the nucleus will be about 30 fm for a 5 GeV pion and about 7 fm for a 5 GeV kaon or a 7 GeV $`\rho `$ meson. These lengths have to be compared with the typical size of nuclear radii, e.g. 2.7 fm for 12C and 7.1 fm for 208Pb. Because it suppresses the final state interactions, the formation time has a big effect on photonuclear production cross sections at high energies. Turning this argument around, exclusive and semi-inclusive photoproduction of mesons on nuclei  offer a possibility to study these formation times if the FSI are well under control.

Usually, the FSI are modeled within Glauber theory  and treated purely absorptive. A more realistic description of the FSI has to take also regeneration of the mesons studied into account. For this we use a coupled channel semi-classical transport model based on the Boltzmann-Uehling-Uhlenbeck (BUU) equation. Originally developed to describe heavy ion collisions  at SIS energies it has been extended in later works to investigate also inclusive particle production in heavy ion collisions up to 200 AGeV and $`\pi `$  and $`p`$ induced as well as photon and electron induced reactions in the resonance region . Inclusive photoproduction of mesons at energies between 1 and 7 GeV has been investigated in . An attractive feature of this model is its capability to describe a large variety of very different reaction types in a consistent way.

Since photon induced reactions are known to be shadowed ($`\sigma _{\gamma A}<A\sigma _{\gamma N}`$) above photon energies of approximately 1 GeV , one needs a way to account for this shadowing effect in photoproduction. This is straightforward within Glauber theory, but it is not clear how to account for this initial state coherence length effect in a semi-classical transport model for the FSI. A first attempt at combining the quantum mechanical coherence in the entrance channel with our incoherent treatment of the FSI has been made in . In Sec. II we present a new, improved way to implement shadowing in our model. In Sec. III the results for exclusive $`\rho ^0`$ photoproduction on 208Pb at 7 GeV are compared with the predictions from simple Glauber theory. In Sec. III we also demonstrate the importance of the coupled channel treatment of the FSI by looking at the nuclear $`K^+`$ photoproduction cross section. We summarize our results in Sec. IV.

## II Model

In this section we describe how we account for shadowing in our model and only sketch the principles of the transport model itself. For a detailed description of the transport model used we refer to . In our model the incoherent reaction of a high energy photon with a nucleus takes place in two steps. In the first step the photon reacts with one nucleon inside the nucleus (impulse approximation) and produces some final state $`X`$. In this process nuclear effects like Fermi motion, binding energies and Pauli blocking of the final state nucleons are taken into account. In the second step the final state $`X`$ is propagated within the transport model. Except for the elastic vector meson and exclusive strangeness production (see ) we use the Lund string model FRITIOF  to describe high energy photoproduction on the nucleon. The same string model is also used to deal with high energy particle collisions in the FSI. Since FRITIOF does not accept photons as incoming particles, $`\rho `$ dominance was used in  and the photon was passed as a massless $`\rho ^0`$. This led to an excellent description of charged particle multiplicities in $`\gamma p`$ collisions.

Here we generalize to vector meson dominance (VMD)  by writing the incoming photon state as

$$|\gamma =\left(1\underset{V=\rho ,\omega ,\varphi }{}\frac{e^2}{2g_V^2}\right)|\gamma _0+\underset{V=\rho ,\omega ,\varphi }{}\frac{e}{g_V}|V,$$
(3)
and pass the photon as a massless $`\rho ^0`$, $`\omega `$ or $`\varphi `$ with a probability corresponding to the strength of the vector meson coupling to the photon times its nucleonic cross section $`\sigma _{VN}`$

$$P(V)=\frac{\left(\frac{e}{g_V}\right)^2\sigma _{VN}}{\sigma _{\gamma N}}.$$
(4)
These probabilities follow from using (3) in the optical theorem:

$$\sigma _{\gamma N}=\left(1\underset{V=\rho ,\omega ,\varphi }{}\frac{e^2}{2g_V^2}\right)^2\sigma _{\gamma _0N}+\underset{V=\rho ,\omega ,\varphi }{}\left(\frac{e}{g_V}\right)^2\sigma _{VN}.$$
(5)
As can be seen from (3) and (5) there is also a finite probability

$$P(\gamma _0)=1\underset{V=\rho ,\omega ,\varphi }{}P(V)$$
(6)
that the ’bare’ photon $`\gamma _0`$ has to be passed to FRITIOF; the component $`\gamma _0`$ does not get shadowed in the nucleus. In a generalized VMD model for example the $`\gamma _0`$ includes contributions from heavy intermediate hadronic states which, as we show later, are not shadowed because of the large momentum transfer $`q_V`$ that is needed to put them on their mass shell. Since FRITIOF does not accept a ’bare’ photon as input we replace it again by a vector meson $`V`$(=$`\rho ^0`$, $`\omega `$ or $`\varphi `$), with the probability

$$P_{\gamma _0}(V)=\frac{\left(\frac{e}{g_V}\right)^2\sigma _{VN}}{_{V^{}=\rho ,\omega ,\varphi }\left(\frac{e}{g_V^{}}\right)^2\sigma _{V^{}N}}.$$
(7)
The particle production in FRITIOF can be decomposed into two parts. First there is a momentum transfer taking place between the two incoming hadrons leaving two excited strings with the quantum numbers of the initial hadrons. After that the two strings fragment into the observed particles. As a formation time we use 0.8 fm/c in the rest frame of each hadron; during this time the hadrons do not interact with the rest of the system.

Up to now we did not take shadowing into account. In Glauber theory this corresponds to using only the left amplitude in Fig. 1, where a photon directly produces some hadron $`X`$ at nucleon $`j`$. The left amplitude alone leads to the unshadowed incoherent photoproduction cross section

$$\sigma _{\gamma AXA^{}}^{\text{unshadowed}}=\sigma _{\gamma NXN}d^2b𝑑z_j\stackrel{~}{n}(\stackrel{}{b},z_j)e^{\sigma _{XN}_{z_j}^{\mathrm{}}𝑑zn(\stackrel{}{b},z)}$$
(8)
where $`n(\stackrel{}{r})`$ denotes the nucleon number density, $`\stackrel{~}{n}(\stackrel{}{r})`$ the number density of nucleons with the correct charge to produce the hadron $`X`$ and $`\sigma _{XN}`$ is the total $`XN`$ cross section. The exponential damping factor in (8) describes the absorption of the particle $`X`$ on its way out of the nucleus.

In Glauber theory shadowing arises from the interference of the left amplitude in Fig. 1 with the second amplitude of order $`\alpha _{em}`$ which is shown on the right hand side. In this process the photon produces a vector meson $`V`$ on nucleon $`i`$ without excitation of the nucleus. This vector meson then scatters at fixed impact parameter $`\stackrel{}{b}`$ (eikonal approximation) through the nucleus to nucleon $`j`$ and produces the final state meson $`X`$, leaving the nucleus in the same excited state as in the direct process. Off-diagonal scattering, where a vector meson $`V`$ scatters into a different vector meson $`V^{}`$, is usually neglected. In view of the upcoming discussion of exclusive $`\rho ^0`$ photoproduction on Pb we state here the Glauber result  for the incoherent vector meson production cross section:

$`\sigma _{\gamma AVA^{}}`$ $`=`$ $`\sigma _{\gamma NVN}{\displaystyle d^2b_{\mathrm{}}^{\mathrm{}}𝑑z_jn(\stackrel{}{b},z_j)e^{\sigma _{VN}_{z_j}^{\mathrm{}}𝑑zn(\stackrel{}{b},z)}}`$ (10)
$`\times \left|1{\displaystyle _{\mathrm{}}^{z_j}}𝑑z_in(\stackrel{}{b},z_i){\displaystyle \frac{\sigma _{VN}}{2}}(1i\alpha _V)e^{iq_Vz_i}\mathrm{exp}\left[{\displaystyle \frac{1}{2}}\sigma _{VN}(1i\alpha _V){\displaystyle _{z_i}^{z_j}}𝑑zn(\stackrel{}{b},z)\right]\right|^2.`$
Again $`n(\stackrel{}{r})`$ denotes the nucleon number density, $`\sigma _{VN}`$ the total $`VN`$ cross section and $`\alpha _V`$ the ratio of real to imaginary part of the $`VN`$ forward scattering amplitude. The momentum transfer

$$q_V\frac{m_V^2}{2E_\gamma }$$
(11)
arises from putting the vector meson on its mass shell. A large momentum transfer $`q_V`$ is suppressed by the elastic nuclear formfactor and leads to less shadowing due to the oscillating factor $`\mathrm{exp}[iq_Vz_i]`$ in the integrand of Eq. (10). Note that $`q_V`$ is just the inverse of the coherence length, i.e. the distance, given by the uncertainty principle, which the photon can travel as a $`V`$ fluctuation. For the quantities in (10) we use the parameterizations of Model I of Ref.  with which one obtains a very good description of the shadowing effect in nuclear photoabsorption down to the onset region .

One clearly sees from (10) how the FSI separate from the ’initial state interactions’ of the photon. We now replace the exponential damping factor $`\mathrm{exp}[\sigma _{VN}_{z_j}^{\mathrm{}}𝑑zn(\stackrel{}{b},z)]`$, which corresponds to purely absorptive FSI, by a transport model. This allows us to incorporate a wider class of FSI. In addition we want to include events where the final vector meson $`V`$ is not produced in the primary reaction but via sidefeeding. This leads to less reduction of the nuclear production cross section even with the same $`\sigma _{VN}`$. We thus also need to know how all the other possible primary reactions are shadowed. We, therefore, start from (3) and use Glauber theory to calculate how the single $`V`$ components of the photon change due to multiple scattering on the way to nucleon $`j`$ where the state $`X`$ is produced :

$$|\gamma (\stackrel{}{r}_j)=\left(1\underset{V=\rho ,\omega ,\varphi }{}\frac{e^2}{2g_V^2}\right)|\gamma _0+\underset{V=\rho ,\omega ,\varphi }{}\frac{e}{g_V}\left(1\overline{\mathrm{\Gamma }_V}(\stackrel{}{r}_j)\right)|V.$$
(12)
Here the (photon energy dependent) nuclear profile functions for the different vector meson components $`V`$ are given as

$$\overline{\mathrm{\Gamma }_V}(\stackrel{}{b},z_j)=_{\mathrm{}}^{z_j}𝑑z_in(\stackrel{}{b},z_i)\frac{\sigma _{VN}}{2}(1i\alpha _V)e^{iq_V(z_iz_j)}\mathrm{exp}\left[\frac{1}{2}\sigma _{VN}(1i\alpha _V)_{z_i}^{z_j}𝑑zn(\stackrel{}{b},z)\right].$$
(13)
Note that the $`\gamma _0`$ component is, by definition, not modified due to the presence of the nuclear medium. The cross section for the photon to react with nucleon $`j`$ at position $`\stackrel{}{r}_j`$ inside the nucleus can be deduced via (12) from the optical theorem:

$$\sigma _{\gamma N}(\stackrel{}{r}_j)=\left(1\underset{V=\rho ,\omega ,\varphi }{}\frac{e^2}{2g_V^2}\right)^2\sigma _{\gamma _0N}+\underset{V=\rho ,\omega ,\varphi }{}\left(\frac{e}{g_V}\right)^2\left|1\overline{\mathrm{\Gamma }_V}(\stackrel{}{r}_j)\right|^2\sigma _{VN}.$$
(14)
Like for the photon in vacuum each term gives the relative weight for the corresponding photon component to be passed to FRITIOF. When integrated over the whole nucleus one gets from Eq. (14) the total incoherent photonuclear cross section

$$\sigma _{\gamma A}^{inc}=d^3r_jn(\stackrel{}{r}_j)\sigma _{\gamma N}(\stackrel{}{r}_j)$$
(15)
which is shown in Fig. 2 by the dashed line together with the total nuclear photoabsorption cross section as derived from the optical theorem in :

$`\sigma _{\gamma A}^{tot}`$ $`=`$ $`A\sigma _{\gamma N}{\displaystyle \underset{V=\rho ,\omega ,\varphi }{}}{\displaystyle \frac{k_V}{2k}}\left({\displaystyle \frac{e}{g_V}}\right)^2\text{Re}\{\sigma _{VN}^2(1\alpha _V)^2{\displaystyle }d^2b{\displaystyle _{\mathrm{}}^{\mathrm{}}}dz_i{\displaystyle _{z_i}^{\mathrm{}}}dz_jn(\stackrel{}{b},z_i)n(\stackrel{}{b},z_j)`$ (17)
$`\times e^{iq_V(z_iz_j)}\mathrm{exp}[{\displaystyle \frac{1}{2}}\sigma _{VN}(1i\alpha _V){\displaystyle _{z_i}^{z_j}}dz^{}n(\stackrel{}{b},z^{})]\}.`$
Here $`k`$ and $`k_V`$ denote the momentum of the photon and the vector meson respectively and we again use the parameterizations of Model I of Ref. . More than 90% of the difference between those two cross sections stems from coherent $`\rho ^0`$ photoproduction.

In Fig. 3 we show how strongly the $`\rho ^0`$ and the $`\varphi `$ components of a real 20 GeV photon are separately shadowed in Pb. We plot the number density of the nucleons reacting with the $`V`$ components of the photon

$$a_{eff}^V(\stackrel{}{r}_j)=n(\stackrel{}{r}_j)\frac{1}{\sigma _{\gamma N}}\left(\frac{e}{g_V}\right)^2\left|1\overline{\mathrm{\Gamma }_V}(\stackrel{}{r}_j)\right|^2\sigma _{VN}$$
(18)
as a function of $`\stackrel{}{r}_j`$. One clearly sees that due to its smaller nucleonic cross section the $`\varphi `$ component is less shadowed than the $`\rho ^0`$ component at the backside of the nucleus. As a consequence strangeness production (e.g. $`K`$ photoproduction), where the primary reaction is preferably triggered by the $`\varphi `$ component of the photon, is less shadowed than, e. g., $`\pi `$ photoproduction. This dependence of the strength of shadowing on the reaction type is new compared to the treatment of shadowing in  and can also be seen directly from the second amplitude in Fig. 1 because of the occurrence of the scattering process $`VNXN`$ at nucleon $`j`$.

As already mentioned above the purely absorptive FSI of the Glauber model are very different from the coupled channel description of a transport model. The transport model we use is based on the BUU equation that describes the time evolution of the phase space density $`f_i(\stackrel{}{r},\stackrel{}{p},t)`$ of particles of type $`i`$ that can interact via binary reactions . In our case these particles are the nucleons of the target nucleus as well as the baryonic resonances and mesons ($`\pi `$, $`\eta `$, $`\rho `$, $`K`$, …) that can either be produced in the primary $`\gamma N`$ reaction or during the FSI. For particles of type $`i`$ the BUU equation looks as follows:

$$\left(\frac{}{t}+\frac{H}{\stackrel{}{r}}\frac{}{\stackrel{}{r}}\frac{H}{\stackrel{}{r}}\frac{}{\stackrel{}{p}}\right)f_i(\stackrel{}{r},\stackrel{}{p},t)=I_{coll}[f_1,\mathrm{}f_i,\mathrm{},f_M].$$
(19)
In the case of baryons the Hamilton function $`H`$ includes a mean field potential which in our model depends on the particle position and momentum. The collision integral on the right hand side accounts for the creation and annihilation of particles of type $`i`$ in a collision as well as elastic scattering from one position in phase space into another. For fermions Pauli blocking is taken into account in $`I_{coll}`$ via blocking factors. For each particle type $`i`$ such a BUU equation exists; all are coupled via the mean field and the collision integral. This leads to a system of coupled differential-integral equations which we solve via a test particle ansatz for the phase space density (for details see ). Since the collision integral also accounts for particle creation in a collision the observed outgoing particle $`X`$ cannot only be produced in the primary reaction, but can also be created by sidefeeding in which a particle $`Y`$ is created first which propagates and then, by FSI, produces $`X`$. In addition the state $`X`$ might get absorbed on its way out of the nucleus but be fed in again in a later interaction. Both cases can a priori not be ignored, but are usually neglected in Glauber models.

## III Results

Exclusive vector meson photo- and electroproduction on nuclei is an ideal tool to study effects of the coherence length, formation time and color transparency. Exclusive $`\rho ^0`$ electroproduction has been investigated in the HERMES experiment  at photon energies between 10 GeV and 20 GeV and virtuality $`Q^2`$ 5 GeV2. The calculations for meson production on nuclei are usually done within Glauber theory . As already mentioned above the FSI in Glauber theory are usually purely absorptive. This means that for the reaction $`\gamma A\rho ^0A^{}`$ the primary reaction has to be $`\gamma N\rho ^0N`$. If one treats the FSI via an absorptive optical potential one gets an exponential damping $`\mathrm{exp}[\sigma _{\rho N}_{z_j}^{\mathrm{}}𝑑zn(\stackrel{}{b},z)]`$ of the nuclear production cross section.

We are presently working at incorporating also incoming virtual photons into the formalism developed in Sec. II and at enlarging the configuration space for the FSI. Here we discuss, therefore, only results obtained with real photons and at a lower energy. In Fig. 4 we show the results of our model for the mass differential cross section of incoherent, exclusive $`\rho ^0`$ photoproduction on 208Pb for $`E_\gamma =7`$ GeV. In this case ’exclusive’ means that the final state consists of a $`\pi ^+\pi ^{}`$ pair and 208 bound nucleons. The solid line represents a calculation with the primary reaction restricted to $`\gamma N\rho ^0N`$. It already includes the effects of shadowing, Fermi motion, Pauli blocking and the nucleon potential, but no FSI. The dotted line shows the effect of the FSI without a formation time of the $`\rho ^0`$ in $`\gamma N\rho ^0N`$. The Glauber model yields quantitatively the same effect of the FSI. This means that in this case (only $`\gamma N\rho ^0N`$ as primary reaction) FSI processes like $`\rho ^0N\pi N`$, $`\pi N\rho ^0N`$, where the primary $`\rho ^0`$ gets absorbed first and is fed into the outgoing channel by a later FSI, are negligible. If one assumes a formation time of $`\tau _f=0.8`$ fm/c for the $`\rho ^0`$, one gets the result indicated by the dash-dotted line. Due to the finite formation time there is considerably less absorption and the nuclear production cross section increases. If the observed spectrum looked like this, one would in Glauber theory be lead to the conclusion of a finite $`\rho ^0`$ formation time .

However, one will get a similar result with $`\tau _f=0`$ if one allows for other primary reaction besides $`\gamma N\rho ^0N`$ and uses a coupled channel model. This can be seen by looking at the dashed line in Fig. 4. We find that about 60% of the additional $`\rho ^0`$ stem from inelastic $`\rho ^0`$ production in the primary reaction, e.g. $`\gamma N\rho ^0\pi N`$ where the $`\pi `$ gets absorbed during the FSI. We now apply an exclusivity measure like the one used in the HERMES experiment

$$2\text{ GeV}<\mathrm{\Delta }E=\frac{p_Y^2m_N^2}{2m_N}<0.6\text{ GeV},$$
(20)
where $`m_N`$ is the nucleon mass and

$$p_Y=p_N+p_\gamma p_\rho $$
(21)
is the 4-momentum of the undetected final state. Here $`p_\gamma `$ and $`p_\rho `$ denote the 4-momenta of the incoming photon and the detected $`\pi ^+\pi ^{}`$ pair and $`p_N`$ is the 4-momentum of the struck nucleon which, for the calculation of $`p_Y`$, is assumed to be at rest. Using (20) leads to a decrease of the cross section (dash-dot-dotted line) because some of the inelastic primary events are excluded. If the exclusivity measure was good enough to single out only the elastic $`\gamma N\rho ^0N`$ reactions from the primary events, the dash-dot-dotted curve would coincide with the dotted line and Glauber theory would be applicable. Since this is not the case, one still extracts a too large formation time when using Glauber theory.

One therefore needs a further constraint in addition to the exclusivity measure (20) which becomes apparent by looking at the differential cross section $`\frac{d\sigma }{dt}`$ in Fig. 5; the meaning of the lines is as before. For $`|t|>0.1`$ GeV2 the full calculation with exclusivity measure (dash-dot-dotted line) gives the same result as the one with the primary reaction $`\gamma N\rho ^0N`$ and FSI (dotted line). In this kinematic regime Glauber theory can therefore be used. For $`|t|<0.1`$ GeV2, however, the exclusivity measure (20) cannot distinguish between elastic $`\rho ^0`$ photoproduction ($`\gamma N\rho ^0N`$) and other primary reactions, e.g. inelastic $`\rho ^0`$ photoproduction ($`\gamma N\rho ^0X`$, $`XN`$). At low values of $`|t|`$ there exist many states $`X`$ with invariant masses which are not excluded by the exclusivity measure (20). In addition we find that about 25% of the finally accepted $`\rho ^0`$ in this $`t`$ region are not produced in the primary reaction but stem from sidefeeding in the FSI. In the HERMES experiment one makes a lower $`|t|`$ cut to get rid of the coherent $`\rho ^0`$ photoproduction contribution. In the case of lead and $`E_\gamma =7`$ GeV the coherent part can be neglected above $`|t|0.05`$ GeV2. For Glauber theory to be reliable one has to increase this threshold to approximately $`|t|=0.1`$ GeV2. Glauber theory can thus be trusted only under certain kinematic constraints.

So far we have discussed a case in which a strongly interacting particle, the $`\rho ^0`$ meson, is produced. In the following we will now discuss the special effects that appear in a coupled channel treatment of the FSI when a weakly interacting particle, such as the $`K^+`$ meson, is considered. In Fig. 6 we therefore show the cross section for the reaction $`\gamma AK^+X`$ in the photon energy range 1-7 GeV for 208Pb. This reaction had already been investigated in ; the new treatment of shadowing and the initial interactions as outlined in Eqs. (3)-(7) lead to an increase of the total yield at 7 GeV by about 20%. The solid curve in Fig. 6 represents the results of the full calculation (including shadowing, FSI, $`\tau _f=0.8`$ fm/c, etc.). By comparison with the calculation without shadowing (dash-dotted line) one sees how important shadowing becomes with increasing photon energy. At 7 GeV it reduces the nuclear production cross section to about 65%.

The importance of a coupled channel treatment of the FSI becomes clear when comparing the full calculation with the one without FSI (dashed line). Since the $`\overline{s}`$ quark cannot be absorbed in medium the FSI can just increase the $`K^+`$ yield via processes like $`\pi NK^+Y`$ ($`Y=\mathrm{\Sigma },\mathrm{\Lambda }`$) for example. With decreasing formation time the primarily produced pions have a greater chance to produce $`K^+`$ in the FSI. As a consequence of this, a shorter formation time will lead to an increase of the nuclear $`K^+`$ photoproduction cross section as can be seen from the dotted line. An enhancement of the $`K^+`$ production cross section due to FSI cannot be explained by purely absorptive FSI as in simple Glauber theory where it would necessarily be interpreted as being due to a longer formation time of the $`K^+`$.

## IV Summary

High energy photoproduction off nuclei offers a promising possibility to study the physics of hadron formation. Necessary for such studies is a reliable model of the FSI to extract the formation time from the production cross sections. Whereas Glauber models allow for a straightforward implementation of the nuclear shadowing effect they usually have the disadvantage of a purely absorptive treatment of the FSI. As we have shown the latter can lead to a wrong estimate of the formation time. A more realistic treatment of the FSI is possible within a coupled channel transport model. We have presented a method to combine such an incoherent treatment of FSI with coherence length effects in the entrance channel which can easily be extended to higher energies and virtual photons. We have shown that in particular the production of mesons with long mean free path will be affected by the coupled channel effects in the FSI.

## V Acknowledgments

This work was supported by DFG.