C
C
C  James O. Jensen
C  jojensen@crdec.apgea.army.mil
C  Chemical Research, Development and Engineering Center
C  SMCCR-RSP-C
C  Aberdeen Proving Ground
C  Maryland   21010-5423
C
C  Janet L. Jensen
C  jljensen@crdec6.apgea.army.mil
C  Chemical Research, Development and Engineering Center
C  SMCCR-DDT
C  Aberdeen Proving Ground
C  Maryland   21010-5423
C
C
      program mover 
C
      implicit double precision(a-h,o-z)
C 
C  The final frequencies are stored in FREQ.
C  The final vectors are stored in VEC.  The vectors are stored
C  in packed form ( ie. no extra zero's between the vectors).
C  BUF is just a scratch array.
C 
C  The final Geometry is stored in the array C. The geometry is
C   also stored with no extra zero's. 
C  The atomic numbers are stored in the array IAN. 
C 
C  This version of the program has BUF and VEC dimensioned to the
C   largest possible case.  A later version of this code should contain
C   a call to MALLOC.
C
      parameter (maxatm=400) 
      parameter (mxatm3=1200) 
      parameter (lenbuf=1440000)
C
      Dimension Buf(LenBuf)
      dimension vec(lenbuf) 
      dimension freq(mxatm3)
      integer getcmdc
      character*140 FilNam, itype*1
      Common /Mol/    NAtoms,ICharg,Multip,NAE,NBE,NE,NBasis,
     $                IAN(401),Atmchg(400),C(1200)
C
C  Size of the MOL common - used in reading the common from the checkfile. 
C
C
      LRwMol = 4*MaxAtm + InToWP(8+MaxAtm)
C
C   RYD is a constant that converts atomic units into wavenumbers.
C
       ryd = 5140.48707686d0
C
C  The program is currently getting the name of the checkfile
C   from the command line.  
C
      FilNam = ' '
      len=getcmdc(1,FilNam)
      if(len.lt.1) then
        write(6,'('' Specify input file name on command line.'')')
        stop
        endif
C
C  OPEN FILE FOR READING  - the checkfile
C
      call fopen(1,4,FilNam(1:len),0,iout)
C
      FilNam = ' '
      len=getcmdc(2,FilNam)
      if(len.lt.1) then
        write(6,'('' Specify output file name on command line.'')')
        stop
        endif
C
      open(unit=16,file=FilNam(1:len) )
C
C
C  NUM IS THE NUMBER OF THE SUBFILE TO READ
C    997 is the /MOL/ common block stored on the checkfile.
C
      num = 997 
C
C  OPEN MOL SUBFILE FOR READING
C
      call fileio(2,-num,0,0,0)
C
C  READ MOL
C    natoms is the first entry in the common.
C
      Call FileIO(2,num,LRwMol,NAtoms,0)
C
      write(16,1001) natoms 
1001  format(3x,i8) 
C
C
      write(16,1002) (ian(ii),ii=1,natoms)
1002  format(3x,8i9)
C
C  CONVERT COORDINATES INTO Z-MATRIX ORIENTAION 
C
      call zcon(c,natoms)
C
        do 40 i=1,natoms
        nndx = (i-1) * 3 
        write(16,1003) (c(nndx+ii),ii=1,3)
1003    format(3x,3f15.8) 
  40    continue
C
C
C  NUM IS THE NUMBER OF THE SUBFILE TO READ
C   585 is the second derivative matrix in cartesian coords.
C   stored on the checkfile.
C
      num= 585  
C
C  ITQRY is a useful command to see how much is left in the subfile.
C
      left = itqry(num)
      lread = left 
C     write(6,*) ' no on buffer =', left
C
C  OPEN SUBFILE FOR READING
C
      call fileio(2,-num,0,0,0)
C
C      if (lenbuf.lt.left) then
C         write(6,*) ' increase size of lenbuf '
C         write(6,*) ' left = ',left
C         write(6,*) ' lenbuf ', lenbuf  
C         stop 
C      endif 
C
C  READ SUBFILE
C
      call fileio(2,num,lread,buf,0)
C
      n3atms = natoms * 3
C
C  OUTPAK makes a nice output for symmetric matrices stored in
C   lower triangular form..
C
C     call outpak(buf,n3atms,1,6)  
C
C  Mass weight the second derivative matrix.
C
      call mswt(buf,ian,natoms)
C
C     call outpak(buf,n3atms,1,6)  
C
C  DIAGONALIZE THE MASS WEIGHTED SECOND DERIVATIVE MATRIX
C   FIND THE EIGENVECTORS AND EIGENVALUES
C 
      call eigen(buf,vec,n3atms,0) 
C
C
C  Put in all of the scaling factors in the frequencies
C    The vibrational frequencies are related to the 
C    square roots of the eigenvalues of the massweighted
C    second derivative matrix. 
C 
       do 50 jj = 1,n3atms 
       nndx = jj * (jj + 1) / 2
       ff = buf(nndx)
       if(ff.ge.0.0d0) then
         freq(jj) = dsqrt(ff) * ryd 
       else
         freq(jj) = - dsqrt(-ff) * ryd 
       endif
   50  continue
C
C
C  See what is there.  Print frequencies and displacement vectors.
C
       do 75 jj = 1,n3atms 
       write(16,1004)  freq(jj)
 1004  format(3x,f30.16)  
       mmdx = n3atms * (jj-1) 
         do 65 j = 1,natoms
         nndx = 3*(j-1) + mmdx 
         write (16,1005) (vec(nndx+ii),ii=1,3)
 1005    format(3x,3f15.8)
   65    continue
   75  continue
C
      stop 
      end
C
C
C
      SUBROUTINE OUTPAK (MATRIX,NROW,NCTL,NOUT)
C...........VERSION = 09/05/73/03                                               
C.......................................................................        
C
C OUTPAK PRINTS A REAL*8 SYMMETRIC MATRIX STORED IN ROW-PACKED LOWER            
C
C TRIANGULAR FORM (SEE DIAGRAM BELOW) IN FORMATTED FORM WITH NUMBERED           
C
C ROWS AND COLUMNS.  THE INPUT IS AS FOLLOWS:                                   
C
C        MATRIX(*)...........PACKED MATRIX                                      
C
C        NROW................NUMBER OF ROWS TO BE OUTPUT                        
C
C        NCTL................CARRIAGE CONTROL FLAG: 1 FOR SINGLE SPACE,         
C                                                   2 FOR DOUBLE SPACE,         
C                                                   3 FOR TRIPLE SPACE.         
C
C        NOUT................UNIT NUMBER FOR OUTPUT
C
C THE MATRIX ELEMENTS ARE ARRANGED IN STORAGE AS FOLLOWS:
C
C        1
C        2    3
C        4    5    6
C        7    8    9   10
C       11   12   13   14   15
C       16   17   18   19   20   21
C       22   23   24   25   26   27   28
C
C AND SO ON.
C
C OUTPAK IS SET UP TO HANDLE 6 COLUMNS/PAGE WITH A 6F20.14 FORMAT               
C
C FOR THE COLUMNS.  IF A DIFFERENT NUMBER OF COLUMNS IS REQUIRED, CHANGE        
C
C FORMATS 1000 AND 2000, AND INITIALIZE KCOL WITH THE NEW NUMBER OF             
C
C COLUMNS.
C
C AUTHOR:  NELSON H.F. BEEBE, QUANTUM THEORY PROJECT, UNIVERSITY OF
C          FLORIDA, GAINESVILLE
C
C.......................................................................        
      REAL*8 MATRIX,COLUMN
      INTEGER BEGIN,ASA,BLANK,CTL
      DIMENSION MATRIX(1),ASA(3)
      DATA KCOL/3/, COLUMN/8HCOLUMN  /, ASA/4H    , 4H0   , 4H-   /,
     &     BLANK/4H    /, ZERO/0.D+00/
      CTL = BLANK
      IF ((NCTL.LE.3).AND.(NCTL.GT.0)) CTL = ASA(NCTL)
C
C LAST IS THE LAST COLUMN NUMBER IN THE ROW CURRENTLY BEING PRINTED
C
      LAST = MIN0(NROW,KCOL)
C
C BEGIN IS THE FIRST COLUMN NUMBER IN THE ROW CURRENTLY BEING PRINTED.
C
      NCOL=1
C.....BEGIN NON STANDARD DO LOOP.
      BEGIN=1
 1050 NCOL = 1
      WRITE (NOUT,1000) (COLUMN,I,I = BEGIN,LAST)
      DO 40 K = BEGIN,NROW
      KTOTAL = (K*(K-1))/2 + BEGIN - 1
      DO 10 I = 1,NCOL
      IF (MATRIX(KTOTAL+I) .NE. ZERO) GO TO 20
   10 CONTINUE
      GO TO 30
   20 WRITE (NOUT,2000) CTL,K,(MATRIX(I+KTOTAL),I=1,NCOL)
   30 IF (K .LT. (BEGIN+KCOL-1)) NCOL = NCOL + 1
   40 CONTINUE
      LAST = MIN0(LAST+KCOL,NROW)
      BEGIN=BEGIN+NCOL
      IF (BEGIN.LE.NROW) GO TO 1050
 1000 FORMAT (/12X,5(5X,A6,I4,5X),(5X,A6,I4))
 2000 FORMAT (A1,4H ROW,I4,2X,6D20.12)
      RETURN
      END
C