C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ C C FCONV - Conversion between internal and symmetry harmonic force C fields and GAMESS -> VIBCA. C (and generation of internal coordinate input data for VIBCA) c C C version 19.IV.2016 ----- Zbigniew KISIEL ----- C __________________________________________________ C | Institute of Physics, Polish Academy of Sciences | C | Al.Lotnikow 32/46, Warszawa, POLAND | C | kisiel@ifpan.edu.pl | C | http://info.ifpan.edu.pl/~kisiel/prospe.htm | C_________________________/-------------------------------------------------- C C C Modification history: c c 17.03.99: Adapted to deal with output from PC-Gamess of Granovsky c 7.05.02: Further patch for PC-GAMESS 5.1 (compatibility with other GAMESS C now probably lost) C 6.06.03: Echo of GAMESS frequencies and intensities to FC.VIB C 29.07.03: Debugging C 27.12.06: NMAX=100 + PCGAMESS7 compatibility C 10.08.07: correction to extraction of intensities and frequencies C 19.04.09: correction for the 'REDUCED MASS' line in PCG7.1 C 19.04.16: some debugging C C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ C C Data set has two forms, depending on required operation. For C conversions between symmetric and internal coordinates it should be C as follows: C C -Comment C -NINT, NSYM, IFORCE C -U matrix (NSYM lines of NINT elements) C -F matrix (either NINTxNINT or NSYMxNSYM - lower triangle only) C -Optionally definitions of internal coordinates (NINT lines), with five C numbers each, terminated by a line with five zeros. Coordinate C definition line can also carry an alphanumeric describing the C nature of the coordinate, this should be enclosed between = signs, C eg. =alpha=. These descriptors are used for annotation of extracted C potential energy contributions. C C IFORCE declares whether force constants are in aJ/A**2 etc (IFORCE=1) C or Nm-1 etc (IFORCE=2). C C----------------------------------------------------------------------------- C C For conversion of GAMESS output to form compatible with VIBCA the C data set should be a GAMESS output generated with options which print C the Hessian in internal coordinates. The first set of internal C coordinate definitions can carry coordinate identifiers in the same form C as used for internal coordinates above (appended at the end of C each coordinate line) C C FCONV will read GAMESS output created on both PC's and from UNIX C machines (the difference being extra space at beginning of each line) C and will output a complete VIBCA deck to FC.VIB. c Only for the case of linear-angle bends C the coordinates in the lines for the dummy atoms have to be reedited C according to the declarations in the $LIBE group in GAMESS. C C----------------------------------------------------------------------------- C IMPLICIT REAL*8(A-H,O-Z) PARAMETER (NMAX=100,NFMAX=10000,NTERMS=10000,nmomax=3*NMAX) COMMON /BLK0/F,U,UP,FF,W,FVALS REAL*8 F(NMAX,NMAX),U(NMAX,NMAX),UP(NMAX,NMAX),FF(NMAX,NMAX), * W(NMAX,NMAX),FVALS(NFMAX),vibfre(nmomax),vibint(nmomax) REAL*8 ATCORD(NMAX,3) INTEGER INTCOR(NMAX,5),IPOINT(NTERMS,3),LLAB(NMAX),ICHARG(NMAX) CHARACTER FILNAM*50,FNAM(2)*8,COMENT*78,units(2)*40,label(nmax)*10 CHARACTER NAMES(36)*4 character line*80 DATA NAMES/' 1H',' ', * ' 6LI',' 9BE',' 10B',' 12C',' 14N',' 16O',' 19F',' ', * ' ',' ',' ','28SI',' 31P',' 32S','35CL',' ', * ' ',' ',' ',' ',' ',' ',' ',' ', * ' ',' ',' ',' ',' ',' ',' ',' ', * '79BR',' '/ C DATA UNITS/'aJ/A**2, aJ/A, aJ', * 'N/m, 10**(-10)N/rad, 10**(-20)J/rad**2'/ DATA FNAM/'internal','symmetry'/ C WRITE(*,1) 1 FORMAT(1X//' ',76(1H_)/' |',T79,'|'/ * ' |',18X,' FCONV - force constant conversions',T79,'|'/ * ' |',76(1H_),'|'/' version 19.IV.2016',T64,'Zbigniew KISIEL'/) C C...input C 3 WRITE(*,2) 2 FORMAT(1X/5x,' 1 = internal to symmetry',/ * 5x,' 2 = symmetry to internal',/ * 5x,' 3 = GAMESS to VIBCA ..... ',$) READ(*,'(I5)',ERR=3)IFLAG IF(IFLAG.LT.1.OR.IFLAG.GT.3)GOTO 3 4 WRITE(*,5) 5 FORMAT(1X/5x,' Name of data file: ',$) READ(*,'(A)',ERR=4)FILNAM OPEN(2,FILE=FILNAM,ERR=4,status='OLD') READ(2,'(A)')COMENT IF(IFLAG.EQ.3)GOTO 400 C C...input for symmetry <-> internal conversions C READ(2,*)NINT,NSYM,IFORCE IF(IFORCE.EQ.1)FMULT=1.D0 IF(IFORCE.EQ.2)FMULT=0.01D0 IF(NINT.GT.NMOMAX.OR.NSYM.GT.NMOMAX)THEN WRITE(*,'(1X/'' ***** Too many coordinates '',A//)')CHAR(7) STOP END IF DO 6 I=1,NSYM READ(2,*)(U(I,J),J=1,NINT) DO 7 J=1,NINT UP(J,I)=U(I,J) 7 CONTINUE 6 CONTINUE C NLIM=NINT IF(IFLAG.EQ.2)NLIM=NSYM JJ=1 DO 8 I=1,NLIM READ(2,*)(F(I,J),J=1,JJ) DO 9 J=1,JJ F(J,I)=F(I,J) 9 CONTINUE JJ=JJ+1 8 CONTINUE GOTO 401 C C . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . C C...INPUT FROM GAMESS OUTPUT C C...establish the number of atoms C 400 READ(2,'(A)',end=500)LINE IF(LINE(1:16).NE.'ATOM ATOMIC'.AND. * LINE(2:17).NE.'ATOM ATOMIC')GOTO 400 C IUNIX=0 IF(LINE(2:17).EQ.'ATOM ATOMIC')IUNIX=1 NATMS=0 READ(2,'(A)',end=502)LINE c 600 READ(2,'(A)',end=502)LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) IF(LINE(14:14).EQ.'.')THEN NATMS=NATMS+1 if(natms.gt.nmax)then write(*,'(1x/'' ***** ERROR: Too many atoms:''// * 10x,'' found ='',i8/ * 10x,''maximum ='',i8//)')natms,nmax endif GOTO 600 ENDIF C write(*,'(1x/i5,'' Atoms identified'')')NATMS if(natms.gt.0)goto 602 write(*,'(1x//'' ***** INPUT ERROR''//)') stop C C...error conditions C 500 write(*,501) 501 format(1x//' ***** END OF FILE REACHED on scan for atoms:'/ * ' string ''ATOM ATOMIC'' not found'//) stop C 502 write(*,503)natms 503 format(1x//' ***** END OF FILE REACHED on scan for atoms:'/ * ' reading atom:',i8//) stop C C...establish the number of internal coordinates C 602 READ(2,'(A)',end=800)LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) IF(LINE(21:40).NE.'INTERNAL COORDINATES')GOTO 602 c DO 603 I=1,5 READ(2,'(A)',end=506)LINE 603 CONTINUE c c pcgamess71 c c 1 STRETCH xxxx2xxxx1xxxxXxxxxXxxxxXxxxxX 2.6948311 1.4260433 c 1 STRETCH 2 1 2.6948311 1.4260433 c 2 STRETCH 3 2 2.6942242 1.4257222 c 3 BEND 3 2 1 2.0359574 116.6517647 c c pcgamess62 c c 7 STRETCH 8 3 1.8267751 0.9666878 c 8 STRETCH 9 5 1.8369443 0.9720691 c 9 BEND 3 2 1 2.0097374 115.1494730 c c wingamess27jun2006 c c 1 STRETCH 2 1 2.6948311 1.4260433 c 2 STRETCH 3 2 2.6942242 1.4257222 c 3 BEND 3 2 1 2.0359574 116.6517647 c c 604 READ(2,'(A)',end=506)LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) IF(LINE(41:41).EQ.'.'.or.LINE(53:53).EQ.'.')THEN <-- PCG7 NINT=NINT+1 if(nint.gt.nmomax)then write(*,'(1x/'' ***** ERROR: Too many internals:''// * 10x,'' found ='',i8/ * 10x,''maximum ='',i8//)')nint,nmomax stop endif GOTO 604 ENDIF write(*,'(i5,'' Internal coordinates identified'')')NINT if(nint.gt.0)goto 510 write(*,'(1x//'' ***** INPUT ERROR''//)') stop C C...error conditions C 506 write(*,507)nint 507 format(1x//' ***** END OF FILE REACHED on scan for internals:'/ * ' reading internal:',i8//) stop c 800 write(*,801) 801 format(1x//' ***** END OF FILE REACHED on scan for internals:'/ * ' Cannot find the string "INTERNAL COORDINATES"' * //) stop C C...determine the comment C 510 REWIND(2) 605 READ(2,'(A)',end=511)LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) IF(LINE(13:17).NE.'$DATA'.AND.LINE(13:17).NE.'$data')GOTO 605 <-- PCG7 c READ(2,'(A)')LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) COMENT=LINE(12:80) goto 650 C C...error conditions C 511 write(*,514) 514 format(1x//' ***** END OF FILE REACHED on scan for comment'/ * ' string ''$DATA'' not found'//) stop C 650 WRITE(*,601)COMENT 601 FORMAT(1X/1X,A//5x, * ' force constant scaling factor = ',$) READ(*,*,ERR=650)FMABI NMODES=3*NATMS C C...read the atomic coordinates C 421 READ(2,'(A)')LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) IF(LINE(1:16).NE.'ATOM ATOMIC')GOTO 421 READ(2,'(A)')LINE DO 422 I=1,NATMS READ(2,'(A)')LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) READ(LINE(11:),'(F5.1,F17.10,2F20.10)') * CHARGE,(ATCORD(I,J),J=1,3) ICHARG(I)=CHARGE c write(*,'(1x,F10.5,3F15.6)')CHARGE,(ATCORD(I,J),J=1,3) IF(ICHARG(I).GT.36)ICHARG(I)=36 DO 423 J=1,3 ATCORD(I,J)=ATCORD(I,J)*0.5291771D0 423 CONTINUE 422 CONTINUE C C...definitions of internal coordinates: C 406 READ(2,'(A)')LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) IF(LINE(21:40).NE.'INTERNAL COORDINATES')GOTO 406 DO 405 I=1,5 READ(2,'(A)')LINE 405 CONTINUE C LASTAT=NATMS NLINBE=0 DO 518 I=1,NINT READ(2,'(A)')LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) c if(line(41:41).eq.'.')then READ(LINE(13:30),'(6I3)')(INTCOR(I,II),II=2,5) IF(LINE(5:11).EQ.'STRETCH' )INTCOR(I,1)=1 IF(LINE(5:8) .EQ.'BEND' )INTCOR(I,1)=2 IF(LINE(5:12).EQ.'PLA.BEND')INTCOR(I,1)=3 IF(LINE(5:12).EQ.'LIN.BEND')INTCOR(I,1)=4 IF(LINE(5:11).EQ.'TORSION' )INTCOR(I,1)=6 ENDIF c if(line(53:53).eq.'.')then READ(LINE(15:44),'(6I5)')(INTCOR(I,II),II=2,5) IF(LINE(7:14).EQ.'STRETCH ')INTCOR(I,1)=1 IF(LINE(7:14).EQ.'BEND ')INTCOR(I,1)=2 IF(LINE(7:14).EQ.'PLA.BEND')INTCOR(I,1)=3 IF(LINE(7:14).EQ.'LIN.BEND')INTCOR(I,1)=4 IF(LINE(7:14).EQ.'TORSION ')INTCOR(I,1)=6 ENDIF C C...convert from ZMAT to VIBCA internal coordinate definitions for torsion, C out-of-plane bend and linear angle bend respectively C C str. bend. tor. plabe. libe. C GAMESS: 1 2 3 4 5 C VIBCA: 1 2 6 3 4 C C...torsion GAMESS: A-B-C-D VIBCA: A-B-C-D (simple torsion) C IF(INTCOR(I,1).EQ.6)THEN GOTO 412 ENDIF C C...planar bend GAMESS: A-B-C VIBCA: A-D-B C \D \C IF(INTCOR(I,1).EQ.3)THEN ITEMP=INTCOR(I,3) INTCOR(I,3)=INTCOR(I,4) INTCOR(I,4)=INTCOR(I,5) INTCOR(I,5)=ITEMP GOTO 412 ENDIF C C...linear bend IF(INTCOR(I,1).EQ.4)THEN LASTAT=LASTAT+1 INTCOR(I,5)=LASTAT NLINBE=NLINBE+1 GOTO 412 ENDIF C C...label field 412 NEQ=0 IFIRST=0 ILAST=0 DO 519 J=65,80 IF(LINE(J:J).EQ.'=')THEN IF(NEQ.EQ.0)THEN NEQ=1 IFIRST=J+1 GOTO 519 ENDIF IF(NEQ.EQ.1)THEN ILAST=J-1 GOTO 520 ENDIF ENDIF 519 CONTINUE 520 IF(IFIRST+ILAST.NE.0)THEN LABEL(I)=LINE(IFIRST:ILAST) LLAB(I)=ILAST-IFIRST+1 ELSE LABEL(I)=' ' LLAB(I)=1 ENDIF 518 CONTINUE C C...force constant matrix (HESSIAN MATRIX): C set NCOLFM to the number of columns in the force constant output C (5 for NT, 10 for UNIX) c 407 READ(2,'(A)')LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) IF(LINE(1:26).NE.'HESSIAN MATRIX IN INTERNAL')GOTO 407 c ncolfm=10 ncolfm=5 C READ(2,'(A)')LINE READ(2,'(A)')LINE READ(2,'(A)')LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) if(LINE(1:16).eq.'TO CONVERT UNITS')THEN READ(2,'(A)')LINE READ(2,'(A)')LINE READ(2,'(A)')LINE ELSE BACKSPACE(2) BACKSPACE(2) ENDIF C write(*,'(1x/ * '' Processing HESSIAN MATRIX IN INTERNAL COORDINATES:''/)') NFIRST=1 NLAST=ncolfm c 408 IF(NLAST.GT.NINT)NLAST=NINT READ(2,'(A)')LINE READ(2,'(A)')LINE write(*,'(1x,a)')line(1:62) READ(2,'(A)')LINE DO 410 I=1,NINT READ(2,*)NN,(F(I,J),J=NFIRST,NLAST) 410 CONTINUE IF(NLAST.LT.NINT)THEN NFIRST=NFIRST+ncolfm NLAST=NLAST+ncolfm GOTO 408 ENDIF C C...Conversion from atomic units of force constants (including multiplication C by a scaling factor C C stretch-stretch: hartree/bohr**2 = 15.5690 aJ/Angstr**2 C stretch-bend: hartree/(bohr*radian) = 8.2388 aJ/(Angstrom*radian) C bend-bend: hartree/radian**2 = 4.3598 aJ/radian**2 C DO 411 I=1,NINT DO 411 J=1,NINT FMULT=FMABI IF(INTCOR(I,1).EQ.1.AND.INTCOR(J,1).EQ.1)THEN FMULT=15.5690D0*FMULT GOTO 420 ENDIF IF(INTCOR(I,1).EQ.1.OR.INTCOR(J,1).EQ.1)THEN FMULT=8.2388D0*FMULT GOTO 420 ENDIF FMULT=4.3598D0*FMULT 420 FF(I,J)=FMULT*F(I,J) 411 CONTINUE OPEN(3,FILE='FC.RES',STATUS='UNKNOWN') WRITE(3,25)COMENT IFORCE=1 C C...write F-matrix to FC.RES C DO 116 I=1,NINT WRITE(3,114)(FF(I,J),J=1,NINT) 116 CONTINUE C C...read calculated IR frequencies and intensities (no longer used???) C 402 READ(2,'(A)',end=403)LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) if(LINE(28:41).NE.'IR INTENSITIES')goto 402 c NM2=0 462 READ(2,'(A)',END=403)LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) if(LINE(6:15).ne.'FREQUENCY:')goto 462 NM1=NM2+1 NM2=NM1+4 if(nm1.gt.nmodes)goto 403 if(nm2.gt.nmodes)nm2=nmodes read(line,461)(vibfre(i),i=nm1,nm2) 461 format(17x,5(F12.5)) c READ(2,'(A)')LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) read(line,461)(vibint(i),i=nm1,nm2) GOTO 462 C C . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . C C...INTERNAL <-> SYMMETRY CONVERSION AND OUTPUT C 401 OPEN(3,FILE='FC.RES',STATUS='UNKNOWN') WRITE(3,25)COMENT 25 FORMAT(1X,78(1H-)//1X,A//1X,78(1H-)/) WRITE(3,10)NINT,NSYM 10 FORMAT(I6,' internal coordinates'/I6,' symmetry coordinates') WRITE(3,1000)UNITS(IFORCE) 1000 FORMAT(1X/' Units of force constants are: 'a/20x, * 'for stretch, stretch-bend, and bend respectively.') WRITE(3,20) 20 FORMAT(1X/' U matrix:'/) DO 21 I=1,NSYM WRITE(3,23)(U(I,J),J=1,NINT) 21 CONTINUE 23 FORMAT(1X,20F8.5) WRITE(3,11)FNAM(IFLAG) 11 FORMAT(1X/' Force field in ',A,' coordinates:'/) DO 12 I=1,NLIM if(iforce.eq.1)WRITE(3,114)(F(I,J),J=1,NLIM) if(iforce.eq.2)WRITE(3,14)(F(I,J),J=1,NLIM) 12 CONTINUE 114 FORMAT(1X,20F8.4) 14 FORMAT(1X,20F8.2) GOTO(100,200),IFLAG C C C...internal to symmetry: F = U f U' C 100 CALL MULT(U,F,W,NSYM,NINT,NINT) CALL MULT(W,UP,FF,NSYM,NINT,NSYM) WRITE(3,11)FNAM(3-IFLAG) C C...output symmetry force constant matrix to FC.RES C DO 15 I=1,NSYM if(iforce.eq.1)WRITE(3,114)(FF(I,J),J=1,NSYM) if(iforce.eq.2)WRITE(3,14)(FF(I,J),J=1,NSYM) 15 CONTINUE GOTO 17 C C C...symmetry to internal: f = U'F U C 200 CALL MULT(UP,F,W,NINT,NSYM,NSYM) CALL MULT(W,U,FF,NINT,NSYM,NINT) WRITE(3,11)FNAM(3-IFLAG) C C...output internal force constant matrix to FC.RES C DO 16 I=1,NINT if(iforce.eq.1)WRITE(3,114)(FF(I,J),J=1,NINT) if(iforce.eq.2)WRITE(3,14)(FF(I,J),J=1,NINT) 16 CONTINUE C C...read data for VIBCA C C...read coordinate definitions DO 18 I=1,NINT READ(2,*,ERR=17,END=17)(INTCOR(I,II),II=1,5) C BACKSPACE(2) READ(2,'(A)')LINE NEQ=0 DO 380 J=1,80 IF(LINE(J:J).EQ.'=')THEN IF(NEQ.EQ.0)THEN NEQ=1 IFIRST=J+1 GOTO 380 ENDIF IF(NEQ.EQ.1)THEN ILAST=J-1 GOTO 381 ENDIF ENDIF 380 CONTINUE 381 LABEL(I)=LINE(IFIRST:ILAST) LLAB(I)=ILAST-IFIRST+1 C 18 CONTINUE C DO 55 I=1,NINT DO 55 J=1,NINT FF(I,J)=FF(I,J)*FMULT 55 CONTINUE C C . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . C C...G E N E R A T I O N O F VIBCA D E C K, File FC.VIB C C...Identification of diagonal potential terms C 403 NFVALS=0 NTERM=0 DO 39 I=1,NINT NTERM=NTERM+1 IF(I.EQ.1)THEN NFVALS=NFVALS+1 FVALS(NFVALS)=FF(I,I) NFV=NFVALS GOTO 50 ENDIF C DO 48 J=1,NFVALS IF(DINT(100000.D0*FVALS(J)).EQ.DINT(100000.D0*FF(I,I))) * THEN NFV=J GOTO 50 ENDIF 48 CONTINUE NFVALS=NFVALS+1 FVALS(NFVALS)=FF(I,I) NFV=NFVALS C 50 IPOINT(NTERM,1)=I IPOINT(NTERM,2)=NFV IPOINT(NTERM,3)=I 39 CONTINUE C C...Identification of off-diagonal potential terms C DO 19 I=1,NINT-1 DO 19 J=I+1,NINT NTERM=NTERM+1 DO 51 K=1,NFVALS IF(INT(10000.D0*ABS(FVALS(K))).EQ. * INT(10000.D0*ABS(FF(I,J))).AND. * FVALS(K)*FF(I,J).GE.0.D0)THEN NFV=K GOTO 52 ENDIF 51 CONTINUE NFVALS=NFVALS+1 FVALS(NFVALS)=FF(I,J) NFV=NFVALS 52 IPOINT(NTERM,1)=I IPOINT(NTERM,2)=NFV IPOINT(NTERM,3)=J 19 CONTINUE C CLOSE(3) OPEN(3,FILE='FC.VIB',STATUS='UNKNOWN') C C...Header lines and atomic coordinates (for GAMESS input only) C IF(IFLAG.EQ.3)THEN WRITE(3,'(A)')COMENT WRITE(3,'(5I5,f7.2)')NATMS,NFVALS,NLINBE/2,0,1,fmabi DO 425 I=1,NATMS WRITE(3,426)NAMES(ICHARG(I)),0,0,(ATCORD(I,J),J=1,3), * 0,0,0.0D0 425 CONTINUE IF(NLINBE.NE.0)THEN DO 526 I=1,NLINBE WRITE(3,426)' *',0,0,0.0,0.0,0.0,0,0,0.0D0 526 CONTINUE ENDIF 426 FORMAT(1X,A4,2I5,3F10.5,2I5,F10.4) ENDIF C C...Force constant values C do 60 II=1,nfvals,10 iii=ii+9 if(iii.gt.nfvals)iii=nfvals WRITE(3,33)(FVALS(iiii),iiii=ii,iii) 60 continue c33 FORMAT(1X,10F11.6) 33 FORMAT(1X,10F13.8) C C...Potential terms C DO 46 I=1,NTERM IF(IPOINT(I,1).EQ.IPOINT(I,3))THEN K=IPOINT(I,1) WRITE(3,31)(INTCOR(K,J),J=1,5),IPOINT(I,2), * 0,0,0,0,0,LABEL(K)(1:LLAB(K)) ELSE K=IPOINT(I,1) KK=IPOINT(I,3) WRITE(3,131)(INTCOR(K,J),J=1,5),IPOINT(I,2), * (INTCOR(KK,J),J=1,5),LABEL(K)(1:LLAB(K)), * LABEL(KK)(1:LLAB(KK)) ENDIF 46 CONTINUE 31 FORMAT(5I5,I5,5I5,2X,A) 131 FORMAT(5I5,I5,5I5,2X,A,'-',A) WRITE(3,131)0 C C...Termination of the main VIBCA deck C WRITE(3,'(A)')'END OF DATA' WRITE(3,'(A)')'END OF DATA' write(3,'(1x///)') C C C...Transfer frequencies and intensities after the end of the VIBCA deck C even if the code at label 402 above spots them C C if(nm2.gt.0)goto 717 REWIND(2) nm2=0 nmodes=0 700 read(2,'(A)',end=717)LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) if(line(1:16).eq.' FREQUENCY:')THEN read(line(18:78),'(5f12.2)')(vibfre(n),n=nmodes+1,nmodes+5) read(2,'(A)',end=17)LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) if(line(1:16).eq.' REDUCED MASS:')THEN read(2,'(A)',end=17)LINE IF(IUNIX.EQ.1)LINE=LINE(2:80) endif read(line(18:78),'(5f12.2)')(vibint(n),n=nmodes+1,nmodes+5) nmodes=nmodes+5 nm2=nmodes endif goto 700 C 717 WRITE(3,'(1x///1x,78(1h-)//1x, * ''Intensities (DEBYE**2/AMU-ANGSTROM**2) and Vibrational '', * ''Wavenumbers (cm**-1):''//T30,''GAMESS'',T56, * ''VIBCA''/ * '' Int Omega Omega''/ * )') C do 470 i=1,nm2 nm=i-6 if(nm.lt.0)nm=0 if(nm.eq.0.or.vibfre(i).ne.0.d0)then write(3,471)nm,VIBint(i),VIBfre(i) endif 470 continue 471 FORMAT(' Mode',I3,':',5x,2F15.6) write(3,'(1x/1x,78(1H-))') C 17 CLOSE(2) CLOSE(3) WRITE(*,'(1X/'' **** Done ****''//)') STOP END C C----------------------------------------------------------------------------- C SUBROUTINE MULT(A,B,C,N,M,L) C C C MULT WILL MULTIPLY TWO GENERAL, DOUBLE DIMENSIONED, DOUBLE PRECISION C ARRAYS AA AND BB, TO GIVE THE RESULT CC C C CC ( N TIMES L ) = AA ( N TIMES M ) * BB ( M TIMES L ) C C - A IS A TWO DIMENSIONAL ARRAY (NMAX x NMAX) CONTAINING C AA IN THE TOP LEFT HAND CORNER C C - B IS A TWO DIMENSIONAL ARRAY (NMAX x NMAX) CONTAINING C BB IN THE TOP LEFT HAND CORNER C C - C IS A TWO DIMENSIONAL ARRAY (NMAX x NMAX) CONTAINING C THE PRODUCT OF AA AND BB ON OUTPUT C C PARAMETER (NMAX=100) REAL*8 A(NMAX,NMAX),B(NMAX,NMAX),C(NMAX,NMAX),PROD C DO 1 I=1,N DO 1 J=1,L PROD=0.0D0 DO 2 K=1,M PROD=PROD+A(I,K)*B(K,J) 2 CONTINUE C(I,J)=PROD 1 CONTINUE C RETURN END C C----------------------------------------------------------------------------- C-----------------------------------------------------------------------------