C***********************************************************

PROGRAM DMATRIX

C***********************************************************

CHARACTER*20 INPUT_FILE, OUTPUT_FILE, RESULT_FILE

INTEGER TIME_SIZE,SPACE_SIZE1,SPACE_SIZE2,SPACE_SIZE3,

& VOLUME

PARAMETER(TIME_SIZE=4)

PARAMETER(SPACE_SIZE1=4,SPACE_SIZE2=4,SPACE_SIZE3=4)

PARAMETER(VOLUME=TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2

& *SPACE_SIZE3)

PARAMETER(INPUT_FILE="DATA_IN.XYZ")

PARAMETER(OUTPUT_FILE="DATA_OUT.XYZ")

PARAMETER(RESULT_FILE="OUT11")

PARAMETER(MAX_INTERATIONS=100,N_TEST=1)

PARAMETER(MAX_LAN_ITER=400)

PARAMETER(NC2=9,ND=4,MU=4)

REAL K_WILSON

COMPLEX U(NC2,VOLUME,MU)

COMPLEX GAMMA(ND,ND,MU)

COMMON /GAMMA_MATRIX/ GAMMA

COMMON /GAUGE_FIELD/ U

COMMON /WILSON/ K_WILSON

OPEN(UNIT=10,FILE="result",STATUS="UNKNOWN")

OPEN(UNIT=11,FILE="gauge_nadia.dat",status="old")

K_WILSON=0.1

IIOPTION=1

CALL GAMMA_MATRICES

CALL GAUGEFIELD(IIOPTION)

CALL D_MATRIX

END

C**********************************************************

SUBROUTINE GAMMA_MATRICES

C**********************************************************

PARAMETER(ND=4,MU=4)

COMPLEX GAMMA(ND,ND,MU)

COMMON /GAMMA_MATRIX/ GAMMA

DATA ((GAMMA(i,j,1),j=1,4),i=1,4)/

& (-1.,0.) ,(0.,0.), (0.,0.), (0.,0.),

& (0.,0.), (-1.,0.), (0.,0.), (0.,0.),

& (0.,0.), (0.,0.), (1.,0.), (0.,0.),

& (0.,0.), (0.,0.), (0.,0.), (1.,0.) /

C

DATA ((GAMMA(i,j,2),j=1,4),i=1,4)/

& (0.,0.), (0.,0.), (0.,0.), (0.,-1.),

& (0.,0.), (0.,0.), (0.,-1.), (0.,0.),

& (0.,0.), (0.,1.), (0.,0.), (0.,0.),

& (0.,1.), (0.,0.), (0.,0.), (0.,0.)/

C

DATA ((GAMMA(i,j,3),j=1,4),i=1,4)/

& (0.,0.), (0.,0.), (0.,0.), (-1.,0.),

& (0.,0.), (0.,0.), (1.,0.), (0.,0.),

& (0.,0.), (1.,0.), (0.,0.), (0.,0.),

& (-1.,0.), (0.,0.), (0.,0.), (0.,0.)/

C

DATA ((GAMMA(i,j,4),j=1,4),i=1,4)/

& (0.,0.), (0.,0.), (0.,-1.), (0.,0.),

& (0.,0.), (0.,0.), (0.,0.), (0.,1.),

& (0.,1.), (0.,0.), (0.,0.), (0.,0.),

& (0.,0.), (0.,-1.), (0.,0.), (0.,0.)/

RETURN

END

C**********************************************************

SUBROUTINE GAUGEFIELD(IIOPTION)

C**********************************************************

INTEGER TIME_SIZE,SPACE_SIZE1,SPACE_SIZE2,SPACE_SIZE3,

& VOLUME

PARAMETER(TIME_SIZE=4)

PARAMETER(SPACE_SIZE1=4,SPACE_SIZE2=4,SPACE_SIZE3=4)

PARAMETER(VOLUME=TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2

& *SPACE_SIZE3)

PARAMETER(NC2=9,MU=4)

COMPLEX U(NC2,VOLUME,MU)

COMMON /GAUGE_FIELD/ U

IF (IIOPTION.EQ.1) THEN

DO IU=1,MU

read(11,*) (( U(IC,IV,IU),ic=1,9),iv=1,volume)

END DO

END IF

IF (IIOPTION.EQ.2) THEN

DO 10 IV=1,VOLUME

