C 17 JAN 06 - MWS - UPDATE ARGS TO TRFMCX
C 14 NOV 05 - DGF - PAD COMMON BLOCKS ENRGYS AND SCFOPT
C 19 SEP 05 - MWS - ADJUSTMENT FOR H OR I FUNCTION RUNS
C  5 JUL 05 - MWS - SELECT NEW ATOM,BASIS,EFP,PCM,DAF DIMENSIONS
C  1 JUN 05 - MWS - PAD THE COMMON TRFOPT
C 13 FEB 05 - MWS - PAD COMMON BLOCK NSHEL
C  5 FEB 05 - MWS - PAD COMMON MCINP
C  7 SEP 04 - MWS - PAD COMMON BLOCK INTFIL
C 16 JAN 04 - MWS - SYNCH ARGUMENTS FOR NEWTON CALL
C  9 DEC 03 - MWS - SYNCH COMMON BLOCK RUNOPT
C  7 AUG 02 - MWS - SYNC CALLING ARGS TO NEWTON
C 22 MAY 02 - MWS - PARTLY ENABLE USE OF DDI INTEGRAL TRANSFORMATION
C  8 OCT 01 - MWS - MCQDU: FIX THE CASE OF USING SPHERICAL HARMONICS
C 19 SEP 01 - MWS - CONVERT MXAOCI PARAMTER TO MXAO
C  6 SEP 01 - MWS - ADD DUMMY ARGUMENTS TO NAMEIO CALL
C 25 JUN 01 - MWS - ALTER COMMON BLOCK SCFOPT
C 26 OCT 00 - MWS - INTRODUCE MXAO PARAMETER
C 11 JUN 00 - MWS - MCQUD: ALWAYS ALLOCATE HONDO INTEGRAL BUFFER
C  1 MAY 00 - MWS - NO SYMMETRY OPTION IN TRFOPT COMMON
C 21 DEC 99 - MWS - CHANGE INTEGRAL SYMMETRY TEST
C 12 NOV 98 - GDF - CHANGE BIT PACKING TO ISHIFT
C 27 OCT 98 - MAF - MCQUD: ALLOW FOR USE OF SPHERICAL HARMONICS
C 13 APR 98 - MWS - CHANGES RELATING TO USING DETERMINANT CI ELSEWHERE
C 28 SEP 97 - MWS - CONVERT PARALLEL CALLS FROM TCGMSG TO DDI
C 21 APR 97 - MWS - MCQUD: CHANGE CALL TO NEWTON TO PASS FROZEN MOS.
C  2 APR 97 - GMC - QUDCUP: POSITION TO CORRECT EXCITED STATE EIGENVEC
C  8 JAN 97 - GMC - NEW MODULE FROM HONDO FOR FULL QUADRATIC MCSCF
C
C*MODULE MCQUD   *DECK MCQUD
      SUBROUTINE MCQUD
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      DOUBLE PRECISION METHOD
C
      PARAMETER (MXSH=5000, MXGTOT=20000, MXGSH=30, MXG2=MXGSH*MXGSH,
     *           MXATM=2000, MXRT=100)
C
      LOGICAL MICCVG,CVGED,CVGING,DMPING,FINAL,MCNOCI
      LOGICAL SOME,GOPARR,DSKWRK,MASWRK,TDSKWK,DIRTRF,PACK2E
      LOGICAL CANONC,FCORE,FORS,NOCI,EKT,LINSER
      LOGICAL NOTRF,INITAL,DOFOCK,DOEXCH,QUD,FOCAS,SOSCF,FULLNR,DROPC
      LOGICAL DDITRF,DOOOOO,DOVOOO,DOVVOO,DOVOVO,DOVVVO,DOVVVV,DOCORE
C
      PARAMETER (MXNORO=250, MXAO=8192)
C
      COMMON /CASOPT/ CASHFT,CASDII,NRMCAS,FOCAS,SOSCF,DROPC
      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 /FUNCT / E,EGRAD(3,MXATM)
      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,IJKO,IJKT,IDAF,NAV,IODA(950)
      COMMON /MCINP / METHOD,CISTEP,ACURCY,ENGTOL,DAMP,MICIT,NWORD,NORB,
     *                NOROT(2,MXNORO),MOFRZ(15),NPFLG(10),
     *                NOFO,MCFMO,CANONC,FCORE,FORS,NOCI,EKT,LINSER
      COMMON /MCPAR / NFZC,NCORBS1,NACT,NORBS1,NORB1,NTOT
      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),MIN(MXSH),MAX(MXSH),NSHELL
      COMMON /NTNOPT/ ENERGY,ENERG0,DEMAX,SQCDF,ITER,MICRO,NOTRF
      COMMON /ORBSET/ NORBMX,NORBS,NCORBS,NLEVS,NNA,NNB,NNC,NSYM,MSYM,
     *                IDOCC,IVAL,IMCC,ISYM(MXAO),JCODE(MXAO),
     *                NLCS(MXAO),LEVPT(MXAO),LEVNR(MXAO),IOUT(MXAO),
     *                NREFS,IEXCT,NFOCI,INTACT
      COMMON /OUTPUT/ NPRINT,ITOL,ICUT,NORMF,NORMP,NOPK
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
      COMMON /QUDMC / QUDTHR,QUD
      COMMON /RESTAR/ TIMLIM,IREST,NREC,INTLOC,IST,JST,KST,LST
      COMMON /RUNOPT/ RUNTYP,EXETYP,NEVALS,NGLEVL,NHLEVL
      COMMON /SCFOPT/ CONVHF,MAXIT,MCONV,NPUNCH,NPREO(4)
      COMMON /TRFOPT/ CUTTRF,NWDTRF,MPTRAN,ITRFAO,NOSYMT,IPURTF,DIRTRF
C
      PARAMETER (CVGTOL=5.0D-03, DMPTOL=0.2D+00, TEN=1.0D+01,
     *           TEN7=1.0D-07, PT2=0.2D+00, TWOPT2=2.2D+00,
     *           ZERO=0.0D+00, TWO=2.0D+00)
C
      SAVE INITAL
C
      DATA DRTINP/8HDRT     /
      DATA CHECK/8HCHECK   /
      DATA OPTMZE/8HOPTIMIZE/
      DATA INITAL/.TRUE./
C
C     ----- MAIN DRIVER FOR FULLY QUADRATIC MCSCF -----
C
      SOME=NPRINT.NE.-5 .AND. MASWRK
      NPSAVE = NPRINT
      IF(SOME) WRITE(IW,9000)
      IF(SOME) WRITE(IW,9010) ENUCR
C
      IREST=2
      EMC =ZERO
      EMC0=ZERO
      ITER=0
      MICRO=0
      DE =ZERO
      DEAVG=ZERO
      EPSLON=TEN7
C
      IF(DAMP.LE.ZERO) DAMP=ZERO
      DMPING=DAMP.GT.ZERO
      CVGING=.FALSE.
      CVGED =.FALSE.
      MICCVG=.FALSE.
      FINAL =.FALSE.
      MCNOCI=.FALSE.
      NOTRF=.FALSE.
      DOFOCK=.TRUE.
      DOEXCH=.FALSE.
      FULLNR=.FALSE.
C
      TOLLAG = ACURCY
      TOLENG = ENGTOL
      IF(RUNTYP.EQ.OPTMZE  .AND.  NEVALS.GT.0) THEN
         GRMS = DDOT(3*NAT,EGRAD,1,EGRAD,1)
         GRMS = SQRT(GRMS/(3*NAT))
                               CFACT =  1.0D+00
         IF(GRMS.GT.0.005D+00) CFACT =  2.5D+00
         IF(GRMS.GT.0.020D+00) CFACT = 10.0D+00
         IF(GRMS.GT.0.100D+00) CFACT = 20.0D+00
         TOLLAG = CFACT*TOLLAG
         TOLENG = CFACT*TOLENG
      END IF
C
C        FROZEN CORE RUNS DO NOT LEAD TO A SYMMETRIC LAGRANGIAN
C
      IF(FCORE) THEN
         TOLNG2=TOLENG*TEN*TEN
         TOLLG2=TEN*TEN
      ELSE
         TOLNG2=TOLENG*TEN*TEN
         TOLLG2=TOLLAG*TWO
      END IF
C
C         NOTES ON PARALLEL EXECUTION (MAY 17, 2002)
C
C         THIS PROGRAM IS PARTIALLY ENABLED FOR PARALLEL EXECUTION,
C         MEANING THAT IT COMPLETES THE NON-QUADRATIC PRELIMINARY
C         ITERATIONS CORRECTLY WITH THE PRESENT -DDITRF- CODE.
C         AFTER BEGINNING TO RUN THE FULLY QUADRATIC ITERATIONS,
C         THE PROGRAM CORRECTLY DOES A DISTRIBUTED TRANSFORMATION,
C         BUT THEN PRINTS FOUR DDI_DESTROY MESSAGES.  THE GUGA
C         SORTING OF THE OOOO INTEGRALS APPEARED TO RUN CORRECTLY,
C         HOWEVER, AND THE HAMILTONIAN GENERATION PREPARED THE
C         CORRECT NUMBER OF LOOPS.  THE ROUTINE -QUDCUP- FAILED
C         TRYING TO READ TRANSFORMED INTEGRALS FROM DISK, IN
C         THE CALL TO -M1LGRR-.  THIS SEEKS VOOO INTEGRALS, AND
C         NEEDS TO HAVE AN OPTION TO PULL THIS CLASS OF INTEGRALS
C         FROM DISTRIBUTED MEMORY.  EXAMPLES OF USING VOOO INTEGRALS
C         CAN BE FOUND IN MCTWO.SRC.  IT LOOKS LIKE -QUDCUP- MAY BE
C         THE ONLY ROUTINE THAT READS MO INTEGRALS, SO FIXING THIS
C         MIGHT VERY WELL ENABLE THE PARALLEL QUADRATIC CONVERGER.
C         HOWEVER, THE DDI_DESTROY MESSAGES ARE A MYSTERY, AS THEY
C         SHOULD NOT BE DESTROYED UNTIL AFTER A FULLY QUADRATIC
C         ITERATION IS FINISHED.
C
C         UNTIL THESE SNAGS ARE SOLVED, THE INPUT READING -MCIN-
C         ROUTINE HAS AN ABORT OF PARALLEL QUAD JOBS.
C
C         MY TEST INPUT WAS THE FOLLOWING, PLUS EXAM09'S $VEC GROUP:
C  !    THE FINAL E= -75.6386218845, 10 ITERS, RMS GRADIENT=0.021601818
C  !
C   $CONTRL SCFTYP=MCSCF RUNTYP=GRADIENT $END
C   $SYSTEM TIMLIM=8 MWORDS=5 MEMDDI=5 $END
C   $BASIS  GBASIS=N21 NGAUSS=3 $END
C   $GUESS  GUESS=MOREAD  NORB=13 $END
C   $MCSCF  QUAD=.T. CISTEP=GUGA $END
C   $DET    NCORE=1 NACT=6 NELS=8 $END
C   $DRT    GROUP=C2V FORS=.T. NMCC=1 NDOC=4 NVAL=2 $END
C
C   $DATA
C  WATER...3-21G BASIS...FORS-MCSCF...EXPERIMENTAL GEOMETRY
C  CNV 2
C
C  OXYGEN     8.0
C  HYDROGEN   1.0   0.0   0.7572157   0.5865358
C   $END
C
      DSKWRK=.TRUE.
      IF (GOPARR) THEN
         DDITRF = .TRUE.
         DOOOOO = .TRUE.
         DOVOOO = .TRUE.
         DOVVOO = .TRUE.
         DOVOVO = .TRUE.
         DOVVVO = .FALSE.
         DOVVVV = .FALSE.
         DOCORE = .TRUE.
      ELSE
         DDITRF = .FALSE.
         DOOOOO = .FALSE.
         DOVOOO = .FALSE.
         DOVVOO = .FALSE.
         DOVOVO = .FALSE.
         DOVVVO = .FALSE.
         DOVVVV = .FALSE.
         DOCORE = .TRUE.
      ENDIF
C
      TIM0=ZERO
      TIM1=ZERO
      CALL TSECND(TIM0)
      TLEFTS = TIMLIM-TIM0
C
C     ----- CONSTRUCT -DRT- TABLE -----
C
      CALL DRTGEN(NPFLG(1),DRTINP)
C
      IF(DROPC) NORBS = NORBS + NCORBS
C
C     CHEAT SHEET FOR THE ORBITAL COUNTERS:
C        NTOT IS THE TOTAL NUMBER OF AO-S
C        NORB IS THE TOTAL NUMBER OF MO-S (NORB.LE.NTOT)
C        NCORBS IS THE NUMBER OF CORE ORBITALS
C        NACT   IS THE NUMBER OF ACTIVE ORBITALS
C        NORBS  IS THE TOTAL NUMBER OF OCCUPIED ORBITALS (NCORBS+NACT)
C
      NTOT = NUM
      NORB = NQMT
      NORB1 = NORB
      NORBS1 = NORBS
      NCORBS1 = NCORBS
      NACT = NORBS - NCORBS
      TDSKWK = DSKWRK
      DSKWRK = .TRUE.
C
C     ----- OPEN DISK FILE 23 FOR ORBITAL HESSIAN AND GRADIENT -----
C
      IF(INITAL) THEN
         CALL SEQOPN(23,'JKFILE','UNKNOWN',.FALSE.,'UNFORMATTED')
         INITAL = .FALSE.
      ELSE
         CALL SEQREW(23)
      END IF
C
C     ----- EXCHANGE INTEGRALS FOR DIRECT TRANSFORMATION -----
C
      IF(DIRTRF) THEN
         NSH2 = (NSHELL*NSHELL+NSHELL)/2
         CALL BASCHK(LMAX)
                       NANGM =  4
         IF(LMAX.EQ.2) NANGM =  6
         IF(LMAX.EQ.3) NANGM = 10
         IF(LMAX.EQ.4) NANGM = 15
         IF(LMAX.EQ.4) NANGM = 21
         IF(LMAX.EQ.4) NANGM = 28
         MAXG = NANGM**4
         CALL VALFM(LOADFM)
         LXINTS = LOADFM + 1
         LGHOND = LXINTS + NSH2
         LDDIJ  = LGHOND + MAXG
         LAST   = LDDIJ  + 49*MXG2
         NEED = LAST - LOADFM - 1
         CALL GETFM(NEED)
         IF(MASWRK) WRITE(IW,9028)
         CALL EXCHNG(X(LXINTS),X(LGHOND),X(LDDIJ),NSH2,MAXG,INTTYP)
         CALL RETFM(NEED)
      END IF
C
C           **********************
C     ***** BEGIN MCSCF ITERATIONS *****
C           **********************
C
      IF(EXETYP.EQ.CHECK) MAXIT=2
  100 CONTINUE
      ITER=ITER+1
C
C     ----- GET CI WAVEFUNCTION AND ENERGY (E IS SET BY GUG2DM) -----
C
      IF (DROPC) THEN
         CALL TRFMCX(NPFLG(2),NCORBS,NORBS,NORB,DOFOCK,DOEXCH,
     *               DDITRF,DOOOOO,DOVOOO,DOVVOO,DOVOVO,
     *               DOVVVO,DOVVVV,DOCORE)
      ELSE
         CALL TRFMCX(NPFLG(2),0,0,0,.FALSE.,.FALSE.,
     *               DDITRF,DOOOOO,DOVOOO,DOVVOO,DOVOVO,
     *               DOVVVO,DOVVVV,DOCORE)
      END IF
      CALL GUGSRT(NPFLG(3),DDITRF)
      CALL GUGAEM(NPFLG(3))
      CALL GUGADG(NPFLG(4))
      CALL GUG2DM(NPFLG(6))
C
      EMC0=EMC
      EMC =E
      ETOT=E
      DE0=DE
      DE =EMC-EMC0
      IF(ITER.EQ.1) DEAVG=ZERO
      IF(ITER.EQ.2) DEAVG= ABS(DE)
      IF(ITER.GE.3) DEAVG=( ABS(DE)+ ABS(DE0)+PT2*DEAVG)/TWOPT2
C
      EPSLON=EPSLON/TEN
      IF(ABS(DE).LE.ACURCY) EPSLON=TEN7
      IF(EPSLON.LT.TEN7) EPSLON=TEN7
      TOLMIC= ABS(EMC)* SQRT(EPSLON)
C
      CVGING=DE.LT.CVGTOL.AND.DAMP.LT.DMPTOL
      IF(DMPING.AND.     CVGING                 ) DAMP=ZERO
      IF(DMPING.AND..NOT.CVGING.AND.DAMP.EQ.ZERO) DAMP=DMPTOL
C
      IF(ITER.GT.2.AND..NOT.CVGING) CALL NTNDMP(DE,DE0,DEAVG,DAMP)
      IF(.NOT.DMPING.AND.DAMP.GT.ZERO) DMPING=.TRUE.
C
C     ***** START MICRO-ITERATION *****
C
      ENERG0=ZERO
      MICRO=0
  500 CONTINUE
      MICRO=MICRO+1
      NOTRF=.FALSE.
C
C     ----- MCSCF ORBITAL IMPROVEMENT -----
C
      CALL NEWTON(CISTEP,NPFLG(8),METHOD,ACURCY,DAMP,EKT,
     *            CANONC,NPUNCH,FCORE,FORS,NOCI,
     *            NWORD,NORB,NOROT,MXNORO,IROTMX,JROTMX,MOFRZ,
     *            FINAL,MCNOCI,DOEXCH,DROPC,FULLNR)
C
      CALL TSECND(TIM1)
      TIM0=TIM1
      IF((ITER.EQ.1) .AND. MASWRK) WRITE(IW,9030)
      IF (MASWRK) WRITE(IW,9040) ITER,EMC,DE,DEMAX,SQCDF,MICRO,DAMP
      CALL FLSHBF(IW)
C
      NPRINT = -5
      IF(NPFLG(1).LE.0) NPFLG(1)=-5
      IF(NPFLG(2).LE.0) NPFLG(2)=-5
      IF(NPFLG(3).LE.0) NPFLG(3)=-5
      IF(NPFLG(4).LE.0) NPFLG(4)=-5
      IF(NPFLG(6).LE.0) NPFLG(6)=-5
      IF(NPFLG(8).LE.0) NPFLG(8)=-5
C
      MICCVG=DEMAX.LT.TOLMIC.OR.ENERGY.GT.ENERG0
      ENERG0=ENERGY
      IF(MICCVG.OR.MICRO.GE.MICIT) GO TO 510
      IF (DROPC) THEN
         CALL TRFMCX(NPFLG(2),NCORBS,NORBS,NORB,DOFOCK,DOEXCH,
     *               DDITRF,DOOOOO,DOVOOO,DOVVOO,DOVOVO,
     *               DOVVVO,DOVVVV,DOCORE)
      ELSE
         CALL TRFMCX(NPFLG(2),0,0,0,.FALSE.,.FALSE.,
     *               DDITRF,DOOOOO,DOVOOO,DOVVOO,DOVOVO,
     *               DOVVVO,DOVVVV,DOCORE)
      END IF
      GO TO 500
C
C     ***** END OF MICRO-ITERATIONS *****
C
  510 CONTINUE
      IF(EXETYP.EQ.CHECK) GO TO 400
C
C     ----- TEST FOR CONVERGENCE -----
C
      CVGED = CVGING .AND. ABS(DE).LT.TOLNG2 .AND. DEMAX.LT.TOLLAG
      IF(CVGED .AND. MASWRK) WRITE(IW,9050)
      IF(CVGED) GO TO 1000
      CVGED = CVGING .AND. ABS(DE).LT.TOLENG .AND. DEMAX.LT.TOLLG2
      IF(CVGED .AND. MASWRK) WRITE(IW,9060)
      IF(CVGED) GO TO 1000
C
C     ----- TEST FOR TIME LIMIT -----
C
      CALL TSECND(TIM1)
      TLEFT = TIMLIM-TIM1
      IF(TLEFT .LT.  1.5D+00*(TLEFTS-TLEFT)/ITER) THEN
         IF (MASWRK) THEN
            WRITE(IW,9070)
            WRITE(IW,9090)
         END IF
         EMC=ZERO
         E  =ZERO
         GO TO 1000
      END IF
C
C     ---- TEST FOR ITERATION LIMIT -----
C
      IF(ITER.GE.MAXIT) THEN
         IF (MASWRK) THEN
            WRITE(IW,9080)
            WRITE(IW,9090)
         END IF
         EMC=ZERO
         E  =ZERO
         GO TO 1000
      END IF
      IF(SQCDF.GT.QUDTHR.OR.DMPING) GO TO 100
C
C     ----- SWITCH TO QUADRATIC PROCEDURE -----
C     FULL PRINTOUT ON THE FIRST SUCH ITERATION
C
      IF(MASWRK) WRITE(IW,9005)
      CALL TIMIT(1)
      DO 395 I=1,10
         IF(NPFLG(I).EQ.-5) NPFLG(I) = 0
  395 CONTINUE
      ITSAVE = ITER
C
C          FULLY QUADRATIC ITERATIONS BEGIN HERE
C
  400 CONTINUE
      ITER=ITER+1
      IF (DROPC) THEN
         CALL TRFMCX(NPFLG(2),NCORBS,NORBS,NORB,DOFOCK,DOEXCH,
     *               DDITRF,DOOOOO,DOVOOO,DOVVOO,DOVOVO,
     *               DOVVVO,DOVVVV,DOCORE)
      ELSE
         CALL TRFMCX(NPFLG(2),0,0,0,.FALSE.,.FALSE.,
     *               DDITRF,DOOOOO,DOVOOO,DOVVOO,DOVOVO,
     *               DOVVVO,DOVVVV,DOCORE)
      END IF
      NOTRF=.TRUE.
      CALL NEWTON(CISTEP,NPFLG(8),METHOD,ACURCY,DAMP,EKT,
     *            CANONC,NPUNCH,FCORE,FORS,NOCI,
     *            NWORD,NORB,NOROT,MXNORO,IROTMX,JROTMX,MOFRZ,
     *            FINAL,MCNOCI,DOEXCH,DROPC,FULLNR)
      CALL GUGSRT(NPFLG(3),DDITRF)
      CALL GUGAEM(NPFLG(3))
      IF (DROPC) NORBS=NORBS+NCORBS
      CALL GUGACI(NPFLG(4))
      CALL GUG2TM(NPFLG(6))
      CALL QUDCUP(NPFLG(8))
      CALL QUDEQS(NPFLG(8))
      CALL GUG2DM(NPFLG(6))
C
      EMC0=EMC
      EMC =ENERGY
      ETOT=ENERGY
      DE0=DE
      DE =EMC-EMC0
C
      CALL TSECND(TIM1)
      TIM0=TIM1
      IF((ITER-ITSAVE).EQ.1  .AND. MASWRK) WRITE(IW,9030)
      IF (MASWRK) WRITE(IW,9040) ITER,EMC,DE,DEMAX,SQCDF,MICRO,DAMP
      CALL FLSHBF(IW)
C
      ENERG0=ENERGY
      CVGING=DE.LT.CVGTOL.AND.DAMP.LT.DMPTOL
      IF(DMPING.AND.     CVGING                 ) DAMP=ZERO
      IF(DMPING.AND..NOT.CVGING.AND.DAMP.EQ.ZERO) DAMP=DMPTOL
C
      IF(ITER.GT.2.AND..NOT.CVGING) CALL NTNDMP(DE,DE0,DEAVG,DAMP)
      IF(.NOT.DMPING.AND.DAMP.GT.ZERO) DMPING=.TRUE.
