text:
# Stochastic Aggregation: Scaling Properties

## I Introduction

In the previous study, we introduced a stochastic aggregation process involving both active and passive clusters. We generalized Smoluchowski’s rate equations and obtained exact results for several kernels. In this study, we apply stochastic aggregation to reaction-diffusion, coarsening, and ballistic agglomeration problems. Our goal is to examine the range of validity the mean-field results, and to determine whether the overall scaling behavior extends to low dimensional systems.

The rate equations approach is mean-field in nature, i.e., it is valid only when spatial correlations are absent. Formally, it is applicable only in infinite spatial dimension, or in the presence of an effective mixing mechanism. This mean-field theory should also be asymptotically exact when the spatial dimension is sufficiently high. In low spatial dimensions, however, significant spatial correlations eventually develop, and the rate equation approach does not apply in the long time limit. We therefore focus on one-dimensional systems where spatial correlations are most pronounced.

We performed numerical simulations of stochastic aggregation processes with both diffusive and ballistic particle transport. The simulations show that the scaling behavior suggested by the mean-field theory is indeed generic, as it extends to one-dimensional systems. We find that two nontrivial model-dependent exponents characterize the survival probabilities of the particles and monomers, respectively. Smoluchowski’s theory provides reasonable estimates for these exponents.

Additionally, we studied the limiting mass distribution of passive clusters. Surprisingly, over a substantial mass range, this distribution depends only weakly on the underlying transport mechanism. Furthermore, mean-field theory provides an excellent approximation for the limiting mass distribution.

The rest of this paper is organized as follows. The general scaling behavior is outlined in Sec. II. Predictions of the mean-field theory are summarized in Sec. III. Numerical simulations of stochastic aggregation processes with diffusive and ballistic transport mechanisms are described in Secs. IV and V, respectively. A discussion of the results is presented in Sec. VI.

## II Scaling Properties

Stochastic aggregation involves two types of clusters: active and passive. Initially, the system consists of active monomers only. When two active clusters merge, the newly-born aggregate remains active with probability $`p`$, or becomes passive (i.e., it never aggregates again) with probability $`q=1p`$. Eventually, all active clusters are depleted and the system consists of passive clusters only. This process can be viewed as an aggregation-annihilation process since it interpolates between aggregation ($`p=1`$) and annihilation ($`p=0`$).

Quantities of interest include $`A_k(t)`$ and $`P_k(t)`$, the distributions of active and passive clusters at time $`t`$, as well as the final distribution of passive clusters, $`P_k(\mathrm{})`$. As shown in , two conservation laws underly this system. The first is mass conservation, $`k[A_k(t)+P_k(t)]=\mathrm{const}`$. The second conservation law reflects the fact that changes in the overall densities are coupled, $`qA(t)+(1+q)P(t)=\mathrm{const}`$, where $`A(t)=A_k(t)`$ and $`P(t)=P_k(t)`$ are the number densities of active and passive clusters, respectively.

Therefore, it is sufficient to study the time evolution of the number density and the mass density of active clusters, $`A(t)`$ and $`M(t)=kA_k(t)`$, respectively. The latter quantity is the survival probability of an active particle, i.e., the probability that it still belongs to an active cluster at time $`t`$. Smoluchowski’s theory suggests that both quantities decay algebraically in the long time limit

$$A(t)t^\nu ,M(t)t^{\nu \psi }.$$
(1)
As will be shown below, this as well as other scaling properties suggested by this theory hold qualitatively even for low dimensional stochastic aggregation processes. While the decay exponent $`\nu `$ is typically robust in that it depends only on the major characteristics of the process such as the spatial dimension or the transport mechanism, the exponent $`\psi \psi (p)`$ is non-universal as it depends on the details of the model, i.e., on the probability $`p`$. In turn, this implies a non-universal growth law for the average mass of an active cluster $`k=M/At^{\nu (1\psi )}`$.

For the system to follow a scaling behavior, the average mass must be the only relevant scale in the long time limit, and conversely, any scale characterizing the initial mass distribution must be “erased” eventually. In other words, the mass distribution is characterized by a single rescaled variable

$$A_k(t)t^{\nu (\psi 2)}F\left(kt^{\nu (\psi 1)}\right),$$
(2)
with the time dependent prefactor fixed by the decay laws (1).

This scaling behavior is similar to that found for deterministic aggregation-annihilation processes and for aggregation-annihilation of domains in coarsening processes. These studies suggest that another independent exponent describes the decay of small clusters. Specifically, the monomer density decays according to

$$A_1(t)t^{\nu \delta },$$
(3)
with a model-dependent exponent $`\delta \delta (p)`$. The monomer density decay reflects the small argument behavior of the scaling function $`F(\xi )\xi ^\sigma `$ with $`\delta 1=(1\psi )(1+\sigma )`$. One of our main results is that the mass distribution of active clusters is described by a set of nontrivial exponents $`(\psi ,\delta )`$. These exponents can be viewed as persistence exponents as they characterize the survival probability of an active particle, and an active monomer .

Several properties of the scaling exponents are general. For instance, the inequalities $`\psi 1\delta `$ hold since $`A_1A_kkA_k`$. The two exponents are equal $`\psi =\delta =1`$ in the annihilation case ($`p=0`$), since $`A_k(t)=A(t)\delta _{k,1}`$. In the aggregation limit ($`p=1`$) the mass density of active clusters is conserved and therefore $`\psi =0`$.

We now turn to the mass distribution of passive clusters. The Smoluchowski theory suggests that the same scaling form underlies both mass distributions

