On the systematic errors of the pseudofermion algorithm

490

Nuclear Physics B (Proc. Suppl.) 9 (1989) 490-494 North-Holland, Amsterdam

On the Systematic Errors of the Pseudofermion Algorithm J. Potvin't, M. Campostrinitt, K.J. Moriarty "x and C. Rebbit tDepartment of Physics, Boston University, Boston, MA 02115; 'Department of Mathematics, Statistics and Computer Science, Dalhousie University, Halifax, Nova Scotia, B3H 3J5, Canada; * INFN Sezione di Piss, Italy; x John yon Neumann National Supercomputer Center and Institute for Advanced Studies, Princeton, NJ 03540 We report on a detailed study of the systematic errors in the pseudofermion algorithm on a l0s x 22 lattice at gauge coupling fl = 5.7 and with three flavors of staggered dynamical quarks of mass m = 0.10,0.05. The Wilson loop and the hadronlc mass spectrum were studied at four values of the acceptance (60, 70, 81, and 89%). We suggest a new parameter which is more accurate in describing the magnitude of the stepsize. The Wilson loop and the hadronic mass spectrum were also studied on a l0s x 24 volume with fl = 5.47,n! = 2, m = 0.05 and acceptance 80, 89%, in order to compare with the results of a calculation by Gottlieb, Liu, Renken, Sugar and Toussalnt with the hybrid algorithm.

1.

INTRODUCTION

The simulation of lattice QCD with dynamical fermions is a most difficult problem involving large amounts of computer time. So far, only approximate but fast algorithms have permitted the implementation of such a calculation [1]. Since the Brookhaven meeting where the first results in hadron spectroscopy in full QCD were presented [2], a good deal of work has been done in exploring further the effects of the errors involved in these algorithms for realistic values of the coupling, quark masses and lattice volumes. Most approximate algorithms currently used nowa days can be divided into two groups. There are those based on the Metropolis approach such as the pseudofermion [3] and the bush factorized [4] algorithms; there are also those based on the iteration of stochastic or deterministic equations of motion such as the Langevin equation [5], the hybrid method [6] and the microcanonical algorithm [7]. Most have been used in calculations of the hadron mass spectrum on lattices larger than 84 [8-10, 13-15]. What is important at this stage, however, is to understand the effects of the systematic errors inherent to these algorithms in spectroscopy calculations, and to assess their accuracy in relation to each other. Some of that pioneering work has already been done for the Langevin and the bush factorized method and for the pseudofermion algorithm at small volumes [11, 18] and with the Gross Neveu model [12]. We present in these proceedings an analysis of the systematic errors in the pseudofermlon algorithm on Wilson loops and hadron masses, as well as the results of a comparison with data obtained with 0920-5632/89/$03.50 © Elsevier Science Publishers B.V. (North-Holland Physics Publishing Division)

the hybrid algorithm [15]. Our study of the systematic errors on loops and hadron masses were done on a 10s x 32 lattice with gauge coupling 3 = s.7 and quark mass ma = 0.I0 and 0.05, with acceptances 50, 60, 70, 81, and 89%. The number of iterations in each runs were 9,000, 8,000, 15,000, 15,000, 41,000 iterations respectively. On the other hand, the comparison with the results of the hybrid algorithm was done on a l0 s x 24 lattice, with ~ = 5.4"/, ma = .05 and acceptance 80 and 92% (16,000, 12,000 iterations respectively), thus using; exactly the same parameters as in the hybrid study of Gottlieb et al. [15]. 2.

THE PSEUDOFERMION ALGORTIHM

The pseudofermion algorithm of Fucito, Marinari, Parisi and Rebbi [3] has been discussed extensively in the literature [3, 16]. Let us summarize here its salient features. The pseudofermion algorithm is a Metropolis procedure [17] in which the link variables U," are tentatively changed into U,"' with

R", = E,p(ipe~Ai)

(1)

the ®i being normal Gaussian random numbers and ~ su(3) Gell Mann matrices. This proposed change is accepted or rejected, depending on the sign of the variation in the full action S. [17]: 6s = sw,"'] - sT."]

J. Potvin et al./Systematic errors of the pseudofermion algorithm = 5s,.(s)+ ~-~T,[I,((D'+m)(D+'+m)-h*(D+m)(D + +m))] =

+

0.584

491

' I r ' ~ ~ ' I il I' W{I,I) 0.244

r

i

1

,

i

i

v-~

W(3,3)

-

m 2 + DD+

)]