C
C     ----- TEST FOR CONVERGENCE -----
C
      CVGED = CVGING .AND. ABS(DE).LT.TOLNG2 .AND. DEMAX.LT.TOLLAG
      IF(CVGED .AND. MASWRK) WRITE(IW,9050)
      IF(CVGED) GO TO 1000
      CVGED = CVGING .AND. ABS(DE).LT.TOLENG .AND. DEMAX.LT.TOLLG2
      IF(CVGED .AND. MASWRK) WRITE(IW,9060)
      IF(CVGED) GO TO 1000
C
C     ----- TEST FOR TIME LIMIT -----
C
      CALL TSECND(TIM1)
      TLEFT = TIMLIM-TIM1
      IF(TLEFT .LT.  1.5D+00*(TLEFTS-TLEFT)/ITER) THEN
         IF (MASWRK) THEN
            WRITE(IW,9070)
            WRITE(IW,9090)
         END IF
         EMC=ZERO
         E  =ZERO
         GO TO 1000
      END IF
C
C     ---- TEST FOR ITERATION LIMIT -----
C
      IF(ITER.GE.MAXIT) THEN
         IF (MASWRK) THEN
            WRITE(IW,9080)
            WRITE(IW,9090)
         END IF
         EMC=ZERO
         E  =ZERO
         GO TO 1000
      END IF
C
      E=ZERO
C
      IF(NPFLG(2).LE.0) NPFLG(2)=-5
      IF(NPFLG(3).LE.0) NPFLG(3)=-5
      IF(NPFLG(4).LE.0) NPFLG(4)=-5
      IF(NPFLG(6).LE.0) NPFLG(6)=-5
      IF(NPFLG(8).LE.0) NPFLG(8)=-5
      GO TO 400
C
C           **********************
C     ***** END OF MCSCF ITERATION *****
C           **********************
C
C     ----- IF -CVGED- DO ONE MORE PASS FOR CONSISTENCY -----
C
 1000 CONTINUE
      IF(CVGED .AND. EXETYP.NE.CHECK) THEN
         IF (DROPC) THEN
            CALL TRFMCX(NPFLG(2),NCORBS,NORBS,NORB,DOFOCK,DOEXCH,
     *                  DDITRF,DOOOOO,DOVOOO,DOVVOO,DOVOVO,
     *                  DOVVVO,DOVVVV,DOCORE)
         ELSE
            CALL TRFMCX(NPFLG(2),0,0,0,.FALSE.,.FALSE.,
     *                  DDITRF,DOOOOO,DOVOOO,DOVVOO,DOVOVO,
     *                  DOVVVO,DOVVVV,DOCORE)
         END IF
         CALL GUGSRT(NPFLG(3),DDITRF)
         CALL GUGAEM(NPFLG(3))
         CALL GUGADG(NPFLG(4))
         CALL GUG2DM(NPFLG(6))
         NOTRF=.TRUE.
         CALL NEWTON(CISTEP,NPFLG(8),METHOD,ACURCY,DAMP,EKT,
     *               CANONC,NPUNCH,FCORE,FORS,NOCI,
     *               NWORD,NORB,NOROT,MXNORO,IROTMX,JROTMX,MOFRZ,
     *               FINAL,MCNOCI,DOEXCH,DROPC,FULLNR)
      END IF
C
C     ----- IF AVERAGE ENERGY CALCULATION, AND USER DESIRES,
C           GET -DM1- AND -DM2- FOR STATE THAT CONTRIBUTES
C           MOST TO THE AVERAGE ENERGY
C
CMWS  IF(AVG  .AND.  UNAVG  .AND. CVGED) THEN
CMWS     AVG=.FALSE.
CMWS     CALL GUG2DM(NPFLG(6))
CMWS  END IF
C
C     ----- PRINT WAVEFUNCTION -----
C
      NPRINT = NPSAVE
      IF (MASWRK) WRITE(IW,9100) EMC,ITER
      FINAL=.TRUE.
C
      CALL NEWTON(CISTEP,NPFLG(8),METHOD,ACURCY,DAMP,EKT,
     *            CANONC,NPUNCH,FCORE,FORS,NOCI,
     *            NWORD,NORB,NOROT,MXNORO,IROTMX,JROTMX,MOFRZ,
     *            FINAL,MCNOCI,DOEXCH,DROPC,FULLNR)
C
C     ----- END OF MCSCF -----
C
      IF(CVGED) IREST=0
      IF (MASWRK) WRITE(IW,9110)
      CALL TIMIT(1)
C
C     ----- DO NATURAL ORBITAL -CI- IF REQUESTED -----
C
      IF(NOCI  .AND.  CVGED) THEN
         IF (MASWRK) WRITE(IW,9120)
         IF (DROPC) THEN
            CALL TRFMCX(0,NCORBS,NORBS,NORB,DOFOCK,DOEXCH,
     *                  DDITRF,DOOOOO,DOVOOO,DOVVOO,DOVOVO,
     *                  DOVVVO,DOVVVV,DOCORE)
         ELSE
            CALL TRFMCX(0,0,0,0,.FALSE.,.FALSE.,
     *                  DDITRF,DOOOOO,DOVOOO,DOVVOO,DOVOVO,
     *                  DOVVVO,DOVVVV,DOCORE)
         END IF
         CALL GUGSRT(0,DDITRF)
         CALL GUGAEM(0)
         CALL GUGADG(0)
         CALL GUG2DM(0)
         MCNOCI=.TRUE.
         CALL NEWTON(CISTEP,NPFLG(8),METHOD,ACURCY,DAMP,EKT,
     *               CANONC,NPUNCH,FCORE,FORS,NOCI,
     *               NWORD,NORB,NOROT,MXNORO,IROTMX,JROTMX,MOFRZ,
     *               FINAL,MCNOCI,DOEXCH,DROPC,FULLNR)
      END IF
C
      DSKWRK = TDSKWK
C
      RETURN
C
 9000 FORMAT(/10X,17(1H-)/10X,'MCSCF CALCULATION'/10X,17(1H-))
 9005 FORMAT(6X,9(1H-),'SWITCHING TO FULLY QUADRATIC PROCEDURE',9(1H-))
 9010 FORMAT(/' ----- NUCLEAR ENERGY ----- = ',F20.10)
 9028 FORMAT(/1X,'EVALUATING EXCHANGE INTEGRALS FOR DIRECT',
     *           ' INTEGRAL TRANSFORMATION')
 9030 FORMAT(/1X,'ITER',4X,'TOTAL ENERGY',10X,'DEL(E)',
     *        4X,'LAG.ASYMM.',2X,'SQCDF',2X,'MICIT',
     *        3X,'DAMP')
 9040 FORMAT(1X,I3,2F18.9,F10.6,1P,E10.3,0P,I3,F9.4)
 9050 FORMAT(/10X,20(1H-)/10X,20HLAGRANGIAN CONVERGED/10X,20(1H-))
 9060 FORMAT(/10X,16(1H-)/10X,16HENERGY CONVERGED    /10X,16(1H-))
 9070 FORMAT(1X,'RUNNING OUT OF CPU TIME...')
 9080 FORMAT(1X,'EXCESSIVE NUMBER OF ITERATIONS...')
 9090 FORMAT(1X,'MCSCF IS NOT CONVERGED!'/
     *       1X,'A $VEC GROUP OF CURRENT MO-S IS IN THE PUNCH FILE'/
     *       1X,'USE THIS WITH GUESS=MOREAD TO RESTART THIS RUN')
 9100 FORMAT(/1X,'FINAL MCSCF ENERGY IS',F20.10,' AFTER',I4
     *       ,' ITERATIONS')
 9110 FORMAT(1X,'.....DONE WITH MCSCF ITERATIONS.....')
 9120 FORMAT(//10X,34(1H-)/10X,'CI BASED ON MCSCF NATURAL ORBITALS'/
     *         10X,34(1H-))
      END
C*MODULE MCQUD   *DECK GUG2TM
      SUBROUTINE GUG2TM(NPRINT)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      PARAMETER (MXATM=2000, MXRT=100, NEXTMX=200, MXAO=8192)
      PARAMETER (NNAM=7)
C
      LOGICAL AVG
      LOGICAL OUT,SOME,GOPARR,DSKWRK,MASWRK
C
      COMMON /CIFILS/ NFT11,NFT12,NFT13,NFT14,NFT15,NFT16,IDAF20,NEMEMX
      COMMON /DM2DM2/ IH00,IH01,NWKS1
      COMMON /EMPRT / ILEVEL,JLEVEL,KLEVEL,LLEVEL,OUT,SOME
      COMMON /ENRGYS/ EN,EELCT,ETOT,SZ,SZZ,EC,ESCF,EERD,E1,E2,
     *                VEN,VEE,EPOT,EKIN,ECI(MXRT),STATN,EDFT(2)
      COMMON /FMCOM / H(1)
      COMMON /FUNCT / E,EG(3,MXATM)
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
      COMMON /LOOPS / X(5),IX(12),IJKLAD,NINTMX,NGRPS
      COMMON /ORBSET/ NORBMX,NORBS,NCORBS,NLEVS,NNA,NNB,NNC,NSYM,MSYM,
     *                IDOCC,IVAL,IMCC,ISYM(MXAO),JCODE(MXAO),
     *                NLCS(MXAO),LEVPT(MXAO),LEVNR(MXAO),IOUT(MXAO),
     *                NREFS,IEXCT,NFOCI,INTACT
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
      COMMON /RUNLAB/ TITLE(10),A(MXATM),B(MXATM),BFLAB(MXAO)
      COMMON /RUNOPT/ RUNTYP,EXETYP,NEVALS,NGLEVL,NHLEVL
      COMMON /TM2AVG/ WAVG(MXRT),IAVG(MXRT),NAVG,AVG
      COMMON /TM2BUG/ ICI
C
      DIMENSION WSTATE(MXRT)
      DIMENSION TITLE1(10)
      DIMENSION QNAM(NNAM),KQNAM(NNAM)
C
      PARAMETER (ZERO=0.0D+00, ONE=1.0D+00)
C
      DATA GUGDM2/8HGUGDM2  /
      DATA QNAM/8HCUTOFF  ,8HNWORD   ,8HWSTATE  ,8HNDAR    ,
     *          8HLDAR    ,8HNBOXMX  ,8HNOMEM   /
      DATA KQNAM/3,1,-3,1,1,1,1/
      DATA CHECK/8HCHECK   /
C
      SOME=NPRINT.NE.-5 .AND. MASWRK
      OUT =NPRINT.GE.1 .AND. MASWRK
      IF(SOME) WRITE(IW,9999)
C
C     ----- READ NAMELIST -$GUGDM2- -----
C
      CUTOFF=1.0D-09
      NWORD=0
      KQNAM(3) = 10*MXRT + 3
      CALL VCLR(WSTATE,1,MXRT)
      WSTATE(1) = ONE
      NDAR=4000
      NBOXMX=200
      CALL RASIZE(LDAR)
      NOMEM=0
C
      CALL NAMEIO(IR,JRET,GUGDM2,NNAM,QNAM,KQNAM,
     *            CUTOFF,NWORD,WSTATE,NDAR,LDAR,NBOXMX,NOMEM,
     *            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,  0,0,0,0,0)
      IF(JRET.EQ.2) THEN
         IF (MASWRK) WRITE (IW,9010)
         CALL ABRT
      END IF
C
C     ----- GET CORE MEMORY -----
C
      CALL GOTFM(NGOT)
      IF(NWORD.GT.0) NGOT=MIN(NWORD,NGOT)
C
C     ----- READ CI STATE WAVEFUNCTION -----
C
      CALL SEQREW(NFT12)
      READ(NFT12) NSTATE,NWKS,TITLE,TITLE1
      NCI = NSTATE
C
C     ----- NORMALIZE THE STATE WEIGHTS -----
C
      NAVG = 0
      WSUM = ZERO
      DO 20 I=1,MXRT
         IF(WSTATE(I).LT.ZERO)THEN
           IF(MASWRK) WRITE(IW,*)' STATE ',I,' HAS WEIGHT ',WSTATE(I)
           CALL ABRT
         END IF
         WSUM = WSUM + WSTATE(I)
         IF(WSTATE(I).EQ.ZERO) GO TO 20
            NAVG = NAVG+1
            IAVG(NAVG) = I
            WAVG(NAVG) = WSTATE(I)
   20 CONTINUE
C
      DO 30 I=1,NAVG
         WAVG(I) = WAVG(I)/WSUM
   30 CONTINUE
C
      IF(NAVG.NE.1) GO TO 9040
C
C     ----- READ CI WAVEFUNCTION DATA -----
C
      CALL VALFM(LOADFM)
      IH00  = LOADFM + 1
      IH01  = IH00   + NWKS*NAVG
      NEED1 = IH01   - LOADFM - 1
      CALL GETFM(NEED1)
C
      JOFF = IH00
      IST = 1
      DO 50 I=1,NSTATE
         IF(I.EQ.IAVG(IST)) THEN
            CALL SQREAD(NFT12,H(JOFF),NWKS)
            JOFF = JOFF + NWKS
            IST = IST+1
         ELSE
            READ (NFT12)
         END IF
   50 CONTINUE
C
C     ----- CALCULATE AVERAGE STATE ENERGY -----
C
      E=ZERO
      DO 100 I=1,NAVG
         ICI=IAVG(I)
         IF(ICI.GT.NCI) GO TO 9020
         E=E+ECI(ICI)*WAVG(I)
  100 CONTINUE
C
C     ----- READ IN -DRT- DATA -----
C
      CALL SEQREW(NFT11)
      NORBMX = 0
      NROWS = 0
      READ (NFT11) NORBMX,NORBS,NSYM,NROWS,NWKS,LEVFRM,NEMEMX,
     *             NREFS,IEXCT,NFOCI,INTACT,NCORBS
      READ (NFT11) TITLE1
C
      IH1=IH01+NWKS
      IH2=IH1+4*NROWS
      IH3=IH2+NROWS
      IH4=IH3+NROWS
      IH5=IH1
      IH6=IH4+NROWS
      IH7=IH6+NROWS
      IH8=IH7+NROWS
      IH11=IH8+NROWS
      IH12=IH11+NROWS*4
      NKL=NSYM*NORBS
      NIJ=(NORBS*(NORBS+1))/2
      NLNG=LEVFRM*NSYM
      IH13=IH12+NIJ
      IH14=IH13+NIJ
      IH15=IH14+NKL
      IH16=IH15+NKL
      IH17=IH16+NEXTMX
      IH18=IH17+NEXTMX
      IH19=IH18+NEXTMX
      IH20=IH19+NLNG
      IH21=IH20+NLNG
      IH22=IH21+NLNG
      IH23=IH22+NLNG
      IH26=IH23+NLNG
      NEED2 = IH26 - IH01
      CALL GETFM(NEED2)
      NEED = NEED1 + NEED2
C
      IJKLAD=IH26
      NROWS4=NROWS*4
      CALL DRTTAP(NFT11,NROWS,NROWS4,NKL,NIJ,H(IH2),H(IH3),H(IH4),
     *            H(IH5),H(IH6),H(IH7),H(IH8),H(IH01),H(IH11),
     *            H(IH12),H(IH13),H(IH14),H(IH15),H(IH16),H(IH17),
     *            H(IH18),NWKS,NEXT,NGRPS,MXNINT)
C
C        LOOK AHEAD AT MEMORY -NEED3- TO BE ALLOCATED IN TM2GUG
C
      NEED3=NINTMX*NWKS
      IF(EXETYP.EQ.CHECK) THEN
         MXNEED = NEED1 + NEED2 + NEED3
         IF(MASWRK) WRITE(IW,9950) MXNEED,NGOT
      END IF
      LEFT=NGOT-NEED
      IF(NEED3.GT.LEFT) GO TO 9000
C
C     ----- WRITE HEADER RECORD ON -TM2- FILE -----
C
      CALL SEQREW(NFT15)
      ITAP=3
      WRITE(NFT15) ITAP,NWKS,NEMEMX,CUTOFF,TITLE,TITLE1
      IF(SOME) WRITE(IW,9992) TITLE,TITLE1,NWKS
      DO 110 I=1,NAVG
         IF(SOME) WRITE(IW,9991) IAVG(I),WAVG(I)
  110 CONTINUE
      IF(SOME) WRITE(IW,9984)
C
C     ----- CALCULATE -CI- 2-PARTICLE TRANSITION MATRIX -----
C
      CALL TM2GUG(NWKS)
C
      IF(OUT) WRITE(IW,9995)
      IF (EXETYP .EQ. CHECK) THEN
         CALL TM2FIN(NWKS)
      ELSE
         CALL LOOPY2(NROWS,H( IH3),H( IH6),H( IH7),H( IH8),H(IH11),
     *               H( IH2),H(IH12),H(IH13),H(IH14),H(IH15),H(IH16),
     *               H(IH17),H(IH18),H( IH5),NEXT,NROWS4,NIJ,NKL,
     *               H(IH19),H(IH20),H(IH21),H(IH22),H(IH23),H( IH4),
     *               NLNG,LEVFRM,NWKS)
      END IF
C
C     ----- RESET CORE MEMORY -----
C
      CALL RETFM(NEED2)
      CALL RETFM(NEED1)
C
C     ----- SORT -TM2- ELEMENTS -----
C
      CALL SRT2TM(NWKS,NPRINT,NWORD,NORBMX,NORBS,NSYM,
     *            NGRPS,MXNINT,NFT11,NFT15,NFT16)
C
      IF(SOME) WRITE(IW,9997)
      IF(SOME) CALL TEXIT(2,2)
      RETURN
 9000 CONTINUE
      IF (MASWRK) WRITE(IW,9994) NEED3,LEFT
      CALL ABRT
      RETURN
 9020 CONTINUE
      IF (MASWRK) WRITE(IW,9989) ICI,NCI
      CALL ABRT
      RETURN
 9040 CONTINUE
      IF (MASWRK) WRITE(IW,9986)
      CALL ABRT
      RETURN
 9999 FORMAT(/,10X,40(1H-),
     1       /,10X,'2-PARTICLE TRANSITION MATRIX CALCULATION',
     2       /,10X,40(1H-))
 9997 FORMAT(' ...... END OF 2-PARTICLE TRANSITION MATRIX ',
     1       'CALCULATION ...... ')
 9995 FORMAT(' CALCULATE 2-PARTICLE TRANSITION MATRIX ')
 9994 FORMAT(' *** ERROR *** MORE MEMORY IS NEEDED IN GUG2TM.'/
     *       ' NEED MEMORY = NXNINT*NO.CSFS = ',I10,
     *       ' REMAINING MEMORY = ',I10/
     *       ' SEE MXNINT IN $DRT, BUT DO NOT ADJUST MXNINT SO LOW'/
     *       ' THAT THE NUMBER OF INTEGRALS GROUPS EXCEEDS 8-10.'/)
 9992 FORMAT(/,' NAME OF ENERGY MATRIX TAPE ... ',10A8,
     1       /,' NAME OF THE -DRT- TABLE    ... ',10A8,
     2       /,' NUMBER OF CONFIGURATIONS = ',I7)
 9991 FORMAT(' STATE NUMBER = ',I5,' WEIGHT = ',F8.5)
 9989 FORMAT(' INCONSISTENCY IN -GUGTM2- .    ICI,NCI = ',2I5)
 9986 FORMAT(' QUADRATIC MCHF PROGRAM CAN NOT HANDLE',
     1       ' STATE AVERAGING. STOP ')
 9984 FORMAT(' -GUGA- LOOPS ARE GENERATED HERE TOO.')
 9010 FORMAT(1X,'**** ERROR READING $GUGDM2.')
 9950 FORMAT(/,' STATISTICS FOR FORMING TRANSITION DENSITY MATRIX:'/
     1         ' NUMBER OF WORDS NEEDED    = ',I10,
     2         ' NUMBER OF WORDS AVAILABLE = ',I10)
      END
C*MODULE MCQUD   *DECK TM2GUG
      SUBROUTINE TM2GUG(NWKS)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      PARAMETER (MXRT=100)
C
      LOGICAL OUT,SOME,AVG
C
      COMMON /CIFILS/ NFT11,NFT12,NFT13,NFT14,NFT15,NFT16,IDAF20,NEMEMX
      COMMON /DM2DM2/ IH00,IH01,NWKS1
      COMMON /EMPRT / ILEVEL,JLEVEL,KLEVEL,LLEVEL,OUT,SOME
      COMMON /FMCOM / GIJKL(1)
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
      COMMON /LOOPS / A,C,D,VAL,CUTOFF,LAD,NUWK,NLWK,IUWK,JUWK,ITRACK,
     *                IX(6),IJKLAD,NINTMX,NGRPS
      COMMON /MACHIN/ ISINGL,MAXFM,MAXSM,LIMFM,LIMSM
      COMMON /SDM2SV/ SQRT2,LGRP,NLEXT,NLINT,NRECO
      COMMON /TM2AVG/ WAVG(MXRT),IAVG(MXRT),NAVG,AVG
C
      DIMENSION CIVEC(1),INDX(1)
      EQUIVALENCE (CIVEC(1),GIJKL(1)),(INDX(1),GIJKL(1))
C
      PARAMETER (PT5=0.5D+00, TWO=2.0D+00)
C
      SQRT2= SQRT(TWO)
      LGRP=1
      NLEXT=0
      NLINT=0
C
C           GET GIJKL BUFFER SPACE
C
      CALL VALFM(LOADFM)
      IJKLAD = LOADFM + 1
      LAST   = IJKLAD + NINTMX*NWKS
      NEED   = LAST   - LOADFM - 1
      CALL GETFM(NEED)
      GO TO 20
C
C      ------------------
       ENTRY TM2GRP(NWKS)
C      ------------------
       IF(LGRP.GT.NGRPS) THEN
          WRITE(IW,*) 'CONFUSION IN ENTRY POINT TM2GRP'
          CALL ABRT
       END IF
       IADDRS=IJKLAD
       DO 10 IWKS=1,NWKS
       CALL SQWRIT(NFT15,GIJKL(IADDRS),NINTMX)
   10  IADDRS=IADDRS+NINTMX
       LGRP=LGRP+1
C
C     ----- INITIALIZE TM2 ELEMENTS IN -LGRP- -----
C
   20  CONTINUE
       CALL VCLR(GIJKL(IJKLAD),1,NINTMX*NWKS)
       RETURN
C
C     ------------
      ENTRY TM2EXT
C     ------------
      NLEXT=NLEXT+1
      IF(NUWK.EQ.0) RETURN
      GO TO 200
C
C     ------------------
      ENTRY TM2OUT(NWKS)
C     ------------------
      NLINT=NLINT+1
C
C     ----- ACCOUNT FOR POSITION OF -INDX- ARRAY -----
C
  200 CONTINUE
      INDOFF=(IH01-1) * ISINGL
      IND=IUWK+INDOFF
      JND=JUWK+INDOFF
      ICIOFF=NWKS*(NAVG-1)
      DO 320 I=1,NLWK
      II=INDX(IND)+ICIOFF
      JJ=INDX(JND)+ICIOFF
      DO 310 J=1,NUWK
      VALI=CIVEC(II + IH00 -1)
      VALJ=CIVEC(JJ + IH00 -1)
C
      NAD=LAD-(IJKLAD-1)
      LAD1= NAD*NWKS + IJKLAD -1
      LAD2= LAD1 + NWKS
      LAD3= LAD2 + NWKS
