C 25 Mar 08 - JHS - NEW MODULE FOR NEO-NONORTHOGONAL-CI
C
C*MODULE NEONCI   *DECK NEONCI_DRV
      SUBROUTINE NEONCI_DRV
C
C     THIS IS THE DRIVER FOR NEO NONORTHOGANAL CI
C     IT DETERMINES IF A REQUEST TO STORE THE ELEC AND NUC
C     MOS WAS MADE, IF A REQUEST TO DO A NEONCI WAS ASKED OR
C     IF NOTHING ADDITIONAL IS REQUIRED. IT DIRECTS THE
C     COURSE OF THE PROGRAM ACCORDING TO THE REQUEST
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      LOGICAL NEONCI
C
      PARAMETER(MXATM=2000)
C
      COMMON /INFOA / NAT,ICH,MUL,NUM,NQMT,NE,NA,NB,
     *                ZAN(MXATM),C(3,MXATM),IAN(MXATM)
      COMMON /IOFILE/ IR,IW,IP,IS,IJKT,IDAF,NAV,IODA(950)
      COMMON /NORTHO/ NEONCI,LOCORB
      COMMON /NUCMOI/ NUNIQN,IUNIQN(20),IUNIQT(20),NQMNUC,IQMNUC(20),
     *                 IQNTYP(20),NUMNB,NUCST,NAUXNB,IAUXNB(20),NUMULT,
     *                 NNA,NNB,NTAUXB
C
C     IF "LOCORB" EXISTS IN INPUT THEN IDENTIFY THE FILE
C     AS A LEFT OR RIGHT  LOCALIZED NEO CALC AND STORE
C     THE CONVERGED HF ORBITALS FOR LATER USE IN THE NEONCI:
C            1 --> LEFT LOCALIZED NEO CALCULATION
C            2 --> RIGHT LOCALIZED NEO CALCULATION
      L1 = NUM
      N1 = NUMNB
C
      IF(LOCORB.EQ.1 .OR. LOCORB.EQ.2) THEN
C
C        ---  PRINT LOCALIZED ORBITALS TO DAF  ---
C
      CALL ORBLOCAL(LOCORB,L1,N1)
      RETURN
C
      ELSE IF(NEONCI) THEN
C
C        --- INPUT FILE IS A NEONCI RUN - CALL NEONCI ROUTINE ---
C
      WRITE(IW,9000)
      CALL NEONCI_SETUP
C
      RETURN
C
      ELSE
C
C        --- DO NOTHING - INPUT FILE IS ORDINARY ---
C
      ENDIF
C
      RETURN
 9000 FORMAT(/5X,42('-')/
     *        10X,' NEO-NONORTHOGONAL CI (NEO-NCI) ',/
     *        9X,'PROGRAM WRITTEN BY JONATHAN H SKONE'/
     *        5X,42('-')/)
      END
C
C*MODULE NEONCI   *DECK ORBLOCAL
      SUBROUTINE ORBLOCAL(LOCORB,L1,N1)
C
C     THIS ROUTINE READS THE LOCALIZED NEO-HF ELECTRONIC AND NUCLEAR
C     ORBITALS AND WRITES THEM TO THE DAF FILE FOR LATER USE IN A
C     NEONCI CALCULATION
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      PARAMETER (MXATM=2000)
      PARAMETER (MXAO=8192)
C
      LOGICAL GOPARR,DSKWRK,MASWRK
C
      DIMENSION TIMSTR(3)
C
      COMMON /FMCOM / X(1)
      COMMON /CONV  / DENTOL,EN,ETOT,EHF,EHF0,DIFF,ITER,ICALP,ICBET
      COMMON /IOFILE/ IR,IW,IP,IS,IJKT,IDAF,NAV,IODA(950)
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
      COMMON /RUNLAB/ TITLE(10),ANAM(MXATM),BNAM(MXATM),BFLAB(MXAO)
C
      MPUNCH = 2
      L3 = L1*L1
      N3 = N1*N1
C
C     --- GET MEMORY TO READ EIGENVECTORS ---
C
      CALL VALFM(LOADFM)
      LEMO  = LOADFM + 1
      LNMO  = LEMO   + L3
      LAST  = LNMO   + N3
      NEED  = LAST - LOADFM -1
      CALL GETFM(NEED)
C
C     --- PUNCH LEFT LOCALIZED ORBITALS TO DAT FILE ---
C
      WRITE(IW,6000)
      CALL DAREAD(IDAF,IODA,X(LEMO),L3,15,0)
      CALL DAREAD(IDAF,IODA,X(LNMO),N3,444,0)
      IF(LOCORB.EQ.1) THEN
        IF(MPUNCH.GT.0 .AND. MASWRK) THEN
          CALL TMDATE(TIMSTR)
          WRITE(IP,7000) TIMSTR,TITLE,ETOT,EN,ITER
          WRITE(IP, FMT='(8H $LOCEL1)')
          CALL PUSQL(X(LEMO),L1,L1,L1)
          WRITE(IP, FMT='(5H $END)')
C
          CALL TMDATE(TIMSTR)
          WRITE(IP,8000) TIMSTR,TITLE,ETOT,EN,ITER
          WRITE(IP, FMT='(8H $LOCPR1)')
          CALL PUSQL(X(LNMO),N1,N1,N1)
          WRITE(IP, FMT='(5H $END)')
        END IF
      ELSE
          CALL TMDATE(TIMSTR)
          WRITE(IP,7100) TIMSTR,TITLE,ETOT,EN,ITER
          WRITE(IP, FMT='(8H $LOCEL2)')
          CALL PUSQL(X(LEMO),L1,L1,L1)
          WRITE(IP, FMT='(5H $END)')
C
          CALL TMDATE(TIMSTR)
          WRITE(IP,8100) TIMSTR,TITLE,ETOT,EN,ITER
          WRITE(IP, FMT='(8H $LOCPR2)')
          CALL PUSQL(X(LNMO),N1,N1,N1)
          WRITE(IP, FMT='(5H $END)')
      ENDIF
      CALL RETFM(NEED)
C
      RETURN
C
 6000 FORMAT(/1X,' THE LOCALIZED ELECTRON AND PROTON MOS WILL BE',/
     *        1X,' WRITTEN TO THE DAF SO THEY MAY BE LATER USED',/
     *        1X,' FOR A NEO-NCI CALCULATION',/)
 7000 FORMAT('--- LEFT LOCALIZED ELEC ORBITALS --- GENERATED AT ',
     *       3A8/10A8/
     *       'E(RHF)=',F20.10,', E(NUC)=',F16.10,',',I5,' ITERS')
 7100 FORMAT('--- RIGHT LOCALIZED ELEC ORBITALS --- GENERATED AT ',
     *       3A8/10A8/
     *       'E(RHF)=',F20.10,', E(NUC)=',F16.10,',',I5,' ITERS')
 8000 FORMAT('--- LEFT LOCALIZED NUC ORBITALS --- GENERATED AT ',
     *       3A8/10A8/
     *       'E(RHF)=',F20.10,', E(NUC)=',F16.10,',',I5,' ITERS')
 8100 FORMAT('--- RIGHT LOCALIZED NUC ORBITALS --- GENERATED AT ',
     *       3A8/10A8/
     *       'E(RHF)=',F20.10,', E(NUC)=',F16.10,',',I5,' ITERS')
      END
C
C*MODULE NEONCI   *DECK BASLOC
      SUBROUTINE BASLOC(INC,NATOM)
C
C     THIS ROUTINE DETERMINES HOW MANY ELECTRONIC BASIS FUNCTIONS ARE
C     ON THE SECOND BASIS FUNCTION CENTER
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
      PARAMETER (MXSH=5000,MXGTOT=20000)
      COMMON /NSHEL / EX(MXGTOT),CS(MXGTOT),CP(MXGTOT),CD(MXGTOT),
     *                CF(MXGTOT),CG(MXGTOT),CH(MXGTOT),CI(MXGTOT),
     *                KSTART(MXSH),KATOM(MXSH),KTYPE(MXSH),KNG(MXSH),
     *                KLOC(MXSH),KMIN(MXSH),KMAX(MXSH),NSHELL
      COMMON /NUCMOI/ NUNIQN,IUNIQN(20),IUNIQT(20),NQMNUC,IQMNUC(20),
     *                IQNTYP(20),NUMNB,NUCST,NAUXNB,IAUXNB(20),NUMULT,
     *                NNA,NNB,NTAUXB
C
      LAT = NATOM
      ICOUNT = 0
      DO I = 1,NSHELL
        IF(KATOM(I).EQ.LAT) THEN
             IFUNC = KTYPE(I)
           IF (IFUNC .EQ. 1) ICOUNT = ICOUNT + 1
           IF (IFUNC .EQ. 2) ICOUNT = ICOUNT + 3
           IF (IFUNC .EQ. 3) ICOUNT = ICOUNT + 9
        ENDIF
      END DO
      INC = ICOUNT*NQMNUC
C
      RETURN
      END
C
C*MODULE NEONCI   *DECK ORBRECALL1
      SUBROUTINE ORBRECALL1(VEC,ORB,NI,N,IORB,IPART)
C
C     THIS ROUTINE READS THE STORED SINGLE CENTER ELECTRONIC OR NUCLEAR
C     MOS FROM FILE AND ASSIGNS THEM TO THE ARRAY ORB
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      PARAMETER(MXATM=2000)
C
      CHARACTER*8 STVEC
C
      LOGICAL GOPARR,DSKWRK,MASWRK
C
      COMMON /FMCOM / X(1)
      COMMON /IOFILE/ IR,IW,IP,IS,IJKT,IDAF,NAV,IODA(950)
      COMMON /INFOA / NAT,ICH,MUL,NUM,NQMT,NE,NA,NB,
     *                ZAN(MXATM),C(3,MXATM),IAN(MXATM)
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
C
      DIMENSION ORB(N,N),VEC(NI,NI)
C
      IZERO = 0
      IONE  = 1
      ITWO  = 2
C
      CALL VCLR(ORB,1,N*N)
C
      IF(IPART.EQ.1) THEN
C
C      ELECTRONIC ORBITALS
C
        IF(IORB.EQ.1) THEN
          STVEC = ' $LOCEL1'
        ELSE
          STVEC = ' $LOCEL2'
        ENDIF
C
C    NOW GET ORTHONORMAL TRANSFORMATION MATRIX FOR ORBITALS
C
         NUM1 = N
         NUM2 = (N*N+N)/2
         NUM3 = N*N
         LOADFM=0
         CALL VALFM(LOADFM)
C
         LS    = LOADFM + 1
         LQ    = LS     + NUM2
         LEIG  = LQ     + NUM3
         LIWRK = LEIG   + NUM1
         LSCR  = LIWRK  + NUM1
         LAST  = LSCR   + NUM1*8
         NEEDEXORB = LAST - LOADFM - 1
         CALL GETFM(NEEDEXORB)
C
         CALL DAREAD(IDAF,IODA,X(LS),NUM2,NFIL,0)
         CALL QMTSYM(X(LS),VEC,X(LQ),X(LEIG),X(LSCR),X(LIWRK),L0,NUM1,
     *               NUM2,NUM3,.FALSE.)
C
      ELSE
C
C      PROTON ORBITALS
C
        IF(IORB.EQ.1) THEN
          STVEC = ' $LOCPR1'
        ELSE
          STVEC = ' $LOCPR2'
        ENDIF
      ENDIF
C
      NSTM=0
      NUM1=NI
C
C     FOR NOW NMOS IS SET TO ALL ORBITALS
C
      NMOS=NI
      IF(NMOS.EQ.0) CALL ABRT
C
C     --- POSITION THE INPUT FILE TO STVEC ---
C
      CALL SEQREW(IR)
      JJ = 0
      ICC = 0
      CALL FNDGRP(IR,STVEC,IEOF)
      IF (IEOF.EQ.0) GO TO 200
C
C     --- STVEC CARD NOT FOUND, KILL JOB ---
C
      IF (MASWRK) WRITE(IW,9020) STVEC
      CALL ABRT
      STOP
C
C     --- STVEC CARD FOUND, READ THE ORBITALS ---
C
  200 CONTINUE
C
      IF (MASWRK) THEN
         DO 280 J = 1,NMOS
            IMAX = 0
            IC = 0
  240       CONTINUE
               IMIN = IMAX+1
               IMAX = IMAX+5
               IC = IC+1
               IF(IMAX .GT. NUM1) IMAX = NUM1
               READ(IR,9040,END=300,ERR=300) JJ,ICC,(VEC(I,J),
     *                                  I=IMIN+NSTM,IMAX+NSTM)
               MODJ=MOD(J,100)
               IF(JJ.EQ.MODJ . AND.  ICC.EQ.IC) GO TO 260
                  WRITE(IW,9060) J,IC,JJ,ICC
                  IF (GOPARR) CALL DDI_BCAST(350,'I',ITWO,1,MASTER)
                  CALL ABRT
                  STOP
  260          CONTINUE
            IF(IMAX .LT. NUM1) GO TO 240
  280    CONTINUE
         IF (GOPARR) CALL DDI_BCAST(350,'I',IZERO,1,MASTER)
C
      ELSE
C        --- SLAVE WORK ---
C
         IF (GOPARR) CALL DDI_BCAST(350,'I',IEND,1,MASTER)
         IF (IEND.EQ.IONE) GO TO 300
         IF (IEND.EQ.ITWO) CALL ABRT
      END IF
C
C       --- NOW COPY SINGLE CENTER ORBS TO TWO CENTER ORBITAL ARRAY ---
C
       DO I = 1,NI
         DO J = 1,NI
           ORB(I,J) = VEC(I,J)
         END DO
       END DO
C
C        --- FOR RIGHT LOCALIZED ELEC ---
C
      IF(IPART.EQ.1 .AND. IORB.EQ.2)  THEN
        CALL BASLOC(INC,NAT)
C   HERE THE NUMBER OF BF FOR ONE BFC (INC) IS DOUBLED
C   AND SUBTRACTED FROM THE TOTAL NUMBER OF B.F. AND
C   1 IS ADDED SO THE STARTING POSITION IS AT THE START OF
C   BFC TWO FOR QM NUC 1.
        K = N-2*INC+1
        DO i = N-INC+1,N
         DO j = 1,N
           ORB(I,J) = ORB(K,J)
           ORB(K,J) = 0.0D+00
         ENDDO
         K = K+1
        ENDDO
C
C        --- FOR RIGHT LOCALIZED PROTON ---
C
      ELSE IF(IPART.EQ.2 .AND. IORB.EQ.2) THEN
         K=1
       DO I = NI+1,N
        DO J = 1,N
          ORB(I,J) = ORB(K,J)
          ORB(K,J) = 0.0D+00
        ENDDO
        K = K+1
       ENDDO
      ENDIF
C
      IF(IPART.EQ.1) THEN
C
C            ORTHONORMALIZE THE ORBITALS
C            BACK-TRANSFORM THE MO'S
C
        TOLZ=1.0D-08
        TOLE=1.0D-05
C
        CALL CLENMO(ORB,NUM1,NUM1,TOLZ,TOLE,IW,.FALSE.)
        CALL ORTHO(X(LQ),X(LS),ORB,X(LSCR),MIN(NUM1,L0),L0,NUM1,
     *             NUM2,NUM1)
        CALL TFSQB(ORB,X(LQ),X(LSCR),L0,NUM1,NUM1)
        CALL RETFM(NEEDEXORB)
      ENDIF
C
      RETURN
C
C     --- PREMATURE END OF $VEC INPUT ENCOUNTERED ---
C
  300 CONTINUE
      IF (MASWRK) THEN
         IF (GOPARR) CALL DDI_BCAST(350,'I',IONE,1,MASTER)
         WRITE(IW,9100) J,IC
      END IF
      CALL ABRT
      STOP
C
 9020 FORMAT(1X,' NO ',A8,' INPUT FOUND - STOP' )
 9040 FORMAT(I2,I3,5E15.8)
 9060 FORMAT(' *** ERROR READING THE MO COEFFICIENTS AT',2I4/
     *       '     THE INPUT VALUES WERE',2I4)
 9100 FORMAT(' *** ERROR: PREMATURE END OF ORBITAL INPUT ENCOUNTERED'/
     *       '            LOOKING FOR ORBITAL',I4,' ELEMENT',I4)
C
      END
C
C*MODULE NEONCI   *DECK ORBRECALL2
      SUBROUTINE ORBRECALL2(VEC,N,IORB,IPART)
C
C     THIS ROUTINE READS THE STORED TWO CENTER ELECTRONIC OR NUCLEAR
C     MOS FROM FILE AND ASSIGNS THEM TO THE ARRAY ORB
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      CHARACTER*8 STVEC
C
      LOGICAL GOPARR,DSKWRK,MASWRK
C
      COMMON /FMCOM / X(1)
      COMMON /IOFILE/ IR,IW,IP,IS,IJKT,IDAF,NAV,IODA(950)
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
C
      DIMENSION VEC(N,N)
C
      IZERO = 0
      IONE  = 1
      ITWO  = 2
C
      IF(IPART.EQ.1) THEN
C
C      ELECTRONIC ORBITALS
C
      NFIL = 12
C
        IF(IORB.EQ.1) THEN
          STVEC = ' $LOCEL1'
        ELSE
          STVEC = ' $LOCEL2'
        ENDIF
C
C    NOW GET ORTHONORMAL TRANSFORMATION MATRIX FOR ELEC ORBITALS
C
         NUM1 = N
         NUM2 = (N*N+N)/2
         NUM3 = N*N
         LOADFM=0
         CALL VALFM(LOADFM)
C
         LS    = LOADFM + 1
         LQ    = LS + NUM2
         LEIG  = LQ + NUM3
         LIWRK = LEIG + NUM1
         LSCR  = LIWRK + NUM1
         LAST  = LSCR + NUM1*8
         NEEDEXORB = LAST - LOADFM - 1
         CALL GETFM(NEEDEXORB)
C
         CALL DAREAD(IDAF,IODA,X(LS),NUM2,NFIL,0)
         CALL QMTSYM(X(LS),VEC,X(LQ),X(LEIG),X(LSCR),X(LIWRK),L0,NUM1,
     *               NUM2,NUM3,.FALSE.)
C
      ELSE
C
C      PROTON ORBITALS
C
      NFIL = 441
C
        IF(IORB.EQ.1) THEN
          STVEC = ' $LOCPR1'
        ELSE
          STVEC = ' $LOCPR2'
        ENDIF
      ENDIF
C
      NSTM=0
      NUM1=N
C
C     FOR NOW SET NMOS TO TOTAL NUMBER OF ORBITALS
C
      NMOS=N
      IF(NMOS.EQ.0) CALL ABRT
C
C     --- POSITION THE INPUT FILE TO STVEC ---
C
      CALL SEQREW(IR)
      JJ = 0
      ICC = 0
      CALL FNDGRP(IR,STVEC,IEOF)
      IF (IEOF.EQ.0) GO TO 200
