Persistent currents in mesoscopic rings: A numerical and renormalization group study V. Meden Institut f¨r Theoretische Physik, Universit¨t G¨ttingen, Bunsenstr. 9, D-37073 G¨ttingen, Germany u a o o arXiv:cond-mat/0209588v1 [cond-mat.str-el] 25 Sep 2002 U. Schollw¨ck o Max-Planck-Institut f¨r Festk¨rperforschung, Heisenbergstr. 1, D-70569 Stuttgart, Germany u o (September 25, 2002) The persistent current in a lattice model of a one-dimensional interacting electron system is systematically studied using a complex version of the density matrix renormalization group algorithm and the functional renormalization group method. We mainly focus on the situation where a single impurity is included in the ring penetrated by a magnetic flux. Due to the interplay of the electron-electron interaction and the impurity the persistent current in a system of N lattice sites vanishes faster then 1/N . Only for very large systems and large impurities our results are consistent with the bosonization prediction obtained for an effective field theory. The results from the density matrix renormalization group and the functional renormalization group agree well for interactions as large as the band width, even though as an approximation in the latter method the flow of the two-particle vertex is neglected. This confirms that the functional renormalization group method is a very powerful tool to investigate correlated electron systems. The method will become very useful for the theoretical description of the electronic properties of small conducting ring molecules. 71.10.Pm, 73.23.Ra, 73.63.-b of the transmission amplitude |T (kF )| of the potential at the Fermi wave vector kF . For a vanishing impurity, i.e. |T (kF )| → 1, I(φ) has a saw tooth like shape, which gets rounded off if |T (kF )| is decreased. In the limit of small |T (kF )|, I(φ) is proportional to |T (kF )| sin φ. For the tight-binding lattice model supplemented by a single weak hopping matrix element the persistent current at half-filling can also be calculated analytically13 and the same characteristics can be found (see Sect. II). Compared to the Fermi liquid behavior of higher dimensional systems a large class of models of homogeneous one-dimensional interacting electrons has a significantly different low-energy physics. These models belong to the Luttinger liquid universality class.15 The low-energy excitations of Luttinger liquids are not given by fermionic quasi-particles, but are of collective, bosonic nature. This leads e.g. to a typical power-law decay of correlation functions. The low-energy physics of Luttinger liquids is characterized by a set of interaction and filling factor dependent parameters.15 For the case of spinless fermions on which we focus any pair of the four parameters vJ (velocity of current excitations), vN (velocity relevant if particles are added), vc (velocity of charge excitations at constant number of particles), and the Luttinger liquid parameter K can be used. Within bosonization and for the impurity free case the persistent current in a Luttinger liquid is periodic in φ and of saw tooth like shape with slope vJ /(πL).15–17 In Sect. III we will compare our numerical results to this prediction. The low-energy physics of Luttinger liquids is strongly affected by the presence of a single impurity.18–23 The problem is usually mapped onto an effective continuum field theory using bosonization, where terms which are I. INTRODUCTION The experimental observation of persistent currents in mesoscopic metallic and semiconducting rings pierced by a magnetic flux1–6 has led to many theoretical investigations focusing on the interplay of electron-electron interaction and disorder in such systems.7 This interplay is considered to be one of the possible reasons for the large current observed in the experiments.7 Despite these studies a quantitative theoretical understanding of the observed amplitude of the currents for the three-dimensional rings is still missing. To gain theoretical insight the simplified situation of one-dimensional rings with interaction and disorder was studied using exact diagonalization for systems of very few lattice sites (up to 16)8–10 and the self-consistent HartreeFock approximation.10,11 Here we consider the further simplified problem of the persistent current in a onedimensional ring of interacting electrons in the presence of a single impurity and penetrated by a magnetic flux. Within an effective continuum field theory it has theoretically been investigated using bosonization12 and conformal field theory.13,14 At first glance this problem seems to be of purely academic interest, but the fast progress in the design and manipulation of conducting ring molecules gives a perspective that such systems might be accessible to experiments in the very near future. For a continuum model of non-interacting onedimensional electrons the leading behavior of the persistent current I(φ) in the system size L can be calculated in the presence of an arbitrary potential scatterer.12 It is a periodic function of the flux, vanishes as 1/L, and its shape and size are determined by the absolute value 1 calculate I(φ) with very high precision and for systems of up to N = 128 lattice sites.29 Additionally the functional renormalization group (RG) method introduced recently into the theory of strongly correlated electrons is used.30,31 It has been applied to two-dimensional correlated electron systems32 and one-dimensional Luttinger liquids.33–35 We have used this method before to study a local observable33,34 - the local spectral weight close an impurity - which is also dominated by the interplay of electron-electron interaction and impurity. In contrast the persistent current is a property of the entire system. Within the RG approach the flow equations are closed by neglecting the flow of the two-particle vertex. Nonetheless DMRG and RG agree quantitatively for interactions of the order of the band width. This confirms that the functional RG is a very powerful tool to investigate strongly correlated electrons. For the single impurity case we indeed find that the persistent current vanishes faster then N −1 . To analyze N I(φ) in more detail we expand the current in a Fourier series, demonstrate that for large impurities and very large system sizes the behavior of the first Fourier coefficient is consistent with a power-law decay with exponent −αB , and that the higher order coefficients decay even faster. Thus in the N → ∞ limit I(φ) is proportional to sin φ with an amplitude which vanishes as N −αB −1 . However, this universal bosonization prediction only holds for very long chains respectively very large impurities. For smaller systems and impurities the asymptotic limit is not reached and the current displays a more complex behavior. It can quantitavely be described using the functional RG method. This makes the method an ideal tool to investigate the electronic properties of small, one-dimensional molecular rings which might be experimentally accessible very soon. expected to be irrelevant in the low-energy limit are neglected.18–22 Within this field theory the leading L dependence of the persistent current was obtained by an additional self-consistent approximation using the analogy to the problem of quantum coherence in a dissipative environment.12 This approach gives a current which for large L vanishes as L−αB −1 , with αB = 1/K − 1, and independent of the bare impurity strength is of purely sinusoidal shape.12 For repulsive interaction one has K < 1 and thus αB > 0. αB is also the exponent of the powerlaw suppression (as a function of energy) of the local spectral weight close to an open boundary and the chemical potential.21,24 This explains the index B which stands for boundary. Many of the bosonization results for observables which are dominated by the interplay of a single impurity and the electron-electron interaction can be understood in terms of a single particle picture, which for Luttinger liquids has of course to be used with caution: For large L the effective transmission amplitude near the Fermi points is suppressed with respect to the non-interacting transmission by a factor of L−αB . Combining this with the behavior of the persistent current in the non-interacting case the above result, obtained without the use of the single particle language, can be derived. For several reasons it is desirable to directly show the above behavior of the persistent current in microscopic lattice models avoiding bosonization. Mapping such models on the field theory involves approximations. Their validity can be questioned and they lead to a loss of information; the scales of the microscopic model and corrections to the expected power-law scaling hidden in the irrelevant terms are lost. Furthermore even within the field theory an additional approximation is necessary to determine the leading (in the system size) behavior of the persistent current. A knowledge of the detailed shape of I(φ) beyond the leading behavior and for microscopic models is of special importance if one wants to compare theoretical results to experiments. Several attempts have been made in this direction using the model of spinless fermions with nearest neighbor hopping and interaction. Here we also focus on this model. The persistent current can be calculated by taking the derivative of the groundstate energy E0 (φ) with respect to the flux φ penetrating the ring. Instead of calculating the full functional form of I(φ) the so called phase sensitivity ∆E0 , which is the difference of the groundstate energy at flux 0 and π has numerically been determined using the density matrix renormalization group (DMRG) method.25,26 The phase sensitivity can be considered a crude measure for the persistent current, but of course does not contain information on the detailed shape of the current as a function of φ. Very recently ∆E0 has again been studied numerically using DMRG and the quantum Monte Carlo method.27 ∆E0 instead of I(φ) (which implies calculating E0 (φ) for several φ and numerically taking the derivative) is calculated because this way the hamiltonian matrix remains real and the numerical effort is reduced considerably.28 Here we use complex DMRG to II. THE MODEL The Hamiltonian for the impurity free ring penetrated by a magnetic flux φ (measured in units of the flux quantum φ0 = hc/e) with nearest neighbor interaction U is given by N c† cj+1 eiφ/N + c† cj e−iφ/N j j+1 H=− j=1 N nj nj+1 , +U (1) j=1 in standard second-quantized notation. The hopping matrix element and the lattice constant are set to one and periodic boundary conditions are used. For U = 0 the groundstate energy of the model can easily be calculated. At temperature T = 0 and for fixed N the persistent current follows from E0 (φ) by taking the derivative with respect to φ 2 I(φ) = − dE0 (φ) dφ obtain I(φ) it is not necessary to use this field theoretical approach since Eq. (7) of Ref. 13 which determines the allowed wave vectors can be solved directly. Written in terms of |T (kF )| the current is given by (2) as can be shown using the Hellman-Feynman theorem. To leading order in 1/N this gives I(φ) = − vF × πN I(φ) = φ for NF odd and − π ≤ φ < π φ − π for NF even and 0 ≤ φ < 2π , which is exactly the expression obtained in Ref. 12 for the continuum model. Later we will be interested in the small |T (kF )| limit. Therefore we expand Eq. (7) up to third order in |T (kF )| where vF denotes the Fermi velocity and NF is the number of particles.36 Both functions have to be continued periodically. Eq. (2) also holds for non-vanishing interaction and if impurity terms are added. The above evenodd effect can also be observed if interaction and impurities are included. We from now on only considere even NF . At φ = 0 the Luttinger liquid parameter K and the velocity of current excitations vJ of the model Eq. (1) can be determined using the Bethe ansatz.37 For halffilling the resulting integral equations have been solved analytically with the results K= 2 arccos (−U/2) π N I(φ) = −1 1 − (U/2)2 2 arccos (−U/2) π − arccos (−U/2) π (3) −1 . (4) For this filling the model is a Luttinger liquid for −2 < U ≤ 2. At |U | = 2 the model shows phase transitions to a charge density wave groundstate (U = 2) and a phase separated state (U = −2). To leading order in 1/N the bosonization prediction for the persistent current in a homogeneous Luttinger liquid is a saw tooth like curve with slope −vJ /(πN ), i.e. compared to the non-interacting case vF is replaced by the velocity determining the current excitations. To compare our data with this result also away from half-filling we numerically solved the Bethe ansatz integral equations following Refs. 37 and 38 and determined vJ (see Sect. III). In our study we will always stay in the Luttinger liquid phase. To H we add a hopping impurity Hh = (1 − ρ) c† c1 eiφ/N + c† cN e−iφ/N 1 N . (5) with ρ between 0 (no hopping between sites N and 1) and 1 (no impurity). In the non-interacting limit, at halffilling, and at wavevector kF the transmission coefficient of the hopping impurity (for N → ∞) is given by39 |T (kF )|2 = 4ρ2 , (1 + ρ2 )2 vF |T (kF )| 1 + |T (kF )|2 /8 sin (φ) 2 vF + |T (kF )|2 sin (2φ) 2π vF + |T (kF )|3 sin (3φ) + O |T (kF )|4 , 16 (8) which at the same time gives an expansion in a Fourier series. The expansion explicitly shows that in the limit |T (kF )| → 0 of a strong impurity the current becomes more and more of sinusoidal shape. A simple approximation which allows to study the effect of the interaction and impurity on the persistent current simultaneously is the Hartree-Fock approximation. For the bulk properties of homogeneous one-dimensional correlated electron systems this approximation does not capture any of the Luttinger liquid features and is thus of very limited usefulness. Furthermore, when applying the self-consistent Hartree-Fock approximation to a model with a single impurity, the self-consistent iterative solution of the Hartree-Fock equations will drive the system into a charge density wave groundstate with a finite single-particle gap which is qualitatively incorrect since a single impurity cannot change bulk properties of the system. Nevertheless self-consistent Hartree-Fock has been used to determine the persistent current in onedimensional, interacting, and disordered rings10,11 and also for the single impurity case.11 In contrast if one is interested in the local properties close to a boundary or impurity non-self-consistent Hartree-Fock provides useful informations.24,33,34 In Sect. V we also present results for the persistent current calculated within the nonself-consistent Hartree-Fock approximation. They were obtained by numerically determining the groundstate of H + Hh for U = 0 and fixed φ and N . From this the expectation values nj 0 and c† cj , which determine j+1 0 the mean-field hamiltonian HMF , can be calculated. HMF can then be diagonalized numerically. There are two possibilities to determine the groundstate energy within this approximation. One can either take the Slater determinant expectation value of H + Hh using the groundstate of HMF or determine the energy via the one-particle propagator40 by the formula and (vJ = vc K) vJ = π vF arccos (|T (kF )| cos [φ − π]) |T (kF )| sin φ , (7) πN 1 − |T (kF )|2 cos2 φ (6) which provides us with a measure for the strength of the impurity. For U = 0, NF = N/2, and to leading order in 1/N the persistent current of the hamiltonian H + Hh has been calculated using conformal field theory.13 To 3 N H + Hh MF = d 1 c† (τ ′ ) cj (τ ) − lim j ′ →τ τ 2 dτ j=1 DMRG (n, U ) extracted as fit parameters. In Fig. 2 the vJ from such fits are compared to the exact vJ (n, U ) obtained from the Bethe ansatz. Data for n = 1/4 (with N = 64) and n = 1/2 (with N = 60) and different U are shown. For both fillings the Bethe ansatz and DMRG results are indistinguishable on the scale of the plot. The relative error of the DMRG data is of the order of 10−4 . This confirms the bosonization prediction, shows that the complex DMRG can be applied successfully, and that the current velocity can to very high precision be extracted from finite size data as small as 60 lattice sites without using finite size scaling. In Ref. 25 data for the phase sensitivity of the translational invariant model at half-filling were obtained using DMRG. Our results for DMRG vJ (n = 1/2, U ) are consistent with the ones presented there. As seen in Fig. 2 the dependence of vJ on U gets weaker if the filling is getting smaller. This happens because for smaller fillings the lattice model is closer to the electron gas model with quadratic dispersion. The latter model is Galilean-invariant which leads to vJ = vF independent of the interaction. MF N eiφ/N c† cj+1 j − j=1 MF +(1 − ρ) eiφ/N c† c1 N + c.c. MF + c.c. (9) where the expectation values are taken using the ground(†) state of HMF . cj (τ ) denotes the annihilation (creation) operator in the imaginary time Heisenberg representation. Transforming the first term into the Matsubara NF frequency representation it can be written as l=1 εMF , l MF with the eigenenergies εl of HMF . Only if the selfconsistent Hartree-Fock approximation is used, which for the reason given above we do not consider in the present context, both possibilities give the same result. Formally both approximations to the energy are correct to leading order in U , but it turned out that the latter method gives better results compared to the high-precision DMRG data. For the results presented in Sec. V we thus used Eq. (9). There it will be shown that the currents calculated using the Hartree-Fock approximation are qualitatively wrong, since in this method correlation effects are neglected which in the present context are of great importance. -28 -28.1 U=0.5 U=1 U=1.5 E0 -28.2 III. COMPLEX DMRG -28.3 Using the DMRG algorithm the groundstate energy of an interacting one-dimensional many fermion system can be calculated to high precision.41 To determine the persistent current using Eq. (2) for the hamiltonian H + Hh [Eqs. (1) and (5)] the DMRG procedure has to be generalized to complex hamiltonian matrices. Calculation time scales up by about a factor of 4, memory usage by a factor of 2. This limits the performance of the method, which is however numerically very stable as hermiticity is conserved. We have kept up to 400 states, ensuring that energies and derived currents are essentially exact. As a test of our program we first studied the impurity free case given by the hamiltonian Eq. (1). We calculated the groundstate energy as a function of φ for 0 ≤ φ ≤ π. Results for quarter-filling n = 1/4, N = 64 and U = 0.5, 1, 1.5 are shown in Fig. 1. Bosonization predicts that to leading order in 1/N the current is given by I(φ) = − vJ (φ − π) πN -28.4 -28.5 0 1 1.5 2 2.5 3 φ FIG. 1. Groundstate energy as a function of the flux φ for quarter-filling, N = 64, and different U . The symbols are DMRG data and the lines are quadratic fits (see the text). 2 vJ 1.8 1.6 (10) and thus the groundstate energy by vJ (φ − π)2 . E0 (φ) = const. + 2πN 0.5 n=1/2, Bethe ansatz n=1/4, Bethe ansatz n=1/2, DMRG n=1/4, DMRG 1.4 (11) 0 0.5 1 1.5 2 U FIG. 2. Current velocity vJ (n, U ) as obtained from the Bethe ansatz and from the DMRG. As shown in Fig. 1 the DMRG data nicely lie on curves of this form with the constant and the current velocity 4 The set of Eqs. (13) and (14) is complemented by a differential equation for the “zero-particle” vertex γ0 (with the Boltzmann constant kB set to one) IV. RG METHOD Besides DMRG we will also use the functional RG method in the version using one-particle irreducible vertex functions.30,31 In collaboration with W. Metzner and K. Sch¨nhammer we have successfully applied this o method in the past to determine the local spectral weight of a Luttinger liquid close to an open boundary and an impurity.33,34 In the method one introduces a cut-off parameter Λ in the non-interacting propagator G0 cutting out degrees of freedom on energy scales less than Λ and derives an exact hierarchy of coupled differential flow equations for the one-particle irreducible vertex functions by differentiating with respect to Λ, where Λ flows from ∞ to 0. For φ = 0 the flow equation for the selfenergy42 of the hamiltonian H + Hh has been given in Ref. 34. As in this reference we here also neglect the flow of the two-particle vertex which closes the set of differential equations; this leads to a energy independent selfenergy. Within the approximation the results obtained are at least correct to leading order in U , but our work presented in Refs. 33 and 34 shows that the Luttinger liquid scaling of the impurity (i.e. the transmission) is included, which makes the RG a promising method also in the present context. If a frequency cutoff G0,Λ (iω) = Θ(|ω| − Λ)G0 (iω) lim T T →0 1 Tr 2π U d Λ Σ = dΛ j,j±1 2π ˆ E0 (φ) − µ N T →0 −1 − ΣΛ (13) (14) −1 . (15) In Eqs. (13) and (14) the flow from Λ = ∞ down to a scale Λ0 much larger then the band width has already been included. The flow is continued from Λ0 downwards with the initial conditions ΣΛ0 = U , 1 ≤ j ≤ N j,j ΣΛ0 = 0 , j,j±1 1