(2)

(We follow the notation of l:tef. [16]). In the pseudofermion algorithm, Eq. 2 is approximated by

i

I

I

I00 90

I

I

I

1

80

70

'

'

jo.2 , t-

I

I

i

I

1

I00 9 0

60

I

80

1~

!

70

I

1

IJ

60

o.le4

6SPI =SS'=(3)+~Tr( D 6D+m = ++6--D-DDD+ D+I)

(3)

(with 6D = D' - D ---D[U,~'] - D[U,~]) i.e. a l~rst order

evaluation of the logarithm in the numerator of Eq. 2. This algorithm uses two additional approximations. First, the same matrix (m=+ DD+)~s used to update all the links in a given sweep through the lattice, or,in other words, that matrix is never changed when the change g : ' is accepted. Second, that matrix (m2+ DD+) -1 is calculated via an inner Monte Carlo simulation of a path integral featuring the SU(3) columns ¢~,) (the so-called "pseudofermions") and the action ¢(DD + +m)¢[3, 16, 17]. The pseudofermion algorithm becomes exact in the limit p ~ 0 (i.e. at acceptance 100%), and for a large number of sweeps through the O's in the inner Monte Carlo process. In the calculations described below, 8 Metropolis hits on the gauge variables were performed during each sweep. Moreover, 24 pseudofermions sweeps in the inner Monte Carlo were used. See Ref. [8] for more details. In most studies so far, the proper choice of the acceptance has been determined by looking for no variations in the Wilson loop for several values of the acceptance. In Fig. 1 we show that such criterion may be misleading since the loop appears to change more between 80 and 90% acceptance than in between 60 and 70%. Moreover, it has been rarely stressed that the algorithm can perform badly even at high acceptances when the quark mass is tuned to small values [18]. In itself, the acceptance proves to be a poor indicator of the pseudofermion algorithm's accuracy. In what follows, we would like to motivate a better parameter, and outline a procedure which extracts the results in the limit p ~ 0 from those data obtained at finite p.

3.

p/(,,~)=

A useful relationship between the exact evaluation of 65 (Eq. 2) and its pseudofermion approximation can be obtained via the expansion of (2) in powers of (D'D '+ - D D + ) / ( m 2 + DD+)[19]. Defining the ra-

W(2,2)

o.I 83

0.030

'

1

W(2,5)

I'