C
C     --- STVEC CARD NOT FOUND, KILL JOB ---
C
      IF (MASWRK) WRITE(IW,9020) STVEC
      CALL ABRT
      STOP
C
C     --- STVEC CARD FOUND, READ THE ORBITALS ---
C
  200 CONTINUE
C
      IF (MASWRK) THEN
         DO 280 J = 1,NMOS
            IMAX = 0
            IC = 0
  240       CONTINUE
               IMIN = IMAX+1
               IMAX = IMAX+5
               IC = IC+1
               IF(IMAX .GT. NUM1) IMAX = NUM1
               READ(IR,9040,END=300,ERR=300) JJ,ICC,(VEC(I,J),
     *                                  I=IMIN+NSTM,IMAX+NSTM)
               MODJ=MOD(J,100)
               IF(JJ.EQ.MODJ . AND.  ICC.EQ.IC) GO TO 260
                  WRITE(IW,9060) J,IC,JJ,ICC
                  IF (GOPARR) CALL DDI_BCAST(350,'I',ITWO,1,MASTER)
                  CALL ABRT
                  STOP
  260          CONTINUE
            IF(IMAX .LT. NUM1) GO TO 240
  280    CONTINUE
         IF (GOPARR) CALL DDI_BCAST(350,'I',IZERO,1,MASTER)
C
      ELSE
C        --- SLAVE WORK ---
C
         IF (GOPARR) CALL DDI_BCAST(350,'I',IEND,1,MASTER)
         IF (IEND.EQ.IONE) GO TO 300
         IF (IEND.EQ.ITWO) CALL ABRT
      END IF
C
      IF(IPART.EQ.1) THEN
C
C            ORTHONORMALIZE THE ORBITALS
C            BACK-TRANSFORM THE MO'S
C
        TOLZ=1.0D-08
        TOLE=1.0D-05
C
        CALL CLENMO(VEC,NUM1,NMOS,TOLZ,TOLE,IW,.FALSE.)
        CALL ORTHO(X(LQ),X(LS),VEC,X(LSCR),MIN(NMOS,L0),L0,NUM1,
     *             NUM2,NUM1)
        CALL TFSQB(VEC,X(LQ),X(LSCR),L0,NUM1,NUM1)
        CALL RETFM(NEEDEXORB)
      ENDIF
C
C
      RETURN
C
C     --- PREMATURE END OF $VEC INPUT ENCOUNTERED ---
C
  300 CONTINUE
      IF (MASWRK) THEN
         IF (GOPARR) CALL DDI_BCAST(350,'I',IONE,1,MASTER)
         WRITE(IW,9100) J,IC
      END IF
      CALL ABRT
      STOP
C
 9020 FORMAT(1X,' NO ',A8,' INPUT FOUND - STOP' )
 9040 FORMAT(I2,I3,5E15.8)
 9060 FORMAT(' *** ERROR READING THE MO COEFFICIENTS AT',2I4/
     *       '     THE INPUT VALUES WERE',2I4)
 9100 FORMAT(' *** ERROR: PREMATURE END OF ORBITAL INPUT ENCOUNTERED'/
     *       '            LOOKING FOR ORBITAL',I4,' ELEMENT',I4)
      END
C
C*MODULE NEONCI   *DECK NEONCI_RUN
      SUBROUTINE NEONCI_RUN(LCI,LCIT,LCI2,NCIRUN,NCISYM,
     *           NUCFCI,CI,CIFIN1D,WRK1,TOTOVLP,SDIAG,
     *           SEIGEN,SWRK,ISWRK,SVEC,SX1,SX2,CEIGEN,
     *           WRK2,IWRK2,VINDIAG)
C
C     THIS ROUTINE CALCULATES THE NEONCI MATRIX ELEMENTS AND THEN SOLVES
C     THE NEONCI SECULAR EQUATION
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      LOGICAL PACK2E,DBUG
      LOGICAL NCISYM,NUCFCI
      LOGICAL ISW
C
      DOUBLE PRECISION NCIRUN
C
      PARAMETER (MXATM=2000,MXRT=100)
      PARAMETER (TWO=2.0D+00)
C
      COMMON /ENRGYS/ ENUCR,EELCT,ETOT,SZ,SZZ,ECORE,ESCF,EERD,E1,E2,
     *                VEN,VEE,EPOT,EKIN,ESTATE(MXRT),STATN,EDFT(2)
      COMMON /FMCOM / X(1)
      COMMON /INFOA / NAT,ICH,MUL,NUM,NQMT,NE,NA,NB,
     *                  ZAN(MXATM),C(3,MXATM),IAN(MXATM)
      COMMON /INTFIL/ NINTMX,NHEX,NTUPL,PACK2E,INTTYP,IGRDTYP
      COMMON /IOFILE/ IR,IW,IP,IS,IJKT,IDAF,NAV,IODA(950)
      COMMON /MEMINF/ MEMSET
      COMMON /NCIFLG/ NUMCNT,NCITRN,NCIPRT,NUMNCI,NUMEXE,NUMEXP,
     *                IEXE(40),IEXP(40),CUTNCI
      COMMON /NUCMOI/ NUNIQN,IUNIQN(20),IUNIQT(20),NQMNUC,IQMNUC(20),
     *                  IQNTYP(20),NUMNB,NUCST,NAUXNB,IAUXNB(20),NUMULT,
     *                  NNA,NNB,NTAUXB
      COMMON /RUNOPT/ RUNTYP,EXETYP,NEVALS,NGLEVL,NHLEVL
C
      DIMENSION CI(LCIT),CIFIN1D(LCIT),WRK1(LCI)
      DIMENSION TOTOVLP(LCIT),SDIAG(LCI,LCI)
      DIMENSION SEIGEN(LCI),SWRK(LCI,8),ISWRK(LCI),SVEC(LCI,LCI),
     *          SX1(LCI,LCI),SX2(LCI,LCI)
      DIMENSION CEIGEN(LCI),WRK2(LCI,8),IWRK2(LCI),VINDIAG(LCI,LCI)
C
      DATA DEBUG/8HDEBUG   /
      DATA CHECK/8HCHECK   /

C
C        --- DEFINE ARRAY LENGTH SIZES ---
C
      L1 = NUM
      L2 = (L1*L1+L1)/2
      L3 = L1 * L1
      N1 = NUMNB
      N3 = N1*N1
      N2 = (N1*N1+N1)/2
C
      DBUG = .FALSE.
      IF(NCIRUN.NE.DEBUG) THEN
         IF(EXETYP.EQ.DEBUG) DBUG = .TRUE.
      ELSE
         DBUG = .TRUE.
      ENDIF
C
C        --- IF NCIRUN=CHECK PRINT OUT MEMORY REQUIREMENTS AND RETURN ---
C
      IF(NCIRUN.EQ.CHECK .OR. EXETYP.EQ.CHECK) THEN
c      IMEMMAT = LCIT*3 +LCI*5 + 16*LCI + LCI2*5 + LCIT*12
      IMEMMAT = LCIT*15 +LCI*21 + LCI2*5
C
      IMEM = 0
      IMEM = IMEM + L3*2
      IMEM = IMEM + N3*2
      IMEM = IMEM + L2*L1*NA
      IMEM = IMEM + NINTMX
      IMEM = IMEM + NINTMX
      IMEM = IMEM + L1*NA
      IMEM = IMEM + L3
      IMEM = IMEM + L3
      IMEM = IMEM + NA*NA*NA
      IMEM = IMEM + NA*NA*NA
      IMEM = IMEM + L2
      IMEM = IMEM + L1*NA
      IMEM = IMEM + IMEMMAT
C
      IMEM2 = 0
      IMEM2 = IMEM2 + L3*2
      IMEM2 = IMEM2 + N3*2
      IMEM2 = IMEM2 + L2*NA*NA
      IMEM2 = IMEM2 + NINTMX
      IMEM2 = IMEM2 + NINTMX
      IMEM2 = IMEM2 + L1*NA
      IMEM2 = IMEM2 + L3
      IMEM2 = IMEM2 + L3
      IMEM2 = IMEM2 + NA*NA*NA
      IMEM2 = IMEM2 + NA*NA*NA
      IMEM2 = IMEM2 + NA*NA
      IMEM2 = IMEM2 + IMEMMAT
C
      NDISK = NA*NA*NA*NA
      WRITE(IW,900) IMEM,IMEM2,NDISK
      RETURN
      ENDIF
C
C        --- GET MEMORY FOR ARRAYS IN ROUTINE ---
C
      CALL GOTFM(NGOTMX)
      MEMSET = NGOTMX
C
      IF(DBUG) THEN
      PMEMVAL = ((MEMSET-NGOTMX)*8.0D+0)/1.0D+6
      WRITE(IW,*)
      WRITE(IW,*) "MEMTAGOPEN                                         "
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) " MEMORY AVAILABLE AT START OF NEONCI:",
     *            (NGOTMX*8.0D+0)/1.0D+6," MB"
      WRITE(IW,*) "                                                    "
      WRITE(IW,*) "================================================="
      WRITE(IW,*) "        **** VALUES FOR DIMENSIONS IN WORDS****  "
      WRITE(IW,*) "-------------------------------------------------"
      WRITE(IW,*) "          ELECTRONIC "
      WRITE(IW,*) " L1 EBF:                 ",L1
      WRITE(IW,*) " L2 (EBF*EBF+EBF)/2:     ",L2
      WRITE(IW,*) " L3 EBF*EBF:             ",L3
      WRITE(IW,*) " L4 ELEC*ELEC:           ",NE*NE
      WRITE(IW,*) "          NUCLEAR    "
      WRITE(IW,*) " N1 NBF:                 ",N1
      WRITE(IW,*) " N2 (NBF*NBF+NBF)/2:     ",N2
      WRITE(IW,*) " N3 NBF*NBF:             ",N3
      WRITE(IW,*) " N4 NUC*NUC:             ",N3
      WRITE(IW,*) "================================================="
      WRITE(IW,*) "MEMTAGEND                                           "
      WRITE(IW,*) "MEMSTAT:  ", PMEMVAL
      WRITE(IW,*)
      CALL FLSHBF(IW)
      ENDIF
C
      LOADFM=0
      CALL VALFM(LOADFM)
      LEMO1     = LOADFM + 1
      LEMO2     = LEMO1  + L3
      NMO1      = LEMO2  + L3
      NMO2      = NMO1   + N3
      LE1       = NMO2   + N3
      LE2       = LE1    + L3
      NP1       = LE2    + L3
      NP2       = NP1    + N3
      LAST      = NP2    + N3
      NEEDA = LAST - LOADFM - 1
      CALL GETFM(NEEDA)
C
      CALL GOTFM(NMEM1)
C
      IF(DBUG) THEN
      PMEMVAL = ((MEMSET-NMEM1)*8.0D+0)/1.0D+6
      WRITE(IW,*) "                                                    "
      WRITE(IW,*) "MEMTAGOPEN                                          "
      WRITE(IW,*) "================================================="
      WRITE(IW,*) "     **** MEMORY NEEDED FOR MOS in WORDS****     "
      WRITE(IW,*) "-------------------------------------------------"
      WRITE(IW,*) "          ELECTRONIC  "
      WRITE(IW,'(A15,I10,A6)') " ELEC MO 1:     ",L3, " WORDS"
      WRITE(IW,'(A15,I10,A6)') " ELEC MO 2:     ",L3, " WORDS"
      WRITE(IW,*) "          NUCLEAR     "
      WRITE(IW,'(A15,I10,A6)') " NUC  MO 1:     ",N3, " WORDS"
      WRITE(IW,'(A15,I10,A6)') " NUC  MO 2:     ",N3, " WORDS"
      WRITE(IW,*) "                                                    "
      WRITE(IW,*) "MEMORY LEFT AFTER CALL FOR MO NEEDS:",
     *            (NMEM1*8.0D+0)/1.0D+6," MB"
      WRITE(IW,*) "================================================="
      WRITE(IW,*) "MEMTAGEND                                           "
      WRITE(IW,*) "MEMSTAT:  ", PMEMVAL
      WRITE(IW,*)
      CALL FLSHBF(IW)
      ENDIF
C
C        --- IF A SINGLE CENTER CALC DETERMINE NUM OF EXTRA B.F. ---
C
      IF(NUMCNT.EQ.1) THEN
         CALL BASLOC(INC,NAT)
c           note by MWS.  The next two were spelled LEX1/NEX1, and
c           L10/N10 in the memory allocations below were undefined.
c           Presumably the next two are really L10/N10.
         L10 = L1-INC ! SIZE OF ELEC ORB VEC FROM SINGLE CENTER CALC
         N10 = N1/2   ! SIZE OF NUCLEAR ORB VEC FROM SINGLE CENTER CALC
C
C        ---  READ IN SINGLE CENTER LOCALIZED ORBITALS ---
C
         LOADFM=0
         CALL VALFM(LOADFM)
C
         LEMOI1    = LOADFM + 1
         LEMOI2    = LEMOI1 + L10*L10
         NMOI1     = LEMOI2 + L10*L10
         NMOI2     = NMOI1  + N10*N10
         LAST      = NMOI2  + N10*N10
         NEEDEXORB = LAST - LOADFM - 1
         CALL GETFM(NEEDEXORB)
C
         IORB = 1
         IPART = 1
         CALL ORBRECALL1(X(LEMOI1),X(LEMO1),L10,L1,IORB,IPART)
         IPART = 2
         CALL ORBRECALL1(X(NMOI1),X(NMO1),N10,N1,IORB,IPART)
         IORB = 2
         IPART = 1
         CALL ORBRECALL1(X(LEMOI2),X(LEMO2),L10,L1,IORB,IPART)
         IPART = 2
         CALL ORBRECALL1(X(NMOI2),X(NMO2),N10,N1,IORB,IPART)
C
      ELSE
C
C        --- READ IN TWO CENTER LOCALIZED ORBITALS  ---
C
         IORB = 1
         IPART = 1
         CALL ORBRECALL2(X(LEMO1),L1,IORB,IPART)
         IPART = 2
         CALL ORBRECALL2(X(NMO1),N1,IORB,IPART)
         IORB = 2
         IPART = 1
         CALL ORBRECALL2(X(LEMO2),L1,IORB,IPART)
         IPART = 2
         CALL ORBRECALL2(X(NMO2),N1,IORB,IPART)
C
      ENDIF
C
      IF (NCIPRT.NE.1) THEN
      WRITE(IW,9001)
      CALL PRSQ(X(LEMO1),L1,L1,L1)
      WRITE(IW,9002)
      CALL PRSQ(X(NMO1),N1,N1,N1)
      WRITE(IW,9003)
      CALL PRSQ(X(LEMO2),L1,L1,L1)
      WRITE(IW,9004)
      CALL PRSQ(X(NMO2),N1,N1,N1)
      ENDIF
C
      WRITE(IW,9309)
C
      CALL VALFM(LOADFM)
      IOCCE = LOADFM + 1
      IVIRE = IOCCE + LCIT
      IOCCP = IVIRE + LCIT
      IVIRP = IOCCP + LCIT
      IROP  = IVIRP + LCIT
      IRVP  = IROP  + LCIT
      ICOP  = IRVP  + LCIT
      ICVP  = ICOP  + LCIT
      IROE  = ICVP  + LCIT
      IRVE  = IROE  + LCIT
      ICOE  = IRVE  + LCIT
      ICVE  = ICOE  + LCIT
      LAST  = ICVE  + LCIT
      NEEDIMAP = LAST - LOADFM -1
      CALL GETFM(NEEDIMAP)

      CALL GETIMAP(X(IOCCE),X(IVIRE),X(IOCCP),X(IVIRP),X(IROP),X(IRVP),
     *             X(ICOP),X(ICVP),X(IROE),X(IRVE),X(ICOE),X(ICVE),
     *             LCI,LCIT,NUCFCI,DBUG)
C
      ILEN = LCI
      ISW  = .TRUE.
      IDEC = 1
      LDX  = 1
C
C       --- GET THE DETERMINANT FOR ELEMENT I ---
C
      DO IDX =1,LCIT
        CALL GETDET(X(LEMO1),X(LEMO2),X(NMO1),X(NMO2),X(LE1),X(LE2),
     *            X(NP1),X(NP2),LCI,LCIT,L1,N1,IDX,ISW,ILEN,IDEC,LDX,
     *            X(IROP),X(IRVP),X(ICOP),X(ICVP),X(IROE),X(IRVE),
     *            X(ICOE),X(ICVE))
C
       IF(NCIRUN.EQ. DEBUG) THEN
C
C       --- PRINT THE ORBITALS USED FOR NCI ELEMENT I ---
C
      WRITE(IW,9200) IDX
      CALL PRSQ(X(LE1),L1,L1,L1)
      WRITE(IW,9210) IDX
      CALL PRSQ(X(NP1),N1,N1,N1)
      WRITE(IW,9220) IDX
      CALL PRSQ(X(LE2),L1,L1,L1)
      WRITE(IW,9230) IDX
      CALL PRSQ(X(NP2),N1,N1,N1)
      ENDIF
C
C        --- IF NCISYM OPTION SELECTED SKIP CALC OF DIAGONAL ELEMENTS ---
C
        IF (NCISYM .AND. LCI .EQ. 2)  THEN
          IF(IDX.NE.2) THEN
            CI(IDX) = ETOT
            CI(IDX) = ETOT
            TOTOVLP(IDX) = 1.0D+0
            TOTOVLP(IDX) = 1.0D+0
            GOTO 700
          ENDIF
        ENDIF
C
C        --- CONSTRUCT NCI ELEMENT IDX ---
C
      CALL NCIELEM(VALUE,X(LE1),X(LE2),X(NP1),X(NP2),TOTOV2,
     *             PROVLP,ELECOVLP,ENUCR,L1,N1,DBUG,NUCFCI)
      CI(IDX) = VALUE
      TOTOVLP(IDX) = TOTOV2
      IF(NCIPRT.EQ.3 .AND. IDX .EQ. 2) THEN
        SE12 = DABS(ELECOVLP)
        SP12 = DABS(PROVLP)
      ENDIF
C
  700 CONTINUE
      ENDDO
C
C        --- MAKE LABEL ARRAY NCI VECTORS ---
C
      CALL VALFM(LOADFM)
      LABEO = LOADFM + 1
      LABEV = LABEO  + LCI
      LABPO = LABEV  + LCI
      LABPV = LABPO  + LCI
      LAST  = LABPV  + LCI
      NEEDLAB = LAST - LOADFM -1
      CALL GETFM(NEEDLAB)
C
      DO I = 1,LCI
         X(LABEO +I -1) = X(IROE +I -1)
         X(LABEV +I -1) = X(IRVE +I -1)
C
         X(LABPO +I -1) = X(IROP +I -1)
         X(LABPV +I -1) = X(IRVP +I -1)
      ENDDO
C
      CALL RETFM(NEEDIMAP)
C
      WRITE(IW,9310)
      CALL TIMIT(1)