DO 20 IU=1,MU

DO 30 IC=1,NC2

U(IC,IV,IU)=0.

30 CONTINUE

U(1,IV,IU)=(1,0)

U(5,IV,IU)=(1,0)

U(9,IV,IU)=(1,0)

20 CONTINUE

10 CONTINUE

END IF

RETURN

END

C**********************************************************

SUBROUTINE D_MATRIX

C**********************************************************

INTEGER TIME_SIZE,SPACE_SIZE1,SPACE_SIZE2,SPACE_SIZE3,

& VOLUME

PARAMETER(TIME_SIZE=4)

PARAMETER(SPACE_SIZE1=4,SPACE_SIZE2=4,SPACE_SIZE3=4)

PARAMETER(VOLUME=TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2

& *SPACE_SIZE3)

PARAMETER(MAX_INTERATIONS=100,N_TEST=1)

PARAMETER(MAX_LAN_ITER=400)

PARAMETER(NCND=12,NC2=9,NC=3,ND=4,MU=4)

PARAMETER (NSS=NC*ND*VOLUME)

COMPLEX GAMMA(ND,ND,MU)

COMPLEX UGAMMA(ND,NC,ND,NC,MU)

COMPLEX UGAMMAADJ(ND,NC,ND,NC,MU)

COMPLEX*16 D(NSS,NSS)

COMPLEX U(NC2,VOLUME,MU)

COMMON /GAMMA_MATRIX/ GAMMA

COMMON /GAUGE_FIELD/ U

CALL INITIAL(D)

DO IT=1,TIME_SIZE

DO IX=1,SPACE_SIZE1

DO IY=1,SPACE_SIZE2

DO IZ=1,SPACE_SIZE3

CALL MULTIPLY(U,IS,UGAMMA,ugammaadj)

DO IC=1,NC

DO JC=1,NC

DO ID=1,ND

DO JD=1,ND

ITT1=IT-1

ITT2=IT+1

IXX1=IX-1

IXX2=IX+1

IYY1=IY-1

IYY2=IY+1

IZZ1=IZ-1

IZZ2=IZ+1

IF(IT.EQ.1) ITT1=TIME_SIZE

IF(IT.EQ.TIME_SIZE) ITT2=1

IF(IX.EQ.1) IXX1=SPACE_SIZE1

IF(IX.EQ.SPACE_SIZE1) IXX2=1

IF(IY.EQ.1) IYY1=SPACE_SIZE2

IF(IY.EQ.SPACE_SIZE2) IYY2=1

IF(IZ.EQ.1) IZZ1=SPACE_SIZE3

IF(IZ.EQ.SPACE_SIZE3) IZZ2=1

ISSS=id+(ic-1)*4+(iT-1)*12+(IX-1)*TIME_SIZE*12+(IY-1)*TIME_SIZE*

& SPACE_SIZE1*12+(IZ-1)*TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2*12

IF(IC.EQ.JC.AND.ID.EQ.JD) THEN

D(Isss,isss)=CMPLX(1,0)

END IF

c**************************************************************************

JSSS=Jd+(Jc-1)*4+(iTT2-1)*12+(IX-1)*TIME_SIZE*12+(IY-1)*TIME_SIZE*

& SPACE_SIZE1*12+(IZ-1)*TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2*12

D(ISSS,JSSS)=UGAMMA(ID,IC,JD,JC,1)

JSSS=Jd+(Jc-1)*4+(iT-1)*12+(IXX2-1)*TIME_SIZE*12

& +(IY-1)*TIME_SIZE*

& SPACE_SIZE1*12+(IZ-1)*TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2*12

D(ISSS,JSSS)=UGAMMA(ID,IC,JD,JC,2)

JSSS=Jd+(Jc-1)*4+(iT-1)*12+(IX-1)*TIME_SIZE*12+

& (IYY2-1)*TIME_SIZE*

& SPACE_SIZE1*12+(IZ-1)*TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2*12

D(ISSS,JSSS)=UGAMMA(ID,IC,JD,JC,3)

JSSS=Jd+(Jc-1)*4+(iT-1)*12+(IX-1)*TIME_SIZE*12+(IY-1)*TIME_SIZE*

& SPACE_SIZE1*12+(IZZ2-1)*TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2*12

D(ISSS,JSSS)=UGAMMA(ID,IC,JD,JC,4)

