Thermoelectric and Seebeck coefficients of granular metals
aa r X i v : . [ c ond - m a t . m e s - h a ll ] F e b Thermoelectric and Seebeck coefficients of granular metals
Andreas Glatz and I. S. Beloborodov Materials Science Division, Argonne National Laboratory, Argonne, Illinois 60439, USA Department of Physics and Astronomy, California State University, Northridge, California 91330, USA (Dated: June 20, 2018)In this work we present a detailed study and derivation of the thermopower and thermoelectriccoefficient of nano-granular metals at large tunneling conductance between the grains, g T ≫ E c ≫ T /g T > δ , where δ is the mean energy level spacing for asingle grain and E c its charging energy. We show that the electron-electron interaction leads to anincrease of the thermopower with decreasing grain size and discuss our results in the light of futuregeneration thermoelectric materials for low temperature applications. The behavior of the figure ofmerit depending on system parameters like grain size, tunneling conductance, and temperature ispresented. PACS numbers: 73.63.-b, 72.15.Jf, 73.23.Hk
I. INTRODUCTION
Thermoelectric materials with high efficiency is a ma-jor research area in condensed matter physics and mate-rials science for several decades now. Due to recent ad-vances in nano-fabrication, these materials promise nextgeneration devices for conversion of thermal energy toelectrical energy and vice versa. A measure for the per-formance or efficiency of a thermoelectric material is thedimensionless figure of merit , usually denoted as ZT ,where T is the temperature. It depends on the ther-mopower or Seebeck coefficient, and the electric and ther-mal conductivities . However, the Wiedemann-Franz law has defeated much progress in increasing theperformance of bulk materials, since it directly relateselectric and thermal conductivity whereas the figure ofmerit is proportional to the quotient of both. The See-beck coefficient of a material measures the magnitude ofan induced thermoelectric voltage in response to a tem-perature difference across that material. If the tempera-ture difference ∆ T between the two ends of a material issmall, then the thermopower, S , of the material is definedas S = − ∆ V / ∆ T , where ∆ V is the voltage differenceacross a sample.In general, thermoelectric devices are used as con-verters for either electrical power into heating/cooling( Peltier effect ) or sources of different temperature intoelectricity (
Seebeck effect ). These devices are usuallymuch simpler, especially without moving parts, than con-ventional devices, e.g. two-phase compressors for cooling,and therefore more reliable. However, for both effectsthe materials need to have good electrical conductivityto minimize ohmic heating and at the same time to bebad thermal conductors to avoid thermal equilibration ofthe temperature gradient. Therefore, the aim is to cre-ate materials which optimize these parameters togetherwith the thermopower. Currently thermoelectric devicesbased on p- and n-type-doped semiconductor junctionsarchive only about 12% of the maximal theoretical effi- ciency (as compared to 60% in conventional cooling sys-tems) .To be competitive compared with conventional refrig-erators, one must develop thermoelectric materials withlarge ZT . The highest figure of merit for bulk ther-moelectric materials is about 1, but in order to matchthe efficiency of mechanical systems, ZT ∼ .However, for ZT & . Although it is possiblein principle to develop homogeneous materials with thatlarge figure of merit, there are no concrete devices on thehorizon. Especially promising for further improvementin efficiency are inhomogeneous/granular thermoelectricmaterials in which one can directly control the systemparameters. In Ref. [11] a figure of merit at 300 K of 2 . ZT of 3 . K for a bulk material withnanoscale inclusions were reported.Overall, recent years have seen a remarkable progressin the design of granular conductors with controllablestructure parameters. Granules can be capped with or-ganic (ligands) or inorganic molecules which connect andregulate the coupling between them. Altering the sizeand shape of granules one can regulate quantum con-finement effects. In particular, tuning microscopic pa-rameters one can vary the granular materials from be-ing relatively good metals to pronounced insulators as afunction of the strength of electron tunneling couplingsbetween conducting grains . This makes granularconductors a perfect exemplary system for studying ther-moelectric and related phenomena.All these experimental achievements and technologicalprospects call for a comprehensive theory able to providequantitative description of not only the electric but alsothermoelectric properties of granular conductors, whichcan in future serve as basis for a clever design of devicesfor a new generation of nano-thermoelectrics.Most theoretical progress so far was archived by nu-merical solution of phenomenological models . How-ever, no analytical results obtained from a microscopicmodel for coupled nanodot/grain systems is available upto now. Thus, the fundamental question that remainsopen is how thermoelectric coefficient and thermopowerbehave in nanogranular thermoelectric materials. Here,we make a step towards answering this question for gran-ular metals at intermediate temperatures by generalizingour approach recently developed for the description ofelectric and heat transport . In particular, we will an-swer the question to what extend quantum and confine-ment effects in nanostructures are important in changing ZT .In this paper we investigate the thermopower S , ther-moelectric coefficient η , and the figure of merit ZT ofgranular samples focusing on the case of large tunnel-ing conductance between the grains, g T ≫
1. With-out Coulomb interaction the granular system would be agood metal in this limit and our task is to include charg-ing effects in the theory. We furthermore restrict ourconsiderations to the case of intermediate temperatures E c ≫ T /g T > δ, (1)where δ is the mean level spacing of a single grain and E c the charging energy. The left inequality means thatthe temperature is not high enough such that Coulombeffects are pronounced. It also allows us to perform allcalculations up to logarithmic accuracy. The right in-equality allows us to consider the electronic motion as co-herent within the grains, however this coherence does notextend to scales larger than the size of a single grain .In Ref. [21] we presented a few major results which weextend here and describe the derivations in much moredetail.The paper is organized as follows: In Sec. II we sum-marize our main results and discuss their range of ap-plicabilities, in Sec. III we introduce the model, and inSec. IV we outline the derivation of thermoelectric co-efficient of granular metals in and without the presenceof interaction which is the main result of this paper. Inthe following section we discuss the behavior of the See-beck coefficient and figure of merit (Sec. V) as a func-tion of sample parameters. Finally, in Sec. VI we discussour findings and present possible further applications ofour method. Important details of our calculations arepresented in several comprehensive Appendixes: in Ap-pendix A we calculate the thermoelectric coefficient ofhomogeneous disorder metals in the absence of interac-tion. In Appendix B we derive the heat and electric cur-rents of granular metals, and in Appendixes C and D wecalculate the thermoelectric coefficient of granular metalswithout and in the presence of interaction, respectively. II. RESULTS AND SUMMARY
In this section we summarize our results and discusstheir range of applicabilities. The main results of ourwork are as follows: (i) We derive the expression for the thermoelectric coefficient η of granular metals that in-cludes corrections due to Coulomb interaction at tem-peratures T > g T δ , where δ is the mean level spacing ofa single grain η = η (0) (cid:18) − g T d ln g T E c T (cid:19) . (2a)Here η (0) = − ( π / eg T a − d ( T /ε F ) , (2b)is the thermoelectric coefficient of granular materials inthe absence of electron-electron interaction with e beingthe electron charge, a the size of a single grain, d =2 , ε F being the Fermienergy, and E c = e /a is the charging energy.The condition for the temperature range of our theoryensures that the argument of the logarithm in Eq. (2a)is much larger than 1, such that all numerical prefactorsunder the logarithm can be neglected. Furthermore, italso defines a critical lower limit for the grain size whenthe charging energy E c becomes of order of the meanenergy level spacing δ .At the temperatures under consideration, the electronmotion is coherent within the grains, but coherence doesnot extend to scales larger than the size a of a singlegrain . Under these conditions, the electric conductivity σ and the electric thermal conductivity κ are given by theexpressions σσ (0) = 1 − ln( g T E c /T ) / (2 πdg T ) , (3a) κκ (0) = 1 − ln[ g T E c /T ]2 πdg T + 12 π g T (cid:26) γ, d = 3ln g T E c T , d = 2 . (3b)where σ (0) = 2 e g T a − d and κ (0) = L σ (0) T, (3c)are the electric (including spin) and thermal conductivi-ties of granular metals in the absence of Coulomb inter-action with L = π / e being the Lorentz number. Wemention that at temperature T > g T δ the correction tothe thermoelectric coefficient, Eq. (2a), has a T ln T de-pendence in both d = 2 , T dependence in all dimensions as well.(ii) Using the above results, we obtain the expressionfor thermopower S of granular metals S = S (0) (cid:18) − π − πg T d ln g T E c T (cid:19) , (4a)where S (0) = − ( π / /e )( T /ε F ) , (4b)is the thermopower of granular metals in the absence ofCoulomb interaction.(iii) Finally, we find the figure of merit to be: ZZ (0) = 1 − π − πg T d ln g T E c T − π g T (cid:26) γ, d = 3ln g T E c T , d = 2 , (5a)where Z (0) T = ( π / T /ε F ) , (5b)is the bare figure of merit of granular materials and γ ≈ .
355 is a numerical coefficient. In Sec. V we present plotsof Z in dependence of various sample parameters. Wefind, that the influence of granularity is most effective forsmall grain sizes and the presence of Coulomb interactiondecreases the figure of merit.At this point we remark that all results are obtainedin the absence of phonons which become relevant only athigher temperatures. At the end of this paper we willbriefly discuss their influence. III. MODEL
We start our considerations with the introduction ofour model. We consider a d − dimensional array of metal-lic grains with Coulomb interaction between electrons.The motion of electrons inside the grains is diffusive, i.e.,the electron’s mean free path ℓ is smaller than the grainsize a , and they tunnel from grain to grain. We assumethat the sample would be a good metal in the absence ofCoulomb interaction. However, we also assume that thetunneling conductance g T is still smaller than the grainconductance g , meaning that the granular structure ispronounced and the resistivity is controlled by tunnelingbetween grains.Each grain is characterized by two energy scales: (i)the mean energy level spacing δ , and (ii) the chargingenergy E c = e /a (for a typical grain size of a ≈ E c is of the order of 2000K) and we assume that thecondition δ ≪ E c is fulfilled.The system of coupled metallic grains is described bythe Hamiltonian ˆ H = P i ˆ H i , where the sum is takenover all grains in the system andˆ H i = X k ξ k ˆ a † ik ˆ a ik + X j = i e ˆ n i ˆ n j C ij + X j,p,q ( t pqij ˆ a † ip ˆ a jq + c.c.) . (6)The first term on the right hand side (r. h. s.) of Eq. (6)describes the i -th isolated disordered grain, ˆ a † i,k (ˆ a i,k ) arethe creation (annihilation) operators for an electron inthe state k and ξ k = k / m − µ with µ being the chemicalpotential.The second term describes the charging energy, C ij isthe capacitance matrix and ˆ n i = P k ˆ a † ik ˆ a ik is the numberoperator for electrons in the i -th grain. The Coulomb in-teraction is long ranged and its off-diagonal componentscannot be neglected. Note that, since metallic grains have an infinite dielectric constant, the effective dielec-tric constant of the whole sample can be considerablylarger than the dielectric constant of its insulating com-ponent. Thus the effective single-grain charging energycan be much less than the electrostatic energy of a singlegrain in vacuum.The last term in Eq. (6) is the tunneling part of theHamiltonian where t ij are the tunnel matrix elementsbetween grains i and j which we consider to be randomGaussian variables defined by their correlators: h t ∗ p q ij t p ′ q ′ ij i = t ij δ p p ′ δ q q ′ , (7) t ij = t = const . The dimensionless tunneling conduc-tance is related to the average matrix elements as g T = 2 π ( t /δ ) . (8)The conductance g T is defined per one spin component,such that, for example, the high temperature (Drude)conductivity of a periodic granular sample is σ (0) =2 e g T a − d . IV. THERMOELECTRIC COEFFICIENT
Having introduced our model in the previous section,we come now to the main methodical part of our work,the derivation of the thermoelectric coefficient η . In gen-eral the three kinetic coefficients are: the electric conduc-tivity σ , the thermoelectric coefficient η , and the thermalconductivity κ . They are related to the Matsubara re-sponse functions L ( αβ ) with α, β ∈ { e, h } as j ( e ) = − ( L ( ee ) / ( e T )) ∇ ( eV ) − ( L ( eh ) / ( eT )) ∇ T, j ( h ) = − ( L ( eh ) / ( eT )) ∇ ( eV ) − ( L ( hh ) /T ) ∇ T. (9a)Here j ( e ) ( j ( h ) ) is the electric (thermal) current and V the electrostatic potential. From Eq. (9a) one finds that σ = L ( ee ) /T, η = L ( eh ) /T , (9b) S = − ∆ V / ∆ T = L ( eh ) / ( T L ( ee ) ) , where the response functions are given by Kubo formulas L ( αβ ) = − ıT ∂a d ∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω → /T Z dτ e ı Ω m τ h T τ j ( α ) ( τ ) j ( β ) (0) i Ω m →− ı Ω+ δ , (10)with T τ being the time ordering operator for the currentswith respect to the imaginary time τ . Thus, to calculatethe thermoelectric coefficient η and thermopower S ofgranular metals one has to know the explicit form of theelectric j ( e ) and thermal j ( h ) currents.The electric current j ( e ) i through grain i is defined as j ( e ) i = X j ˆ ( e ) ij = e d ˆ n i /dt = ıe [ˆ n i , ˆ H ] . (11a) FIG. 1: (color online) Vertices corresponding to the ther-mal current operator, Eqs. (12b): vertex (a) corresponds to b ( h, ij and (b) to b ( h, ij . The solid lines denote the propaga-tor of electrons, the thick wavy line describes Coulomb in-teraction, the tunneling vertices are described by the circles, ω n = πT (2 n +1) and Ω m = 2 πmT are Fermionic and BosonicMatsubara frequencies respectively ( n, m ∈ Z ). Straightforward calculations (see Appendix B) lead toˆ ( e ) ij = ıe X k,q ( t kqij ˆ a † ik ˆ a jq − t qkji ˆ a † jq ˆ a ik ) . (11b)For granular metals the thermal current operator j ( h ) i = X j ˆ ( h ) ij , (12a)can be obtained as follows. The energy content of eachgrain changes as a function of time, such that d ˆ H i /dt = i [ ˆ H i , ˆ H ]. Energy conservation requires that this energyflows to other grains in the system, d ˆ H i /dt ≡ P j ˆ ( h ) ij .Calculating the commutator [ ˆ H i , ˆ H ], we obtain (for de-tails see Appendix B)ˆ ( h ) ij = ˆ ( h, ij + ˆ ( h, ij , (12b)ˆ ( h, ij = ı X k,q ξ k + ξ q h t kqij ˆ a † ik ˆ a jq − t qkji ˆ a † jq ˆ a ik i , (12c)ˆ ( h, ij = − e X m " { ˆ n i ; ˆ ( e ) jm } + C im − { ˆ n j ; ˆ ( e ) im } + C jm , (12d)where { ˆ A ; ˆ B } + denotes the anti-commutator. The con-tribution b ( h, ij is the heat current in the absence ofelectron-electron interaction, while the second term ˆ ( h, ij appears due to Coulomb interaction. Equation (12b) im-plies that the thermal current operator must be associ-ated with two different vertices in diagram representa-tion, Fig. 1. We remark that Eq. (6) suggests also afinite contribution to j ( h ) i proportional to t – which in-deed exists. However, it vanishes when summed over thesample (for details see Appendix B).For large tunneling conductance, the Matsubara ther-mal current - electric current correlator can be analyzedperturbatively in 1 /g T , using the diagrammatic tech-nique discussed in Ref. [18] that we briefly outline below.The self-energy of the averaged single electron Green’sfunction has two contributions: The first one correspondsto scattering by impurities inside a single grain while the FIG. 2: (color online) Diagrams describing the thermoelectriccoefficient of granular metals at temperatures
T > g T δ : dia-gram (a) corresponds to η in Eq. (2a). Diagrams (b)-(d) de-scribe first order corrections to the thermoelectric coefficientof granular metals η (1) in Eq. (25) due to electron-electron in-teraction. The solid lines denote the propagator of electrons,the wavy lines describe effective screened electron-electronpropagator, and the (red) triangles describe the elastic in-teraction of electrons with impurities. The tunneling verticesare described by the circles. The sum of the diagrams (b)-(d)results in the thermoelectric coefficient correction η (1) givenin Eq. (33). second is due to the process of scattering between thegrains. The former results only in a small renormaliza-tion of the relaxation time which depends in general onthe electron energy ω as τ − ω = τ − [1 + ( d/ − ω/ε F ] , (13)which is a result of the renormalization of the density ofstates at the Fermi surface [see Eqs. (A8) and (A9)].In the following we outline the calculation of thermo-electric coefficient η in the non-interacting case and itscorrection due to Coulomb interaction. A detailed deriva-tion of both can be found in Appendices C and D, respec-tively. A. non-interacting case
First, we consider the thermoelectric coefficient η (0) of granular metals in the absence of interaction. Theexpression for the thermoelectric coefficient in linear re-sponse theory is η (0) = ı ∂a d T ∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 Q (0) . (14)Here a is the grain size and Q (0) the correlator of theheat current, −→ ( h, (see Fig. 1a), and electric current, −→ ( e ) , shown in Fig. 2a. Notice, that there is an impor-tant difference between calculations of the thermoelectriccoefficient η and thermopower S and the calculations ofthe electric σ and thermal κ conductivities in Eqs. (3a,3b). Indeed, to calculate σ and κ it was sufficient to ap-proximate the tunneling matrix element t pq as a constant t which is evaluated at the Fermi surface and neglect vari-ations of t pq with energy which occur on the scale T /ε F .However, this approximation is insufficient for calculatingthe thermoelectric coefficient η and thermopower S sincethe dominant contribution to these quantities vanishesdue to particle-hole symmetry such that both quantitiesare proportional to the small parameter T /ε F . Since itis necessary to take into account terms of order of T /ε F in order to obtain a nonzero result for η and S , the cor-responding expansions must be carried out to this orderfor all quantities which depend on energy.For granular metals the important element of the dia-gram is the tunneling matrix elements t kqij describing thecoupling between grains i and j . Therefore we derive anexpression for t kqij in the following, assuming that i and jare nearest neighbor grains and t kqij is independent of theposition in the sample. In order to calculate the energydependence of these elements we assume the tunnelingbarrier between grains to be a delta potential. For theone-particle Hamiltonian b H = − ~ m d dx + λδ ( x ), with λ being the strength of the potential, the transmission ratefor a single particle with energy ε p = ε F + ξ p is T p = (cid:18) mλ ~ ε p (cid:19) − ≃ (cid:18) mλ ~ ε F (1 − ξ p /ε F ) (cid:19) − (15)= T (cid:0) T − − ξ p /ε F (cid:1) − . Here T = mλ ~ ε F is the transmission rate at the Fermilevel and we use the fact that ξ p ≪ ε F . Next, weconsider the case of large barriers, in this regime T p ≃ T (1 + ξ p /ε F ). In granular systems we have many chan-nels and have to consider tunneling processes with en-ergy ξ i in grain i and with ξ j in neighboring grain j : t ∝ N (cid:16) T p i + T p j (cid:17) . Thus, the final expression for tun-neling matrix element is t ( ξ i , ξ j ) = t (cid:18) ξ i + ξ j ε F (cid:19) , (16)where t is the constant tunneling matrix element evalu-ated at the Fermi surface. For convenience we use i = 1and j = 2 in the following. Using Eq. (16) we obtainthe following expression for the correlation function inEq. (14) Q (0) = − set T h−→ n (0) e · −→ n (0) h i X ω n a d +2 (17a) × Z d d p (2 π ) d d d p (2 π ) d (cid:20) ξ + ξ (cid:21) (cid:20) ξ + ξ ε F (cid:21) × G ( p , ω n ) G ( p , ω n + Ω m ) , where −→ n (0) α is the unit vector in direction of the current α ∈ { e, h } , h−→ n (0) e · −→ n (0) h i = 1 /d is the result of averaging over angles, the summation goes over Fermionic Matsub-ara frequencies ω n = 2 πT ( n + 1 / G ( p, ω n ) is theMatsubara Green’s function G ( p, ω n ) = [ ıω n − ξ p ± ı/ (2 τ ω )] − , (17b)with ξ p = ε p − ε F being the electron energy with respectto the Fermi energy [ ε p = p / (2 m )] and the (energy) re-laxation time τ ω is defined in Eq. (13). To shorten thenotation, we neglect the momentum argument in the fol-lowing and attach the grain index to the Green’s function G . The two ξ -factors under the integration in (17a) arisefrom the heat current, Eq. (12b), and the energy correc-tion of the tunneling element, Eq. (16). The momentumintegrals in Eq. (17a) are transformed into energy inte-grals taking into account first order corrections in ξ/ε F of the Jacobian [see Eq. (A7) of Appendix A].We first perform the analytical continuation over theFermionic Matsubara frequencies ω n = 2 πT ( n + 1 /
2) inEq. (17a). In order to accomplish that, the analyticalstructure of diagram (a) in Fig. 2 needs to be analyzed,which gives rise to three different regions for the Matsub-ara summations: I = ] −∞ ; − Ω m ] , I = ] − Ω m ; 0[ , I =[0; ∞ [, in which we can determine whether the Green’sfunction is retarded, G − ( ω ), or advanced, G + ( ω ): S = X n ∈ I G − ( ω n ) G − ( ω n + Ω m ) (18a)= Z − dω πıT tanh (cid:16) ω T (cid:17) G − ( − ıω + ı Ω) G − ( − ıω ) ,S = X n ∈ I G +1 ( ω n ) G − ( ω n + Ω m ) (18b)= Z − dω πıT tanh (cid:16) ω T (cid:17) (cid:2) G − ( − ıω ) G +2 ( − ıω − ı Ω) − G − ( − ıω + ı Ω) G +2 ( − ıω ) (cid:3) ,S = X n ∈ I G +1 ( ω n ) G +2 ( ω n + Ω m ) (18c)= Z − dω πıT tanh (cid:16) ω T (cid:17) G +1 ( − ıω ) G +2 ( − ıω − ı Ω) . To calculate η (0) in Eq. (14) we now consider the deriva-tive of S + S + S with respect to Bosonic frequencies ∂∂ Ω (cid:12)(cid:12) Ω=0 . For brevity we omit arguments − ıω of Green’sfunctions, leading to ∂∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 ( S + S + S ) = (19) Z − dω πıT tanh (cid:16) ω T (cid:17) ∂∂ω (cid:20) τ ω G − G +1 G − G +2 (cid:21) . Next, one can perform the integration over variables ξ and ξ using Eqs. (13) and (20) and the residuum theo-rem Z dξ dξ g ( ξ , ξ ) G − G +1 G − G +2 (20a)= 4 π τ (cid:20) ω + ı/τ + d ε F (2 ω + ı/τ ) (cid:21) , with g ( ξ , ξ ) = ξ + ξ + d ξ + ξ ) ε F . (20b)As a result we obtain the following expression for thederivative of correlation function ∂∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 Q (0) = − πs ıet a d +2 (cid:16) ν (0) d (cid:17) T ε F (21) × Z dω [2 + ( d/ − ı/ ( τ ε F )] ω cosh ( ω/ (2 T ))= − π s ı et a d +2 (cid:16) ν (0) d (cid:17) T ε F . Here all contributions of order 1 /ε F or smaller are ne-glected in the final expression. Substituting this result,Eq. (21), into Eq. (14) we finally obtain the following ex-pression for the non-interacting thermoelectric coefficientof granular metals η (0) = − sπ et a d +2 (cid:16) ν (0) d (cid:17) Tε F . (22)One can re-write this expression using the relations: ν (0) d D d = g T a − d ; ν (0) d = ( δa d ) − ; t = g T δ / (2 π ),where D d is the diffusion constant, g T the tunneling con-ductance, and δ the mean level spacing, giving η (0) = − sπ eg T a − d ( T /ε F ) . (23) B. correction due to Coulomb interaction
Now, we consider the correction η (1) to the thermoelec-tric coefficient of granular metals due to electron-electroninteraction η = η (0) + η (1) , (24)where η (0) is given by Eq. (22). Analogously to η (0) , η (1) can be obtained from η (1) = ı ∂a d T ∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 (cid:16) Q (1) + Q (2) + Q (3) (cid:17) , (25)where the diagrams Q (1) , Q (2) , Q (3) contributing to η (1) are shown in Fig. 2 (b, c, d). Detailed calculations of η (1) are presented in Appendix D. However, here we outlinethe main steps of this derivation. These three diagramsinclude the effect of elastic scattering of electron at im-purities described by diffusons D − = τ ω ( | Ω i | + ǫ q δ ) , (26a)with ǫ q = 2 g T " d − ′ X a cos( −→ q · −→ a ) , (26b) where P ′ a stands for summation over all directions andorientations n ± a −→ e (0) j o , and the effect of the dynamicallyscreened Coulomb potential e V ( q, Ω i ) = D V ( q, Ω i ) D e V ( q, Ω i ) = 2 E c ( q ) τ ω [ | Ω i | + 4 E c ( q ) ǫ q ] [ | Ω i | + ǫ q δ ] , (27a) V ( q, Ω i ) = (cid:18) E c ( q ) + 2 ǫ q | Ω i | + ǫ q δ (cid:19) − , where we use the notation E c ( q ) = e a d − ln( qa ) , d = 1 π/q, d = 22 π/q , d = 3 . (27b)Each diagram in Fig. 2 has also two types of renormalizedinteraction vertices: (i) the inter-grain vertexΦ (1) ω (Ω i ) = Z a d d −→ q (2 π ) d E c ( q ) P ′ a cos( −→ q · −→ a ) τ ω [ | Ω i | + 4 E c ( q ) ǫ q ] [ | Ω i | + ǫ q δ ] , (28)and (ii) the intra-grain vertexΦ (2) ω (Ω i ) = Z a d d −→ q (2 π ) d E c ( q ) 2 dτ ω [ | Ω i | + 4 E c ( q ) ǫ q ] [ | Ω i | + ǫ q δ ] . (29)Explicitly, the contribution Q (1) in Eq. (25) [diagram (b)in Fig. 2] is given by Q (1) = − s d et T a d +2 (cid:16) ν (0) d (cid:17) (30a) × X ω n , Ω i Z dξ dξ g F ( s s s s )1 Φ (1) ω (Ω i ) , where the function g is defined in Eq. (20b) and we usethe notation F ( s s s s )1 = G s ( ω n + Ω m ) G s ( ω n + Ω m + Ω i ) × G s ( ω n ) G s ( ω n + Ω i ) , (30b)with s i = ± denote the analytic structure of the Green’sfunctions implying restrictions on the frequency summa-tion.For the contribution Q (2) in Eq. (25) [diagram (c) inFig. 2] we have the following expression Q (2) = − s d et T a d +2 (cid:16) ν (0) d (cid:17) (31a) × X ω n , Ω i Z dξ dξ g F ( s s s )2 Φ (2) ω (Ω i ) , where we use the notation F ( s s s )2 = [ G s ( ω n + Ω m )] G s ( ω n + Ω m + Ω i ) G s ( ω n ) . (31b)The diagram Q (3) , shown in Fig. 2 (d), describes thecontribution of the correlation function with the interac-tion part of the heat current operator, b ( h, ij [second term FIG. 3: (color online) Plots of the dimensionless figure ofmerit
Z/Z (0) vs. grain size a (in nm) - for different valuesof dimensionless tunneling conductance g T (see legend): theupper panel is for the three-dimensional (3D) case and thelower for the two-dimensional (2D). All curves are plotted for T = 100K. At this temperature, the dimensionless bare figureof merit for granular metals is Z (0) T ≈ − . in the right hand side of Eq. (12b)], and has therefore adifferent structure in comparison with contributions Q (1) and Q (2) : Q (3) = − s d et T a d +2 (cid:16) ν (0) d (cid:17) (32a) × X ω n , Ω i Z dξ dξ g F ( s s s )3 Φ (Ω i , q ) , where F ( s s s )3 = G s ( ω n + Ω m + Ω i ) G s ( ω n + Ω m ) G s ( ω n ) ,g = 2 (cid:18) d ε F ( ξ + ξ ) (cid:19) . (32b)The main contribution to η (1) from diagram Q (3) is ofthe order of ( T /ε F ) , whereas Q (1) and Q (2) have 1 /ε F contributions. Therefore we will not consider diagram Q (3) any further, but keep contributions of order T /ε F only.Thus, the first order interaction corrections to the ther-moelectric coefficient are only generated by diagrams (b)and (c) in Fig. 2. Substituting Eqs. (30a) and (31a)into Eq. (25) after summation over Fermionic, ω n , andBosonic, Ω i , frequencies and analytical continuation weobtained (see Appendix D for details) η (1) = − η (0) πg T (cid:16) a π (cid:17) d Z d d q ln (cid:20) E c ( q ) ǫ q T (cid:21) , (33)where the q -integration goes over the d -dimensionalsphere with radius π/a . Integrating over q in Eq. (33)we obtain the following expressions, neglecting all con-stants under the logarithm, in two (2 D ) and three (3 D )
8 10 Z / Z ( ) g T a=5nma=10nma=25nma=50nm Z / Z ( ) FIG. 4: (color online) Plots of the dimensionless figure ofmerit
Z/Z (0) vs. tunneling conductance g T - for differentvalues of grain sizes a (see legend): the upper panel is forthe three-dimensional (3D) case and the lower for the two-dimensional (2D). All curves are plotted for T = 100K. dimensions: η (1)2 D = − η (0) g T ln E c g T T , η (1)3 D = − η (0) g T ln E c g T T , (34)which lead to Eq. (2a).
V. THERMOPOWER AND FIGURE OF MERIT
Now we have the expressions for all three kinetic coef-ficients σ , κ , and η for granular metals to order T /ε F .Based on these, we can derive other thermodynamicquantities – in particular we discuss the thermopower(Seebeck coefficient) and figure of merit in this section.Both quantities are relevant parameters for thermocou-ples, the thermopower is a measure the change of voltagedue to a temperature gradient and the figure of merit ameasure for the performance of the device.The thermopower is related to the kinetic coefficientsas (see Eq. (9b)) S = η/σ . (35a)Again, we only consider terms up to order T /ε F andobtain the expression S = η (0) σ (0) (cid:18) η (1) η (0) − σ (1) σ (0) (cid:19) , (35b)which results in Eq. (4a).The dimensionless figure of merit is related to the ki-netic coefficients as ZT = T η / ( σκ ) = S σT /κ, (36a)
0 20 40 60 80 100
T [K] Z / Z ( ) FIG. 5: (color online) Plots of the dimensionless figure ofmerit
Z/Z (0) vs. temperature T for g T = 5 (upper twographs) and g T = 2 (lower two graphs) and grain size a = 5 , d = 2 ,
3. Legends are next to the corre-sponding curves. giving ZT = ( S (0) ) σ (0) Tκ (0) (cid:18) − κ (1) κ (0) + σ (1) σ (0) + 2 S (1) S (0) (cid:19) , (36b)resulting in Eq. (5a), which has lowest order ( T /ε F ) [Eq. 5b]. In Eq. (36b) κ (1) , σ (1) , and S (1) are correctionsto the thermal conductivity, electrical conductivity, andthe Seebeck coefficient due to Coulomb interaction, re-spectively. The numerical coefficient two in front of thelast term reflects the fact that the Seebeck coefficient ap-pears squared in the definition of ZT . Using Eqs. (4a)and (36b) one can see that the second term of the right-hand-side of Eq. (5a) originates due to correction to theSeebeck coefficient, Eq. (4a). The origin of the third termin Eq. (5a) is due to correction to the thermal conductiv-ity, κ , last term on the right hand side of Eq. (3b). Wemention that the second term on the right hand side ofEq. (3b) cancels with the correction to the electrical con-ductivity, Eq. (3a), after the substitution into Eq. (36b).In Fig. 3 the dependence of Z on the grain size a andin Fig. 4 on the tunneling conductance g T for two- andthree-dimensional samples are shown. Fig. 5 shows thetemperature dependence of the figure of merit. Theseplots show that the correction term to ZT is most effec-tive for small grain at not very high tunneling conduc-tance and at low temperatures. VI. DISCUSSIONS
In the presence of interaction effects and not very lowtemperatures
T > g T δ , granular metals behave differ-ently from homogeneous disorder metals. However, in theabsence of interactions the result for η (0) below Eq. (2a)[or Eq. C15] coincides with the thermoelectric coefficientof homogeneous disordered metals (see Appendix A), η (0) hom = − (2 / ep F ( τ T ) , (37) with p F being the Fermi momentum. One can expectthat at low temperatures, T < g T δ , even in the presenceof Coulomb interaction the behavior of thermoelectric co-efficient and thermopower of granular metals is similar tothe behavior of η hom and S hom , however this temperaturerange is beyond the scope of the present paper. Our re-sults for thermopower (4a) and figure of merit (5a) showthat the influence of Coulomb interaction is most effectivefor small grains. The thermopower S decreases with thegrain size which is a result of the delicate competition ofthe corrections of thermoelectric coefficient (2a) and theelectric conductivity (3a). In particular, if the numericalprefactor of the correction to η would be slightly smaller,the sign of the correction to S would change.Above we only considered the electron contribution tothe figure of merit. At higher temperatures T > T ∗ ,where T ∗ ∼ q g T c ph /l ph a, (38)is a characteristic temperature with l ph and c ph beingthe phonon scattering length and phonon velocity respec-tively , phonons will provide an independent, additionalcontribution to thermal transport, κ ph = T l ph /c ph . (39)However, the phonon contribution can be neglected fortemperatures g T δ < T < T ∗ . (40)A detailed study of the influence of phonons at high tem-peratures, including room temperature, will be subject ofa forthcoming work.So far, we ignored the fact that electron-electron in-teractions also renormalize the chemical potential µ . Ingeneral, this renormalization may affect the kinetic coef-ficients: the thermal current vertex, Fig 1, as well as theelectron Green’s functions depend on µ . In particular oneneeds to replace ∇ ( eV ) → ∇ ( eV + µ ) in Eq. (9a). To firstorder in the interactions, the renormalization of µ onlyleads to corrections to diagram (a) in Fig. 2. As it canbe easily shown, for this diagram the renormalization ofthe heat and electric current vertices is exactly canceledby the renormalization of the two electron propagators.Therefore, the renormalization of the chemical potentialdoes not affect our results in the leading order.Finally, we remark that the bare figure of merit Z (0) T for granular metals at g T > − . Therefore these materials are not suitablefor solid-state refrigerators , but should be replaced bygranular semiconductors with g T <
1. However, the caseof granular metals is still relevant for low temperatureapplications in, e.g., thermocouples. Therefore we con-clude this paragraph discussing the dimensionless figureof merit ZT of granular materials at weak coupling be-tween the grains, g T ≪
1. In this regime the electronic
FIG. 6: (color online) Lowest order diagram for the heat-electric current correlator for homogeneous disordered metals.The external Bosonic frequency is denoted by Ω (wavy lines)and the internal Fermionic frequency by ω (straight lines).The electric and heat current vertexes are −→ e and −→ h , respec-tively. contribution to thermal conductivity κ e of granular met-als was recently investigated in Ref. 24, where it wasshown that κ e ∼ g T T /E c . (41)In this regime the electric conductivity of granular metalsobeys the law σ ∼ g T exp( − E c /T ) . (42)However, an expression for the thermoelectric coefficientin this region is not available yet, but recently it hasbeen proposed, based on experiment, that nanostruc-tured thermoelectric materials in the low coupling region(AgPb m SbT e m , Bi Te /Sb Te , or CoSb ) can have higher figures of merit than their bulk counter-parts.In conclusion, we have investigated the thermoelectriccoefficient and thermopower of granular nanomaterialsin the limit of large tunneling conductance between thegrains and temperatures T > g T δ . We have shown towhat extend quantum and confinement effects in granu-lar metals are important in changing ZT depending onsystem parameters. We also presented the details of ourcalculations. Acknowledgements
We thank Frank Hekking, NickKioussis, and Gang Lu for useful discussions. A. G. wassupported by the U.S. Department of Energy Office ofScience under the Contract No. DE-AC02-06CH11357.
APPENDIX A: THERMOELECTRICCOEFFICIENT OF HOMOGENEOUSDISORDERED METALS IN THE ABSENCE OFINTERACTION
In order to demonstrate important steps of our cal-culations, we present a derivation of the thermoelectriccoefficient for homogeneous disordered metals in the ab-sence of interaction in this appendix. In linear responsetheory the thermoelectric coefficient can be written as η (0) = ı ∂L d T ∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 A (0) , (A1) where the diagrammatic representation of the correlator A (0) is shown in Fig. 6 with electric and heat currentvertexes −→ e = em X p −→ p b a † p b a p , −→ h = 1 m X p −→ p ξ p b a † p b a p , (A2)respectively. The sums over momentum aretransformed into integrals P p → (cid:0) L π (cid:1) d R d d p bywhich we can transform the sum over the mo-mentum product to R d d p e R d d p h −→ p e · −→ p h = R d d p e R d d p h |−→ p e | h−→ n e · −→ n h i δ ( −→ p e − −→ p h ) (cid:0) πL (cid:1) d ,where −→ n α is the unit vector in direction of the current α ∈ { e, h } . Averaging over angles gives h−→ n e · −→ n h i = 1 /d ,and therefore: A (0) ( k, Ω m ) = − sd em T X ω n L d Z d d p (2 π ) d |−→ p | × ξ p G ( p + k, ω n + Ω m ) G ( p, ω n ) , (A3)where s is the spin degeneracy factor, Ω m an externalbosonic Matsubara frequency, and G ( p, ω n ) the momen-tum dependent Matsubara Green’s function. In the fol-lowing we only consider the case of zero external momen-tum, k = 0. For the advanced (in C + ) and retarded (in C − ) Green’s function, we use G ± ( p, ω n ) = [ ıω n − ξ p ± ı/ (2 τ ω )] − , (A4)where ξ p = ε p − ε F is the electron energy with respectto the Fermi energy [ ε p = p / (2 m )] and the (energy)relaxation time τ ω depends on the (real) frequency ω .The momentum integral in Eq. (A3) is transformedinto an energy integral as follows Z d d p (2 π ) d f ( |−→ p | ) = Ω d d − m d (2 π ) d Z ∞ f ( √ mε ) ε d − dε, (A5)where Ω d is the value of the angular integral (Ω , , = { , π, π } ). Since the Green’s function depends on ξ ,we need to rewrite this integral using the dε = dξ and ε d/ − = ε d/ − F (1 + ξ/ε F ) d/ − Z ∞ ε d/ − dε = ε d/ − F Z ∞− ε F (1 + ξ/ε F ) d/ − dξ (A6) ≈ ε d/ − F Z ∞−∞ [1 + ( d/ − ξ/ε F ] dξ . Combining Eqs. (A5) and (A6) we obtain Z d d p (2 π ) d f ( |−→ p | ) ≈ ν (0) d Z ∞−∞ (cid:20) d − ξε F (cid:21) (A7) × f (cid:16)p m ( ξ − ε F ) (cid:17) dξ, with ν (0) d = d/ − Ω( d )(2 π ) d m d/ ε d/ − F being the density ofstates (DOS) at the Fermi surface.0The DOS (and τ ω ) depends on ω via ν d ( ω ) = 1 π Z d d p (2 π ) d ℑ G − (A8)= ν (0) d π (cid:20) ωε F (cid:21) d − Z (1 + x/ e x ) d/ − x + 1 dx , with e x = 2 τ ω ε F (1 + ω/ε F ) and the symbol ℑ stands forimaginary part. However, the term x/ e x is neglected inthe integral, hence [ R R (cid:0) x + 1 (cid:1) − dx = π ] τ − ω = τ − ν d ( ω ) ν (0) d ≈ τ − [1 + ( d/ − ω/ε F ] . (A9)For conveniens we drop the momentum argument of G in the following, as well as the − ıω argument.Using all the above equations we finally obtain η (0) = − sd ıem ν (0) d ∂∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 X ω n Z dξ (A10) × [1 + ( d/ − ξ/ε F ] ξ ( ξ − ε F ) × ı ( ω n + Ω m ) − ξ + ı sgn( ω n + Ω m ) / (2 τ ) × ıω n − ξ + ı sgn( ω n ) / (2 τ ) . Here the derivative by the real external frequency, Ω,requires the analytic continuation of the Matsubara ex-pression. For the following calculation of the ω n sum andenergy integral, it is convenient to introduce the notation g ξ ≡ [1 + ( d/ − ξ/ε F ] ξ ( ξ − ε F ) (A11)= − ξε F + (2 − d/ ξ + O (cid:2) ξ /ε F (cid:3) . Next, we split the sum over the Fermionic Matsubarafrequencies ω n = 2 πT ( n + 1 /
2) into three intervals I =] −∞ ; − Ω m ] , I = ] − Ω m ; 0[ , I = [0; ∞ [ such that we canspecify the analytical structure of the Green’s function,i.e. η (0) = − sd ıem ν (0) d ∂∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 Z dξ g ξ × " X n ∈ I G − ( ω n + Ω m ) G − ( ω n )+ X n ∈ I G + ( ω n + Ω m ) G − ( ω n )+ X n ∈ I G + ( ω n + Ω m ) G + ( ω n ) . Now, we perform the analytical continuation of the sumsto real frequencies and at the same time of the externalfrequency: S = X n ∈ I G − ( ω n + Ω m ) G − ( ω n ) (A12)= Z R − dω πıT tanh( ω T ) G − ( − ıω ) G − ( − ıω + ı Ω) ,S = X n ∈ I G + ( ω n + Ω m ) G − ( ω n ) = Z R − dω πıT tanh( ω T )(A13) × (cid:2) G + ( − ıω − ı Ω) G − ( − ıω ) − G + ( − ıω ) G − ( − ıω + ı Ω) (cid:3) ,S = X n ∈ I G + ( ω n + Ω m ) G + ( ω n ) (A14)= Z R dω πıT tanh( ω T ) G + ( − ıω − ı Ω) G + ( − ıω ) . Next, we rewrite the Ω-derivative of the integrant interms of ω -derivatives, (arguments − ıω are omitted): ∂∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 G ( − ıω ∓ ı Ω) = ± ∂∂ω G , (A15)which is valid for both, advanced and retarded functions,and therefore we can simplify1 ∂∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 (cid:2) − G + ( − ıω − ı Ω) G − + G + G − ( − ıω + ı Ω) (cid:3) = (cid:0) G + (cid:1) G − + G + (cid:0) G − (cid:1) = − ∂∂ω (cid:0) G + G − (cid:1) , (A16) ∂∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 (cid:2) − G − G − ( − ıω + ı Ω) (cid:3) = − G − (cid:0) G − (cid:1) = 12 ∂∂ω (cid:0) G − (cid:1) ,∂∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 (cid:2) G + ( − ıω − ı Ω) G + (cid:3) = − (cid:0) G + (cid:1) G + = 12 ∂∂ω (cid:0) G + (cid:1) . Using Eqs. (A15) and (A16) we can write the following ∂∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 [ S + S + S ] = Z R dω πıT tanh ω T ∂∂ω (cid:8) ( G + − G − ) (cid:9) . (A17)We notice that at this point the ω -integration should notbe done by parts, since the boundary term is important.We can now do the ξ -integration using the fact that G + − G − = − ( ı/τ ) G + G − Ξ( ω ) ≡ Z dξ g ξ (cid:18) ıG + G − τ (cid:19) (A18)= − Z g ξ dξ ( τ [ ω − ξ + ı/ (2 τ )] [ ω − ξ − ı/ (2 τ )]) = − π (cid:20) τ ω (4 − d ) + 4 − d τ − ωτ ε F (cid:21) . Here, the higher order terms in g ξ are neglected. InEq. (A18) we need to keep the terms proportional to ω only, therefore after ω -expansion of τ = τ ω we obtain:Ξ( ω ) = − π (cid:2) τ ω (4 − d ) + ω τ ( d − (cid:3) = − πτ ω . (A19)And, finally η (0) = − sd ıeν (0) d m Z R dω πıT tanh( ω T ) ∂∂ω Ξ( ω )= − s d em ν (0) d τ T Z R ω dω cosh ( ω/ T )= − π s d em ν (0) d ( τ T ) . (A20)To perform the last integration we used the integral R x dx cosh ( x ) = π /
6. In d = 3 the density of states hasthe form ν (0) d =3 = mp F π , leading to η (0)3 D = − s ep F ( τ T ) . (A21) APPENDIX B: HEAT AND ELECTRICCURRENT OPERATORS OF GRANULARMETALS
In this Appendix we derive an expression for the heatand electric current operators of granular metals in the presence of Coulomb interaction. The Hamiltonian forthe granular system is b H = P i b ǫ i where b ǫ i = X k ξ k b a † i,k b a i,k + e X j b n i C − ij b n j (B1)+ 12 X j,k,q h t kqij b a † i,k b a j,q + t qkji b a † j,q b a i,k i ≡ b ǫ ( e ) i + b ǫ ( c ) i + b ǫ ( t ) i . Here we introduce the notation b n i ≡ P k b a † i,k b a i,k for thenumber of electron within a grain i . The creation (annihi-lation) operators b a † i,k ( b a i,k ) satisfy the anti-commutationrelations nb a † i,k ; b a j,q o + = δ ij δ kq and nb a ( † ) i,k ; b a ( † ) j,q o + = 0.The electric current through grain i is − ı d b n i dt = hb n i ; b H i ≡ − ıe X j b ( e ) ij . (B2)And the heat current through grain i is − ı d b ǫ i dt = hb ǫ i ; b H i ≡ − ı X j b ( h ) ij . (B3)Accordingly the total current operators are given by −→ e,h = −→ n (0) e,h X i,j b ( e,h ) ij , (B4)where −→ n (0) e,h are the unit directions of the current (electricfield, E , gradient) or heat (temperature, T , gradient)flows.
1. electric current operator
First, we calculate the electric current. Since the num-ber operator of electrons b n i commutes with the first twoterms of b H in Eq. (B2), i.e. hb n i ; b H i = hb n i ; b H ( t ) i , we onlyneed to calculate2 hb n i ; b H ( t ) i = 12 X k,i ′ ,j,k ′ ,q hb a † i,k b a i,k (cid:16) t k ′ qi ′ j b a † i ′ ,k ′ b a j,q + t qk ′ ji ′ b a † j,q b a i ′ ,k ′ (cid:17) − (cid:16) t k ′ qi ′ j b a † i ′ ,k ′ b a j,q + t qk ′ ji ′ b a † j,q b a i ′ ,k ′ (cid:17) b a † i,k b a i,k i = X j,k,q (cid:16) t kqij b a † i,k b a j,q − t qkji b a † j,q b a i,k (cid:17) ≡ b C ( nt ) i (B5)As a result we get the following expression for the electriccurrent operator b ( e ) ij = ıe X k,q (cid:16) t kqij b a † i,k b a j,q − t qkji b a † j,q b a i,k (cid:17) . (B6)
2. heat current operator
Second we turn to the heat current, for which we needto calculate nine commutators (see Eq. (B3)) b C ( αβ ) i ≡ hb ǫ ( α ) i ; b H ( β ) i , (B7)where α ∈ { e, c, t } and operator b ǫ ( α ) i was defined inEq. (B1).Since operators b n i commute with each other, four ofthe commutators in Eq. (B7) vanish: b C ( ee ) i = b C ( cc ) i = b C ( ec ) i = b C ( ce ) i = 0 . (B8) The heat current operator can be conveniently writtenas a sum of two contributions b ( h ) ij = b ( h, ij + b ( h, ij , (B9)where, the non-interacting part b ( h, ij of the heat currentoriginates from the sum of the commutators b C ( et ) i + b C ( te ) i and the interacting part b ( h, ij from the sum of b C ( ct ) i + b C ( tc ) i . a. non-interacting part of heat current operator We first calculate the non-interacting part of the heatcurrent operator, b ( h, ij . That is, we consider the sum ofcommutators b C ( et ) i + b C ( te ) i . A straightforward calculationleads to b C ( et ) i = 12 X i ′ ,j,k ′ ,q,k ξ k (cid:16)b a † i,k b a i,k h t k ′ qi ′ j b a † i ′ ,k ′ b a j,q + t qk ′ ji ′ b a † j,q b a i ′ ,k ′ i − h t k ′ qi ′ j b a † i ′ ,k ′ b a j,q + t qk ′ ji ′ b a † j,q b a i ′ ,k ′ i b a † i,k b a i,k (cid:17) = X j,k,q ξ k (cid:16) t kqij b a † i,k b a j,q − t qkji b a † j,q b a i,k (cid:17) , (B10) b C ( te ) i = 12 X j,k,q,i ′ ,k ′ ξ k ′ (cid:16)h t kqij b a † i,k b a j,q + t qkji b a † j,q b a i,k i b a † i ′ ,k ′ b a i ′ ,k ′ − b a † i ′ ,k ′ b a i ′ ,k ′ h t kqij b a † i,k b a j,q + t qkji b a † j,q b a i,k i(cid:17) = 12 X j,k,q (cid:16) ξ q t kqij b a † i,k b a j,q + ξ k t qkji b a † j,q b a i,k − ξ k t kqij b a † i,k b a j,q − ξ q t qkji b a † j,q b a i,k (cid:17) . (B11)Using Eqs. (B10) and (B11) we obtain for the sum b C ( et ) i + b C ( te ) i = 12 X j,k,q (cid:16) ( ξ q + ξ k ) t kqij b a † i,k b a j,q − ( ξ k + ξ q ) t qkji b a † j,q b a i,k (cid:17) . (B12) As a result the non-interacting part of the heat currentoperator of granular metals has the form b ( h, ij = ı X k,q ξ k + ξ q h t kqij b a † i,k b a j,q − t qkji b a † j,q b a i,k i . (B13)3 b. interacting part of heat current operator To obtain the expression for the interacting part ofthe heat current operator b ( h, ij we need to considerthe commutators of the form hb ǫ ( t ) i ′ j ; b n i i , where b ǫ ( t ) ij = P k,q h t kqij b a † i,k b a j,q + t qkji b a † j,q b a i,k i . Using the calculation for b C ( nt ) i [Eq. (B5)] we get2 hb ǫ ( t ) i ′ j ; b n i i = X k,k ′ ,q (cid:16) t k ′ qi ′ j (cid:16) δ ij δ kq b a † i ′ ,k ′ b a i,k − δ ii ′ δ kk ′ b a † i,k b a j,q (cid:17) + t qk ′ ji ′ (cid:16) δ ii ′ δ kk ′ b a † j,q b a i,k − δ ij δ kq b a † i,k b a i ′ ,k ′ (cid:17)(cid:17) (B14)= [ δ ij − δ ii ′ ] X k,q (cid:16) t kqi ′ j b a † i ′ ,k b a j,q − t kqji ′ b a † j,k b a i ′ ,q (cid:17) = − ıe [ δ ij − δ ii ′ ] b ( e ) i ′ j . For the following steps of the calculation of b C ( ct ) i and b C ( tc ) i , we need the commutator hb n m ; b ( e ) ij i − ıe hb n m ; b ( e ) ij i = X k,q hb n m (cid:16) t kqij b a † i,k b a j,q − t qkji b a † j,q b a i,k (cid:17) − (cid:16) t kqij b a † i,k b a j,q − t qkji b a † j,q b a i,k (cid:17) b n m i = − ( δ mj − δ mi ) X k,q (cid:16) t kqij b a † i,k b a j,q + t qkji b a † j,q b a i,k (cid:17) = 2 ( δ mi − δ mj ) b ǫ ( t ) ij . (B15)To calculate the interacting part of the heat current operator b ( h, ij we need the following commutators4 e b C ( ct ) i = 2 X j,i ′ ,j ′ (cid:16)b n i C − ij b n j b ǫ ( t ) i ′ j ′ − b ǫ ( t ) i ′ j ′ b n i C − ij b n j (cid:17) = ıe X j,i ′ ,j ′ C − ij (cid:16) [ δ ij ′ − δ ii ′ ] b ( e ) i ′ j ′ b n j + [ δ jj ′ − δ ji ′ ] b n i b ( e ) i ′ j ′ (cid:17) (B16)= ıe X j C − ij X i ′ b ( e ) i ′ i b n j − X j ′ b ( e ) ij ′ b n j + X i ′ b n i b ( e ) i ′ j − X j ′ b n i b ( e ) jj ′ = ıe X j,m C − ij h(cid:16)b ( e ) mi − b ( e ) im (cid:17) b n j + b n i (cid:16)b ( e ) mj − b ( e ) jm (cid:17)i , e b C ( tc ) i = 2 X j,i ′ ,j ′ (cid:16)b ǫ ( t ) ij b n i ′ C − i ′ j ′ b n j ′ − b n i ′ C − i ′ j ′ b n j ′ b ǫ ( t ) ij (cid:17) = − ıe X j,i ′ ,j ′ C − i ′ j ′ (cid:16) [ δ j ′ j − δ j ′ i ] b n i ′ b ( e ) ij + [ δ i ′ j − δ i ′ i ] b ( e ) ij b n j ′ (cid:17) (B17)= − ıe X j X i ′ C − i ′ j b n i ′ b ( e ) ij + X j ′ C − jj ′ b ( e ) ij b n j ′ − X i ′ C − i ′ i b n i ′ b ( e ) ij − X j ′ C − ij ′ b ( e ) ij b n j ′ = − ıe X j,m h(cid:0) C − mj − C − mi (cid:1) b n m b ( e ) ij + (cid:0) C − jm − C − im (cid:1) b ( e ) ij b n m i . b ( e ) mi = − b ( e ) im and C − mj = C − jm , and hb n m ; b ( e ) ij i = 2 ıe ( δ mi − δ mj ) b ǫ ( t ) ij [Eq. B15] we obtain b C ( ct ) i + b C ( tc ) i = ıe X j,m (cid:16) C − ij h(cid:16)b ( e ) mi − b ( e ) im (cid:17) b n j + b n i (cid:16)b ( e ) mj − b ( e ) jm (cid:17)i − h(cid:0) C − mj − C − mi (cid:1) b n m b ( e ) ij + (cid:0) C − jm − C − im (cid:1) b ( e ) ij b n m i(cid:17) = ıe X j,m (cid:18) C − ij hb ( e ) mi b n j + b n i b ( e ) mj i − (cid:20) C − mj nb n m ; b ( e ) ij o + − C − mi nb ( e ) ij ; b n m o + (cid:21)(cid:19) = ıe X j,m C − ij hb ( e ) mi b n j + b n j b ( e ) mi − ıe ( δ jm − δ ji ) b ǫ ( t ) mi + b n i b ( e ) mj + b ( e ) mj b n i + 2 ıe ( δ im − δ ij ) b ǫ ( t ) mj i − (cid:20) C − mj nb n m ; b ( e ) ij o + − C − mi nb ( e ) ij ; b n m o + (cid:21) = − e X j (cid:16) − b ǫ ( t ) ji + b ǫ ( t ) ji + b ǫ ( t ) ij − b ǫ ( t ) ji (cid:17)| {z } =0 + ıe X j,m C − ij (cid:20)nb ( e ) mi ; b n j o + + nb n i ; b ( e ) mj o + (cid:21) − (cid:20) C − mj nb n m ; b ( e ) ij o + − C − ji nb ( e ) im ; b n j o + (cid:21) = ıe X j,m (cid:18) C − ij nb n i ; b ( e ) mj o + − C − mj nb n m ; b ( e ) ij o + (cid:19) . (B18)As a result for the interaction part of the heat current operator of granular metals we obtain b ( h, ij = − e X m (cid:18) C − im nb n i ; b ( e ) jm o + − C − jm nb n j ; b ( e ) im o + (cid:19) . (B19)So far we have omitted the last commutator b C ( tt ) i in Eq. (B7), which would be an additional contribution to thenon-interacting part of the heat current operator, b ( h, ij in Eq. (B13). However, if this term is summed over i itvanishes and does not contribute to b ( h, ij . However, for completeness we present the calculation of b C ( tt ) i here as well.We define the commutator b C ( tt ) i as b C ( tt ) i = X j,i ′ ,j ′ (cid:16)b ǫ ( t ) ij b ǫ ( t ) i ′ j ′ − b ǫ ( t ) i ′ j ′ b ǫ ( t ) ij (cid:17) ≡ b A − b B, (B20)and simplify first the expression for operator 4 b B .4 b ǫ ( t ) i ′ j ′ b ǫ ( t ) ij = X k,q,k ′ ,q ′ (cid:16) t k ′ q ′ i ′ j ′ b a † i ′ ,k ′ b a j ′ ,q ′ + t q ′ k ′ j ′ i ′ b a † j ′ ,q ′ b a i ′ ,k ′ (cid:17) (cid:16) t kqij b a † i,k b a j,q + t qkji b a † j,q b a i,k (cid:17) (B21)= X k,q,k ′ ,q ′ (cid:16) t k ′ q ′ i ′ j ′ t kqij h δ j ′ i δ q ′ k b a † i ′ ,k ′ b a j,q − δ i ′ j δ k ′ q b a † i,k b a j ′ ,q ′ i + t q ′ k ′ j ′ i ′ t kqij h δ i ′ i δ k ′ k b a † j ′ ,q ′ b a j,q − δ j ′ j δ q ′ q b a † i,k b a i ′ ,k ′ i ++ t k ′ q ′ i ′ j ′ t qkji h δ j ′ j δ q ′ q b a † i ′ ,k ′ b a i,k − δ i ′ i δ k ′ k b a † j,q b a j ′ ,q ′ i + t q ′ k ′ j ′ i ′ t qkji h δ i ′ j δ k ′ q b a † j ′ ,q ′ b a i,k − δ j ′ i δ q ′ k b a † j,q b a i ′ ,k ′ i(cid:17) + 4 b A . b C ( tt ) i = − X j,k,q X i ′ ,k ′ (cid:16) t k ′ ki ′ i t kqij b a † i ′ ,k ′ b a j,q − t qk ′ ji ′ t kqij b a † i,k b a i ′ ,k ′ (cid:17) − X j ′ ,q ′ (cid:16) t qq ′ jj ′ t kqij b a † i,k b a j ′ ,q ′ − t q ′ kj ′ i t kqij b a † j ′ ,q ′ b a j,q (cid:17) (B22)+ X i ′ ,k ′ (cid:16) t k ′ qi ′ j t qkji b a † i ′ ,k ′ b a i,k − t kk ′ ii ′ t qkji b a † j,q b a i ′ ,k ′ (cid:17) − X j ′ ,q ′ (cid:16) t kq ′ ij ′ t qkji b a † j,q b a j ′ ,q ′ − t q ′ qj ′ j t qkji b a † j ′ ,q ′ b a i,k (cid:17) = − X j,m,k,q,p h t pkmi t kqij b a † m,p b a j,q + t pkmi t kqij b a † m,p b a j,q − t qpjm t kqij b a † i,k b a m,p − t qpjm t kqij b a † i,k b a m,p + t pqmj t qkji b a † m,p b a i,k + t pqmj t qkji b a † m,p b a i,k − t kpim t qkji b a † j,q b a m,p − t kpim t qkji b a † j,q b a m,p i = − X j,m,k,q,p h t pkmi t kqij b a † m,p b a j,q − t qpjm t kqij b a † i,k b a m,p + t pqmj t qkji b a † m,p b a i,k − t kpim t qkji b a † j,q b a m,p i = 12 X j,m,k,q,p h t qpjm t kqij b a † i,k b a m,p − t pqmj t qkji b a † m,p b a i,k i = X j,m,k,q,p t qpjm t kqij hb a † i,k b a m,p − b a † m,p b a i,k i . The underlined terms cancel since t pkmi = t kpim , which isalso used in the last step. If this expression is summedover i , we can exchange indices m and i (and k and q ) inthe second summand. As a result we obtain X i b C ( tt ) i = 0 , (B23)i.e. there is no additional contribution to the full heatcurrent operator −→ h in Eq. (B4) from this term.
3. summary: electric and heat current operators
To summarize this appendix we explicitly write the ex-pressions for electric b ( e ) ij and heat ˆ ( h ) ij current operators b ( e ) ij = ıe X k,q (cid:16) t kqij b a † i,k b a j,q − t qkji b a † j,q b a i,k (cid:17) . (B24)ˆ ( h ) ij = ˆ ( h, ij + ˆ ( h, ij , (B25)ˆ ( h, ij = ı X k,q ξ k + ξ q h t kqij ˆ a † ik ˆ a jq − t qkji ˆ a † jq ˆ a ik i , ˆ ( h, ij = − e X m " { ˆ n i ; ˆ ( e ) jm } + C im − { ˆ n j ; ˆ ( e ) im } + C jm . APPENDIX C: THERMOELECTRICCOEFFICIENT OF GRANULAR METALS INTHE ABSENCE OF INTERACTION
In this appendix we consider the thermoelectric coeffi-cient η (0) of granular metals in the absence of interaction FIG. 7: (color online) Lowest order diagram for the heat-electric current correlator of granular metals. The externalBosonic frequency is denoted by Ω (wavy lines) and the inter-nal Fermionic frequency by ω . The electric and heat currentvertexes are −→ e and −→ h , respectively. The ± denote the pos-sible analytical structure of the Green’s functions (straightlines). in analogy to Appendix A. The expression for the ther-moelectric coefficient in linear response theory is η (0) = ı ∂a d T ∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 Q (0) . (C1)Here a is the grain size and Q (0) the correlator of the heatcurrent, −→ (0) h (see also Fig. 1a), and electric current, −→ e ,shown in Fig. 7.For granular metals the important element of the di-agram is the tunneling matrix elements t kqij describingthe coupling between grains i and j . Therefore we de-rive an expression for t kqij in the following, assuming thati and j are nearest neighbor grains and t kqij is indepen-dent of the position in the sample. In order to calculatethe energy dependence of these elements we assume thetunneling barrier between grains to be a delta potential.For the one-particle Hamiltonian b H = − ~ m d dx + λδ ( x )the transmission rate for a single particle with energy6 ε p = ε F + ξ p is T p = (cid:18) mλ ~ ε p (cid:19) − ≃ (cid:18) mλ ~ ε F (1 − ξ p /ε F ) (cid:19) − (C2)= T (cid:18) ~ ε F mλ − ξ p /ε F (cid:19) − . Here T (cid:16) mλ ~ ε F (cid:17) − is the transmission rate at ε p = ε F and we use the fact that ξ p ≪ ε F . Next, we con-sider the case of large barriers, in this regime T p ≃ T (1 + ξ p /ε F ). In granular systems we have many chan-nels and have to consider tunneling processes with en-ergy ξ in grain i = 1 and with ξ in grain j = 2: t ∝ N (cid:0) T p + T p (cid:1) . So, we obtain t ( ξ , ξ ) ≃ t (cid:18) ξ + ξ ε F (cid:19) . (C3)Therefore we have the following expression for correlationfunction in Eq. (C1) Q (0) = − set T D −→ n (0) e · −→ n (0) h E X ω n a d +2 Z d d p (2 π ) d Z d d p (2 π ) d (cid:18) ξ + ξ (cid:19) (cid:18) ξ + ξ ε F (cid:19) G ( p , ω n ) G ( p , ω n + Ω m ) ≈ − s d et T a d +2 (cid:16) ν (0) d (cid:17) X ω n Z dξ dξ g ( ξ , ξ ) G ( ξ , ω n ) G ( ξ , ω n + Ω m ) , (C4)where −→ n (0) α is the unit vector in direction of the current α ∈ { e, h } , h−→ n (0) e · −→ n (0) h i = 1 /d is the result of averagingover angles, the summation goes over Fermionic Matsub-ara frequencies ω n = 2 πT ( n + 1 / G ( p, ω n ) is theGreen’s function defined in Eq. (A4) of Appendix A withmomenta/energies p i / ξ i of grain i . To shorter the nota-tion in the following we neglect the momentum argumentand attach the grain index to G . In Eq. (C4) we intro-duce the notation g ( ξ , ξ ) = ( ξ + ξ ) (cid:18) ξ + ξ ε F (cid:19) (C5) × (cid:18) (cid:20) d − (cid:21) ξ ε F (cid:19) (cid:18) (cid:20) d − (cid:21) ξ ε F (cid:19) . The factors in this order arise from: the heat currentvertex Eq. (12c), the energy correction to tunneling ele-ments Eq. (16), and corrections to the DOS due to finiteFermi energy Eq. (A7). In the linear order in ξ/ε F weobtain g ( ξ , ξ ) = ξ + ξ + d ξ + ξ ) ε F . (C6)We first perform the analytical continuation inEq. (C4) (for convenience the grain index is written asindex to the Green’s functions) Z dξ dξ g ( ξ , ξ ) X n ∈ I G − ( ω n ) G − ( ω n + Ω m ) | {z } S + X n ∈ I G +1 ( ω n ) G − ( ω n + Ω m ) | {z } S + X n ∈ I G +1 ( ω n ) G +2 ( ω n + Ω m ) | {z } S . (C7)After analytical continuation we obtain S = − Z dω πıT tanh (cid:16) ω T (cid:17) G − ( − ıω + ı Ω) G − ( − ıω ) , (C8) S = − Z dω πıT tanh (cid:16) ω T (cid:17) (cid:2) G − ( − ıω ) G +2 ( − ıω − ı Ω) − G − ( − ıω + ı Ω) G +2 ( − ıω ) (cid:3) , (C9) S = − Z dω πıT tanh (cid:16) ω T (cid:17) G +1 ( − ıω ) G +2 ( − ıω − ı Ω) . (C10)Now we consider the derivative of the sum of ( S + S + S ) with respect to Bosonic frequency ∂∂ Ω (cid:12)(cid:12) Ω=0 . For brevity7we omit arguments − ıω of Green’s functions ∂∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 ( S + S + S ) = Z dω πıT tanh (cid:16) ω T (cid:17) (cid:20)(cid:18) ∂∂ω G − (cid:19) G − − G − (cid:18) ∂∂ω G +2 (cid:19) − (cid:18) ∂∂ω G − (cid:19) G +2 + G +1 (cid:18) ∂∂ω G +2 (cid:19)(cid:21) = Z dω πıT tanh (cid:16) ω T (cid:17) (cid:18) ∂∂ω G − (cid:19) (cid:0) G − − G +2 (cid:1) + (cid:18) ∂∂ω G +2 (cid:19) (cid:0) G +1 − G − (cid:1)| {z } exchange indices = Z dω πıT tanh (cid:16) ω T (cid:17) (cid:0) G − − G +2 (cid:1) ∂∂ω (cid:0) G − − G +1 (cid:1)| {z } [( G − − G +2 ) ∂∂ω ( G − − G +1 ) + ( G − − G +1 ) ∂∂ω ( G − − G +2 )]= Z dω πıT tanh (cid:16) ω T (cid:17) ∂∂ω (cid:2)(cid:0) G − − G +1 (cid:1) (cid:0) G − − G +2 (cid:1)(cid:3) = Z − dω πıT tanh (cid:16) ω T (cid:17) ∂∂ω (cid:20) τ ω G − G +1 G − G +2 (cid:21) . (C11)Now one can perform the integration over variables ξ and ξ using Eqs. (C7), (C11), and the residuum theorem Z dξ dξ g ( ξ , ξ ) G − G +1 G − G +2 (C12)= Z dξ dξ ξ + ξ + d ξ + ξ ) ε F ( ω − ξ − ı/ (2 τ )) ( ω − ξ + ı/ (2 τ )) | {z } ξ , = ω + ı/ (2 τ ) ( ω − ξ − ı/ (2 τ )) ( ω − ξ + ı/ (2 τ ))= 2 πı Z dξ ω + ı/ (2 τ ) + ξ + d ε F ( ω + ı/ (2 τ ) + ξ ) ( − ı/τ ) ( ω − ξ − ı/ (2 τ )) ( ω − ξ + ı/ (2 τ )) | {z } ξ , = ω + ı/ (2 τ ) = 4 π τ (cid:20) ω + ı/τ + d ε F (2 ω + ı/τ ) (cid:21) , where the ξ i, -values below the braces denote the poles in the complex plane used to perform the integration. As aresult we obtain the following expression for the derivative of correlation function ∂∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 Q (0) = s d et a d +2 (cid:16) ν (0) d (cid:17) Z dω πı tanh (cid:16) ω T (cid:17) ∂∂ω (cid:20) π (cid:20) ω + ı/τ + d ε F (2 ω + ı/τ ) (cid:21)(cid:21) (C13)= πs ıd et a d +2 (cid:16) ν (0) d (cid:17) Z dω tanh (cid:16) ω T (cid:17) ∂∂ω (cid:20) d ε F (2 ω + ı/τ [( d/ − ω/ε F ]) (cid:21) = − πs ıd et a d +2 (cid:16) ν (0) d (cid:17) (2 T ) − d ε F Z dω (2 + ( d/ − ı/ ( τ ε F )) ω cosh ( ω/ (2 T ))= − πs ıd et a d +2 (cid:16) ν (0) d (cid:17) (2 T ) d ε F π . In the second line in Eq. (C13) the derivative is takeninto account (removing the boundary terms of the par-tial integration) and the contributions of order 1 /ε F orsmaller are neglected in the last line.Finally, we obtain the following expression for non-interacting thermoelectric coefficient of granular metals η (0) = − sπ et a d +2 (cid:16) ν (0) d (cid:17) Tε F . (C14)One can re-write this expression using the relations: ν (0) d D d = ga − d ; ν (0) d = ( δa d ) − ; t = gδ / (2 π ), where D d is the diffusion constant, g the tunneling conductance, and δ the mean level spacing, giving η (0) = − sπ eg T a − d ( T /ε F ) . (C15) APPENDIX D: THERMOELECTRICCOEFFICIENT OF GRANULAR METALS WITHINTERACTION
In this appendix we consider the correction η (1) tothe thermoelectric coefficient of granular metals due toelectron-electron interaction, i.e., η = η (0) + η (1) , (D1)8 FIG. 8: (color online) Diagram describing the correction tothe thermoelectric coefficient due to electron-electron interac-tion corresponding to the term Q (1) in Eq. (D2). The externalBosonic frequency is denoted by Ω (wavy lines) and the inter-nal Fermionic frequency by ω (straight lines). The electric andheat current vertexes are −→ e and −→ h (without Coulomb contri-bution), respectively. The (red) triangles denote the diffusons D and the thick wavy line the screen Coulomb interaction. where η (0) was calculated in Appendix C, Eq. (C14). Thestructure of the diagrams Q (1) , Q (2) , Q (3) contributing to η (1) are shown in Fig. 2 and we can write η (1) = ı ∂a d T ∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 (cid:16) Q (1) + Q (2) + Q (3) (cid:17) . (D2)These diagrams include the effect of elastic scatteringof electron at impurities described by diffusons, D − = τ ω ( | Ω i | + ǫ q δ ), and the effect of the dynamically screenedCoulomb potential e V ( q, Ω i ) = D V ( q, Ω i ) D : e V ( q, Ω i ) = 2 E c ( q ) τ ω [ | Ω i | + 4 E c ( q ) ǫ q ] [ | Ω i | + ǫ q δ ] , (D3) V ( q, Ω i ) = (cid:18) E c ( q ) + 2 ǫ q | Ω i | + ǫ q δ (cid:19) − . The renormalized interaction vertices are: (i) inter-grainΦ (1) ω (Ω i ) = a d Z d −→ q (2 π ) d e V ( q, Ω i ) ′ X a cos( −→ q · −→ a ) (D4)= Z a d d −→ q (2 π ) d τ − ω E c ( q ) P ′ a cos( −→ q · −→ a )[ | Ω i | + 4 E c ( q ) ǫ q ] [ | Ω i | + ǫ q δ ] , and (ii) intra-grain:Φ (2) ω (Ω i ) = a d d Z d −→ q (2 π ) d e V ( q, Ω i ) (D5)= Z d d d −→ q (2 π ) d d τ − ω E c ( q )[ | Ω i | + 4 E c ( q ) ǫ q ] [ | Ω i | + ǫ q δ ] , with ǫ q = 2 g T (cid:0) d − P ′ a cos( −→ q · −→ a ) (cid:1) where P ′ a standsfor summation over all directions and orientations n ± a −→ e (0) j o , and E c ( q ) = e C ( q ) = e a d − ln( qa ) , d = 1 π/q, d = 22 π/q , d = 3 . (D6)Explicitly, the contribution Q (1) in Eq. (D2) is given by FIG. 9: (color online) Diagram describing correction to thethermoelectric coefficient due to electron-electron interactioncorresponding to term Q (2) in Eq. (D2). All notations are thesame as in Fig. 8 Q (1) = − s d et T a d +2 (cid:16) ν (0) d (cid:17) (D7) × X ω n , Ω i Z dξ dξ g F ( s s s s )1 Φ (1) ω (Ω i )= Q (1 , ±± ) + Q (1 , ∓∓ ) , with g ≡ g ( ξ , ξ ) defined in Eq. (C6) of Appendix C.Here we introduce the function F ( s s s s )1 = G s ( ω n + Ω m ) G s ( ω n + Ω m + Ω i ) × G s ( ω n ) G s ( ω n + Ω i ) , (D8)with s i denote the analytic structure of the Green’s func-tions implying restrictions on the frequency summation– in principle there are 16 different combinations of the s i ’s, see Fig. 8. However, only the two diagrams Q (1 , ±± ) and Q (1 , ∓∓ ) contribute to the correction η (1) . The otheranalytical structures do not have either valid frequencydomains or poles are located in only one half plane of C .For the contribution Q (2) in Eq. (D2), Fig. 9, we havethe following expression Q (2) = − s d et T a d +2 (cid:16) ν (0) d (cid:17) X ω n , Ω i Z dξ dξ g (D9) × F ( s s s )2 Φ (2) ω (Ω i ) = 2 h Q (2 , ± +) + Q (2 , ∓− ) + Q (2 , ±− ) i , where F ( s s s )2 = [ G s ( ω n + Ω m )] G s ( ω n + Ω m + Ω i ) G s ( ω n ) . (D10)Due to symmetry all three contributing diagrams in theright hand side of Eq. (D9) have a factor of 2. Again, outof the eight possible combinations for s , s , and s onlythree combinations – Q (2 , ± +) , Q (2 , ±− ) , Q (2 , ∓− ) – have avalid or non-zero analytical structure.The diagram Q (3) , shown in Fig. 10, describes the con-tribution of the correlation function with the interactionpart of the heat current operator, b ( h, ij (see Appendix B),and has therefore a different structure in comparison withcontributions Q (1) and Q (2) : Q (3) = − s d et T a d +2 (cid:16) ν (0) d (cid:17) (D11) × X ω n , Ω i Z dξ dξ g F ( s s s )3 Φ (Ω i , q ) , FIG. 10: (color online) Diagram describing correction to thethermoelectric coefficient due to electron-electron interactioncorresponding to term Q (3) in Eq. (D2). All notations are thesame as in Fig. 8, but here the heat vertex corresponds to b h, i and the diagram has only one diffuson. with F ( s s s )3 = G s ( ω n + Ω m + Ω i ) G s ( ω n + Ω m ) G s ( ω n ) , (D12)and g ( ξ , ξ ) = 2 (cid:18) ξ + ξ ε F (cid:19) (D13) × (cid:18) (cid:18) d − (cid:19) ξ ε F (cid:19) (cid:18) (cid:18) d − (cid:19) ξ ε F (cid:19) = 2 (cid:18) d ε F ( ξ + ξ ) (cid:19) + O ( ξ /ε F ) . Since the linear part of function g ( ξ , ξ ) has a factor ε − F , the main contribution to η (1) from the diagram Q (3) is of the order of T /ε F , whereas Q (1) and Q (2) have1 /ε F contributions, which we are considering here only.Therefore we will not consider diagram Q (3) any further.In the following we discuss the five diagrams contribut-ing to η (1) in details, especially their analytical structureand the resulting restrictions on the frequency summa-tions: Q (1 , ±± ) , Q (1 , ∓∓ ) , Q (2 , ± +) , Q (2 , ±− ) , and Q (2 , ∓− ) .
1. Calculation of contribution Q (1 , ±± ) in Eq. (D7). Here we discuss the contribution Q (1 , ±± ) introducedin Eq. (D7). The analytical structure of this diagram(Fig. 8), defined by indexes s to s in Eq. (D7), demands ω n + Ω m > , ω n + Ω m + Ω i < ω n > , ω n + Ω i < . These inequalities define the limits of frequency summa-tions: 0 < ω n < − Ω m − Ω i and Ω i < − Ω m . (D15)First, we calculate the ξ -integrals in Eq. (D7) using theresidue theorem: Z dξ dξ ξ + ξ + d ξ + ξ ) ε F ( ı ( ω n + Ω m ) − ξ + ı/ (2 τ )) ( ı ( ω n + Ω m + Ω i ) − ξ − ı/ (2 τ )) | {z } ξ , = ı ( ω n +Ω m +Ω i ) − ı/ (2 τ ) ∈ C − ( ıω n − ξ + ı/ (2 τ )) | {z } ξ , = ıω n + ı/ (2 τ ) ∈ C + ( ı ( ω n + Ω i ) − ξ − ı/ (2 τ ))= − (2 πı ) ı (2 ω n + Ω m + Ω i ) − d ε F (2 ω n + Ω m + Ω i ) ( − ı Ω i + ı/τ ) ( ı Ω i − ı/τ ) = 4 π (Ω i − /τ ) (cid:20) ı (2 ω n + Ω m + Ω i ) − d ε F (2 ω n + Ω m + Ω i ) (cid:21) = 4 π (Ω i − /τ ) g ( ı (2 ω n + Ω m + Ω i )) ≈ π τ g ( ı (2 ω n + Ω m + Ω i )) , (D16)where in the last line we introduce the function g ( z ) = z + [ d/ (2 ε F )] z . (D17)For the final approximation we used the fact that | Ω i | τ ≪
1. The poles for the residual are written be-low the underbraces.Using the result of integration over ξ and ξ inEq. (D16), we can simplify the expression for this dia-gram to Q (1 , ±± ) = − π sd et T τ a d +2 (cid:16) ν (0) d (cid:17) (D18) × X <ω n < − Ω m − Ω i Ω i< − Ω m Φ (1) ω (Ω i ) g ( ı (2 ω n + Ω m + Ω i )) .
2. Calculation of contribution Q (1 , ∓∓ ) in Eq. (D7) Here we discuss the contribution Q (1 , ∓∓ ) defind inEq. (D7). The analytical structure of this diagram de-mands ω n + Ω m < , ω n + Ω m + Ω i > ω n < , ω n + Ω i > , which defines the limits of the frequency summations: − Ω i < ω n < − Ω m and Ω i > Ω m . (D20)We first perform the ξ -integrals in Eq. (D7):0 Z dξ dξ (cid:16) ξ + ξ + d ξ + ξ ) ε F (cid:17) [ ı ( ω n + Ω m ) − ξ − ı/ (2 τ )] [ ı ( ω n + Ω m + Ω i ) − ξ + ı/ (2 τ )] | {z } ξ , = ı ( ω n +Ω m +Ω i )+ ı/ (2 τ ) ∈ C + [ ıω n − ξ − ı/ (2 τ )] | {z } ξ , = ıω n − ı/ (2 τ ) ∈ C − [ ı ( ω n + Ω i ) − ξ + ı/ (2 τ )]= 4 π ı (2 ω n + Ω m + Ω i ) − d ε F (2 ω n + Ω m + Ω i ) ( − ı Ω i − ı/τ ) ( ı Ω i + ı/τ ) = 4 π g ( ı (2 ω n + Ω m + Ω i ))(Ω i + 1 /τ ) ≈ π τ g ( ı (2 ω n + Ω m + Ω i )) . (D21)Substituting the result of Eq. (D21) back into Eq. (D7)we obtain Q (1 , ∓∓ ) = − π sd et T τ a d +2 (cid:16) ν (0) d (cid:17) (D22) × X − Ω i <ω j vvn< − Ω m Ω i> Ω m Φ (1) ω (Ω i ) g ( ı (2 ω n + Ω m + Ω i )) .
3. Calculation of contribution Q (2 , ± +) in Eq. (D9) Here we discuss the contribution Q (2 , ± +) introducedin Eq. (D9). The analytical structure of this diagramdemands ω n + Ω m > , ω n + Ω m + Ω i < ω n > < ω n < − Ω m − Ω i and Ω i < − Ω m . (D24) For symmetry reasons we write both versions of the in-tegrals (diffusons on grain 1 & 2) and therefore alreadytake the factor 2 in front of the Q (2) subdiagrams into ac-count in the right hand side of Eq. (D9). Furthermore weintroduce the short notation for the Green’s functions: a ± ( ξ ) ≡ ı ( ω n + Ω m ) − ξ ± ı/ (2 τ ) , (D25) b ± ( ξ ) ≡ ı ( ω n + Ω m + Ω i ) − ξ ± ı/ (2 τ ) ,c ± ( ξ ) ≡ ıω n − ξ ± ı/ (2 τ ) . In the following we need only the poles of functions b ± ( ξ )and c ± ( ξ ): ξ ( b ) ± = ı ( ω n + Ω m + Ω i ) ± ı/ (2 τ ) ∈ C ± and ξ ( c ) ± = ıω n ± ı/ (2 τ ) ∈ C ± . We also use the function g introduced in Eq. (D9).Therefore the ξ -integrals in Eq. (D9) can be writtenas: Z dξ dξ g ( ξ + ξ ) h(cid:0) a ( ξ ) b − ( ξ ) c + ( ξ ) (cid:1) − + (cid:0) a ( ξ ) b − ( ξ ) c + ( ξ ) (cid:1) − i (D26)= Z dξ dξ g ( ξ + ξ ) a ( ξ ) b − ( ξ ) c + ( ξ ) + a ( ξ ) b − ( ξ ) c + ( ξ ) a ( ξ ) b − ( ξ ) c + ( ξ ) a ( ξ ) b − ( ξ ) c + ( ξ )= 2 πı Z dξ a ( ξ ) b − ( ξ ) c + ( ξ ) g (cid:16) ξ + ξ ( c )+ (cid:17) c + ( ξ )2 − g (cid:16) ξ + ξ ( b ) − (cid:17) a ( ξ ) b − ( ξ ) a (cid:16) ξ ( b ) − (cid:17) = 2 πı Z dξ g (cid:16) ξ + ξ ( c )+ (cid:17) a ( ξ ) b − ( ξ ) − g (cid:16) ξ + ξ ( b ) − (cid:17) a (cid:16) ξ ( b ) − (cid:17) c + ( ξ ) . Here we executed the ξ -integral and used the fact thatthe first term in the second line has only one pole in C + (factor πı ). We are left with two term where the secondone has only a single ξ -pole in C + which gives another factor πı . However, both terms give the same result4 π g (cid:16) ξ ( c )+ + ξ ( b ) − (cid:17) a (cid:16) ξ ( b ) − (cid:17) = 4 π g ( ı (2 ω n + Ω m + Ω i ))( − ı Ω i + ı/τ ) ≈ − π τ g ( ı (2 ω n + Ω m + Ω i )) . (D27)1Substituting the result of Eq. (D27) into Eq. (D9) weobtain 2 Q (2 , ± +) = 2 π s d et T τ a d +2 (cid:16) ν (0) d (cid:17) (D28) × X <ω n < − Ω m − Ω i Ω i< − Ω m Φ (2) ω (Ω i ) g ( ı (2 ω n + Ω m + Ω i )) . The notation introduced in this Appendix allows us towrite down the Q (2 , ∓− ) and Q (2 , ±− ) contributions inEq. (D9) by just changing the +/- indices.
4. Calculation of contribution Q (2 , ∓− ) in Eq. (D9) Here we discuss the contribution Q (2 , ∓− ) introducedin Eq. (D9). The analytical structure of this diagram demands ω n + Ω m < , ω n + Ω m + Ω i > , (D29) ω n < , which defines the limits of frequency summations: − Ω m − Ω i < ω n < − Ω m and Ω i > . (D30)Using notations introduced in Eq. (D25) the ξ -integralsin Eq. (D9) can be calculated as: Z dξ dξ g ( ξ + ξ ) h(cid:0) a − ( ξ ) b + ( ξ ) c − ( ξ ) (cid:1) − + (cid:0) a − ( ξ ) b + ( ξ ) c − ( ξ ) (cid:1) − i (D31)= Z dξ dξ g ( ξ + ξ ) a − ( ξ ) b + ( ξ ) c − ( ξ ) + a − ( ξ ) b + ( ξ ) c − ( ξ ) a − ( ξ ) b + ( ξ ) c − ( ξ ) a − ( ξ ) b + ( ξ ) c − ( ξ )= 2 πı Z dξ a − ( ξ ) b + ( ξ ) c − ( ξ ) − g (cid:16) ξ + ξ ( c ) − (cid:17) c − ( ξ )2 + g (cid:16) ξ + ξ ( b )+ (cid:17) a − ( ξ ) b + ( ξ ) a − (cid:16) ξ ( b )+ (cid:17) = 2 πı Z dξ − g (cid:16) ξ + ξ ( c ) − (cid:17) a − ( ξ ) b + ( ξ ) + g (cid:16) ξ + ξ ( b )+ (cid:17) a − (cid:16) ξ ( b )+ (cid:17) c − ( ξ ) = 4 π g (cid:16) ξ ( c ) − + ξ ( b )+ (cid:17) a − (cid:16) ξ ( b )+ (cid:17) ≈ − π τ g ( ı (2 ω n + Ω m + Ω i )) . Substituting the result of Eq. (D31) into Eq. (D9) weobtain 2 Q (2 , ∓− ) = 2 π sd et T τ a d +2 (cid:16) ν (0) d (cid:17) (D32) × X − Ω m − Ω i <ω n < − Ω m Ω i> Φ (2) ω (Ω i ) g ( ı (2 ω n + Ω m + Ω i )) .
5. Calculation of contribution Q (2 , ±− ) in Eq. (D9) Here we discuss the contribution Q (2 , ±− ) introducedin Eq. (D9). The analytical structure of this diagram demands ω n + Ω m > , ω n + Ω m + Ω i < , (D33) ω n < , which defines the limits of frequency summations: − Ω m < ω n < i < − Ω m , (D34)and the disjunct region − Ω m < ω n < − Ω m − Ω i and − Ω m < Ω i < . (D35)Using notations introduced in Eq. (D25) the ξ -integralsin Eq. (D9) for this diagram can be calculated as:2 Z dξ dξ g ( ξ + ξ ) h(cid:0) a ( ξ ) b − ( ξ ) c − ( ξ ) (cid:1) − + (cid:0) a ( ξ ) b − ( ξ ) c − ( ξ ) (cid:1) − i (D36)= Z dξ dξ g ( ξ + ξ ) a ( ξ ) b − ( ξ ) c − ( ξ ) + a ( ξ ) b − ( ξ ) c − ( ξ ) a ( ξ ) b − ( ξ ) c − ( ξ ) a ( ξ ) b − ( ξ ) c − ( ξ )= 2 πı Z dξ a ( ξ ) b − ( ξ ) c − ( ξ ) − g (cid:16) ξ + ξ ( c ) − (cid:17) c − ( ξ )2 + g (cid:16) ξ + ξ ( b ) − (cid:17) a ( ξ ) b − ( ξ ) a (cid:16) ξ ( b ) − (cid:17) = 2 πı Z dξ − g (cid:16) ξ + ξ ( c ) − (cid:17) a ( ξ ) b − ( ξ ) − g (cid:16) ξ + ξ ( b ) − (cid:17) a (cid:16) ξ ( b ) − (cid:17) c − ( ξ ) = − π g (cid:16) ξ ( c ) − + ξ ( b ) − (cid:17) a (cid:16) ξ ( b ) − (cid:17) ≈ π τ g ( ı (2 ω n + Ω m + Ω i − /τ )) . Substituting the result of Eq. (D36) into Eq. (D9) taking into account the two disjunct regions for the frequencysummations, we obtain2 Q (2 , ±− ) = 2 π sd et T τ a d +2 (cid:16) ν (0) d (cid:17) X − Ω m <ω n < Ω i< − Ω m + X − Ω m <ω n < − Ω m − Ω i − Ω m< Ω i< Φ (2) ω (Ω i ) g ( ı (2 ω n + Ω m + Ω i − /τ )) . (D37)
6. Analytical continuation
Here we combine all five contributions to calculate correction η (1) : Q (1 , ±± ) , Q (1 , ∓∓ ) , Q (2 , ± +) , Q (2 , ∓− ) , and Q (2 , ±− ) introduced in the right hand sides of Eqs. (D7), (D9) and discussed in the previous subsections (D 1) through (D 5).We focus in particular on the analytical continuation of the Matsubara to real frequencies, both for the Fermionicand Bosonic frequencies. Using Eqs. (D18), (D22), (D28), (D32), and (D37) one can write Q (1) + Q (2) = Q (1 , ±± ) + Q (1 , ∓∓ ) + 2 (cid:16) Q (2 , ± +) + Q (2 , ∓− ) + Q (2 , ±− ) (cid:17) (D38)= λT − X <ω n < − Ω m − Ω i Ω i< − Ω m τ Φ (1) ω (Ω i ) g ( ı (2 ω n + Ω m + Ω i )) − X − Ω i <ω n < − Ω m Ω i> Ω m τ Φ (1) ω (Ω i ) g ( ı (2 ω n + Ω m + Ω i ))+ X <ω n < − Ω m − Ω i Ω i< − Ω m τ Φ (2) ω (Ω i ) g ( ı (2 ω n + Ω m + Ω i )) + X − Ω m − Ω i <ω n < − Ω m Ω i> τ Φ (2) ω (Ω i ) g ( ı (2 ω n + Ω m + Ω i )) − X − Ω m <ω n < Ω i< − Ω m τ Φ (2) ω (Ω i ) g ( ı (2 ω n + Ω m + Ω i − /τ )) − X − Ω m <ω n < − Ω m − Ω i − Ω m< Ω i< τ Φ (2) ω (Ω i ) g ( ı (2 ω n + Ω m + Ω i − /τ )) , where we introduced the notation λ = π sd et T a d +2 (cid:16) ν (0) d (cid:17) and the functions Φ (1) ω (Ω i ) and Φ (2) ω (Ω i ) were definedin Eqs. (D4) and (D5) respectively.Next, we perform the summation over Fermionic frequencies ω n in Eq. (D38) by shifts and the corresponding3analytical continuation X <ω n < − Ω m − Ω i Ω i< − Ω m f ( ω n , Ω i ) = X ω n > Ω i< − Ω m − X ω n > − Ω m − Ω i Ω i< − Ω m f ( ω n , Ω i ) = X ω n > Ω i< − Ω m [ f ( ω n , Ω i ) − f ( ω n − Ω m − Ω i , Ω i )] , (D39) X − Ω i <ω n < − Ω m Ω i> Ω m f ( ω n , Ω i ) = X ω n < − Ω m Ω i> Ω m − X ω n < − Ω i Ω i> Ω m f ( ω n , Ω i ) = X ω n < Ω i> Ω m [ f ( ω n − Ω m , Ω i ) − f ( ω n − Ω i , Ω i )] , X − Ω m − Ω i <ω n < − Ω m Ω i> f ( ω n , Ω i ) = X ω n < − Ω m Ω i> − X ω n < − Ω m − Ω i Ω i> f ( ω n , Ω i ) = X ω n < Ω i> [ f ( ω n − Ω m , Ω i ) − f ( ω n − Ω i − Ω m , Ω i )] , X − Ω m <ω n < Ω i< − Ω m f ( ω n , Ω i ) = X ω n < Ω i< − Ω m − X ω n < − Ω m Ω i< − Ω m f ( ω n , Ω i ) = X ω n < Ω i< − Ω m [ f ( ω n , Ω i ) − f ( ω n − Ω m , Ω i )] , X − Ω m <ω n < − Ω m − Ω i − Ω m< Ω i< f ( ω n , Ω i ) = X ω n < − Ω m − Ω i − Ω m< Ω i< − X ω n < − Ω m − Ω m< Ω i< f ( ω n , Ω i )= X ω n < − Ω m< Ω i< [ f ( ω n − Ω m − Ω i , Ω i ) − f ( ω n − Ω m , Ω i )]= X ω n < Ω i< [ f ( ω n − Ω m − Ω i , Ω i ) − f ( ω n − Ω m , Ω i )] − X ω n < Ω i< − Ω m [ f ( ω n − Ω m − Ω i , Ω i ) − f ( ω n − Ω m , Ω i )] . Here the function f ( ω n , Ω i ) is the product of the func-tions τ Φ ( α ) ω (Ω i ) and g . For the analytic continuation weneed to consider the ω -dependence of function Φ ( α ) ω (Ω i ),in particular we can use the fact that (which follows di-rectly from the definition of functions Φ ( α )0 in Eqs. (D4)and (D5)) Φ ( α )0 τ = Φ ( α ) ω τ ω , (D40) where Φ ( α )0 is the same as Φ ( α ) ω but with τ ω replaced by τ and is therefore ω -independent. Using Eq. (D40), weobtain4 Q (1) + Q (2) = Q (1 , ±± ) + Q (1 , ∓∓ ) + 2 (cid:16) Q (2 , ± +) + Q (2 , ∓− ) + Q (2 , ±− ) (cid:17) (D41)= λT πıT " − X Ω i < − Ω m h Φ (1)0 (Ω i ) − Φ (2)0 (Ω i ) i τ Z dω tanh (cid:16) ω T (cid:17) [ g (2 ω + ı (Ω m + Ω i )) − g (2 ω − ı (Ω m + Ω i ))]+ X Ω i > Ω m Φ (1)0 (Ω i ) τ Z dω tanh (cid:16) ω T (cid:17) [ g (2 ω + ı ( − Ω m + Ω i )) − g (2 ω + ı (Ω m − Ω i ))] − X Ω i > Φ (2)0 (Ω i ) τ Z dω tanh (cid:16) ω T (cid:17) [ g (2 ω + ı ( − Ω m + Ω i )) − g (2 ω − ı (Ω m + Ω i ))]+ X Ω i < − Ω m Φ (2)0 (Ω i ) τ Z dω tanh (cid:16) ω T (cid:17) [ g (2 ω + ı (Ω m + Ω i − /τ )) − g (2 ω + ı ( − Ω m + Ω i − /τ ))]+ X Ω i < Φ (2)0 (Ω i ) τ Z dω tanh (cid:16) ω T (cid:17) [ g (2 ω − ı (Ω m + Ω i + 1 /τ )) − g (2 ω + ı ( − Ω m + Ω i − /τ ))] − X Ω i < − Ω m Φ (2)0 (Ω i ) τ Z dω tanh (cid:16) ω T (cid:17) [ g (2 ω − ı (Ω m + Ω i + 1 /τ )) − g (2 ω + ı ( − Ω m + Ω i − /τ ))] . Next, we consider the integrands in Eq. (D41) and introduce the short hand notations a ≡ g (2 ω + ı (Ω m + Ω i )) − g (2 ω − ı (Ω m + Ω i )) = 2 ı (Ω m + Ω i ) (1 + 2 dω/ε F ) , (D42) a ≡ g (2 ω + ı ( − Ω m + Ω i )) − g (2 ω + ı (Ω m − Ω i )) = 2 ı (Ω i − Ω m ) (1 + 2 dω/ε F ) ,a ≡ g (2 ω + ı ( − Ω m + Ω i )) − g (2 ω − ı (Ω m + Ω i )) = 2 ı Ω i + 2Ω i d (2 ωı + Ω m ) /ε F ,a ≡ g (2 ω + ı (Ω m + Ω i − /τ )) − g (2 ω + ı ( − Ω m + Ω i − /τ )) = 2Ω m ı + 2Ω m d/ε F (1 /τ ω + 2 ıω − Ω i ) ,a ≡ g (2 ω − ı (Ω m + Ω i + 1 /τ )) − g (2 ω + ı ( − Ω m + Ω i − /τ )) = − ı Ω i − d Ω i /ε F (1 /τ ω + 2 ıω − Ω m ) . Now, we extract only terms which are linear in ω inEq. (D42) and of order 1 /ε F (thus the τ ω in a and a does not give a contribution). Therefore we obtain a ≃ ωı (Ω m + Ω i ) 2 d/ε F (D43) a ≃ ωı (Ω i − Ω m ) 2 d/ε F a ≃ ωı Ω i d/ε F a ≃ ωı Ω m d/ε F a ≃ − ωı Ω i d/ε F . Substituting the result of Eqs. (D43) back into Eq. (D41)we obtain5 Q (1) + Q (2) = Q (1 , ±± ) + Q (1 , ∓∓ ) + 2 (cid:16) Q (2 , ± +) + Q (2 , ∓− ) + Q (2 , ±− ) (cid:17) (D44)= 2 τ ı dλT πıT ε F " − X Ω i < − Ω m h Φ (1)0 (Ω i ) − Φ (2)0 (Ω i ) i Z dω tanh( ω/ T ) ω (Ω m + Ω i )+ X Ω i > Ω m Φ (1)0 (Ω i ) Z dω tanh( ω/ T ) ω (Ω i − Ω m ) − X Ω i > Φ (2) ω (Ω i ) Z dω tanh( ω/ T ) ω Ω i + X Ω i < − Ω m Φ (2)0 (Ω i ) Z dω tanh( ω/ T ) ω Ω m − X Ω i < Φ (2)0 (Ω i ) Z dω tanh( ω/ T ) ω Ω i + X Ω i < − Ω m Φ (2)0 (Ω i ) Z dω tanh( ω/ T ) ω Ω i = e λ " − X Ω i < − Ω m h Φ (1)0 (Ω i ) − Φ (2)0 (Ω i ) i (Ω m + Ω i ) + X Ω i > Ω m Φ (1)0 (Ω i ) (Ω i − Ω m ) − X Ω i > Φ (2)0 (Ω i ) Ω i + X Ω i < − Ω m Φ (2)0 (Ω i ) Ω m − X Ω i < Φ (2)0 (Ω i ) Ω i + X Ω i < − Ω m Φ (2)0 (Ω i ) Ω i . In the last two lines of Eq. (D44) we introduced the no-tation e λ = τ dλπε F R dωω tanh( ω/ T ) where we use the factthat R dωω tanh( ω/ T ) −→ = − ( πT ) , neglecting the in-finite boundary terms of the partial integration.Using the fact that the functions Φ (2)0 (Ω i ) in Eq. (D44)depend only on the absolute value of the Bosonic fre-quency Ω, we finally obtain Q (1) + Q (2) (D45)= − λ X Ω i < − Ω m h Φ (1)0 (Ω i ) − Φ (2)0 (Ω i ) i (Ω m + Ω i ) , where we use the notation λ = − π s e (cid:16) τ t ν (0) d a d +1 (cid:17) T ε F . The internal frequencysummation is also done by analytical continuation, butfor Bosonic frequencies: Ω i −→ − ı e Ω + η . Before thefinal integration we should take the external frequencyderivative and finally calculate the q -integrals of theΦ -functions. Using Eqs. (D4) and (D5) we haveΦ (1)0 (Ω i ) − Φ (2)0 (Ω i ) = a d Z d q (2 π ) d (D46) × E c ( q ) h d − P ′ a cos( q · a ) i τ ( | Ω i | + 4 E c ( q ) ǫ q ) ( | Ω i | + ǫ q δ ) . Using Eq. (D46) one can calculate the sum over internalbosonic frequencies Ω i in Eq. (D45) S Ω = X Ω i < − Ω m Ω m + Ω i ( | Ω i | + 4 E c ( q ) ǫ q ) ( | Ω i | + ǫ q δ ) (D47)= X Ω i < − Ω m Ω m + Ω i ( − Ω i + 4 E c ( q ) ǫ q ) ( − Ω i + ǫ q δ )= X Ω i < Ω i ( − Ω i + Ω m + 4 E c ( q ) ǫ q ) ( − Ω i + Ω m + ǫ q δ )= Z [ − πıT ] − ( − ı e Ω) coth( e Ω / T ) d e Ω[ ı e Ω + Ω m + 4 E c ( q ) ǫ q ][ ı e Ω + Ω m + ǫ q δ ]= Z [ − πT ] − e Ω coth( e Ω / T ) d e Ω[ e Ω − ı (Ω m + 4 E c ( q ) ǫ q )][ e Ω − ı (Ω m + ǫ q δ )]= Z [ − πT ] − x coth( x ) dx [ x − ı T (Ω m + 4 E c ( q ) ǫ q )][ x − ı T (Ω m + ǫ q δ )] . Only the terms proportional to the external frequencyΩ m in Eq. (D47) contribute to the correction to the ther-moelectric coefficient η (1) in Eq. (D2). Taking the deriva-tive of both sides of Eq. (D47) we obtain6 ∂∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 S Ω = − πT Z ∂∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 x coth( x ) dx (cid:0) x − Ω2 T − ıT E c ( q ) ǫ q (cid:1) (cid:0) x − Ω2 T − ı T ǫ q δ (cid:1) (D48)= − πT Z x coth( x ) dx " x − ıa ) ( x − ıb ) + 1( x − ıa ) ( x − ıb ) = − πT Z x coth( x ) dx [2 x − ı ( a + b )]( x − ıa ) ( x − ıb ) ≡ − πT I a,b , where a = 2 E c ( q ) ǫ q /T and b = ǫ q δ/ (2 T ). Finally, we usethe approximation: x coth( x ) ≈ (cid:26) , | x | < | x | , | x | ≥ , (D49)and obtain for the integral I a,b in Eq. (D48) the followingresult I a,b = ı ln (cid:2)(cid:0) a (cid:1) / (cid:0) b (cid:1)(cid:3) ( a − b ) . (D50)For a ≫ ≫ bI a,b = ı aa = ı T ln [2 E c ( q ) ǫ q /T ] E c ( q ) ǫ q . (D51)Thus, for Eq. (D48) we obtain ∂∂ Ω (cid:12)(cid:12)(cid:12)(cid:12) Ω=0 S Ω = − ı πT ln [2 E c ( q ) ǫ q /T ] E c ( q ) ǫ q . (D52)Substituting Eq. (D52) into Eq. (D46) we obtain for thecorrection to the thermoelectric coefficient the followingresult η (1) = − e λ πT (2 π ) d τ (D53) × Z d q h d − P ′ a cos( q · a ) i ln [2 E c ( q ) ǫ q /T ] ǫ q .
7. Final integration and results
In Eq. (D53) we are left with the final integration overinternal momenta q . Therefore we need the functionaldependence of ǫ q and E c ( q ) on q , which were introducedaround Eq. (D6) as ǫ q = 2 g T d − ′ X a cos( q · a ) ! (D54) E c ( q ) = e C ( q ) = e a d − ln( qa ) , d = 1 π/q, d = 22 π/q , d = 3 . (D55) The final q -integral can therefore be written as η (1) = − e λ πT (2 π ) d g T τ Z d q ln [2 E c ( q ) ǫ q /T ] . (D56)The q -integral is cutoff at q = π/a . In d = 3 the ln -argument is finite at q = 0. In d = 2 the volume elementmakes the integrand at q = 0 finite. The coefficient inEq. (D56) can be simplified to − e λ πT (2 π ) d g T τ = π s π ) d et a d +2 g T (cid:16) ν (0) d (cid:17) Tε F = − η (0) πg T (cid:16) a π (cid:17) d , (D57)such that in d = 3 we obtain: η (1)3 D = − η (0) πg T (cid:16) a π (cid:17) Z d q (D58) × ln (cid:20) πe g T aT − cos( aq x ) − cos( aq y ) − cos( aq z )( aq ) (cid:21) = − η (0) πg T (cid:16) a π (cid:17) π a (cid:20) π (cid:18) e a g T π/T (cid:19) + c (cid:21) . Here c = R d q ln h (3 − cos( πq x ) − cos( πq y ) − cos( πq z )) / ( πq ) i and integration is over the unit sphere. However, thisnumerical constant (and all other inside the logarithm)can be neglected, since E c /T ≫ g T ≫ d = 2 we obtain η (1)2 D = − η (0) πg T (cid:16) a π (cid:17) Z d q (D59) × ln (cid:20) e g T πT a (3 − cos( aq x ) − cos( aq y )) / ( aq ) (cid:21) = − η (0) πg T (cid:16) a π (cid:17) π a (cid:20) ln 8 e g T πT a + c (cid:21) , with c = R d q π ln h − cos( πq x ) − cos( πq y πq i . a. Final results Now, we can summarize our final results in d = 2 , η (0) = − sπ et a d +2 (cid:16) ν (0) d (cid:17) Tε F , (D60) η (1)3 D = sπ ea g T t (cid:16) ν (0)3 (cid:17) Tε F ln (cid:18) e g T T a (cid:19) (D61)= − η (0) g T ln ( E c g T /T ) ,η (1)2 D = sπ ea g T (cid:16) t ν (0)2 (cid:17) Tε F ln (cid:18) e g T T a (cid:19) (D62)= − η (0) g T ln ( E c g T /T ) , or combined in a compact way (valid for d = 2 ,
3) weobtain the following final result for the thermoelectriccoefficient of granular metals η = η (0) (cid:18) − πg T ln g T E c T (cid:19) . (D63) D. M. Rowe and C. M. Bhandari, Modern thermoelectrics(Reston, Reston, VA, 1983). A. A. Abrikosov, Fundamentals of the theory of metals(North Holland, Amsterdam, 1988). G. D. Mahan, Many-particle phyics, (Kluwer, New York,2000) G. D. Mahan, J. Appl. Phys. , 4551 (1991). L. E. Bell, Science , 1457 (2008). F. DiSalve, Science , 703 (1999). B. Sales, Science , 1248 (2002). T. Tritt and M.A. Subramanian, MRS Bull. , 188 (2006). G.J. Snyder and E. Toberer, Nat. Mater. , 105 (2008). A. Majumdar, Science , 777 (2004). R. Venkatasubramanian, E. Siivola, T. Colpitts, and BO’Quinn, Nature , 597 (2001). T. C. Harman, P. J. Taylor, M. P. Walsh, B. E. LaForge,Science , 2229 (2002). C. Collier, T. Vossmeyer, and J. Heath, Ann. Rev. of Phys.Chem. , 371 (1998). C. Murray, C. Kagan, and M. Bawendi, Ann. Rev. ofMaterials Science, , 545 (2000). T. Tran, I. S. Beloborodov, X. M. Lin, V. M. Vinokur, andH. M. Jaeger, Phys. Rev. Lett. , 076806 (2005). T. E. Humphrey and H. Linke, Phys. Rev. Lett. , 096601(2005). D. A. Broido and T. L. Reinecke, Phys. Rev. B , 13797 (1995). I. S. Beloborodov, A. V. Lopatin, V. M. Vinokur, and K.B. Efetov, Rev. Mod. Phys. , 469 (2007). I. S. Beloborodov, K. B. Efetov, A. V. Lopatin, andV. M. Vinokur Phys. Rev. Lett. , 246801 (2003). I. S. Beloborodov, A. V. Lopatin, F. W. J. Hekking,R. Fazio, and V. M. Vinokur, Europhys. Lett. , 435(2005). A. Glatz and I. S. Beloborodov, Phys. Rev. B ,041404(R) (2009). K. B. Efetov and A. Tschersich, Europhys. Lett. , 114,(2002); Phys. Rev. B , 174205 (2003). M. Jonson and G. D. Mahan, Phys. Rev. B , 4223 (1980). V. Tripathi and Y. L. Loh, Phys. Rev. Lett. , 046805(2006). K. F. Hsu, S. Loo, F. Guo, W. Chen, J. S. Dyck, C. Uher,T. Hogan, E. K. Polychroniadis, M. G. Kanatzidis, Science , 818 (2004). B. Poudel, Q. Hao, Y. Ma, Y. Lan, A. Minnich, B. Yu,X. Yan, D. Wang, A. Muto, D. Vashaee, X. Chen, J. Liu,M. S. Dresselhaus, G. Chen, Z. Ren, Science , 634(2008). J. L. Mi, X. B. Zhao, T. J. Zhu, and J. P. Tu, Appl. Phys.Lett.91