C
      IF(II.EQ.JJ) VALI=VALI*PT5
      IF(II.EQ.JJ) VALJ=VALJ*PT5
C
C     ----- CALCULATE TRANSITION DENSITY MATRIX ELEMENTS -----
C
      GO TO (261,262,263,264,265,266,267,268,269,270,
     1       271,272,273,274,275,276,277,278),ITRACK
  261 GIJKL(II +LAD1)=GIJKL(II +LAD1) + VALJ * A
      GIJKL(JJ +LAD1)=GIJKL(JJ +LAD1) + VALI * A
      GO TO 300
  262 GIJKL(II +LAD2)=GIJKL(II +LAD2) + VALJ * A
      GIJKL(JJ +LAD2)=GIJKL(JJ +LAD2) + VALI * A
      GO TO 300
  263 GIJKL(II +LAD3)=GIJKL(II +LAD3) + VALJ * A
      GIJKL(JJ +LAD3)=GIJKL(JJ +LAD3) + VALI * A
      GO TO 300
  264 GIJKL(II +LAD1)=GIJKL(II +LAD1) + VALJ
      GIJKL(II +LAD2)=GIJKL(II +LAD2) + VALJ
      GIJKL(JJ +LAD1)=GIJKL(JJ +LAD1) + VALI
      GIJKL(JJ +LAD2)=GIJKL(JJ +LAD2) + VALI
      GO TO 300
  265 GIJKL(II +LAD1)=GIJKL(II +LAD1) + VALJ
      GIJKL(II +LAD2)=GIJKL(II +LAD2) +(VALJ + VALJ)
      GIJKL(JJ +LAD1)=GIJKL(JJ +LAD1) + VALI
      GIJKL(JJ +LAD2)=GIJKL(JJ +LAD2) +(VALI + VALI)
      GO TO 300
  266 GIJKL(II +LAD2)=GIJKL(II +LAD2) + VALJ
      GIJKL(II +LAD3)=GIJKL(II +LAD3) + VALJ
      GIJKL(JJ +LAD2)=GIJKL(JJ +LAD2) + VALI
      GIJKL(JJ +LAD3)=GIJKL(JJ +LAD3) + VALI
      GO TO 300
  267 GIJKL(II +LAD1)=GIJKL(II +LAD1) +(VALJ * A)
      GIJKL(II +LAD3)=GIJKL(II +LAD3) +(VALJ * A)
      GIJKL(JJ +LAD1)=GIJKL(JJ +LAD1) +(VALI * A)
      GIJKL(JJ +LAD3)=GIJKL(JJ +LAD3) +(VALI * A)
      GO TO 300
  268 GIJKL(II +LAD1)=GIJKL(II +LAD1) +(VALJ * A)
      GIJKL(II +LAD2)=GIJKL(II +LAD2) +(VALJ * A)
      GIJKL(II +LAD3)=GIJKL(II +LAD3) +(VALJ * A)
      GIJKL(JJ +LAD1)=GIJKL(JJ +LAD1) +(VALI * A)
      GIJKL(JJ +LAD2)=GIJKL(JJ +LAD2) +(VALI * A)
      GIJKL(JJ +LAD3)=GIJKL(JJ +LAD3) +(VALI * A)
      GO TO 300
  269 GIJKL(II +LAD1)=GIJKL(II +LAD1) +(VALJ * A)
      GIJKL(II +LAD2)=GIJKL(II +LAD2) +(VALJ * A)* D
      GIJKL(JJ +LAD1)=GIJKL(JJ +LAD1) +(VALI * A)
      GIJKL(JJ +LAD2)=GIJKL(JJ +LAD2) +(VALI * A)* D
      GO TO 300
  270 GIJKL(II +LAD1)=GIJKL(II +LAD1) +(VALJ * A)* D
      GIJKL(II +LAD2)=GIJKL(II +LAD2) +(VALJ * A)
      GIJKL(JJ +LAD1)=GIJKL(JJ +LAD1) +(VALI * A)* D
      GIJKL(JJ +LAD2)=GIJKL(JJ +LAD2) +(VALI * A)
      GO TO 300
  271 GIJKL(II +LAD1)=GIJKL(II +LAD1) +(VALJ * A)
      GIJKL(II +LAD3)=GIJKL(II +LAD3) +(VALJ * A)* D
      GIJKL(JJ +LAD1)=GIJKL(JJ +LAD1) +(VALI * A)
      GIJKL(JJ +LAD3)=GIJKL(JJ +LAD3) +(VALI * A)* D
      GO TO 300
  272 GIJKL(II +LAD2)=GIJKL(II +LAD2) +(VALJ * A)* D
      GIJKL(II +LAD3)=GIJKL(II +LAD3) +(VALJ * A)
      GIJKL(JJ +LAD2)=GIJKL(JJ +LAD2) +(VALI * A)* D
      GIJKL(JJ +LAD3)=GIJKL(JJ +LAD3) +(VALI * A)
      GO TO 300
  273 GIJKL(II +LAD2)=GIJKL(II +LAD2) +(VALJ * A)
      GIJKL(II +LAD3)=GIJKL(II +LAD3) +(VALJ * A)
      GIJKL(JJ +LAD2)=GIJKL(JJ +LAD2) +(VALI * A)
      GIJKL(JJ +LAD3)=GIJKL(JJ +LAD3) +(VALI * A)
      GO TO 300
  274 GIJKL(II +LAD1)=GIJKL(II +LAD1) + VALJ
      GIJKL(II +LAD3)=GIJKL(II +LAD3) + VALJ
      GIJKL(JJ +LAD1)=GIJKL(JJ +LAD1) + VALI
      GIJKL(JJ +LAD3)=GIJKL(JJ +LAD3) + VALI
      GO TO 300
  275 GIJKL(II +LAD2)=GIJKL(II +LAD2) + VALJ
      GIJKL(II +LAD3)=GIJKL(II +LAD3) - VALJ
      GIJKL(JJ +LAD2)=GIJKL(JJ +LAD2) + VALI
      GIJKL(JJ +LAD3)=GIJKL(JJ +LAD3) - VALI
      GO TO 300
  276 GIJKL(II +LAD1)=GIJKL(II +LAD1) + VALJ * SQRT2
      GIJKL(JJ +LAD1)=GIJKL(JJ +LAD1) + VALI * SQRT2
      GO TO 300
  277 GIJKL(II +LAD1)=GIJKL(II +LAD1) +(VALJ * A)
      GIJKL(II +LAD2)=GIJKL(II +LAD2) +(VALJ * A)
      GIJKL(JJ +LAD1)=GIJKL(JJ +LAD1) +(VALI * A)
      GIJKL(JJ +LAD2)=GIJKL(JJ +LAD2) +(VALI * A)
      GO TO 300
  278 GIJKL(II +LAD2)=GIJKL(II +LAD2) +(VALJ * A)
      GIJKL(II +LAD3)=GIJKL(II +LAD3) -(VALJ * A)
      GIJKL(JJ +LAD2)=GIJKL(JJ +LAD2) +(VALI * A)
      GIJKL(JJ +LAD3)=GIJKL(JJ +LAD3) -(VALI * A)
  300 CONTINUE
       IF(OUT) WRITE(IW,9999) NLINT,NLEXT,ILEVEL,JLEVEL,KLEVEL,LLEVEL,
     * NAD,NUWK,NLWK,IUWK,JUWK,ITRACK,A,D,VALJ,VALI,
     * II ,GIJKL(II +LAD1),GIJKL(II +LAD2),GIJKL(II +LAD3),
     * JJ ,GIJKL(JJ +LAD1),GIJKL(JJ +LAD2),GIJKL(JJ +LAD3)
C
      II=II+1
  310 JJ=JJ+1
      IND=IND+1
  320 JND=JND+1
      RETURN
C
C     ------------------
      ENTRY TM2FIN(NWKS)
C     ------------------
      CALL SQWRIT(NFT15,GIJKL(IJKLAD),NINTMX*NWKS)
C
      NLTOT=NLINT+NLEXT
C
      IF(OUT) THEN
         WRITE(IW,404) NLTOT
         WRITE(IW,405) NLINT
         WRITE(IW,406) NLEXT
      END IF
C
      NEED = NINTMX*NWKS
      CALL RETFM(NEED)
      RETURN
C
  404 FORMAT(/,I12,' IS THE TOTAL NUMBER OF GENERATED LOOPS')
  405 FORMAT(I12,' WERE CREATED BY THE LOOP-DRIVEN ALGORITHM')
  406 FORMAT(I12,' WERE CREATED IMPLICITLY')
 9999 FORMAT(1X,2I10,4I4,5I6,I3,2F12.6,2E16.8,
     * /,I10,3E16.8,I10,3E16.8)
      END
C*MODULE MCQUD   *DECK SRT2TM
      SUBROUTINE SRT2TM(NWKS,NPRINT,NWORD,NORBMX,NORBS,NSYM,
     *                  NGRPS,MXNINT,NFT11,NFT15,NFT16)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      LOGICAL OUT
      LOGICAL GOPARR,DSKWRK,MASWRK
C
      PARAMETER (MXAO=8192)
C
      COMMON /FMCOM / X(1)
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
C
      DIMENSION TITLE(10),ISYM(MXAO),MAP(MXAO)
C
      OUT =NPRINT.GE.3 .AND. MASWRK
      IF(OUT) WRITE(IW,9999)
C
C     ----- GET CORE MEMORY -----
C
      CALL GOTFM(NGOT)
      IF(NWORD.GT.0) NGOT=MIN(NWORD,NGOT)
C
C     ----- READ IN -DRT- DATA -----
C
      CALL SEQREW(NFT11)
      READ (NFT11)
      READ (NFT11) TITLE
      READ (NFT11) MAP
      READ (NFT11) ISYM
      DO 10 I=1,12
   10 READ (NFT11)
      M1= NORBS
      M2=(M1*(M1+1))/2
      MN= M1*NSYM
C
      CALL VALFM(LOADFM)
      IJADD = LOADFM + 1
      IJGRP = IJADD + M2
      KADD  = IJGRP + M2
      LADD  = KADD  + MN
      LAST  = LADD  + MN
      NEED1 = LAST - LOADFM - 1
      CALL GETFM(NEED1)
C
      CALL SRTDRT(NFT11,X(IJADD),X(IJGRP),X(KADD),X(LADD),M2,MN)
      IF(OUT) WRITE(IW,9991) TITLE
C
C     ----- SORT THE -TM2- ELEMENTS -----
C
      I20   = LAST
      LAST  = I20  + NWKS*MXNINT
      NEED2 = LAST - I20
C
      NEED = NEED1 + NEED2
      IF(NEED.GT.NGOT) THEN
         MORE=NEED-NGOT
         IF (MASWRK) WRITE(IW,9995) NGOT,MORE
      END IF
      CALL GETFM(NEED2)
C
      NFTI=NFT15
      NFTO=NFT16
      CALL GETTM2(X(I20),M1,NFTI,NFTO,
     *            X(IJADD),X(IJGRP),X(KADD),X(LADD),
     *            ISYM,MAP,NSYM,NORBMX,NORBS,NGRPS,NWKS,MXNINT,OUT)
      CALL SEQREW(NFT15)
      CALL SEQREW(NFT16)
      IF(OUT) WRITE(IW,9992)
C
C     ----- RESET CORE MEMORY -----
C
      CALL RETFM(NEED2)
      CALL RETFM(NEED1)
      RETURN
C
 9999 FORMAT(/,10X,37(1H-),/,10X,'2-PARTICLE TRANSITION MATRIX SORTING',
     1       /,10X,37(1H-))
 9995 FORMAT(' SRT2TM: NOT ENOUGH CORE AVAILABLE. LEFT = ',I10,
     1       ' NEED ',I10,' MORE.')
 9992 FORMAT(' ...... END OF TRANSITION MATRIX SORTING ...... ')
 9991 FORMAT(/,' NAME OF THE DISTINCT ROW TABLE @ ',10A8)
      END
C*MODULE MCQUD   *DECK SRTDRT
      SUBROUTINE SRTDRT(NFT11,IJADD,IJGRP,KADD,LADD,M2,MN)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      DIMENSION IJADD(M2),IJGRP(M2),KADD(MN),LADD(MN)
C
      READ(NFT11) IJADD
      READ(NFT11) IJGRP
      READ(NFT11) KADD
      READ(NFT11) LADD
C
      RETURN
      END
C*MODULE MCQUD   *DECK GETTM2
      SUBROUTINE GETTM2(GIJKL,L1,NFTI,NFTO,
     * IJADD,IJGRP,KADD,LADD,ISYM,MAP,NSYM,NORBMX,NORBS,NGRPS,NWKS,
     * MXNINT,OUT)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      LOGICAL OUT
      PARAMETER (MXAO=8192)
C
      COMMON /IJPAIR/ IA(MXAO)
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
C
      DIMENSION ISYM(1),IJADD(1),IJGRP(1),KADD(1),LADD(1)
      DIMENSION MAP(1)
      DIMENSION GIJKL(NWKS,MXNINT)
      DIMENSION MULT8(8),LKUPSM(64),KLOFF(8)
C
      DATA PT5 /0.5D+00/
      DATA MULT8  /0,8,16,24,32,40,48,56/
      DATA LKUPSM /1,2,3,4,5,6,7,8,
     1             2,1,4,3,6,5,8,7,
     2             3,4,1,2,7,8,5,6,
     3             4,3,2,1,8,7,6,5,
     4             5,6,7,8,1,2,3,4,
     5             6,5,8,7,2,1,4,3,
     6             7,8,5,6,3,4,1,2,
     7             8,7,6,5,4,3,2,1/
C
C     ----- SORT -TM2- ELEMENTS -----
C
      NFZC=0
      DO 10 I=1,NORBMX
      IF(MAP(I).LT.0) NFZC=NFZC+1
   10 CONTINUE
      DO 20 I=1,NSYM
   20 KLOFF(I)=(I-1)*NORBS
C
C     ----- SKIP HEADER RECORD ON -TM2- FILE -----
C
      CALL SEQREW(NFTO)
C
C     ----- START NEW LOOPS -----
C
      CALL SEQREW(NFTI)
      READ(NFTI)
      IGRP=0
  100 IGRP=IGRP+1
      IF(IGRP.GT.NGRPS) GO TO 400
C
      LEN = MXNINT*NWKS
      CALL SQREAD(NFTI,GIJKL,LEN)
C
C     ----- PROCESS 1E-DM FIRST -----
C
      DO 130 K0=1,L1
      DO 120 L0=1,K0
      KMAP=MAP(K0+NFZC)
      LMAP=MAP(L0+NFZC)
      K=MAX0(KMAP,LMAP)
      L=MIN0(KMAP,LMAP)
      MB=3
      IF(K.EQ.L) MB=2
      KL=IA(K)+K
      KKSYM=ISYM(K)
      IF(KKSYM.NE.ISYM(L)) GO TO 120
      LSMPT=L + KLOFF(KKSYM)
      LAD =IJADD(KL)+KADD(K)+LADD(LSMPT)+MB
      LGRP=IJGRP(KL)
      IF(LGRP.NE.IGRP) GO TO 120
      DO 110 IWKS = 1,NWKS
      VAL=GIJKL(IWKS,LAD)
      IF(L0.NE.K0) VAL=VAL*PT5
      GIJKL(IWKS,LAD)=VAL
      IF(OUT) WRITE(IW,9998) K0,L0,K,L,IWKS,LAD,LGRP,VAL
  110 CONTINUE
      CALL QUDWTM(NFTO,1,0,0,K0,L0,GIJKL(1,LAD),NWKS)
  120 CONTINUE
  130 CONTINUE
C
C     ----- PROCESS 2E-DM NOW -----
C
      DO 340 I0=1,L1
      DO 330 J0=1,I0
      DO 320 K0=1,I0
      L0MAX=K0
      IF(K0.EQ.I0) L0MAX=J0
      DO 310 L0=1,L0MAX
C
      IMAP=MAP(I0+NFZC)
      JMAP=MAP(J0+NFZC)
      KMAP=MAP(K0+NFZC)
      LMAP=MAP(L0+NFZC)
      I=MAX0(IMAP,JMAP)
      J=MIN0(IMAP,JMAP)
      K=MAX0(KMAP,LMAP)
      L=MIN0(KMAP,LMAP)
      IF(I-K) 210,220,230
  210 N=K
      K=I
      I=N
      N=L
      L=J
      J=N
      GO TO 230
  220 IF(J-L) 210,230,230
  230 CONTINUE
      IF(J-K) 240,250,260
  240 II=I
      JJ=K
      KK=J
      LL=L
      IF(KK.GE.LL) GO TO 280
      KK=L
      LL=J
      IF(JJ.GT.KK) GO TO 270
      MB=2
      GO TO 290
  250 II=I
      JJ=J
      KK=K
      LL=L
      IF(JJ.GT.LL) GO TO 280
      IF(II.EQ.JJ) GO TO 280
      JJ=I
      KK=I
      MB=2
      GO TO 290
  260 II=I
      JJ=J
      KK=K
      LL=L
      MB=2
      GO TO 290
  270 MB=3
      GO TO 290
  280 MB=1
  290 CONTINUE
      IJ=IA(II)+JJ
      IISYM=ISYM(II)
      IJMUL=MULT8(IISYM)+ISYM(JJ)
      IJSYM=LKUPSM(IJMUL)
      IJKMUL=MULT8(IJSYM)+ISYM(KK)
      IJKSYM=LKUPSM(IJKMUL)
      LSYM=ISYM(LL)
      IF(LSYM.NE.IJKSYM) GO TO 310
      KSMPT=KK+KLOFF(IJSYM)
      LSMPT=LL+KLOFF(IJKSYM)
      LAD =IJADD(IJ)+KADD(KSMPT)+LADD(LSMPT)+MB
      LGRP=IJGRP(IJ)
      IF(LGRP.NE.IGRP) GO TO 310
C
      DO 300 IWKS=1,NWKS
      VAL=GIJKL(IWKS,LAD)
      VAL=VAL+VAL
      IF(J0.NE.I0            ) VAL=VAL*PT5
      IF(L0.NE.K0            ) VAL=VAL*PT5
      IF(K0.NE.I0.OR.L0.NE.J0) VAL=VAL*PT5
      GIJKL(IWKS,LAD)=VAL
      IF(OUT) WRITE(IW,9999) I0,J0,K0,L0,I,J,K,L,II,JJ,KK,LL,IWKS,
     *                       LAD,LGRP,VAL
  300 CONTINUE
C
      CALL QUDWTM(NFTO,2,I0,J0,K0,L0,GIJKL(1,LAD),NWKS)
  310 CONTINUE
  320 CONTINUE
  330 CONTINUE
  340 CONTINUE
      GO TO 100
C
C     ----- FINISH
C
  400 CONTINUE
      RETURN
 9999 FORMAT(12I5,2I10,I5,E20.12)
 9998 FORMAT(4I5,2I10,I5,E20.12)
      END
C*MODULE MCQUD   *DECK LOOPY2
      SUBROUTINE LOOPY2(NROWS,NABC,NLWKS,NUWKS,PUWK,IWGHT,NELEC,
     *                  IJADD,IJGRP,KADD,LADD,INEXT,JMNNXT,JMXNXT,
     *                  IARC,NEXT,NROWS4,NIJ,NKL,NUWKY,NLWKX,PUWKY,
     *                  NLWKW,NLWKY,NABCS,NLNG,LEVFRM,NWKS)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      INTEGER PUWK,PUWKY
      LOGICAL OUT,SOME,GOPARR,DSKWRK,MASWRK
C
      PARAMETER (MXAO=8192)
C
      COMMON /EMPRT / ILEVEL,JLEVEL,KLEVEL,LLEVEL,OUT,SOME
      COMMON /IOFILE/ IR,IW,IP,IJK,IPK,IDAF,NAV,IODA(950)
      COMMON /LOOPS / ACF,CCF,D,VAL,CUTOFF,LAD,NUWK,NLWK,IUWK,JUWK,
     *                ITRAK,LEVI,JMIN,JMAX,IAD,IISYMT,KSBT,INTSAD,
     *                NINTMX,NGRPS
      COMMON /ORBSET/ NORBMX,NORBS,NCORBS,NLEVS,NNA,NNB,NNC,NSYM,MSYM,
     *                IDOCC,IVAL,IMCC,ISYM(MXAO),JCODE(MXAO),
     *                NLCS(MXAO),LEVPT(MXAO),LEVNR(MXAO),IOUT(MXAO),
     *                NREFS,IEXCT,NFOCI,INTACT
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
C
      DIMENSION ISEGM(128),JSEGM(128),IMAIN(128),ISUB(128)
      DIMENSION IUWKMN(128),IUWKSB(128),ITRACK(128),ISHIFT(4)
      DIMENSION ACOEF(128),BCOEF(128)
      DIMENSION ISMOFF(8)
      DIMENSION LMIN(128)
      DIMENSION NABC(NROWS),NLWKS(NROWS),NUWKS(NROWS),IARC(NROWS4),
     * PUWK(NROWS),IWGHT(NROWS4),IJADD(NIJ),IJGRP(NIJ),KADD(NKL),
     * LADD(NKL),INEXT(1),JMNNXT(1),JMXNXT(1)
      DIMENSION NELEC(NROWS),NLCSMN(22)
      DIMENSION LCOND(8)
      DIMENSION JSEGNR(22),JSEGPT(22),IARCMN(228),IARCSB(228)
      DIMENSION ITRK(228),JCOND(228),KCOND(228),NXTSEG(228)
      DIMENSION COEFFS(20,21),CFS(420)
      DIMENSION JSEGPX(3)
      DIMENSION MULTS(8),LKUPSM(64)
      DIMENSION NUWKY(1),NLWKX(1),PUWKY(1),NLWKW(1),NLWKY(1),NABCS(1)
C
      EQUIVALENCE (COEFFS(1,1),CFS(1))