C*************************************************************************

JSSS=Jd+(Jc-1)*4+(iTT1-1)*12+(IX-1)*TIME_SIZE*12+(IY-1)*TIME_SIZE*

& SPACE_SIZE1*12+(IZ-1)*TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2*12

D(ISSS,JSSS) =UGAMMAadj(ID,IC,JD,JC,1)

JSSS=Jd+(Jc-1)*4+(iT-1)*12+(IXX1-1)*TIME_SIZE*12+(IY-1)*TIME_SIZE*

& SPACE_SIZE1*12+(IZ-1)*TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2*12

D(ISSS,JSSS)=UGAMMAadj(ID,IC,JD,JC,2)

JSSS=Jd+(Jc-1)*4+(iT-1)*12+(IX-1)*TIME_SIZE*12+(IYY1-1)*TIME_SIZE*

& SPACE_SIZE1*12+(IZ-1)*TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2*12

D(ISSS,JSSS)=UGAMMAadj(ID,IC,JD,JC,3)

JSSS=Jd+(Jc-1)*4+(iT-1)*12+(IX-1)*TIME_SIZE*12+(IY-1)*TIME_SIZE*

& SPACE_SIZE1*12+(IZZ1-1)*TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2*12

D(ISSS,JSSS)=UGAMMAadj(ID,IC,JD,JC,4)

END DO

END DO

END DO

END DO

END DO

END DO

END DO

END DO

CALL EIGENVALUES(D)

RETURN

END

C**********************************************************

SUBROUTINE INITIAL(D)

C**********************************************************

INTEGER TIME_SIZE,SPACE_SIZE1,SPACE_SIZE2,SPACE_SIZE3,

& VOLUME

PARAMETER(TIME_SIZE=4)

PARAMETER(SPACE_SIZE1=4,SPACE_SIZE2=4,SPACE_SIZE3=4)

PARAMETER(VOLUME=TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2

& *SPACE_SIZE3)

PARAMETER(NC=3,ND=4)

PARAMETER(NSS=nC*ND*VOLUME)

complex*16 D(NSS,NSS)

DO IT=1,TIME_SIZE

DO IX=1,SPACE_SIZE1

DO IY=1,SPACE_SIZE2

DO IZ=1,SPACE_SIZE3

DO JT=1,TIME_SIZE

DO JX=1,SPACE_SIZE1

DO JY=1,SPACE_SIZE2

DO JZ=1,SPACE_SIZE3

DO IC=1,NC

DO JC=1,NC

DO ID=1,ND

DO JD=1,ND

IS=IT+(IX-1)*TIME_SIZE+(IY-1)*TIME_SIZE*SPACE_SIZE1+

& (IZ-1)*TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2

JS=JT+(JX-1)*TIME_SIZE+(JY-1)*TIME_SIZE*SPACE_SIZE1+

& (JZ-1)*TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2

ISSS=ID+(IC-1)*4+(IS-1)*12

JSSS=JD+(JC-1)*4+(JS-1)*12

D(ISSS,JSSS)=cmplx(0,0)

END DO

END DO

END DO

END DO

END DO

END DO

END DO

END DO

END DO

END DO

END DO

END DO

RETURN

END

C**********************************************************

SUBROUTINE EIGENVALUES(D)

C**********************************************************

INTEGER TIME_SIZE,SPACE_SIZE1,SPACE_SIZE2,SPACE_SIZE3,

& VOLUME

PARAMETER(TIME_SIZE=4)

PARAMETER(SPACE_SIZE1=4,SPACE_SIZE2=4,SPACE_SIZE3=4)

PARAMETER(VOLUME=TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2

& *SPACE_SIZE3)

PARAMETER(NC=3,ND=4)

PARAMETER(NSS=NC*ND*VOLUME)

PARAMETER(NSSTEST=NC*ND*TIME_SIZE)

COMPLEX*16 D(NSS,NSS)

INTEGER NMAX,LDA,LDV,LWORK

PARAMETER(NMAX=NSS,LDA=NMAX,LDV=NMAX,LWORK=64*NMAX)

CHARACTER*1 JOB

INTEGER IFAIL,N

COMPLEX*16 V(LDV,NMAX),W(NMAX),WORK(LWORK)

REAL*8 RWORK(2*NMAX)

IFAIL=0