C
C       --- CONVERT FROM ROW STORED LOWER TRIANGULAR ---
C       --- TO COLUMN STORED LOWER TRIANGULAR MATRIX ---
C
      CALL LTRI_RTOC(CI,LCIT,LCI)
      CALL LTRI_RTOC(TOTOVLP,LCIT,LCI)

C
C        --- NOW SOLVE NCI EQUATION ---
C
      WRITE(IW,*)
      WRITE(IW,9410)
      CALL PRTRI(CI,LCI)
C
      WRITE(IW,9510)
      CALL PRTRI(TOTOVLP,LCI)
C
C        --- DIAGONALIZE OVERLAP MATRIX ---
C
       CALL VCLR(SDIAG,1,LCI2)
       CALL GLDIAG(LCI,LCI,LCI,TOTOVLP,SWRK,SEIGEN,SVEC,
     *             IERR,ISWRK)
C
C        --- FORM S-1/2 FROM DIAGONALIZED S ---
C
       DO I = 1,LCI
          J = I
                 SDIAG(I,J) = 1.00D+0/DSQRT(SEIGEN(I))
       ENDDO
C
       CALL MRARBR(SVEC,LCI,LCI,LCI,SDIAG,LCI,LCI,SX1,LCI)
       CALL TRPOSQ(SVEC,LCI)
       CALL MRARBR(SX1,LCI,LCI,LCI,SVEC,LCI,LCI,SX2,LCI)
C
C        --- FORM THE ORTHOGONALIZED NCI MATRIX ---
C
       CALL TFTRI(CIFIN1D,CI,SX2,WRK1,LCI,LCI,LCI)
C
      WRITE(IW,9610)
      CALL PRTRI(CIFIN1D,LCI)
C
C        --- DIAGONALIZE NCI MATRIX ---
C
      CALL GLDIAG(LCI,LCI,LCI,CIFIN1D,WRK2,CEIGEN,VINDIAG,
     *             IERR,IWRK2)
C
      WRITE(IW,9615)
      DO I = 1,NUMNCI
         WRITE(IW,9700) I,CEIGEN(I)
         CALL PREIGV(VINDIAG,LCI,LCI,I,X(LABEO),X(LABEV),X(LABPO),
     *               X(LABPV))
      ENDDO
      CALL RETFM(NEEDLAB)
C
      DELTA_E = (CEIGEN(2)-CEIGEN(1))*2.1947463D+5
      WRITE(IW,9710) DELTA_E
      WRITE(IW,9910) DELTA_E/TWO
      IF(NCIPRT.EQ.3) THEN
         WRITE(IW,9997) CEIGEN(1), CEIGEN(2),CI(1),CI(2),TOTOV2
         WRITE(IW,9998) SE12,SP12
      ENDIF
      WRITE(IW,9999)
C
      CALL RETFM(NEEDA)
C
      RETURN
C
  900 FORMAT(/1X,' MEMORY REQUIREMENTS FOR NEO-NCI CALCULATION'/
     *        2X,43('-')/
     *        1X,' MEMORY NEEDED FOR NEO-NCI CALC',
     *        ' USING 3-STEP TRANSFORMATION: ',I12, ' WORDS'/
     *        1X,' MEMORY NEEDED FOR NEO-NCI CALC',
     *        ' USING 2-STEP TRANSFORMATION: ',I12,' WORDS'/
     *        1X,' AMOUNT OF DISK SPACE NEEDED FOR MO TRANSFORMED',
     *        ' INTEGRALS:   ',I12,' WORDS'/)
 9001 FORMAT(//10X,28('-')/
     *        10X,'LEFT LOCALIZED ELECTRON ORBS'/
     *        10X,28('-'))
 9002 FORMAT(//10X,26('-')/
     *        10X,'LEFT LOCALIZED PROTON ORBS'/
     *        10X,26('-'))
 9003 FORMAT(//10X,29('-')/
     *        10X,'RIGHT LOCALIZED ELECTRON ORBS'/
     *        10X,29('-'))
 9004 FORMAT(//10X,27('-')/
     *        10X,'RIGHT LOCALIZED PROTON ORBS'/
     *        10X,27('-'))
 9200 FORMAT(//11X,18('-')/
     *        10X,'< ELECTRON VEC ',I4,' |'/
     *        11X,18('-'))
 9210 FORMAT(//11X,16('-')/
     *        10X,'< PROTON VEC ',I4,' |'/
     *        11X,16('-'))
 9220 FORMAT(//11X,18('-')/
     *        10X,'| ELECTRON VEC ',I4,' >'/
     *        11X,18('-'))
 9230 FORMAT(//11X,16('-')/
     *        10X,'| PROTON VEC ',I4,' >'/
     *        11X,16('-'))
 9309 FORMAT(/1X,'BEGIN FORMING THE NEO-NCI ',
     *           'HAMILTONIAN MATRIX ....')
 9310 FORMAT(/1X,' .... FINISHED FORMING THE NEO-NCI ',
     *           'HAMILTONIAN MATRIX '/)
 9410 FORMAT(/5X,15('-'),
     *       /5X,' NEO-NCI MATRIX',
     *       /5X,15('-'))
 9510 FORMAT(/5X,19('-'),
     *       /5X,' OVERLAP (S) MATRIX',
     *       /5X,19('-'))
 9610 FORMAT(/5X,30('-'),
     *       /5X,' ORTHOGONALIZED NEO-NCI MATRIX',
     *       /5X,30('-'))
 9615 FORMAT(//1X,' NEO-NCI EIGENVALUES AND EIGENVECTORS',
     *       /1X,37('-'))
 9700 FORMAT(/1X,'STATE = ',I2,' ENERGY =',F18.10)
 9710 FORMAT(//1X,' VIBRONIC SPLITTING =',F18.10,' cm-1')
 9910 FORMAT(/1X,' VIBRONIC COUPLING  =',F18.10,' cm-1')
 9997 FORMAT(/1X,' THE FOLLOWING OUTPUT IS FOR EXTRACTION BY A PERL',
     *        ' SCRIPT'//1X,' ENERG OUT',6F18.10)
 9998 FORMAT(/1X,' SE12 AND SP12 ',2F18.10)
 9999 FORMAT(/1X,' .... END OF NEO-NCI CALCULATION ....'/)
      END
C
C*MODULE NEONCI   *DECK GETDET
      SUBROUTINE GETDET(VECIN1,VECIN2,NVECIN1,NVECIN2,VEC1,VEC2,NVEC1,
     *                  NVEC2,NMT,NMT2,L1,N1,INDX,ISW,ILEN,IDEC,LDX,
     *                  IROP,IRVP,ICOP,ICVP,IROE,IRVE,ICOE,ICVE)
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      LOGICAL ISW
C
      DIMENSION VECIN1(L1,L1),VECIN2(L1,L1),NVECIN1(N1,N1),
     *          NVECIN2(N1,N1),VEC1(L1,L1),VEC2(L1,L1),NVEC1(N1,N1),
     *          NVEC2(N1,N1),IROP(NMT2),IRVP(NMT2),ICOP(NMT2),
     *          ICVP(NMT2),IROE(NMT2),IRVE(NMT2),ICOE(NMT2),ICVE(NMT2)
C
C     THIS ROUTINE GETS THE ROW AND COLUMN ELECTRONIC AND NUCLEAR
C     DETERMINANTS FOR NCI ELEMENT IDX BASED ON THE PROTON AND
C     ELECTRON ROW AND COLUMN MAPS FOR IDX THAT WERE PREVIOUSLY
C     GENERATED IN THE ROUTINE GETIMAP
C
      IEVE = MOD(LDX,2)
C
C     INDICIES OF SINGLE EXCITATION
C
      IROEX = IROE(INDX)
      IRVEX = IRVE(INDX)
      IROPX = IROP(INDX)
      IRVPX = IRVP(INDX)
      ICOEX = ICOE(INDX)
      ICVEX = ICVE(INDX)
      ICOPX = ICOP(INDX)
      ICVPX = ICVP(INDX)
C
C     ELECTRON DETERMINANT
C
C     <VEC| ROW
C
      IF(IEVE.NE.0) THEN
        CALL CPYVEC(VEC1,VECIN1,L1,IROEX,IRVEX)
        CALL CPYVEC(NVEC1,NVECIN1,N1,IROPX,IRVPX)
      ELSE
        CALL CPYVEC(VEC1,VECIN2,L1,IROEX,IRVEX)
        CALL CPYVEC(NVEC1,NVECIN2,N1,IROPX,IRVPX)
      ENDIF
        LDX = LDX +1
C
C     |VEC> COLUMN
C
      IF(ISW) THEN   ! THIS MAY BE INCORRECT FOR NMT GT 2
        CALL CPYVEC(VEC2,VECIN1,L1,ICOEX,ICVEX)
        CALL CPYVEC(NVEC2,NVECIN1,N1,ICOPX,ICVPX)
      ELSE
        CALL CPYVEC(VEC2,VECIN2,L1,ICOEX,ICVEX)
        CALL CPYVEC(NVEC2,NVECIN2,N1,ICOPX,ICVPX)
      ENDIF
C
      IF(INDX.EQ.ILEN) THEN
          ISW  = .NOT.ISW
          ILEN = ILEN + NMT -IDEC
          IDEC = IDEC + 1
          LDX  = IDEC
      ENDIF
C
      RETURN
      END
C
C*MODULE NEONCI   *DECK CPYVEC
      SUBROUTINE CPYVEC(V1,VTMP,L,IO,IV)
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      DIMENSION V1(L,L),VTMP(L,L)
C
      L3 = L*L
      CALL DCOPY(L3,VTMP,1,V1,1)
      DO I = 1,L
         V1(I,IO) = VTMP(I,IV)
         V1(I,IV) = VTMP(I,IO)
      ENDDO
C
      RETURN
      END
C*MODULE NEONCI   *DECK LTRI_RTOC
      SUBROUTINE LTRI_RTOC(A,NT,N)
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      COMMON /FMCOM / X(1)
C
      DIMENSION A(NT)
C
C     THIS ROUTINE INTERCHANGES THE 1D ARRAY OF A LOWER DIAGONAL MATRIX WHICH
C     IS STORED IN ORDER OF ROWS TO A LOWER DIAGONAL MATRIX WHICH IS STORED
C     IN ORDER OF COLUMNS
C     IN:  A     --> 1D ARRAY FOR A LOWER DIAGONAL MATRIX OF SIZE NT
C
      CALL VALFM(LOADFM)
      IA  = LOADFM +1
      LAST = IA + NT
      NEED = LAST - LOADFM -1
      CALL GETFM(NEED)
C
      CALL DCOPY(NT,A,1,X(IA),1)
C
      M = N
      L = 0
      K = 1
C
      DO I = 1, N
        DO J = 1,I
          A(K) = X(IA+I+L-1)
          M = M - 1
          L = L + M
          K = K+1
        ENDDO
          M = N
          L = 0
      ENDDO
C
      CALL RETFM(NEED)
C
      RETURN
      END
C
C*MODULE NEONCI   *DECK NCIELEM
      SUBROUTINE NCIELEM(VAL,EMO1,EMO2,PMO1,PMO2,TDET,PRODET,ELECDET,
     *                  ENUCR,L1,N1,DBUG,NUCFCI)
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      LOGICAL PACK2E,DBUG
      LOGICAL NUCFCI
C
      PARAMETER (MXATM=2000,MXAO=8192)
C
      COMMON /FMCOM / X(1)
      COMMON /IJPAIR/ IA(MXAO)
      COMMON /INFOA / NAT,ICH,MUL,NUM,NQMT,NE,NA,NB,
     *                ZAN(MXATM),C(3,MXATM),IAN(MXATM)
      COMMON /INTFIL/ NINTMX,NHEX,NTUPL,PACK2E,INTTYP,IGRDTYP
      COMMON /IOFILE/ IR,IW,IP,IS,IJKT,IDAF,NAV,IODA(950)
      COMMON /MEMINF/ MEMSET
      COMMON /NCIFLG/ NUMCNT,NCITRN,NCIPRT,NUMNCI,NUMEXE,NUMEXP,
     *                IEXE(40),IEXP(40),CUTNCI
      COMMON /NUCMOI/ NUNIQN,IUNIQN(20),IUNIQT(20),NQMNUC,IQMNUC(20),
     *                IQNTYP(20),NUMNB,NUCST,NAUXNB,IAUXNB(20),NUMULT,
     *                NNA,NNB,NTAUXB
C
      DIMENSION EMO1(L1,L1),EMO2(L1,L1),PMO1(N1,N1),PMO2(N1,N1)
C
      L2 = (L1*L1+L1)/2
      L3 = L1 * L1
      L4 = NE*NE
      N2 = (N1*N1+N1)/2
      N3 = N1*N1
      N4 = NQMNUC*NQMNUC
C
C        --- SET UP IA ARRAY FOR READING INTEGRALS  ---
C
      II = 0
      ILIM = MAX(N1,L1)
      DO 100 I = 1,ILIM
         IA(I) = II
         II = II + I
  100 CONTINUE
C
C        --- CALL AO TO MO TRANSFORMATION FOR ELEC AND NUC ---
C
      CALL TSECND(TIM1)
      CALL H2E(EMO1,EMO2,L1,0,NA,.TRUE.,DBUG,NUCFCI)
      CALL TSECND(TIM2)
      TPASS = TIM2 - TIM1
      WRITE(IW,9140) TPASS
C
      IF(NQMNUC.GT.1) THEN
      CALL TSECND(TIM1)
      CALL H2E(PMO1,PMO2,N1,1,NNA,.FALSE.,DBUG,NUCFCI)
      CALL TSECND(TIM2)
      TPASS = TIM2 - TIM1
      WRITE(IW,9150) TPASS
      ENDIF
C
C        --- CALL AO TO MO TRANSFORMATION FOR MIXED ELEC NUC ---
C
      CALL TSECND(TIM1)
      CALL H2EN(EMO1,EMO2,PMO1,PMO2,N1,L1,NNA,NA,DBUG,NUCFCI)
      CALL TSECND(TIM2)
      TPASS = TIM2 - TIM1
      WRITE(IW,9160) TPASS
C
C        --- GET MEMORY FOR OVERLAP AND ONE PARTICLE TERMS ---
C
      CALL VALFM(LOADFM)
      LMIN1     = LOADFM + 1
      LMOVE     = LMIN1  + L4
      L1E       = LMOVE  + L4    ! L1E is triangular 1 elec hamiltonian
      L1E2      = L1E    + L2    ! L1E2 is 1 elec square hamiltonian
      LH1       = L1E2   + L3
      LIPVTS    = LH1    + L3
      LWORKS    = LIPVTS + NE
      NMIN1     = LWORKS + NE
      NMOVE     = NMIN1  + N4
      N1E       = NMOVE  + N4
      N1E2      = N1E    + N2
      NH1       = N1E2   + N3
      NIPVTS    = NH1    + N3
      NWORKS    = NIPVTS + NQMNUC
      LAST      = NWORKS + NQMNUC
C
      NEEDMT = LAST - LOADFM - 1
      CALL GETFM(NEEDMT)
      CALL GOTFM(NGOT3)
C
      IF(DBUG) THEN
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) "MEMTAGOPEN                                         "
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) "  MEMORY NEEDED FOR OVERLAP AND ONE PARTICLE TERMS "
      WRITE(IW,*) "---------------------------------------------------"
      WRITE(IW,'(A25,I15,A6)') " LMIN1 (ELEC MINOR): ",L4," WORDS"
      WRITE(IW,'(A25,I15,A6)') " LMOVE               ",L4," WORDS"
      WRITE(IW,'(A25,I15,A6)') " L1E:                ",L2," WORDS"
      WRITE(IW,'(A25,I15,A6)') " L1E2:               ",L3," WORDS"
      WRITE(IW,'(A25,I15,A6)') " LH1:                ",L3," WORDS"
      WRITE(IW,'(A25,I15,A6)') " LIPVTS:             ",NE," WORDS"
      WRITE(IW,'(A25,I15,A6)') " LWORKS:             ",NE," WORDS"
      WRITE(IW,'(A25,I15,A6)') " NMIN1 (NUC MINOR):  ",N4," WORDS"
      WRITE(IW,'(A25,I15,A6)') " NMOVE               ",N4," WORDS"
      WRITE(IW,'(A25,I15,A6)') " N1E:                ",N2," WORDS"
      WRITE(IW,'(A25,I15,A6)') " N1E2:               ",N3," WORDS"
      WRITE(IW,'(A25,I15,A6)') " NH1:                ",N3," WORDS"
      WRITE(IW,'(A25,I15,A6)') " NIPVTS:             ",NQMNUC," WORDS"
      WRITE(IW,'(A25,I15,A6)') " NWORKS:             ",NQMNUC," WORDS"
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) " TOTAL:               ",(NEEDMT*8.0D+0)/1.0D+6," MB"
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) "MEMORY LEFT AFTER CALL OVERLAP AND ONE PARTICLE INT",
     *            (NGOT3*8.0D+0)/1.0D+6," MB"
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) "MEMTAGEND                                          "
                 PMEMVAL = ((MEMSET-NGOT3)*8.0D+0)/1.0D+6
      WRITE(IW,*) "MEMSTAT:  ", PMEMVAL
      WRITE(IW,*) "                                                   "

      WRITE(IW,*) "MEMTAGOPEN                                         "
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) "  MEMORY HELD TILL END OF CI ELEMENT CONSTRUCTION  "
      WRITE(IW,*) "---------------------------------------------------"
      WRITE(IW,*) "GRAND TOTAL:  ",NEEDMT*(8.0D+0/1.0D+6)," MB"
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) "MEMTAGEND                                          "
      WRITE(IW,*) "                                                   "
      CALL FLSHBF(IW)
      ENDIF
C
C        --- GET TEMPORARY MEMORY FOR OVERLAP AND DET MINOR CALLS ---
C
      CALL VALFM(LOADFM1)
      LIAA      =  LOADFM1 + 1
      NIAA      =  LIAA    + L1*2
      LSPIN     =  NIAA    + N1*2
      NSPIN     =  LSPIN   + L1*2
      LS        =  NSPIN   + N1*2
      NS        =  LS      + L2
      LS2D      =  NS      + N2
      NS2D      =  LS2D    + L3
      LAST1     =  NS2D    + N3
C
      NEEDMOV   =  LAST1 - LOADFM1 -1
      CALL GETFM(NEEDMOV)
      CALL GOTFM(NGOT4)
