optimization/7563: g77 3.1.1 internal error: Segmentation Fault at -O2 -funroll-loops -m64
David G Hough
david.hough@sun.com
Fri Aug 9 16:56:00 GMT 2002
>Number: 7563
>Category: optimization
>Synopsis: g77 3.1.1 internal error: Segmentation Fault at -O2 -funroll-loops -m64
>Confidential: no
>Severity: serious
>Priority: medium
>Responsible: unassigned
>State: open
>Class: ice-on-legal-code
>Submitter-Id: net
>Arrival-Date: Fri Aug 09 16:56:00 PDT 2002
>Closed-Date:
>Last-Modified:
>Originator: David G Hough
>Release: 3.1.1
>Organization:
Sun Microsystems
>Environment:
System: SunOS dgh 5.8 Generic_108528-15 sun4u sparc SUNW,Ultra-60
Architecture: sun4
host: sparc-sun-solaris2.8
build: sparc-sun-solaris2.8
target: sparc-sun-solaris2.8
configured with: /net/validator/export/set/validator18/hwvalid/compdir/gcc-3.1.1/configure --prefix=/net/validator/export/set/validator18/hwvalid/compdir/g311 --enable-languages=c,f77
>Description:
g77 -c qcd.f -O2 -funroll-loops -m64
qcd.f: In program `qcd2':
qcd.f:290: internal error: Segmentation Fault
>How-To-Repeat:
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C THIS SOFTWARE COPYRIGHTED BY CALTECH, 1983, 1986, 1987
C
C PURE GAUGE SU(3) MONTE CARLO
C THIS CODE HAS A LONG HISTORY:
C EARLY VERSIONS OF THIS DEVELOPED FOR THE MARK I
C HYPERCUBE BY STEVE OTTO AND PAUL STOLORZ (19
C JON FLOWER AND STEVE OTTO EXTENDED AND IMPROVED THE
C CODE FOR THE MARK II MACHINES. (1984-86)
C THE CURRENT INCARNATION IS A "PORT" TO THE NCUBE
C DONE BY CLIVE BAILLIE, ED FELTEN, STEVE OTTO,
C AND DAVID WALKER (1987)
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C PURE GAUGE SU(3) THEORY IN 3+1 DIMENSIONS
C LINK VARIABLE UPDATING DONE VIA CABIBBO-MARINARI HEAT BA
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE
C AND DAVID WALKER, SUMMER 1987
C
PROGRAM qcd2
C
INTEGER dims
PARAMETER (dims=4)
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mleg=18*l1*l2*l3)
INTEGER umax
PARAMETER (umax=72*l1*l2*l3*l4)
PARAMETER (mbar=18*l1*l2*l3*l4)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mleg)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
PARAMETER (lismax=(l1*l2*l3*l4)/2+1)
COMMON /colcom/ lisred(lismax) , lisblk(lismax) , nred , nblack
C
DOUBLEPRECISION u1(umax)
COMMON /comu1 / u1
C
DOUBLEPRECISION ubar(mbar) , sbar(mbar)
COMMON /comubr/ ubar , sbar
C
INTEGER numstg
COMMON /comfil/ numstg
C
INTEGER i , j , k , numoff , avgnum , flg
INTEGER r , t , numbar , dir , rmin , rmax , tmin(13) , tmax(13)
INTEGER r1max , r2max , r1 , r2 , r3 , tmx , tmn
CHARACTER c
INTEGER sweeps , pinval
C
DOUBLEPRECISION fc(2)
c isvalid, unit 8 changed to unit 0, etc. for /valid
CHARACTER*7 isvalid
C
write(6,*) ' BEGIN perfect.qcd2.F,1.12,99/03/12 '
C
call validflush
call startvalidtime
C
CALL prnset(42)
C
c WRITE(8,201)
c201 FORMAT(1X,'OUTPUT FOR PERFECT CLUB BENCHMARK: QCD'/
c * 1X,'$Revision: 1.12 $ $Author: lpointer $'/)
c WRITE(8,202)
c202 FORMAT(1X,'VALIDATION PARAMETERS:'/)
C
WRITE (6,*) 'PURE GAUGE QCD USING PSEUDO HEAT BATH'
WRITE (6,*) '-------------------------------------'
WRITE (6,*) ' '
WRITE (6,*) 'ENTER COMMAND (? FOR HELP)'
WRITE (6,*) ' '
100 CONTINUE
READ (5,99008) c
C
IF ( c.EQ.'R' ) THEN
CALL setsed()
ELSEIF ( c.EQ.'D' ) THEN
CALL dmpsed()
ELSEIF ( c.EQ.'S' ) THEN
CALL setpar()
ELSEIF ( c.EQ.'P' ) THEN
CALL projec()
ELSEIF ( c.EQ.'U' ) THEN
WRITE (6,*)
& 'ENTER NUMBER OF TEMPORAL LINK HITS TO AVERAGE OVER'
READ (5,*) avgnum
WRITE (6,*) avgnum
C CLOCK = SECNDS(0)
CALL linkbr(avgnum,3,ubar)
C WRITE(6,101) SECNDS(0)-CLOCK
WRITE (6,*) 'UBARS MADE.'
ELSEIF ( c.EQ.'M' ) THEN
WRITE (6,*) 'ENTER DIR (0<->X,1<->Y,2<->Z), RMIN, RMAX'
READ (5,*) dir , rmin , rmax
WRITE (6,*) dir , rmin , rmax
DO r = rmin , rmax
WRITE (6,*) 'FOR R = ' , r , ' ENTER TMIN, TMAX '
READ (5,*) tmin(r+1) , tmax(r+1)
WRITE (6,*) tmin(r+1) , tmax(r+1)
ENDDO
C
DO r = 1 , 13
DO t = 1 , 13
which(r,t) = 0
ENDDO
ENDDO
C
DO r = rmin , rmax
DO t = 1 , 13
IF ( (t.GE.tmin(r+1)) .AND. (t.LE.tmax(r+1)) )
& which(r+1,t+1) = 1
ENDDO
ENDDO
C
WRITE (6,*) 'ENTER NUMBER OF SPATIAL LINK HITS TO AVERAGE OVER'
READ (5,*) avgnum
WRITE (6,*) avgnum
WRITE (6,*) 'MEASURING WILSON (ON-AXIS) LOOPS IN DIRECTION ' ,
& dir
C CLOCK = SECNDS(0)
CALL measur(dir,rmin,rmax,avgnum)
C WRITE(6,102) SECNDS(0)-CLOCK
fc(1) = numstg
DO r = rmin , rmax
DO t = tmin(r+1) , tmax(r+1)
fc(2) = wloop(r+1,t+1)
fc(2) = fc(2)/fc(1)
WRITE (6,1) r , t , fc(2)
1 format ('W' , i3 , 'X' , i3 , ':' , e19.9)
c WRITE(8,*) 'W',R,'X',T,':',FC(2)
C
ENDDO
ENDDO
C
ival = 0
IF ( abs((fc(2)-.633259)/.633259).GE.1E-2 ) ival = 1
C
ELSEIF ( c.EQ.'I' ) THEN
numoff = 0
ELSEIF ( c.EQ.'A' ) THEN
WRITE (6,*) 'ENTER R1,R2,R3,TMN,TMX'
READ (5,*) r1 , r2 , r3 , tmn , tmx
WRITE (6,*) r1 , r2 , r3 , tmn , tmx
WRITE (6,*) 'MEASURING BARYON LOOPS ...'
C CLOCK = SECNDS(0)
DO t = tmn , tmx
CALL qqqmea(r1,r2,r3,t)
ENDDO
C WRITE(6,103) SECNDS(0)-CLOCK
C
write(6,11) fc(2),0.633259d0
11 format(' fc(2) ',f10.6,' vs ', f10.6)
IF ( ival.EQ.0 ) THEN
c WRITE(8,203)
99001 FORMAT (//1X,'RESULTS FOR THIS RUN ARE: VALID')
isvalid = 'VALID'
write(6,*) 'qcd2=LG results are VALID'
ELSE
c WRITE(8,204)
isvalid = 'INVALID'
99002 FORMAT (//1X,'RESULTS FOR THIS RUN ARE: INVALID')
write(6,*) 'qcd2=LG results are INVALID'
call validexitbad
ENDIF
C
ELSEIF ( c.EQ.'E' ) THEN
WRITE (6,*) 'ENTER DIR (0<->X,1<->Y,2<->Z)'
READ (5,*) dir
WRITE (6,*) dir
WRITE (6,*) 'ENTER NUMBER OF SPATIAL LINK HITS TO AVERAGE OVER'
READ (5,*) numbar
WRITE (6,*) numbar
DO i = 1 , 13
DO j = 1 , 13
DO k = 1 , 6
res(i,j,k) = 0.0
ENDDO
ENDDO
ENDDO
WRITE (6,*) 'ENTER R1MAX AND R2MAX :'
READ (5,*) r1max , r2max
WRITE (6,*) r1max , r2max
WRITE (6,*) 'NOW ENTER THE QUADRANT TABLE OF TMINS ;'
WRITE (6,*) r1max , 'ROWS &' , r2max , 'COLS (0 ==> OMIT)'
READ (5,*) ((which(i,j),j=1,r2max),i=1,r1max)
WRITE (6,*) 'MEASURING LOOPS IN DIRECTION ' , dir
C
C CLOCK = SECNDS(0)
CALL rotmea(dir,numbar,r1max,r2max)
C WRITE(6,105) SECNDS(0)-CLOCK
DO i = 1 , r1max
DO j = 1 , r2max
IF ( which(i,j).NE.0 ) THEN
tmx = min(12,which(i,j)+5)
DO t = which(i,j) , tmx - 1
WRITE (6,2) i + 1 , j + 1 ,
& t , res(i,j,t-which(i,j)+1)
& /real(numstg)
2 format ('R (' , i3 , ',' , i3 , ') X' ,
& i3 , ':' , e19.9)
c WRITE(8,*) 'R (',I+1,',',J+1,') X',T,':',
c . RES(I,J,T-WHICH(I,J)+1)/DOUBLEPRECISION(NUMSTG)
ENDDO
ENDIF
ENDDO
ENDDO
C
IF ( abs((res(1,1,1-which(1,1)+1)/real(numstg)-.434777)
& /.434777).GE.1E-2 ) ival = 1
IF ( abs((res(1,1,2-which(1,1)+1)/real(numstg)-.299319)
& /.299319).GE.1E-2 ) ival = 1
IF ( abs((res(1,1,3-which(1,1)+1)/real(numstg)-.205511)
& /.205511).GE.1E-2 ) ival = 1
IF ( abs((res(1,1,4-which(1,1)+1)/real(numstg)-.136702)
& /.136702).GE.1E-2 ) ival = 1
IF ( abs((res(1,1,5-which(1,1)+1)/real(numstg)-9.14717E-02)
& /9.14717E-02).GE.1E-2 ) ival = 1
C
ELSEIF ( c.EQ.'C' ) THEN
WRITE (6,*) '0 <-> COLD START, 1 <-> RESTART'
READ (5,*) flg
WRITE (6,*) flg
CALL init(flg)
WRITE (6,*) 'INITIALIZED SUCCESSFULLY.'
ELSEIF ( c.EQ.'B' ) THEN
CALL lattot()
WRITE (6,*) 'BACKED UP SUCCESSFULLY'
ELSEIF ( c.EQ.'F' ) THEN
WRITE (6,*) 'INPUT # OF SWEEPS AND PROJECT INTERVAL'
READ (5,*) sweeps , pinval
WRITE (6,*) sweeps , pinval
WRITE (6,*) 'RUNNING...'
C CLOCK = SECNDS(0)
CALL update(sweeps,pinval)
C WRITE(6,104) SECNDS(0)-CLOCK
WRITE (6,*) 'LATTICE UPDATED'
ELSEIF ( c.EQ.'?' ) THEN
WRITE (6,*) 'COMMANDS (AND THEIR PARAMETERS) ARE:'
WRITE (6,*) '===================================='
WRITE (6,*) 'C : INITIALIZE (COLD OR READ)'
WRITE (6,*) 'D : UPDATES (SWEEPS, PINVAL)'
WRITE (6,*) 'B : BACK UP LATTICE (FILENAME)'
WRITE (6,*) 'S : SET PARAMETERS'
WRITE (6,*) 'R : RESET RANDOM NUMBER SEEDS'
WRITE (6,*) 'D : DUMP CURRENT SEED TO FILE'
WRITE (6,*) 'P : PROJECT MATRICES ONTO SU(3)'
WRITE (6,*) 'M : MEASURE ON AXIS LOOPS'
WRITE (6,*) 'E : MEASURE OFF AXIS LOOP'
WRITE (6,*) 'A : MEASURE BARYON LOOPS'
WRITE (6,*) 'U : MAKE U BARS'
WRITE (6,*) 'Q : QUIT'
ELSEIF ( c.EQ.'Q' ) THEN
WRITE (6,*) 'FINISHED.'
C
call printvalidtime(1.0d0,"perfect.qcd2.F,1.12,99/03/12")
C
c WRITE(8,*)
OPEN (UNIT=1,FILE='SEED',ACCESS='SEQUENTIAL',STATUS='UNKNOWN')
CLOSE (1,STATUS='delete')
STOP
ELSE
WRITE (6,*) 'UNKNOWN COMMAND (? FOR HELP)'
ENDIF
GOTO 100
C
99003 FORMAT (/5X,'Time for LINKBR = ',F11.4,' seconds')
99004 FORMAT (/5X,'Time for MEASUR = ',F11.4,' seconds')
99005 FORMAT (/5X,'Time for QQQMEA = ',F11.4,' seconds')
99006 FORMAT (/5X,'Time for UPDATE = ',F11.4,' seconds')
99007 FORMAT (/5X,'Time for ROTMEA = ',F11.4,' seconds')
99008 FORMAT (A1)
END
**==cfill.f
C
SUBROUTINE cfill(a,b)
DOUBLEPRECISION a(*) , b(*)
C
DO i = 1 , 12
a(i) = b(i)
ENDDO
a(13) = a(3)*a(11) - a(4)*a(12) - a(5)*a(9) + a(6)*a(10)
a(14) = -a(4)*a(11) - a(3)*a(12) + a(6)*a(9) + a(5)*a(10)
a(15) = a(5)*a(7) - a(6)*a(8) - a(1)*a(11) + a(2)*a(12)
a(16) = -a(6)*a(7) - a(5)*a(8) + a(2)*a(11) + a(1)*a(12)
a(17) = a(1)*a(9) - a(2)*a(10) - a(3)*a(7) + a(4)*a(8)
a(18) = -a(2)*a(9) - a(1)*a(10) + a(4)*a(7) + a(3)*a(8)
C
RETURN
END
**==choos.f
C HEAT BATH - USES CABIBBO-MARINARI ALGORITHM INSTEAD OF
C METROPOLIS
C 3 SUBGROUP VERSION
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE
C AND DAVID WALKER, SUMMER 1987
C
SUBROUTINE choos(sst,pu)
DOUBLEPRECISION sst(*)
DOUBLEPRECISION pu(*)
C
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mleg=18*l1*l2*l3)
INTEGER dims
PARAMETER (dims=4)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mleg)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
C
DOUBLEPRECISION r0 , r1 , r2 , r3 , r
DOUBLEPRECISION a0 , a1 , a2 , a3 , xx , efac
DOUBLEPRECISION rfac , norm
DOUBLEPRECISION part(18) , plaq(18) , mat(18)
DOUBLEPRECISION dis , b0 , b1 , b2 , b3
DOUBLEPRECISION temp(18)
INTEGER i , j
C
DO i = 1 , 18
part(i) = 0.0
mat(i) = pu(i)
ENDDO
C
DO j = 1 , 6
DO i = 1 , 18
part(i) = part(i) + sst(i+18*(j-1))
ENDDO
ENDDO
C
DO i = 1 , 3
CALL mult(mat,part,plaq)
IF ( i.EQ.1 ) THEN
r0 = 0.5*(plaq(1)+plaq(9))
r1 = 0.5*(plaq(2)-plaq(10))
r2 = 0.5*(plaq(3)-plaq(7))
r3 = 0.5*(plaq(4)+plaq(8))
ELSEIF ( i.EQ.2 ) THEN
r0 = 0.5*(plaq(9)+plaq(17))
r1 = 0.5*(plaq(10)-plaq(18))
r2 = 0.5*(plaq(11)-plaq(15))
r3 = 0.5*(plaq(12)+plaq(16))
ELSE
r0 = 0.5*(plaq(1)+plaq(17))
r1 = 0.5*(plaq(2)-plaq(18))
r2 = 0.5*(plaq(5)-plaq(13))
r3 = 0.5*(plaq(6)+plaq(14))
ENDIF
r = r0*r0 + r1*r1 + r2*r2 + r3*r3
r = sqrt(r)
r0 = r0/r
r1 = r1/r
r2 = r2/r
r3 = r3/r
xx = beta*r*2.0/3.0
efac = exp(-2.0*xx)
rfac = 1.0 - efac
C
50 CONTINUE
a0 = 1.0 + log(efac+rfac*pranf(0))/xx
xtemp = pranf(0)
IF ( (xtemp*xtemp).GT.(1.0-a0*a0) ) GOTO 50
norm = sqrt(1.0-a0*a0)
C
100 CONTINUE
a1 = 2.0*pranf(0) - 1.0
a2 = 2.0*pranf(0) - 1.0
a3 = 2.0*pranf(0) - 1.0
dis = a1*a1 + a2*a2 + a3*a3
IF ( dis.GT.1.0 ) GOTO 100
C
norm = norm/sqrt(dis)
a1 = a1*norm
a2 = a2*norm
a3 = a3*norm
b0 = a0*r0 + a1*r1 + a2*r2 + a3*r3
b1 = -a0*r1 + a1*r0 + a2*r3 - a3*r2
b2 = a2*r0 + a3*r1 - a0*r2 - a1*r3
b3 = -a2*r1 + a3*r0 - a0*r3 + a1*r2
IF ( i.EQ.1 ) THEN
temp(1) = b0*mat(1) - b1*mat(2) - b2*mat(3) - b3*mat(4)
temp(2) = b1*mat(1) + b0*mat(2) + b3*mat(3) - b2*mat(4)
temp(3) = b2*mat(1) - b3*mat(2) + b0*mat(3) + b1*mat(4)
temp(4) = b3*mat(1) + b2*mat(2) - b1*mat(3) + b0*mat(4)
temp(7) = b0*mat(7) - b1*mat(8) - b2*mat(9) - b3*mat(10)
temp(8) = b1*mat(7) + b0*mat(8) + b3*mat(9) - b2*mat(10)
temp(9) = b2*mat(7) - b3*mat(8) + b0*mat(9) + b1*mat(10)
temp(10) = b3*mat(7) + b2*mat(8) - b1*mat(9) + b0*mat(10)
temp(13) = b0*mat(13) - b1*mat(14) - b2*mat(15) - b3*mat(16)
temp(14) = b1*mat(13) + b0*mat(14) + b3*mat(15) - b2*mat(16)
temp(15) = b2*mat(13) - b3*mat(14) + b0*mat(15) + b1*mat(16)
temp(16) = b3*mat(13) + b2*mat(14) - b1*mat(15) + b0*mat(16)
mat(1) = temp(1)
mat(2) = temp(2)
mat(3) = temp(3)
mat(4) = temp(4)
mat(7) = temp(7)
mat(8) = temp(8)
mat(9) = temp(9)
mat(10) = temp(10)
mat(13) = temp(13)
mat(14) = temp(14)
mat(15) = temp(15)
mat(16) = temp(16)
ELSEIF ( i.EQ.2 ) THEN
temp(3) = b0*mat(3) - b1*mat(4) - b2*mat(5) - b3*mat(6)
temp(4) = b1*mat(3) + b0*mat(4) + b3*mat(5) - b2*mat(6)
temp(5) = b2*mat(3) - b3*mat(4) + b0*mat(5) + b1*mat(6)
temp(6) = b3*mat(3) + b2*mat(4) - b1*mat(5) + b0*mat(6)
temp(9) = b0*mat(9) - b1*mat(10) - b2*mat(11) - b3*mat(12)
temp(10) = b1*mat(9) + b0*mat(10) + b3*mat(11) - b2*mat(12)
temp(11) = b2*mat(9) - b3*mat(10) + b0*mat(11) + b1*mat(12)
temp(12) = b3*mat(9) + b2*mat(10) - b1*mat(11) + b0*mat(12)
temp(15) = b0*mat(15) - b1*mat(16) - b2*mat(17) - b3*mat(18)
temp(16) = b1*mat(15) + b0*mat(16) + b3*mat(17) - b2*mat(18)
temp(17) = b2*mat(15) - b3*mat(16) + b0*mat(17) + b1*mat(18)
temp(18) = b3*mat(15) + b2*mat(16) - b1*mat(17) + b0*mat(18)
mat(3) = temp(3)
mat(4) = temp(4)
mat(5) = temp(5)
mat(6) = temp(6)
mat(9) = temp(9)
mat(10) = temp(10)
mat(11) = temp(11)
mat(12) = temp(12)
mat(15) = temp(15)
mat(16) = temp(16)
mat(17) = temp(17)
mat(18) = temp(18)
ELSE
temp(1) = b0*mat(1) - b1*mat(2) - b2*mat(5) - b3*mat(6)
temp(2) = b1*mat(1) + b0*mat(2) + b3*mat(5) - b2*mat(6)
temp(5) = b2*mat(1) - b3*mat(2) + b0*mat(5) + b1*mat(6)
temp(6) = b3*mat(1) + b2*mat(2) - b1*mat(5) + b0*mat(6)
temp(7) = b0*mat(7) - b1*mat(8) - b2*mat(11) - b3*mat(12)
temp(8) = b1*mat(7) + b0*mat(8) + b3*mat(11) - b2*mat(12)
temp(11) = b2*mat(7) - b3*mat(8) + b0*mat(11) + b1*mat(12)
temp(12) = b3*mat(7) + b2*mat(8) - b1*mat(11) + b0*mat(12)
mat(1) = temp(1)
mat(2) = temp(2)
mat(5) = temp(5)
mat(6) = temp(6)
mat(7) = temp(7)
mat(8) = temp(8)
mat(11) = temp(11)
mat(12) = temp(12)
ENDIF
ENDDO
C
CALL cfill(pu,mat)
C
RETURN
END
**==cpymat.f
SUBROUTINE cpymat(outmat,inmat,n)
DOUBLEPRECISION inmat(*) , outmat(*)
INTEGER n
C
DO i = 1 , n
outmat(i) = inmat(i)
ENDDO
C
RETURN
END
**==icpymt.f
SUBROUTINE icpymt(outmat,inmat,n)
INTEGER inmat(*) , outmat(*) , n
C ADDED BY LPOINTER 4-3-89
DO i = 1 , n
outmat(i) = inmat(i)
ENDDO
RETURN
END
**==dmpsed.f
SUBROUTINE dmpsed()
C
OPEN (UNIT=1,FILE='SEED',ACCESS='SEQUENTIAL',STATUS='UNKNOWN')
REWIND (1)
ino = igetsd(0)
WRITE (1,99001) ino
99001 FORMAT (I12)
WRITE (6,*) 'O.K. ...... DONE IT !'
CLOSE (1)
END
**==getpos.f
C
CCCCCCCCCCCCCCCCCCCCCCCCC END OF ROUTINE DMPSED CCCCCCCCCCCCCCCCCCCC
C
SUBROUTINE getpos(latt,coord,index)
INTEGER coord(4) , latt(4) , index
C
j = index
DO i = 1 , 4
coord(i) = mod(j,latt(i))
j = j/latt(i)
ENDDO
C
RETURN
END
**==hit.f
SUBROUTINE hit(coord,origin,dir)
INTEGER coord(4) , dir , origin
C
INTEGER dims
PARAMETER (dims=4)
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mbar=18*l1*l2*l3)
INTEGER umax
PARAMETER (umax=72*l1*l2*l3*l4)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mbar)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
C
DOUBLEPRECISION u1(umax)
COMMON /comu1 / u1
C
INTEGER othend , site , cnt
INTEGER pu
DOUBLEPRECISION st(18,6)
INTEGER ptemp(24)
C
INTEGER pstr(24,4)
SAVE pstr
C
DATA pstr/2 , 5 , 6 , 14 , 6 , 5 , 2 , 14 , 3 , 5 , 7 , 14 , 7 ,
& 5 , 3 , 14 , 4 , 5 , 8 , 14 , 8 , 5 , 4 , 14 , 1 , 6 , 5 ,
& 14 , 5 , 6 , 1 , 14 , 3 , 6 , 7 , 14 , 7 , 6 , 3 , 14 , 4 ,
& 6 , 8 , 14 , 8 , 6 , 4 , 14 , 1 , 7 , 5 , 14 , 5 , 7 , 1 ,
& 14 , 2 , 7 , 6 , 14 , 6 , 7 , 2 , 14 , 4 , 7 , 8 , 14 , 8 ,
& 7 , 4 , 14 , 1 , 8 , 5 , 14 , 5 , 8 , 1 , 14 , 2 , 8 , 6 ,
& 14 , 6 , 8 , 2 , 14 , 3 , 8 , 7 , 14 , 7 , 8 , 3 , 14/
C
C CALL CPYMAT(PTEMP,PSTR(1,DIR+1),24) LPOINTER 4-3-89
CALL icpymt(ptemp,pstr(1,dir+1),24)
C
IF ( coord(dir+1).EQ.latt1(dir+1) ) THEN
othend = 0
ELSE
othend = coord(dir+1) + 1
ENDIF
C
coord(dir+1) = othend
site = coord(1)*mov(1) + coord(2)*mov(2) + coord(3)*mov(3)
& + coord(4)*mov(4)
cnt = 1
DO j = 1 , 6
coord(dir+1) = othend
CALL syslop(coord,site,ptemp(cnt),st(1,j),nchar)
cnt = nchar + cnt
ENDDO
C
pu = origin + rot(dir+1)
CALL choos(st(1,1),u1(pu+1))
C
RETURN
END
**==init.f
SUBROUTINE init(flag)
C
INTEGER numstg
COMMON /comfil/ numstg
C
INTEGER dims
PARAMETER (dims=4)
INTEGER umax
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mleg=18*l1*l2*l3,umax=72*l1*l2*l3*l4)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mleg)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
C
DOUBLEPRECISION u1(umax)
COMMON /comu1 / u1
C
DOUBLEPRECISION facfi , facif
COMMON /comfac/ facfi , facif
C
INTEGER flag , bignum , movx , buff , rem , bufnum
INTEGER pu
CHARACTER*5 bname
C
movx = 72
buff = 256
bignum = 4*latt(1)*latt(2)*latt(3)*latt(4)
C
IF ( flag.EQ.1 ) THEN
WRITE (6,*) 'BACKUP FILE NAME ? : '
READ (5,99001) bname
OPEN (UNIT=1,FILE=bname,ACCESS='SEQUENTIAL',STATUS='UNKNOWN',
& FORM='UNFORMATTED')
REWIND (1)
bignum = movx*numstg
bufnum = bignum/buff
rem = mod(bignum,buff)
WRITE (6,*) 'BACKUP FILE NAME = ' , bname , 'BUFNUM = ' ,
& bufnum , 'REM = ' , rem
pu = 1
DO i = 0 , bufnum - 1
READ (1) (u1(j),j=pu,pu+buff)
pu = pu + buff
WRITE (6,*) 'B'
ENDDO
WRITE (6,*) 'E'
CLOSE (1)
ELSE
pu = 1
DO i = 0 , bignum - 1
u1(pu) = 1.0
u1(pu+1) = 0.0
u1(pu+2) = 0.0
u1(pu+3) = 0.0
u1(pu+4) = 0.0
u1(pu+5) = 0.0
u1(pu+6) = 0.0
u1(pu+7) = 0.0
u1(pu+8) = 1.0
u1(pu+9) = 0.0
u1(pu+10) = 0.0
u1(pu+11) = 0.0
u1(pu+12) = 0.0
u1(pu+13) = 0.0
u1(pu+14) = 0.0
u1(pu+15) = 0.0
u1(pu+16) = 1.0
u1(pu+17) = 0.0
pu = pu + rot(2)
ENDDO
ENDIF
C
RETURN
99001 FORMAT (A)
END
**==lattot.f
C
C THE ROUTINES USED TO INITIALIZE AND DUMP THE LATTICE.
C
C CONVERTED FORM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, SUMMER 1987
C
SUBROUTINE lattot()
C
INTEGER umax
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (umax=72*l1*l2*l3*l4)
DOUBLEPRECISION u1(umax)
COMMON /comu1 / u1
C
INTEGER numstg
COMMON /comfil/ numstg
C
INTEGER bignum
INTEGER bufnum , rem
INTEGER pu
INTEGER buff , movx
CHARACTER*5 bname
C
WRITE (6,*) 'ENTER 5-LETTER BACKUP FILE NAME ? : '
READ (5,99002) bname
WRITE (6,*) bname
OPEN (UNIT=1,FILE=bname,ACCESS='SEQUENTIAL',STATUS='UNKNOWN')
C OPEN(UNIT=1,FILE = BNAME,ACCESS = 'SEQUENTIAL',STATUS='UNKNOWN',
C .FORM='UNFORMATTED')
REWIND (1)
movx = 72
buff = 256
bignum = movx*numstg
bufnum = bignum/buff
rem = mod(bignum,buff)
WRITE (6,*) 'BACKUP FILE NAME = ' , bname , 'BUFNUM = ' , bufnum ,
& 'REM = '
pu = 1
DO i = 0 , bufnum - 1
C WRITE (1) (U1(J),J=PU,PU+BUFF)
WRITE (1,99001) (u1(j),j=pu,pu+buff)
99001 FORMAT (9F8.3)
pu = pu + buff
WRITE (6,*) 'B'
ENDDO
WRITE (6,*) 'E'
CLOSE (1)
RETURN
99002 FORMAT (A)
END
**==legmak.f
C
C THIS ROUTINE MAKES BARRED LEGS ONLY IN THE ARRAY LEGS
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, 1987
C
SUBROUTINE legmak(t,tt)
INTEGER tt , t
C
INTEGER dims
PARAMETER (dims=4)
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mleg=18*l1*l2*l3)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mleg)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
C
INTEGER strb(26)
C
INTEGER site , coord(4) , plegs
C
DO i = 1 , t
strb(i) = 18
ENDDO
C
DO i = t + 1 , 2*t + 1
strb(i) = 30
ENDDO
C
strb(2*t+1) = 14
C
coord(4) = tt
site = tt*mov(4)
plegs = 1
DO i = 0 , latt(3) - 1
coord(3) = i
DO j = 0 , latt(2) - 1
coord(2) = j
DO k = 0 , latt(1) - 1
coord(1) = k
CALL syslop(coord,site,strb,legs(plegs),nchar)
plegs = plegs + 18
site = site + mov(1)
ENDDO
ENDDO
ENDDO
C
RETURN
END
**==linkbr.f
C LINKBR MULTI-HITS AND AVERAGES U MATRICES IN ANY DIRECTION
C AND PUTS THE RESULTS IN THE ARRAY POINTED TO BY "WHERE".
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, 1987
C
SUBROUTINE linkbr(n,dir,where)
DOUBLEPRECISION where(*)
INTEGER n , dir
C
INTEGER dims
PARAMETER (dims=4)
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mleg=18*l1*l2*l3)
INTEGER umax
PARAMETER (umax=72*l1*l2*l3*l4)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mleg)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
C
DOUBLEPRECISION u1(umax)
COMMON /comu1 / u1
C
INTEGER coord(4) , othend , site , endsit , pt , cnt
INTEGER pu
INTEGER pstr(24)
DOUBLEPRECISION avmat(18) , mat2(18) , st(18,6)
C
INTEGER env(24,4)
SAVE env
C
DATA env/2 , 5 , 6 , 14 , 6 , 5 , 2 , 14 , 3 , 5 , 7 , 14 , 7 ,
& 5 , 3 , 14 , 4 , 5 , 8 , 14 , 8 , 5 , 4 , 14 , 1 , 6 , 5 ,
& 14 , 5 , 6 , 1 , 14 , 3 , 6 , 7 , 14 , 7 , 6 , 3 , 14 , 4 ,
& 6 , 8 , 14 , 8 , 6 , 4 , 14 , 1 , 7 , 5 , 14 , 5 , 7 , 1 ,
& 14 , 2 , 7 , 6 , 14 , 6 , 7 , 2 , 14 , 4 , 7 , 8 , 14 , 8 ,
& 7 , 4 , 14 , 1 , 8 , 5 , 14 , 5 , 8 , 1 , 14 , 2 , 8 , 6 ,
& 14 , 6 , 8 , 2 , 14 , 3 , 8 , 7 , 14 , 7 , 8 , 3 , 14/
C
pt = 1
DO icord4 = 0 , latt(4) - 1
coord(4) = icord4
DO icord3 = 0 , latt(3) - 1
coord(3) = icord3
DO icord2 = 0 , latt(2) - 1
coord(2) = icord2
DO icord1 = 0 , latt(1) - 1
coord(1) = icord1
cnt = 1
site = icord1*mov(1) + icord2*mov(2) + icord3*mov(3)
& + icord4*mov(4)
pu = site + rot(dir+1)
CALL cpymat(avmat,u1(pu+1),18)
IF ( n.GT.1 ) THEN
DO j = 1 , 18
mat2(j) = u1(pu+j)
ENDDO
IF ( coord(dir+1).EQ.latt1(dir+1) ) THEN
othend = 0
endsit = site - latt1(dir+1)*mov(dir+1)
ELSE
othend = coord(dir+1) + 1
endsit = site + mov(dir+1)
ENDIF
C CALL CPYMAT(PSTR,ENV(1,DIR+1),24) LPOINTER 4-3-89
CALL icpymt(pstr,env(1,dir+1),24)
C
DO j = 1 , 6
coord(dir+1) = othend
CALL syslop(coord,endsit,pstr(cnt),st(1,j),
& nchar)
cnt = cnt + nchar
ENDDO
DO j = 1 , n - 1
CALL choos(st(1,1),mat2)
DO i = 1 , 18
avmat(i) = avmat(i) + mat2(i)
ENDDO
ENDDO
DO i = 1 , 18
avmat(i) = avmat(i)/n
ENDDO
ENDIF
DO k = 1 , 18
where(pt+k-1) = avmat(k)
ENDDO
pt = pt + 18
ENDDO
ENDDO
ENDDO
ENDDO
C
RETURN
END
**==lpstrm.f
C
C S. W. OTTO AND J. W. FLOWER, SUMMER 1985.
C
C CONSTRUCTS THE STRINGS FOR MEASURING PLANAR RxT WILSON LOOPS.
C IF R=1 THEN INSTEAD OF FETCHING A BARRED LEG FOR THE T DIRECTION W
C CONSTRUCT IT EXPLICITLY USING THE REGULAR LINKS.
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, 1987
C
SUBROUTINE lpstrm(pstrt,pstru,pstrb,r,t,dir)
INTEGER pstrt(60) , pstru(60) , pstrb(60)
INTEGER r , t , dir
C
INTEGER rr , tt , cnt , temp
C
C FIRST GET THE TOP OF THE LOOP BARRING AS MANY LINKS AS POSSIBLE
C
cnt = 1
DO tt = 0 , t - 1
pstrt(cnt) = 26
cnt = cnt + 1
ENDDO
C
temp = 23 + dir
DO rr = 0 , r - 1
pstrt(cnt) = temp
cnt = cnt + 1
ENDDO
C
temp = 5 + dir
pstrt(cnt) = temp
cnt = cnt + 1
C
temp = 19 + dir
DO rr = 1 , r - 2
pstrt(cnt) = temp
cnt = cnt + 1
ENDDO
C
IF ( r.GE.2 ) THEN
temp = 5 + dir
pstrt(cnt) = temp
cnt = cnt + 1
ENDIF
C
DO tt = 0 , t - 1
pstrt(cnt) = 30
cnt = cnt + 1
ENDDO
pstrt(cnt) = 14
C
C NOW CONSTRUCT THE TIMELIKE LEG. IF R=1 MAKE EXPLICITLY
C
cnt = 1
temp = 23 + dir
DO rr = 0 , r - 1
pstru(cnt) = temp
cnt = cnt + 1
ENDDO
C
IF ( r.GT.1 ) THEN
pstru(cnt) = 36
cnt = cnt + 1
ELSE
DO tt = 0 , t - 1
pstru(cnt) = 4
cnt = cnt + 1
ENDDO
DO tt = 0 , t - 1
pstru(cnt) = 30
cnt = cnt + 1
ENDDO
ENDIF
C
temp = 27 + dir
DO rr = 0 , r - 1
pstru(cnt) = temp
cnt = cnt + 1
ENDDO
C
pstru(cnt) = 14
cnt = cnt + 1
C
C FINALLY MAKE THE BOTTOM OF THE LOOP
C
cnt = 1
temp = 1 + dir
pstrb(cnt) = temp
C
cnt = cnt + 1
temp = 15 + dir
DO rr = 1 , r - 2
pstrb(cnt) = temp
cnt = cnt + 1
ENDDO
C
IF ( r.GE.2 ) THEN
temp = 1 + dir
pstrb(cnt) = temp
cnt = cnt + 1
ENDIF
C
temp = 27 + dir
DO rr = 0 , r - 1
pstrb(cnt) = temp
cnt = cnt + 1
ENDDO
C
pstrb(cnt) = 14
C
RETURN
END
**==make.f
C
C STEVE OTTO (AND JON FLOWER) JULY 1984
C
C SETS UP INITIAL MOVEMENT TABLES.
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, SUMMER 1987
C
SUBROUTINE make
INTEGER dims
PARAMETER (dims=4)
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mleg=18*l1*l2*l3)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mleg)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
C
DO i = 1 , 4
rot(i) = 18*(i-1)
ENDDO
C
mov(1) = 72
mov(2) = mov(1)*latt(1)
mov(3) = mov(2)*latt(2)
mov(4) = mov(3)*latt(3)
C
DO i = 1 , 4
latt1(i) = latt(i) - 1
ENDDO
C
RETURN
END
**==measur.f
C
C JON FLOWER, STEVE OTTO MARCH 1985. MOD: JUNE 85
C
C MEASURES WILSON LOOPS UP TO RxT BARRING EVERYTHING IN SIGHT....
C ALSO MAKES LEGS.
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, 1987
C
SUBROUTINE measur(dir,rmin,rmax,numhit)
INTEGER dir , rmin , rmax , numhit
C
INTEGER dims
PARAMETER (dims=4)
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mleg=18*l1*l2*l3)
INTEGER umax
PARAMETER (umax=72*l1*l2*l3*l4)
PARAMETER (mbar=18*l1*l2*l3*l4)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mleg)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
C
DOUBLEPRECISION ubar(mbar) , sbar(mbar)
COMMON /comubr/ ubar , sbar
C
INTEGER r , nchar , coord(4) , tflag , site , t , plegs
INTEGER strt(60) , stru(60) , strb(60)
DOUBLEPRECISION up(18) , down(18) , top(18)
DOUBLEPRECISION bottom(18) , stemp(18,2)
C
DO r = 0 , 12
DO t = 0 , 12
wloop(r+1,t+1) = 0.0
ENDDO
ENDDO
C
CALL linkbr(numhit,dir,sbar)
C
C ONE TIME SLICE AT A TIME
C
DO icord4 = 0 , latt(4) - 1
coord(4) = icord4
DO t = 1 , 12
tflag = 0
DO r = rmin , rmax
IF ( which(r+1,t+1).EQ.1 ) tflag = 1
ENDDO
C
C TFLAG=1 MEANS WE NEED LEGS OF THIS T
C
IF ( tflag.EQ.1 ) CALL legmak(t,icord4)
C
DO r = rmin , rmax
IF ( (r.LE.t) .AND. (which(r+1,t+1).EQ.1) ) THEN
CALL lpstrm(strt,stru,strb,r,t,dir)
site = icord4*mov(4)
plegs = 1
DO icord3 = 0 , latt(3) - 1
coord(3) = icord3
DO icord2 = 0 , latt(2) - 1
coord(2) = icord2
DO icord1 = 0 , latt(1) - 1
coord(1) = icord1
C
C FETCH AND DAGGER THE DOWN LEG
C
CALL udag(legs(plegs),down)
C
C FIND THE TOP
C
CALL syslop(coord,site,strt,top,nchar)
C
C FIND THE UP
C
CALL syslop(coord,site,stru,up,nchar)
C
C FIND THE BOTTOM
C
CALL syslop(coord,site,strb,bottom,nchar)
C
C PUT THEM ALL TOGETHER
C
CALL mult(bottom,up,stemp(1,1))
CALL mult(stemp(1,1),top,stemp(1,2))
CALL mult(stemp(1,2),down,stemp(1,1))
wloop(r+1,t+1) = wloop(r+1,t+1)
& + (stemp(1,1)+stemp(9,1)+stemp(17,1))
& /3.0
plegs = plegs + 18
site = site + mov(1)
ENDDO
ENDDO
ENDDO
ENDIF
ENDDO
ENDDO
ENDDO
C
RETURN
END
**==mult.f
SUBROUTINE mult(a,b,c)
DOUBLEPRECISION a(18) , b(18) , c(18)
C
c(1) = a(1)*b(1) - a(2)*b(2) + a(7)*b(3) - a(8)*b(4) + a(13)*b(5)
& - a(14)*b(6)
c(2) = a(1)*b(2) + a(2)*b(1) + a(7)*b(4) + a(8)*b(3) + a(13)*b(6)
& + a(14)*b(5)
c(3) = a(3)*b(1) - a(4)*b(2) + a(9)*b(3) - a(10)*b(4) + a(15)*b(5)
& - a(16)*b(6)
c(4) = a(3)*b(2) + a(4)*b(1) + a(9)*b(4) + a(10)*b(3) + a(15)*b(6)
& + a(16)*b(5)
c(5) = a(5)*b(1) - a(6)*b(2) + a(11)*b(3) - a(12)*b(4) + a(17)
& *b(5) - a(18)*b(6)
c(6) = a(5)*b(2) + a(6)*b(1) + a(11)*b(4) + a(12)*b(3) + a(17)
& *b(6) + a(18)*b(5)
c(7) = a(1)*b(7) - a(2)*b(8) + a(7)*b(9) - a(8)*b(10) + a(13)
& *b(11) - a(14)*b(12)
c(8) = a(1)*b(8) + a(2)*b(7) + a(7)*b(10) + a(8)*b(9) + a(13)
& *b(12) + a(14)*b(11)
c(9) = a(3)*b(7) - a(4)*b(8) + a(9)*b(9) - a(10)*b(10) + a(15)
& *b(11) - a(16)*b(12)
c(10) = a(3)*b(8) + a(4)*b(7) + a(9)*b(10) + a(10)*b(9) + a(15)
& *b(12) + a(16)*b(11)
c(11) = a(5)*b(7) - a(6)*b(8) + a(11)*b(9) - a(12)*b(10) + a(17)
& *b(11) - a(18)*b(12)
c(12) = a(5)*b(8) + a(6)*b(7) + a(11)*b(10) + a(12)*b(9) + a(17)
& *b(12) + a(18)*b(11)
c(13) = a(1)*b(13) - a(2)*b(14) + a(7)*b(15) - a(8)*b(16) + a(13)
& *b(17) - a(14)*b(18)
c(14) = a(1)*b(14) + a(2)*b(13) + a(7)*b(16) + a(8)*b(15) + a(13)
& *b(18) + a(14)*b(17)
c(15) = a(3)*b(13) - a(4)*b(14) + a(9)*b(15) - a(10)*b(16) + a(15)
& *b(17) - a(16)*b(18)
c(16) = a(3)*b(14) + a(4)*b(13) + a(9)*b(16) + a(10)*b(15) + a(15)
& *b(18) + a(16)*b(17)
c(17) = a(5)*b(13) - a(6)*b(14) + a(11)*b(15) - a(12)*b(16)
& + a(17)*b(17) - a(18)*b(18)
c(18) = a(5)*b(14) + a(6)*b(13) + a(11)*b(16) + a(12)*b(15)
& + a(17)*b(18) + a(18)*b(17)
C
RETURN
END
**==prnse2.f
C
SUBROUTINE prnse2(iseed)
PARAMETER (mult=1103515245,iadd=12345)
INTEGER aaa , bbb , randx
COMMON /compra/ aaa , bbb , randx
SAVE /compra/
aaa = mult
bbb = iadd
randx = iseed
RETURN
END
**==observ.f
C
C
C J. W. FLOWER AUG 1985
C
C CONSTRUCTS THE THREE QUARK OBSERVABLE :
C
C (1) (2) (3)
C EPSILON EPSILON U U U
C IJK PQR IP JQ KR
C
C THIS CODE IS CERTAINLY NON-OPTIMAL BUT SACRIFICES
C SPEED FOR CODING SIMPLICITY
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, 1987
C
SUBROUTINE observ(tot,u1,u2,u3)
DOUBLEPRECISION tot(2) , u1(2,9) , u2(2,9)
DOUBLEPRECISION u3(2,9)
C
INTEGER p , q , r , epsilo(3,3,3)
DOUBLEPRECISION fac
C
DO i = 1 , 3
DO j = 1 , 3
DO k = 1 , 3
epsilo(i,j,k) = 0
ENDDO
ENDDO
ENDDO
C
epsilo(1,2,3) = 1
epsilo(2,3,1) = 1
epsilo(3,1,2) = 1
epsilo(2,1,3) = -1
epsilo(1,3,2) = -1
epsilo(3,2,1) = -1
C
DO i = 0 , 2
DO p = 0 , 2
DO j = 0 , 2
DO q = 0 , 2
DO k = 0 , 2
IF ( epsilo(i+1,j+1,k+1).NE.0 ) THEN
DO r = 0 , 2
IF ( epsilo(p+1,q+1,r+1).NE.0 ) THEN
fac = epsilo(i+1,j+1,k+1)
& *epsilo(p+1,q+1,r+1)
tot(1) = tot(1) + fac*u1(1,3*i+p+1)
& *u2(1,3*j+q+1)*u3(1,3*k+r+1)
tot(1) = tot(1) - fac*u1(2,3*i+p+1)
& *u2(2,3*j+q+1)*u3(1,3*k+r+1)
tot(1) = tot(1) - fac*u1(1,3*i+p+1)
& *u2(2,3*j+q+1)*u3(2,3*k+r+1)
tot(1) = tot(1) - fac*u1(2,3*i+p+1)
& *u2(1,3*j+q+1)*u3(2,3*k+r+1)
tot(2) = tot(2) + fac*u1(1,3*i+p+1)
& *u2(1,3*j+q+1)*u3(2,3*k+r+1)
tot(2) = tot(2) + fac*u1(1,3*i+p+1)
& *u2(2,3*j+q+1)*u3(1,3*k+r+1)
tot(2) = tot(2) + fac*u1(2,3*i+p+1)
& *u2(1,3*j+q+1)*u3(1,3*k+r+1)
tot(2) = tot(2) - fac*u1(2,3*i+p+1)
& *u2(2,3*j+q+1)*u3(2,3*k+r+1)
ENDIF
ENDDO
ENDIF
ENDDO
ENDDO
ENDDO
ENDDO
ENDDO
C
RETURN
END
**==prnset.f
C
SUBROUTINE prnset(iseed)
COMMON /ran / irand0 , irand1 , iab0 , iab1 , icb0 , icb1
SAVE /ran/
C
PARAMETER (mult0=20077,mult1=16838,iadd0=12345,iadd1=0)
PARAMETER (ipow16=2**16)
C
iseed0 = mod(iseed,ipow16)
iseed1 = iseed/ipow16
iab0 = mult0
iab1 = mult1
icb0 = iadd0
icb1 = iadd1
CALL lmult(j0,j1,iseed0,iseed1,iab0,iab1)
CALL ladd(irand0,irand1,j0,j1,icb0,icb1)
RETURN
END
**==pranf.f
C
CCCCCCCCCCCCCCCCCCCCC END OF ROUTINE PRNSET CCCCCCCCCCCCCCCCCCCCCCC C
C
FUNCTION pranf(idummy)
COMMON /ran / irand0 , irand1 , iab0 , iab1 , icb0 , icb1
SAVE /ran/
C
DOUBLE PRECISION randb , divfac , pow16
PARAMETER (divfac=2.147483648D9,pow16=6.5536D4)
C
randb = (dble(irand0)+pow16*dble(irand1))/divfac
pranf = real(randb)
CALL lmult(j0,j1,irand0,irand1,iab0,iab1)
CALL ladd(irand0,irand1,j0,j1,icb0,icb1)
C
RETURN
END
**==ladd.f
C
CCCCCCCCCCCCCCCCCCCCC END OF ROUTINE PRANF CCCCCCCCCCCCCCCCCCCCCCC C
C
SUBROUTINE ladd(i0,i1,j0,j1,k0,k1)
C
PARAMETER (ipow16=2**16,ipow15=2**15)
C
i0 = mod(k0+j0,ipow16)
i1 = mod((k0+j0)/ipow16+k1+j1,ipow15)
C
RETURN
END
**==lmult.f
C
CCCCCCCCCCCCCCCCCCCCC END OF ROUTINE LADD CCCCCCCCCCCCCCCCCCCCCCC C
C
SUBROUTINE lmult(i0,i1,j0,j1,k0,k1)
C
PARAMETER (ipow30=2**30,ipow16=2**16,ipow15=2**15)
C
k = k0*j0
IF ( k.LT.0 ) k = (k+ipow30) + ipow30
i0 = mod(k,ipow16)
l = k0*j1 + k1*j0
IF ( l.LT.0 ) l = (l+ipow30) + ipow30
k = k/ipow16 + l
IF ( k.LT.0 ) k = (k+ipow30) + ipow30
i1 = mod(k,ipow15)
C
RETURN
END
**==igetsd.f
C
CCCCCCCCCCCCCCCCCCCCCCCCCCC END OF ROUTINE LMULT CCCCCCCCCCCCCCCCCCCCCCC
C
FUNCTION igetsd(idummy)
C
COMMON /ran / irand0 , irand1 , iab0 , iab1 , icb0 , icb1
SAVE /ran/
C
igetsd = irand0 + irand1*2**16
C
RETURN
END
**==projec.f
C
CCCCCCCCCCCCCCCCCCCCCCCCCCC END OF ROUTINE IGETSD CCCCCCCCCCCCCCCCCCCC
C
C
C STEVE OTTO AND JON FLOWER, JULY 1984.
C
C DUE TO THE STRANGE SHAPE OF THE SU(3) GROUP MANIFOLD (WHATEVER IT
C THE LINK MATRICES TEND TO DRIFT OFF AND MUST BE PROJECTED BACK
C PERIODICALLY. THE SITUATION IS OBVIOUSLY NOT HELPED BY THE LIMITED
C ACCURACY OF THE FIXED POINT FORMAT. */
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, 1987
C
SUBROUTINE projec
C
INTEGER dims
PARAMETER (dims=4)
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mleg=18*l1*l2*l3)
INTEGER umax
PARAMETER (umax=72*l1*l2*l3*l4)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mleg)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
C
DOUBLEPRECISION u1(umax)
COMMON /comu1 / u1
C
DOUBLEPRECISION x(12) , z(6) , fac
INTEGER site
INTEGER pu
C
numlin = 4*latt(1)*latt(2)*latt(3)*latt(4)
C
DO i = 0 , numlin - 1
site = rot(2)*i
pu = site
DO j = 1 , 12
x(j) = u1(pu+j)
ENDDO
z(1) = x(3)*x(11) - x(4)*x(12) - x(5)*x(9) + x(6)*x(10)
z(2) = -x(4)*x(11) - x(3)*x(12) + x(6)*x(9) + x(5)*x(10)
z(3) = x(5)*x(7) - x(6)*x(8) - x(1)*x(11) + x(2)*x(12)
z(4) = -x(6)*x(7) - x(5)*x(8) + x(2)*x(11) + x(1)*x(12)
z(5) = x(1)*x(9) - x(2)*x(10) - x(3)*x(7) + x(4)*x(8)
z(6) = -x(2)*x(9) - x(1)*x(10) + x(4)*x(7) + x(3)*x(8)
fac = 0.0
DO j = 1 , 6
fac = fac + x(j)*x(j)
ENDDO
fac = 1.0/sqrt(fac)
DO j = 1 , 6
x(j) = x(j)*fac
ENDDO
x(7) = z(3)*x(5) - z(4)*x(6) - z(5)*x(3) + z(6)*x(4)
x(8) = -z(3)*x(6) - z(4)*x(5) + z(5)*x(4) + z(6)*x(3)
x(9) = z(5)*x(1) - z(6)*x(2) - z(1)*x(5) + z(2)*x(6)
x(10) = -z(5)*x(2) - z(6)*x(1) + z(1)*x(6) + z(2)*x(5)
x(11) = z(1)*x(3) - z(2)*x(4) - z(3)*x(1) + z(4)*x(2)
x(12) = -z(1)*x(4) - z(2)*x(3) + z(3)*x(2) + z(4)*x(1)
fac = 0.0
DO j = 1 , 6
fac = fac + x(6+j)*x(6+j)
ENDDO
fac = 1.0/sqrt(fac)
DO j = 1 , 6
x(6+j) = x(6+j)*fac
ENDDO
CALL cfill(u1(pu+1),x)
ENDDO
C
RETURN
END
**==qqqlpc.f
C
C J. W. FLOWER & STEVE OTTO, SEPTEMBER 1984.
C
C USED TO ACTUALLY CALCULATE THE THREE QUARK OBSERVABLES. */
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, 1987
C
SUBROUTINE qqqlpc(r1,r2,r3,t,lpstr,tot)
INTEGER r1 , r2 , r3 , t
INTEGER lpstr(512)
DOUBLEPRECISION tot(2)
C
INTEGER dims
PARAMETER (dims=4)
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mleg=18*l1*l2*l3)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mleg)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
C
INTEGER coord(4) , odir , site , cnt
INTEGER pstr(512)
DOUBLEPRECISION u1(18) , u2(18) , u3(18)
DOUBLEPRECISION bits(18,3) , temp(18)
C
site = t*mov(4)
C
coord(4) = t
DO icord3 = 0 , latt(3) - 1
coord(3) = icord3
DO icord2 = 0 , latt(2) - 1
coord(2) = icord2
DO icord1 = 0 , latt(1) - 1
coord(1) = icord1
C
C FIRST GET THE MATRIX CORRESPONDING TO DIRECTION R1.
C
C CALL CPYMAT(PSTR,LPSTR,512) LPOINTER 4-3-89
CALL icpymt(pstr,lpstr,512)
cnt = 1
DO i = 1 , 3
CALL syslop(coord,site,pstr(cnt),bits(1,i),nchar)
cnt = nchar + cnt
ENDDO
IF ( r1.GT.0 ) THEN
CALL mult(bits(1,1),bits(1,2),temp)
CALL mult(temp,bits(1,3),u1)
ELSE
DO i = 1 , 18
u1(i) = bits(i,2)
ENDDO
ENDIF
C
C NEXT GET THE MATRIX CORRESPONDING TO DIRECTION R2.
C
DO i = 1 , 3
CALL syslop(coord,site,pstr(cnt),bits(1,i),nchar)
cnt = cnt + nchar
ENDDO
IF ( r2.GT.0 ) THEN
CALL mult(bits(1,1),bits(1,2),temp)
CALL mult(temp,bits(1,3),u2)
ELSE
DO i = 1 , 18
u2(i) = bits(i,2)
ENDDO
ENDIF
C
C NOW LOOP OVER THE TWO (CURRENTLY) DIRECTIONS FOR THE R3 LEG
C CALCULATING THE OBSERVABLES AS WE GO.
C
DO odir = 1 , 2
DO i = 1 , 3
CALL syslop(coord,site,pstr(cnt),bits(1,i),nchar)
cnt = cnt + nchar
ENDDO
IF ( r3.GT.0 ) THEN
CALL mult(bits(1,1),bits(1,2),temp)
CALL mult(temp,bits(1,3),u3)
ELSE
DO i = 1 , 18
u3(i) = bits(i,2)
ENDDO
ENDIF
C
C CALCULATE OBSERVABLES
C
CALL observ(tot,u1,u2,u3)
C
ENDDO
C
C MOVE TO NEXT SITE IN LATTICE.
C
site = site + mov(1)
ENDDO
ENDDO
ENDDO
C
RETURN
END
**==qqqlps.f
C
C J .W .FLOWER AUG 1985
C
C
C FROM ABOVE THE THREE QUARK ORIENTATION IS SPECIFIED AS FOLLOWS :
C
C
C * ^
C |
C | R3
C |
C * ------------|-------------------- *
C <---- R1 ---->|<-------- R2 -------->
C
C THE LIMB WITH R1 & R2 IS TAKEN TO LIE IN DIRECTION "DIR" AND THE
C OTHER LIMB IS TAKEN IN THE TWO TRANSVERSE DIRECTIONS SUCCESSIVELY.
C NOTE : SINCE WE HAVE ACCESS TO ONLY BARRED 'LEGS' WE NEED TO MAKE
C SOME TIMELIKE PIECES ON THE FLY.
C */
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, 1987
C
SUBROUTINE qqqlps(r1,r2,r3,t,dir,str)
INTEGER r1 , r2 , r3 , t , dir
INTEGER str(*)
INTEGER ctr , odir , temp
C
C DO R1 LIMB - GO BACK IN "DIR" FIRST.
C
ctr = 1
temp = 5 + dir
DO i = 0 , r1 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
C
temp = 23 + dir
DO i = 0 , r1 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
C
str(ctr) = 14
ctr = ctr + 1
C
temp = 27 + dir
DO i = 0 , r1 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
C
IF ( r1.NE.0 ) THEN
str(ctr) = 36
ctr = ctr + 1
ELSE
DO i = 0 , t - 1
str(ctr) = 4
ctr = ctr + 1
ENDDO
DO i = 0 , t - 1
str(ctr) = 30
ctr = ctr + 1
ENDDO
ENDIF
C
temp = 23 + dir
DO i = 0 , r1 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
C
str(ctr) = 14
ctr = ctr + 1
C
temp = 27 + dir
DO i = 0 , r1 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
C
DO i = 0 , t - 1
str(ctr) = 26
ctr = ctr + 1
ENDDO
C
temp = 1 + dir
DO i = 0 , r1 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
C
DO i = 0 , t - 1
str(ctr) = 30
ctr = ctr + 1
ENDDO
C
str(ctr) = 14
ctr = ctr + 1
C
C DO R2 LIMB - GO FORWARD IN "DIR" FIRST.
C
temp = 1 + dir
DO i = 0 , r2 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
C
temp = 27 + dir
DO i = 0 , r2 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
C
str(ctr) = 14
ctr = ctr + 1
C
temp = 23 + dir
DO i = 0 , r2 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
C
IF ( r2.NE.0 ) THEN
str(ctr) = 36
ctr = ctr + 1
ELSE
DO i = 0 , t - 1
str(ctr) = 4
ctr = ctr + 1
ENDDO
DO i = 0 , t - 1
str(ctr) = 30
ctr = ctr + 1
ENDDO
ENDIF
C
temp = 27 + dir
DO i = 0 , r2 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
C
str(ctr) = 14
ctr = ctr + 1
C
temp = 23 + dir
DO i = 0 , r2 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
C
DO i = 0 , t - 1
str(ctr) = 26
ctr = ctr + 1
ENDDO
C
temp = 5 + dir
DO i = 0 , r2 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
C
DO i = 0 , t - 1
str(ctr) = 30
ctr = ctr + 1
ENDDO
C
str(ctr) = 14
ctr = ctr + 1
C
C NOW GET OTHER LEG ... LOOP OVER TRANSVERSE DIRECTIONS
C .... AND ALSO GO FORWARD AND (EVENTUALLY) BACK TO R3.
C
DO ntrans = 1 , 2
odir = mod((dir+ntrans),3)
temp = 1 + odir
DO i = 0 , r3 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
temp = 27 + odir
DO i = 0 , r3 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
str(ctr) = 14
ctr = ctr + 1
C
temp = 23 + odir
DO i = 0 , r3 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
IF ( r3.NE.0 ) THEN
str(ctr) = 36
ctr = ctr + 1
ELSE
DO i = 0 , t - 1
str(ctr) = 4
ctr = ctr + 1
ENDDO
DO i = 0 , t - 1
str(ctr) = 30
ctr = ctr + 1
ENDDO
ENDIF
temp = 27 + odir
DO i = 0 , r3 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
str(ctr) = 14
ctr = ctr + 1
C
temp = 23 + odir
DO i = 0 , r3 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
DO i = 0 , t - 1
str(ctr) = 26
ctr = ctr + 1
ENDDO
temp = 5 + odir
DO i = 0 , r3 - 1
str(ctr) = temp
ctr = ctr + 1
ENDDO
DO i = 0 , t - 1
str(ctr) = 30
ctr = ctr + 1
ENDDO
str(ctr) = 14
C
ctr = ctr + 1
ENDDO
C
RETURN
END
**==qqqmea.f
C
C J. W. FLOWER & STEVE OTTO SEPTEMBER 1984
C
C MEASURES THREE QUARK OBSERVABLES USING UBARS.
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, 1987
C
SUBROUTINE qqqmea(r1,r2,r3,tt)
INTEGER r1 , r2 , r3 , tt
C
INTEGER dims
PARAMETER (dims=4)
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mleg=18*l1*l2*l3)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mleg)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
C
INTEGER numstg
COMMON /comfil/ numstg
C
INTEGER dir , t
CMWB 6.22.88 CHARACTER*512 LPSTR
INTEGER lpstr(512)
DOUBLEPRECISION tot(2) , vals(2,15)
C
vals(1,tt) = 0.0
vals(2,tt) = 0.0
C
DO t = 0 , latt(4) - 1
CALL legmak(tt,t)
tot(1) = 0.0
tot(2) = 0.0
DO dir = 0 , 2
CALL qqqlps(r1,r2,r3,tt,dir,lpstr)
CALL qqqlpc(r1,r2,r3,t,lpstr,tot)
ENDDO
vals(1,tt) = vals(1,tt) + tot(1)
vals(2,tt) = vals(2,tt) + tot(2)
ENDDO
C
vals(1,tt) = vals(1,tt)/real(numstg*3*2)
vals(2,tt) = vals(2,tt)/real(numstg*3*2)
WRITE (6,1) r1 , r2 , r3 , tt , vals(1,tt) , vals(2,tt)
1 format ('Q (' , 3i3, ')' , ' X ' , i3 , ' : ' ,
& 2e19.9)
C
RETURN
END
**==rotlps.f
C MAKE THE STRINGS USED TO MEASURE OFF-AXIS WILSON LOOPS.
C STRING COMES IN THREE PARTS.
C (i) OUTWARD MOVES AND FAKE BACK
C (ii) FETCH LEG AND FAKE BACK.
C (iii GET TOP PART AND FAKE BACK.
C BARS ALL T LINKS EXCEPT IF R1+R2=1.
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, 1987
C
SUBROUTINE rotlps(pstr,r1,r2,t,dir1,dir2)
INTEGER pstr(256)
INTEGER r1 , r2 , t , dir1 , dir2
C
INTEGER r , dir , cnt , temp
C
C GO OUT IN SPATIAL DIRECTION.
C
num1 = 0
num2 = 0
cnt = 1
100 CONTINUE
IF ( (num1.LT.r1) .AND. (num2.LT.r2) ) THEN
temp = 1 + dir2
pstr(cnt) = temp
cnt = cnt + 1
num2 = num2 + 1
IF ( (num1+num2).LT.(r1+r2-1) ) THEN
temp = 15 + dir1
ELSE
temp = 1 + dir1
ENDIF
pstr(cnt) = temp
cnt = cnt + 1
num1 = num1 + 1
GOTO 100
ENDIF
C
IF ( num1.LT.r1 ) THEN
dir = dir1
ELSE
dir = dir2
ENDIF
C
DO r = (num1+num2) , (r1+r2-1)
IF ( (dir.EQ.dir1) .AND. (r.LT.r1+r2-1) .AND. (r.GT.0) ) THEN
temp = 15 + dir
ELSE
temp = 1 + dir
ENDIF
pstr(cnt) = temp
cnt = cnt + 1
ENDDO
C
C NOW FAKE BACK... THIS PART IS EASY !
C
temp = 27 + dir1
DO i = 0 , r1 - 1
pstr(cnt) = temp
cnt = cnt + 1
ENDDO
temp = 27 + dir2
DO i = 0 , r2 - 1
pstr(cnt) = temp
cnt = cnt + 1
ENDDO
pstr(cnt) = 14
cnt = cnt + 1
C
C NOW GO FETCH THE LEG. IF R1+R2 = 1 THEN WE HAVE TO DO
C THE LEG ON THE FLY
C
temp = 23 + dir1
DO i = 0 , r1 - 1
pstr(cnt) = temp
cnt = cnt + 1
ENDDO
temp = 23 + dir2
DO i = 0 , r2 - 1
pstr(cnt) = temp
cnt = cnt + 1
ENDDO
IF ( (r1+r2).NE.1 ) THEN
pstr(cnt) = 36
cnt = cnt + 1
ELSE
DO i = 0 , t - 1
pstr(cnt) = 4
cnt = cnt + 1
ENDDO
DO i = 0 , t - 1
pstr(cnt) = 30
cnt = cnt + 1
ENDDO
ENDIF
temp = 27 + dir1
DO i = 0 , r1 - 1
pstr(cnt) = temp
cnt = cnt + 1
ENDDO
temp = 27 + dir2
DO i = 0 , r2 - 1
pstr(cnt) = temp
cnt = cnt + 1
ENDDO
pstr(cnt) = 14
cnt = cnt + 1
C
C NOW FETCH THE TOP PART OF THE LOOP.
C
temp = 23 + dir1
DO i = 0 , r1 - 1
pstr(cnt) = temp
cnt = cnt + 1
ENDDO
temp = 23 + dir2
DO i = 0 , r2 - 1
pstr(cnt) = temp
cnt = cnt + 1
ENDDO
DO i = 0 , t - 1
pstr(cnt) = 26
cnt = cnt + 1
ENDDO
C
num1 = 0
num2 = 0
200 CONTINUE
IF ( (num1.LT.r1) .AND. (num2.LT.r2) ) THEN
temp = 5 + dir2
pstr(cnt) = temp
cnt = cnt + 1
num2 = num2 + 1
IF ( (num1+num2).LT.(r1+r2-1) ) THEN
temp = 19 + dir1
ELSE
temp = 5 + dir1
ENDIF
pstr(cnt) = temp
cnt = cnt + 1
num1 = num1 + 1
GOTO 200
ENDIF
C
IF ( num1.LT.r1 ) THEN
dir = dir1
ELSE
dir = dir2
ENDIF
DO r = (num1+num2) , (r1+r2-1)
IF ( (dir.EQ.dir1) .AND. (r.LT.(r1+r2-1)) .AND. (r.GT.0) ) THEN
temp = 19 + dir
ELSE
temp = 5 + dir
ENDIF
pstr(cnt) = temp
cnt = cnt + 1
ENDDO
C
DO i = 0 , t - 1
pstr(cnt) = 30
cnt = cnt + 1
ENDDO
C
pstr(cnt) = 14
C
RETURN
END
**==rotmea.f
C
C J. W. FLOWER OCT 1985
C
C THIS CODE MEASURES THE OFF AXIS WILSON LOOPS.
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, 1987
C
SUBROUTINE rotmea(dir,numhit,r1max,r2max)
INTEGER dir , numhit , r1max , r2max
C
INTEGER dims
PARAMETER (dims=4)
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mleg=18*l1*l2*l3)
PARAMETER (mbar=18*l1*l2*l3*l4)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mleg)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
C
DOUBLEPRECISION ubar(mbar) , sbar(mbar)
COMMON /comubr/ ubar , sbar
C
INTEGER tmax , t , plegs , cnt
INTEGER coord(4) , nchar , tflag , site , tmp
DOUBLEPRECISION bot(18) , top(18) , up(18)
DOUBLEPRECISION down(18) , stemp(18,2)
DOUBLEPRECISION temp
INTEGER pstr(256) , lpstr(256)
C
CALL linkbr(numhit,dir,sbar)
C
DO icord4 = 0 , latt(4) - 1
coord(4) = icord4
DO t = 1 , 12
tflag = 0
DO i = 1 , r1max
DO j = 1 , r2max
IF ( which(i,j).NE.0 ) THEN
tmax = min(12,which(i,j)+5)
IF ( (which(i,j).LE.t) .AND. (t.LE.tmax) )
& tflag = 1
ENDIF
ENDDO
ENDDO
IF ( tflag.EQ.1 ) THEN
CALL legmak(t,icord4)
C
C NOW LOOP THROUGH THE REQUIRED LOOPS.
C
DO i = 1 , r1max
DO j = 1 , r2max
IF ( which(i,j).NE.0 ) THEN
wloop(i,j) = 0.
tmax = min(12,which(i,j)+5)
IF ( (t.GE.which(i,j)) .AND. (t.LE.tmax) ) THEN
CALL rotlps(lpstr,i,j,t,dir,mod((dir+1),3))
site = coord(4)*mov(4)
plegs = 1
DO icord3 = 0 , latt(3) - 1
coord(3) = icord3
DO icord2 = 0 , latt(2) - 1
coord(2) = icord2
DO icord1 = 0 , latt(1) - 1
coord(1) = icord1
CALL udag(legs(plegs),down)
C CALL CPYMAT(PSTR,LPSTR,256) LPOINTER 4-3-89
CALL icpymt(pstr,lpstr,256)
cnt = 1
CALL syslop(coord,site,pstr(cnt),
& bot,nchar)
cnt = cnt + nchar
CALL syslop(coord,site,pstr(cnt),up,
& nchar)
cnt = cnt + nchar
CALL syslop(coord,site,pstr(cnt),
& top,nchar)
CALL mult(bot,up,stemp(1,1))
CALL mult(stemp(1,1),top,stemp(1,2))
CALL mult(stemp(1,2),down,stemp(1,1)
& )
wloop(i,j) = wloop(i,j)
& + (stemp(1,1)+stemp(9,1)
& +stemp(17,1))/3.0
plegs = plegs + 18
site = site + mov(1)
ENDDO
ENDDO
ENDDO
ENDIF
ENDIF
ENDDO
ENDDO
C
DO i = 1 , r1max
DO j = 1 , r2max
IF ( which(i,j).NE.0 ) THEN
temp = wloop(i,j)
tmp = t - which(i,j) + 1
res(i,j,tmp) = res(i,j,tmp) + temp
ENDIF
ENDDO
ENDDO
ENDIF
ENDDO
ENDDO
C
RETURN
END
**==setcol.f
SUBROUTINE setcol(latt)
INTEGER latt(4)
C
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (lismax=(l1*l2*l3*l4)/2+1)
COMMON /colcom/ lisred(lismax) , lisblk(lismax) , nred , nblack
C
nred = 0
nblack = 0
iptr = 0
index4 = 0
DO i = 1 , latt(4)
index3 = index4
DO j = 1 , latt(3)
index2 = index3
DO k = 1 , latt(2)
index1 = index2
DO l = 1 , latt(1)
IF ( index1.EQ.0 ) THEN
nred = nred + 1
lisred(nred) = iptr
ELSE
nblack = nblack + 1
lisblk(nblack) = iptr
ENDIF
iptr = iptr + 72
index1 = 1 - index1
ENDDO
index2 = 1 - index2
ENDDO
index3 = 1 - index3
ENDDO
index4 = 1 - index4
ENDDO
C
RETURN
END
**==setpar.f
C
C STEVE OTTO AND JON FLOWER JULY 1984.
C
C RECEIVES PARAMETER INPUT FOR THE SU3 ROUTINES.
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, 1987
C
SUBROUTINE setpar
C
INTEGER dims
PARAMETER (dims=4)
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mleg=18*l1*l2*l3)
INTEGER umax
PARAMETER (umax=72*l1*l2*l3*l4)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mleg)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
C
INTEGER numstg
COMMON /comfil/ numstg
C
100 CONTINUE
WRITE (6,*) 'ENTER X,Y,Z,T LATTICE SIZE ?'
READ (5,*) latt
WRITE (6,*) latt
num = latt(1)*latt(2)*latt(3)*latt(4)
C
IF ( num*72.GT.umax ) THEN
WRITE (6,*) 'NOT ENOUGH ROOM FOR GAUGE FIELD'
WRITE (6,*) '*******************************'
WRITE (6,*) 'TRY AGAIN'
GOTO 100
ENDIF
C
numsit = latt(1)*latt(2)*latt(3)*latt(4)
numstg = latt(1)*latt(2)*latt(3)*latt(4)
WRITE (6,*) 'ENTER VALUE OF BETA '
READ (5,*) beta
WRITE (6,*) beta
C
CALL make()
C
CALL setcol(latt)
C
WRITE (6,*) 'PARAMETERS NOW SET.'
WRITE (6,*) 'NOW NEED TO INITIALIZE..........'
C
RETURN
END
**==setsed.f
SUBROUTINE setsed()
C
INTEGER seed , mode
C
WRITE (6,*) 'SELECT MODE OF SEED SETTING: '
WRITE (6,*) '0 <-> BY HAND'
WRITE (6,*) '1 <-> READ IN FROM STANDARD FILE, SEED'
READ (5,*) mode
IF ( mode.EQ.0 ) THEN
WRITE (6,*) 'INPUT SEED'
READ (5,*) seed
CALL prnset(seed)
ELSE
OPEN (UNIT=1,FILE='SEED',ACCESS='SEQUENTIAL',STATUS='UNKNOWN')
REWIND (1)
READ (1,99001) seed
99001 FORMAT (I12)
CALL prnset(seed)
CLOSE (1)
ENDIF
WRITE (6,*) 'SEED SET'
RETURN
END
**==syslop.f
C
CCCCCCCCCCCCCCCCCCCCCCCCC END OF ROUTINE SETSED CCCCCCCCCCCCCCCCCCCC
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
C THIS SOFTWARE COPYRIGHTED BY CALTECH, 1983, 1986, 1987 C
C C
C PURE GAUGE SU(3) MONTE CARLO C
C C
C THIS CODE HAS A LONG HISTORY: C
C EARLY VERSIONS OF THIS DEVELOPED FOR THE MARK I C
C HYPERCUBE BY STEVE OTTO AND PAUL STOLORZ (1981-83) C
C JON FLOWER AND STEVE OTTO EXTENDED AND IMPROVED THE C
C CODE FOR THE MARK II MACHINES. (1984-86) C
C THE CURRENT INCARNATION IS A "PORT" TO THE NCUBE C
C DONE BY CLIVE BAILLIE, ED FELTEN, STEVE OTTO, C
C AND DAVID WALKER (1987) C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C SYSLOP SYSTEMATICALLY CALCULATES WILSON LOOPS
C FOR SU3 THEORY IN 3+1 DIM WITH THE INT STRINGS
C AS INPUT, ALSO KEEPS TRACK OF THE LENGTH OF THE STRING.
C
C CODING OF LATTICE PATHS:
C 1, 2, 3, 4 +X,+Y,+Z,+T
C 5, 6, 7, 8 -X,-Y,-Z,-T
C 15,16,17,18 X,+Y,+Z,+T USING BARS.
C 19,20,21,22 -X,-Y,-Z,-T USING BARS.
C 23,24,25,26 +X,+Y,+Z,+T FAKE STEPS
C 27,28,29,30 -X,-Y,-Z,-T FAKE STEPS
C 36 RETURNS LEG
C
C CONVERTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, 1987
C
SUBROUTINE syslop(coord,site,ptr,ploop,nn)
INTEGER coord(4) , nn , site , sites
INTEGER ptr(*)
DOUBLEPRECISION ploop(18)
C
INTEGER dims
PARAMETER (dims=4)
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mleg=18*l1*l2*l3)
INTEGER umax
PARAMETER (umax=72*l1*l2*l3*l4)
PARAMETER (mbar=18*l1*l2*l3*l4)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mleg)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
C
DOUBLEPRECISION u1(umax)
COMMON /comu1 / u1
C
DOUBLEPRECISION ubar(mbar) , sbar(mbar)
COMMON /comubr/ ubar , sbar
C
DOUBLEPRECISION ftemp(18,2) , filmat(18)
DOUBLEPRECISION extra(18)
INTEGER pu , pbar
INTEGER plegs , tindex , setflg
C
sites = site
index = 0
setflg = 0
nn = 0
C
100 CONTINUE
nn = nn + 1
ic = ptr(nn)
IF ( ic.EQ.14 ) THEN
C
C
DO i = 1 , 18
ploop(i) = ftemp(i,index+1)
ENDDO
C
RETURN
ELSE
C
IF ( ic.LE.4 ) THEN
ind = ic
pu = sites + rot(ind)
IF ( setflg.EQ.0 ) THEN
CALL cpymat(ftemp(1,index+1),u1(pu+1),18)
setflg = 1
ELSE
CALL cpymat(filmat,u1(pu+1),18)
tindex = 1 - index
CALL mult(ftemp(1,index+1),filmat,ftemp(1,tindex+1))
index = tindex
ENDIF
coord(ind) = mod(coord(ind)+1,latt1(ind)+1)
IF ( coord(ind).EQ.0 ) THEN
sites = sites - mov(ind)*latt1(ind)
ELSE
sites = sites + mov(ind)
ENDIF
C
ELSEIF ( ic.LE.8 ) THEN
ind = ic - 4
IF ( coord(ind).EQ.0 ) THEN
coord(ind) = latt1(ind)
sites = sites + mov(ind)*latt1(ind)
ELSE
coord(ind) = coord(ind) - 1
sites = sites - mov(ind)
ENDIF
pu = sites + rot(ind)
IF ( setflg.EQ.0 ) THEN
CALL cpymat(filmat,u1(pu+1),18)
CALL udag(filmat,ftemp(1,index+1))
setflg = 1
ELSE
CALL cpymat(extra,u1(pu+1),18)
CALL udag(extra,filmat)
tindex = 1 - index
CALL mult(ftemp(1,index+1),filmat,ftemp(1,tindex+1))
index = tindex
ENDIF
C
ELSEIF ( ic.EQ.36 ) THEN
lsite = sites - coord(4)*mov(4)
plegs = (lsite/72)*18
DO i = 1 , 18
ftemp(i,index+1) = legs(plegs+i)
ENDDO
setflg = 1
C
ELSEIF ( ic.LE.18 ) THEN
ind = ic - 14
pbar = (sites/72)*18
IF ( ind.EQ.4 ) THEN
CALL cpymat(filmat,ubar(pbar+1),18)
ELSE
CALL cpymat(filmat,sbar(pbar+1),18)
ENDIF
IF ( setflg.EQ.0 ) THEN
DO j = 1 , 18
ftemp(j,index+1) = filmat(j)
ENDDO
setflg = 1
ELSE
tindex = 1 - index
CALL mult(ftemp(1,index+1),filmat,ftemp(1,tindex+1))
index = tindex
ENDIF
coord(ind) = mod(coord(ind)+1,latt1(ind)+1)
IF ( coord(ind).EQ.0 ) THEN
sites = sites - mov(ind)*latt1(ind)
ELSE
sites = sites + mov(ind)
ENDIF
C
ELSEIF ( ic.LE.22 ) THEN
ind = ic - 18
IF ( coord(ind).EQ.0 ) THEN
coord(ind) = latt1(ind)
sites = sites + mov(ind)*latt1(ind)
ELSE
coord(ind) = coord(ind) - 1
sites = sites - mov(ind)
ENDIF
pbar = (sites/72)*18
IF ( ind.EQ.4 ) THEN
CALL cpymat(extra,ubar(pbar+1),18)
ELSE
CALL cpymat(extra,sbar(pbar+1),18)
ENDIF
IF ( setflg.EQ.0 ) THEN
CALL udag(extra,ftemp(1,index+1))
setflg = 1
ELSE
CALL udag(extra,filmat)
tindex = 1 - index
CALL mult(ftemp(1,index+1),filmat,ftemp(1,tindex+1))
index = tindex
ENDIF
C
ELSEIF ( ic.LE.26 ) THEN
ind = ic - 22
coord(ind) = mod(coord(ind)+1,latt1(ind)+1)
IF ( coord(ind).EQ.0 ) THEN
sites = sites - mov(ind)*latt1(ind)
ELSE
sites = sites + mov(ind)
ENDIF
C
ELSEIF ( ic.LE.30 ) THEN
ind = ic - 26
IF ( coord(ind).EQ.0 ) THEN
coord(ind) = latt1(ind)
sites = sites + mov(ind)*latt1(ind)
ELSE
coord(ind) = coord(ind) - 1
sites = sites - mov(ind)
ENDIF
C
ENDIF
GOTO 100
ENDIF
END
**==udag.f
C
C UDAG EVALUATES THE COMPLEX CONJUGATE OF THE INPUT 3x3 MATRIX
C
SUBROUTINE udag(pin,pout)
DOUBLEPRECISION pin(18) , pout(18)
C
pout(1) = pin(1)
pout(7) = pin(3)
pout(13) = pin(5)
pout(2) = -pin(2)
pout(8) = -pin(4)
pout(14) = -pin(6)
pout(3) = pin(7)
pout(9) = pin(9)
pout(15) = pin(11)
pout(4) = -pin(8)
pout(10) = -pin(10)
pout(16) = -pin(12)
pout(5) = pin(13)
pout(11) = pin(15)
pout(17) = pin(17)
pout(6) = -pin(14)
pout(12) = -pin(16)
pout(18) = -pin(18)
C
RETURN
END
**==update.f
C
C J. W. FLOWER AND STEVE OTTO JULY 1984.
C
C UPDATES THE LATICE "SWEEPS" TIMES WITH PROJECTIONS EVERY "PINVAL"
C
C CONVERRTED FROM C TO FORTRAN BY KEN BARISH, CLIVE BAILLIE,
C AND DAVID WALKER, 1987
C
SUBROUTINE update(sweeps,pinval)
INTEGER sweeps , pinval
C
INTEGER dims
PARAMETER (dims=4)
PARAMETER (l1=12,l2=12,l3=12,l4=12)
PARAMETER (mleg=18*l1*l2*l3)
INTEGER umax
PARAMETER (umax=72*l1*l2*l3*l4)
INTEGER latt(dims) , latt1(dims) , mov(dims) , rot(dims) ,
& which(13,13)
DOUBLEPRECISION wloop(13,13) , res(13,13,6)
DOUBLEPRECISION beta , legs(mleg)
COMMON /comsys/ beta , wloop , res , legs , latt , latt1 , mov ,
& rot , which
PARAMETER (lismax=(l1*l2*l3*l4)/2+1)
COMMON /colcom/ lisred(lismax) , lisblk(lismax) , nred , nblack
C
C
DOUBLEPRECISION u1(umax)
COMMON /comu1 / u1
C
INTEGER coord(4) , dir , origin
C
DO loop = 0 , sweeps - 1
C
C FIRST DO RED POINTS
C
DO i = 1 , nred
origin = lisred(i)
index = origin/72
CALL getpos(latt,coord,index)
DO dir = 0 , 3
CALL hit(coord,origin,dir)
ENDDO
ENDDO
C
C NOW DO BLACK POINTS
C
DO i = 1 , nblack
origin = lisblk(i)
index = origin/72
CALL getpos(latt,coord,index)
DO dir = 0 , 3
CALL hit(coord,origin,dir)
ENDDO
ENDDO
C
IF ( mod(loop,pinval).EQ.(pinval-1) ) CALL projec()
ENDDO
C
RETURN
END
>Fix:
>Release-Note:
>Audit-Trail:
>Unformatted:
More information about the Gcc-bugs
mailing list