$$P_k(t)t^{\nu (\psi 2)}G\left(kt^{\nu (\psi 1)}\right).$$
(4)
In contrast with the active cluster distribution, the passive cluster distribution approaches a nontrivial final distribution $`P_k(\mathrm{})`$. Such a time independent final distribution is consistent with the above scaling form only when the scaling function diverges, $`F(\xi )\xi ^\gamma `$ in the limit $`\xi 0`$, with $`\gamma =(2\psi )/(1\psi )`$. As a result, the final mass distribution of passive clusters decays algebraically in the large mass limit

$$P_k(\mathrm{})k^\gamma \mathrm{with}\gamma =\frac{2\psi }{1\psi }.$$
(5)
At a given time $`t`$, this decay is realized for clusters whose mass $`k`$ does not exceed the characteristic mass $`k^{}t^{\nu (1\psi )}`$. Note also that $`0<\psi <1`$ implies $`2<\gamma <\mathrm{}`$. Generally, the mass conservation restricts the large mass decay exponent to $`\gamma >2`$. Since the $`\psi `$ exponent varies between $`0`$ and $`1`$, we see that the entire range of acceptable exponents is realized by tuning the probability $`p`$.

## III Mean-Field Theory

It is well established that spatial correlations can be safely neglected only in spatial dimensions larger than some upper critical dimension, $`d>d_c`$. For example, for irreversible aggregation with mass-independent diffusion and reaction rates, one has $`d_c=2`$; for a general aggregation process, however, the upper critical dimension may be arbitrarily large. Below the upper critical dimension, substantial spatial correlations develop, and the most important features including the scaling exponents and the scaling functions are changed. Generally, the lower the spatial dimension, the larger the difference with the mean-field predictions.

Although the Smoluchowski rate equations approach does not apply in low spatial dimensions, it can still serve as a useful approximation after an appropriate modification. This can be accomplished by replacing the overall reaction rate with an effective density-dependent reaction rate $`rr(A)`$

$`{\displaystyle \frac{dA_k}{dt}}`$ $`=`$ $`r\left({\displaystyle \frac{p}{2}}{\displaystyle \underset{i+j=k}{}}A_iA_jA_kA\right),`$ (6)
$`{\displaystyle \frac{dP_k}{dt}}`$ $`=`$ $`r\left({\displaystyle \frac{q}{2}}{\displaystyle \underset{i+j=k}{}}A_iA_j\right).`$ (7)
We are primarily interested in situations where aggregation is independent of the mass, and therefore we use a mass independent rate kernel. The reaction rate $`r(A)`$ is model dependent. In reaction-diffusion processes, the reaction rate decays algebraically with the density (see, e.g., Refs.). Assuming $`r(A)A^\alpha `$ yields $`\frac{dA}{dt}A^{\alpha +2}`$, and consequently, the density decay exponent is found

$`\nu ={\displaystyle \frac{1}{1+\alpha }}.`$ (8)
In general, a reduction in the reaction rate, i.e., $`\alpha >0`$, leads to a slowing down in the density decay rate, $`\nu <1`$. Apart from the change in $`\nu `$, all other aspects of this approximation are identical to the Smoluchowski theory with a constant rate kernel. Indeed, the above rate equations reduce to the Smoluchowski’s rate equations with a redefined time variable, $`t\tau =_0^t𝑑t^{}r(t^{})`$. In particular, the scaling exponents $`\psi `$ and $`\delta `$ are independent of $`\alpha `$:

$$\psi =2\frac{1p}{2p},\delta =\frac{2}{2p}.$$
(9)
One can verify the expected limiting behaviors $`\psi (1)=0`$, and $`\psi (0)=\delta (0)=1`$. Furthermore, the scaling functions are as in the constant kernel solution , and for example, $`F(\xi )`$ is purely exponential. The corresponding small argument exponents $`\gamma =2/p`$ and $`\sigma =0`$ follow from $`\psi `$ and $`\delta `$ using the aforementioned scaling relations. The final mass distribution of passive clusters is independent of the reaction rate $`r`$

$$P_k(\mathrm{})=\frac{q}{p}\frac{\mathrm{\Gamma }(1+2/p)\mathrm{\Gamma }(k)}{\mathrm{\Gamma }(k+2/p)}.$$
(10)
Below, we compare these mean-field predictions with simulation results for one-dimensional stochastic aggregation where spatial correlations are most pronounced. We also examine the role of the aggregates’ transport mechanism by considering both diffusive and ballistic transport.

## IV Diffusive Transport

In diffusive stochastic aggregation, identical particles are placed onto a $`d`$-dimensional lattice. All particles perform independent random walks, i.e., they hop to a randomly chosen nearest-neighbor site with a constant rate. If this site is occupied, the two particles coalesce irreversibly, and with probability $`p`$ the resulting aggregate remains active, while with probability $`q=1p`$ it becomes passive. Effectively, passive particles are removed from the system.

In the case of single-species reaction diffusion processes, the effective reaction rate can be obtained from dimensional analysis. Eq. (6) implies $`[r]=[L]^d[T]^1`$, and since the reaction rate can only be a function of the diffusion coefficient $`[D]=[L]^2[T]^1`$ and the density $`[A]=[L]^d`$, one finds $`rDA^{(2d)/d}`$. Hence, $`\alpha =(2d)/d`$ and Eq. (8) yields the correct decay exponents $`\nu =d/2`$ below the upper critical dimension $`d_c=2`$.

