      SUBROUTINE DERI1(C,NORBS,COORD,NUMBER,WORK,GRAD
     1                 ,F,MINEAR,FD,WMAT,HMAT,FMAT)
      implicit double precision (A-H,O-Z)
      include 'SIZES'
      include 'MOP6FILE'
      dimension WMAT(MPACK),HMAT(MPACK*2),FMAT(MPACK*2)
*********************************************************************
*
*     DERI1 COMPUTE THE NON-RELAXED DERIVATIVE OF THE NON-VARIATIONALLY
*     OPTIMIZED WAVEFUNCTION ENERGY WITH RESPECT TO ONE CARTESIAN
*     COORDINATE AT A TIME
*                             AND
*     COMPUTE THE NON-RELAXED FOCK MATRIX DERIVATIVE IN M.O BASIS AS
*     REQUIRED IN THE RELAXATION SECTION (ROUTINE 'DERI2').
*
*   INPUT
*     C(NORBS,NORBS) : M.O. COEFFICIENTS.
*     COORD  : CARTESIAN COORDINATES ARRAY.
*     NUMBER : LOCATION OF THE REQUIRED VARIABLE IN COORD.
*     WORK   : WORK ARRAY OF SIZE N*N.
*     WMAT     : WORK ARRAYS FOR d<PQ|RS> (2-CENTERS  A.O)
*   OUTPUT
*     C,COORD,NUMBER : NOT MODIFIED.
*     GRAD   : DERIVATIVE OF THE HEAT OF FORMATION WITH RESPECT TO
*              COORD(NUMBER), WITHOUT RELAXATION CORRECTION.
*     F(MINEAR) : NON-RELAXED FOCK MATRIX DERIVATIVE WITH RESPECT TO
*              COORD(NUMBER), EXPRESSED IN M.O BASIS, SCALED AND
*              PACKED, OFF-DIAGONAL BLOCKS ONLY.
*     FD     : IDEM BUT UNSCALED, DIAGONAL BLOCKS, C.I-ACTIVE ONLY.
*
************************************************************************
      common /MOLKST/ NUMAT,NAT(NUMATM),NFIRST(NUMATM),NMIDLE(NUMATM)
     1               ,NLAST(NUMATM), NDUMY1, NELECS,NALPHA,NBETA
     2               ,NCLOSE,NOPEN,NDUMY,FRACT
     3       /VECTOR/ CDUM(MORB2),EIGS(MAXORB),WDUM(MORB2),EIGB(MAXORB)
     4       /DENSTY/ P(MPACK), PA(MPACK), PB(MPACK)
      common /CIBITS/ NMOS,LAB,NELEC,NBO(3)
     1       /HMATRX/ H(MPACK)
     2       /XYIJKL/ XY(NMECI,NMECI,NMECI,NMECI)
     3       /CIVECT/ CONF(NMECI**4+NMECI**2)
      common /FOKMAT/ FDUMY(MPACK), SCALAR(MPACK)
      common /NVOMAT/ DIAG(MPACK/2)
     1       /KEYWRD/ KEYWRD
      common /NUMCAL/ NUMCAL
      dimension COORD(*),C(NORBS,NORBS),WORK(NORBS,NORBS),F(*),FD(*)
      character KEYWRD*241
      logical DEBUG
      data ICALCN /0/
C
      if (ICALCN.NE.NUMCAL) then
         DEBUG=INDEX(KEYWRD,'DERI1').NE.0
         IPRT=MIOB
         LINEAR=NORBS*(NORBS+1)/2
         ICALCN=NUMCAL
      end if
      if (DEBUG) call TIMER('BEFORE DERI1')
      STEP=1.D-3
C
C     2 POINTS FINITE DIFFERENCE TO GET THE INTEGRAL DERIVATIVES
C     ----------------------------------------------------------
C     STORED IN HMAT AND WMAT, WITHOUT DIVIDING BY THE STEP SIZE.
C
      NATI=(NUMBER-1)/3+1
      NATX=NUMBER-3*(NATI-1)
      call DHCORE (COORD,HMAT,WMAT,ENUCL2,NATI,NATX,STEP)