C
      DATA MULTS  /0,8,16,24,32,40,48,56/
      DATA LKUPSM /1,2,3,4,5,6,7,8,
     1             2,1,4,3,6,5,8,7,
     2             3,4,1,2,7,8,5,6,
     3             4,3,2,1,8,7,6,5,
     4             5,6,7,8,1,2,3,4,
     5             6,5,8,7,2,1,4,3,
     6             7,8,5,6,3,4,1,2,
     7             8,7,6,5,4,3,2,1/
      DATA JSEGNR /16,34,52,63,75,92,102,118,128,137,148,155,162,172,
     1 179,186,193,200,207,214,221,228/
      DATA JCOND  /12*1,4*-1,13*1,5*-1,13*1,5*-1,176*0/
      DATA KCOND  /7*1,9*0,5*1,13*0,5*1,13*0,7*1,4*0,7*1,5*0,11*1,6*0,
     1 5*1,5*0,10*1,6*0,5*1,5*0,5*1,4*0,5*1,6*0,1,1,5*0,1,1,5*0,4*1,
     2 62*0/
      DATA ITRK  /10,2,2,3,7,3,7,1,3*9,5*1,0,2*10,2,10,4*11,2*9,3,3,6*0,
     1  10,10,2,10,4*11,2*9,2*3,27*0,1,3*0,4*12,8*0,2*3,1,2*0,2*1,5*0,
     2  2*1,4*0,4*1,6*0,1,0,0,1,1,14*0,1,35*0,8,6*0,8,6*0,13,6*0,13,33*0
     3 /
      DATA NXTSEG /3*0,17,15,18,16,10,5,5,4,7,2,2,3,3,20,3*21,22,11,11,
     1 12,12,6,6,7,7,4*2,3,19,3*22,21,11,11,13,13,6,6,9,9,4*3,2,3*0,
     2 21,21,22,22,4*4,3*0,21,21,22,22,4*5,7,3*0,21,21,22,22,19,19,20,
     3 20,4*6,7,9,0,21,21,20,20,4*7,8,0,0,21,21,22,22,19,19,20,20,4*8,
     4 7,9,0,22,22,19,19,4*9,8,0,21,21,22,22,4*10,0,21,21,22,22,4*11,
     5 12,13,21,21,3*12,12,14,22,22,4*13,14,21,21,22,22,4*14,12,13,
     6 0,0,4*15,16,0,0,4*16,15,0,0,4*17     ,18,0,0,4*18,17,0,0,4*19,
     7 20,0,0,20,3*20,19,0,0,4*21,22,0,0,4*22,21/
      DATA IARCMN /4,3,2,3,4,2,4,4,2,3,4,3,3,4,2,4,1,3,4,3*2,4,3,4,1,2,
     1 1,3,1,2,3,4,2,1,2,4,3*3,4,2,4,1,3,1,2,1,2,3,4,3,2,3,4,3,4,2,4,
     2 1,2,3,4,2,3,4,3,4,2,4,1,2,3,4,3,2,3,4,3,4,2,4,1,2,1,3,1,2,3,4,
     3 3,2,2,2,4,1,2,1,2,3,4,2,2,3,3,4,2,4,1,2,1,3,1,2,3,4,3,2,3,3,4,
     4 1,3,1,2,3,4,3,1,1,3,1,2,1,2,3,4,1,1,3,1,2,1,2,3,4,3,2,1,2,1,2,
     5 3,4,2,1,3,1,2,3,4,3,1,3,1,2,1,2,3,4,3,2,1,2,1,2,3,4,2,1,3,1,2,
     6 3,4,3,1,2,1,2,3,4,2,1,3,1,2,3,4,3,3,4,1,2,3,4,3,2,4,1,2,3,4,2,
     7 1,2,1,2,3,4,2,1,3,1,2,3,4,3/
      DATA IARCSB/4,3,2,1,2,1,3,1,2,3,4,2,1,2,1,3,4,3,4,2,3,1,3,1,2,3,4,
     1 2,4,1,2,3,4,3,4,2,4,3,2,1,2,1,3,2,4,3,4,1,2,3,4,2,2,3,4,1,2,1,3,
     2 1,2,3,4,2,3,4,1,2,1,3,1,2,3,4,2,2,3,4,1,2,1,3,3,4,2,4,1,2,3,4,
     3 2,3,3,1,3,3,4,1,2,3,4,3,2,3,1,2,1,3,3,4,2,4,1,2,3,4,2,3,2,1,2,
     4 2,4,1,2,3,4,2,4,2,4,3,4,1,2,3,4,4,2,4,3,4,1,2,3,4,2,3,3,4,1,2,
     5 3,4,3,2,4,1,2,3,4,2,2,4,3,4,1,2,3,4,2,3,3,4,1,2,3,4,3,2,4,1,2,
     6 3,4,2,3,4,1,2,3,4,3,2,4,1,2,3,4,2,1,2,1,2,3,4,2,1,3,1,2,3,4,3,
     7 3,4,1,2,3,4,3,2,4,1,2,3,4,2/
      DATA JSEGPX /12,29,47/
      DATA NLCSMN /6*1,0,1,14*0/
      DATA ZERO,ONE,TWO,TEN10 /0.0D+00,1.0D+00,2.0D+00,1.0D-10/
      DATA PT5 /0.5D+00/
C
       CRITE = TEN10
       ROOT2 = SQRT(TWO)
       ROOTN2 = -ROOT2
       TOOR2 = ONE / ROOT2
       TOORN2 = -TOOR2
       JSEGPT(1)=0
       DO 130 I=1,21
 130   JSEGPT(I+1)=JSEGNR(I)
       DO 135 I=1,2
       DO 135 J=1,21
 135   COEFFS(I,J)=ZERO
       DO 137 I=1,NIJ
 137   IJADD(I)=IJADD(I)+(INTSAD-1)
       DO 140 I=3,20
       A = I-2
       COEFFS(I, 1) =  SQRT(A/(A+ONE))
       COEFFS(I, 2) = -COEFFS(I,1)
       COEFFS(I, 3) =  COEFFS(I,1)/SQRT(TWO)
       COEFFS(I, 4) = -COEFFS(I,3)
       COEFFS(I, 5) =  SQRT((A+ONE)/A)
       COEFFS(I, 6) = -COEFFS(I,5)
       COEFFS(I, 7) =  COEFFS(I,5)/SQRT(TWO)
       COEFFS(I, 8) = -COEFFS(I,7)
       COEFFS(I, 9) =  SQRT((A+TWO)/(A*TWO))
       COEFFS(I,10) = -COEFFS(I,9)
       COEFFS(I,11) =  SQRT(A/(TWO*(A+TWO)))
       COEFFS(I,12) = -COEFFS(I,11)
       COEFFS(I,13) =  SQRT(TWO/(A*(A+ONE)))
       COEFFS(I,14) =  SQRT(A*(A+TWO))/(A+ONE)
       COEFFS(I,15) = -COEFFS(I,14)
       COEFFS(I,16) =  SQRT((A-ONE)*(A+TWO)/(A*(A+ONE)))
       COEFFS(I,17) = -COEFFS(I,16)
       COEFFS(I,18) = -SQRT(TWO/(A*(A+TWO)))
       COEFFS(I,19) =  ONE/A
       COEFFS(I,20) = -COEFFS(I,19)
       COEFFS(I,21) = -SQRT(TWO)/A
 140   CONTINUE
       DO 150 K=1,4
       ISHIFT(K)=(K-1)*NROWS
 150   CONTINUE
       DO 155 I=1,NROWS
       NELEC(I)=NELEC(I)+NELEC(I)+NABC(I)
 155   CONTINUE
       DO 160 I=1,NSYM
       ISMOFF(I)=(I-1)*NORBS
 160   LCOND(I)=0
       I=ISYM(1)
       LCOND(I)=1
       LCOND(1)=1
       NSM=0
       DO 170 IORB=2,NORBS
       DO 165 I=1,NSYM
       IF(LCOND(I).EQ.0) GO TO 165
       LKUP=MULTS(I)+ISYM(IORB)
       J=LKUPSM(LKUP)
       IF(LCOND(J).GT.0) GO TO 165
       LCOND(J)=IORB
       NSM=NSM+1
       IF(NSM.EQ.NSYM) GO TO 175
 165   CONTINUE
 170   CONTINUE
 175   CONTINUE
       DO 180 I=1,NSYM
       IF(LCOND(I).EQ.0) LCOND(I)=NORBS+1
 180   CONTINUE
C
C     ----- LOOP OVER GROUPS OF INTEGRALS -----
C
       INXT=1
       DO 500 LGRP=1,NGRPS
       IF(OUT) WRITE(IW,187) LGRP
 187   FORMAT(' PROCESS INTEGRALS FROM GROUP',I4)
 190   CONTINUE
       I=INEXT(INXT)
       LEVI=I+1
       JMAX=JMXNXT(INXT)
       JMIN=JMNNXT(INXT)
       ILEVEL=I
       IAD=(I*(I-1))/2
       IJ=IAD+JMAX
       IF(LGRP.NE.IJGRP(IJ)) GO TO 490
C----
      NGO=0
      IF(NGO.EQ.0) GO TO 195
C----
       IF(LEVI.GT.LEVFRM) GO TO 195
       CALL LOPEXT(NABC,NLWKS,NUWKS,PUWK,NUWKY,NLWKX,PUWKY,NLWKW,NLWKY,
     *  IJADD,KADD,LADD,NELEC,NABCS,NROWS,NIJ,NKL,NLNG,LEVFRM)
       GO TO 485
 195   CONTINUE
       IISYM=ISYM(I)
       LEV=LEVI
       LEVM=LEV-1
       NR=LEVNR(LEV)
       NPT=LEVPT(LEV)
C
       DO 480 IROW=1,NR
       NPT=NPT+1
       ISEGM(LEV)=1
       ISEG=1
       IMN=NPT
       ISB=NPT
       KSEG=0
       KSEGMX=JSEGNR(ISEG)
       LMIN(LEV)=LCOND(IISYM)
       IUWKMN(LEV)=PUWK(NPT)
       IUWKSB(LEV)=PUWK(NPT)
       IMAIN(LEV)=NPT
       ISUB(LEV)=NPT
       NUWK=NUWKS(NPT)
       ACOEF(LEV)=ONE
C
C     ----- TEST NEXT SEGMENT OF GROUP -----
C
200   KSEG=KSEG+1
      IF(KSEG.GT.KSEGMX) GO TO 440
      KMN=IARCMN(KSEG)
      IARPT=IMN+ISHIFT(KMN)
      KMN=IARC(IARPT)
      IF(KMN.EQ.0) GO TO 200
      KSB=IARCSB(KSEG)
      JARPT=ISB+ISHIFT(KSB)
      KSB=IARC(JARPT)
      IF(KSB.EQ.0) GO TO 200
      JSEGM(LEV)=KSEG
      IUWKMN(LEVM)=IUWKMN(LEV)+IWGHT(IARPT)
      IUWKSB(LEVM)=IUWKSB(LEV)+IWGHT(JARPT)
      LMIN(LEVM)=LMIN(LEV)
      IF(JCOND(KSEG))220,240,230
220   CONTINUE
      IF(LEVM.LE.JMIN) GO TO 440
      GO TO 240
230   CONTINUE
      IF(LEVM.GT.JMAX) GO TO 420
      JLEVEL=LEVM
      IJ=IAD+LEVM
      JAD=IJADD(IJ)
      IJKUP=MULTS(IISYM)+ISYM(LEVM)
      IJSYM=LKUPSM(IJKUP)
      LMIN(LEVM)=LCOND(IJSYM)
      KSMPTX=ISMOFF(IJSYM)
240   CONTINUE
      IF(KCOND(KSEG).EQ.0) GO TO 260
      KLEVEL=LEVM
      KSMPT=LEVM+KSMPTX
      KAD=JAD+KADD(KSMPT)
      IJKKUP=MULTS(IJSYM)+ISYM(LEVM)
      IJKSYM=LKUPSM(IJKKUP)
      LMIN(LEVM)=LCOND(IJKSYM)
      LSMPTX=ISMOFF(IJKSYM)
260   CONTINUE
      IF(ITRK(KSEG))280,280,270
270   ITRACK(LEVM)=ITRK(KSEG)
      GO TO 290
280   ITRACK(LEVM)=ITRACK(LEV)
290   CONTINUE
      GO TO(50,1,1,1,3,1,4,1,44,45,51,1,1,3,1,4,5,52,53,6,54,40,41,1,7,
     A 46,47,1,8,1,6,2,2,9,10,52,53,11,55,42,43,1,12,48,49,1,13,1,2,11,
     B 2,36,77,77,78,77,79,77,80,1,1,1,1,87,68,82,68,69,67,70,71,75,
     C 76,71,83,87,68,82,68,69,67,70,68,69,67,70,71,75,76,71,83,83,
     D 6,6,74,6,8,1,16,16,1,17,18,19,19,20,18,21,19,20,18,21,1,22,23,
     E 1,24,24,11,11,81,11,13,1,27,27,1,28,1,3,2,4,2,1,2,2,1,29,63,64,
     F 65,66,71,72,73,71,84,85,56,57,1,30,1,1,86,56,58,1,1,31,1,32,
     G 59,60,61,62,1,33,88,1,34,35,1,5,1,6,2,2,9,1,10,1,2,11,2,36,
     H 1,5,1,6,2,2,9,1,10,1,2,11,2,36,1,4,1,37,2,2,38,1,3,1,2,37,2,39,
     I 1,5,1,6,2,2,9,1,10,1,2,11,2,36),KSEG
   1   ACOEF(LEVM)=ACOEF(LEV)
       GO TO 120
   2   ACOEF(LEVM)=-ACOEF(LEV)
       GO TO 120
   3   IA = NABC(IMN) + 2
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
   4   IA = NABC(IMN) + 83
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
   5   IA = NABC(IMN) + 82
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
   6   IA = NABC(IMN) + 261
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
   7   IA = NABC(IMN) + 1
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
   8   IA = NABC(IMN) + 102
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
   9   IA = NABC(IMN) + 362
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  10   IA = NABC(IMN) + 3
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  11   IA = NABC(IMN) + 263
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  12   IA = NABC(IMN) + 84
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  13   IA = NABC(IMN) + 23
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  16   IA = NABC(IMN) + 281
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  17   IA = NABC(IMN) + 402
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  18   IA = NABC(IMN) + 162
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  19   IA = NABC(IMN) + 222
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  20   IA = NABC(IMN) + 143
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  21   IA = NABC(IMN) + 42
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  22   IA = NABC(IMN) + 302
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  23   IA = NABC(IMN) + 303
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  24   IA = NABC(IMN) + 342
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  27   IA = NABC(IMN) + 283
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  28   IA = NABC(IMN) + 404
       ACOEF(LEVM)=ACOEF(LEV)*CFS(IA)
       GO TO 120
  29   ACOEF(LEVM) = ACOEF(LEV) * ROOT2
       GO TO 120
  30   IA = NABC(IMN) + 301
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  31   IA = NABC(IMN) + 304
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  32   IA = NABC(IMN) + 244
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  33   IA = NABC(IMN) + 322
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  34   IA = NABC(IMN) + 243
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  35   IA = NABC(IMN) + 242
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  36   IA = NABC(IMN) + 384
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  37   IA = NABC(IMN) + 262
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  38   IA = NABC(IMN) + 363
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  39   IA = NABC(IMN) + 383
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  40   IA = NABC(IMN) + 122
       IB = IA - 61
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       BCOEF(LEVM) = ACOEF(LEV) * CFS(IB)
       GO TO 120
  41   IB = NABC(IMN) + 162
       ACOEF(LEVM) = ACOEF(LEV) * TOORN2
       BCOEF(LEVM) = ACOEF(LEV) * CFS(IB)
       GO TO 120
  42   IA = NABC(IMN) + 43
       IB = IA + 81
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       BCOEF(LEVM) = ACOEF(LEV) * CFS(IB)
       GO TO 120
  43   IB = NABC(IMN) + 222
       ACOEF(LEVM) = ACOEF(LEV) * TOORN2
       BCOEF(LEVM) = ACOEF(LEV) * CFS(IB)
       GO TO 120
  44   IB=NABC(IMN)+221
       ACOEF(LEVM) = ACOEF(LEV) * TOOR2
       BCOEF(LEVM) = ACOEF(LEV) * CFS(IB)
       GO TO 120
  45   IB = NABC(IMN) + 163
       ACOEF(LEVM) = ACOEF(LEV) * TOOR2
       BCOEF(LEVM) = ACOEF(LEV) * CFS(IB)
       GO TO 120
  46   IB = NABC(IMN) + 162
       ACOEF(LEVM) = ACOEF(LEV) * TOOR2
       BCOEF(LEVM) = ACOEF(LEV) * CFS(IB)
       GO TO 120
  47   IA = NABC(IMN) + 122
       IB = IA - 81
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       BCOEF(LEVM) = ACOEF(LEV) * CFS(IB)
       GO TO 120
  48   IB = NABC(IMN) + 222
       ACOEF(LEVM) = ACOEF(LEV) * TOOR2
       BCOEF(LEVM) = ACOEF(LEV) * CFS(IB)
       GO TO 120
  49   IA = NABC(IMN) + 43
       IB = IA + 101
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       BCOEF(LEVM) = ACOEF(LEV) * CFS(IB)
       GO TO 120
  50   ACOEF(LEVM) = ACOEF(LEV) + ACOEF(LEV)
       D=PT5
       GO TO 120
  51   ACOEF(LEVM)=ACOEF(LEV)*ROOT2
       GO TO 120
  52   ACOEF(LEVM) = -ACOEF(LEV)
       D= -ONE
       GO TO 120
  53   ACOEF(LEVM) = -ACOEF(LEV) - ACOEF(LEV)
       D=-PT5
       GO TO 120
  54   IA = NABC(IMN) + 362
       D=ONE/CFS(IA)
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  55   IA = NABC(IMN) + 384
       D=ONE/CFS(IA)
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  56   ACOEF(LEVM) = ACOEF(LEV)
       D = -ONE
       GO TO 120
  57   IA = NABC(IMN) + 82
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       D=-ONE
       GO TO 120
  58   IA = NABC(IMN) + 3
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       D=-ONE
       GO TO 120
  59   IA = NABC(IMN) + 123
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       D=-ONE
       GO TO 120
  60   IA = NABC(IMN) + 222
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       D=-ONE
       GO TO 120
  61   IA = NABC(IMN) + 62
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       D=-ONE
       GO TO 120
  62   IA = NABC(IMN) + 162
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       D=-ONE
       GO TO 120
  63   IA = NABC(IMN) + 42
       IB = IA + 81
       ACOF = ACOEF(LEV) * CFS(IA)
       BCOF = BCOEF(LEV) * CFS(IB)
       D = ACOF + BCOF
       IF( ABS(D).LT.CRITE) GO TO 110
       ACOEF(LEVM) = D
       D = (ACOF - BCOF) / D
       GO TO 120
  64   IB = NABC(IMN) + 222
       ACOF = ACOEF(LEV) * TOORN2
       BCOF = BCOEF(LEV) * CFS(IB)
       D = ACOF + BCOF
       IF( ABS(D).LT.CRITE) GO TO 110
       ACOEF(LEVM) = D
       D = (ACOF - BCOF) / D
       GO TO 120
  65   IA = NABC(IMN) + 123
       IB = IA - 61
       ACOF = ACOEF(LEV) * CFS(IA)
       BCOF = BCOEF(LEV) * CFS(IB)
       D = ACOF + BCOF
       IF( ABS(D).LT.CRITE) GO TO 110
       ACOEF(LEVM) = D
       D = (ACOF - BCOF) / D
       GO TO 120
  66   IB = NABC(IMN) + 162
       ACOF = ACOEF(LEV) * TOORN2
       BCOF = BCOEF(LEV) * CFS(IB)
       D = ACOF + BCOF
       IF( ABS(D).LT.CRITE) GO TO 110
       ACOEF(LEVM) = D
       D = (ACOF - BCOF) / D
       GO TO 120
  67   IB = NABC(IMN) + 162
       DX=ACOEF(LEV)*TOORN2
       D=DX+BCOEF(LEV)*CFS(IB)
       IF( ABS(D).LT.CRITE) GO TO 111
       ACOEF(LEVM) = D
       D=-(DX+DX)/D
       GO TO 120
  68   IB = NABC(IMN) + 222
       DX=ACOEF(LEV)*TOORN2
       D=DX+BCOEF(LEV)*CFS(IB)
       IF( ABS(D).LT.CRITE) GO TO 111
       ACOEF(LEVM) = D
       D=-(DX+DX)/D
       GO TO 120
  69   IA = NABC(IMN) + 62
       IB = IA + 81
       DX=ACOEF(LEV)*CFS(IA)
       D=DX+BCOEF(LEV)*CFS(IB)
       IF( ABS(D).LT.CRITE) GO TO 111
       ACOEF(LEVM) = D
       D=-(DX+DX)/D
       GO TO 120
  70   IA = NABC(IMN) + 143
       IB = IA - 101
       DX=ACOEF(LEV)*CFS(IA)
       D=DX+BCOEF(LEV)*CFS(IB)
       IF( ABS(D).LT.CRITE) GO TO 111
       ACOEF(LEVM) = D
       D=-(DX+DX)/D
       GO TO 120
  71   ACOEF(LEVM) = ACOEF(LEV)
       BCOEF(LEVM) = BCOEF(LEV)
       GO TO 120
  72   IB = NABC(IMN) + 322
       ACOEF(LEVM) = -ACOEF(LEV)
       BCOEF(LEVM) = BCOEF(LEV) * CFS(IB)
       GO TO 120
  73   IB = NABC(IMN) + 323
       ACOEF(LEVM) = -ACOEF(LEV)
       BCOEF(LEVM) = BCOEF(LEV) * CFS(IB)
       GO TO 120
  74   IA = NABC(IMN) + 21
       ACOEF(LEVM)=ACOEF(LEV)*CFS(IA)
       GO TO 120
  75   IB = NABC(IMN) + 302
       ACOEF(LEVM) = ACOEF(LEV)
       BCOEF(LEVM) = BCOEF(LEV) * CFS(IB)
       GO TO 120
  76   IB = NABC(IMN) + 303
       ACOEF(LEVM) = ACOEF(LEV)
       BCOEF(LEVM) = BCOEF(LEV) * CFS(IB)
       GO TO 120
  77   ACOEF(LEVM)=ACOEF(LEV)*TOORN2
       D=-TWO
       GO TO 120
  78   ACOEF(LEVM)=ACOEF(LEV)*ROOTN2
       D=-TWO
       GO TO 120
  79   IA = NABC(IMN) + 62
       ACOEF(LEVM)=ACOEF(LEV)*CFS(IA)
       D=-TWO
       GO TO 120
  80   IA = NABC(IMN) + 143
       ACOEF(LEVM)=ACOEF(LEV)*CFS(IA)
       D=-TWO
       GO TO 120
  81   IA = NABC(IMN) + 104
       ACOEF(LEVM)=ACOEF(LEV)*CFS(IA)
       GO TO 120
  82   ACOEF(LEVM) = ACOEF(LEV) * ROOTN2
       D=-TWO
       GO TO 120
  83   IA = NABC(IMN) + 342
       ACOEF(LEVM) = BCOEF(LEV) * CFS(IA)
       GO TO 120
  84   IA = NABC(IMN) + 243
       ACOEF(LEVM) = BCOEF(LEV) * CFS(IA)
       GO TO 120
  85   IA = NABC(IMN) + 242
       ACOEF(LEVM) = BCOEF(LEV) * CFS(IA)
       GO TO 120
  86   IA = NABC(IMN) + 241
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
  87   IB = NABC(IMN) + 162
       DX=ACOEF(LEV)*TOORN2
       D=DX+BCOEF(LEV)*CFS(IB)
       IF( ABS(D).LT.CRITE) GO TO 111
       ACOEF(LEVM)=D
       D=-(DX+DX)/D
       GO TO 120
  88   IA = NABC(IMN) + 323
       ACOEF(LEVM) = ACOEF(LEV) * CFS(IA)
       GO TO 120
 110   ITRACK(LEVM)=3
       ACOEF(LEVM)=ACOF-BCOF
       GO TO 120
 111   ITRACK(LEVM) = 2
       ACOEF(LEVM)=-(DX+DX)
 120   CONTINUE
       IF(NXTSEG(KSEG).GT.0) GO TO 400
       IF(ISYM(LEVM).NE.IJKSYM) GO TO 200
       LLEVEL=LEVM
       LSMPT=LEVM+LSMPTX
       LAD=KAD+LADD(LSMPT)
       IF(KMN-KSB) 300,380,300
 300   LEVL=LEVM
       KSEGMX=4
 310   LEV=LEVM
       LEVM=LEV-1
       IF(LEVM.GT.0) GO TO 315
       IF (MASWRK) WRITE(IW,313)
 313   FORMAT(' PROBLEMS WITH PARTIAL SPACE')
       CALL ABRT
       RETURN
 315   CONTINUE
       KSEG=0
       IMAIN(LEV)=KMN
       IMN=KMN
       ISUB(LEV)=KSB
       ISB=KSB
 320   KSEG=KSEG+1
       IF(KSEG.GT.KSEGMX) GO TO 360
       IARPT=IMN+ISHIFT(KSEG)
       KMN=IARC(IARPT)
       IF(KMN.LE.0) GO TO 320
       JARPT=ISB+ISHIFT(KSEG)
       KSB=IARC(JARPT)
       IF(KSB.LE.0) GO TO 320
       JSEGM(LEV)=KSEG
       IUWKMN(LEVM)=IUWKMN(LEV)+IWGHT(IARPT)
       IUWKSB(LEVM)=IUWKSB(LEV)+IWGHT(JARPT)
       IF(KMN-KSB) 310,340,310
 340   NLWK = NLWKS(KMN)
       IUWK = IUWKMN(LEVM)
       JUWK = IUWKSB(LEVM)
       ITRAK = ITRACK(LEVL)
       ACF = ACOEF(LEVL)
       CALL TM2OUT(NWKS)
       GO TO 320
 360   IF(LEV.EQ.LEVL) GO TO 440
       LEVM=LEV
       LEV=LEVM+1
       IMN=IMAIN(LEV)
       ISB=ISUB(LEV)
       KSEG=JSEGM(LEV)
       GO TO 320
 380   NLWK = NLWKS(KMN)
       IUWK = IUWKMN(LEVM)
       JUWK = IUWKSB(LEVM)
       ITRAK = ITRACK(LEVM)
       ACF = ACOEF(LEVM)
       CALL TM2OUT(NWKS)
       GO TO 200
 400   CONTINUE
       IF(LEVM.LE.LMIN(LEVM)) GO TO 200
       ISEG=NXTSEG(KSEG)
       IF(NLCSMN(ISEG).GT.0.AND.NELEC(KMN).EQ.0) GO TO 200