C
      IF(DBUG) THEN
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) "MEMTAGOPEN                                         "
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) "       TEMPORARY MEMORY NEEDED OVERLAP INTS        "
      WRITE(IW,*) "---------------------------------------------------"
      WRITE(IW,'(A25,I15,A6)') " LIAA:               ",L1*2," WORDS"
      WRITE(IW,'(A25,I15,A6)') " NIAA:               ",N1*2," WORDS"
      WRITE(IW,'(A25,I15,A6)') " LPSIN:              ",L1*2," WORDS"
      WRITE(IW,'(A25,I15,A6)') " NSPIN:              ",N1*2," WORDS"
      WRITE(IW,'(A25,I15,A6)') " LS:                 ",L2," WORDS"
      WRITE(IW,'(A25,I15,A6)') " NS:                 ",N2," WORDS"
      WRITE(IW,'(A25,I15,A6)') " LS2D:               ",L3," WORDS"
      WRITE(IW,'(A25,I15,A6)') " NS2D:               ",N3," WORDS"
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) " TOTAL:       ",(NEEDMOV*8.0D+0)/1.0D+6," MB"
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) "MEMORY LEFT AFTER CALL FOR TEMP OVERLAP ARRAYS ",
     *            (NGOT4*8.0D+0)/1.0D+6," MB"
      WRITE(IW,*) "MEMTAGEND                                          "
                 PMEMVAL = ((MEMSET-NGOT4)*8.0D+0)/1.0D+6
      WRITE(IW,*) "MEMSTAT:  ", PMEMVAL
      WRITE(IW,*) "                                                   "
      CALL FLSHBF(IW)
      ENDIF
C
C        ---  GET MOLECULAR ORBITAL OVERLAP ---
C
      CALL MOVLP(EMO1,EMO2,X(LS),X(LS2D),X(LMOVE),L1,NE,12,
     *           X(LIAA),X(LSPIN),0,DBUG)
      CALL MOVLP(PMO1,PMO2,X(NS),X(NS2D),X(NMOVE),N1,NQMNUC,441,
     *           X(NIAA),X(NSPIN),1,DBUG)

      CALL RETFM(NEEDMOV)
C
C        --- GET 1ST ORDER MINOR OF MO OVERLAP ---
C
      CALL DETMINOR(X(LMOVE),X(LMIN1),ELECDET,X(LIPVTS),X(LWORKS),NE)
      CALL DETMINOR(X(NMOVE),X(NMIN1),PRODET,X(NIPVTS),X(NWORKS),NQMNUC)
C
C        --- CALL AO TO MO TRANSFORMATION FOR 1 ELEC AND 1 PROTON ---
C
      CALL H1E(X(LH1),X(L1E),X(L1E2),EMO1,EMO2,L1,11,DBUG)
      CALL H1E(X(NH1),X(N1E),X(N1E2),PMO1,PMO2,N1,440,DBUG)
C
C        --- GET SUM OF ELECTRONIC ONE AND TWO PARTICLE ENERGIES ---
C
      CALL VALFM(LOADFM)
      LISPIN    =  LOADFM + 1
      LA        =  LISPIN + L1*2
      LAST      =  LA     + NE
      NEED2E1   =  LAST   - LOADFM - 1
      CALL GETFM(NEED2E1)
C
      CALL VALFM(LOADFM2)
      LH2       = LOADFM2 + 1
      LAST2     =  LH2    + NA*NA*NA
      NEED2E2   = LAST2   - LOADFM2 - 1
      CALL GETFM(NEED2E2)
      CALL GOTFM(NGOT1)
C
      IF(DBUG) THEN
      WRITE(IW,*)
      WRITE(IW,*) "MEMTAGOPEN                                         "
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) " MEMORY NEEDED FOR ELEC MO 2 PARTICLE INTS"
      WRITE(IW,*) "---------------------------------------------------"
      WRITE(IW,'(A25,I10,A6)') " LH2 (2 ELEC MO INTS): ",
     * NA*NA*NA," WORDS"
      WRITE(IW,'(A25,I10,A6)') " LISPIN              : ",L1*2," WORDS"
      WRITE(IW,'(A25,I10,A6)') " LA                  : ",NE," WORDS"
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) "MEMORY LEFT AFTER CALL FOR READING 2ELEC MO INTS :",
     *            (NGOT1*8.0D+0)/1.0D+6," MB"
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) "MEMTAGEND                                          "
                 PMEMVAL = ((MEMSET-NGOT1)*8.0D+0)/1.0D+6
      WRITE(IW,*) "MEMSTAT:  ", PMEMVAL
      WRITE(IW,*)
      CALL FLSHBF(IW)
      ENDIF
C
      CALL VALCON(VAL1,VAL2,X(LH1),X(LH2),X(LMIN1),X(LMOVE),X(LA),
     *             X(LISPIN),L1,NE,PRODET,NA,ELECDET,DBUG)
      CALL RETFM(NEED2E2)
C
C        --- GET SUM OF NUCLEAR ONE AND TWO PARTICLE ENERGIES ---
C
      CALL VALFM(LOADFMN2)
C
      NH2       =  LOADFMN2 + 1
      LASTN2    =  NH2      + NNA*NNA*NNA*NNA
      NEEDN2    =  LASTN2   - LOADFMN2 - 1
      CALL GETFM(NEEDN2)
      CALL GOTFM(NGOTN2)
C
      IF(DBUG) THEN
      WRITE(IW,*)
      WRITE(IW,*) "MEMTAGOPEN                                         "
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) " MEMORY NEEDED FOR NUC MO 2 PARTICLE INTS"
      WRITE(IW,*) "---------------------------------------------------"
      WRITE(IW,'(A25,I10,A6)') " NH2 (2 NUC MO INTS): ",
     * NNA*NNA*NNA*NNA," WORDS"
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) "MEMORY LEFT AFTER CALL FOR READING 2NUC MO INTS :",
     *            (NGOTN2*8.0D+0)/1.0D+6," MB"
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) "MEMTAGEND                                          "
                 PMEMVAL = ((MEMSET-NGOTN2)*8.0D+0)/1.0D+6
      WRITE(IW,*) "MEMSTAT:  ", PMEMVAL
      WRITE(IW,*)
      CALL FLSHBF(IW)
      ENDIF
C
      CALL VALCONII(VAL3,VAL4,X(NH1),X(NH2),X(NMIN1),X(NMOVE),
     *             N1,NQMNUC,ELECDET,NNA,NINTMX,PRODET,DBUG)
      CALL RETFM(NEEDN2)
C
C        --- GET SUM OF MIXED NUCLEAR ELECTRONIC TWO PARTICLE ENERGIES ---
C
      CALL VALFM(LOADFM)
      LNH2   =  LOADFM + 1
      LAST   =  LNH2    + NNA*NNA*NA*NA
      NEEDNE =  LAST - LOADFM -1
      CALL GETFM(NEEDNE)
      CALL GOTFM(NGOT2)
C
      IF(DBUG) THEN
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) "MEMTAGOPEN                                         "
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) " MEMORY NEEDED TO READ ELEC-NUC AO TO MO INTS  "
      WRITE(IW,*) "---------------------------------------------------"
      WRITE(IW,'(A25,I10,A6)') " LNH2 (2 ELEC-NUC MO INTS): ",
     *             NNA*NNA*NA*NA," WORDS"
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) "MEMORY LEFT AFTER CALL FOR ELEC-NUC MO INTS        ",
     *            (NGOT2*8.0D+0)/1.0D+6," MB"
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) "MEMTAGEND                                          "
                 PMEMVAL = ((MEMSET-NGOT2)*8.0D+0)/1.0D+6
      WRITE(IW,*) "MEMSTAT:  ", PMEMVAL
      WRITE(IW,*) "                                                   "
      CALL FLSHBF(IW)
      ENDIF
C
      CALL VALCONIII(VAL5,X(LNH2),X(LMIN1),X(LA),X(LISPIN),L1,NE,
     *               NQMNUC,X(NMIN1),NA,NNA,NINTMX,DBUG)
      CALL RETFM(NEEDNE)
      CALL RETFM(NEED2E1)
C
C        --- SUM ALL QM ENERGY TERMS ---
C
      EVAL= VAL1 + VAL2 + VAL3 + VAL4 + VAL5
      IF(NCIPRT.EQ.3) WRITE(IW,9610) EVAL
C
C        --- EVALUATE THE TOTAL OVERLAP BETWEEN STATES (TDET) ---
C
      TDET = PRODET*ELECDET
C
C        --- SUM THE TOTAL QM ENERGY WITH THE CLASSICAL NUC ENERGY ---
C
       ENUCROV = TDET*ENUCR
       VAL = EVAL + ENUCROV
C
C        --- WRITE ALL OVERLAPS TO FILE AND THE TOTAL ENERGY ---
C
      IF(NCIPRT.EQ.3) THEN
        WRITE(IW,9710) PRODET,ELECDET,TDET
        WRITE(IW,9810) ENUCROV,VAL
      ENDIF
C
      CALL RETFM(NEEDMT)
C
      RETURN
C
 9140 FORMAT(/1X,' TWO ELECTRON AO TO MO TRANSFORMATION ',
     *       'TOOK ',F12.2,' SECONDS')
 9150 FORMAT(/1X,' TWO NUCLEI AO TO MO TRANSFORMATION ',
     *       'TOOK   ',F12.2,' SECONDS')
 9160 FORMAT(/1X,' ELEC-NUC AO TO MO TRANSFORMATION ',
     *       'TOOK     ',F12.2,' SECONDS')
 9610 FORMAT(/2X,3('-'),' TOTAL QUANTUM NUC AND ELEC ENERGIES ',
     *        3('-'),/
     *       /1X,'     TOTAL OF QM TERMS =',F18.10)
 9710 FORMAT(/2X,10('-'),' DETERMINANT OVERLAPS ',11('-')/
     *       /1X,'          PROT OVERLAP =',F18.10,/
     *       /1X,'          ELEC OVERLAP =',F18.10,/
     *       /1X,'         TOTAL OVERLAP =',F18.10)
 9810 FORMAT(/2X,8('-'),' CLASSICAL NUCLEAR ENERGY  ',8('-'),/
     *       /1X,'         NUC REPULSION =',F18.10,/
     *       /2X,14('-'),' TOTAL ENERGY  ',14('-'),/
     *       /1X,'                 E TOT =',F18.10)
      END
C
C*MODULE NEONCI   *DECK MOVLP
      SUBROUTINE MOVLP(ORBI,ORBII,S,S2D,SMO,N,NPART,LFIL,IA,ISPIN,IFLAG,
     *                 DBUG)
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      LOGICAL DBUG
C
      COMMON /IOFILE/ IR,IW,IP,IS,IJKT,IDAF,NAV,IODA(950)
C
      DIMENSION SMO(NPART,*),ORBI(N,*),ORBII(N,*),S2D(N,*),S(*)
      DIMENSION IA(N*2),ISPIN(N*2)
C
      L1=N
      L2=(N*N+N)/2
C
C        --- READ IN AO OVERLAP INTS ---
C
      CALL DAREAD(IDAF,IODA,S,L2,LFIL,0)
      IF(DBUG) THEN
         WRITE(IW,900)
         CALL PRTRI(S,L1)
         CALL FLSHBF(IW)
      ENDIF
      CALL CPYTSQ(S,S2D,L1,1)
C
C        --- FORM MO OVERLAP MATRIX ---
C
      CALL VCLR(SMO,1,NPART*NPART)
C
C        --- SETUP SPIN MOLECULAR ORBITAL OVERLAP MATRIX ---
C
      IF(IFLAG.EQ.0) THEN
C
C        --- SETUP SPIN TO SPATIAL INDEX ---
C
      L = 1
      J = 1
      DO  I = 1,N*2
        IA(I) = J
        IF(L.EQ.2) THEN
          J = J +1
          L = 1
        ELSE
          L = L + 1
        ENDIF
      ENDDO
C
C        --- SET UP SPIN INDEX ARRAY ---
C
      IDEX = 1
C
      DO I = 1,N*2
         ISPIN(I) = IDEX
         IF (IDEX.EQ.1) THEN
             IDEX = 2
         ELSE
             IDEX = 1
         ENDIF
      ENDDO
C
C        --- ELECTRON MO OVERLAP MATRIX ---
C
      DO K = 1,NPART
         KK = IA(K)
         KS = ISPIN(K)
        DO L = 1,NPART
           LL = IA(L)
           LS = ISPIN(L)
          IF(KS.EQ.LS) THEN
            DO I = 1,N
              DO J = 1,N
         SMO(K,L) = SMO(K,L) + ORBI(I,KK)*ORBII(J,LL)*S2D(I,J)
              ENDDO
            ENDDO
          ENDIF
        ENDDO
      ENDDO
C
        IF(DBUG) THEN
          WRITE(IW,901)
          CALL PRSQ(SMO,NPART,NPART,NPART)
          CALL FLSHBF(IW)
        ENDIF
C
      ELSE
C
C        --- NUCLEI MO OVERLAP MATRIX ---
C        ---  ALL HIGH SPIN ---
C
      DO K = 1,NPART
         DO L = 1,NPART
           DO I = 1,N
             DO J = 1,N
              SMO(K,L) = SMO(K,L) + ORBI(I,K)*ORBII(J,L)*S2D(I,J)
             ENDDO
           ENDDO
         ENDDO
      ENDDO
        IF(DBUG) THEN
          WRITE(IW,902)
          CALL PRSQ(SMO,NPART,NPART,NPART)
          CALL FLSHBF(IW)
        ENDIF
      ENDIF
C
      RETURN
  900 FORMAT(/5X,22('-')/
     *       5X,'ATOMIC ORBITAL OVERLAP'/
     *       5X,22('-')/)
  901 FORMAT(/5X,26('-')/
     *       5X,'ELECTRONIC MO SPIN OVERLAP'/
     *       5X,26('-')/)
  902 FORMAT(/5X,22('-')/
     *       5X,'PROTON MO SPIN OVERLAP'/
     *       5X,22('-')/)
      END
C
C*MODULE NEONCI   *DECK DETMINOR
      SUBROUTINE DETMINOR(SMO,SMIN1,SDET,IPVTS,WORKS,NPART)
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      DIMENSION SMO(NPART,NPART),SMIN1(NPART,NPART)
      DIMENSION DET(2),IPVTS(NPART),WORKS(NPART)
C
C  --- CALCULATE DETERMINANT OF FULL SPIN MOLECULAR ORB OVERLAP MATRIX ---
C
      IF(NPART.EQ.1) THEN
        SDET = SMO(1,1)
        SMIN1(1,1) = 1.0D+0
      ELSE
        CALL DGEFA(SMO,NPART,NPART,IPVTS,INFOS)
        CALL DGEDI(SMO,NPART,NPART,IPVTS,DET,WORKS,11)
        SDET = DET(1) * ((10.0D+00)**DET(2))
C
C  --- CALCULATE FIRST DEGREE MINOR MATRIX ---
C
       DO I = 1,NPART
        DO J = 1, NPART
         SMIN1(I,J) = SDET*SMO(J,I)
        ENDDO
       ENDDO
      ENDIF
C
      RETURN
      END
C
C*MODULE NEONCI   *DECK H1E
      SUBROUTINE H1E(H1,HAO1,HAO2,CA,CB,N,LFIL,DBUG)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      LOGICAL DBUG
C
      COMMON /IOFILE/ IR,IW,IP,IS,IJKT,IDAF,NAV,IODA(950)
C
      DIMENSION H1(N,*),HAO1(*),HAO2(N,*),CA(N,*),CB(N,*)
C
C        --- ONE PARTICLE INTEGRAL AO TO MO TRANSFORMATION ---
C
      CALL VCLR(H1,1,N*N)
      CALL DAREAD(IDAF,IODA,HAO1,(N*N+N)/2,LFIL,0)
C
      IF(DBUG) THEN
        WRITE(IW,900)
        CALL PRTRIL(HAO1,N)
        CALL FLSHBF(IW)
      ENDIF
C
      CALL CPYTSQ(HAO1,HAO2,N,1)
C
      DO K =1,N
        DO L =1,N
          DO I=1,N
            DO J=1,N
             H1(K,L)= H1(K,L)+CA(I,K)*CB(J,L)*HAO2(I,J)
            ENDDO
          ENDDO
        ENDDO
      ENDDO
C
      IF(DBUG) THEN
        WRITE(IW,901)
        CALL PRSQ(H1,N,N,N)
        CALL FLSHBF(IW)
      ENDIF
C
      RETURN
  900 FORMAT(/5X,25('-')/
     *       5X,'ONE PARTICLE AO INTEGRALS'/
     *       5X,25('-')/)
  901 FORMAT(/5X,25('-')/
     *       5X,'ONE PARTICLE MO INTEGRALS'/
     *       5X,25('-')/)
      END
C
C*MODULE NEONCI   *DECK H2E
      SUBROUTINE H2E(CA,CB,N,IFLAG,INCOR,IWRT,DBUG,NUCFCI)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      LOGICAL IWRT,PACK2E,DBUG,NUCFCI
C
      COMMON /FMCOM/ X(1)
      COMMON /INTFIL/ NINTMX,NHEX,NTUPL,PACK2E,INTTYP,IGRDTYP
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
      COMMON /MEMINF/ MEMSET
C
      DIMENSION CA(N,*),CB(N,*)
C
      L2 = (N*N+N)/2
C
      IJKO = 8 ! ELEC-ELEC AO DISK FILE
      NFT  = 9 ! ELEC-ELEC MO DISK FILE
C
      IF(IFLAG.EQ.1) THEN  ! SET TO READ WRITE NUC-NUC
      IJKO = 68  ! NUC-NUC AO INT DISK FILE
      NFT  = 69  ! NUC-NUC MO INT DISK FILE
      CALL SEQOPN(NFT,'NUMOIN','UNKNOWN',.FALSE.,'UNFORMATTED')
      ELSE
      CALL SEQOPN(NFT,'MOINTS','UNKNOWN',.FALSE.,'UNFORMATTED')
      ENDIF
C
C     THERE ARE TWO AO TO MO TRANSFORMATIONS FOR NONORTHOGONAL
C     ORBITALS CURRENTLY AVAILABLE. THE PROGRAM WILL CHOOSE THE
C     BEST METHOD BASED ON THE AVAILABLE RESOURCES.
C
C     THE FIRST METHOD IS A THREE STEP TRANSFORMATION.
C     IT REQUIRES A LARGE AMOUNT OF MEMORY BUT IS THE
C     FASTEST TRANSFORMATION METHOD PRESENTLY AVAILABLE.
C     THE SECOND METHOD PERFORMS THE TRANSFORMATION IN TWO
C     STEPS AND IS LESS MEMORY INTENSIVE THAN THE FIRST
C     METHOD BUT IS SIGNIFICANLY SLOWER. BOTH METHODS ARE
C     PARTIALLY SEGMENTED. UNFORTUNATELY A FULLY SEGMENTED
C     TRANSFORMATION IS CURRENTLY UNAVAILABLE. THE NUMBER OF
C     TWO PARTICLE MOLECULAR ORBITAL INTEGRALS TO BE WRITTEN AND
C     READ TO DISK HAS BEEN CHANGED FOR ELECTRONS ONLY FROM THE
C     FIXED SIZE (NINTMX=15000) TO A DYNAMIC SIZE (NINTMX=INCOR**3)
C     WHEREBY IT NOW READS AND WRITES THE MO INTEGRALS IN CHUNKS
C     OF THE CORE ORBITALS CUBED. THIS DOES TWO THINGS: IT MAKES
C     THE IO MORE EFFICIENT AND IT ALLOWS THE TWO PARTICLE MO
C     INTEGRALS TO BE READ FROM DISK IN AN ORDERED FORMAT WHICH
C     MAKES ASSIGNING THE INDICIES AS THE ORBITALS ARE READ BACK
C     INTO MEMORY VERY STRAIGHTFORWARD. THE TWO PARTICLE
C     NONORTHOGONAL MO INTEGRAL ARRAY IS OF THE SIZE
C     (OCCUPIED ORBITALS)**4 WHICH IS THE AMOUNT OF DISK SPACE
C     NEEDED TO TRANSFORM THE TWO PARTICLE INTEGRALS. THIS ARRAY
C     COULD BE FURTHER REDUCED IN SIZE BY 1/2 DUE TO HERMICITY
C     BUT PRESENTLY IS NOT.
C
      NDISK = INCOR*INCOR*INCOR*INCOR
      NINTMX2 = NINTMX
      IF (IWRT) NINTMX2 = INCOR*INCOR*INCOR  ! FOR ELEC WRITE LOOP IS INCOR^3
