PROGRAM GGBEST implicit real*8(a-h,o-z) implicit integer*4 (i-n) parameter (nat0=1000,nag0=40) DIMENSION r(3,nat0),iz(nat0),ri(3,nat0) CHARACTER*50 FIL CHARACTER*80 FN LOGICAL lzmat,auto CHARACTER*1 OK,DUM real*8 XX(3),YY(3),A(3),P(3),E0(3),SSE0E0,SSPP,AD, 1SSAE0,SSPE0,SSXXXX,SSYYYY,SSXXYY,APX,APY,APZ,SIGN, 2CB,CA,BANGLE,ANGLE,PI,SSE0AP dimension itypep(nag0),L0p(nag0),L1p(nag0),L2p(nag0),L3p(nag0) logical lex,lef CHARACTER*10 fnp(nag0) nat=0 natt=0 ld=0 c read also dummy atoms: DUM='N' inquire(file='AUTO',exist=auto) if(auto)then FIL='G98.OUT' lzmat=.true. else if(iargc().eq.0)then WRITE(*,*)' Full filename of the Gaussian output:' READ(*,'(A)')FIL WRITE(*,*)' Use the Z-matrix or standard orientation (Z/S) ?' READ(*,'(A)')OK else call getarg(1,FIL) call getarg(2,OK) if(iargc().gt.2)call getarg(3,DUM) if(DUM.eq.'y')DUM='Y' endif lzmat=.true. if(ok.eq.'s'.or.ok.eq.'S')lzmat=.false. endif inquire(file='LIST.PAR',exist=lex) if(lex)then write(6,*)'LIST.PAR found' open(45,file='LIST.PAR') read(45,*)nag if(nag.gt.nag0)then write(6,*)'Too many parameters' stop endif do 21 iag=1,nag read(45,4500)fnp(iag) write(6,4500)fnp(iag) 4500 format(a10) read(45,*)itypep(iag),L0p(iag),L1p(iag),L2p(iag),L3p(iag) 21 write(6,6005)itypep(iag),L0p(iag),L1p(iag),L2p(iag),L3p(iag) 6005 format(5i6) close(45) endif NG=0 nb=0 ES=99999.0d0 ESMIN=99999.0d0 nl=0 write(6,6009) 6009 format( 1' Reading output and retaining only the lowest-energy structure') inquire(file=FIL,exist=lef) if(.not.lef)goto 8989 OPEN(2,FILE=FIL) 1 READ(2,2000,END=1000,ERR=1000)FN 2000 FORMAT(A80) nl=nl+1 c c Read Energy: if(FN(1:9).eq.' SCF Done')then is=1 40 is=is+1 if(FN(is-1:is-1).ne.'='.and.is.lt.80)goto 40 it=is 44 it=it+1 if(FN(it:it).eq.' '.or.FN(it+1:it+1).ne.' ')goto 44 read(FN(is:it+1),*,err=921,end=921)ES 921 continue endif if(FN(28:34).eq.'EUMP2 =')then ES=0.0d0 read(FN(35:80),*,err=923,end=923)ES 923 continue endif if(FN(1:9).eq.' Energy= ')then ES=0.0d0 read(FN(9:26),*,err=922,end=922)ES 922 continue endif c c Read Geometry: IF((lzmat.and.(FN(19:39).EQ.'Z-Matrix orientation:'.OR. 1 FN(26:46).EQ.'Z-Matrix orientation:'.OR. 1 FN(20:37).EQ.'Input orientation:'.OR. 1 FN(27:44).EQ.'Input orientation:')) 1 .OR. 1 ((.not.lzmat).and. 2 (FN(20:40).EQ.'Standard orientation:'.OR. 2 FN(26:46).EQ.'Standard orientation:')))THEN NG=NG+1 ig98=0 if(FN(26:46).EQ.'Z-Matrix orientation:'.OR. 1 FN(27:44).EQ.'Input orientation:'.OR. 1 FN(26:46).EQ.'Standard orientation:')ig98=1 DO 4 I=1,4 nl=nl+1 4 READ(2,*) l=0 ld=0 lt=0 5 READ(2,2000)FN nl=nl+1 IF(FN(2:4).NE.'---')THEN l=l+1 lt=lt+1 BACKSPACE 2 nl=nl-1 if(ig98.eq.0)then READ(2,*)KA,KA,(r(i,l),i=1,3) else READ(2,*)KA,KA,r(1,l),(r(i,l),i=1,3) endif nl=nl+1 iz(l)=KA IF(KA.EQ.-1.and.DUM.ne.'Y')l=l-1 IF(KA.EQ.-1)ld=ld+1 GOTO 5 ENDIF nat=l natt=lt c write(6,*)'NL',nl c write(6,*)'NG',NG c write(6,*)'ES',ES c write(6,*)'ESMIN',ESMIN if(NG.eq.1.or.ES.lt.ESMIN)THEN ESMIN=ES nb=NG do 34 l=1,nat do 34 ix=1,3 34 ri(ix,l)=r(ix,l) ENDIF ENDIF GOTO 1 1000 CLOSE(2) write(6,*)nl,' lines' c 8989 if(NG.eq.0)then if(lex)then open(51,file='LIST.TXT',access='append') write(51,*)'No structure found' close(51) endif else write(6,*)NG,' structures' if(ld.ne.0)then write(*,*)natt,' atoms,',ld,'dummy,',nat,' taken' else write(*,*)nat,' atoms' endif write(6,*)'Best: ',nb write(6,*)'Appending FILE.X' OPEN(4,FILE='FILE.X',access='append') 990 write(4,4400)FIL(1:40),ESMIN write(6,4400)FIL(1:40),ESMIN 4400 format(A40,f20.8) write(4,*)nat do 35 l=1,nat if(iz(l).eq.-1)iz(l)=0 35 write(4,4000)iz(l),(ri(i,l),i=1,3),(0,i=1,7),0.0d0 4000 format(I3,3F12.6,7(1x,i1),f4.1) CLOSE(4) if(lex)then write(6,*)'Appending LIST.TXT' open(51,file='LIST.TXT',access='append') do 221 iag=1,nag itype=itypep(iag) L0=L0p(iag) L1=L1p(iag) L2=L2p(iag) L3=L3p(iag) PI=3.141592653589D0 A(1)=ri(1,L2)-ri(1,L3) A(2)=ri(2,L2)-ri(2,L3) A(3)=ri(3,L2)-ri(3,L3) P(1)=ri(1,L0)-ri(1,L1) P(2)=ri(2,L0)-ri(2,L1) P(3)=ri(3,L0)-ri(3,L1) E0(1)=ri(1,L1)-ri(1,L2) E0(2)=ri(2,L1)-ri(2,L2) E0(3)=ri(3,L1)-ri(3,L2) SSE0E0=E0(1)*E0(1)+E0(2)*E0(2)+E0(3)*E0(3) DO 8 KK=1,3 8 E0(KK)=E0(KK)/DSQRT(SSE0E0) SSPP=P(1)*P(1)+P(2)*P(2)+P(3)*P(3) AD=DSQRT(SSPP) SSAE0=A(1)*E0(1)+A(2)*E0(2)+A(3)*E0(3) SSPE0=P(1)*E0(1)+P(2)*E0(2)+P(3)*E0(3) DO 9 KK=1,3 XX(KK)=A(KK)-E0(KK)*SSAE0 9 YY(KK)=P(KK)-E0(KK)*SSPE0 SSXXXX=XX(1)*XX(1)+XX(2)*XX(2)+XX(3)*XX(3) SSYYYY=YY(1)*YY(1)+YY(2)*YY(2)+YY(3)*YY(3) SSXXYY=XX(1)*YY(1)+XX(2)*YY(2)+XX(3)*YY(3) APX=XX(2)*YY(3)-XX(3)*YY(2) APY=XX(3)*YY(1)-XX(1)*YY(3) APZ=XX(1)*YY(2)-XX(2)*YY(1) SSE0AP=APX*E0(1)+APY*E0(2)+APZ*E0(3) SIGN=1.0D0 IF (SSE0AP.NE.0.0D0)SIGN=-SSE0AP/ABS(SSE0AP) CB=0.0D0 IF (AD.NE.0.0D0)CB=-SSPE0/AD CA=0.0D0 IF (SSXXXX.NE.0.0D0.AND.SSYYYY.NE.0.0D0) 1 CA=-SSXXYY/DSQRT(SSXXXX)/DSQRT(SSYYYY) BANGLE=0.0D0 ANGLE=0.0D0 IF (CB.EQ.-1.0D0)BANGLE=PI IF (ABS(CB).LT.1.0D0)BANGLE=ACOS(CB) IF (CA.LE.-1.0D0) ANGLE=PI IF (ABS(CA).LT.1.0D0) ANGLE=ACOS(CA) ANGLE=180.0D0* ANGLE/PI BANGLE=180.0D0*BANGLE/PI ANGLE=SIGN*ANGLE if(ANGLE.lt.-180.0d0)ANGLE=360.0d0+ANGLE c torsion, distance, band angle if(itype.eq.2)WRITE(51,322)AD 322 format(f15.6,$) if(itype.eq.3)WRITE(51,322)BANGLE 221 if(itype.eq.4)WRITE(51,322)ANGLE WRITE(51,3221)ESMIN,NG,NB 3221 format(f20.8,2i4) close(51) endif endif STOP END