To examine the above scaling picture we performed numerical simulations of diffusive stochastic aggregation processes in one dimension. Unless noted otherwise, the data was obtained from an average over 10 independent realizations in a system of size $`L=10^7`$ with periodic boundary conditions. Initially, all sites were occupied. First, we verified that the number density, the mass density, and the monomer density indeed decay algebraically in the long time limit, in accord with Eqs. (1) and (3). The case $`p=1/2`$ is shown in Fig. 1, and the corresponding decay exponents were found: $`\nu =0.500(1)`$, $`\psi =0.6193(3)`$, and $`\delta =1.460(2)`$. Mean-field theory correctly predicts $`\nu =1/2`$. Furthermore, the predictions $`\psi =2/3`$ and $`\delta =4/3`$ provide a reasonable approximation. One can compare with the case of disordered (Sinai) diffusion where a real-space decimation procedure was used to determine exact values of these exponents. Remarkably, the disorder case exponent $`\psi =0.61937`$ is in excellent agreement with the simulation value. There is a small discrepancy with the second exponent $`\delta =1.47041`$. Additionally, we verified that the densities of active and passive clusters follow the scaling forms of Eqs. (2) and (4), respectively (see Fig. 2). In agreement with the mean-field theory, the scaling functions decay exponentially for large masses.

We also studied how the exponents vary with the probability $`p`$, as shown in Figs. 3 and 4. The exact exponents found for the disordered case by Le Doussal and Monthus provide an excellent approximation (within 0.1%) for $`\psi `$. In the case of $`\delta `$, a different behavior emerges in the aggregation limit, $`p1`$, where the exact value is $`\delta =3`$ , and the disagreement with both mean-field theory and the disordered case are most pronounced.

The above scaling arguments suggest that the limiting mass distribution of passive clusters decays algebraically with the exponent $`\gamma =(2\psi )/(1\psi )`$. For $`p=1/2`$, one therefore expects $`\gamma 3.627`$ (compare with $`\gamma =3.62722`$ and $`\gamma =4`$, predicted by the disordered case and the mean-field theory). This corresponds to a very strong suppression of large masses, and therefore, it is much more difficult to confirm this behavior numerically. Nevertheless, our simulations (Fig. 5) are consistent with a power law decay with an exponent $`\gamma 3.6`$.

In one dimension, the diffusion-controlled stochastic aggregation is equivalent to the Potts model with zero-temperature Glauber dynamics. For the $`Q`$-state Potts model with spatially uncorrelated initial conditions, aggregation of domain walls occurs with probability $`p=\frac{Q2}{Q1}`$, and annihilation occurs with probability $`q=\frac{1}{Q1}`$. Therefore, the above can be reformulated in terms of domain walls rather than aggregates. In the coarsening context, the domain wall mass (or number) distribution is dual to domain number distribution .

## V Ballistic Transport

The situation when particles move ballistically involves several complications. First, while the annihilation limit is uniquely defined , the aggregation limit has various realizations. In traffic flows, the velocity of a newly-born cluster is the smaller of the two velocities , while in application to astrophysics and granular gases the velocity follows from momentum conservation . Second, numerical results for the annihilation case  and analytical results for the traffic case  show that the initial conditions are remembered forever, in contrast with the diffusive case. Specifically, the small velocity characteristics of the initial velocity distribution influence the long time asymptotic behavior, including the scaling exponents.

We consider the momentum conserving case, also known as “ballistic aggregation” or “sticky gas” . The initial velocities are assigned according to the distribution $`P_0(v)`$. The mass (momentum) of a newly-born cluster is equal to the sum of masses (momenta) of the two colliding clusters. After an agglomeration event, the newborn particle remains active with probability $`p`$, or becomes passive with probability $`q=1p`$.

To apply the Smoluchowski rate equations approach, we again use dimensional analysis to calculate the decay exponent $`\nu `$. The collision rate is $`rva^{d1}`$, where $`v`$ is the typical velocity and $`a`$ is the typical radius of an aggregate. A particle of radius $`a`$ contains of the order $`a^d`$ monomers whose initial momenta are uncorrelated. Momentum conservation therefore implies $`va^{d/2}`$. Using $`a^dM/AA^{\psi 1}`$ gives the collision rate $`ra^{(d2)/2}A^{(d2)(\psi 1)/2d}`$. From Eq. (8) one finds

$$\nu =\frac{2d}{d+2+\psi (d2)},$$
(11)
with $`\psi `$ given by Eq. (9). Apart from the exponent $`\nu `$, features such as the exponential mass distribution and the exponents $`\psi `$ and $`\delta `$ are given by the mean-field theory outlined above. In two dimensions, the collision rate does not depend on $`\psi `$ and hence, the asymptotic behavior $`At^1`$ agrees with that found for deterministic ballistic agglomeration. For $`d2`$, stochastic and deterministic asymptotics differ: stochasticity enhances decay of the number density $`A`$ for $`d<2`$ and weakens it for $`d>2`$. A more detailed mean-field theory can be carried. It yields a factorizing joint mass-velocity distribution, with an exponential mass distribution, and a Gaussian velocity distribution .

In the aggregation case, $`\psi =0`$ and therefore the correct scaling exponent $`\nu =2d/(d+2)`$  is recovered from Eq. (11). For the annihilation case, however, initial conditions are “remembered” forever and therefore the above dimensional arguments no longer hold. The predicted exponent in the annihilation case is always mean-field $`\nu =1`$, while one-dimensional numerical simulations yield an exponent continuously varying from $`0`$ to $`1`$ depending on the initial velocity distribution $`P_0(v)`$, e.g., $`\nu 0.8`$ for uniform initial distributions.