C----  IF(LVFRM.EQ.LEVM.AND.ISEG.LT.4) GO TO 460
       LEV=LEVM
       LEVM=LEV-1
       ISEGM(LEV)=ISEG
       KSEG=JSEGPT(ISEG)
       IMN=KMN
       IMAIN(LEV)=KMN
       ISB=KSB
       ISUB(LEV)=KSB
       KSEGMX=JSEGNR(ISEG)
       GO TO 200
 420   CONTINUE
       KSEG=JSEGPX(ISEG)
       GO TO 200
 440   CONTINUE
       IF(LEV.EQ.LEVI) GO TO 480
       LEVM=LEV
       LEV=LEVM+1
       ISEG=ISEGM(LEV)
       IMN=IMAIN(LEV)
       ISB=ISUB(LEV)
       KSEG=JSEGM(LEV)
       KSEGMX=JSEGNR(ISEG)
       GO TO 200
C--- 460   CONTINUE
C---       IUWK=IUWKMN(LEVM)
C---       JUWK=IUWKSB(LEVM)
C---       ACF=ACOEF(LEVM)
C---       IISYMT=IISYM
C---       KSBT=KSB
C---       CALL THREXT(NABC,NLWKS,NUWKS,PUWK,NUWKY,NLWKX,PUWKY,
C---                   NLWKW,NLWKY,IJADD,KADD,LADD,NELEC,NABCS,
C---                   NROWS,NIJ,NKL,NLNG)
C---       GO TO 200
 480   CONTINUE
C
C     ----- END OF LOOP OVER ROWS OF -LEVI- -----
C
 485   CONTINUE
       INXT=INXT+1
       IF(INXT.LE.NEXT) GO TO 190
       GO TO 500
C
C     ----- READ IN NEW GROUP OF INTEGRALS -----
C
 490   CONTINUE
       CALL TM2GRP(NWKS)
 500   CONTINUE
       CALL TM2FIN(NWKS)
       RETURN
       END
C*MODULE MCQUD   *DECK LOOPXT
      SUBROUTINE LOOPXT(NABCB,NLWKS,NUWKS,PUWK,NUWKY,NLWKX,PUWKY,NLWKW,
     *                  NLWKY,IJADD,KADD,LADD,NABCA,NABCS,NROWS,NIJ,
     *                  NKL,NLNG,LEVFRM)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      INTEGER PUWK,PUWKX,PUWKW,PUWKY
      INTEGER PUWK3W,PUWK3X,PUWK1W,PUWK1X,PUWK2W,PUWK2X
      LOGICAL SOME,OUT,GOPARR,DSKWRK,MASWRK
C
      PARAMETER (MXAO=8192)
C
      COMMON /EMPRT / ILEVEL,JLEVEL,KLEVEL,LLEVEL,OUT,SOME
      COMMON /IOFILE/ IR,IW,IP,IJK,IPK,IDAF,NAV,IODA(950)
      COMMON /LOOPS / ACOEF,CCOEF,DCOEF,VAL,CUTOFF,LAD,NUWK,NLWK,
     *                IUWK,JUWK,ITRACK,LEVI,JMIN,JMAX,IADT,JSMT,KSB,
     *                INTSAD,NINTMX,NGRPS
      COMMON /ORBSET/ NORBMX,NORBS,NCORBS,NLEVS,NNA,NNB,NNC,NSYM,MSYM,
     *                IDOCC,IVAL,IMCC,ISYM(MXAO),JCODE(MXAO),
     *                NLCS(MXAO),LEVPT(MXAO),LEVNR(MXAO),IOUT(MXAO),
     *                NREFS,IEXCT,NFOCI,INTACT
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
C
      DIMENSION NUWKW(8),NUWKX(8),ISMOFF(8),PUWKW(8),PUWKX(8)
      DIMENSION MSMOFF(8),MFRMPT(25)
      DIMENSION NLWKW(NLNG),NLWKX(NLNG),NLWKY(NLNG),NUWKY(NLNG),
     *          PUWKY(NLNG)
      DIMENSION NABCB(NROWS),NLWKS(NROWS),NUWKS(NROWS),PUWK(NROWS),
     *          IJADD(NIJ),KADD(NKL),LADD(NKL),NABCA(NROWS),NABCS(NROWS)
      DIMENSION MULTS(8),LKUPSM(64)
C
      DATA MULTS  /0,8,16,24,32,40,48,56/
      DATA LKUPSM /1,2,3,4,5,6,7,8,
     1             2,1,4,3,6,5,8,7,
     2             3,4,1,2,7,8,5,6,
     3             4,3,2,1,8,7,6,5,
     4             5,6,7,8,1,2,3,4,
     5             6,5,8,7,2,1,4,3,
     6             7,8,5,6,3,4,1,2,
     7             8,7,6,5,4,3,2,1/
      DATA ONE,TWO /1.0D+00,2.0D+00/
C
      SQRT2= SQRT(TWO)
      DO 20 ISM=1,NSYM
      ISMOFF(ISM)=(ISM-1)*LEVFRM
      MSMOFF(ISM)=(ISM-1)*NORBS
      MFRM=LEVPT(LEVFRM)
      DO 10 I=1,LEVFRM
      IPT=ISMOFF(ISM)+I
      NLWKW(IPT)=0
      NLWKX(IPT)=0
      NLWKY(IPT)=0
      NUWKY(IPT)=0
      PUWKY(IPT)=0
10    CONTINUE
      NUWKW(ISM)=0
      NUWKX(ISM)=0
      PUWKW(ISM)=0
      PUWKX(ISM)=0
20    CONTINUE
      DO 30 IM=1,25
30    MFRMPT(IM)=0
      DO 100 LEV=2,LEVFRM
      NPT=LEVPT(LEV)
      NR=LEVNR(LEV)
      DO 90 IR0=1,NR
      NPT=NPT+1
      IA=NABCA(NPT)
      IF(IA.GE.2) GO TO 930
      IB=NABCB(NPT)
      ISM=NABCS(NPT)
      ISMPT=ISMOFF(ISM)+LEV
      IF(IA.EQ.1) GO TO 70
      IF(IB.EQ.2) GO TO 80
      IF(IB.EQ.0) GO TO 90
C
C     ----- THIS IS A Y POINT -----
C
      NUWKY(ISMPT)=NUWKS(NPT)
      NLWKY(ISMPT)=NLWKS(NPT)
      PUWKY(ISMPT)=PUWK(NPT)
      IF(LEV.LT.LEVFRM) GO TO 90
      GO TO 90
70    CONTINUE
C
C     ----- THIS IS A W POINT -----
C
      NLWKW(ISMPT)=NLWKS(NPT)
      IF(LEV.LT.LEVFRM) GO TO 90
      NUWKW(ISM)=NUWKS(NPT)
      PUWKW(ISM)=PUWK(NPT)
      MFRMPT(IR0)=ISM
      GO TO 90
80    CONTINUE
C
C     ----- THIS IS A X POINT -----
C
      NLWKX(ISMPT)=NLWKS(NPT)
      IF(LEV.LT.LEVFRM) GO TO 90
      NUWKX(ISM)=NUWKS(NPT)
      PUWKX(ISM)=PUWK(NPT)
      MFRMPT(IR0)=-ISM
90    CONTINUE
100   CONTINUE
      RETURN
C
C     ----- ENTRY POINT -LOPEXT- -----
C
      ENTRY LOPEXT(NABCB,NLWKS,NUWKS,PUWK,NUWKY,NLWKX,
     * PUWKY,NLWKW,NLWKY,IJADD,KADD,LADD,NABCA,NABCS,NROWS,NIJ,
     * NKL,NLNG,LEVFRM)
      I=LEVI-1
      NLWK=1
      ILEVEL=I
      IAD=(I*(I-1))/2
      JSM=ISYM(I)
      LKUPJ=MULTS(JSM)
      JKIND=1
      CCOEF=ONE
      DCOEF=-ONE
      NUWK5W=NUWKW(1)
      JUWK5W=PUWKW(1)
      IUWK5W=JUWK5W+NLWKW(I)
      DO 800 J=JMIN,JMAX
      JLEVEL=J
      IJ=IAD+J
      JAD=IJADD(IJ)
      ISJ=ISYM(J)
      LKUPSJ=MULTS(ISJ)
      LKUP=LKUPJ+ISJ
      IHSM3=LKUPSM(LKUP)
      LKUPK=MULTS(IHSM3)
      KSMPTX=MSMOFF(IHSM3)
      JSOFF=ISMOFF(ISJ)
      IH3OFF=ISMOFF(IHSM3)
      NUWK3W=NUWKW(IHSM3)
      NUWK3X=NUWKX(IHSM3)
      IPT1=IH3OFF+I
      IPT2=JSOFF+J
      PUWK3W=PUWKW(IHSM3)
      PUWK3X=PUWKX(IHSM3)
      IUWK3W=PUWK3W+NLWKW(IPT1)+NLWKY(IPT2)
      IUWK3X=PUWK3X+NLWKX(IPT1)+NLWKY(IPT2)
      IF(I.NE.J) GO TO 160
C
C     ----- CASES WHERE I=J -----
C
      JKIND=2
      IUWK9W=JUWK5W+NLWKW(I+1)-1
160   CONTINUE
      KKIND=JKIND
      DO 750 K=1,J
      KLEVEL=K
      KSMPT=K+KSMPTX
      KAD=JAD+KADD(KSMPT)
      ISK=ISYM(K)
      LKUP=LKUPK+ISK
      LSM=LKUPSM(LKUP)
      LKUPL=MULTS(LSM)
      LSOFF=ISMOFF(LSM)
      LSMPTX=MSMOFF(LSM)
      KSOFF=ISMOFF(ISK)
      ACOEF=ONE
      IPT1=IH3OFF+K
      JUWK3W=PUWK3W+NLWKW(IPT1)
      JUWK3X=PUWK3X+NLWKX(IPT1)
      IF(K.EQ.J) GO TO 170
      LKUP=LKUPJ+ISK
      IHSM1=LKUPSM(LKUP)
      IH1OFF=ISMOFF(IHSM1)
      NUWK1W=NUWKW(IHSM1)
      NUWK1X=NUWKX(IHSM1)
      PUWK1W=PUWKW(IHSM1)
      PUWK1X=PUWKX(IHSM1)
      IPT1=IH1OFF+I
      IPT2=KSOFF+K
      IUWK1W=PUWK1W+NLWKW(IPT1)+NLWKY(IPT2)
      IUWK1X=PUWK1X+NLWKX(IPT1)+NLWKY(IPT2)
      IPT1=IH1OFF+J
      JUWK1W=PUWK1W+NLWKW(IPT1)
      JUWK1X=PUWK1X+NLWKX(IPT1)
      LKUP=LKUPSJ+ISK
      IHSM2=LKUPSM(LKUP)
      IH2OFF=ISMOFF(IHSM2)
      NUWK2W=NUWKW(IHSM2)
      NUWK2X=NUWKX(IHSM2)
      PUWK2W=PUWKW(IHSM2)
      PUWK2X=PUWKX(IHSM2)
      IPT1=IH2OFF+I
      IUWK2W=PUWK2W+NLWKW(IPT1)
      IUWK2X=PUWK2X+NLWKX(IPT1)
      IPT1=IH2OFF+J
      JUWK2W=PUWK2W+NLWKW(IPT1)+NLWKY(IPT2)
      JUWK2X=PUWK2X+NLWKX(IPT1)+NLWKY(IPT2)
      GO TO 180
C
C     ----- CASES WHERE J=K -----
C
170   CONTINUE
      KKIND=KKIND+2
      JUWK7W=JUWK5W+NLWKW(K+1)-1
      JUWK9W=JUWK5W+NLWKW(K)
180   CONTINUE
      LKIND=KKIND
      DO 700 L=1,K
      IF(ISYM(L).NE.LSM) GO TO 700
      LLEVEL=L
      LSMPT=L+LSMPTX
      LAD=KAD+LADD(LSMPT)
      IHL=LSOFF+L
      NLWKYL=NLWKY(IHL)
      IF(L.EQ.K)LKIND=LKIND+4
      GO TO (250,300,350,400,500,550,700,600),LKIND
250   CONTINUE
C
C     ----- CASES WHERE ALL INDICES ARE DIFFERENT -----
C
C  TYPE 1W
      NUWK=NUWK1W
      IUWK=IUWK1W
      JUWK=JUWK1W+NLWKYL
      ITRACK=6
      CALL TM2EXT
C  TYPE 1X
      NUWK=NUWK1X
      IUWK=IUWK1X
      JUWK=JUWK1X+NLWKYL
      ITRACK =15
      CALL TM2EXT
C  TYPE 2W
      NUWK=NUWK2W
      IUWK=IUWK2W+NLWKYL
      JUWK=JUWK2W
      ITRACK=4
      CALL TM2EXT
C  TYPE 2X
      NUWK=NUWK2X
      IUWK=IUWK2X+NLWKYL
      JUWK=JUWK2X
      ITRACK=10
      CALL TM2EXT
C  TYPE 3W
      NUWK=NUWK3W
      IUWK=IUWK3W
      JUWK=JUWK3W+NLWKYL
      ITRACK=14
      CALL TM2EXT
C  TYPE 3X
      NUWK=NUWK3X
      IUWK=IUWK3X
      JUWK=JUWK3X+NLWKYL
      ITRACK=11
      CALL TM2EXT
      GO TO 700
300   CONTINUE
C
C     ----- CASES WHERE I=J -----
C
C  TYPE 8W
      NUWK=NUWK1W
      IUWK=IUWK1W
      JUWK=JUWK1W+NLWKYL
      ITRACK=4
      CALL TM2EXT
C  TYPE 8X
      NUWK=NUWK1X
      IUWK=IUWK1X
      JUWK=JUWK1X+NLWKYL
      ITRACK=10
      CALL TM2EXT
C  TYPE 10W
      NUWK=NUWK5W
      IUWK=IUWK9W
      JUWK=JUWK3W+NLWKYL
      ITRACK=16
      CALL TM2EXT
      GO TO 700
350   CONTINUE
C
C     ----- CASES WHERE J=K (TYPES 6,7) -----
C
C  TYPE 6W
      NUWK=NUWK3W
      IUWK=IUWK3W
      JUWK=JUWK3W+NLWKYL
      ITRACK=4
      CALL TM2EXT
C  TYPE 6X
      NUWK=NUWK3X
      IUWK=IUWK3X
      JUWK=JUWK3X+NLWKYL
      ITRACK=9
      CALL TM2EXT
C  TYPE 7W
      NUWK=NUWK5W
      IUWK=IUWK5W+NLWKYL
      JUWK=JUWK7W
      ITRACK=16
      CALL TM2EXT
      GO TO 700
400   CONTINUE
C
C     ----- CASES WHERE I=J=K (TYPE 12) -----
C
C  TYPE 12(A+C)W
      NUWK=NUWK5W
      IUWK=IUWK5W+NLWKYL
      JUWK=JUWK5W+NLWKW(L+1)-1
      ACOEF=SQRT2
      ITRACK=13
      CALL TM2EXT
C  TYPE 12(B+C)W
      IUWK=IUWK9W
      JUWK=JUWK9W+NLWKYL
      ITRACK=7
      CALL TM2EXT
C  TYPE 12WY AND 12XY
      ACOEF=ONE
      ITRACK=3
      IPT1=L
      DO 410 IZ=1,NSYM
      NLWK=NLWKY(IPT1)
      IF(NLWK.EQ.0) GO TO 410
      LKUP=LKUPL+IZ
      IHSYM=LKUPSM(LKUP)
      KPT1=ISMOFF(IHSYM)+I
      KPT2=ISMOFF(IHSYM)+L
      NUWK=NUWKW(IHSYM)
      IUWK=PUWKW(IHSYM)+NLWKW(KPT1)
      JUWK=PUWKW(IHSYM)+NLWKW(KPT2)
      CALL TM2EXT
      NUWK=NUWKX(IHSYM)
      IUWK=PUWKX(IHSYM)+NLWKX(KPT1)
      JUWK=PUWKX(IHSYM)+NLWKX(KPT2)
      CALL TM2EXT
410   IPT1=IPT1+LEVFRM
C  TYPE 12YZ
      NLWK=1
      IPT1=LSOFF+I
      NUWK=NUWKY(IPT1+1)
      IUWK=PUWKY(IPT1+1)+NLWKY(IPT1)
      JUWK=PUWKY(IPT1+1)+NLWKYL
      CALL TM2EXT
C  TYPES 12WZ AND 12XZ
      IZMAX=I-1
      IZMIN=L+1
      IF(IZMAX.LT.IZMIN) GO TO 145
      DO 440 IZ=IZMIN,IZMAX
      IZSYM=ISYM(IZ)
      LKUP=LKUPL+IZSYM
      IHSYM=LKUPSM(LKUP)
      KPT1=ISMOFF(IHSYM)+I
      KPT2=ISMOFF(IHSYM)+IZ
      IPT1=ISMOFF(IZSYM)+IZ
      NUWK=NUWKW(IHSYM)
      IUWK=PUWKW(IHSYM)+NLWKW(KPT1)+NLWKY(IPT1)
      JUWK=PUWKW(IHSYM)+NLWKW(KPT2)+NLWKYL
      ACOEF=ONE
      CALL TM2EXT
      NUWK=NUWKX(IHSYM)
      IUWK=PUWKX(IHSYM)+NLWKX(KPT1)+NLWKY(IPT1)
      JUWK=PUWKX(IHSYM)+NLWKX(KPT2)+NLWKYL
      ACOEF=-ONE
      CALL TM2EXT
440   CONTINUE
145   CONTINUE
      ACOEF=ONE
      GO TO 700
500   CONTINUE
C
C     ----- CASES WHERE L=K (TYPES 4,5) -----
C
C  TYPE 4W
      NUWK=NUWK1W
      IUWK=IUWK1W
      JUWK=JUWK1W+NLWKYL
      ITRACK=4
      CALL TM2EXT
C  TYPE 4X
      NUWK=NUWK1X
      IUWK=IUWK1X
      JUWK=JUWK1X+NLWKYL
      ITRACK=10
      CALL TM2EXT
C  TYPE 5W
      NUWK=NUWK5W
      IUWK=IUWK3W
      JUWK=JUWK5W+NLWKW(L+1)-1
      ITRACK=16
      CALL TM2EXT
      GO TO 700
550   CONTINUE
C
C     ----- CASES WHERE I=J AND K=L (TYPES 11,13) -----
C
C  TYPE 11W
      NUWK=NUWK1W
      IUWK=IUWK1W
      JUWK=JUWK1W+NLWKYL
      ITRACK=4
      CALL TM2EXT
C  TYPE 11X
      NUWK=NUWK1X
      IUWK=IUWK1X
      JUWK=JUWK1X+NLWKYL
      ITRACK=10
      CALL TM2EXT
C  TYPE 13W
      NUWK=NUWK5W
      IUWK=IUWK9W
      JUWK=JUWK5W+NLWKW(L+1)-1
      ITRACK=1
      CALL TM2EXT
      GO TO 700
600   CONTINUE
C
C     ----- CASES WHERE I=J=K=L (TYPE 14) -----
C
C  TYPE 14(A+B)W
      NUWK=NUWK5W
      IUWK=IUWK9W
      JUWK=IUWK
      ITRACK=5
      CALL TM2EXT
C  TYPES 14WY AND 14XY
      ITRACK=2
      IPT1=I
      DO 630 IZ=1,NSYM
      NLWK=NLWKY(IPT1)
      IF(NLWK.EQ.0) GO TO 630
      LKUP=LKUPL+IZ
      IHSYM=LKUPSM(LKUP)
      KPT1=ISMOFF(IHSYM)+I
      NUWK=NUWKW(IHSYM)
      IUWK=PUWKW(IHSYM)+NLWKW(KPT1)
      JUWK=IUWK
      CALL TM2EXT
      NUWK=NUWKX(IHSYM)
      IUWK=PUWKX(IHSYM)+NLWKX(KPT1)
      JUWK=IUWK
      CALL TM2EXT
630   IPT1=IPT1+LEVFRM
C  TYPE 14YZ
      NLWK=1
      NUWK=NUWKY(IHL+1)
      IUWK=PUWKY(IHL+1)+NLWKYL
      JUWK=IUWK
      CALL TM2EXT
700   CONTINUE
750   CONTINUE
800   CONTINUE
      RETURN
C
C     ----- ENTRY POINT -THREXT- -----
C
C     ----- THREXT DOES THE THREE EXTERNAL LOOPS -----
C
      ENTRY THREXT(NABCB,NLWKS,NUWKS,PUWK,NUWKY,NLWKX,
     * PUWKY,NLWKW,NLWKY,IJADD,KADD,LADD,NABCA,NABCS,NROWS,NIJ,
     * NKL,NLNG)
      CCOEF=ONE
      DCOEF=-ONE
      ACOEFT=ACOEF
      IUWKT=IUWK
      JUWKT=JUWK
      NLWK=1
      JMN=JMIN
      IF(JMIN.LT.2)JMN=2
      JMX=JMAX
      IF(JMAX.GE.LEVFRM)JMX=LEVFRM-1
      KSBL=KSB-MFRM
      JFSYM=MFRMPT(KSBL)
      IF(JFSYM)2010,940,1000