'''

o.I 82

..

•

o.I 81 0.180

, i IOO 90

i

'' ' ''

I

,

80

I~

I

TO

I

lOO2T

I

60

,

I

I00 90

,

I

80

,

!

TO

l

I

I

60

ACCEPTANCE W~

Fig. 1 The Wilson loop versus the acceptance.

tios of the fermionic determinants (S, is the fermion action [16]) as -

E , V (-6Sp)

~,~ -

~=pC-6Sf,')

and

and assuming p ~ i and p ,~ (m)', one can write

= Ape + ~-T~

\ r~=+ D D + ]

+0

L\ -~ ~-D-D-I (4)

Roughly speaking, each of these terms contribute like: Ape ~ 1 + p/m= 6D6D + rn2 + DD+

p2 m2

D'D'+-DD+)' ~

:~ -D-K;

p2 p3 p" ~ .~-X ' . . , ' ~ ,

Notice that the second approximation of the pseudofermion algorithm i.e. the use of the same (m' + DD) -1 for all the updates in a given sweep will contribute additional terms to (4) in powers of [19] p2 pS pS

~-~, ~-~, ~--~... Equation (4) defines a more natural measure of the typical step size involved in the update. To lowest order, the step size is proportional to p/m= and is more natural than "p" or the acceptance as a measure of

492

J. Potvin et al./Systematic errors of the pseudofermion algorithm

w[L~)

.s.-}\ \ .s86-

m 3. 3 :

\ \

B=5.70

\ .584.582- " ~ , ~ % m a .580 I

~ I 2 4 92

81

,'

= .IO ,' ~

6 8 10 70

89

•

I 30

I

20 81

'

I

70

! 50

! 40 .....

I

60

a/(m=)"

I $0

,I ACC

(%)

Fig. 2 The Wilson loop versus a/m= and the acceptance, for two v~lues of the quark mass.

the stepsize. A similar parameter has been suggested in the past in the context of the Langevin algorithm [5JAlso to lowest order in p/m 2, equation (4) shows how to extrapolate to the exact results (i.e. the limit p - . 0) from those obtained at finite p. Fig. 2 displays the Wilson loop W(1,1) versus p/m= for two quark masses (m = .I0, .05), and the possible extrapolations needed to calculate the exact value of W(1,1). Among other things, it illustrates how dangerous linear extrapolations can be for large values of p/m z. The range p/m = <_ i0 is the interval commonly used with other approximate algorithms such as the Langevln equation (where p " ~ the Langevin time step). Our formalism also stresses the fact that at fixed p (or acceptance), the effective step size in the update increases when going to small quark mass "m." One should never consider extrapolating to ra - . 0 at fixed p in order to probe QCD in realistic situations. An extrapolation in a should be done first. Figure 3 shows the dependence of various h~dron masses on Mm = (the techniques to extract hadron masses have been presented in Ref. [8]). In the case m = .10, the general trend is a decrease with p/m = in agreement with the fact that configurations U# become more ordered (i.e. W(1,1)-~ unity) with smaller p/m ~ as a result of a more accurate contribution of fermion effects. A similar trend can he noticed in the case m = .05, hut with greater difficulty due to statistical accuracy.

We can get an idea of the typical length scale in our simulation by comparing the data at similar step sizes. In particular, we can use the hadron masses at m --- .05 and .1O with acceptance of 89% and 70% respectively. We find the pion to display the proper Goldstone behavior. Moreover, using the physical value of the mass of the rho, we obtain an inverse lattice spacing of 1327(200) Mev. 4. COMPARISON WITH THE HYBRID ALGORITHM We have run a second set of simulations with the purpose of comparing our results with those of the hybrid algorithm as performed by Gottlieb, Liu, Renken, Sugar and Toussaint [15]. Using the same lattice size (lo s × 24), gauge coupling (fl = 5.4T), quark flavor number (n! = 2, staggered fermions) and quark mass (m = .05), we have run at acceptance 81 and 92%, and calculated the Wilson loops and the hadron mass spectrum (with 81% acc. for the latter only). The differences between those two sets of results illustrate nicely the possible disagreements between the two algorithms when applied to realistic situations (i.e. large volumes, and small quark masses). Let us mention that the masses of Gottlleb et al. shown in Fig. 4 were obtained by analyzing their hadron propagators with our own fitting programs. The hybrid algorithm of Gottlieb et al. is a scheme where finite differences are approximated with a second order approximation, in contrast to Eq. 3, which

J. Potvin et M, / Systematic errors of the pseudofermion algorithm

m a : ~ a = 0.05

mha

ma=~a-0.10

493

1.30.

mh8

B 1,50

1.20

1.40

1.10.

1,30

1.00

~

1,20

It

a1(1++}

0.80.

1.10

} P' '},

foCO++)

1,00 0.90 0.80 0.70

I

I

I

I

I

I

I

I

1

2

3

4

5

8

7

8

L

p/cm,,)2

is first order. Their typical step size is parametrlzed by Ar which has 0.04 here. Figure 4 summarizes the results of our comparison. Their Wilson loop averages to .5418(1) in con-

pc.~

0.70. 0.60-

I

0.40

P

=c.)

I

10

I

20

30

• t)/rrna~ 2

Fig. 3b The hadron mass spectrum versus p/m 2.

lO~ 24

B =5.47

ma =0.05

2.0 1.8

trast to ours, .S3Sl(5) (92% acc.). In genera, the

1.6

hadron masses appear to agree within error; however, our masses appear systematically higher, in accord with the fact that the hybrid algorithm is more exact, having a more ordered gauge configurations.

1,4

P

m~"

!

W(U) t~S~4

---t

1.2 I.O O.S

m~

!

0.6

CONCLUSIONS

/rW[x)

|

0.50-

Fig. 3a The hadron mass spectrum versus p/rn,2.

5.

aIC1+÷)

0.90.

l

mz

l

,.p

I

1*'4

0.4

Tm

0.2

Like in the Langevin or the hybrid methods, the results obtained from the pseudofermion algorithm have to be carefully monitored with respect to the step size. For small enough values of p/m 2, it should be possible to extrapolate to zero - p/ra 2, thus obtaining error free results. On the other hand, a detailed comparison of our hadron masses with those of the hybrid algorithm show general agreement between the two methods.

O,O

+S30 A A O

n

A

B

A

O

S2O

CnrnpoStdnl, Morlmrty, Potvln, and Rel~i Gottlleb, Liu, flenken, Sugar, and Tounaint

Fig. 4 The hadron mass spectrum, as calculated by the pseudofermion algorithm (Campostrini et al.) and the hybrid algorithm (Gottlieb et al.)

494

J. Potvin et al./ Systematic errors of the pseudofermion algorithm

ACKNOWLEDGEMENTS We thank S. Gottlieb and D. Toussalnt for useful discussions and for sending us their hadron propagators before publication. We would also like to thank Lloyd M. Thorndyke and Carl S. Ledbetter of ETA Systems, Inc. and Robert M. Price, Tom Roberts and Gil Williams of Control Data Corporation for their continued interest, support and encouragement. We also gratefully acknowledge financial or computer support from the following source: U.S. Department of Energy (Contract DE-AC02-86ER 40284), the John Von Neumann Super computer Center (Grant Nos. 110128, 171812, 171813, 551701-551705), the Control Data Corporation (Pacer Fellowship 85-88PCR01), ETA systems (Grant Nos. 304 658 and 1312963), the Natural Sciences and Engineering Research Council of Canada (Grant Nos. A8420 and A9030), Dalhousie University, Technical University of Nova Scotia, the Canada/Nova Scotia Technology transfer and Industrial Innovation Agreement, and Energy, Mines and ~.esources~ Government of Canada. REFERENCES [1.] There has been a recent proposal by Kennedy and collaborators (Nucl. Phys. Suppl. B4 (1988) 576) about an algorithm which simulates exactly the effects of the dynamical quarks, at a reasonable cost in computer time. Results of the first QCD applications of this so-called "Hybrid Monte Carlo algorithm" have been presented at this conference. See the contribution by Gupta, Ukawa and Weingarten in these proceedings. [2.] See the contributions by J. Potvin, M. Campostrini K. Moriarty and C. Rebbi, and by M. Fukugita in Proc. of the NATO ARW 641/86 Workshop on Lattice Gauge Theory, Brookhaven National Lab., Sept. 1986; It. Satz, I. Harrity and J. Potvin, eds; Plenum, 1987. Also see E. Laermann, F. Langhammer, I. Schmitt and P.M. Zerwas, Phys. Left. 173B (1986) 437 and 443.

[3.] F. Fucito, E. Marinari, G. Parisi and C. Rebbi Nucl. Phys B180 (1981) 369; H.W. Hamber, E. Marinari, G. Parlsi and C. Rebbi, Phys. Lett. 124B (1983) 99. [4.] Ph. deForcrand and I.O. Stamatescu Nucl. Phys. B261 (1985) 613, and I.O. Stamatescu in these proceedings [5.] G. Batrouni, G.R. Katz, A.S. Kronfeld, G.P. Lepage, B. Svetitsky and K.G. Wilson, Phys. Rev. D32 (1985) 2736; A. Ukawa and M. Fukugita, Phys. Rev. Lett. 55 (1985) 1854. [6.] S. Duane, Nucl. Phys. B257 [FS14] (1985) 652; J. B. Kogut, Nucl. Phys. B270 [FS16] (1986) 169. [7.] D. Callaway and A. Rahman, Phys. Rev. Left. 49 (1982) 613. [8.] M. Campostrini, K.J.M. Moriarty, J. Potvin and C. Rebbi, Phys. Left. 193B (1987) 78. [9.] E. Laermann in these proceedings [1O.] Y. Koike, M. Fukugita and A. Ukawa, Kyoto preprir RIFP-761, July, 1988. [11.] A. Billoire, Ph. de Forcrand and E. Marinari, Nucl. Phys. B270 [FS16] (1986) 33. [12.] M. Campostrini, G. Curci and P. Rossi, Nucl. Phys. B (in press). [13.] H. Hamber, Phys. Lett. 193B (1987) 292. [14.] M.P. Grady, D.K. Sinclair and J.B. Kogut, Phys. Lett. B200 (1988) 149. [15.] S. Gottlieb, W. Liu, R.L. Renken, R.L. Sugar and D. Toussaint~ 1988, Fermilab preprint, Fermilab PUB-88/62-T. [16.] M. Campostrini, K.J.M. Moriarty, J. Potvin and C. Rebbi, Comput. Phys. Commun. 50 (1988) 395. [17.] M. Creutz, L. Jacobs and C. Rebbi, Phys. Rep. 93 (1983) 201. [18.] R.V. Gavai, A. Gocksch and U.M. Heller, Nucl. Phys. B283 (1987) 381. [19.] M. Campostrini, K.J.M. Moriarty, J. Potvin and C. Rebbi, in preparation.