C
C HMAT HOLDS THE ONE-ELECTRON DERIVATIVES OF ATOM NATI FOR DIRECTION
C      NATX W.R.T. ALL OTHER ATOMS
C WMAT HOLDS THE TWO-ELECTRON DERIVATIVES OF ATOM NATI FOR DIRECTION
C      NATX W.R.T. ALL OTHER ATOMS
      STEP=0.5D0/STEP
C
C     NON-RELAXED FOCK MATRIX DERIVATIVE IN A.O BASIS.
C     ------------------------------------------------
C     STORED IN FMAT, DIVIDED BY STEP.
C
      call SCOPY (LINEAR,HMAT,1,FMAT,1)
      call DFOCK2 (FMAT,P,PA,WMAT,NUMAT,NFIRST,NMIDLE,NLAST,NATI)
C
C  FMAT HOLDS THE ONE PLUS TWO - ELECTRON DERIVATIVES OF ATOM NATI FOR
C       DIRECTION NATX W.R.T. ALL OTHER ATOMS
C
C       DERIVATIVE OF THE SCF-ONLY ENERGY (I.E BEFORE C.I CORRECTION)
C
      GRAD=(HELECT(NORBS,P,HMAT,FMAT)+ENUCL2)*STEP
C     TAKE STEP INTO ACCOUNT IN FMAT
      do 10 I=1,LINEAR
   10 FMAT(I)=FMAT(I)*STEP
C
C     RIGHT-HAND SIDE SUPER-VECTOR F = C' FMAT C USED IN RELAXATION
C     -----------------------------------------------------------
C     STORED IN NON-STANDARD PACKED FORM IN F(MINEAR) AND FD.
C     THE SUPERVECTOR IS THE NON-RELAXED FOCK MATRIX DERIVATIVE IN
C     M.O BASIS: F(IJ)= ( (C' * FOCK * C) (I,J) )   WITH I.GT.J .
C     F IS SCALED AND PACKED IN SUPERVECTOR FORM WITH
C                THE CONSECUTIVE FOLLOWING OFF-DIAGONAL BLOCKS:
C             1) OPEN-CLOSED  I.E. F(IJ)=F(I,J) WITH I OPEN & J CLOSED
C                                  AND I RUNNING FASTER THAN J,
C             2) VIRTUAL-CLOSED SAME RULE OF ORDERING,
C             3) VIRTUAL-OPEN   SAME RULE OF ORDERING.
C     FD IS PACKED OVER THE C.I-ACTIVE M.O WITH
C                THE CONSECUTIVE DIAGONAL BLOCKS:
C             1) CLOSED-CLOSED   IN CANONICAL ORDER, WITHOUT THE
C                                DIAGONAL ELEMENTS,
C             2) OPEN-OPEN       SAME RULE OF ORDERING,
C             3) VIRTUAL-VIRTUAL SAME RULE OF ORDERING.
C
C     PART 1 : WORK(N,N) = FMAT(N,N) * C(N,N)
      do 20 I=1,NORBS
   20 call SUPDOT (WORK(1,I),FMAT,C(1,I),NORBS,1)