1000  CONTINUE
      IH1OFF=ISMOFF(JFSYM)
      LKUPJ=MULTS(JSMT)
      LKUP=LKUPJ+JFSYM
      IFSYM=LKUPSM(LKUP)
      IF(JMN.LT.JMX) RETURN
      DO 2000 J=JMN,JMX
      KKIND=1
      JLEVEL=J
      IJ=IADT+J
      JAD=IJADD(IJ)
      ISJ=ISYM(J)
      LKUPSJ=MULTS(ISJ)
      LKUP=LKUPJ+ISJ
      KSM=LKUPSM(LKUP)
      LKUPK=MULTS(KSM)
      IPT1=IH1OFF+J
      JUWK1W=JUWKT+NLWKW(IPT1)
      KSMPTX=MSMOFF(KSM)
      JSOFF=ISMOFF(ISJ)
      IPT2=JSOFF+J
      IUWK3W=IUWKT+NLWKY(IPT2)
      DO 1900 K=1,J
      IHSM3=1
      IF(IFSYM.EQ.ISJ)IHSM3=0
      KLEVEL=K
      KSMPT=K+KSMPTX
      KAD=JAD+KADD(KSMPT)
      ISK=ISYM(K)
      LKUP=LKUPK+ISK
      LSM=LKUPSM(LKUP)
      LSOFF=ISMOFF(LSM)
      LSMPTX=MSMOFF(LSM)
      KSOFF=ISMOFF(ISK)
      IPT1=IH1OFF+K
      JUWK3W=JUWKT+NLWKW(IPT1)
      IHSM1=1
      IF(IFSYM.EQ.ISK)IHSM1=0
      LKUP=LKUPSJ+ISK
      IHSYM2=LKUPSM(LKUP)
      IHSM2=1
      IF(IHSYM2.EQ.JFSYM)IHSM2=0
      IF(J.EQ.K) GO TO 1200
      IPT2=KSOFF+K
      IUWK1W=IUWKT+NLWKY(IPT2)
      JUWK2W=JUWK1W+NLWKY(IPT2)
      GO TO 1250
1200  CONTINUE
      KKIND=KKIND+1
      JUWK7W=JUWKT+NLWKW(K+1)-1
1250  CONTINUE
      LKIND=KKIND
      DO 1800 L=1,K
      IF(ISYM(L).NE.LSM) GO TO 1800
      LLEVEL=L
      LSMPT=L+LSMPTX
      LAD=KAD+LADD(LSMPT)
      IHL=LSOFF+L
      NLWKYL=NLWKY(IHL)
      IF(L.EQ.K)LKIND=LKIND+2
      GO TO (1300,1400,1500,1800),LKIND
1300  CONTINUE
C
C     ----- CASES WHERE NO INDICES ARE THE SAME (1-3) -----
C
C TYPE 1W
      IF(IHSM1.NE.0) GO TO 1310
      IUWK=IUWK1W
      JUWK=JUWK1W+NLWKYL
      ITRACK=13
      CALL TM2EXT
C TYPE 2W
1310  CONTINUE
      IF(IHSM2.NE.0) GO TO 1320
      IUWK=IUWKT+NLWKYL
      JUWK=JUWK2W
      ITRACK=17
      CALL TM2EXT
C TYPE 3W
1320  CONTINUE
      IF(IHSM3.NE.0) GO TO 1800
      IUWK=IUWK3W
      JUWK=JUWK3W+NLWKYL
      ITRACK=7
      CALL TM2EXT
      GO TO 1800
1400  CONTINUE
C
C     ----- CASES WHERE J=K (6,7) -----
C
C TYPE 6W
      IF(IHSM3.NE.0) GO TO 1410
      IUWK=IUWK3W
      JUWK=JUWK3W+NLWKYL
      ITRACK=17
      CALL TM2EXT
C TYPE 7W
1410  CONTINUE
      IF(JFSYM.NE.1) GO TO 1800
      IUWK=IUWKT+NLWKYL
      JUWK=JUWK7W
      ITRACK=1
      ACOEF=ACOEF*SQRT2
      CALL TM2EXT
      ACOEF=ACOEFT
      GO TO 1800
1500  CONTINUE
C
C     ----- CASES WHERE L=K (4,5) -----
C
C TYPE 4W
      IF(IHSM1.NE.0) GO TO 1510
      IUWK=IUWK1W
      JUWK=JUWK1W+NLWKYL
      ITRACK=17
      CALL TM2EXT
C TYPE 5W
1510  CONTINUE
      IF(JFSYM.NE.1) GO TO 1800
      IUWK=IUWK3W
      JUWK=JUWKT+NLWKW(L+1)-1
      ITRACK=1
      ACOEF=ACOEF*SQRT2
      CALL TM2EXT
      ACOEF=ACOEFT
1800  CONTINUE
1900  CONTINUE
2000  CONTINUE
      RETURN
2010  CONTINUE
      JFSYM=-JFSYM
      IH1OFF=ISMOFF(JFSYM)
      LKUPJ=MULTS(JSMT)
      LKUP=LKUPJ+JFSYM
      IFSYM=LKUPSM(LKUP)
      IF(JMN.LT.JMX) RETURN
      DO 4000 J=JMN,JMX
      KKIND=1
      JLEVEL=J
      IJ=IADT+J
      JAD=IJADD(IJ)
      ISJ=ISYM(J)
      LKUPSJ=MULTS(ISJ)
      LKUP=LKUPJ+ISJ
      KSM=LKUPSM(LKUP)
      LKUPK=MULTS(KSM)
      KSMPTX=MSMOFF(KSM)
      IPT1=IH1OFF+J
      JUWK1X=JUWKT+NLWKX(IPT1)
      JSOFF=ISMOFF(ISJ)
      IPT2=JSOFF+J
      IUWK3X=IUWKT+NLWKY(IPT2)
      DO 3900 K=1,J
      IHSM3=1
      IF(IFSYM.EQ.ISJ)IHSM3=0
      KLEVEL=K
      KSMPT=K+KSMPTX
      KAD=JAD+KADD(KSMPT)
      ISK=ISYM(K)
      LKUP=LKUPK+ISK
      LSM=LKUPSM(LKUP)
      LSOFF=ISMOFF(LSM)
      LSMPTX=MSMOFF(LSM)
      KSOFF=ISMOFF(ISK)
      IPT1=IH1OFF+K
      JUWK3X=JUWKT+NLWKX(IPT1)
      IHSM1=1
      IF(IFSYM.EQ.ISK)IHSM1=0
      LKUP=LKUPSJ+ISK
      IHSYM2=LKUPSM(LKUP)
      IHSM2=1
      IF(IHSYM2.EQ.JFSYM)IHSM2=0
      IF(J.EQ.K) GO TO 3200
      IPT2=KSOFF+K
      IUWK1X=IUWKT+NLWKY(IPT2)
      JUWK2X=JUWK1X+NLWKY(IPT2)
      GO TO 3250
3200  CONTINUE
      KKIND=KKIND+1
3250  CONTINUE
      LKIND=KKIND
      DO 3800 L=1,K
      IF(ISYM(L).NE.LSM) GO TO 3800
      LLEVEL=L
      LSMPT=L+LSMPTX
      LAD=KAD+LADD(LSMPT)
      IHL=LSOFF+L
      NLWKYL=NLWKY(IHL)
      IF(L.EQ.K)LKIND=LKIND+2
      GO TO (3300,3400,3500,3800),LKIND
3300  CONTINUE
C
C     ----- CASES WHERE NO INDICES ARE THE SAME (1-3) -----
C
C TYPE 1X
      IF(IHSM1.NE.0) GO TO 3310
      IUWK=IUWK1X
      JUWK=JUWK1X+NLWKYL
      ITRACK=18
      CALL TM2EXT
C TYPE 2X
3310  CONTINUE
      IF(IHSM2.NE.0) GO TO 3320
      IUWK=IUWKT+NLWKYL
      JUWK=JUWK2X
      ITRACK=10
      CALL TM2EXT
C TYPE 3X
3320  CONTINUE
      IF(IHSM3.NE.0) GO TO 3800
      IUWK=IUWK3X
      JUWK=JUWK3X+NLWKYL
      ITRACK=11
      CALL TM2EXT
      GO TO 3800
3400  CONTINUE
C
C     ----- CASES WHERE J=K (6,7) -----
C
C TYPE 6X
      IF(IHSM3.NE.0) GO TO 3800
      IUWK=IUWK3X
      JUWK=JUWK3X+NLWKYL
      ITRACK=9
      CALL TM2EXT
      GO TO 3800
3500  CONTINUE
C
C     ----- CASES WHERE L=K (4,5) -----
C
C TYPE 4X
      IF(IHSM1.NE.0) GO TO 3800
      IUWK=IUWK1X
      JUWK=JUWK1X+NLWKYL
      ITRACK=10
      CALL TM2EXT
3800  CONTINUE
3900  CONTINUE
4000  CONTINUE
      RETURN
C
  930 CONTINUE
      IF(MASWRK) WRITE(IW,*) 'LOOPXT: PROBLEM WITH IA=',IA
      CALL ABRT
      RETURN