We have simulated the stochastic aggregation process on a one-dimensional ring with $`10^6`$ particles. The initial velocity distribution was uniform in $`[1,1]`$. We measured the scaling exponent $`\psi `$ via the scaling relation $`MA^\psi `$, rather than directly versus time, since the exponent $`\nu (p)`$ is not known analytically. We have found that the mean-field prediction, $`\psi =(22p)/(2p)`$, provides a reasonable approximation for the exponent $`\psi `$, as shown in Fig. 6. Furthermore, this approximation should improve in higher dimensions.

We compared the mean-field prediction for the mass distribution of passive clusters, Eq. (10), with the numerically obtained distributions in both ballistic and diffusive cases. Interestingly, the rate equations provide an excellent approximation for small and moderate masses (see Fig. 7). Given that the discrepancy in $`\psi `$ is maximal for the case $`p=1/2`$, one may expect an even better approximation for other values of $`p`$. Noting the strong decay of this distribution, the contribution of very large masses is extremely small; for example, $`P_{100}(\mathrm{})2.4\times 10^7`$ for $`p=1/2`$. Hence, the most pronounced part of the distribution is well approximated by the rate equations theory. Surprisingly, the transport mechanism does not play an important role as far as the final mass distribution of passive clusters is concerned.

## VI Discussion

We have investigated diffusion- and ballistic-controlled stochastic aggregation in one dimension. We have seen that the rate equations approach captures the overall scaling behavior and additionally it provides reasonable estimates for the decay exponents. In general, the mass distribution is characterized by two nontrivial model-dependent decay exponents.

In the diffusion-controlled case, the exponent $`\psi `$ underlying the survival probability of a particle is in excellent agreement with the exact results from the disordered case. In fact, one cannot dismiss the possibility that the disordered and the pure values are identical, based on numerics alone. However, there is an evident discrepancy in the exponent $`\delta `$ as the disordered case exponent diverges logarithmically in the aggregation limit. Stochastic aggregation is equivalent to domain coarsening in the zero-temperature Potts-Glauber model. The above exponents $`(\psi ,\delta )`$ characterize the domain wall number distribution in analogy with $`(\psi _D,\delta _D)`$ for the domain number distribution . In the latter case as well, exact values calculated for the disordered case provide an excellent approximation for the domain exponents. In general, the particle survival probability exponent $`\psi `$ is robust, while the monomer survival probability exponent $`\delta `$ is very sensitive to the details of the process.

In the ballistic-controlled case, we have shown that even in the absence of a consistent mean-field theory, some characteristics such as the exponent $`\psi `$ are well approximated by the rate equations. Understanding of reaction processes with an underlying ballistic transport remains largely incomplete. The asymptotic behavior is highly sensitive to the initial conditions, and the critical dimension is apparently infinite. In fact, exact analytical results are available mostly in the aggregation limit .

The most intriguing property of the stochastic aggregation is the profound lack of universality. Indeed, the weak dependence on the transport mechanism is in contrast with the strong dependence on the parameter $`p`$. For example, our numerical results show that the final distribution of passive clusters is very close in diffusion- and ballistic-controlled situations. Another very impressive manifestation of this is the excellent agreement between the values of the exponent $`\psi (p)`$ in the disordered and pure cases.

This research was supported by the DOE (W-7405-ENG-36), NSF (DMR9632059), and ARO (DAAH04-96-1-0114).


format: text

text:
# From localization to delocalization in the quantum Coulomb glass

## I Introduction

The discovery of a metal-insulator transition (MIT) in the two-dimensional electron gas in Si-MOSFETs  has induced renewed attention to the transport properties of disordered electrons. This MIT is in conflict with the theory of localization for non-interacting electrons which predicts that all states are localized in 2D. The electron density in the Si-MOSFETs is very low which makes the electron-electron interaction particularly important. Thus it is generally assumed that some type of interaction effect is responsible for this MIT. One of the most remarkable findings about the MIT in Si-MOSFETs is that an in-plane magnetic field (which does not couple to the orbital motion of the electrons) strongly suppresses the conducting phase . This suggests that the spin degrees of freedom play an important role for the transition. A complete understanding has, however, not yet been obtained. There have been attempts to explain the experiments based on the perturbative renormalization group , non-perturbative effects , or the transition being a superconductor-insulator transition rather than a MIT .

In order to attack the problem of disordered interacting electrons numerically we have developed  an efficient method, the Hartree-Fock based diagonalization (HFD) which is related to the quantum-chemical configuration interaction approach. We have used this method to study the influence of interactions on the conductance in one , two , and three  dimensions. We found a delocalizing tendency of the interactions for strong disorder but a localizing one for weak disorder. Similar results have been obtained by means of the density-matrix renormalization group  in one dimension and exact diagonalization in two dimensions . Since in most of the numerical work in the literature spinless electrons were considered, there are not many results about the importance of the spin degrees of freedom.

In this work we address this question by generalizing the HFD method to spinful electrons. We then use it study the influence of the spin degrees of freedom on the Kubo-Greenwood conductance.

## II Model and method

The generic model for spinless interacting disordered electrons is the quantum Coulomb glass . In this paper we use a straight-forward generalization of the quantum Coulomb glass to spinful electrons. It is defined on a regular hypercubic lattice with $`g=L^d`$ ($`d`$ is the spatial dimensionality) sites occupied by $`N=N_{}+N_{}=2Kg`$ electrons ($`0<K<1`$). To ensure charge neutrality each lattice site carries a compensating positive charge of $`2Ke`$. The Hamiltonian is given by