CALL F02GBF('N',NMAX,D,LDA,W,V,LDV,RWORK,WORK,LWORK,IFAIL)

WRITE(10,14) ('(',REAL(W(I)),',', IMAG(W(I)),')',I=1,NMAX)

14 format((3x,4(A,F7.4,A,F7.4,A,:)))

write(10,*) " "

write(10,*) " "

write(10,*) " "

RETURN

END

C***********************************************************

SUBROUTINE MULTIPLY(U,IS,UGAMMA,ugammaadj)

C***********************************************************

INTEGER TIME_SIZE,SPACE_SIZE1,SPACE_SIZE2,SPACE_SIZE3,

& VOLUME

PARAMETER(TIME_SIZE=4)

PARAMETER(SPACE_SIZE1=4,SPACE_SIZE2=4,SPACE_SIZE3=4)

PARAMETER(VOLUME=TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2

& *SPACE_SIZE3)

PARAMETER(NC=3,NC2=9,ND=4,MU=4)

COMPLEX IUNIT(ND,ND)

REAL K_WILSON

COMPLEX GAMMA(ND,ND,MU)

complex ugammaadj(nd,nc,nd,nc,mu)

COMPLEX Utemp(NC,NC),Ugamma(ND,NC,ND,NC,MU),Uadj(NC,NC)

COMPLEX U(NC2,VOLUME,MU)

COMMON /GAMMA_MATRIX/ GAMMA

COMMON /WILSON/ K_WILSON

common/unit/IUNIT

CALL UNITGAMMA

DO IMU=1,4

CALL ORDERING(U,IS,IMU,Utemp)

END DO

DO IMU=1,4

DO I1=1,4

DO I2=1,3

DO J1=1,4

DO J2=1,3

UGAMMA(I1,I2,J1,J2,IMU)=-K_WILSON*(IUNIT(I1,J1)-

& GAMMA(I1,J1,IMU))*Utemp(I2,J2)

END DO

END DO

END DO

END DO

END DO

CALL ADJLINK(Utemp,Uadj)

DO IMU=1,4

DO I1=1,4

DO I2=1,3

DO J1=1,4

DO J2=1,3

UGAMMAadj(I1,I2,J1,J2,IMU)=-K_WILSON*(IUNIT(I1,J1)+

& GAMMA(I1,J1,IMU))*UADJ(I2,J2)

END DO

END DO

END DO

END DO

END DO

RETURN

END

C**********************************************************

SUBROUTINE UNITGAMMA

C**********************************************************

PARAMETER (ND=4)

COMPLEX IUNIT(ND,ND)

COMMON /UNIT/IUNIT

DO ID=1,ND

DO JD=1,ND

IUNIT(ID,JD)=CMPLX(0.0,0.0)

IF (ID.EQ.JD) THEN

IUNIT(ID,JD)=CMPLX(1,0)

END IF

END DO

END DO

RETURN

END

C**********************************************************

SUBROUTINE ORDERING(Uin,IS,IMU,Uout)

C**********************************************************

INTEGER TIME_SIZE,SPACE_SIZE1,SPACE_SIZE2,SPACE_SIZE3,

& VOLUME

PARAMETER(TIME_SIZE=4)

PARAMETER(SPACE_SIZE1=4,SPACE_SIZE2=4,SPACE_SIZE3=4)

PARAMETER(VOLUME=TIME_SIZE*SPACE_SIZE1*SPACE_SIZE2

& *SPACE_SIZE3)

PARAMETER(NC2=9,MU=4)

COMPLEX Uin(NC2,VOLUME,MU)

COMPLEX Uout(NC2)

DO I=1,NC2

Uout(I)=Uin(I,IS,IMU)

END DO

RETURN

END

C***********************************************

SUBROUTINE ADJLINK(X,Y)

C***********************************************

DIMENSION X(18),Y(18)

Y(1)=X(1)

Y(2)=-X(2)

Y(3)=X(7)

Y(4)=-X(8)

Y(5)=X(13)

Y(6)=-X(14)

Y(7)=X(3)

Y(8)=-X(4)

Y(9)=X(9)

Y(10)=-X(10)

Y(11)=X(15)

Y(12)=-X(16)

Y(13)=X(5)

Y(14)=-X(6)

Y(15)=X(11)

Y(16)=-X(12)

Y(17)=X(17)

Y(18)=-X(18)

RETURN

END