C
C $Id: fmix.F,v 1.2 1998/07/16 16:39:46 jjv5 Exp arjan $
C
C------------------------------------------------------------------------
      SUBROUTINE FMIX(ITER,ICNVRG)
C
C     FOCK MATRIX MIXING SCHEME TO ACCELERATE SCF CONVERGENCE.
C
C     CALL AFTER FOCK MATRIX IS CONSTRUCTED, AFTER ENERGY EVALUATION,
C     BEFORE CALL TO MOSUB.
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
#include "divcon.dim"
#include "divcon.h"
      LOGICAL MIXED,BIGSTP
C
      SAVE MIXED,DE1,DP1
C
      ICNVRG = 0
      DCNVRG = 1.0D-5
      PCNVRG = 1.0D-5
      IIMAX = IIMAT(NATOMS+1)-1
      IJMAX = IJMAT(IP1(NATOMS+1))-1
      IF(ITER.EQ.1)THEN
        FNEW = 1.0D0
        DE2 = 1.0D0
        DP2 = 1.0D0
        MIXED = .FALSE.
      ENDIF
      IF(ITER.GT.1)THEN
        CALL DECALC(DE2,DP2)
C-RDC c       write(0,*) ' DE=',DE2,'  DP=',DP2
C-RDC c       write(8,*) ' DE=',DE2,'  DP=',DP2
      ENDIF
      IF(ITER.GT.2)THEN
        IF(ABS(DE2).LT.DCNVRG.AND.ABS(DE1).LT.DCNVRG.AND.
     .    DP2.LT.PCNVRG.AND.DP1.LT.PCNVRG)THEN
          ICNVRG = 1
          RETURN
        ENDIF
        IF(.NOT.MIXED)THEN
          MIXED = DP2.LT.0.002D0.AND.DP1.LT.0.004D0
        ENDIF
        IF(MIXED)THEN
          FMAX = 1.5D0
          SSIGN = DE1*DE2
          IF(SSIGN.GT.0.0D0.AND.ABS(DE2).GT.ABS(DE1))THEN
            FNEW = FMAX
            FOLD = 1.0D0 - FMAX
          ELSEIF(SSIGN.LT.0.0D0)THEN
            FNEW = 1.0D0
            FOLD = 0.0D0
          ELSE
            DENOM = DE1 - DE2
            DMIN = ABS(DE1)*1.0D-3
            IF(ABS(DENOM).LT.DMIN) DENOM = SIGN(DMIN,DENOM)
            FNEW = DE1/DENOM
            IF(FNEW.GT.FMAX) FNEW = FMAX
            IF(FNEW.LT.0.5D0) FNEW = 0.5D0
            FOLD = 1.0D0 - FNEW
          ENDIF
C-RDC c         write(0,*) ' FOLD=',FOLD,'  FNEW=',FNEW
C-RDC c         write(8,*) ' FOLD=',FOLD,'  FNEW=',FNEW
          DO 40 II=1,IIMAX
            FDIAG(II) = FOLD*FIIOLD(II) + FNEW*FDIAG(II)
 40       CONTINUE
          DO 50 IJ=1,IJMAX
            FDIAT(IJ) = FOLD*FIJOLD(IJ) + FNEW*FDIAT(IJ)
 50       CONTINUE
        ENDIF
      ENDIF
C
      DO 90 II=1,IIMAX
        FIIOLD(II) = FDIAG(II)
 90   CONTINUE
      DO 100 IJ=1,IJMAX
        FIJOLD(IJ) = FDIAT(IJ)
 100  CONTINUE
      DE1 = DE2
      DP1 = DP2
      RETURN
      END
C
C
C
      SUBROUTINE DECALC(DE,DP)
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
#include "divcon.dim"
#include "divcon.h"
C
      IIMAX = IIMAT(NATOMS+1)-1
      IJMAX = IJMAT(IP1(NATOMS+1))-1
C
      DE = 0.0D0
      DP = 0.0D0
      DO 100 II=1,IIMAX
        DPDIAG = PDIAG(II) - PIIOLD(II)
        DE = DE + FDIAG(II)*DPDIAG
        DP = DP + DPDIAG**2
 100  CONTINUE
C
      DEI = 0.0D0
      DPI = 0.0D0
      II = 0
      DO 200 I=1,NATOMS
        NORBSI = NATORB(IATNUM(I))
        IF(NORBSI.EQ.0) GO TO 200
        DO 150 IORB=1,NORBSI
          II = II + IORB
          DPDIAG = PDIAG(II) - PIIOLD(II)
          DEI = DEI + FDIAG(II)*DPDIAG
          DPI = DPI + DPDIAG**2
 150    CONTINUE
 200  CONTINUE
      DE = DE - 0.5D0*DEI
      DP = DP - 0.5D0*DPI
C
      DO 300 IJ=1,IJMAX
        DPDIAT = PDIAT(IJ) - PIJOLD(IJ)
        DE = DE + FDIAT(IJ)*DPDIAT
        DP = DP + DPDIAT**2
 300  CONTINUE
      DP = DSQRT(DP/(IIMAX+IJMAX))
      RETURN
      END