$$H=t\underset{ij,\sigma }{}(c_{i\sigma }^{}c_{j\sigma }+h.c.)+\underset{i,\sigma }{}\phi _in_{i\sigma }+\frac{1}{2}\underset{ij,\sigma ,\sigma ^{}}{}(n_{i\sigma }K)(n_{j\sigma ^{}}K)U_{ij}+U_H\underset{i}{}n_in_i$$
(1)
where $`c_{i\sigma }^{}`$ and $`c_{i\sigma }`$ are the creation and annihilation operators at site $`i`$ and spin $`\sigma `$, and $`ij`$ denotes all pairs of nearest-neighbor sites. $`t`$ is the strength of the hopping term, i.e., the kinetic energy, and $`n_{i\sigma }`$ is the occupation number of spin state $`\sigma `$ at site $`i`$. We parametrize the interaction $`U_{ij}=e^2/r_{ij}`$ by its value $`U`$ between nearest-neighbor sites. The Coulomb repulsion between two electrons at the same site is described by the Hubbard interaction $`U_H`$ The random potential values $`\phi _i`$ are chosen independently from a box distribution of width $`2W_0`$ and zero mean. The boundary conditions are periodic and the Coulomb interaction is treated in the minimum image convention (which implies a cut-off at a distance of $`L/2`$).

A numerically exact solution of a quantum many-particle system requires the diagonalization of a matrix whose dimension increases exponentially with system size. This severely limits the possible sample sizes. In order to overcome this problem we have developed the HFD method. The basic idea is to work in a truncated Hilbert space consisting of the corresponding Hartree-Fock (Slater) ground state and the low-lying excited Slater states. For each disorder configuration three steps have to be performed: (i) find the Hartree-Fock solution of the problem, (ii) determine the $`B`$ Slater states with the lowest energies, and (iii) calculate and diagonalize the Hamilton matrix in the subspace spanned by these states. The number $`B`$ of new basis states determines the quality of the approximation, reasonable values have to be found empirically.

## III Results

The conductance of a quantum many-particle system can be obtained from linear-response theory. It is essentially determined by the current-current correlation function of the ground state. The real (dissipative) part of the conductance (in units of $`e^2/h`$) is given by the Kubo-Greenwood formula ,

$$\mathrm{Re}G^{xx}(\omega )=\frac{2\pi ^2}{\omega }\underset{\nu }{}|0|j^x|\nu |^2\delta (\omega +E_0E_\nu )$$
(2)
where $`\omega `$ denotes the frequency. $`j^x`$ is the $`x`$ component of the current operator and $`\nu `$ denotes the eigenstates of the Hamiltonian. Eq. (2) describes an isolated system while in a real d.c. transport experiment the sample is connected to contacts and leads. This results in a finite life time $`\tau `$ of the eigenstates leading to an inhomogeneous broadening $`\gamma =\tau ^1`$ of the $`\delta `$ functions in (2) . To suppress the discreteness of the spectrum of a finite system, $`\gamma `$ should be at least of the order of the single-particle level spacing. For our systems this requires a comparatively large $`\gamma 0.05`$. We tested different $`\gamma `$ and found that the conductance values depend on $`\gamma `$ but the qualitative results do not.

The main results of this paper are summarized in Fig. 1 which shows the disorder and interaction dependence of the typical conductance for both spinless and spinful electrons on a two-dimensional lattice of $`4\times 4`$ sites.

The qualitative behavior in both cases is similar: In the strongly localized regime (small $`t`$) a moderate interaction delocalizes the electrons while a sufficiently strong interaction always strongly suppress the conductance. This is the precursor of a Wigner crystal or Wigner glass. With decreasing disorder (increasing $`t`$) the interaction-induced enhancement of the conductance also decreases and eventually vanishes. The behavior of the conductance can be attributed to the competition of two effects: First, the interactions destroy the phase of the electrons and thus the interference necessary for localization. This is particularly effective if the localization length is small to begin with. Second, the interactions introduce an additional source of randomness which tends to increase the localization.

A comparison of the cases of spinless and spinful electrons shows that the interaction induced delocalization is significantly larger for spinful electrons. Moreover, the enhancement seems to vanish at a larger kinetic energy (which we did not reach in the simulations). A systematic investigation of the dependence of the conductance on $`U`$ and $`U_H`$ will be published elsewhere.

In summary, we have studied the influence of electron-electron interactions on Anderson localization for spinless and spinful electrons in two dimensions. For strong disorder moderate interactions significantly enhance the transport. This enhancement is much stronger for spinful than for spinless electrons. Identifying a real phase transition and thus establishing a connection between these findings and the experiments on Si-MOSFETs requires a finite-size scaling analysis of the conductance. This remains a task for the future.

This work was supported in part by the NSF under grant no. DMR–98–70597 and by the DFG under grant no. SFB 393/C2. T.V. thanks the Aspen Center for Physics and the University of Oregon for hospitality during the completion of this paper.


text:
# Interaction effects on the conductance peak height statistics in quantum dots

## Abstract

A random interaction matrix model is used to study the statistics of conductance peak heights in Coulomb blockade quantum dots. When the single-particle dynamics conserves time-reversal symmetry, the peak height statistics is insensitive to the interaction strength. But when the single-particle dynamics breaks time-reversal symmetry, the peak height statistics exhibits a crossover from unitary to orthogonal symmetry as the interaction strength increases. This crossover is driven by the time-reversal symmetry of the interaction. Our random interaction matrix model describes features of both the measured peak height and peak spacing statistics.

Random matrix theory (RMT) provides a useful tool for describing the universal statistical fluctuations of the spectrum and eigenfunctions of a quantum system whose associated classical dynamics is chaotic. RMT has been successfully applied to the study of mesoscopic phenomena in quantum dots – submicron 2D devices where electrons are confined by electrostatic potentials . In dots with irregular shapes the single-electron dynamics is mostly chaotic, and the use of RMT is justified within a single-particle framework. In dots that are strongly coupled to leads, or open dots, the approximation of non-interacting quasi-particles is reasonable, and interactions can be considered indirectly through their effect on the electron coherence time. When a finite dephasing time is included, RMT can describe quantitatively the observed conductance statistics in open dots .