C
C     THREE STEP TRANSFORMATION (PQ RS) --> (PQ RL) --> (PQ KL) --> (IJ KL)
C
      LOADFM=0
      CALL VALFM(LOADFM)
C
      IXPQRJ = LOADFM + 1
      IXX    = IXPQRJ + L2*N*INCOR
      IIX    = IXX    + NINTMX
      ICT    = IIX    + NINTMX
      IX1    = ICT    + N*INCOR
      IY1    = IX1    + N*N
      IXBUF  = IY1    + N*N
      IIXBUF = IXBUF  + NINTMX2
      IXXINT = IIXBUF + NINTMX2
      IXTEMP = IXXINT + L2
      LAST   = IXTEMP + N*INCOR
      NEEDINT1 = LAST - LOADFM - 1
      CALL GOTFM(NGOTINT0)
C
      IF(NUCFCI) GOTO 700
      IF(IWRT) WRITE(IW,9008)
      IF(.NOT.IWRT)  WRITE(IW,9009)
      WRITE(IW,9010) N,INCOR
      WRITE(IW,9011) NEEDINT1,NGOTINT0,NDISK
  700 CONTINUE
C
C        --- IF ENOUGH MEMORY IS AVAILABLE PROCEED WITH 3 STEP TRANSFORMATION ---
C
      IF(NEEDINT1.LT.NGOTINT0)  THEN
C
      IF(.NOT. NUCFCI) WRITE(IW,9012)
      CALL GETFM(NEEDINT1)
      CALL GOTFM(NGOTINT1)
C
      IF(DBUG) THEN
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) "MEMTAGOPEN                                         "
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) " MEMORY NEEDED FOR ARRAYS IN TWOINT1 SUBROUTINE     "
      WRITE(IW,*) "---------------------------------------------------"
      WRITE(IW,'(A15,I10,A6)') " IXPQRJ:   ",L2*N*INCOR," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IXX:      ",NINTMX," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IIX:      ",NINTMX," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IXXINT:   ",L2," WORDS"
      WRITE(IW,'(A15,I10,A6)') " ICT:      ",N*INCOR," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IX1:      ",N*N," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IY1:      ",N*N," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IXBUF:   ",NINTMX2," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IIXBUF:  ",NINTMX2," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IXTEMP:  ",INCOR*N," WORDS"
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) " MEMORY LEFT AFTER CALL FOR TWOINT1 ROUTINE ARRAYS ",
     * (NGOTINT1*8.0D+0)/1.0D+6," MB"
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) "MEMTAGEND                                          "
               PMEMVAL = ((MEMSET-NGOTINT1)*8.0D+0)/1.0D+6
      WRITE(IW,*) "MEMSTAT:  ", PMEMVAL
      WRITE(IW,*) "                                                   "
      CALL FLSHBF(IW)
      ENDIF
C
      CALL   TWOINT1(X(IXPQRJ),X(IXX),X(IIX),CA,X(ICT),
     *               N,L2,NINTMX,X(IX1),X(IY1),X(IXBUF),X(IIXBUF),CB,
     *               INCOR,X(IXXINT),X(IXTEMP),NFT,NINTMX2,DBUG)
C
      CALL RETFM(NEEDINT1)
C
      ELSE
      WRITE(IW,9015)
C
C     TWO STEP TRANSFORMATION (PQ RS) --> (PQ KL) --> (IJ KL)
C
      LOADFM = 0
      CALL VALFM(LOADFM)
C
      IXPQKL = LOADFM + 1
      IXX    = IXPQKL + L2*INCOR*INCOR
      IIX    = IXX    + NINTMX
      ICT    = IIX    + NINTMX
      IX1    = ICT    + N*INCOR
      IY1    = IX1    + N*N
      IXBUF  = IY1    + N*N
      IIXBUF = IXBUF  + NINTMX2
      IXXINT = IIXBUF + NINTMX2
      LAST   = IXXINT + INCOR*INCOR
      NEEDINT2 = LAST - LOADFM - 1
      IF(.NOT. NUCFCI) WRITE(IW,9016) NEEDINT2,NGOTINT0
C
        IF(NEEDINT2.LT.NGOTINT0) THEN
        IF(.NOT. NUCFCI) WRITE(IW,9018)
C
C
      CALL GETFM(NEEDINT2)
      CALL GOTFM(NGOTINT2)
C
      IF(DBUG) THEN
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) "MEMTAGOPEN                                         "
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) " MEMORY NEEDED FOR ARRAYS IN TWOINT2 SUBROUTINE     "
      WRITE(IW,*) "---------------------------------------------------"
      WRITE(IW,'(A15,I10,A6)') " IXPQKL:   ",L2*INCOR*INCOR," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IXX:      ",NINTMX," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IIX:      ",NINTMX," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IXXINT:   ",INCOR*INCOR," WORDS"
      WRITE(IW,'(A15,I10,A6)') " ICT:      ",N*INCOR," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IX1:      ",N*N," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IY1:      ",N*N," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IXBUF:    ",NINTMX2," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IIXBUF:   ",NINTMX2," WORDS"
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) " MEMORY LEFT AFTER CALL FOR TWOINT2 ROUTINE ARRAYS ",
     * (NGOTINT2*8.0D+0)/1.0D+6," MB"
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) "MEMTAGEND                                          "
               PMEMVAL = ((MEMSET-NGOTINT2)*8.0D+0)/1.0D+6
      WRITE(IW,*) "MEMSTAT:  ", PMEMVAL
      WRITE(IW,*) "                                                   "
      CALL FLSHBF(IW)
      ENDIF
C
      CALL   TWOINT2(X(IXPQKL),X(IXX),X(IIX),CA,X(ICT),
     *               N,L2,NINTMX,X(IX1),X(IY1),X(IXBUF),X(IIXBUF),
     *               CB,INCOR,X(IXXINT),NFT,NINTMX2,DBUG)
C
      CALL RETFM(NEEDINT2)
C
        ELSE
          WRITE(IW,9019) NEEDINT2-NGOTINT0
          CALL ABRT
        ENDIF
      ENDIF
C
C        --- RESET IJKO TO READ ELEC-ELEC AO INTS ---
      IJKO=8
C
      RETURN
C
 9008 FORMAT(//5X,36('-')/
     *        5X,'TWO ELECTRON INTEGRAL TRANSFORMATION'/
     *        5X,36('-')/)
 9009 FORMAT(//5X,35('-')/
     *        5X,'TWO NUCLEAR INTEGRAL TRANSFORMATION'/
     *        5X,35('-')/)
 9010 FORMAT(1X,' TOTAL NUMBER OF ATOMIC ORBITALS             = ',I5/
     *       1X,' TOTAL NUMBER OF OCCUPIED MOLECULAR ORBITALS = ',I5/)
 9011 FORMAT(1X,' MEMORY NEEDED FOR THREE STEP TRANSFORMATION = ',
     *        I12,' WORDS'/
     *       1X,' MEMORY AVAILABLE FOR TRANSFORMATION         = ',
     *        I12,' WORDS'/
     *       1X,' DISK SPACE NEEDED FOR STORING INTEGRALS     = ',
     *        I12,' WORDS'/)
 9012 FORMAT(1X,' CHOOSING THREE STEP TRANSFORMATION...'/)
 9015 FORMAT(1X,' *** NOT ENOUGH MEMORY FOR THREE STEP TRANSFORMATION ',
     *       '! ***'/
     *       1X,' *** WILL TRY THE SLOWER TWO STEP TRANSFORMATION ***'/)
 9016 FORMAT(1X,' MEMORY NEEDED FOR TWO STEP TRANSFORMATION   = ',
     *        I12,' WORDS'/
     *        1X,' MEMORY AVAILABLE FOR TRANSFORMATION         = ',
     *        I12,' WORDS')
 9018 FORMAT(/1X,' CHOOSING THE TWO STEP TRANSFORMATION...'/)
 9019 FORMAT(/1X,' *** NOT ENOUGH MEMORY FOR TWO STEP TRANSFORMATION !',
     *      ' ***'/
     *        1X,' *** INCREASE MEMORY BY AT LEAST',I12,' WORDS ***'/)
      END
C
C*MODULE NEONCI   *DECK TWOINT2
      SUBROUTINE TWOINT2(XPQKL,XX,IX,CA,CT,NUM,L2,NINTMX,X1,Y1,
     *                   XBUF,IXBUF,CB,INCOR,XXINT,NFT,NINTMX2,DBUG)
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      LOGICAL DBUG
C
      DIMENSION XPQKL(L2,INCOR*INCOR),XX(NINTMX),
     *          IX(NINTMX), CA(NUM,NUM),CT(INCOR,NUM),X1(NUM,NUM),
     *          Y1(NUM,NUM), XBUF(NINTMX2),IXBUF(NINTMX2),CB(NUM,NUM),
     *          XXINT(INCOR*INCOR)
C
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
C
      CALL VCLR(XPQKL,1,L2*INCOR*INCOR)
      CALL SEQREW(IJKO)
   77 CALL PREAD(IJKO,XX,IX,NXX,NINTMX)
C
C        --- 1ST AND 2ND INDEX TO BE TRANSFORMED ---
      CALL TWOTRFNO(XPQKL,CA,CB,L2,NUM,INCOR,XX,IX,NXX,0)
      IF(NXX.GT.0) GO TO 77
C
      IF(DBUG) THEN
        WRITE(IW,9010)
        WRITE(IW,900)
        CALL FLSHBF(IW)
      ENDIF
C
C        --- 3RD AND 4TH INDEX TO BE TRANSFORMED ---
C
      NUMX=0
      NREC=0
      MIJ=0
      DO 10 MJ=1,INCOR
      DO 10 MP=1,NUM
C
   10 CT(MJ,MP)=CB(MP,MJ)
      MI=1
   60 CONTINUE
      DO MJ=1,INCOR
C
         MAJ=0
         MIJ=MIJ+1
         MPQ=0
         DO 20 MP=1,NUM
         DO 20 MQ=1,MP
            MPQ=MPQ+1
            X1(MP,MQ)=XPQKL(MPQ,MIJ)
   20       X1(MQ,MP)=X1(MP,MQ)
C
         CALL MRARBR(X1,NUM,NUM,NUM,CA,NUM,INCOR,Y1,NUM)
         CALL MRARBR(CT,INCOR,INCOR,NUM,Y1,NUM,INCOR,X1,NUM)
C
         DO 30 MK=1,INCOR ! NOW LOOPS OVER ALL K values
           DO 30 ML=1,INCOR
C
            MAJ=MAJ+1
            DUM=X1(ML,MK)
            XXINT(MAJ)=DUM
   30      CONTINUE
C
      MJX = MJ
      CALL TRFWTNO(NFT,XXINT,NUMX,XBUF,NREC,
     *             DBUG,MI,MJX,IXBUF,INCOR,NINTMX2)
C
      ENDDO
C
      MI=MI+1
      IF(MI.LE.INCOR) GO TO 60 ! ONLY LOOPING OVER OCCUPIED ORBS
       NUMX = -NUMX
       NREC = NREC +1
       CALL PWRIT(NFT,XBUF,IXBUF,NUMX,NINTMX2) ! WRITE LAST BIT OF INTEGRALS
C
      IF(DBUG) WRITE(IW,9011)
      CALL FLSHBF(IW)
C
      RETURN
C
  900 FORMAT(/5X,'PRINTING OUT THE MO TRANSFORMED 2 PARTICLE TERMS')
 9010 FORMAT(1X,' FIRST AND SECOND INDICIES ARE TRANSFORMED')
 9011 FORMAT(1X,' THIRD AND FOURTH INDICIES ARE TRANSFORMED')
      END
C
C*MODULE NEONCI   *DECK TWOINT1
      SUBROUTINE TWOINT1(XPQRJ,XX,IX,CA,CT,NUM,L2,NINTMX,X1,Y1,
     *                  XBUF,IXBUF,CB,INCOR,XXINT,XTEMP,
     *                  NFT,NINTMX2,DBUG)
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      LOGICAL DBUG
C
      DIMENSION XPQRJ(INCOR,L2,NUM),XX(NINTMX),
     *          IX(NINTMX), CA(NUM,NUM),CT(INCOR,NUM),X1(NUM,NUM),
     *          Y1(NUM,NUM), XBUF(NINTMX2),IXBUF(NINTMX2),CB(NUM,NUM),
     *          XXINT(L2),XTEMP(NUM,INCOR)
C
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
C
      NXPQRJ = INCOR*L2*NUM
      CALL VCLR(XPQRJ,1,NXPQRJ)
C
      CALL SEQREW(IJKO)
   77 CALL PREAD(IJKO,XX,IX,NXX,NINTMX)
C
C        --- READ IN AO INTEGRALS AND TRANSFORM FIRST INDEX ---
C
      CALL ONETRF(XPQRJ,CA,L2,NUM,INCOR,XX,IX,NXX,0)
      IF(NXX.GT.0) GO TO 77
C
C        --- TRANSFORM SECOND INDEX ---
C
      CALL TR2MCNO(XPQRJ,CB,XTEMP,NUM,L2,INCOR,INCOR)
C
C
      IF(DBUG) THEN
       WRITE(IW,9010)
       WRITE(IW,900)
       CALL FLSHBF(IW)
      ENDIF
C
C        --- 3RD AND 4TH INDEX TO BE TRANSFORMED ---
C
      NUMX=0
      NREC=0
      DO 10 MJ=1,INCOR
      DO 10 MP=1,NUM
C
   10 CT(MJ,MP)=CB(MP,MJ)
      MI=1
   60 CONTINUE
      DO MJ=1,INCOR
C
         MAJ=0
         MPQ=0
         DO 20 MP=1,NUM
         DO 20 MQ=1,MP
            MPQ=MPQ+1
            X1(MP,MQ)=XPQRJ(MI,MPQ,MJ)
   20       X1(MQ,MP)=X1(MP,MQ)
C
         CALL MRARBR(X1,NUM,NUM,NUM,CA,NUM,INCOR,Y1,NUM)
         CALL MRARBR(CT,INCOR,INCOR,NUM,Y1,NUM,INCOR,X1,NUM)
C
         DO 30 MK=1,INCOR ! NOW LOOPS OVER ALL K VALUES
           DO 30 ML=1,INCOR
C
            MAJ=MAJ+1
            DUM=X1(ML,MK)
            XXINT(MAJ)=DUM
   30      CONTINUE
      MJX = MJ
      CALL TRFWTNO(NFT,XXINT,NUMX,XBUF,NREC,
     *             DBUG,MI,MJX,IXBUF,INCOR,NINTMX2)
C
      ENDDO
C
      MI=MI+1
      IF(MI.LE.INCOR) GO TO 60 ! LOOP ONLY OVER OCCUPIED ORBS
       NUMX = -NUMX
       NREC = NREC +1
       CALL PWRIT(NFT,XBUF,IXBUF,NUMX,NINTMX2) ! write in chunks of nintmx
C
      IF(DBUG) WRITE(IW,9011)
C
      RETURN
C
  900 FORMAT(/5X,'PRINTING OUT THE MO TRANSFORMED 2 PARTICLE TERMS')
 9010 FORMAT(1X,' FIRST AND SECOND INDICIES ARE TRANSFORMED')
 9011 FORMAT(1X,' THIRD AND FOURTH INDICIES ARE TRANSFORMED')
      END
C
C*MODULE NEONCI   *DECK H2EN
      SUBROUTINE H2EN(CA,CB,CNA,CNB,N1,L1,INCOR_NUC,INCOR,DBUG,NUCFCI)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      LOGICAL PACK2E,DBUG,NUCFCI
C
      COMMON /FMCOM/ X(1)
      COMMON /INTFIL/ NINTMX,NHEX,NTUPL,PACK2E,INTTYP,IGRDTYP
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
      COMMON /MEMINF/ MEMSET
C
      DIMENSION CA(L1,*),CB(L1,*),CNA(N1,*),CNB(N1,*)
C
      IF(.NOT. NUCFCI) WRITE(IW,9108)
C
      N2         = (N1*N1+N1)/2
      INCOR3     = INCOR*INCOR
      INCOR_NUC3 = INCOR_NUC*INCOR_NUC
C
      IJKO = 67
      NFT  = 71
      CALL SEQOPN(NFT,'NUELMO','UNKNOWN',.FALSE.,'UNFORMATTED')
C
      CALL GOTFM(NGOTMIX0)
      LOADFM=0
      CALL VALFM(LOADFM)
C
      IXPQRJ = LOADFM + 1
      IXX    = IXPQRJ + N2*L1*INCOR
      IIX    = IXX    + NINTMX
      ICT    = IIX    + NINTMX
      IX1    = ICT    + N1*INCOR_NUC
      IY1    = IX1    + N1*N1
      IXBUF  = IY1    + N1*N1
      IXXINT = IXBUF  + NINTMX
      IXTEMP = IXXINT + INCOR*INCOR
      IIXBUF = IXTEMP + INCOR*L1
      LAST   = IIXBUF + NINTMX
      NEEDMIX = LAST - LOADFM - 1
C
      IF(.NOT. NUCFCI) WRITE(IW,9109) NEEDMIX,NGOTMIX0,INCOR3*INCOR_NUC3
C
      IF(NEEDMIX.LT.NGOTMIX0)  THEN
      CALL GETFM(NEEDMIX)