C
C     PART 2 : F(IJ) =  (C' * WORK) (I,J) ... OFF-DIAGONAL BLOCKS.
      L=1
      if (NBO(2).NE.0 .AND. NBO(1).NE.0) then
C        OPEN-CLOSED
         call MTXM (C(1,NBO(1)+1),NBO(2),WORK,NORBS,F(L),NBO(1))
         L=L+NBO(2)*NBO(1)
      end if
      if (NBO(3).NE.0 .AND. NBO(1).NE.0) then
C        VIRTUAL-CLOSED
         call MTXM (C(1,NOPEN+1),NBO(3),WORK,NORBS,F(L),NBO(1))
         L=L+NBO(3)*NBO(1)
      end if
      if (NBO(3).NE.0 .AND. NBO(2).NE.0) then
C        VIRTUAL-OPEN
         call MTXM (C(1,NOPEN+1),NBO(3),WORK(1,NBO(1)+1),NORBS,F(L),NBO(
     12))
      end if
C     SCALE F ACCORDING TO THE DIAGONAL METRIC TENSOR 'SCALAR '.
      do 30 I=1,MINEAR
   30 F(I)=F(I)*SCALAR(I)
      if (DEBUG) then
         write(MIOB,*) ' F IN DERI1'
         J=MIN(20,MINEAR)
         write(MIOB,'(5F12.6)') (F(I),I=1,J)
      end if
C
C     PART 3 : SUPER-VECTOR FD, C.I-ACTIVE DIAGONAL BLOCKS, UNSCALED.
      L=1
      Nend=0
      do 50 LOOP=1,3
         NINIT=Nend+1
         Nend =NEND+NBO(LOOP)
         N1=MAX(NINIT,NELEC+1   )
         N2=MIN(Nend ,NELEC+NMOS)
         if (N2.LT.N1) go to 50
         do 40 I=N1,N2
            if (I.GT.NINIT) then
               call MXM (C(1,I),1,WORK(1,NINIT),NORBS,FD(L),I-NINIT)
               L=L+I-NINIT
            end if
   40    continue
   50 continue
C
C     NON-RELAXED C.I CORRECTION TO THE ENERGY DERIVATIVE.
C     ----------------------------------------------------
C
C     C.I-ACTIVE FOCK EIGENVALUES DERIVATIVES, STORED IN FD(continueD).
      LCUT=L
      do 60 I=NELEC+1,NELEC+NMOS
         FD(L)=DOT(C(1,I),WORK(1,I),NORBS)
   60 L=L+1
C
C     C.I-ACTIVE 2-ELECTRONS INTEGRALS DERIVATIVES. STORED IN XY.
C   FMAT IS USED HERE AS SCRATCH SPACE
C
      call DIJKL1 (C(1,NELEC+1),NORBS,NATI,WMAT,FMAT,HMAT,FMAT)
      do 70 I=1,NMOS
         do 70 J=1,NMOS
            do 70 K=1,NMOS
               do 70 L=1,NMOS
   70 XY(I,J,K,L)=XY(I,J,K,L)*STEP
C
C     BUILD THE C.I MATRIX DERIVATIVE, STORED IN WMAT.
      call MECID (FD(LCUT-NELEC),GSE,EIGB,WORK)
      if (DEBUG) then
         write(MIOB,*) ' GSE:',GSE
C#      write(MIOB,*) ' EIGB:',(EIGB(I),I=1,10)
C#      write(MIOB,*) ' WORK:',(WORK(I,1),I=1,10)
      end if
      call MECIH (WORK,WMAT,NMOS,LAB)
C
C     NON-RELAXED C.I CONTRIBUTION TO THE ENERGY DERIVATIVE.
      call SUPDOT (WORK,WMAT,CONF,LAB,1)
      GRAD=(GRAD+DOT(CONF,WORK,LAB))*23.061D0
      if (DEBUG) then
         write(IPRT,'('' * * * GRADIENT COMPONENT NUMBER'',I4)') NUMBER
         write(IPRT,'('' NON-RELAXED C.I-ACTIVE FOCK EIGENVALUES '',
     1                ''DERIVATIVES (E.V.)'')')
         write(IPRT,'(8F10.4)') (FD(LCUT-1+I),I=1,NMOS)
         write(IPRT,'('' NON-RELAXED 2-ELECTRONS DERIVATIVES (E.V.)''/
     1''    I    J    K    L       d<I(1)J(1)|K(2)L(2)>'')')
         do 80 I=1,NMOS
            do 80 J=1,I
               do 80 K=1,I
                  LL=K
                  if (K.EQ.I) LL=J
                  do 80 L=1,LL
   80    write(IPRT,'(4I5,F20.10)')
     1              NELEC+I,NELEC+J,NELEC+K,NELEC+L,XY(I,J,K,L)
         write(IPRT,'('' NON-RELAXED GRADIENT COMPONENT'',F10.4,
     1'' KCAL/MOLE'')')GRAD
         call TIMER('AFTER DERI1')
      end if
      return
      end