As the dot-leads coupling is made weaker, the charge on the dot becomes quantized, and electron-electron interactions cannot be ignored. In such almost closed dots, the conductance displays sharp peaks versus the gate voltage. The measured peak height distributions  were found to agree with RMT predictions . However, it is not clear why RMT should describe the peak height statistics in dots with strong electron-electron interactions. For example, the spacings between peaks were observed to have a Gaussian-like distribution  and not the Wigner-Dyson distribution that is expected in the constant interaction model plus single-particle RMT. Numerical simulations of an Anderson model of a small disordered dot ($`10`$ electrons) with Coulomb interactions also found Gaussian peak spacing distributions , while the peak height distributions showed only a weak dependence on interactions .

RMT is not limited to non-interacting systems. Indeed it was first introduced to explain the neutron resonance data in the compound nucleus. But the neutron resonances are measured at finite excitations, while the linear conductance experiments in quantum dots probe the ground state of the system with a varying number of electrons. To understand the statistics that emerges from the interplay between one-body chaos and interactions in these dots, it is necessary to construct a random matrix model that includes interactions explicitly and reduces to one-body RMT in the absence of interactions . We can break the time-reversal symmetry of the one-body Hamiltonian (e.g., with a magnetic field), but the time-reversal symmetry of the two-body interaction should be preserved. An important question is whether such a random matrix model can reproduce both the measured peak height and peak spacing statistics.

A random interaction matrix model (RIMM) was introduced in Ref.  to study the peak spacing distribution in chaotic dots. The model combines a random two-body interaction  with a random one-body Hamiltonian, and describes a crossover of the peak spacing statistics from a Wigner-Dyson distribution to a Gaussian distribution in terms of a parameter that measures the fluctuations of the interaction matrix elements. Here, we use the RIMM to study the effects of interactions on the conductance peak height statistics. For a GOE one-body statistics, we find that the partial width (to decay into one of the leads) has a GOE Porter-Thomas distribution independently of the interaction strength. However, for a GUE one-body statistics, we find a crossover from a GUE to a GOE Porter-Thomas distribution as the interaction strength increases. It is well known that, in the absence of interactions, a complete crossover from GOE to GUE statistics can be induced by a magnetic field. Our results indicate that this crossover is not complete when strong interactions are turned on because of a competition between an asymptotic symmetry of the one-body Hamiltonian and the GOE symmetry of the time-reversal–conserving two-body interaction. Finally, for the case of GUE one-body statistics, we find that a Gaussian-like peak spacing distribution develops at interaction strengths for which the peak height statistics is still close to GUE. This explains the observed peak spacing statistics  and peak height statistics  within a single random matrix model.

The $`n`$-th conductance peak height $`G_n`$ corresponds to the tunneling of an electron into a dot with $`n1`$ electrons to form a dot with $`n`$ electrons. At low temperatures, $`G_n=(e^2/h)(\pi \overline{\mathrm{\Gamma }}/4kT)g_n`$, where

$`g_n={\displaystyle \frac{1}{\overline{\mathrm{\Gamma }}}}{\displaystyle \frac{\mathrm{\Gamma }_n^l\mathrm{\Gamma }_n^r}{\mathrm{\Gamma }_n^l+\mathrm{\Gamma }_n^r}}.`$ (1)
$`\mathrm{\Gamma }_n^{l(r)}`$ is the partial width of the ground state of the $`n`$-electron dot to decay into an electron in the left (right) lead and the ground state of the dot with $`n1`$ electrons:

$$\mathrm{\Gamma }_n\left|\mathrm{\Phi }_{\mathrm{g}.\mathrm{s}.}(n)|\psi ^{}(𝒓)|\mathrm{\Phi }_{\mathrm{g}.\mathrm{s}.}(n1)\right|^2.$$
(2)
$`\psi ^{}(𝒓)`$ is the creation operator of an electron at the point $`𝒓`$ \[$`𝒓=𝒓_{l(r)}`$ for the left (right) point contact\], and $`\mathrm{\Phi }_{\mathrm{g}.\mathrm{s}.}(n)`$ is the ground state wavefunction of the $`n`$-electron dot. For non-interacting electrons, $`\mathrm{\Phi }_{\mathrm{g}.\mathrm{s}.}(n)`$ is a Slater determinant of the $`n`$ lowest single-particle eigenfunctions in the dot, and Eq. (2) reduces to $`\mathrm{\Gamma }_n|\varphi _n(𝒓)|^2`$, where $`\varphi _n`$ is the $`n`$-th single-particle wavefunction. If the single-particle dynamics is chaotic, $`|\varphi _n(𝒓)|^2`$ satisfies Porter-Thomas statistics, leading to universal distributions of the conductance peak heights (1) that are sensitive only to the underlying symmetry class .

To determine how interactions might modify the Porter-Thomas statistics of the partial widths, we use the RIMM of spinless interacting fermions

$`H={\displaystyle \underset{ij}{}}h_{ij}a_i^{}a_j+{\displaystyle \frac{1}{4}}{\displaystyle \underset{ijkl}{}}\overline{u}_{ijkl}a_i^{}a_j^{}a_la_k.`$ (3)
The one-body matrix elements $`h_{ij}`$ are chosen from the appropriate Gaussian random matrix ensemble, while the anti-symmetrized two-body matrix elements $`\overline{u}_{ij;kl}=u_{ij;kl}u_{ij;lk}`$ form a GOE in the two-particle space