C
      IF(DBUG) THEN
      CALL GOTFM(NGOTMIX1)
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) "MEMTAGOPEN                                         "
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) " MEMORY NEEDED FOR ARRAYS IN TWOMIX1 SUBROUTINE     "
      WRITE(IW,*) "---------------------------------------------------"
      WRITE(IW,'(A15,I10,A6)') " IXPQRJ:  ",N2*L1*INCOR," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IXX:     ",NINTMX," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IIX:     ",NINTMX," WORDS"
      WRITE(IW,'(A15,I10,A6)') " ICT:     ",N1*INCOR_NUC," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IX1:     ",N1*N1," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IY1:     ",N1*N1," WORDS"
      WRITE(IW,'(A15,I10,A6)') " XBUF:   ",NINTMX," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IXBUF:   ",NINTMX," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IXXINT:  ",INCOR*INCOR," WORDS"
      WRITE(IW,'(A15,I10,A6)') " IXTEMP:  ",L1*INCOR," WORDS"
      WRITE(IW,*) "                                                   "
      WRITE(IW,*) "MEMORY LEFT AFTER CALL FOR TWOMIX1 ROUTINE ARRAYS  ",
     *            (NGOTMIX1*8.0D+0)/1.0D+6," MB"
      WRITE(IW,*) "==================================================="
      WRITE(IW,*) "MEMTAGEND                                          "
               PMEMVAL = ((MEMSET-NGOTMIX1)*8.0D+0)/1.0D+6
      WRITE(IW,*) "MEMSTAT:  ", PMEMVAL
      WRITE(IW,*) "                                                   "
      CALL FLSHBF(IW)
      ENDIF
C
      CALL   HEN_TRAN(X(IXPQRJ),X(IXX),X(IIX),CA,X(ICT),
     *               L1,NINTMX,X(IX1),X(IY1),X(IXBUF),CB,
     *               N1,N2,CNA,CNB,INCOR,INCOR_NUC,
     *               X(IXXINT),X(IXTEMP),NFT,X(IIXBUF),DBUG)

C
         CALL RETFM(NEEDMIX)
      ELSE
         WRITE(IW,9010) NEEDMIX-NGOTMIX0
         CALL ABRT
      ENDIF
C
C        --- RESET IJKO TO 2 ELEC AO INTS ---
C
      IJKO = 8
C
      RETURN
C
 9108 FORMAT(//,5X,40('-'),
     *       /,5X,'ELECTRON-NUCLEAR INTEGRAL TRANSFORMATION',
     *       /,5X,40('-'))
 9109 FORMAT(/1X,' MEMORY NEEDED FOR ELEC-NUC TRANSFORMATION = ',
     *        I12,' WORDS'/
     *       1X,' MEMORY AVAILABLE FOR TRANSFORMATION       = ',
     *        I12,' WORDS'/
     *       1X,' DISK SPACE NEEDED FOR STORING INTEGRALS   = ',
     *        I12,' WORDS'/)
 9010 FORMAT(/1X,' *** NOT ENOUGH MEMORY FOR ELEC-NUC TRANSFORMATION !',
     *      ' ***'/
     *        1X,' *** INCREASE MEMORY BY AT LEAST',I12,' WORDS ***'/)
      END
C
C*MODULE NEONCI   *DECK HEN_TRAN
      SUBROUTINE HEN_TRAN(XPQRJ,XX,IX,CA,CT,NUM,NINTMX,X1,Y1,
     *                    XBUF,CB,N1,N2,CNA,CNB,INCOR,INCOR_NUC,
     *                    XXINT,XTEMP,NFT,IXBUF,DBUG)
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      LOGICAL DBUG
C
      DIMENSION XPQRJ(INCOR,N2,NUM),
     *          XX(NINTMX),IX(NINTMX),
     *          CA(NUM,NUM),CT(INCOR_NUC,N1),X1(N1,N1),Y1(N1,N1),
     *          XBUF(*),CB(NUM,NUM),IXBUF(NINTMX),
     *          CNA(N1,N1),CNB(N1,N1),
     *          XXINT(INCOR*INCOR),XTEMP(NUM,*)
C
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
C
      CALL VCLR(XPQRJ,1,INCOR*N2*NUM)
      CALL SEQREW(IJKO)
   77 CALL PREAD(IJKO,XX,IX,NXX,NINTMX)
C
C        --- READ IN ELEC-NUC AO INTS AND TRANSFORM FIRST ELEC INDEX ---
C
      CALL NE_ONETRF(XPQRJ,CA,N2,NUM,INCOR,XX,IX,NXX)
      IF(NXX.GT.0) GO TO 77
C
C        --- TRANSFORM SECOND ELEC INDEX ---
C
      CALL TR2MCNO(XPQRJ,CB,XTEMP,NUM,N2,INCOR,INCOR)
C
      IF(DBUG) THEN
         WRITE(IW,900)
         WRITE(IW,9110)
      ENDIF
C
C        --- 3RD AND 4TH INDEX TO BE TRANSFORMED ---
C
      NUMX=0
      NREC=0
      DO 10 MJ=1,INCOR_NUC
      DO 10 MP=1,N1
   10 CT(MJ,MP)=CNB(MP,MJ)
C
      MI=1
   60 CONTINUE
      DO MJ=1,INCOR
C
         MAJ=0
         MPQ=0
         DO 20 MP=1,N1
         DO 20 MQ=1,MP
            MPQ=MPQ+1
            X1(MP,MQ)=XPQRJ(MI,MPQ,MJ)
   20       X1(MQ,MP)=X1(MP,MQ)
C
         CALL MRARBR(X1,N1,N1,N1,CNA,N1,INCOR_NUC,Y1,N1)
         CALL MRARBR(CT,INCOR_NUC,INCOR_NUC,N1,Y1,N1,INCOR_NUC,X1,N1)
C
C        --- MK AND ML RUN OVER THE PROTON INDICIES FOR WRITING INTS TO DISK ---
C        --- MI AND MJ RUN OVER THE ELECTRON INDICIES FOR WRITING INTS TO DISK---
C
         DO 30 MK=1,INCOR_NUC
           DO 30 ML=1,INCOR_NUC
            MAJ=MAJ+1
            DUM=X1(ML,MK)
            XXINT(MAJ)=-1.0D+0*DUM
   30      CONTINUE
C
      MJX=MJ
      CALL TRFWTNO(NFT,XXINT,NUMX,XBUF,NREC,
     *             DBUG,MI,MJX,IXBUF,INCOR_NUC,NINTMX)
      ENDDO
C
      MI=MI+1
      IF(MI.LE.INCOR) GO TO 60 ! LOOP OVER OCCUPIED ELEC ORBS ONLY
      NUMX = -NUMX
      NREC = NREC +1
      CALL PWRIT(NFT,XBUF,IXBUF,NUMX,NINTMX) ! WRITE LAST BIT OF INTS TO DISK
C
      IF(DBUG) WRITE(IW,9111)
      CALL FLSHBF(IW)
C
      RETURN
C
  900 FORMAT(/5X,'PRINTING OUT THE MO TRANSFORMED MIXED PARTICLE TERMS')
 9110 FORMAT(1X,' FIRST AND SECOND ELEC INDICIES ARE TRANSFORMED')
 9111 FORMAT(1X,' THIRD AND FOURTH PROT INDICIES TRANSFORMED')
      END
C
C*MODULE NEONCI   *DECK VALCON
      SUBROUTINE VALCON(VAL1,VAL2,H1,H2,DMINR1,SSMO,IA,ISPIN,
     *                   L1,NPART,PRODET,INCOR,ELECDET,DBUG)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      LOGICAL DBUG
C
      PARAMETER(MXATM=2000)
C
      COMMON /FMCOM/ X(1)
      COMMON /INFOA / NAT,ICH,MUL,NUM,NQMT,NE,NA,NB,
     *                ZAN(MXATM),C(3,MXATM),IAN(MXATM)
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
      COMMON /NCIFLG/ NUMCNT,NCITRN,NCIPRT,NUMNCI,NUMEXE,NUMEXP,
     *                IEXE(40),IEXP(40),CUTNCI
C
      DIMENSION H1(L1,L1), H2(INCOR,INCOR,INCOR),
     *           DMINR1(NPART,NPART)
      DIMENSION SSMO(NPART,NPART),IA(NPART)
      DIMENSION ISPIN(L1*2)
C
      INCOR2 = INCOR*INCOR
      INCOR3 = INCOR2*INCOR
      NFT = 9
      VAL1=0.0D+00
      VAL2=0.0D+00
C
C        --- SETUP INDEX ARRAY ---
C
      KK= 1
      JJ= 1
      LL= 1
      DO I = 1,NE
         IA(I) = JJ
         IF(LL.EQ.2) THEN
           JJ=JJ+1
           LL=1
         ELSE
           LL=LL+1
         ENDIF
      ENDDO
C
C        --- SET UP SPIN INDEX ARRAY ---
C
      IDEX = 1
      DO I = 1,L1*2
         ISPIN(I) = IDEX
         IF (IDEX.EQ.1) THEN
            IDEX = 2
         ELSE
            IDEX = 1
         ENDIF
      ENDDO
C
C        --- ONE ELEC TERM ---
      IF(DBUG) WRITE(IW,900)
      DO I = 1,NE
        II=IA(I)
        DO J = 1,NE
           JJ=IA(J)
           VAL1= VAL1 + H1(II,JJ)*DMINR1(I,J)*PRODET
         IF(DBUG)  WRITE(IW,901) H1(II,JJ),II,JJ,DMINR1(I,J),I,J
        ENDDO
      ENDDO
C
C        --- TWO ELEC TERM ---
C
      CALL VALFM(LOADFM)
      LXX     = LOADFM + 1
      LIX     = LXX    + INCOR3
      LIA     = LIX    + INCOR3
      LAST    = LIA    + INCOR
      NEED2   = LAST - LOADFM - 1
      CALL GETFM(NEED2)
C
C  READ 2 ELECTRON INTEGRALS FROM DAF FILE
C
      IF(DBUG) WRITE(IW,902)
      I1 = 1
      I2 = 2
      CALL SEQREW(NFT)
      DO I = 1, INCOR
       II=IA(I)
       CALL RDNOCI2E(NFT,H2,INCOR,INCOR,INCOR3,X(LIA),
     *                X(LXX),X(LIX),INCOR3)
C
       IS1 = ISPIN(I1)
       DO J = 1,NE
         JS1 = ISPIN(J)
         JJ=IA(J)
        DO K = 1,NE
        KS1=ISPIN(K)
        KK=IA(K)
         DO L = 1,NE
          LS1=ISPIN(L)
          LL=IA(L)
         IF (IS1.NE.JS1.OR.KS1.NE.LS1) GOTO 400
C
         DMIN2= (SSMO(J,I1)*SSMO(L,K)-
     *                 (SSMO(L,I1)*SSMO(J,K)))*ELECDET
         VAL2 = VAL2 + H2(JJ,KK,LL)*DMIN2*PRODET
         IF(DBUG)  WRITE(IW,903) H2(JJ,KK,LL),II,JJ,KK,LL,
     *             DMIN2,I1,J,K,L
  400 CONTINUE
         IF(I2.GT.NE) GOTO 500
         IS2 = ISPIN(I2)
         IF (IS2.NE.JS1.OR.KS1.NE.LS1) GOTO 500
C
      DMIN2= (SSMO(J,I2)*SSMO(L,K)-
     *                 (SSMO(L,I2)*SSMO(J,K)))*ELECDET
         VAL2 = VAL2 + H2(JJ,KK,LL)*DMIN2*PRODET
         IF(DBUG)  WRITE(IW,903) H2(JJ,KK,LL),II,JJ,KK,LL,
     *             DMIN2,I2,J,K,L
  500 CONTINUE
         ENDDO
        ENDDO
       ENDDO
       I1 = I1 + 2
       I2 = I2 + 2
      ENDDO
       CALL SEQREW(NFT)
       CALL RETFM(NEED2)
C
      VAL2 = VAL2/2.0D+00
C
      IF (NCIPRT.EQ.3 .OR. DBUG) WRITE(IW,9010)  VAL1,VAL2
C
      RETURN
C
  900 FORMAT(1X,' ONE ELECTRON TERM CONSTRUCTION ')
  901 FORMAT(1X,' H1(I,J)',F14.8,I5,I5,
     *       '  MINOR1(ISPIN,JSPIN)',F14.8,I5,I5)
  902 FORMAT(/1X,' TWO ELECTRON TERM CONSTRUCTION ')
  903 FORMAT(1X,' H2(I,J,K,L)',F14.8,I5,I5,I5,I5,
     *       '  MINOR2(ISPIN,JSPIN,KSPIN,LSPIN)',F14.8,I5,I5,I5,I5)
 9010 FORMAT(//5X,32('-')/
     *        5X,' ENERGY COMPONENTS AND OVERLAPS ',/
     *        5X,32('-')/
     *       /2X,13('-'),' ELECTRON TERMS ',13('-'),/
     *       /1X,'       ONE ELEC ENERGY =',F18.10,
     *       /1X,'       TWO ELEC ENERGY =',F18.10)
      END
C
C*MODULE NEONCI   *DECK VALCONII
      SUBROUTINE VALCONII(VAL1,VAL2,H1,H2,DMINR1,SSMO,
     *                    L1,NPART,ELECDET,INCOR,NINTMX,PRODET,DBUG)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      LOGICAL DBUG
C
      COMMON /FMCOM/ X(1)
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
      COMMON /NCIFLG/ NUMCNT,NCITRN,NCIPRT,NUMNCI,NUMEXE,NUMEXP,
     *                IEXE(40),IEXP(40),CUTNCI
      COMMON /NUCMOI/ NUNIQN,IUNIQN(20),IUNIQT(20),NQMNUC,IQMNUC(20),
     *                IQNTYP(20),NUMNB,NUCST,NAUXNB,IAUXNB(20),NUMULT,
     *                NNA,NNB,NTAUXB
C
      DIMENSION H1(L1,L1), H2(INCOR,INCOR,INCOR,INCOR),
     *           DMINR1(NPART,NPART)
      DIMENSION SSMO(NPART,NPART)
C
      INCOR2 = INCOR*INCOR
      INCOR4 = INCOR2*INCOR2
      VAL1=0.0D+00
      VAL2=0.0D+00
      NFT = 69
C
C        --- ONE NUC TERM ---
C
      IF(DBUG) WRITE(IW,900)
      DO I = 1,NQMNUC
        DO J = 1,NQMNUC
        VAL1= VAL1 + H1(I,J)*DMINR1(I,J)*ELECDET
         IF(DBUG)  WRITE(IW,901) H1(I,J),I,J,DMINR1(I,J),I,J
        ENDDO
      ENDDO
C
      IF(NQMNUC.GT.1) THEN
C
C        --- TWO NUC TERM ---
C
      IF(DBUG) WRITE(IW,902)
C
      CALL VALFM(LOADFM)
      LXX     = LOADFM + 1
      LIX     = LXX    + NINTMX
      LIA     = LIX    + NINTMX
      LAST    = LIA    + INCOR
      NEED2N   = LAST - LOADFM - 1
      CALL GETFM(NEED2N)
C
      CALL RDNOCI2N_NE(NFT,H2,INCOR,INCOR,INCOR4,X(LIA),X(LXX),
     *             X(LIX),NINTMX)
C
      DO I = 1,NQMNUC
        DO J = 1,NQMNUC
         DO K = 1,NQMNUC
          DO L = 1,NQMNUC
            DMIN2= (SSMO(J,I)*SSMO(L,K)-(SSMO(L,I)*SSMO(J,K)))*PRODET
            VAL2 = VAL2 + H2(I,J,K,L)*DMIN2*ELECDET
         IF(DBUG)  WRITE(IW,903) H2(I,J,K,L),I,J,K,L,
     *             DMIN2,I,J,K,L
         ENDDO
        ENDDO
       ENDDO
      ENDDO
      VAL2 = VAL2/2.0D+00
      CALL RETFM(NEED2N)
      ENDIF
C
      IF (NCIPRT.EQ.3 .OR. DBUG) WRITE(IW,9020) VAL1,VAL2
C
      RETURN
C
  900 FORMAT(1X,' ONE QM NUCLEUS TERM CONSTRUCTION ')
  901 FORMAT(1X,' H1(I,J)',F14.8,I5,I5,
     *       '  MINOR1(I,J)',F14.8,I5,I5)
  902 FORMAT(/1X,' TWO QM NUCLEI TERM CONSTRUCTION ')
  903 FORMAT(1X,' H2(I,J,K,L)',F14.8,I5,I5,I5,I5,
     *       '  MINOR2(I,J,K,L)',F14.8,I5,I5,I5,I5)
 9020 FORMAT(/2X,14('-'),' PROTON TERMS ',15('-'),/
     *       /1X,'        ONE NUC ENERGY =',F18.10,
     *       /1X,'        TWO NUC ENERGY =',F18.10)
      END
C
C*MODULE NEONCI   *DECK VALCONIII
      SUBROUTINE VALCONIII(VAL2,H2,DMINR1,IA,ISPIN,L1,NPART,NUCP,DNM,
     *                      INCOR,INCOR_NUC,NINTMX,DBUG)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      PARAMETER(MXATM=2000)
C
      LOGICAL DBUG
C
      COMMON /FMCOM/ X(1)
      COMMON /INFOA / NAT,ICH,MUL,NUM,NQMT,NE,NA,NB,
     *                ZAN(MXATM),C(3,MXATM),IAN(MXATM)
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
      COMMON /NCIFLG/ NUMCNT,NCITRN,NCIPRT,NUMNCI,NUMEXE,NUMEXP,
     *                IEXE(40),IEXP(40),CUTNCI
      COMMON /NUCMOI/ NUNIQN,IUNIQN(20),IUNIQT(20),NQMNUC,IQMNUC(20),
     *                IQNTYP(20),NUMNB,NUCST,NAUXNB,IAUXNB(20),NUMULT,
     *                NNA,NNB,NTAUXB
C
      DIMENSION H2(INCOR_NUC,INCOR_NUC,INCOR,INCOR), DMINR1(NPART,NPART)
      DIMENSION DNM(NUCP,NUCP)
      DIMENSION ISPIN(L1*2),IA(NPART)
C
      INCOR2 = INCOR*INCOR_NUC
      INCOR4 = INCOR2*INCOR2
      NFT=71
      VAL2=0.0D+00
C
      CALL VCLR(H2,1,INCOR4)
C
C        --- SET UP INDEX ARRAY ---
C
      JJ= 1
      LL= 1
      DO I = 1,NE
        IA(I) = JJ
        IF(LL.EQ.2) THEN
          JJ=JJ+1
          LL=1
        ELSE
          LL=LL+1
        ENDIF
      ENDDO
C
C        --- SET UP SPIN INDEX ARRAY ---
C
      IDEX = 1
      DO I = 1,L1*2
         ISPIN(I) = IDEX
         IF (IDEX.EQ.1) THEN
             IDEX = 2
         ELSE
             IDEX = 1
         ENDIF
      ENDDO
C
C        --- MIXED ELEC NUC TERM ---
C
      CALL VALFM(LOADFM)
      LXX     = LOADFM + 1
      LIX     = LXX    + NINTMX
      LIA     = LIX    + NINTMX
      LAST    = LIA    + INCOR_NUC
      NEED2   = LAST - LOADFM - 1
      CALL GETFM(NEED2)
