A N O T E ON WICK'S T H E O R E M

WILLI-H. STEEB

Institute of Theoretical Physics, University of Kid, West Germany

ABSTRACT. In analogy to Gaudin, but in a more complicated case, Wick's theorem in statistical mechanics is proved by using the commutation rules. As a special case we obtain the result of Gaudin. An application to the Hubbard model is given.

The simplest Hamiltonian for discussing strong correlations in narrow energy bands is the so-called Hubbard model. [1 ] This model describes a lattice model for an electron gas with the Coulomb repulsion completely screened out except between electrons on the same site. The Hamiltonian can be written in its simplest form in Wannier representation as H=-t

~ ~ii+Ao.Cio+~. ~/if~/i~,, iAa t

§

where l'/ia." = CRrCitr. The index i labels the lattice points and A the lattice vectors to the nearest

neighbours. Thus the model contains both itinerant and localized aspects. Various many-body techniques have been tried on the problem of approximating electron-electron correlations. The high temperature perturbation expansion has been considered by Beni et al. [2]. Here the kinetic part is to be treated as a perturbation, and the interacting part is considered as the exactly soluble part. This means that we can calculate exactly the grand thermodynamic potential. Now the aim of the present paper is to develop a Wick's theorem for this case. This theorem is also useful for calculation of upper bounds for the grand potential [3] of the Hubbard model using the inequality ~2 ~< Tr(H -/aNe)W t + (1//3)TrWt s

where Wt, the so-called trial density matrix, is given by

Wt = exp(-k~ ~r/i~,r~i~ - k~ Z. r/i,r/=, kB ~.p r/ij.)/Tr , ex ( . . . . . . . . . . . .

).

It should be noted that, although there are investigations of this density matrix [4], a theorem for actual calculations is lacking. In order to obtain a Wick's theorem for this case we use a method which has been developed by Gaudin [5]. He proved Wick's theorem, where the grand canonical density matrix has been given by W = exp (fl(~2o - K0)) with K0 = ~ E k c*kck. The quantity ~2o is the thermodynamic potential for the Hamiltonian Ko. The chemical potential ~s included in Ek, and c~ and c k are the creation and annihilation operators for the particles of the system. The c~ and c k may have commutation relations appropriate to either Bose-Einstein or Fermi-Dirac statistics. k

Letters in Mathematical Physics 1 (1976) 135-139.All Rights Reserved. Copyright 9 1976 by D. Reidel Publishing Company, Dordrecht-Holland.

.

.

135

In similiar fashion to Gaudin we consider here a more complicated case. Clearly, since we are considering electron-electron correlations we are dealing only with fermions in this paper. Our grand canonical density matrix is given by W = exp(/3(~20 - Ko)) with Ko = 2`! ~.nit ni, + 2`2 ~.nit + 1 2`3~.ni,. It should be emphasized that in this paper and in the Gaudin one, only the sta~tic case, 1 , . concerning normal products, xs considered. +

+

Thus we want to evaluate Tr(al a2 a3 ..... an W) where aiC-(cjr, cD', c j4, cj , ; j = 1, 2 .... N). A straighforward calculation shows that [Ko, Nnio] = 0. Therefore, since Ko commutes with ~i ni~, and ~ni, , the trace vanishes unless the set (a~ a2 ...an) contains an equal number of creation and annitulation operators with spin up and also an equal number of creation and annihilation operators with spin down. This means, for instance, that terms like Tr(c~t cj+ W) vanish. The anticommutators are either 1 or 0 and may therefore be taken outside the trace. Now we introduce a useful abbreviation, namely aij: = ala2a3 ..... ai.l ai+ 1 ...aj_ 1aj+ 1...an_ 1a n .

(1)

Thus we can write the following identity n

Tr(al a~ a3...anW) = ~

( - 1 ) k [al, a k ] Tr(al k w)

k=2

(2)

+

Tr(a2 a3...anal W).

-

This means that we commute al successively to the right. Now we must evaluate the last term on the right hand side of Equation (2). A simple calculation shows that e~ko cT~.e-gko = e r

(c~t + c~,ni~(e#X, - 1 ) ) ;

(3)

eXko c;~e #k~ = eX#~ (c;~ + c;, n i, (e~X1-1)). If we consider the annihilation operators we must replace ?q, 2`2, and 2,3 by -X2, -2`2, and -2`3. From now on Xi(i = 1,2,3) contains the quantity/3. With no loss of generality, let a~ be a creation operator. Clearly, using the A.C.R., we can always obtain this form. Using the cyclic property of the trace we obtain by means of Equation [3] Tr(az a3 a4 ...anal W) -= eXTr(a~ a2 a3 ...anW) +eX(e x, -1)Tr(bl b*~al .a2 ...anW) , (4) where +

x:

h2 for a~ = c i t +

9

(5)