$`P(h)e^{\frac{\beta }{2a^2}\mathrm{Tr}h^2};P(\overline{u})e^{\mathrm{Tr}\overline{u}^2/2U^2}.`$ (4)
The states $`|i=a_i^{}|0`$ describe a fixed basis of $`m`$ single-particle states. $`h`$ is an $`m\times m`$ GOE (GUE) matrix when the single-particle dynamics conserves (breaks) time-reversal symmetry. The two-body interaction is assumed to preserve time-reversal symmetry and forms a GOE, irrespective of the symmetry of the one-body Hamiltonian.

To calculate the statistics of the partial widths $`\mathrm{\Gamma }_n`$, we expand $`\psi ^{}(𝒓)=_i\psi _i(𝒓)a_i^{}`$, where $`\psi _i(𝒓)𝒓|i`$ is the wavefunction of the fixed state $`|i`$. It follows from the orthogonal invariance of the ensemble (4) that the statistics of $`\mathrm{\Gamma }_n`$ is identical to the statistics of

$$\mathrm{\Gamma }_n^i\left|\mathrm{\Phi }_{\mathrm{g}.\mathrm{s}.}(n)|a_i^{}|\mathrm{\Phi }_{\mathrm{g}.\mathrm{s}.}(n1)\right|^2$$
(5)
for any $`i`$.

For each realization $`H`$ of the ensemble (3), we calculate the ground states for $`n1`$ and $`n`$ electrons, and compute $`\mathrm{\Gamma }_n^i`$ using (5). The distributions of the normalized width $`\widehat{\mathrm{\Gamma }}=\mathrm{\Gamma }/\overline{\mathrm{\Gamma }}`$ for a GOE single-particle statistics are shown in Fig. 1 for several values of $`U/\mathrm{\Delta }`$. We show $`P(\mathrm{ln}\widehat{\mathrm{\Gamma }})`$ rather than $`P(\widehat{\mathrm{\Gamma }})`$ in order to display more clearly the small values of $`\widehat{\mathrm{\Gamma }}`$. The distributions are independent of $`U/\mathrm{\Delta }`$ and are well described by the GOE Porter-Thomas distribution $`P_{\mathrm{GOE}}(\mathrm{ln}\widehat{\mathrm{\Gamma }})=(\widehat{\mathrm{\Gamma }}/2)^{1/2}\mathrm{exp}(\widehat{\mathrm{\Gamma }}/2)`$ (solid line). As a reference we also show the GUE Porter-Thomas distribution $`P_{\mathrm{GUE}}(\mathrm{ln}\widehat{\mathrm{\Gamma }})=\widehat{\mathrm{\Gamma }}\mathrm{exp}(\widehat{\mathrm{\Gamma }})`$ (dashed line). The conductance peak heights are calculated from (1) using two uncorrelated “sites” $`i`$ and $`j`$ for the left and right leads. The peak height distributions $`P(g)`$ (for a GOE $`h`$) are shown in the insets of Fig. 1. They are all in good agreement with the GOE distribution $`P_{\mathrm{GOE}}(g)=\sqrt{2/\pi g}e^{2g}`$ (solid lines) irrespective of $`U/\mathrm{\Delta }`$. The dashed lines describe the GUE distribution $`P(g)=4g[K_0(2g)+K_1(2g)]\text{e}^{2g}`$ ($`K_0`$ and $`K_1`$ are Bessel functions).

The results for a GUE one-body statistics are shown in Fig. 2 for the same values of $`U/\mathrm{\Delta }`$ as in Fig. 1. When the interaction strength increases, the width and peak height distributions make a crossover from the corresponding GUE distribution (at $`U=0`$) to the GOE distribution (at large $`U`$). Equivalently, the transition from GOE to GUE statistics due to a time-reversal symmetry-breaking one-body field is not complete because of the competing GOE symmetry of the two-body interaction. The crossover distributions are compared with distributions obtained from an RMT ensemble that describes the crossover between the orthogonal and unitary symmetries. This ensemble is $`H=S+i\alpha A`$ , where $`S`$ and $`A`$ are $`N\times N`$ symmetric and antisymmetric uncorrelated Gaussian matrices and $`\alpha `$ is a parameter . Its wavefunction statistics depends on the parameter $`\lambda \alpha \sqrt{N}/\pi `$. In particular, the width distribution is given in closed form by

$`P_\lambda (\widehat{\mathrm{\Gamma }})={\displaystyle _0^1}dtP_\lambda (t){\displaystyle \frac{1+t^2}{t}}e^{\left(\frac{1+t^2}{t}\right)^2\widehat{\mathrm{\Gamma }}}I_0\left({\displaystyle \frac{1t^4}{t^2}}\widehat{\mathrm{\Gamma }}\right),`$ (6)
where $`I_0`$ is the modified Bessel function of order zero. $`t`$ in (6) is a “shape” parameter that fluctuates in the interval $`[0,1]`$ according to a known distribution

$`P_\lambda (t)=\pi ^2`$ $`(1/t^3t)\lambda ^2e^{\frac{\pi ^2}{2}\lambda ^2\left(t1/t\right)^2}\{\varphi _1(\lambda )+`$ (8)
$`[(t+t^1)^2/41/(2\pi ^2\lambda ^2)][1\varphi _1(\lambda )]\},`$
where $`\varphi _1(\lambda )=_0^1dye^{2\pi ^2\lambda ^2(1y^2)}`$. Eq. (6) reduces to the GOE and GUE Porter-Thomas distributions for $`\lambda =0`$ and $`\lambda 1`$, respectively. In practice, the crossover from the GOE to the GUE already occurs for $`\lambda 1`$.