C
C  READ ELECTRON NUCLEAR MIXED INTEGRALS FROM DAF FILE
C
      CALL RDNOCI2N_NE(NFT,H2,INCOR_NUC,INCOR,INCOR4,X(LIA),X(LXX),
     *             X(LIX),NINTMX)
      CALL RETFM(NEED2)
C
      IF(DBUG) WRITE(IW,902)
C
      DO I = 1,NE
         IS1 = ISPIN(I)
         II=IA(I)
       DO J = 1,NE
         JS1 = ISPIN(J)
         JJ=IA(J)
        DO K = 1,NQMNUC
         DO L = 1,NQMNUC
         IF (IS1.NE.JS1) THEN
         ELSE
         VAL2 = VAL2 + H2(K,L,II,JJ)*DMINR1(I,J)*DNM(K,L)
         IF(DBUG) WRITE(IW,903) H2(K,L,II,JJ),K,L,II,JJ,DMINR1(I,J),
     *            I,J,DNM(K,L),K,L
         ENDIF
         ENDDO
        ENDDO
       ENDDO
      ENDDO
C
      IF (NCIPRT.EQ.3 .OR. DBUG) WRITE(IW,9030)  VAL2
C
      RETURN
C
  902 FORMAT(/1X,' MIXED ELECTRON-NUCLEAR TERM CONSTRUCTION ')
  903 FORMAT(1X,' H2(I,J,K,L)',F14.8,I5,I5,I5,I5,
     *       '  NUC MINOR1(I,J)',F14.8,I5,I5,
     *       '  ELEC MINOR1(KSPIN,LSPIN)',F14.8,I5,I5)
 9030 FORMAT(/2X,14('-'),' ELEC-NUC TERM ',14('-'),/
     *       /1X,'       ELEC-NUC ENERGY =',F18.10)
      END
C
C*MODULE NEONCI   *DECK TWOTRFNO
      SUBROUTINE TWOTRFNO(XPQKL,CA,CB,NUM2,NUM,NOC,XX,IX,NXX,NPACKIC)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      PARAMETER (MXAO=8192)
C
      DIMENSION XPQKL(NUM2,*),CA(NUM,*),CB(NUM,*),XX(*),IX(*)
C
      COMMON /IJPAIR/ IA(MXAO)
      COMMON /PCKLAB/ LABSIZ
C
      NINT = ABS(NXX)
C
      IF(LABSIZ.GE.2) THEN
C
C     ----- TWO BYTES PER LABEL -----
C
         DO 20 M=1,NINT
                       NPACK = NPACKIC+M
            VAL1 = XX(NPACK)
            VAL3 = VAL1
C
*I32                     LABEL = IX( 2*NPACK - 1 )
*I32                     IPACK = ISHFT( LABEL, -16 )
*I32                     JPACK = IAND( LABEL, 65535 )
*I32                     LABEL = IX( 2*NPACK     )
*I32                     KPACK = ISHFT( LABEL, -16 )
*I32                     LPACK = IAND( LABEL, 65535 )
*I64                     LABEL = IX(NPACK)
*I64                     IPACK = ISHFT( LABEL, -48 )
*I64                     JPACK = IAND( ISHFT( LABEL, -32 ), 65535 )
*I64                     KPACK = IAND( ISHFT( LABEL, -16 ), 65535 )
*I64                     LPACK = IAND( LABEL, 65535 )
                       MP = IPACK
                       MQ = JPACK
                       MR = KPACK
                       MS = LPACK
C
            MPQ= IA(MP)+MQ
            MRS= IA(MR)+MS
            IF(MP.EQ.MQ) VAL1 = VAL1+VAL1
            IF(MR.EQ.MS) VAL3 = VAL3+VAL3
            MKL=0
            DO 10 MK=1,NOC
            DO 10 ML=1,NOC
               MKL = MKL+1
               XPQKL(MPQ,MKL) = XPQKL(MPQ,MKL) +
     *               VAL1*(CA(MS,MK)*CB(MR,ML)+CB(MS,ML)*CA(MR,MK))
               XPQKL(MRS,MKL) = XPQKL(MRS,MKL) +
     *               VAL3*(CA(MQ,MK)*CB(MP,ML)+CB(MQ,ML)*CA(MP,MK))
   10       CONTINUE
   20    CONTINUE
C
      ELSE
C     ----- ONE BYTE PER LABEL -----
C
*VDIR NODIVLOOP
         DO 40 M=1,NINT
                       NPACK = NPACKIC+M
            VAL1 = XX(NPACK)
            VAL3 = VAL1
C
*I32                     LABEL = IX(NPACK)
*I32                     IPACK = ISHFT( LABEL, -24 )
*I32                     JPACK = IAND( ISHFT( LABEL, -16 ), 255 )
*I32                     KPACK = IAND( ISHFT( LABEL,  -8 ), 255 )
*I32                     LPACK = IAND( LABEL, 255 )
*I64                     IF ( MOD(NPACK,2) .EQ. 0 ) THEN
*I64                       LABEL = IX( NPACK/2 )
*I64                       IPACK = IAND( ISHFT( LABEL, -24 ), 255 )
*I64                       JPACK = IAND( ISHFT( LABEL, -16 ), 255 )
*I64                       KPACK = IAND( ISHFT( LABEL,  -8 ), 255 )
*I64                       LPACK = IAND( LABEL, 255 )
*I64                     ELSE
*I64                       LABEL = IX( (NPACK/2)+1 )
*I64                       IPACK = ISHFT( LABEL, -56 )
*I64                       JPACK = IAND( ISHFT( LABEL, -48 ), 255 )
*I64                       KPACK = IAND( ISHFT( LABEL, -40 ), 255 )
*I64                       LPACK = IAND( ISHFT( LABEL, -32 ), 255 )
*I64                     END IF
                       MP = IPACK
                       MQ = JPACK
                       MR = KPACK
                       MS = LPACK
C
            MPQ= IA(MP)+MQ
            MRS= IA(MR)+MS
            IF(MP.EQ.MQ) VAL1 = VAL1+VAL1
            IF(MR.EQ.MS) VAL3 = VAL3+VAL3
            MKL=0
            DO 30 MK=1,NOC
            DO 30 ML=1,NOC
               MKL = MKL+1
               XPQKL(MPQ,MKL) = XPQKL(MPQ,MKL) +
     *               VAL1*(CA(MS,MK)*CB(MR,ML)+CB(MS,ML)*CA(MR,MK))
               XPQKL(MRS,MKL) = XPQKL(MRS,MKL) +
     *               VAL3*(CA(MQ,MK)*CB(MP,ML)+CB(MQ,ML)*CA(MP,MK))
   30       CONTINUE
   40     CONTINUE
C
      END IF
C
      RETURN
      END
C
C*MODULE NEONCI   *DECK TRFWTNO
      SUBROUTINE TRFWTNO(NFT,X,NX,XBUF,NREC,
     *                   OUT,MI,MJ,IXBUF,NOC,NINTMX)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      LOGICAL OUT
C
      DIMENSION X(*),XBUF(*),IXBUF(*)
C
      COMMON /PCKLAB/ LABSIZ
C
C        --- FOR PRESENT IJ BEGIN LOOPS OVER K AND L ---
C
      MO = 1
      DO 30 K = 1,NOC
        DO 20 L = 1, NOC
          MK = K
          ML = L
          NX = NX + 1
C
               NPACK = NX
               IPACK = MI
               JPACK = MJ
               KPACK = MK
               LPACK = ML
               IF (LABSIZ .EQ. 2) THEN
*I32             IXBUF( 2*NPACK-1 ) = ISHFT( IPACK, 16 ) + JPACK
*I32             IXBUF( 2*NPACK   ) = ISHFT( KPACK, 16 ) + LPACK
*I64             LABEL = ISHFT( IPACK, 48 ) + ISHFT( JPACK, 32 ) +
*I64 *                   ISHFT( KPACK, 16 ) + LPACK
*I64             IXBUF(NPACK) = LABEL
               ELSE IF (LABSIZ .EQ. 1) THEN
*I32             LABEL = ISHFT( IPACK, 24 ) + ISHFT( JPACK, 16 ) +
*I32 *                   ISHFT( KPACK,  8 ) + LPACK
*I32             IXBUF(NPACK) = LABEL
*I64             IF ( MOD(NPACK,2) .EQ. 0 ) THEN
*I64               LABEL = ISHFT( IPACK, 24 ) + ISHFT( JPACK, 16 ) +
*I64 *                     ISHFT( KPACK,  8 ) + LPACK
*I64               IXBUF( NPACK/2 ) = IXBUF( NPACK/2 ) + LABEL
*I64             ELSE
*I64               LABEL = ISHFT( IPACK, 56 ) + ISHFT( JPACK, 48 ) +
*I64 *                     ISHFT( KPACK, 40 ) + ISHFT( LPACK, 32 )
*I64               IXBUF( (NPACK/2)+1 ) = LABEL
*I64             END IF
               END IF
C
          XBUF(NX)  = X(MO)
          MO = MO + 1
          IF (OUT) CALL TR2OUT(MI,MJ,MK,ML,NX,XBUF(NX),0)
c                  WRITE TO DISK IN CHUNKS OF NINTMX
          IF (NX .EQ. NINTMX) THEN
             CALL PWRIT(NFT,XBUF,IXBUF,NINTMX,NINTMX)
             NREC = NREC+1
             NX = 0
          END IF
   20   CONTINUE
   30 CONTINUE
      RETURN
      END
C
C*MODULE NEONCI  *DECK RDNOCI2N_NE
      SUBROUTINE RDNOCI2N_NE(NFT,X2,NOC1,NOC2,NOC4,IA,XX,IX,NINTMX)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      LOGICAL GOPARR,DSKWRK,MASWRK
C
      DIMENSION X2(NOC1,NOC1,NOC2,NOC2),IA(NOC1),
     *          XX(NINTMX),IX(NINTMX)
C
      COMMON /IOFILE/ IR,IW,IP,IS,IPK,IDAF,NAV,IODA(950)
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
      COMMON /PCKLAB/ LABSIZ
C
      IROW = 0
      DO 110 I=1,NOC1
         IA(I) = IROW
         IROW = IROW+I
  110 CONTINUE
C
      CALL VCLR(X2,1,NOC4)
C
         CALL SEQREW(NFT)
C
C        --- READ TRANSFORMED MIXED NUC ELEC MO INTEGRAL FILE ---
C        --- OR READ TRANSFORMED NUC NUC MO INTEGRAL FILE ---
C
  200    CONTINUE
         CALL PREAD(NFT,XX,IX,NX,NINTMX)
         IF (NX.EQ.0) GO TO 240
         MX = ABS(NX)
         IF (MX.GT.NINTMX) THEN
            IF(MASWRK) WRITE(IW,*) 'INTEGRAL CONFUSION IN -RDCI12-'
            CALL ABRT
         END IF
         DO 220 M = 1,MX
            VAL = XX(M)
            NPACK = M
            IF(LABSIZ .EQ. 2) THEN
*I32           LABEL1 = IX( 2*NPACK - 1 )
*I32           LABEL2 = IX( 2*NPACK     )
*I32           IPACK = ISHFT( LABEL1, -16 )
*I32           JPACK = IAND( LABEL1, 65535 )
*I32           KPACK = ISHFT( LABEL2, -16 )
*I32           LPACK = IAND( LABEL2, 65535 )
*I64           LABEL = IX(NPACK)
*I64           IPACK = ISHFT( LABEL, -48 )
*I64           JPACK = IAND( ISHFT( LABEL, -32 ), 65535 )
*I64           KPACK = IAND( ISHFT( LABEL, -16 ), 65535 )
*I64           LPACK = IAND( LABEL, 65535 )
            ELSE IF (LABSIZ .EQ. 1) THEN
*I32           LABEL = IX(NPACK)
*I32           IPACK = ISHFT( LABEL, -24 )
*I32           JPACK = IAND( ISHFT( LABEL, -16 ), 255 )
*I32           KPACK = IAND( ISHFT( LABEL,  -8 ), 255 )
*I32           LPACK = IAND( LABEL, 255 )
*I64           IF ( MOD(NPACK,2) .EQ. 0 ) THEN
*I64              LABEL = IX( NPACK/2 )
*I64              IPACK = IAND( ISHFT( LABEL, -24 ), 255 )
*I64              JPACK = IAND( ISHFT( LABEL, -16 ), 255 )
*I64              KPACK = IAND( ISHFT( LABEL,  -8 ), 255 )
*I64              LPACK = IAND( LABEL, 255 )
*I64           ELSE
*I64              LABEL = IX( (NPACK/2)+1 )
*I64              IPACK = ISHFT( LABEL, -56 )
*I64              JPACK = IAND( ISHFT( LABEL, -48 ), 255 )
*I64              KPACK = IAND( ISHFT( LABEL, -40 ), 255 )
*I64              LPACK = IAND( ISHFT( LABEL, -32 ), 255 )
*I64           END IF
            END IF
C
C        --- INDICIES READ IN REVERSE CANONICAL ORDER ---
C
            K = IPACK
            L = JPACK
            I = KPACK
            J = LPACK
C
           X2(I,J,K,L) = VAL
  220    CONTINUE
         IF(NX.GT.0) GO TO 200
C
  240    CONTINUE
         CALL SEQREW(NFT)
      RETURN
      END
C
C*MODULE NEONCI   *DECK TR2MCNO
      SUBROUTINE TR2MCNO(XPQRJ,C,X,NUM,NDIM,NORB,NORBS)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      DIMENSION XPQRJ(NORBS,NDIM,NUM),C(NUM,*),X(NUM,*)
C
      PARAMETER (ZERO=0.0D+00, ONE=1.0D+00)
C
      DO 20 MPQ=1,NDIM
         DO 10 MJ=1,NORBS
         DO 10 MR=1,NUM
            X(MR,MJ) = XPQRJ(MJ,MPQ,MR)
   10    CONTINUE
         DO 20 MY=1,NORB
            CALL DGEMV('T',NUM,NORBS,ONE,X(1,1),NUM,
     *                 C(1,MY),1,ZERO,XPQRJ(1,MPQ,MY),1)
   20 CONTINUE
      RETURN
      END

C*MODULE NEONCI  *DECK RDNOCI2E
      SUBROUTINE RDNOCI2E(NFT,X2,NOC1,NOC2,NOC3,IA,XX,IX,NINTMX)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      DIMENSION X2(NOC1,NOC2,NOC2),IA(NOC1),
     *          XX(NINTMX),IX(NINTMX)
C
      COMMON /PCKLAB/ LABSIZ
C
C        --- TWO E- TRANSFORMED INTEGRALS INTO REPLICATED MEMORY ---
C
      IROW = 0
      DO 110 I=1,NOC1
         IA(I) = IROW
         IROW = IROW+I
  110 CONTINUE
C
      CALL VCLR(X2,1,NOC3)
C
C        --- READ TRANSFORMED 2E- INTEGRAL FILE ---
C
         ICOUNT = 1
         CALL PREAD(NFT,XX,IX,NX,NINTMX)
         IF (NX.EQ.0) GO TO 240
         MX = ABS(NX)
         DO 220 M = 1,MX
            VAL = XX(M)
            NPACK = M
            IF(LABSIZ .EQ. 2) THEN
*I32           LABEL1 = IX( 2*NPACK - 1 )
*I32           LABEL2 = IX( 2*NPACK     )
*I32           IPACK = ISHFT( LABEL1, -16 )
*I32           JPACK = IAND( LABEL1, 65535 )
*I32           KPACK = ISHFT( LABEL2, -16 )
*I32           LPACK = IAND( LABEL2, 65535 )
*I64           LABEL = IX(NPACK)
*I64           IPACK = ISHFT( LABEL, -48 )
*I64           JPACK = IAND( ISHFT( LABEL, -32 ), 65535 )
*I64           KPACK = IAND( ISHFT( LABEL, -16 ), 65535 )
*I64           LPACK = IAND( LABEL, 65535 )
            ELSE IF (LABSIZ .EQ. 1) THEN
*I32           LABEL = IX(NPACK)
*I32           IPACK = ISHFT( LABEL, -24 )
*I32           JPACK = IAND( ISHFT( LABEL, -16 ), 255 )
*I32           KPACK = IAND( ISHFT( LABEL,  -8 ), 255 )
*I32           LPACK = IAND( LABEL, 255 )
*I64           IF ( MOD(NPACK,2) .EQ. 0 ) THEN
*I64              LABEL = IX( NPACK/2 )
*I64              IPACK = IAND( ISHFT( LABEL, -24 ), 255 )
*I64              JPACK = IAND( ISHFT( LABEL, -16 ), 255 )
*I64              KPACK = IAND( ISHFT( LABEL,  -8 ), 255 )
*I64              LPACK = IAND( LABEL, 255 )
*I64           ELSE
*I64              LABEL = IX( (NPACK/2)+1 )
*I64              IPACK = ISHFT( LABEL, -56 )
*I64              JPACK = IAND( ISHFT( LABEL, -48 ), 255 )
*I64              KPACK = IAND( ISHFT( LABEL, -40 ), 255 )
*I64              LPACK = IAND( ISHFT( LABEL, -32 ), 255 )
*I64           END IF
            END IF
C
C        --- INDICIES ARE READ IN CANONICAL ORDER ---
C
            I = IPACK
            J = JPACK
            K = KPACK
            L = LPACK
C
           X2(J,K,L) = VAL
          ICOUNT = ICOUNT + 1
  220    CONTINUE
C
  240    CONTINUE
      RETURN
      END
C
C*MODULE NEONCI  *DECK PREIGV
      SUBROUTINE PREIGV(V,M,N,MEIG,LABEO,LABEV,LABPO,LABPV)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      CHARACTER*2 LOC
C
      DIMENSION V(N,M),LABEO(N),LABEV(N),LABPO(N),LABPV(N)
C
      COMMON /IOFILE/ IR,IW,IP,IS,IPK,IDAF,NAV,IODA(950)
      COMMON /NCIFLG/ NUMCNT,NCITRN,NCIPRT,NUMNCI,NUMEXE,NUMEXP,
     *                IEXE(40),IEXP(40),CUTNCI