C
  940 CONTINUE
      IF (MASWRK) WRITE(IW,950)
  950 FORMAT(//,10X,'ENTERED A Y POINT')
      CALL ABRT
      RETURN
      END
C*MODULE MCQUD   *DECK GUGACI
      SUBROUTINE GUGACI(NPRINT)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      LOGICAL SOME,OUT,GOPARR,DSKWRK,MASWRK
C
      COMMON /CIFILS/ NFT11,NFT12,NFT13,NFT14,NFT15,NFT16,IDAF20,NEMEMX
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
      COMMON /FMCOM / X(1)
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
      COMMON /RUNOPT/ RUNTYP,EXETYP,NEVALS,NGLEVL,NHLEVL
C
      DIMENSION TITLE(10),TITLE1(10)
      DIMENSION TITLA(10),TITLA1(10)
C
      DATA NWORD  /0/
      DATA CHECK/8HCHECK   /
C
      SOME=NPRINT.NE.-5 .AND. MASWRK
      OUT =NPRINT.GT.1 .AND. MASWRK
      IF(SOME) WRITE(IW,9999)
C
C     ----- GET CORE MEMORY -----
C
      CALL GOTFM(NGOT)
      IF(NWORD.GT.0) NGOT=MIN(NWORD,NGOT)
C
      CALL SEQREW(NFT15)
      CALL SEQREW(NFT16)
      READ(NFT15) ISET1,NWKS,NEMEMX,CUTOFF,TITLA,TITLA1
      READ(NFT16) ISET2,NWKS,NEMEMX,CUTOFF,TITLE,TITLE1
C
      IF(ISET1.NE.1.OR.ISET2.NE.2) GO TO 9000
      DO 20 I=1,10
      IF(TITLE (I).NE.TITLA (I)) GO TO 9000
      IF(TITLE1(I).NE.TITLA1(I)) GO TO 9000
   20 CONTINUE
      IF(OUT) WRITE(IW,9997) TITLE,TITLE1,NWKS,CUTOFF
C
C     ----- IN CORE -CI- CONSTRUCTION -----
C
      CALL VALFM(LOADFM)
      I10 = LOADFM + 1
      I20=I10+(NWKS*(NWKS+1))/2
      I30=I20+ NWKS
      I40=I30+ NWKS
      I50=I40+ NWKS
      I60=I50+ NEMEMX
      I70=I60+ NEMEMX
      LAST=I70+ NEMEMX
      NEED = LAST - LOADFM - 1
      IF(MASWRK  .AND.  EXETYP.EQ.CHECK) WRITE(IW,9950) NEED,NGOT
      CALL GETFM(NEED)
      IF(EXETYP.EQ.CHECK) GO TO 700
C
      CALL DGHII(NFT15,NWKS,X(I20),X(I30),NEMEMX,X(I50),X(I60))
      CALL DGHIJ(NFT16,NWKS,X(I10),X(I30),NEMEMX,X(I50),X(I60),X(I70))
C
C     ----- SAVE -CI- MATRIX -----
C
      DO 40 I=1,NWKS
         II=(I*(I+1))/2
         X(II+I10-1)=X(I+I20-1)
   40 CONTINUE
C
C     ----- WRITE -CI- MATRIX AS RECORD -31- OF -IDAF- -----
C
      NDUM=(NWKS*(NWKS+1))/2
      CALL DAWRIT(IDAF,IODA,X(I10),NDUM,323,0)
      IF(OUT) CALL PRTRI(X(I10),NWKS)
C
  700 CONTINUE
      IF(SOME) WRITE(IW,9995)
      CALL RETFM(NEED)
      RETURN
C
 9000 CONTINUE
      IF (MASWRK) WRITE(IW,9994) ISET1,TITLA,TITLA1,
     *                           ISET2,TITLE,TITLE1
      CALL ABRT
      RETURN
C
 9999 FORMAT(/,10X,9(1H-),/,10X,'CI-MATRIX',/,10X,9(1H-))
 9997 FORMAT(/,' NAME OF THE ENERGY MATRIX ... ',10A8,
     *       /,' NAME OF THE -DRT- TABLE   ... ',10A8,
     *       /,' NUMBER OF CONFIGURATIONS  ... ',I10,
     *       /,1X,'INTEGRAL CUTOFF=',F20.10)
 9995 FORMAT(' ...... END OF CI-MATRIX CONSTRUCTION ......')
 9994 FORMAT(' -NFT15- AND -NFT16- DO NOT MATCH. ',/,
     *       ' NFT15 HEADER RECORD @ ISET1 = ',I5,2X,10A8,/
     *       ' NFT16 HEADER RECORD @ ISET2 = ',I5,2X,10A8)
 9950 FORMAT(/,' STATISTICS FOR FORMING HAMILTONIAN MATRIX:'/
     1         ' NUMBER OF WORDS NEEDED    = ',I10,
     2         ' NUMBER OF WORDS AVAILABLE = ',I10)
      END
C*MODULE MCQUD   *DECK DGHII
      SUBROUTINE DGHII(NFT,NWKS,H,INDX,NEMEMX,XX1,IX)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      LOGICAL GOPARR,DSKWRK,MASWRK
C
      PARAMETER (MXRT=100)
C
      DIMENSION H(*),INDX(NWKS),XX1(*),IX(*)
C
      COMMON /ENRGYS/ EN,EELCT,ETOT,SZ,SZZ,EC,ESCF,EERD,E1,E2,
     *                VEN,VEE,EPOT,EKIN,ECI(MXRT),STATN,EDFT(2)
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
C
      ENC = EC + EN
      IF (GOPARR) ENC = ENC/NPROC
      DO 100 I = 1,NWKS
         H(I) = ENC
  100 CONTINUE
      READ(NFT) INDX
      NUWK=0
      MX=0
*I64  MASK32 = 2**32 - 1
C
  120 CONTINUE
      CALL PXREAD(NFT,XX1,IX,MX,NEMEMX)
      NX=ABS(MX)
      IF(NX.EQ.0) GO TO 220
C
      DO 200 N = 1,NX
         VAL = XX1(N)
C
         NPACK = N
*I32     IPACK = IX( 2*NPACK - 1 )
*I32     JPACK = IX( 2*NPACK     )
*I64     LABEL = IX(NPACK)
*I64     IPACK = ISHFT( LABEL, -32 )
*I64     JPACK = IAND( LABEL, MASK32 )
         IUWK = IPACK
         NLWK = JPACK
C
         IF (IUWK.EQ.0) THEN
            NUWK = INT(VAL)
            GO TO 200
         END IF
C
         DO 180 I = 1,NLWK
            II = INDX(IUWK)
            DO 160 J = 1,NUWK
               H(II) = H(II)+VAL
               II = II+1
  160       CONTINUE
            IUWK = IUWK+1
  180    CONTINUE
  200 CONTINUE
      IF (MX .GT. 0) GO TO 120
C
  220 CONTINUE
C
      IF (GOPARR) CALL DDI_GSUMF(1850,H,NWKS)
C
      RETURN
      END
C*MODULE MCQUD   *DECK DGHIJ
      SUBROUTINE DGHIJ(NFT,NWKS,H,INDX,NEMEMX,XX1,IX,IXX)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      LOGICAL GOPARR,DSKWRK,MASWRK
C
      PARAMETER (MXAO=8192)
C
      COMMON /IJPAIR/ IA(MXAO)
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
C
      PARAMETER (ZERO=0.0D+00)
C
      DIMENSION H(*),INDX(*),XX1(NEMEMX),IX(NEMEMX),IXX(NEMEMX)
C     DIMENSION NLIST(*),IDLIST(*)
C
      DO 10 I=1,NWKS
   10 IA(I)=(I*(I-1))/2
      NH=NWKS*(NWKS+1)/2
      DO 30 I=1,NWKS
      MAX=I-1
      IF(MAX.EQ.0) GO TO 30
      DO 20 J=1,MAX
      IJ=IA(I)+J
   20 H(IJ)=ZERO
   30 CONTINUE
C
C     ----- OFF-DIAGONAL -CI- ELEMENTS -----
C
      CALL SEQREW(NFT)
      READ(NFT)
*I64  MASK32 = 2**32 - 1
C
  120 CONTINUE
         CALL PKREAD(NFT,XX1,IX,IXX,MX,NEMEMX)
         NX = ABS(MX)
         IF(NX.EQ.0) GO TO 190
C
         DO 180 N=1,NX
            VAL = XX1(N)
C
            NPACK = N
*I32        IPACK = IX( 2*NPACK - 1 )
*I32        JPACK = IX( 2*NPACK     )
*I64        LABEL = IX(NPACK)
*I64        IPACK = ISHFT( LABEL, -32 )
*I64        JPACK = IAND( LABEL, MASK32 )
            IUWK = IPACK
            JUWK = JPACK
C
            IF(IUWK.EQ.0) THEN
               NUWK = INT(VAL)
               GO TO 180
            END IF
C
*I32        NLWK = IXX(N)
*I64        IF (MOD(N,2) .EQ. 0) THEN
*I64          LABEL = IXX( N/2 )
*I64          NLWK  = ISHFT( LABEL, -32 )
*I64        ELSE
*I64          LABEL = IXX( N/2+1 )
*I64          NLWK  = IAND( LABEL, MASK32 )
*I64        END IF
C
            DO 160 I = 1,NLWK
               II = INDX(IUWK)
               JJ = INDX(JUWK)
               DO 150 J = 1,NUWK
                  IF(II.GT.NWKS) GO TO 155
                  IF(JJ.GT.NWKS) GO TO 155
                     NI = II
                     NJ = JJ
                     IF(NI.EQ.0) GO TO 140
                     IF(NJ.EQ.0) GO TO 140
                        IF(II.GE.JJ) GO TO 135
                           NI = NJ
                           NJ = II
  135                   CONTINUE
                        NIJ = (NI*NI-NI)/2 + NJ
                        H(NIJ) = H(NIJ) + VAL
  140                CONTINUE
                     II = II+1
                     JJ = JJ+1
  150          CONTINUE
  155          CONTINUE
               IUWK = IUWK+1
               JUWK = JUWK+1
  160       CONTINUE
  180    CONTINUE
      IF(MX.GT.0) GO TO 120
C
  190 CONTINUE
      IF (GOPARR) CALL DDI_GSUMF(1855,H,NH)
      RETURN
      END
C*MODULE MCQUD   *DECK QUDWTM
      SUBROUTINE QUDWTM(NFT,N,I,J,K,L,X,NX)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      DIMENSION X(NX)
C
      WRITE(NFT) N,I,J,K,L,X
      RETURN
      END
C*MODULE MCQUD   *DECK QUDHSS
      SUBROUTINE QUDHSS(IROT,IC,E,A,B,IA,H,
     *                  NROT,NROT1,NROT2,NIA,NORB,NORBS,NFTT,OUT)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      LOGICAL OUT
C
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
C
      DIMENSION E(NORBS,1),H(1),IROT(NORBS,1),IC(NORB,1),A(1),IA(1)
      DIMENSION B(1)
C
      DATA ZERO /0.0D+00/
C
      DO 10 I=1,NIA
   10 IA(I)=(I*(I-1))/2
      DO 20 I=1,NROT1
   20 B(I)=ZERO
      DO 30 I=1,NROT2
   30 A(I)=ZERO
C
C     ----- LAGRANGIAN CONTRIBUTION TO -A- -----
C
      DO 130 IX=1,NORB
      DO 130 IG=1,NORBS
      IB=IROT(IG,IX)
      IF(IB.EQ.0) GO TO 110
      VAL=E(IG,IX)
      IF(IG.GT.IX) VAL=-VAL
      B(IB)=B(IB)+(VAL+VAL)
  110 CONTINUE
C
      DO 120 K=1,NORB
      IF(IX.GT.NORBS.AND.K.GT.NORBS) GO TO 120
      VAL=E(IG,IX)
      IF(K.GT.IX) VAL=-VAL
      IF(K.LT.IG) VAL=-VAL
      IF(IX.GT.NORBS) IB=IROT(K,IX)
      IF(IX.LE.NORBS) IB=IROT(IX,K)
      JB=IROT(IG,K)
      IF(IB.EQ.0.OR.JB.EQ.0) GO TO 120
      I=MAX0(IB,JB)
      J=MIN0(IB,JB)
      IJ=IA(I)+J
      A(IJ)=A(IJ)+VAL
  120 CONTINUE
  130 CONTINUE
C
C     ----- OTHER CONTRIBUTIONS TO -A- -----
C
      CALL SEQREW(NFTT)
      DO 145 IX=2,NORB
      MAXIG=MIN0(IX-1,NORBS)
      DO 145 IG=1,MAXIG
      IIC=IC(IX,IG)
      IF(IIC.EQ.0) GO TO 145
      READ(NFTT) (H(IJH),IJH=1,IIC)
      IJH=0
      IIA=IROT(IG,IX)
      DO 140 JX=2,IX
      MAXJG=MIN0(JX-1,NORBS)
      IF(JX.EQ.IX) MAXJG=IG
      DO 140 JG=1,MAXJG
      JJC=IC(JX,JG)
      IF(JJC.EQ.0) GO TO 140
      IJH=IJH+1
      JJA=IROT(JG,JX)
      IF(IIA.EQ.0.OR.JJA.EQ.0) GO TO 140
      IJA=IA(IIA)+JJA
      A(IJA)=A(IJA)+H(IJH)
  140 CONTINUE
  145 CONTINUE
C
C     ----- MULTIPLY DIAGONAL ELEMENTS BY 2 -----
C
      DO 150 I=1,NROT
      II=IA(I)+I
  150 A(II)=A(II)+A(II)
C
C     ----- PRINT -A- AND -B- -----
C
      IF(.NOT.OUT) RETURN
      WRITE(IW,9999)
      DO 200 I=1,NROT
         N=IA(I)
         WRITE(IW,9998) (A(N+J),J=1,I)
         WRITE(IW,9997)  B(I)
  200 CONTINUE
      RETURN
 9999 FORMAT(40X,'   HESSIAN MATRIX   ',40X,' GRADIENT ',/)
 9998 FORMAT(5E20.12)
 9997 FORMAT(100X,E20.12,/)
      END
C*MODULE MCQUD   *DECK QUDCUP
      SUBROUTINE QUDCUP(NPRINT)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      DOUBLE PRECISION METHOD
      LOGICAL OUT,SOME,PACK2E,AVG
      LOGICAL GOPARR,DSKWRK,MASWRK
      LOGICAL CANONC,FCORE,FORS,NOCI,EKT,LINSER
C
      PARAMETER (MXNORO=250, MXRT=100)
C
      COMMON /CIFILS/ NFT11,NFT12,NFT13,NFT14,NFT15,NFT16,IDAF20,NEMEMX
      COMMON /FMCOM / X(1)
      COMMON /INTFIL/ NINTMX,NHEX,NTUPL,PACK2E,INTTYP,IGRDTYP
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
      COMMON /MCINP / METHOD,CISTEP,ACURCY,ENGTOL,DAMP,MICIT,NWORD,NORB,
     *                NOROT(2,MXNORO),MOFRZ(15),NPFLG(10),
     *                NOFO,MCFMO,CANONC,FCORE,FORS,NOCI,EKT,LINSER
      COMMON /MCPAR / NFZC,NCORBS,NCI,NORBS,NORB1,NUM
      COMMON /ORDOPT/ IDUM(6),NOSQUR
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
      COMMON /PCKLAB/ LABSIZ
      COMMON /RUNOPT/ RUNTYP,EXETYP,NEVALS,NGLEVL,NHLEVL
      COMMON /TM2AVG/ WAVG(MXRT),IAVG(MXRT),NAVG,AVG
C
      DATA CHECK/8HCHECK   /
C
      SOME=NPRINT.NE.-5 .AND. MASWRK
      OUT =NPRINT.GE.1 .AND. MASWRK
      IF(SOME) WRITE(IW,9999)
C
C     ----- GET CORE MEMORY -----
C
      CALL GOTFM(NGOT)
      IF(NWORD.GT.0) NGOT=MIN(NWORD,NGOT)
C
C     ----- READ IN DRT DATA -----
C
      CALL SEQREW(NFT12)
      READ(NFT12) NROTCI,NWKS
      ISTAT=IAVG(NAVG)
      IF(ISTAT.GT.1) THEN
         DO 30 I=1,ISTAT-1
            CALL SEQADV(NFT12)
   30    CONTINUE
      END IF
C
      IF(OUT) WRITE(IW,9987) NFZC,NCORBS,NCI,NORBS,NORB,NUM,NWKS
      IF(OUT) WRITE(IW,9985) NROTCI,NWKS
C
      N1= NUM
      N2= N1*(N1+1)/2
      N3= N1* N1
      NLGR= NORB*NORBS
      NLG2= NLGR*NWKS
C
      MI    = NORBS
      MIJ   = MI*(MI+1)/2
      NIA   = MAX0(NORB,MIJ)
      MIJKL = MIJ*(MIJ+1)/2
      MA    = NORB-NORBS
      MAJ   = MA*MI
      MAJKL = MAJ*MIJ
      MAB   = MA*(MA+1)/2
      MABKL = MAB*MIJ
      MAJBL = MAJ*(MAJ+1)/2
      MAXINT= MAX0(MIJKL,MAJKL,MABKL,MAJBL)
C
C     ----- MEMORY POINTERS  -------
C
C       ADRESS    LENGTH   COMMENT
C
C       AT I00    NORB    PAIR ROTATION INDECES
C       AT I10    NLG2    LAGRANGIAN MATRIX
C       AT I20    N2      CORE FOCK MATRIX OVER MO'S
C       AT I30    N2      CORE FOCK MATRIX OVER AO'S
C       AT I40    N3      MOLECULAR ORBITALS
C       AT I50    N1      TEMPORARY ARRAY FOR TFTRI
C
      CALL VALFM(LOADFM)
      I00= LOADFM + 1
      I10= I00 + NORB
      I20= I10 + NLG2
      I30= I20 + N2
      NEED1=I30-I00
      I40= I30 + N2
      I50= I40 + N3
      I60= I50 + N1
      NEED2=I60-I30
      CALL GETFM(NEED1)
      CALL GETFM(NEED2)
C
C     ------- INITIALIZE LAGRANGIAN  -------
C
      CALL VCLR(X(I10),1,NLG2)
      IF(EXETYP.EQ.CHECK) GO TO 200
C
C     ----- SET UP ORBITAL PAIR ROTATION INDICES -----
C
      CALL QUDPAR(X(I00),NORB,NORBS)
C
C     ----- READ ORBITALS AND CORE FOCK MATRIX IN ATOMIC BASIS -----
C
       CALL DAREAD(IDAF,IODA,X(I40),N3,15,0)
       CALL DAREAD(IDAF,IODA,X(I30),N2,35,0)
C
C     ----- TRANSFORM CORE FOCK MATRIX TO MOLECULAR BASIS -----
C
      IF (MASWRK) CALL TFTRI(X(I20),X(I30),X(I40),X(I50),NORB,NUM,NUM)
C
  200 CONTINUE
      CALL RETFM(NEED2)
C
C     ----- FORM 2E-CONTRIBUTION TO LAGRANGIAN MATRIX -----
C
      NFTG=NFT16
      NFTX=IJKT
C
C     ----- MEMORY POINTERS -------
C
C     ADRESS    LENGTH    COMMENT
C
C     AT I30    MAXINT    ORDERED INT2
C     AT I40    NWKS      TM1 AND TM2
C     AT I50    NINTMX    INT2 BUFFER
C     AT I60    NINTMX    INT2 LABELS
C     AT I70    NIA       PAIR INDEX ARRAY
C
      I40=I30+MAXINT
      I50=I40+NWKS
      I60=I50+NINTMX
      I70=I60+NINTMX*LABSIZ
      LAST=I70+NIA
      NEED2=LAST-I30
      IF(MASWRK .AND. EXETYP.EQ.CHECK) WRITE(IW,9992) NEED1+NEED2,NGOT
      CALL GETFM(NEED2)
      IF(EXETYP.EQ.CHECK) GO TO 400
C
C      ----- LAGRANGIAN CONTRIBUTION FROM CORE FOCK -----
C
      CALL QUDCOR(X(I00),X(I10),X(I20),X(I40),
     *            NORB,NCORBS,NWKS,NFT12,X(I70),NIA)
C
C     -----  LAGRANGIAN CONTRIBUTION FROM (IJ//KL) -----
C
      CALL SEQREW(NFTX)
      IF (MASWRK) READ(NFTX)
      CALL VCLR(X(I30),1,MAXINT)
      CALL M1EGHR(NFTX,X(I30),X(I50),X(I60),X(I70),NINTMX,MI,NOSQUR)
      CALL QUDEGH(X(I00),X(I10),X(I20),X(I30),X(I40),X(I70),
     *            MI,NCORBS,NWKS,NFTG)
C
      IF(MA.EQ.0) GO TO 400
C
C     -----  LAGRANGIAN CONTRIBUTION FROM (AJ//KL) -----
C
      CALL SEQREW(NFTX)
      IF (MASWRK) READ(NFTX)
      CALL M1LGRR(NFTX,X(I30),X(I50),X(I60),X(I70),
     *            NINTMX,MI,MA,MIJ,NOSQUR)
      CALL QUDLGR(X(I00),X(I10),X(I20),X(I30),X(I40),X(I70),
     *            MI,MA,MIJ,NCORBS,NWKS,NFTG)
C
  400 CONTINUE
C
C     ----- COMPLETE THE LAGRANGIAN MATRIX -----
C
      IF (GOPARR) CALL DDI_GSUMF(2000,X(I10),NLG2)
C
C     ----- WRITE OUT LAGRANGIAN MATRIX -----
C
      IF(NROTCI.GT.ISTAT) THEN
         DO 40 I=ISTAT+1,NROTCI
   40       CALL SEQADV(NFT12)
      END IF
      IF(SOME) WRITE(IW,9996)
      CALL QUDWRC(X(I10),NORB,NORBS,NWKS,X(I20),NFT12,OUT)
C
      IF(SOME) WRITE(IW,9986)
      IF(SOME) CALL TEXIT(2,2)
C
C     ----- RESET CORE MEMORY -----
C
      CALL RETFM(NEED2)
      CALL RETFM(NEED1)
      RETURN
C
 9999 FORMAT(/,10X,22(1H-),/,10X,'MO-CI COUPLING VECTORS',
     1       /,10X,22(1H-))
 9996 FORMAT(/,10X,' MO-CI COUPLING VECTOR')
 9992 FORMAT(/,' STATISTICS FOR FORMING LAGRANGIAN MATRIX IN CORE:'/
     1         ' NUMBER OF WORDS NEEDED    = ',I10,
     2         ' NUMBER OF WORDS AVAILABLE = ',I10)
 9987 FORMAT(/,' IN -QUDCUP- , NFZC,NCORBS,NCI,NORBS,NORB,NUM,NWKS =',
     1       6I4,I8)
 9986 FORMAT(' ...... END OF MO-CI COUPLING CALCULATION ..... ')
 9985 FORMAT(/,' IN -QUDCUP- , NROTCI,NWKS = ',I3,I8)
      END
C*MODULE MCQUD   *DECK QUDPAR
      SUBROUTINE QUDPAR(IB,NORB,NORBS)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      DIMENSION IB(NORB)
C
      DO 10 I =1,NORB
   10    IB(I)=(I-1)*NORBS
      RETURN
      END
C*MODULE MCQUD   *DECK QUDCOR
      SUBROUTINE QUDCOR(IB,EG,F,XCI,NORB,NCORBS,NWKS,NFT12,IA,NIA)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      LOGICAL GOPARR,DSKWRK,MASWRK
C
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
C
      DIMENSION IB(NORB),EG(NWKS,1),F(1),IA(1),XCI(NWKS)
C
      DO 10 I=1,NIA
   10    IA(I)=(I*(I-1))/2
      READ(NFT12) XCI
      IF(NCORBS.EQ.0) RETURN
      IF(.NOT.MASWRK) RETURN
C
      DO 20 M=1,NORB
      DO 20 I=1,NCORBS
         IBMI=IB(M)+ I
         MM = M
         II = I
         MI = IA(MAX0(MM,II)) + MIN0(MM,II)
         FCONT = F(MI)+F(MI)
         CALL DAXPY(NWKS,FCONT,XCI,1,EG(1,IBMI),1)
   20    CONTINUE
      RETURN
      END
C*MODULE MCQUD   *DECK QUDEGH
      SUBROUTINE QUDEGH(IB,EG,F,X,GIJKL,IA,NORBS,NCORBS,NWKS,NFT16)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      LOGICAL IANDJ,KANDL,SAME
      LOGICAL GOPARR,DSKWRK,MASWRK
C
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
C
      DIMENSION IA(1),IB(1),EG(NWKS,1)
      DIMENSION X(1),GIJKL(NWKS),F(1)
C
      DATA ZERO /0.0D+00/
      DATA TWO  /2.0D+00/
C
      NIA = NORBS*(NORBS+1)/2
      DO 10 I=1,NIA
   10    IA(I)=(I*(I-1))/2
C
      CALL SEQREW(NFT16)
C
  100 READ(NFT16,END=2000) ITM,II,JJ,KK,LL,GIJKL
C
C     ----- IF ITM = 1 WE HAVE -TM1- , IF ITM=2 WE HAVE -TM2- -----
C
      IF(ITM.EQ.2) GO TO 300
C
C     ----- PROCESS -TM1- FIRST -----
C
      KG = KK + NCORBS
      LG = LL + NCORBS
      KL = IA(MAX0(KG,LG)) + MIN0(KG,LG)
C
      DO 220 M=1,NORBS
      IF(NCORBS.EQ.0) GO TO 210
C
      DO 200 I=1,NCORBS
         MI=IA(MAX0(M,I))+MIN0(M,I)
         MIKL=IA(MAX0(MI,KL))+MIN0(MI,KL)
         MK=IA(MAX0(M,KG))+MIN0(M,KG)
         IL=IA(MAX0(I,LG))+MIN0(I,LG)
         MKIL=IA(MAX0(MK,IL))+MIN0(MK,IL)
         ML=IA(MAX0(M,LG))+MIN0(M,LG)
         IK=IA(MAX0(I,KG))+MIN0(I,KG)
         MLIK=IA(MAX0(ML,IK))+MIN0(ML,IK)
         XCONT = (X(MIKL)+X(MIKL))+ (X(MIKL)+X(MIKL))
     *                            - X(MKIL) - X(MLIK)
         IF(KK.EQ.LL) XCONT = XCONT / TWO
         IBMI=IB(M)+I
         CALL DAXPY(NWKS,XCONT,GIJKL,1,EG(1,IBMI),1)
  200 CONTINUE
C
  210 CONTINUE
C
      IF (.NOT.MASWRK) GO TO 220
      MM = M
      K1 = KG
      L1 = LG
      IBMK=IB(M) + KG
      ML = IA(MAX0(MM,L1)) + MIN0(MM,L1)
      FCONT=F(ML)
      CALL DAXPY(NWKS,FCONT,GIJKL,1,EG(1,IBMK),1)
      IF(K1.EQ.L1) GO TO 220
      IBML=IB(M) + LG
      MK = IA(MAX0(MM,K1)) + MIN0(MM,K1)
      FCONT=F(MK)
      CALL DAXPY(NWKS,FCONT,GIJKL,1,EG(1,IBML),1)
C
  220 CONTINUE
C
      GO TO 100
C
C     ----- PROCESS -TM2- NOW -----
C
  300 CONTINUE
C
      IGIJKL = II + NCORBS
      JGIJKL = JJ + NCORBS
      KGIJKL = KK + NCORBS
      LGIJKL = LL + NCORBS
C
      IANDJ=II.EQ.JJ
      KANDL=KK.EQ.LL
      SAME =II.EQ.KK.AND.JJ.EQ.LL
      NGAB=1
      NGCD=1
C
      IF(.NOT.IANDJ) NGAB=2
      NAB=IA(IGIJKL)+JGIJKL
      IF(.NOT.KANDL) NGCD=2
      NCD=IA(KGIJKL)+LGIJKL
C
      IG=IGIJKL
      JG=JGIJKL
      KG=KGIJKL
      LG=LGIJKL
      KL=NCD
      NGIJKL=NGCD
      NPASS=1
C
 1000 CONTINUE
      DO 1100 M=1,NORBS
         MM=MAX0(M,JG)
         JJ=MIN0(M,JG)
         MJ=IA(MM)+JJ
         MJX=MAX0(MJ,KL)
         KLX=MIN0(MJ,KL)
         MJKLX=IA(MJX)+KLX
         XMJKL=X(MJKLX)
         IF(XMJKL.EQ.ZERO) GO TO 1100
         MILGR=IB(M)+IG
         IF(NGIJKL.EQ.2)XMJKL=XMJKL+XMJKL
         CALL DAXPY(NWKS,XMJKL,GIJKL,1,EG(1,MILGR),1)
 1100 CONTINUE
C
 1200 NPASS=NPASS+1
      GO TO (1300,1210,1220,1230,1300),NPASS
 1210 IF(IANDJ) GO TO 1200
      IG=JGIJKL
      JG=IGIJKL
      GO TO 1000
 1220 IF(SAME) GO TO 1300
      IG=KGIJKL
      JG=LGIJKL
      KG=IGIJKL
      LG=JGIJKL
      KL=NAB
      NGIJKL=NGAB
      GO TO 1000
 1230 IF(KANDL) GO TO 1300
      IG=LGIJKL
      JG=KGIJKL
      GO TO 1000
C
 1300 CONTINUE
C
      GO TO 100
C
 2000 CONTINUE
C
      RETURN
      END
C*MODULE MCQUD   *DECK QUDLGR
      SUBROUTINE QUDLGR(IB,EG,F,X,GIJKL,IA,MI,MA,MIJ,NCORBS,NWKS,NFT16)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      LOGICAL IANDJ,KANDL,SAME
      LOGICAL GOPARR,DSKWRK,MASWRK
C
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
C
      DIMENSION EG(NWKS,1),F(1)
      DIMENSION IA(1),IB(1)
      DIMENSION X(MA,1),GIJKL(NWKS)
C
      DATA ZERO /0.0D+00/
      DATA TWO  /2.0D+00/
C
      NORB= MA+MI
      NIA = MAX0(MIJ,NORB)
      DO 10 I=1,NIA
   10    IA(I)=(I*(I-1))/2
C
      CALL SEQREW(NFT16)
C
  100 READ(NFT16,END=2000) ITM,II,JJ,KK,LL,GIJKL
C
C     ----- IF ITM = 1 WE HAVE -TM1- , IF ITM=2 WE HAVE -TM2- -----
C
      IF(ITM.EQ.2) GO TO 300
C
C     ----- PROCESS -TM1- FIRST -----
C
      KG = KK + NCORBS
      LG = LL + NCORBS
      KL = IA(MAX0(KG,LG)) + MIN0(KG,LG)
C
      DO 220 M=1,MA
      IF(NCORBS.EQ.0) GO TO 210
C
      DO 200 I=1,NCORBS
      IKL=(I-1)*MIJ +KL
      IL=IA(LG)+ I
      KIL=(KG-1)*MIJ +IL
      IK=IA(KG)+ I
      LIK=(LG-1)*MIJ +IK
      XCONT = (X(M,IKL)+X(M,IKL))+(X(M,IKL)+X(M,IKL))-X(M,KIL)-X(M,LIK)
      IF(KK.EQ.LL) XCONT = XCONT / TWO
      IBMI=IB(M+MI)+I
      CALL DAXPY(NWKS,XCONT,GIJKL,1,EG(1,IBMI),1)
  200 CONTINUE
C
  210 CONTINUE
C
      IF (.NOT.MASWRK) GO TO 220
      MM = M + MI
      K1 = KG
      L1 = LG
      IBMK=IB(M+MI) + KG
      ML = IA(MAX0(MM,L1)) + MIN0(MM,L1)
      FCONT=F(ML)
      CALL DAXPY(NWKS,FCONT,GIJKL,1,EG(1,IBMK),1)
      IF(K1.EQ.L1) GO TO 220
      IBML=IB(M+MI) + LG
      MK = IA(MAX0(MM,K1)) + MIN0(MM,K1)
      FCONT=F(MK)
      CALL DAXPY(NWKS,FCONT,GIJKL,1,EG(1,IBML),1)
C
  220 CONTINUE
C
      GO TO 100
C
C     ----- PROCESS -TM2- NOW -----
C
  300 CONTINUE
C
      IGIJKL = II + NCORBS
      JGIJKL = JJ + NCORBS
      KGIJKL = KK + NCORBS
      LGIJKL = LL + NCORBS
C
      IANDJ=II.EQ.JJ
      KANDL=KK.EQ.LL
      SAME =II.EQ.KK.AND.JJ.EQ.LL
      NGAB=1
      NGCD=1
C
      IF(.NOT.IANDJ) NGAB=2
      NAB=IA(IGIJKL)+JGIJKL
      IF(.NOT.KANDL) NGCD=2
      NCD=IA(KGIJKL)+LGIJKL
C
      IG=IGIJKL
      JG=JGIJKL
      KG=KGIJKL
      LG=LGIJKL
      KL=NCD
      NGIJKL=NGCD
      NPASS=1
C
 1000 CONTINUE
      DO 1100 M=1,MA
         JKLX=(JG-1)*MIJ+KL
         XMJKL=X(M,JKLX)
         IF(XMJKL.EQ.ZERO) GO TO 1100
         MILGR=IB(M+MI)+IG
         IF(NGIJKL.EQ.2)XMJKL=XMJKL+XMJKL
         CALL DAXPY(NWKS,XMJKL,GIJKL,1,EG(1,MILGR),1)
 1100 CONTINUE
C
 1200 NPASS=NPASS+1
      GO TO (1300,1210,1220,1230,1300),NPASS
 1210 IF(IANDJ) GO TO 1200
      IG=JGIJKL
      JG=IGIJKL
      GO TO 1000
 1220 IF(SAME) GO TO 1300
      IG=KGIJKL
      JG=LGIJKL
      KG=IGIJKL
      LG=JGIJKL
      KL=NAB
      NGIJKL=NGAB
      GO TO 1000
 1230 IF(KANDL) GO TO 1300
      IG=LGIJKL
      JG=KGIJKL
      GO TO 1000
 1300 CONTINUE
C
      GO TO 100
C
 2000 CONTINUE
C
      RETURN
      END
C*MODULE MCQUD   *DECK QUDWRC
      SUBROUTINE QUDWRC(EG,NORB,NORBS,NWKS,X,NFT12,OUT)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      LOGICAL OUT
C
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
C
      DIMENSION EG(NWKS,NORBS,NORB)
      DIMENSION X(NORBS,NORB)
C
      DO 20 IWKS=1,NWKS
         DO 10 IG=1,NORBS
         DO 10 IX=1,NORB
            X(IG,IX) = EG(IWKS,IG,IX)
   10    CONTINUE
         WRITE(NFT12) X
         IF(OUT) THEN
            WRITE(IW,9999) IWKS
            CALL LGROUT(X,NORBS,NORB)
         END IF
   20    CONTINUE
      RETURN
 9999 FORMAT(/,' IN -QUDWRC-, WRITE MO-CI COUPLING FOR -IWKS- =',I8)
      END
C*MODULE MCQUD   *DECK QUDEQS
      SUBROUTINE QUDEQS(NPRINT)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C
      DOUBLE PRECISION METHOD
      LOGICAL SOME,OUT,AVG
      LOGICAL GOPARR,DSKWRK,MASWRK
      LOGICAL CANONC,FCORE,FORS,NOCI,EKT,LINSER,NOTRF
C
      PARAMETER (MXNORO=250, MXRT=100)
C
      COMMON /CIFILS/ NFT11,NFT12,NFT13,NFT14,NFT15,NFT16,IDAF20,NEMEMX
      COMMON /FMCOM / X(1)
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
      COMMON /MACHIN/ NWDVAR,MAXFM,MAXSM,LIMFM,LIMSM
      COMMON /MCINP / METHOD,CISTEP,ACURCY,ENGTOL,DAMP,MICIT,NWORD,NORB,
     *                NOROT(2,MXNORO),MOFRZ(15),NPFLG(10),
     *                NOFO,MCFMO,CANONC,FCORE,FORS,NOCI,EKT,LINSER
      COMMON /MCPAR / NFZC,NCORBS,NCI,NORBS,NORB1,NUM
      COMMON /NTNOPT/ ENERGY,ENERG0,DEMAX,SQCDF,ITER,MICRO,NOTRF
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
      COMMON /RUNOPT/ RUNTYP,EXETYP,NEVALS,NGLEVL,NHLEVL
      COMMON /TM2AVG/ WAVG(MXRT),IAVG(MXRT),NAVG,AVG
C
      EXTERNAL SPARAX1,SPARXA1
C
      DIMENSION IX(1)
      EQUIVALENCE(X(1),IX(1))
C
      DATA ZERO   /0.0D+00/
      DATA CHECK/8HCHECK   /
C
      SOME=NPRINT.NE.-5 .AND. MASWRK
      OUT =NPRINT.GE.1 .AND. MASWRK
      IF(SOME) WRITE(IW,9999)
C
C     ----- GET CORE MEMORY -----
C
      CALL GOTFM(NGOT)
      IF(NWORD.GT.0) NGOT=MIN(NWORD,NGOT)
C
      M1=NORB
      M3=M1 * M1
      N1=NUM
      N3=N1 * N1
      NLGR=NORB*NORBS
C
C     ----- READ FIRST RECORD ON FT11 AND FT12 -----
C
      CALL SEQREW(NFT11)
      READ(NFT11)
C
      CALL SEQREW(NFT12)
      READ(NFT12) NSTAT,NWKS
C
      IF(OUT) WRITE(IW,9981) NFZC,NCORBS,NCI,NORBS,NORB,NUM,NWKS
      IF(OUT) WRITE(IW,9983)
C
C     ----- READ ORBITAL ROTATION TABLE FROM -IDAF- -----
C
      CALL VALFM(LOADFM)
      I10= LOADFM + 1
      LAST = I10 + NLGR
      NEED = LAST -LOADFM - 1
      CALL GETFM(NEED)
C
      CALL DAREAD(IDAF,IODA,X(I10),NLGR,322,1)
C
      NROT=0
      DO 10 I=1,NLGR
         IROT=IX((I10-1)*NWDVAR+I)
         IF(IROT.GT.NROT) NROT=IROT
   10    CONTINUE
      IF(OUT) WRITE(IW,9987)
C
C     ----- CHECK FOR ENOUGH MEMORY -----
C
      NEQS  = NROT+NWKS
      NEQS2 = (NEQS*(NEQS+1))/2
      NWKS2 = (NWKS*(NWKS+1))/2
      NROT2 = (NROT*(NROT+1))/2
      LAST1 = NLGR+NEQS+NEQS2+NWKS+NWKS2+NEQS
      LAST2 = NLGR+NEQS+NEQS2+6*NEQS
      LAST3 = NLGR+NEQS+NLGR+M3+N3+NORB
      LAST  = MAX0(LAST1,LAST2,LAST3)
      IF(LAST.GT.NGOT) GO TO 9010
      IF(SOME) WRITE(IW,9992) LAST,NGOT
      IF(MASWRK .AND. EXETYP.EQ.CHECK) WRITE(IW,9992) LAST,NGOT
C
      I20=I10+ NLGR
      I30=I20+ NEQS
      I40=I30+ NEQS2
      NEED1=I30-I20
      CALL GETFM(NEED1)
      CALL GETFM(NEQS2)
      I50=I40+ NWKS
      I60=I50+ NWKS2
      I70=I60+ NEQS
      NEED2=I70-I40
      CALL GETFM(NEED2)
C
      IF(EXETYP.EQ.CHECK) GO TO 50
C
C     ----- READ -MO- PART OF NEWTON-RAPHSON EQUATIONS -----
C
      DO 20 I=1,NEQS2
   20    X(I+I30-1)=ZERO
C
      CALL SEQREW(23)
      CALL SQREAD(23,X(I20),NROT)
      CALL SQREAD(23,X(I30),NROT2)
      IF(OUT) WRITE(IW,9986)
C
C     ----- READ -CI- PART OF NEWTON-RAPHSON EQUATIONS -----
C
      IF(NAVG.NE.1) GO TO 9000
C
      ISTAT=IAVG(NAVG)
      IF(OUT) WRITE(IW,9980) ISTAT
      IF(ISTAT.GT.1) THEN
         DO 30 I=1,ISTAT-1
   30       CALL SEQADV(NFT12)
      END IF
      READ(NFT12) (X(I+I40-1),I=1,NWKS)
      IF(NSTAT.GT.ISTAT) THEN
         DO 40 I=ISTAT+1,NSTAT
   40       CALL SEQADV(NFT12)
      END IF
      IF(OUT) WRITE(IW,9982) NSTAT
C
C     ----- FORM MO-CI COUPLING TERMS -----
C
      CALL QUDCIC(NORB,NORBS,NEQS,NROT,NWKS,
     *            X(I10),X(I50),X(I30),X(I40),NFT12,OUT)
      IF(OUT) WRITE(IW,9985)
C
C     ----- COMPLETE NEWTON-RAPHSON EQUATIONS -----
C
      CALL DAREAD(IDAF,IODA,X(I50),NWKS2,323,0)
      CALL QUDABC(NEQS,NROT,NWKS,X(I30),X(I20),X(I50),X(I40),X(I60),
     *            ENERGY,OUT)
      IF(OUT) WRITE(IW,9984)
C
   50 CONTINUE
      CALL RETFM(NEED2)
C
C     ----- SOLVE LINEAR EQUATIONS -----
C
      I50=I40+NEQS
      I60=I50+NEQS
      I70=I60+NEQS
      I80=I70+NEQS
      I90=I80+NEQS
      LAST=I90+NEQS
      NEED2=LAST-I40
      CALL GETFM(NEED2)
      IF(EXETYP.EQ.CHECK) GO TO 60
C
      CALL SPMRSQ(X(I30),X(I20),X(I40),NEQS,
     *            X(I50),X(I60),X(I70),X(I80),IW,X(I90))
      IF(OUT) THEN
         WRITE(IW,9997)
         WRITE(IW,9996) (X(I+I20-1),I=1,NEQS)
      END IF
C
C     ----- GET NEW -CI- VECTOR -----
C
      CALL QUDVCI(NFT12,ISTAT,NSTAT,NWKS,
     *            X(I30),X(I20+NROT),ENERGY,OUT)
C
   60 CONTINUE
      CALL RETFM(NEED2)
      CALL RETFM(NEQS2)
C
C     ----- GET -SQCDF- AND NEW MO'S -----
C
      I40=I30+NLGR
      I50=I40+N3
      I60=I50+M3
      I70=I60+NORB
      NEED2=I70-I30
      CALL GETFM(NEED2)
      IF(EXETYP.EQ.CHECK) GO TO 70
C
      CALL QUDEXP(NORB,NORBS,NROT,SQCDF,X(I10),X(I20),X(I30),OUT)
      CALL NTNTRF(NUM,NORB,NORBS,NCORBS,FCORE,DAMP,
     *             X(I30),X(I40),X(I50),X(I60))
      IF(OUT) THEN
         WRITE(IW,9995)
         CALL PRSQL(X(I40),NORBS,NUM,NUM)
      END IF
C
      IF(SOME) THEN
         WRITE(IW,9994)
         CALL TEXIT(2,2)
      END IF
C
C     ----- RESET CORE MEMORY -----
C
   70 CONTINUE
      CALL RETFM(NEED2)
      CALL RETFM(NEED1)
      CALL RETFM(NEED)
      RETURN
C
 9000 CONTINUE
      IF (MASWRK) WRITE(IW,9998)
      CALL ABRT
      RETURN
 9010 CONTINUE
      MORE=LAST-NGOT
      IF (MASWRK) WRITE(IW,9988) MORE,NGOT,LAST
      CALL ABRT
      RETURN
C
 9999 FORMAT(/,10X,30(1H-),/,10X,'QUADRATIC ORBITAL OPTIMIZATION',
     1       /,10X,30(1H-))
 9998 FORMAT(' QUADRATIC -MCHF- PROGRAM CAN NOT HANDLE STATE',
     1       ' AVERAGING. STOP ')
 9997 FORMAT(/,' SOLUTION OF NEWTON-RAPHSON EQUATIONS ',/)
 9996 FORMAT(F20.12)
 9995 FORMAT(/,' NEW MOLECULAR ORBITALS')
 9994 FORMAT(' ...... END OF QUADRATIC ORBITAL OPTIMIZATION ...... ')
 9992 FORMAT(/,' STATISTICS FOR QUADRATIC ORBITAL OPTIMIZATION: '/
     1        ' NUMBER OF WORDS NEEDED    = ',I10,
     2        ' NUMBER OF WORDS AVAILABLE = ',I10)
 9988 FORMAT(' NEED MORE CORE TO SOLVE QUADRATIC EQUATIONS.',/,
     1       ' INCREASE CORE BY AT LEAST ',I10,
     2       ' WORDS. NEED,NGOT = ',2I10)
 9987 FORMAT(' ORBITAL ROTATION TABLE READ FROM -IDAF- (328)')
 9986 FORMAT(' -MO- PART OF N-R EQUATIONS READ FROM SEQUENTIAL FILE 23')
 9985 FORMAT(' -MO-CI- COUPLING TERMS READ FROM -NFT12- ')
 9984 FORMAT(' N-R EQUATIONS NOW COMPLETE ')
 9983 FORMAT(' IN -QUDEQS-   1 RECORD(S) READ FROM -NFT12- ')
 9982 FORMAT(' IN -QUDEQS-',I3,' RECORD(S) READ FROM -NFT12- ')
 9981 FORMAT(' IN -QUDEQS- , NFZC,NCORBS,NCI,NORBS,NORB,NUM,NWKS =',
     1       6I4,I8)
 9980 FORMAT(' IN -QUDEQS- , ISTAT = ',I4)
      END
C*MODULE MCQUD   *DECK QUDCIC
      SUBROUTINE QUDCIC(NORB,NORBS,NEQS,NROT,NWKS,
     *                  IROT,G,C,V,NFT,OUT)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      LOGICAL OUT
C
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
C
      DIMENSION IROT(NORBS,NORB),G(NORBS,NORB),C(1),V(1)
C
      DATA ZERO /0.0D+00/
C
      IF(OUT) THEN
         WRITE(IW,9997) NWKS,NORBS,NORB
         CALL PRTRI(C,NEQS)
      END IF
      DO 50 IWKS=1,NWKS
         IEQS=NROT+IWKS
         IC=(IEQS*(IEQS-1))/2
         READ(NFT) G
         IF(OUT) THEN
            CALL LGROUT(G,NORBS,NORB)
         END IF
         DO 40 IX=1,NORB
            DO 30 IG=1,NORBS
               IB=IROT(IG,IX)
               IF(IB.EQ.0) GO TO 30
               VAL=G(IG,IX)
               IF(IG.GT.IX) VAL=-VAL
               C(IC+IB)=C(IC+IB)+(VAL+VAL)
   30          CONTINUE
   40       CONTINUE
   50    CONTINUE
C
      IF(OUT) THEN
         WRITE(IW,9999)
         CALL PRTRI(C,NEQS)
         WRITE(IW,9996)
         DO 60 IX=1,NORB
   60       WRITE(IW,9995) (IROT(IG,IX),IG=1,NORBS)
         DO 80 IB=1,NROT
            DUM=ZERO
            DO 70 IWKS=1,NWKS
               IEQS=NROT+IWKS
               IC=(IEQS*(IEQS-1))/2
   70          DUM=DUM+C(IC+IB)*V(IWKS)
   80       WRITE(IW,9998) IB,DUM
      END IF
C
      RETURN
 9999 FORMAT(/,' PARTIAL HESSIAN MATRIX ')
 9998 FORMAT(' IN -QUDCIC- CHECK -IROT- AND -G- = ',I10,F20.12)
 9997 FORMAT(' IN -QUDCIC- READ MO-CI COUPLING,',
     1       ' FOR IWKS,NORBS,NORB =',I8,2I4)
 9996 FORMAT(' IN -QUDCIC- THE ORBITAL ROTATION TABLE IS = ')
 9995 FORMAT(25I4)
      END
C*MODULE MCQUD   *DECK QUDABC
      SUBROUTINE QUDABC(NEQS,NROT,NWKS,HQD,GQD,HCI,VCI,IA,E0,OUT)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      LOGICAL OUT
C
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
C
      DIMENSION HQD(1),GQD(1),HCI(1),VCI(1),IA(1)
C
      DATA ZERO /0.0D+00/
C
      DO 10 I=1,NEQS
   10 IA(I)=(I*(I-1))/2
C
C     ----- CALCULATE -E0- -----
C
      E0=ZERO
      DO 30 I=1,NWKS
      DO 20 J=1,I
      IJ=IA(I)+J
      HIJ=HCI(IJ)
      DUM=VCI(I)*HIJ*VCI(J)
   20 E0=E0+(DUM+DUM)
   30 E0=E0- DUM
      IF(OUT) WRITE(IW,9999)
      IF(OUT) WRITE(IW,9998) (VCI(I),I=1,NWKS)
      IF(OUT) WRITE(IW,9997)
      IF(OUT) CALL PRTRI(HCI,NWKS)
      IF(OUT) WRITE(IW,9996) E0
C
C     ----- CREATE - ( H - E0 ) -    -----
C
      DO 40 I=1,NWKS
      II=IA(I)+I
   40 HCI(II)=HCI(II)-E0
      IF(OUT) WRITE(IW,9995)
      IF(OUT) CALL PRTRI(HCI,NWKS)
C
C     ----- CALCULATE - ( H - E0 ) * C0 -    -----
C           IT IS THE -CI- PART OF THE GRADIENT
C
      DO 60 I=1,NWKS
      DUM=ZERO
      DO 50 J=1,NWKS
      II=MAX0(I,J)
      JJ=MIN0(I,J)
      IJ=IA(II)+JJ
   50 DUM=DUM+HCI(IJ)*VCI(J)
   60 GQD(NROT+I)=DUM
      IF(OUT) WRITE(IW,9994)
      IF(OUT) WRITE(IW,9998) (GQD(NROT+I),I=1,NWKS)
C
C     ----- USE IT TO CONSTRUCT THE -CI- PART OF HESSIAN -----
C           ( 1 - P ) * ( H - E0 ) * ( 1 - P ) =
C           ( H - E0) - P * ( H - E0 ) - ( H - E0 ) * P
C     NOTE THAT THE LAST TERM  P * ( H - E0 ) * P IS ZERO
C     FINALLY ADD -P- CONTRIBUTION WHERE P = C0 * C0 -----
C
      DO 70 I=1,NWKS
      DO 70 J=1,I
      IJQD=IA(NROT+I)+(NROT+J)
C
      IJCI=IA(I)+J
      HQD(IJQD)=HCI(IJCI)
C
      DUM=GQD(NROT+I)*VCI(J)
      HQD(IJQD)=HQD(IJQD)-DUM
C
      DUM=VCI(I)*GQD(NROT+J)
      HQD(IJQD)=HQD(IJQD)-DUM
C
      HQD(IJQD)=HQD(IJQD)+HQD(IJQD)
   70 CONTINUE
C
C     ----- COMPLETE MO-CI COUPLING TERMS -----
C
      DO 80 I=1,NWKS
      DO 80 J=1,NROT
      IJQD=IA(NROT+I)+J
   80 HQD(IJQD)=HQD(IJQD)+HQD(IJQD)
      DO 90 J=1,NROT
      DO 90 I=1,NWKS
      IJQD=IA(NROT+I)+J
      DUM=GQD(J)*VCI(I)
      HQD(IJQD)=HQD(IJQD)-(DUM+DUM)
   90 CONTINUE
      IF(OUT) WRITE(IW,9993)
      IF(OUT) CALL PRTRI(HQD,NEQS)
C
C     ----- COMPLETE THE GRADIENT -----
C
      DO 100 I=1,NEQS
  100 GQD(I)=-GQD(I)
      DO 110 IWKS=1,NWKS
  110 GQD(NROT+IWKS)=GQD(NROT+IWKS)+GQD(NROT+IWKS)
      IF(OUT) WRITE(IW,9992)
      IF(OUT) WRITE(IW,9998) (GQD(I),I=1,NEQS)
      RETURN
 9999 FORMAT(/,' -C0- VECTOR ',/)
 9998 FORMAT(F20.12)
 9997 FORMAT(/,' -CI- MATRIX ')
 9996 FORMAT(' ENERGY FOR THE -C0- VECTOR, E0 = ',F20.12)
 9995 FORMAT(/,' ( H - E0 ) MATRIX ')
 9994 FORMAT(/,25H ( H - E0 ) ' C0 GRADIENT,/)
 9993 FORMAT(/,' HESSIAN MATRIX ')
 9992 FORMAT(/,' GRADIENT VECTOR ',/)
      END
C*MODULE MCQUD   *DECK QUDVCI
      SUBROUTINE QUDVCI(NFT12,ISTAT,NSTAT,NWKS,V,X,E,OUT)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      LOGICAL OUT
C
      PARAMETER (MXRT=100)
C
      COMMON /ENRGYS/ EN,EELCT,ETOT,SZ,SZZ,EC,ESCF,EERD,E1,E2,
     *                VEN,VEE,EPOT,EKIN,ECI(MXRT),STATN,EDFT(2)
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
C
      DIMENSION V(NWKS,1),X(1)
C
      DATA ZERO,PT5,ONE /0.0D+00,0.5D+00,1.0D+00/
C
      CALL SEQREW(NFT12)
      READ(NFT12)
      DO 10 MSTAT=1,NSTAT
   10 READ(NFT12) (V(I,MSTAT),I=1,NWKS)
C
      DUM=ZERO
      DO 20 I=1,NWKS
   20 DUM=DUM+V(I,ISTAT)*X(I)
      DO 30 I=1,NWKS
   30 X(I)=X(I)-V(I,ISTAT)*DUM
      DELSQ=ZERO
      DO 40 I=1,NWKS
   40 DELSQ=DELSQ+X(I)*X(I)
      DELSQ=DELSQ*PT5
C
      DUM=ZERO
      DO 50 I=1,NWKS
      V(I,ISTAT)=V(I,ISTAT)*(ONE-DELSQ)+X(I)
      DUM=DUM+V(I,ISTAT)*V(I,ISTAT)
   50 CONTINUE
C
C     ----- NORMALIZE NEW -CI- VECTOR -----
C
      DUM=ONE/ SQRT(DUM)
      DO 60 I=1,NWKS
   60 V(I,ISTAT)=V(I,ISTAT)*DUM
C
      DO 80 MSTAT=1,NSTAT
      IF(MSTAT.EQ.ISTAT) GO TO 80
      ECI(MSTAT)=ZERO
      DO 70 IWKS=1,NWKS
   70 V(IWKS,MSTAT)=ZERO
   80 CONTINUE
      ECI(ISTAT)=E
C
      CALL SEQREW(NFT12)
      READ(NFT12)
      DO 90 MSTAT=1,NSTAT
   90 WRITE(NFT12) (V(I,MSTAT),I=1,NWKS)
C
      IF(OUT) WRITE(IW,9999)
      IF(OUT) WRITE(IW,9998) (V(I,ISTAT),I=1,NWKS)
      RETURN
 9999 FORMAT(/,' NEW -CI- VECTOR ',/)
 9998 FORMAT(6F20.12)
      END
C*MODULE MCQUD   *DECK QUDEXP
      SUBROUTINE QUDEXP(NORB,NORBS,NROT,SQCDF,IROT,X,S,OUT)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      LOGICAL OUT
C
      COMMON /IOFILE/ IR,IW,IP,IJKO,IJKT,IDAF,NAV,IODA(950)
C
      DIMENSION IROT(NORBS,1),S(NORBS,1),X(1)
C
      DATA ZERO /0.0D+00/
C
      SQCDF=ZERO
      DO 10 I=1,NROT
   10 SQCDF=SQCDF+X(I)*X(I)
      IF(OUT) WRITE(IW,9999) SQCDF
C
      DO 20 IX=1,NORB
      DO 20 IG=1,NORBS
      S(IG,IX)=ZERO
      IB=IROT(IG,IX)
      IF(IB.EQ.0) GO TO 20
      VAL=X(IB)
      IF(IG.GT.IX) VAL=-VAL
      S(IG,IX)=VAL
   20 CONTINUE
      RETURN
 9999 FORMAT(/,' SQCDF = ',E20.7)
      END
C*MODULE MCQUD   *DECK SPMRSQ
      SUBROUTINE SPMRSQ(A,X,B,N,G,H,XI,XJ,IW,IA)
C
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      LOGICAL GOPARR,DSKWRK,MASWRK
C
      COMMON /PAR   / ME,MASTER,NPROC,IBTYP,IPTIM,GOPARR,DSKWRK,MASWRK
C
C     ----- SOLVE A * X = B , USING AN ITERATIVE PROCEDURE       -----
C
C     ----- NUMERICAL RECIPES (P.70), CAMBRIDGE UNIVERSITY PRESS -----
C          W.H.PRESS, B.P.FLANNERY, S.A.TEUKOLSKY, W.T.VETTERLING
C
      DIMENSION A(1),X(N),B(N),IA(1)
      DIMENSION G(N),H(N),XI(N),XJ(N)
      DATA ZERO   /0.0D+00/
      DATA EPS    /1.0D-05/
C
      DO 5 I=1,N
      B(I)=X(I)
    5 IA(I)=(I*(I-1))/2
C
      EPS2=N*EPS*EPS
      IRST=0
   10 IRST=IRST+1
      CALL SPMRAX(A,X,XI,N,IA)
      RP=ZERO
      BSQ=ZERO
      DO 11 J=1,N
         BSQ=BSQ+B(J)*B(J)
         XI(J)=XI(J)-B(J)
         RP=RP+XI(J)*XI(J)
   11 CONTINUE
      CALL SPMRXA(A,XI,G,N,IA)
      DO 12 J=1,N
      G(J)=-G(J)
      H(J)= G(J)
   12 CONTINUE
      MXITER = 10*N
      DO 19 ITER=1,MXITER
         CALL SPMRAX(A,H,XI,N,IA)
         ANUM=ZERO
         ADEN=ZERO
         DO 13 J=1,N
            ANUM=ANUM+G(J)*H(J)
            ADEN=ADEN+XI(J)*XI(J)
   13    CONTINUE
C
         IF(ADEN.EQ.ZERO) THEN
            IF (MASWRK) WRITE(IW,9998)
            CALL ABRT
            RETURN
         END IF
C
         ANUM=ANUM/ADEN
         DO 14 J=1,N
            XI(J)=X(J)
            X(J)=X(J)+ANUM*H(J)
   14    CONTINUE
         CALL SPMRAX(A,X,XJ,N,IA)
         RSQ=ZERO
         DO 15 J=1,N
            XJ(J)=XJ(J)-B(J)
            RSQ=RSQ+XJ(J)*XJ(J)
   15    CONTINUE
         IF(RSQ.EQ.RP.OR.RSQ.LE.BSQ*EPS2) RETURN
         IF(RSQ.GT.RP) THEN
            DO 16 J=1,N
               X(J)=XI(J)
   16       CONTINUE
            IF(IRST.GE.3) RETURN
            GO TO 10
         END IF
         RP=RSQ
         CALL SPMRXA(A,XJ,XI,N,IA)
         GG=ZERO
         DGG=ZERO
         DO 17 J=1,N
            GG=GG+G(J)*G(J)
            DGG=DGG+(XI(J)+G(J))*XI(J)
   17    CONTINUE
         IF(GG.EQ.ZERO) RETURN
         GAM=DGG/GG
         DO 18 J=1,N
            G(J)=-XI(J)
            H(J)=G(J)+GAM*H(J)
   18    CONTINUE
   19 CONTINUE
      IF (MASWRK) WRITE(IW,9999)
 9999 FORMAT(' IN SPMRSQ : TOO MANY ITERATIONS')
 9998 FORMAT(/,' STOP IN SPMRSQ : VERY SINGULAR MATRIX')
      RETURN
      END
C*MODULE MCQUD   *DECK SPMRAX
      SUBROUTINE SPMRAX(A,X,AX,N,IA)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      DIMENSION A(1),X(N),AX(N),IA(1)
C
C     ----- AX = X * A -----
C
      DO 30 I=1,N
         IPM=IA(I)
         AX(I)=A(IPM+I)*X(I)
         IF(I.EQ.1) GO TO 30
            DO 10 J=1,I-1
   10          AX(I)=AX(I) + A(IPM+J)*X(J)
            DO 20 J=1,I-1
   20          AX(J)=AX(J) + A(IPM+J)*X(I)
   30    CONTINUE
      RETURN
      END
C*MODULE MCQUD   *DECK SPMRXA
      SUBROUTINE SPMRXA(A,X,XA,N,IA)
      IMPLICIT DOUBLE PRECISION(A-H,O-Z)
      DIMENSION A(1),X(N),XA(N),IA(1)
C
C     ----- XA = A * X+ -----
C
      DO 30 I=1,N
         IPM=IA(I)
         XA(I)=A(IPM+I)*X(I)
         IF(I.EQ.1) GO TO 30
            DO 10 J=1,I-1
   10          XA(I)=XA(I) + A(IPM+J)*X(J)
            DO 20 J=1,I-1
   20          XA(J)=XA(J) + A(IPM+J)*X(I)
   30    CONTINUE
      RETURN
      END