For each $`U/\mathrm{\Delta }`$, we find $`\lambda `$ by fitting the distribution (6) to the computed width distribution. The distributions (6), shown by the short-dashed lines in Fig. 2, accurately describe the width distributions of the RIMM with a GUE $`h`$. Good agreement with closed RMT expressions  is also obtained for the peak height distributions $`P_\lambda (g)`$ shown in the insets of Fig. 2.

An important issue is the universality associated with the RIMM (3) and (4). The model depends on three parameters: $`m`$, $`n`$, and $`U/\mathrm{\Delta }`$. The top panel of Fig. 3 shows, for a GUE one-body $`h`$, the crossover parameter $`\lambda `$ as a function of $`U/\mathrm{\Delta }`$ for $`m=10`$ and $`n=4,5,6,7`$ (symbols), and for a “reference” case $`m=12`$ and $`n=4`$ (solid line). The curves depend on both $`m`$ and $`n`$, but can all be scaled on the reference curve after scaling the interaction strength by a constant, $`U_{\mathrm{eff}}f(m,n)U/\mathrm{\Delta }`$ (see bottom panel of Fig. 3). The values of the scaling factor $`f(m,n)`$, shown in the bottom inset of Fig. 3, are essentially the same as those needed to obtain a universal peak spacing statistics . We conclude that the peak spacing statistics as well as the peak height statistics in the RIMM depend only on a single parameter $`U_{\mathrm{eff}}`$.

In the range $`U_{\mathrm{eff}}1`$, $`\lambda 1`$ and the statistics is essentially GUE. For $`U_{\mathrm{eff}}11.5`$, the peak height statistics is close to GUE while the peak spacings already follow a Gaussian distribution. This explains the measured RMT-like peak height distributions  and the Gaussian-like shape of the peak spacing distributions  within a single random matrix model. We remark that the small deviations from GUE statistics observed in the experiment of Ref.  and in the calculations of Ref.  in the presence of a magentic field are consistent with a crossover from GUE to GOE.

The average width $`\overline{\mathrm{\Gamma }}`$ is a monotonically decreasing function of $`U/\mathrm{\Delta }`$ and saturates at large values of $`U`$. The top inset of Fig. 3 shows $`\overline{\mathrm{\Gamma }}(U)/\overline{\mathrm{\Gamma }}(0)`$ as a function of $`U/\mathrm{\Delta }`$ for several $`m`$ and $`n`$. This dependence is well described by $`\overline{\mathrm{\Gamma }}(U)\overline{\mathrm{\Gamma }}(\mathrm{})=(\overline{\mathrm{\Gamma }}(0)\overline{\mathrm{\Gamma }}(\mathrm{}))/(1+bU^2/\mathrm{\Delta }^2)`$. The parameter $`\overline{\mathrm{\Gamma }}(\mathrm{})`$ depends on $`m`$ and $`n`$, while $`b`$ is found to be independent of $`m`$ and to depend only weakly on $`n`$.

Next we compare the predictions of the RIMM with those of a model of a quantum dot. We studied a 2D Anderson model with on-site disorder parameter $`W`$ and hopping matrix element $`V=1`$. The electrons are interacting with a Coulomb interaction whose strength over one lattice spacing $`a`$ is $`U_c=e^2/a`$ . We choose periodic boundary conditions in both directions so that the average width and width statistics are independent of the specific sites to which the leads are attached. With these boundary conditions, a background charge (that was included in the calculations of Ref. ) only shifts the energy levels by an $`n`$-dependent constant but does not affect the wavefunctions.

The left panel of Fig. 4 shows the width statistics in the absence of magnetic flux for a $`4\times 5`$ lattice with disorder parameter $`W=5`$ and for $`n=4`$ electrons. The distributions are approximately described by the GOE Porter-Thomas distribution (solid line) irrespective of the value of $`U_c`$ and in agreement with the RIMM. The dependence of the average width $`\overline{\mathrm{\Gamma }}`$ on the interaction strength (not shown) is similar to that observed in the RIMM. The right panel shows width distributions in the presence of magnetic flux $`\mathrm{\Phi }=0.14\mathrm{\Phi }_0`$ and for $`W=3`$. The $`U_c=0`$ distribution (circles) agrees with the GUE Porter-Thomas distribution (dashed line), while the $`U_c=12`$ distribution (triangles) is described by (6) with $`\lambda =0.17`$. This is the distribution obtained in the RIMM for $`U_{\mathrm{eff}}4`$. For weaker interactions the calculated width distributions exhibit some deviations from (6). We note however that, for the small lattices used, it is difficult to find a disorder strength for which the model is in the metallic diffusive regime and displays universal RMT statistics. In particular, in the presence of a magnetic flux we could not find values of $`W`$ for which both the spectral and wavefunction statistics are GUE.

In conclusion, we have investigated the width and peak height statistics in Coulomb-blockade quantum dots using a random interaction matrix model with interactions that preserve time-reversal symmetry. For a GOE one-body symmetry the statistics is insensitive to the interaction strength. However, at strong interactions, a time-reversal symmetry-breaking field leads only to a partial crossover from GOE to GUE statistics. Our random interaction matrix model can reproduce both the observed Gaussian-like shape of the peak spacing distribution and the RMT statistics of the peak height distributions.

We thank Y. Gefen for useful discussions. This work was supported in part by the U.S. DOE grant No. DE-FG-0291-ER-40608. A.W. acknowledges a fellowship from the Studienstiftung des Deutschen Volkes.