C
C     ----- PRINT OUT EIGENVECTOR M -----
C     -V- IS -N- ROWS BY -1- COLUMN
C
      WRITE (IW,9008)
      WRITE (IW,9028) MEIG
      WRITE (IW,9029)
      WRITE (IW,9009)
      DO  J = 1,N
          VAL = DABS(V(J,MEIG))
          IF(VAL.LT.CUTNCI) GOTO 500
          ITMP = MOD(J,2)
        IF(ITMP.EQ.0) THEN
          LOC = ' R'
        ELSE
          LOC = ' L'
        ENDIF
        IF(LABEO(J).EQ.LABEV(J) .AND. LABPO(J).EQ.LABPV(J)) THEN
          WRITE (IW,9045) LOC,V(J,MEIG)
        ELSE IF(LABEO(J).EQ.LABEV(J)) THEN
          WRITE (IW,9046) LOC,LABPO(J),LABPV(J),V(J,MEIG)
        ELSE IF(LABPO(J).EQ.LABPV(J)) THEN
          WRITE (IW,9047) LOC,LABEO(J),LABEV(J),V(J,MEIG)
        ELSE
          WRITE (IW,9048) LOC,LABEO(J),LABEV(J),LABPO(J),
     *                LABPV(J),V(J,MEIG)
        ENDIF
 500  CONTINUE
      ENDDO
      RETURN
C
 9008 FORMAT(1X)
 9009 FORMAT(1X,23('-'))
 9028 FORMAT(5X,'  ELEC |     PROT |',5X,10(1X,I4,4X))
 9029 FORMAT(4X,'occ',1X,'exc | ',1X,'occ',1X,'exc |')
 9045 FORMAT(A2,9X,' | ',8X,' |',2X,10F12.7)
 9046 FORMAT(A2,9X,' | ',2I4,' |',2X,10F12.7)
 9047 FORMAT(A2,1X,2I4,' | ',8X,' |',2X,10F12.7)
 9048 FORMAT(A2,1X,2I4,' | ',2I4,' |',2X,10F12.7)
      END
C
C*MODULE NEONCI   *DECK NEONCI_SETUP
      SUBROUTINE NEONCI_SETUP
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      DOUBLE PRECISION NCIINP,NCIRUN
      LOGICAL NCISYM,NUCFCI
C
      COMMON /FMCOM / X(1)
C
      DATA NCIINP/8HNCIINP  /
C
C        --- READ THE $NCIINP INPUT DATA ---
C
      CALL NCIIN(NCIINP,NCIRUN,NCISYM,NUCFCI,LCI)
      LCIT = (LCI*LCI+LCI)/2
      LCI2 = LCI*LCI
C
C        --- GET MEMORY FOR NCI AND OVERLAP MATRICES ---
C
      CALL VALFM(LOADFM)
      ICI        = LOADFM     + 1
      ICIFIN1D   = ICI        + LCIT
      IWRK1      = ICIFIN1D   + LCIT
      ITOTOVLP   = IWRK1      + LCI
      ISDIAG     = ITOTOVLP   + LCIT
      ISEIGEN    = ISDIAG     + LCI2
      ISWRK      = ISEIGEN    + LCI
      IISWRK     = ISWRK      + LCI*8
      ISVEC      = IISWRK     + LCI
      ISX1       = ISVEC      + LCI2
      ISX2       = ISX1       + LCI2
      ICEIGEN    = ISX2       + LCI2
      IWRK2      = ICEIGEN    + LCI
      IIWRK2     = IWRK2      + LCI*8
      IVINDIAG   = IIWRK2     + LCI
      LAST       = IVINDIAG   + LCI2
      NEED       = LAST - LOADFM -1
      CALL GETFM(NEED)
C
      CALL NEONCI_RUN(LCI,LCIT,LCI2,NCIRUN,NCISYM,NUCFCI,
     *                X(ICI),X(ICIFIN1D),X(IWRK1),
     *                X(ITOTOVLP),X(ISDIAG),X(ISEIGEN),X(ISWRK),
     *                X(IISWRK),X(ISVEC),X(ISX1),X(ISX2),X(ICEIGEN),
     *                X(IWRK2),X(IIWRK2),X(IVINDIAG))
C
      CALL RETFM(NEED)
C
      RETURN
      END
C
C*MODULE NEONCI   *DECK NCIIN
      SUBROUTINE NCIIN(NCIINP,NCIRUN,NCISYM,NUCFCI,NUMDET)
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      DOUBLE PRECISION NCIINP,NCIRUN
      LOGICAL GOPARR,DSKWRK,MASWRK,OK
      LOGICAL NCISYM,NUCFCI
C
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
      COMMON /NCIFLG/ NUMCNT,NCITRN,NCIPRT,NUMNCI,NUMEXE,NUMEXP,
     *                IEXE(40),IEXP(40),CUTNCI
      COMMON /NUCMOI/ NUNIQN,IUNIQN(20),IUNIQT(20),NQMNUC,IQMNUC(20),
     *                 IQNTYP(20),NUMNB,NUCST,NAUXNB,IAUXNB(20),NUMULT,
     *                 NNA,NNB,NTAUXB
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
      PARAMETER (NNAM=12)
      DIMENSION QNAM(NNAM),KQNAM(NNAM)
      DATA BLANK,CALC,CHECK,DEBUG/8H        ,8HCALC    ,8HCHECK   ,
     *     8HDEBUG   /
      DATA QNAM/8HNUMCNT  ,8HNCIRUN  ,8HNCITRN  ,8HNCIPRT  ,8HNUMNCI  ,
     *          8HNUMEXE  ,8HNUMEXP  ,8HNCISYM  ,8HIEXE    ,8HIEXP    ,
     *          8HNUCFCI  ,8HCUTNCI  /
      DATA KQNAM/1,5,1,1,1,1,1,0,401,401,0,3/
C
C        --- SETUP INPUT TO CONTROL THE NEONCI CALCULATION ---
C
      NUMCNT = 2
      NCIRUN = BLANK
      NCITRN = 1
      NCIPRT = 2
      NUMNCI = 2
      NUMEXE = 0
      NUMEXP = 0
      NCISYM = .FALSE.
      NUCFCI = .FALSE.
      CUTNCI = 1.0D-2
C
C        --- READ THE NEONCI CONTROL INPUT ---
C
      JRET=0
      CALL NAMEIO(IR,JRET,NCIINP,NNAM,QNAM,KQNAM,
     *            NUMCNT,NCIRUN,NCITRN,NCIPRT,NUMNCI,NUMEXE,NUMEXP,
     *            NCISYM,IEXE,IEXP,NUCFCI,CUTNCI,
     *   0,0,
     *   0,0,0,0,0,  0,0,0,0,0,
     *   0,0,0,0,0,  0,0,0,0,0,  0,0,0,0,0,  0,0,0,0,0,
     *   0,0,0,0,0,  0,0,0,0,0,  0,0,0,0,0,  0,0,0,0,0)
C
      IF(JRET .EQ. 2) THEN
         IF (MASWRK) WRITE (IW,9008)
         CALL ABRT
      END IF
C
C     ---- SET THE TOTAL NUMBER OF DETERMINANTS IN CI ----
C
      IF(NUCFCI) NUMEXP = NUMNB-1
        NUMDET = (NUMEXE+1)*(NUMEXP+1)*2
C
C     ---- CHECK SPELLINGS, AND ALSO SET THE OPTION FLAGS ----
C
      NERR = 0
      OK = .FALSE.
      IF(NCIRUN.EQ.BLANK)  NCIRUN = CALC
      IF(NCIRUN.EQ.CALC)   OK = .TRUE.
      IF(NCIRUN.EQ.DEBUG)  OK = .TRUE.
      IF(NCIRUN.EQ.CHECK)  OK = .TRUE.
      IF(.NOT.OK) THEN
         IF (MASWRK) WRITE(IW,9010) 'NCIRUN',NCIRUN
         NERR = NERR+1
      END IF
      IF(NERR.GT.0) THEN
         IF (MASWRK) THEN
         WRITE(IW,*) 'YOUR $NCIINP HAS AT LEAST ONE SPELLING MISTAKE.'
         WRITE(IW,*) 'THE JOB IS STOPPING SO YOU CAN FIX IT.'
         END IF
         CALL ABRT
         STOP
      END IF
C
C        --- CHECK THE NUMBER OF NUCLEAR CENTERS ---
C
      IF(NUMCNT.LT.1 .OR. NUMCNT.GT.2) THEN
         WRITE(IW,9011)
         CALL ABRT
      END IF
C
C        --- CHECK THE TRANSFORMATION METHOD ---
      IF(NCITRN.LT.1 .OR. NCITRN.GT.2) THEN
         WRITE(IW,9012)
         CALL ABRT
      END IF
C        --- CHECK THE PRINT OPTION ---
      IF(NCIPRT.LT.1 .OR. NCIPRT.GT.3) THEN
         WRITE(IW,9014)
         CALL ABRT
      END IF
C        --- CHECK THAT NUMNCI DOES NOT EXCEED NUMDET ---
      IF(NUMNCI.GT.NUMDET) THEN
         WRITE(IW,9015)
         CALL ABRT
      END IF
C        --- CHECK THAT NUMDET NOT > 2 IF NCISYM IS TRUE ---
      IF(NCISYM .AND. NUMDET.GT.2) THEN
         WRITE(IW,9016)
         CALL ABRT
      END IF
C        --- CHECK THAT IF NUCFCI IS TRUE THAT NQMNUC = 1 ---
      IF(NUCFCI .AND. NQMNUC.GT.1) THEN
         WRITE(IW,9017)
         CALL ABRT
      END IF
C
      WRITE(IW,9018) NUMCNT,NCIRUN,NCITRN,NCIPRT,NUMNCI,NUMEXE,
     *               NUMEXP,NCISYM,NUCFCI,CUTNCI

      RETURN
 9008 FORMAT(40H ERROR IN $NCIINP NAMELIST INPUT - STOP )
 9010 FORMAT(1X,'*** ERROR! KEYWORD ',A6,
     *     ' WAS GIVEN AN ILLEGAL VALUE ',A8,'.')
 9011 FORMAT(1X,' *** CURRENTLY ONLY ONE AND TWO CENTER CALCULATIONS',
     *       'ARE AVALIABLE ***'/
     *       1X,' FIX NUMCNT IN NCIINP GROUP TO BE EITHER 1 OR 2')
 9012 FORMAT(1X,' *** CURRENTLY ONLY ONE STEP AND TWO STEP ',
     *       'NONORTHOGONAL AO TO MO TRANSFORMATIONS ARE AVALIABLE FOR',
     *       ' ELECTRONS ***'/
     *       1X,' FIX NCITRN IN NCIINP GROUP TO BE EITHER 1 OR 2')
 9014 FORMAT(1X,' *** PRINT OPTION SELECTED IS NOT VALID ***',/
     *       1X,' FIX NCIPRT IN NCIINP GROUP' )
 9015 FORMAT(1X,' *** NUMBER OF CI STATES TO FIND EXCEEDS THE NUMBER ',
     *      'OF DETERMINANTS USED TO CONSTRUCT THE NEO-NCI MATRIX ***',/
     *       1X,' FIX EITHER NUMNCI OR NCIDET IN THE NCIINP GROUP' )
 9016 FORMAT(1X,' *** NCISYM = .TRUE. AND NUMEXE OR NUMEXP ARE NOT ',
     *      'EQUAL TO ZERO ***',/
     *       1X,' THE NCISYM OPTION CAN NOT BE USED WITH EXCITATIONS')
 9017 FORMAT(1X,' *** NUCFCI = .TRUE. AND NQMNUC > 1 ***',/
     *       1X,' *** YOU CAN ONLY USE THE OPTION NUCFCI FOR ONE',
     *          ' QUANTUM PROTON ***')
 9018 FORMAT(/5X,'$NCIINP OPTIONS'/5X,15(1H-)/
     *  5X,7HNUMCNT=,I8,5X,7HNCIRUN=,A8,5X,7HNCITRN=,I8/
     *  5X,7HNCIPRT=,I8,5X,7HNUMNCI=,I8,5X,7HNUMEXE=,I8/
     *  5X,7HNUMEXP=,I8,5X,7HNCISYM=,L8,5X,7HNUCFCI=,L8/
     *  5X,7HCUTNCI=,E8.1)
      END
C*MODULE NEONCI   *DECK GETIMAP
      SUBROUTINE GETIMAP(IOCCE,IVIRE,IOCCP,IVIRP,IROP,IRVP,ICOP,ICVP,
     *                   IROE,IRVE,ICOE,ICVE,NMT,NMT2,NUCFCI,DBUG)
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      LOGICAL DBUG,NUCFCI
      integer dec
C
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
      COMMON /NCIFLG/ NUMCNT,NCITRN,NCIPRT,NUMNCI,NUMEXE,NUMEXP,
     *                IEXE(40),IEXP(40),CUTNCI
C
      DIMENSION IOCCE(NMT2),IVIRE(NMT2),
     *          IOCCP(NMT2),IVIRP(NMT2),
     *          IROP(NMT2),IRVP(NMT2),ICOP(NMT2),ICVP(NMT2),
     *          IROE(NMT2),IRVE(NMT2),ICOE(NMT2),ICVE(NMT2)
C
C     CREATE MAP OF VECTOR INDICIES TO CONSTRUCT NCI MATRIX
C     THE ORDER IS THAT PROTON EXCITATIONS ARE IN THE INNER
C     LOOP AND ELECTRON EXCITATIONS ARE THE OUTER LOOP SO
C     THAT ALL PROTON EXCITATIONS ARE DONE FOR THE GROUND
C     ELECTRONIC STATE THEN THE ALL PROTON EXCITATIONS ARE
C     DONE FOR THE FIRST EXCITED ELECTRONIC STATE AND SO ON
C
           IOCCE(1) = 1
           IVIRE(1) = 1
        IF(NUMEXE.GT.0) THEN
             II = 2
              J = 1
          DO I = 1, NUMEXE
            IOCCE(II) = IEXE(J)
            IVIRE(II) = IEXE(J+1)
             J = J + 2
            II = II+1
          ENDDO
        ENDIF
C
           IOCCP(1) = 1
           IVIRP(1) = 1
        IF(NUMEXP.GT.0) THEN
          IF(NUCFCI) THEN
            DO I = 1,NUMEXP+1
              IOCCP(I) = 1
              IVIRP(I) = I
            ENDDO
          ELSE
               II = 2
                J = 1
            DO I = 1, NUMEXP
              IOCCP(II) = IEXP(J)
              IVIRP(II) = IEXP(J+1)
                J = J + 2
              II = II+1
            ENDDO
          ENDIF
        ENDIF
C
C     --- MAP FOR PROT INDICIES <PVEC|---
C
      NINC1 =  NUMEXP+1
C
      L  = 1
      K  = 1
      KK = 1
      II = 1
      NCOUNT = NMT
C
      DO I = 1,NMT   ! LOOP OVER COLUMNS
        DO J = I,NMT ! LOOP OVER ROWS
           IROP(II) =  IOCCP(K)
           IRVP(II) =  IVIRP(K)
         IF(L.EQ.2) THEN ! L DETERMINES WHETHER OR NOT
           L = 1         ! TO UPDATE K
           K = K + 1
           IF(K.GT.NINC1) K=1 ! END OF PROT EXC FOR THIS COLUMN
         ELSE
           L = L + 1
         ENDIF
        II = II +1
        ENDDO
C   DETERMINE WHERE TO START K FOR NEXT COLUMN
          NOLD  = NCOUNT
          NCOUNT = NCOUNT - 1
          L = 1 + MOD(NCOUNT,2)
          KK = KK + MOD(NOLD,2)
          IF(KK.GT.NINC1) KK = 1
          K = KK
      ENDDO
C
C     --- MAP FOR PROT INDICIES |PVEC>
C
      K = 1
      II = 1
      L  = 2
      NLEN = NMT
      DO I = 1, NMT
        DO J = I, NMT
           ICOP(II) =  IOCCP(K)
           ICVP(II) =  IVIRP(K)
          II = II + 1
        ENDDO
              K = K + MOD(NLEN,2)
              NLEN = NLEN - 1
        IF(K.GT.NINC1) K = 1
      ENDDO
C
C     --- MAP FOR ELEC INDICIES <EVEC|---
C
      IINC1 =  NUMEXE+1
      IINC2 =  NMT/IINC1

      K   = 1
      KK  = 1
      L   = 1
      II  = 1
      INC = 1
      DEC = 1
C
      DO I = 1,NMT   ! LOOP OVER COLUMNS
        DO J = I,NMT ! LOOP OVER ROWS
           IROE(II) =  IOCCE(K)
           IRVE(II) =  IVIRE(K)
           IF(L.EQ.IINC2) THEN ! NEXT ELEC EXC IN THIS COLUMN
              K = K + 1
              L = 1
           ELSE
              L = L + 1
           ENDIF
        II = II + 1
        ENDDO
C   DETERMINE WHERE TO START K FOR NEXT COLUMN
          IF(INC.EQ.IINC2) THEN
             KK = KK + 1
             INC = 1
          ELSE
             INC = INC + 1
          ENDIF
             K = KK
             L = DEC+1
             IF(L.EQ.IINC2) THEN
               DEC = 0
             ELSE
               DEC = DEC + 1
             ENDIF
      ENDDO
C
C     --- MAP FOR ELEC INDICIES |EVEC>
C
      IINC1 =  NUMEXE+1
      IINC2 =  NMT/IINC1
      K = 1
      II = 1
      INC = 1
      DO I = 1, NMT
        DO J = I, NMT
           ICOE(II) =  IOCCE(K) ! ICOE IS COLUMN OCCUPIED ELEC ORB
           ICVE(II) =  IVIRE(K) ! ICVE IS COLUMN VIRTUAL ELEC ORB
           II = II + 1
        ENDDO
        IF (INC.EQ.IINC2) THEN
            K = K + 1
            INC = 1
        ELSE
              INC = INC + 1
        ENDIF
      ENDDO
C
      IF(DBUG) THEN
        WRITE(IW,*) " ******************** "
        WRITE(IW,*) "     PROTON MAPS      "
        WRITE(IW,*) " ******************** "
        DO I =1,NMT2
          WRITE(IW,*) "IROP I", IROP(I), I
        ENDDO
        WRITE(IW,*) "  "
        DO I =1,NMT2
          WRITE(IW,*) "IRVP I", IRVP(I), I
        ENDDO
        WRITE(IW,*) "  "
        DO I =1,NMT2
          WRITE(IW,*) "ICOP I", ICOP(I), I
        ENDDO
        WRITE(IW,*) "  "
        DO I =1,NMT2
          WRITE(IW,*) "ICVP I", ICVP(I), I
        ENDDO
        WRITE(IW,*) "  "
        WRITE(IW,*) " ******************** "
        WRITE(IW,*) "    ELECTRON MAPS     "
        WRITE(IW,*) " ******************** "
        DO I =1,NMT2
          WRITE(IW,*) "IROE I", IROE(I), I
        ENDDO
        WRITE(IW,*) "  "
        DO I =1,NMT2
          WRITE(IW,*) "IRVE I", IRVE(I), I
        ENDDO
        WRITE(IW,*) "  "
        DO I =1,NMT2
          WRITE(IW,*) "ICOE I", ICOE(I), I
        ENDDO
        WRITE(IW,*) "  "
        DO I =1,NMT2
          WRITE(IW,*) "ICVE I", ICVE(I), I
        ENDDO
        WRITE(IW,*) "  "
      ENDIF
C
      RETURN
      END