2`3 for al = ci, If al = c~o, then b~ = c~_ o and b'~ is the adjoint ofb~. Moreover, we have used the fact that al b lb'~ = b~ b~ia~. Therefore we must investigate the last term on the right hand side o f Equation (4). Again we write down the identity (due to Equation (2)). 136

Tr(b ~b] a~ a~...anW ) =

Tr(a~ a~ as... anW)

+ [bt, a~] Tr(b+~ata3...an W) - [bt, as]Tr(b]ala~a,~...anW) +

+

+ [b~, an]+ "Tr(b]ata2...an_lW ) - Tr(b]a,a:...anb,W ))

(6)

Because o f b ~ b ] = ni_o, we have Ibm, b] ] = 1. Clearly [bt, a~ ] = 0, and the term on the right hand _ + ++ side vanishes. Using a~ a~ b~ - bt al a~, [a~, a~ ] = 1, and a~ = 0, the last term on the right hand side +

of Equation (6) can be written by means of Equation (3) as

(7)

Tr(b~i at a= ...anb ~W) = eM + x, + x, - x Tr(b i b] at a~ ...anW ) . Note that a] is the adjoint of at and bt is a creation operator. Up to now we have eliminated bt inside the trace. As the next step we must eliminate b] (see

Equation (6)). First we compute the second term on the right hand side of Equation (6). For this purpose we use again the identity which is given by Equation (2). Using [at, a] ] = 1 and a~ = 0, + the last term on the right hand side in this identity can be written, as described above, as

(8)

Tr(at asa,~...anb+xW) = e -(x~ + x, + xa - x ) Tr(b] a~ aa...anW) 9 Note that b~ is an annihilation operator.~

In the same manner we calculate the other terms. Note that [a~, b~ ] = 0. Therefore we have com. + pletely eliminated both bt and b~ inside the trace. Thus we can write down our result. Introducing the notation Zo : = e x, + M + xa + eX~ § M + eM + x~ + 1 and using the abbreviation (1), we have

Tr(axa2as...anW) =

eXl x~ +x~- x § 1 n 9 Y~ ( - 1 ) k [at, ak]" Tr(a~ ~W) Z0

ex(e x, - l )

k=2

+

1

e-(Xt +M +ha - x ) + 1 "Z---~. { [ b , , a 2 ] ~k=3(--1)k[b]' ak]" Tr(a2kW) i-1

+[b,ai]'(Z +

n

(--1)~+ltb],a~l'Tr(ai~W)+~; k=2

+

"-'

+[b~,an]'E + k=2

(--1) k+l [b%,ak].Tr(ank W)

(--l)k[b~i, aklTr(aikW))

k=i+l

}

.

+

(9)

+

The upper (lower) signs refer to i even (odd). (2 < i < n). Clearly, if we put Xl equal to 0, then we obtain the result of Gaudin. A helpful relation can be deduced from Equation (9). Let the set M = (ni t, ni 4, ni~"ni4 ; i = 1 ,..., N) be given. If A,B,C,..., Z ~ M and the subscripts, i, of these members are pairwise distinct, then we have Tr(AB ........ ZW) = Tr(AVO" Tr(BW) ....... Tr(ZW)

(lO) 137

In other words, the average of any product of operators factorizes into the product of averages of the operators related to each site. Furthermore, we can see an important special case of Equation (9). If all ai have the same spin our formula reduces to

Tr(al a~ ...anW) --

r

+X2+~., -x + I .~ (_l)k[a,,ak]Tr(a,kW)" Zo k=2 +

(11)

Let us now consider severalexamples. First we investigate three basic cases, namely Tr(ni,W), Tr(ni,W), and Tr(ni~.ni,W). At once we obtain Tr(ni,W) = (e x, +x3 + 1)/Zo; Tr(ni,W) = (e x, +x~ + l)/Zo. In the third case we have Tr(nit ni~W) = 1/Zo. Clearly these results can also be obtained from the grand potential [20. For K the grand potential per lattice site is given by

~'0 N

-(1/3)In(1 +exp(-3~u) + exp(-3),3) + exp(-3(2,, + X~ +),3)))

We get the average values, using the well-known identity

~ exp(-3Ko) = - l d r exp((r-3)Ko ~OKo exp(-rKo). bX---~, Note that in this case we have not included/3 in ~'i" Finally we want to consider a more complicated example. In the high temperature expansion of the Hubbard model, we have to calculate in fourth order terms as

Z ( c.+~a+A~ tCi, ~'C+i~+A~Cia ~'C+ia+A~ ~Ci~Ci~ + +A 4 $Ci~~), ilAli2~2 i 3 A~ i,~ A,~

where i k (k = 1,2,3,4) goes over all lattice sites, and Ak (k = 1,2,3,4) is summed over nearestneighbour vectors. Application of Equation (9) yields N~Z~((n, )(n, ) - ( n t ) 2. ( n , ) - (n~)(n,) 2 + (n,)U(n,) 2) +N Z2(-( nt)(n,) + 2( nt)2 n~) ~

02) -2( n t n~)(n~,) - 2(n~ n~)(n~) + 4( n t n~)( n~ )(n~) + ( n t n~)) +N Z ( 2 ( ( n ~ ) ( n & ) - ( n ~ n , ) ) u ) , where we made use of the A.C.R and the fact that 6o.a= 0. Moreover, the following identities have been taken into account: (exp(y) + 1) / (exp(-y) + 1) = exp(y); (n t)(n,) - (n, n , ) = 138

(exp(2k~ + k2 + k3) --exp ( ~ + ?~2 + k3))/Zg" N and z are respectively the number oflattice sites and number of nearest neighbours. It should be noted that, if we consider the cumulant expansion [6], terms on the right hand side of (12) which are proportional to N: vanish.

REFERENCES 1. 2. 3. 4. 5. 6.

Hubbard, J., Proc. Roy. Soc. A276, 238 (1963). Beni, G., Pincus, P., and Hone, D.,Phys. Rev. B8, 3389 (1973). Steeb, W.-H., and Marsch, E., Phys. Stat. Sol. (b)65,403 (t974). Westwanski, B., Phys. Lett. A45,449 (1973). Gaudin, M.,Nucl. Phys. 15, 89 (1960). Kubo, R.,J. Phys. Soc. Japan 17, 1101 (1962).

139