C***********************************************************************GMSDRV C GMSDRV C GMSDRV C G M S N B O GMSDRV C GMSDRV C GMSDRV C GAMESS VERSION OF NBO PROGRAM GMSDRV C (FOR GAMESS VERSION 27-JUN-05) GMSDRV C GMSDRV C GMSDRV C DRIVER ROUTINES: GMSDRV C GMSDRV C SUBROUTINE RUNNBO GMSDRV C SUBROUTINE FEAOIN(CORE,ICORE,MEMORY,NBOOPT) GMSDRV C SUBROUTINE DELSCF GMSDRV C SUBROUTINE CHKNBO(T) GMSDRV C SUBROUTINE INTS1E(ESNRG,FENRG,IFLG,IFLD) GMSDRV C SUBROUTINE GETEHF(EDEF,ECP,ITS) GMSDRV C SUBROUTINE DSKESX(ESNRG,PLNRG,EXNRG,SNRG,DMCA,DMCB,DMDA,DMDB,XX, GMSDRV C + IX) GMSDRV C SUBROUTINE DIRESX(ESNRG,PLNRG,EXNRG,SNRG,DMCA,DMCB,DMDA,DMDB, GMSDRV C + GHONDO,MAXG,DDIJ) GMSDRV C SUBROUTINE ESANDX(ESNRG,PLNRG,EXNRG,SNRG,VAL,I,J,K,L,DMCA,DMCB, GMSDRV C + DMDA,DMDB) GMSDRV C SUBROUTINE NXTINT(VAL,I,J,K,L,XX,IX,NINTMX,ICNT,IOFF,IREC,IS) GMSDRV C GMSDRV C***********************************************************************GMSDRV * SUBROUTINE RUNNBO GMSDRV C***********************************************************************GMSDRV C 11-Feb-93 EDG Checkpoint NBOs GMSDRV C-----------------------------------------------------------------------GMSDRV * IMPLICIT REAL*8 (A-H,O-Z) GMSDRV * DIMENSION NBOOPT(10) GMSDRV C GMSDRV * COMMON/NBIO/LFNIN,LFNPR,LFNAO,LFNPNA,LFNNAO,LFNPNH,LFNNHO,LFNPNB, GMSDRV * + LFNNBO,LFNPNL,LFNNLM,LFNMO,LFNDM,LFNNAB,LFNPPA,LFNARC, GMSDRV * + LFNDAF,LFNLBL,LFNDEF,LFNBRK(100) GMSDRV C GMSDRV C GAMESS Common Block: GMSDRV C GMSDRV * COMMON /IOFILE/ IR,IW,IP,IS,IPK,IDAF,NAV,IODA(950) GMSDRV * COMMON /FMCOM/ CORE(1) GMSDRV C GMSDRV * LFNIN = IR GMSDRV * LFNPR = IW GMSDRV C GMSDRV C Determine the amount of available memory for the NBO analysis. GMSDRV C GMSDRV * CALL VALFM(ICUR) GMSDRV * CALL GOTFM(MEMORY) GMSDRV C GMSDRV C Set NBO options. GMSDRV C GMSDRV * NBOOPT(1) = 0 GMSDRV * NBOOPT(2) = 0 GMSDRV * NBOOPT(3) = 0 GMSDRV * NBOOPT(4) = 0 GMSDRV * NBOOPT(5) = 0 GMSDRV * NBOOPT(6) = 0 GMSDRV * NBOOPT(7) = 1 GMSDRV * NBOOPT(8) = 0 GMSDRV * NBOOPT(9) = 0 GMSDRV * NBOOPT(10) = 6 GMSDRV C GMSDRV C Perform the NPA/NBO/NLMO analyses. GMSDRV C GMSDRV * CALL NBO(CORE(ICUR+1),MEMORY,NBOOPT) GMSDRV * IF(NBOOPT(10).LT.0) RETURN GMSDRV C GMSDRV C Perform the energetic analysis. GMSDRV C GMSDRV * NBOOPT(7) = 0 GMSDRV * 10 NBOOPT(1) = 2 GMSDRV * CALL NBOEAN(CORE(ICUR+1),MEMORY,NBOOPT,IDONE) GMSDRV * IF(IDONE.NE.0) GOTO 20 GMSDRV * CALL DELSCF GMSDRV * NBOOPT(1) = 3 GMSDRV * CALL NBOEAN(CORE(ICUR+1),MEMORY,NBOOPT,IDONE) GMSDRV * GOTO 10 GMSDRV C GMSDRV * 20 RETURN GMSDRV * END GMSDRV C***********************************************************************GMSDRV * SUBROUTINE FEAOIN(CORE,ICORE,MEMORY,NBOOPT) GMSDRV C***********************************************************************GMSDRV C 12-Jun-95 EDG Updated NBO interface for 10-Mar-95 version GMSDRV C 14-Jan-93 EDG Retrieve T and V integrals from dictionary file GMSDRV C-----------------------------------------------------------------------GMSDRV * IMPLICIT REAL*8 (A-H,O-Z) GMSDRV * DIMENSION CORE(MEMORY),ICORE(MEMORY),NBOOPT(10) GMSDRV C GMSDRV C ----------------------------------------------------------------------GMSDRV C GMSDRV C This routine fetchs basis set information from the GAMESS common GMSDRV C blocks and stores it in the NBO common blocks and direct access file GMSDRV C (DAF) for use by the NBO analysis. GMSDRV C GMSDRV C ----------------------------------------------------------------------GMSDRV C GMSDRV C Routine FEAOIN accesses the following records of the dictionary file:GMSDRV C GMSDRV C 11 --- Core Hamiltonian integrals GMSDRV C 12 --- AO overlap matrix GMSDRV C 13 --- Kinetic energy integrals GMSDRV C 14 --- AO Fock matrix (alpha) GMSDRV C 15 --- AO to MO transformation matrix (alpha) GMSDRV C 16 --- AO density matrix (bond order matrix) (alpha) GMSDRV C 18 --- AO Fock matrix (beta) GMSDRV C 19 --- AO to MO transformation matrix (beta) GMSDRV C 20 --- AO density matrix (bond order matrix) (beta) GMSDRV C 95 --- X dipole integrals GMSDRV C 96 --- Y dipole integrals GMSDRV C 97 --- Z dipole integrals GMSDRV C GMSDRV C ----------------------------------------------------------------------GMSDRV C GMSDRV C NBO Common blocks GMSDRV C GMSDRV C GMSDRV * PARAMETER(MAXATM = 200,MAXBAS = 2000) GMSDRV * COMMON/NBFLAG/ROHF,UHF,CI,OPEN,COMPLX,ALPHA,BETA,MCSCF,AUHF,ORTHO GMSDRV * LOGICAL ROHF,UHF,CI,OPEN,COMPLX,ALPHA,BETA,MCSCF,AUHF,ORTHO GMSDRV * COMMON/NBINFO/ISPIN,NATOMS,NDIM,NBAS,NNAO,MXBO,MXAO,MXAOLM,MUNIT GMSDRV * COMMON/NBOPT/IWDM,IW3C,IWAPOL,IWHYBS,IWPNAO,IWTNAO,IWTNAB, GMSDRV * + IWTNBO,IWFOCK,IWCUBF,IPSEUD,KOPT,IPRINT,IWDETL,IWMULP,ICHOOS, GMSDRV * + IWNBBP,IWMSP,IWFIXDM,IW3CHB,IWNJC,JCORE,JPRINT(100) GMSDRV * COMMON/NBAO/LCTR(MAXBAS),LANG(MAXBAS) GMSDRV * COMMON/NBATOM/IATNO(MAXATM),INO(MAXATM),NORBS(MAXATM),LL(MAXATM), GMSDRV * + LU(MAXATM),IZNUC(MAXATM),IATCR(MAXATM) GMSDRV * COMMON/NBIO/LFNIN,LFNPR,LFNAO,LFNPNA,LFNNAO,LFNPNH,LFNNHO,LFNPNB, GMSDRV * + LFNNBO,LFNPNL,LFNNLM,LFNMO,LFNDM,LFNNAB,LFNPPA,LFNARC, GMSDRV * + LFNDAF,LFNLBL,LFNDEF,LFNBRK(100) GMSDRV C GMSDRV * DIMENSION LABELS(35),WFNS(6),NCORE(12) GMSDRV * LOGICAL WSTATE(6,6) GMSDRV C GMSDRV C GAMESS Common blocks GMSDRV C GMSDRV * PARAMETER (MXGTOT=20000, MXSH=5000, MXATM=2000, MXRT=100) GMSDRV * COMMON /ECP2 / CLP(400),ZLP(400),NLP(400),KFIRST(MXATM,6), GMSDRV * * KLAST(MXATM,6),LMAX(MXATM),LPSKIP(MXATM), GMSDRV * * IZCORE(MXATM) GMSDRV * COMMON /ENRGYS/ ENUCR,EELCT,ETOT,SZ,SZZ,ECORE,ESCF,EERD,E1,E2, GMSDRV * * VEN,VEE,EPOT,EKIN,ESTATE(MXRT),STATN,EDFT(2) GMSDRV * COMMON /INFOA / NAT,ICH,MUL,NUM,NX,NE,NA,NB, GMSDRV * * ZAN(MXATM),C(3,MXATM),IAN(MXATM) GMSDRV * COMMON /IOFILE/ IR,IW,IP,IS,IPK,IDAF,NAV,IODA(950) GMSDRV * COMMON /NSHEL / EX(MXGTOT),CS(MXGTOT),CP(MXGTOT),CD(MXGTOT), GMSDRV * * CF(MXGTOT),CG(MXGTOT),CH(MXGTOT),CII(MXGTOT), GMSDRV * * KSTART(MXSH),KATOM(MXSH),KTYPE(MXSH),KNG(MXSH), GMSDRV * * KLOC(MXSH),KMIN(MXSH),KMAX(MXSH),NSHELL GMSDRV * COMMON /RUNLAB/ TITLE(10),A(MXATM),B(MXATM),BFLAB(8192) GMSDRV * COMMON /SCFOPT/ CONVHF,MAXIT,MCONV,NPUNCH,NPREO(4) GMSDRV * COMMON /WFNOPT/ SCFTYP,CITYP,DFTYPE,CCTYP,MPLEVL,MPCTYP GMSDRV * COMMON /XYZPRP/ X(3),PAD(35) GMSDRV C GMSDRV C Obtain the following information: GMSDRV C GMSDRV C ROHF =.TRUE. If RHF open shell wavefunction GMSDRV C =.FALSE. otherwise GMSDRV C GMSDRV C UHF =.TRUE. If UHF wavefunction GMSDRV C =.FALSE. otherwise GMSDRV C GMSDRV C AUHF =.TRUE. If spin-annihilated UHF wavefunction GMSDRV C =.FALSE. otherwise GMSDRV C GMSDRV C CI =.TRUE. If CI wavefunction GMSDRV C =.FALSE. otherwise GMSDRV C GMSDRV C OPEN =.TRUE. If open shell wavefunction GMSDRV C =.FALSE. otherwise GMSDRV C GMSDRV C COMPLX =.TRUE. If complex wavefunction GMSDRV C =.FALSE. otherwise GMSDRV C (Note: The program is not capable of handling this.) GMSDRV C GMSDRV C NATOMS Number of atomic centers GMSDRV C GMSDRV C NDIM Dimension of matrices (overlap and density) GMSDRV C GMSDRV C NBAS Number of basis functions (.le.NDIM) GMSDRV C GMSDRV C IPSEUD Set to one if pseudopotentials are used. GMSDRV C GMSDRV C IWCUBF This pertains only basis sets with F functions. GMSDRV C GMSDRV C If cartesian F functions are input, set IWCUBF to: GMSDRV C 0, if these are to be transformed to the GMSDRV C standard set of pure F functions GMSDRV C 1, if these are to be transformed to the GMSDRV C cubic set of pure F functions GMSDRV C GMSDRV C If pure F functions are input, set to IWCUBF to: GMSDRV C 0, if these are standard F functions GMSDRV C 1, if these are cubic F functions GMSDRV C GMSDRV C IATNO(I),I=1,NATOMS GMSDRV C List of atomic numbers GMSDRV C GMSDRV C LCTR(I),I=1,NBAS GMSDRV C List of atomic centers of the basis functions GMSDRV C (LCTR(3)=2 if basis function 3 is on atom 2) GMSDRV C GMSDRV C LANG(I),I=1,NBAS GMSDRV C List of angular symmetry information for the basis GMSDRV C functions GMSDRV C GMSDRV C LABELS array contains NBO labels for the atomic orbitals GMSDRV C GMSDRV * DATA LABELS / GMSDRV C GMSDRV C s GMSDRV C --- GMSDRV * + 1, GMSDRV C GMSDRV C px py pz GMSDRV C --- --- --- GMSDRV * + 101, 102, 103, GMSDRV C GMSDRV C dxx dyy dzz dxy dxz dyz GMSDRV C --- --- --- --- --- --- GMSDRV * + 201, 204, 206, 202, 203, 205, GMSDRV C GMSDRV C fxxx fyyy fzzz fxxy fxxz fxyy fyyz fxzz fyzz fxyz GMSDRV C ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- GMSDRV * + 301, 307, 310, 302, 303, 304, 308, 306, 309, 305, GMSDRV C GMSDRV C gxxxx gyyyy gzzzz gxxxy gxxxz gxyyy gyyyz gxzzz gyzzz gxxyy GMSDRV C ----- ----- ----- ----- ----- ----- ----- ----- ----- ----- GMSDRV * + 401, 411, 415, 402, 403, 407, 412, 410, 414, 404, GMSDRV C GMSDRV C gxxzz gyyzz gxxyz gxyyz gxyzz GMSDRV C ----- ----- ----- ----- ----- GMSDRV * + 406, 413, 405, 408, 409 / GMSDRV C GMSDRV C GMSDRV C WSTATE array contains the values which should be set in the NBO GMSDRV C common block /NBFLAG/ depending on wavefunction. GMSDRV C GMSDRV * DATA WSTATE / GMSDRV C logical variable in common NBFLAG GMSDRV C ROHF UHF CI OPEN MCSCF AUHF GMSDRV C ------- ------- ------ ------ ------ ------ GMSDRV C Wavefunction GMSDRV C RHF GMSDRV * + .FALSE., .FALSE., .FALSE., .FALSE., .FALSE., .FALSE., GMSDRV C UHF GMSDRV * + .FALSE., .TRUE. , .FALSE., .TRUE. , .FALSE., .FALSE., GMSDRV C ROHF GMSDRV * + .TRUE. , .FALSE., .FALSE., .TRUE. , .FALSE., .FALSE., GMSDRV C GVB GMSDRV * + .TRUE., .FALSE., .FALSE., .TRUE. , .FALSE., .FALSE., GMSDRV C MCSCF GMSDRV * + .FALSE., .FALSE., .FALSE., .FALSE., .TRUE. , .FALSE., GMSDRV C CI GMSDRV * + .FALSE., .FALSE., .TRUE. , .FALSE., .FALSE., .FALSE./ GMSDRV C GMSDRV C GMSDRV C Wavefunction types: GMSDRV C GMSDRV * DATA WFNS /8HRHF , GMSDRV * + 8HUHF , GMSDRV * + 8HROHF , GMSDRV * + 8HGVB , GMSDRV * + 8HMCSCF , GMSDRV * + 8HCI / GMSDRV C GMSDRV * DATA NONE/4HNONE/ GMSDRV * DATA NCORE/2,10,18,28,36,46,54,68,78,86,100,110/ GMSDRV * DATA ZERO/0.0D0/ GMSDRV * DATA TOANG/0.529177249D0/ GMSDRV C GMSDRV C Store job title on NBODAF: GMSDRV C GMSDRV * DO 5 I = 1,10 GMSDRV * CORE(I) = TITLE(I) GMSDRV * 5 CONTINUE GMSDRV * NFILE = 2 GMSDRV * CALL NBWRIT(CORE,10,NFILE) GMSDRV C GMSDRV C Get the number of atoms from NAT and store the atomic numbers in GMSDRV C IATNO and nuclear charges in IZNUC. (NOTE: atomic numbers and GMSDRV C nuclear charges may not be equivalent if effective core potentials GMSDRV C (ECP) are used.) GMSDRV C GMSDRV * NATOMS = NAT GMSDRV * DO 10 I = 1,NAT GMSDRV * IATNO(I) = ZAN(I) + IZCORE(I) GMSDRV * IZNUC(I) = ZAN(I) GMSDRV * IF(IZCORE(I).NE.0) IPSEUD = 1 GMSDRV * 10 CONTINUE GMSDRV C GMSDRV C If this is a MOPAC wavefunction, reduce IZNUC by the number of core GMSDRV C electrons and turn pseudopotentials on: GMSDRV C GMSDRV * IF(MPCTYP.NE.NONE) THEN GMSDRV * DO 15 I = 1,NAT GMSDRV * DO 14 J = 1,12 GMSDRV * ITMP = IATNO(I) - NCORE(J) GMSDRV * IF(ITMP.GT.0) IZNUC(I) = ITMP GMSDRV * IF(ITMP.LE.0) GOTO 15 GMSDRV * IPSEUD = 1 GMSDRV * 14 CONTINUE GMSDRV * 15 CONTINUE GMSDRV * END IF GMSDRV C GMSDRV C KATOM array contains which atom the shell is on, KMIN and KMAX GMSDRV C determine the components in the shell by pointing to a range in the GMSDRV C LABELS array: GMSDRV C GMSDRV * II = 0 GMSDRV * DO 30 I = 1,NSHELL GMSDRV * IATOM = KATOM(I) GMSDRV * MIN = KMIN(I) GMSDRV * MAX = KMAX(I) GMSDRV * DO 20 J = MIN,MAX GMSDRV * II = II + 1 GMSDRV * LCTR(II) = IATOM GMSDRV * LANG(II) = LABELS(J) GMSDRV * 20 CONTINUE GMSDRV * 30 CONTINUE GMSDRV C GMSDRV * NBAS = II GMSDRV * NDIM = NBAS GMSDRV C GMSDRV C Inititialize various NBO options depending upon the wavefunction GMSDRV C type and basis set type. GMSDRV C GMSDRV C First, turn off the complex orbitals, indicate that the pure set GMSDRV C of F functions is desired when transforming from the cartesian set. GMSDRV C GMSDRV * COMPLX = .FALSE. GMSDRV * IWCUBF = 0 GMSDRV * ORTHO = MPCTYP.NE.NONE GMSDRV C GMSDRV C Next set up the wavefunction flags. GMSDRV C GMSDRV * DO 50 I = 1,6 GMSDRV * ISTATE = I GMSDRV * IF (SCFTYP.EQ.WFNS(I)) GOTO 60 GMSDRV * 50 CONTINUE GMSDRV * CALL NBHALT('FEAOIN: Unknown WFNTYP.') GMSDRV C GMSDRV * 60 ROHF = WSTATE(1,ISTATE) GMSDRV * UHF = WSTATE(2,ISTATE) GMSDRV * CI = WSTATE(3,ISTATE) GMSDRV * OPEN = WSTATE(4,ISTATE) GMSDRV * MCSCF = WSTATE(5,ISTATE) GMSDRV * AUHF = WSTATE(6,ISTATE) GMSDRV C GMSDRV C No Fock matrices for ROHF, MCSCF, or CI wavefunctions: GMSDRV C GMSDRV * IF (ROHF.OR.MCSCF.OR.CI) IWFOCK = 0 GMSDRV C GMSDRV C Expectation values of the Fock operator are in atomic units: GMSDRV C GMSDRV * MUNIT = 0 GMSDRV C GMSDRV C Store NATOMS, NDIM, NBAS, MUNIT, wavefunction flags, ISWEAN: GMSDRV C GMSDRV * ICORE(1) = NATOMS GMSDRV * ICORE(2) = NDIM GMSDRV * ICORE(3) = NBAS GMSDRV * ICORE(4) = MUNIT GMSDRV * ICORE(5) = 0 GMSDRV * IF(ROHF) ICORE(5) = 1 GMSDRV * ICORE(6) = 0 GMSDRV * IF(UHF) ICORE(6) = 1 GMSDRV * ICORE(7) = 0 GMSDRV * IF(CI) ICORE(7) = 1 GMSDRV * ICORE(8) = 0 GMSDRV * IF(OPEN) ICORE(8) = 1 GMSDRV * ICORE(9) = 0 GMSDRV * IF(MCSCF) ICORE(9) = 1 GMSDRV * ICORE(10) = 0 GMSDRV * IF(AUHF) ICORE(10) = 1 GMSDRV * ICORE(11) = 0 GMSDRV * IF(ORTHO) ICORE(11) = 1 GMSDRV * ICORE(12) = 1 GMSDRV * ICORE(13) = 0 GMSDRV * NFILE = 3 GMSDRV * CALL NBWRIT(ICORE,13,NFILE) GMSDRV C GMSDRV C Store IATNO, IZNUC, LCTR, and LANG on NBO DAF: GMSDRV C GMSDRV * II = 0 GMSDRV * DO 70 I = 1,NATOMS GMSDRV * II = II + 1 GMSDRV * ICORE(II) = IATNO(I) GMSDRV * 70 CONTINUE GMSDRV * DO 80 I = 1,NATOMS GMSDRV * II = II + 1 GMSDRV * ICORE(II) = IZNUC(I) GMSDRV * 80 CONTINUE GMSDRV * DO 90 I = 1,NBAS GMSDRV * II = II + 1 GMSDRV * ICORE(II) = LCTR(I) GMSDRV * 90 CONTINUE GMSDRV * DO 95 I = 1,NBAS GMSDRV * II = II + 1 GMSDRV * ICORE(II) = LANG(I) GMSDRV * 95 CONTINUE GMSDRV * NFILE = 4 GMSDRV * CALL NBWRIT(ICORE,2*NATOMS+2*NBAS,NFILE) GMSDRV C GMSDRV C Store the total energy on the NBO DAF: GMSDRV C GMSDRV * CORE(1) = ETOT GMSDRV * CORE(2) = ETOT GMSDRV * NFILE = 8 GMSDRV * CALL NBWRIT(CORE,2,NFILE) GMSDRV C GMSDRV C Store the atomic coordinates on the NBO DAF: (Note that these GMSDRV C coordinates are used in the calculation of dipole moments. GAMESS GMSDRV C requires the Cartesian origin to be at the center of mass!!) GMSDRV C GMSDRV * I = 0 GMSDRV * DO 110 IAT = 1,NATOMS GMSDRV * DO 100 K = 1,3 GMSDRV * I = I + 1 GMSDRV * CORE(I) = (C(K,IAT) - X(K)) * TOANG GMSDRV * 100 CONTINUE GMSDRV * 110 CONTINUE GMSDRV * NFILE = 9 GMSDRV * CALL NBWRIT(CORE,3*NATOMS,NFILE) GMSDRV C GMSDRV C Store the overlap matrix on the NBODAF: GMSDRV C GMSDRV * NFILE = 12 GMSDRV * L2 = NDIM*(NDIM+1)/2 GMSDRV * CALL DAREAD(IDAF,IODA,CORE,L2,NFILE,0) GMSDRV * NFILE = 10 GMSDRV * CALL NBWRIT(CORE,L2,NFILE) GMSDRV C GMSDRV C Store the density matrices on the NBODAF: GMSDRV C GMSDRV * NFILE = 16 GMSDRV * L2 = NDIM*(NDIM+1)/2 GMSDRV * CALL DAREAD(IDAF,IODA,CORE,L2,NFILE,0) GMSDRV * NFILE = 20 GMSDRV * CALL NBWRIT(CORE,L2,NFILE) GMSDRV C GMSDRV * IF(OPEN) THEN GMSDRV * NFILE = 20 GMSDRV * CALL DAREAD(IDAF,IODA,CORE,L2,NFILE,0) GMSDRV * NFILE = 21 GMSDRV * CALL NBWRIT(CORE,L2,NFILE) GMSDRV * END IF GMSDRV C GMSDRV C Store the Fock matrices on the NBODAF: GMSDRV C GMSDRV * IF(IWFOCK.NE.0) THEN GMSDRV * NFILE = 14 GMSDRV * L2 = NDIM*(NDIM+1)/2 GMSDRV * IF(IODA(NFILE).GT.0) THEN GMSDRV * CALL DAREAD(IDAF,IODA,CORE,L2,NFILE,0) GMSDRV * NFILE = 30 GMSDRV * CALL NBWRIT(CORE,L2,NFILE) GMSDRV * END IF GMSDRV C GMSDRV * IF(OPEN) THEN GMSDRV * NFILE = 18 GMSDRV * IF(IODA(NFILE).GT.0) THEN GMSDRV * CALL DAREAD(IDAF,IODA,CORE,L2,NFILE,0) GMSDRV * NFILE = 31 GMSDRV * CALL NBWRIT(CORE,L2,NFILE) GMSDRV * END IF GMSDRV * END IF GMSDRV * END IF GMSDRV C GMSDRV C Store the AO to MO transformation matrices on the NBODAF: GMSDRV C GMSDRV * IF(IODA(15).GT.0) THEN GMSDRV * NFILE = 15 GMSDRV * L3 = NDIM*NDIM GMSDRV * CALL DAREAD(IDAF,IODA,CORE,L3,NFILE,0) GMSDRV * NFILE = 40 GMSDRV * CALL NBWRIT(CORE,L3,NFILE) GMSDRV * IF(OPEN) THEN GMSDRV * NFILE = 19 GMSDRV * CALL DAREAD(IDAF,IODA,CORE,L3,NFILE,0) GMSDRV * NFILE = 41 GMSDRV * CALL NBWRIT(CORE,L3,NFILE) GMSDRV * END IF GMSDRV * END IF GMSDRV C GMSDRV C Store the kinetic energy integrals on the NBODAF: GMSDRV C GMSDRV * IF(IODA(13).GT.0) THEN GMSDRV * NFILE = 13 GMSDRV * L2 = NDIM*(NDIM+1)/2 GMSDRV * CALL DAREAD(IDAF,IODA,CORE,L2,NFILE,0) GMSDRV * NFILE = 18 GMSDRV * CALL NBWRIT(CORE,L2,NFILE) GMSDRV C GMSDRV C Store the nuclear attraction integrals on the NBODAF: GMSDRV C GMSDRV * NFILE = 11 GMSDRV * L2 = NDIM*(NDIM+1)/2 GMSDRV * I1 = 1 GMSDRV * I2 = I1 + L2 GMSDRV * CALL DAREAD(IDAF,IODA,CORE(I2),L2,NFILE,0) GMSDRV * DO 115 I = 0,L2-1 GMSDRV * CORE(I2+I) = CORE(I2+I) - CORE(I1+I) GMSDRV * 115 CONTINUE GMSDRV * NFILE = 19 GMSDRV * CALL NBWRIT(CORE(I2),L2,NFILE) GMSDRV * END IF GMSDRV C GMSDRV C Store the x,y,z dipole integrals on the NBODAF: GMSDRV C GMSDRV * IF(IODA(95).GT.0.AND.IODA(96).GT.0.AND.IODA(97).GT.0) THEN GMSDRV * L2 = NDIM*(NDIM+1)/2 GMSDRV * NFILE = 95 GMSDRV * CALL DAREAD(IDAF,IODA,CORE,L2,NFILE,0) GMSDRV * DO 120 I = 1,L2 GMSDRV * CORE(I) = CORE(I) * TOANG GMSDRV * 120 CONTINUE GMSDRV * NFILE = 50 GMSDRV * CALL NBWRIT(CORE,L2,NFILE) GMSDRV * NFILE = 96 GMSDRV * CALL DAREAD(IDAF,IODA,CORE,L2,NFILE,0) GMSDRV * DO 130 I = 1,L2 GMSDRV * CORE(I) = CORE(I) * TOANG GMSDRV * 130 CONTINUE GMSDRV * NFILE = 51 GMSDRV * CALL NBWRIT(CORE,L2,NFILE) GMSDRV * NFILE = 97 GMSDRV * CALL DAREAD(IDAF,IODA,CORE,L2,NFILE,0) GMSDRV * DO 140 I = 1,L2 GMSDRV * CORE(I) = CORE(I) * TOANG GMSDRV * 140 CONTINUE GMSDRV * NFILE = 52 GMSDRV * CALL NBWRIT(CORE,L2,NFILE) GMSDRV * END IF GMSDRV C GMSDRV C Store the AO basis set info on the NBO DAF: (Note that two integers GMSDRV C and three integer arrays are stored first. Also remember that ICORE GMSDRV C and CORE occupy the same memory.) GMSDRV C GMSDRV * NEXP = 0 GMSDRV * DO 150 I = 1,MXGTOT GMSDRV * IF(EX(I).EQ.ZERO) GOTO 150 GMSDRV * NEXP = I GMSDRV * 150 CONTINUE GMSDRV * DO 160 I = 1,2+3*NSHELL+6*NEXP GMSDRV * CORE(I) = ZERO GMSDRV * 160 CONTINUE GMSDRV * ICORE(1) = NSHELL GMSDRV * ICORE(2) = NEXP GMSDRV C GMSDRV C NCOMP(I) -- the number of components in the Ith shell: GMSDRV C GMSDRV * II = 2 GMSDRV * DO 170 I = 1,NSHELL GMSDRV * II = II + 1 GMSDRV * ICORE(II) = KMAX(I) - KMIN(I) + 1 GMSDRV * 170 CONTINUE GMSDRV C GMSDRV C NPRIM(I) -- the number of gaussian primitives in the Ith shell: GMSDRV C GMSDRV * DO 180 I = 1,NSHELL GMSDRV * II = II + 1 GMSDRV * ICORE(II) = KNG(I) GMSDRV * 180 CONTINUE GMSDRV C GMSDRV C NPTR(I) -- pointer for the Ith shell into the gaussian parameters, GMSDRV C EXP, CS, CP, etc.: GMSDRV C GMSDRV * DO 190 I = 1,NSHELL GMSDRV * II = II + 1 GMSDRV * ICORE(II) = KSTART(I) GMSDRV * 190 CONTINUE GMSDRV C GMSDRV C EXP(I) -- orbital exponents indexed by NPTR: GMSDRV C GMSDRV * DO 200 I = 1,NEXP GMSDRV * II = II + 1 GMSDRV * CORE(II) = EX(I) GMSDRV * 200 CONTINUE GMSDRV C GMSDRV C CS,CP,CD,CF,CG -- orbital coefficients: GMSDRV C GMSDRV * DO 210 I = 1,NEXP GMSDRV * II = II + 1 GMSDRV * CORE(II) = CS(I) GMSDRV * 210 CONTINUE GMSDRV * DO 220 I = 1,NEXP GMSDRV * II = II + 1 GMSDRV * CORE(II) = CP(I) GMSDRV * 220 CONTINUE GMSDRV * DO 230 I = 1,NEXP GMSDRV * II = II + 1 GMSDRV * CORE(II) = CD(I) GMSDRV * 230 CONTINUE GMSDRV * DO 240 I = 1,NEXP GMSDRV * II = II + 1 GMSDRV * CORE(II) = CF(I) GMSDRV * 240 CONTINUE GMSDRV * DO 250 I = 1,NEXP GMSDRV * II = II + 1 GMSDRV * CORE(II) = CG(I) GMSDRV * 250 CONTINUE GMSDRV * NFILE = 5 GMSDRV * CALL NBWRIT(CORE,II,NFILE) GMSDRV * RETURN GMSDRV * END GMSDRV C***********************************************************************GMSDRV * SUBROUTINE DELSCF GMSDRV C***********************************************************************GMSDRV C 3-Jun-05 EDG Add DFT switch to record 70 of DAF GMSDRV C 13-Jun-97 EDG Subroutine overhaul GMSDRV C 1-Aug-96 EDG Add electric fields to NEDA GMSDRV C 13-Jun-96 EDG Updated NBO interface for 22-Nov-95 version GMSDRV C 12-Jun-95 EDG Updated NBO interface for 10-Mar-95 version GMSDRV C 2-Jun-95 EDG Add evaluation of self energies GMSDRV C 27-Apr-95 EDG Add direct calculation of ES and EX GMSDRV C 17-Jan-95 EDG Calculate classical ES contribution in NEDA GMSDRV C 24-Jul-94 EDG Calculate EX contribution to ES in NEDA GMSDRV C 18-Mar-94 EDG Fixed direct SCF NEDA GMSDRV C 3-Mar-94 EDG Add UHF capabilities for NEDA GMSDRV C 5-Nov-93 EDG Add direct energy minimization to accelerate SCF GMSDRV C convergence GMSDRV C 5-Nov-93 EDG Add DAMPing to accelerate SCF convergence GMSDRV C 16-Jun-93 EDG Add ECPs to NEDA and cleanup GMSDRV C 3-Jun-93 EDG Add DIIS to accelerate SCF convergence GMSDRV C 30-May-93 EDG Shift origin of dipole moment to fragment center of GMSDRV C mass GMSDRV C 18-Jan-93 EDG Compute dipole moments for DEF and CP fragment wave GMSDRV C functions GMSDRV C 12-Jan-93 EDG Added COMMON/FMCOM/ and memory allocation routines to GMSDRV C correct problems with CRAY calls to SR TFTRI GMSDRV C-----------------------------------------------------------------------GMSDRV * IMPLICIT REAL*8 (A-H,O-Z) GMSDRV * EXTERNAL PACK GMSDRV C GMSDRV * PARAMETER(MAXATM = 200,MAXBAS = 2000) GMSDRV C GMSDRV * LOGICAL NEW,ERROR,DEL GMSDRV * DIMENSION DDEF(3),DFLD(3),DCP(3),DIND(3),COM(3) GMSDRV * DIMENSION SNRG(MAXATM),FENRG(MAXATM),FPNRG(MAXATM),FSNRG(MAXATM) GMSDRV C GMSDRV C NBO common blocks: GMSDRV C GMSDRV * COMMON/NBEDA/KFRG,NFRG,IFRG(MAXBAS),NFE(MAXBAS),NFEB(MAXBAS) GMSDRV * COMMON/NBBAS/LABEL(MAXBAS,6),NBOUNI(MAXBAS),NBOTYP(MAXBAS), GMSDRV * + IATNO(MAXBAS),IBXM(MAXBAS),LOCC(2*MAXBAS),ISCR(2*MAXBAS) GMSDRV * COMMON/NBINFO/ISPIN,NATOMS,NDIM,NBAS,NNAO,MXBO,MXAO,MXAOLM,MUNIT GMSDRV * COMMON/NBFLAG/ROHF,UHF,CI,OPEN,COMPLX,ALPHA,BETA,MCSCF,AUHF,ORTHO GMSDRV * LOGICAL ROHF,UHF,CI,OPEN,COMPLX,ALPHA,BETA,MCSCF,AUHF,ORTHO GMSDRV * COMMON/NBIO/LFNIN,LFNPR,LFNAO,LFNPNA,LFNNAO,LFNPNH,LFNNHO,LFNPNB, GMSDRV * + LFNNBO,LFNPNL,LFNNLM,LFNMO,LFNDM,LFNNAB,LFNPPA,LFNARC, GMSDRV * + LFNDAF,LFNLBL,LFNDEF,LFNBRK(100) GMSDRV C GMSDRV C GAMESS common blocks: GMSDRV C GMSDRV * PARAMETER (MXGTOT=20000, MXSH=5000, MXGSH=30, MXG2=MXGSH*MXGSH, GMSDRV * + MXATM=2000, MXRT=100, MXFRG=1050, MXDFG=5, GMSDRV * + MXDPPT=MXFRG*MXDFG*12) GMSDRV * COMMON /ATHEAT/ ATHEAT GMSDRV * COMMON /CONV / DENTOL,EN,ETOT,EHF,EHF0,DIFF,ITER,ICALP,ICBET GMSDRV * COMMON /DEMOPT/ DEMCUT,IDEM GMSDRV * COMMON /DFGRID/ DFTTHR,DFTGTHR,SWOFF,SW0,NDFTFG, GMSDRV * + NRAD,NTHE,NPHI,NRAD0,NTHE0,NPHI0 GMSDRV * COMMON /DMPING/ SHIFT0,SHIFTV,DMPCUT,SWDIIS,DIRTHR GMSDRV * COMMON /ECP2 / CLP(400),ZLP(400),NLP(400),KFIRST(MXATM,6), GMSDRV * + KLAST(MXATM,6),LMAX(MXATM),LPSKIP(MXATM), GMSDRV * + IZCORE(MXATM) GMSDRV * COMMON /EFLDC / EVEC(3),EFLDL GMSDRV * COMMON /ENRGYS/ ENUCR,EELCT,ET,SZ,SZZ,ECORE,ESCF,EERD,E1,E2, GMSDRV * + VEN,VEE,EPOT,EKIN,ESTATE(MXRT),STATN,EDFT(2) GMSDRV * COMMON /FMCOM / A(1) GMSDRV * COMMON /FRGINF/ NMPTS(MXFRG),NMTTPT,IEFC,IEFD,IEFQ,IEFO, GMSDRV * * NPPTS(MXFRG),NPTTPT,IEFP, GMSDRV * * NRPTS(MXFRG),NRTTPT,IREP,ICHGP,NFRAG, GMSDRV * * NDPPTS(MXDPPT),NDPTTPT,IEFDP GMSDRV * COMMON /FUNCT / E,EG(3*MXATM) GMSDRV * COMMON /IJPAIR/ IA(8192) GMSDRV * COMMON /INFOA / NAT,ICH,MUL,NUM,NX,NE,NA,NB, GMSDRV * + ZAN(MXATM),C(3,MXATM),IAN(MXATM) GMSDRV * COMMON /INTFIL/ NINTMX,NHEX,NTUPL,PACK2E,INTTYP,IGRDTYP GMSDRV * COMMON /INTOPT/ ISCHWZ,IECP,NECP,IEFLD GMSDRV * COMMON /IOFILE/ IR,IW,IP,IS,IPK,IDAF,NAV,IODA(950) GMSDRV * COMMON /MACHIN/ NWDVAR,MAXFM,MAXSM,LIMFM,LIMSM GMSDRV * COMMON /MASSES/ ZMASS(MXATM) GMSDRV * COMMON /N2ELCT/ N2EL GMSDRV * COMMON /NSHEL / EX(MXGTOT),CS(MXGTOT),CP(MXGTOT),CD(MXGTOT), GMSDRV * + CF(MXGTOT),CG(MXGTOT),CH(MXGTOT),CII(MXGTOT), GMSDRV * + KSTART(MXSH),KATOM(MXSH),KTYPE(MXSH),KNG(MXSH), GMSDRV * + KLOC(MXSH),KMIN(MXSH),KMAX(MXSH),NSHELL GMSDRV * COMMON /OPTSCF/ DIRSCF,FDIFF GMSDRV * COMMON /OUTPUT/ NPRINT,ITOL,ICUT,NORMF,NORMP,NOPK GMSDRV * COMMON /RUNOPT/ RUNTYP,EXETYP,NEVALS,NGLEVL,NHLEVL GMSDRV * COMMON /SCFOPT/ CONVHF,MAXIT,MCONV,NPUNCH,NPREO(4) GMSDRV * COMMON /SYMMOL/ GROUP,COMPLEX,IGROUP,NAXIS,ILABMO GMSDRV * COMMON /WFNOPT/ SCFTYP,CITYP,DFTYPE,CCTYP,MPLEVL,MPCTYP GMSDRV C GMSDRV * LOGICAL PACK2E,DIRSCF,FDIFF,EFLDL GMSDRV C GMSDRV * DATA ZERO,ONE,TWO/0.0D0,1.0D0,2.0D0/ GMSDRV * DATA DEBYE/4.803242D0/ GMSDRV * DATA ANONE/8HNONE / GMSDRV * SAVE ZERO,ONE,TWO GMSDRV * SAVE DEBYE,ESNRG,FENRG GMSDRV * SAVE ANONE GMSDRV C GMSDRV C-----------------------------------------------------------------------GMSDRV C GMSDRV C Reduce GAMESS print level: GMSDRV C GMSDRV * MPRINT = NPRINT GMSDRV * NPRINT = -5 GMSDRV C GMSDRV C Open the NBO DAF: GMSDRV C GMSDRV * NEW = .FALSE. GMSDRV * CALL NBOPEN(NEW,ERROR) GMSDRV * IF(ERROR) CALL NBHALT('Error opening NBO DAF in SR DELSCF.') GMSDRV C GMSDRV C Appropriate GAMESS parameters for NEDA? GMSDRV C GMSDRV * IF(KFRG.GT.0.AND.(IGROUP.NE.1.OR.NOPK.NE.1)) THEN GMSDRV * WRITE(LFNPR,900) GMSDRV * CALL NBHALT(' ') GMSDRV * END IF GMSDRV C GMSDRV C Begin NEDA output: GMSDRV C GMSDRV * IF(KFRG.GT.0) WRITE(LFNPR,910) KFRG GMSDRV C GMSDRV C Allocate memory: GMSDRV C GMSDRV * IFUHF = 0 GMSDRV * IF(UHF) IFUHF = 1 GMSDRV * IFFLD = 0 GMSDRV * ISUM = IEFC + IEFD + IEFQ + IEFO + IEFP + NFRAG GMSDRV * IF(EFLDL) IFFLD = 1 GMSDRV * IF(ISUM.NE.0) IFFLD = 1 GMSDRV * L1 = NDIM GMSDRV * L2 = NDIM * (NDIM + 1) / 2 GMSDRV * L3 = NDIM * NDIM GMSDRV * IOFF = (KFRG - 1) * L2 GMSDRV * CALL BASCHK(LMX) GMSDRV * NANGM = 4 GMSDRV * IF(LMX.EQ.2) NANGM = 6 GMSDRV * IF(LMX.EQ.3) NANGM = 10 GMSDRV * IF(LMX.EQ.4) NANGM = 15 GMSDRV * MAXG = NANGM**4 GMSDRV * CALL VALFM(LOADFM) GMSDRV * LSCR = LOADFM + 1 GMSDRV * LDMDA = LSCR + L3 GMSDRV * LDMDB = LDMDA + NFRG*L2 GMSDRV * LDMFA = LDMDB + NFRG*L2 * IFUHF GMSDRV * LDMFB = LDMFA + NFRG*L2 * IFFLD GMSDRV * LDMCA = LDMFB + NFRG*L2 * IFUHF * IFFLD GMSDRV * LDMCB = LDMCA + NFRG*L2 GMSDRV * LVNUC = LDMCB + NFRG*L2 * IFUHF GMSDRV * LXINT = LVNUC + NFRG*L2 GMSDRV * LYINT = LXINT + L2 GMSDRV * LZINT = LYINT + L2 GMSDRV * LAST = LZINT + L2 GMSDRV * IF(DIRSCF) THEN GMSDRV * LBUF = LAST GMSDRV * LIBUF = LAST GMSDRV * LGHOND = LAST GMSDRV * LDDIJ = LGHOND + MAXG GMSDRV * LAST = LDDIJ + 16*MXG2 GMSDRV * ELSE GMSDRV * LGHOND = LAST GMSDRV * LDDIJ = LAST GMSDRV * LBUF = LAST GMSDRV * LIBUF = LBUF + NINTMX GMSDRV * LAST = LIBUF + NINTMX GMSDRV * END IF GMSDRV * NEED = LAST - LOADFM GMSDRV * CALL GETFM(NEED) GMSDRV C GMSDRV C Write the perturbed density and MOs to the dictionary file: GMSDRV C GMSDRV * IF(UHF) THEN GMSDRV * ALPHA = .TRUE. GMSDRV * BETA = .FALSE. GMSDRV * CALL FETEDA(A(LSCR)) GMSDRV * CALL DAWRIT(IDAF,IODA,A(LSCR),L3,15,0) GMSDRV * CALL FENEWD(A(LSCR)) GMSDRV * CALL DAWRIT(IDAF,IODA,A(LSCR),L2,16,0) GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .TRUE. GMSDRV * CALL FETEDA(A(LSCR)) GMSDRV * CALL DAWRIT(IDAF,IODA,A(LSCR),L3,19,0) GMSDRV * CALL FENEWD(A(LSCR)) GMSDRV * CALL DAWRIT(IDAF,IODA,A(LSCR),L2,20,0) GMSDRV * ELSE GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .FALSE. GMSDRV * CALL FETEDA(A(LSCR)) GMSDRV * CALL DAWRIT(IDAF,IODA,A(LSCR),L3,15,0) GMSDRV * CALL FENEWD(A(LSCR)) GMSDRV * CALL DAWRIT(IDAF,IODA,A(LSCR),L2,16,0) GMSDRV * END IF GMSDRV C GMSDRV C Pull the perturbed density into memory: GMSDRV C GMSDRV * IF(KFRG.NE.0) THEN GMSDRV * IF(UHF) THEN GMSDRV * ALPHA = .TRUE. GMSDRV * BETA = .FALSE. GMSDRV * CALL FEDDEF(A(LDMDA),NFRG) GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .TRUE. GMSDRV * CALL FEDDEF(A(LDMDB),NFRG) GMSDRV * ELSE GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .FALSE. GMSDRV * CALL FEDDEF(A(LDMDA),NFRG) GMSDRV * END IF GMSDRV * END IF GMSDRV C GMSDRV C Modify one-electron integrals and nuclear repulsion: GMSDRV C GMSDRV * CALL INTS1E(ESNRG,FENRG,0,1) GMSDRV C GMSDRV C Evaluate the energy: GMSDRV C GMSDRV * IF(KFRG.EQ.0) THEN GMSDRV * ITS = 1 GMSDRV * CALL GETEHF(ELOC,DUMMY,ITS) GMSDRV * ELSE GMSDRV * CALL GETEHF(EDEF,ECP,MAXIT) GMSDRV * EFLD = ECP GMSDRV * EDFF = ECP GMSDRV C GMSDRV C An electric field? GMSDRV C GMSDRV * IF(IFFLD.NE.0) THEN GMSDRV C GMSDRV C Save the optimized density on the DAF: GMSDRV C GMSDRV * IF(KFRG.EQ.1) THEN GMSDRV * IF(UHF) THEN GMSDRV * ALPHA = .TRUE. GMSDRV * BETA = .FALSE. GMSDRV * CALL DAREAD(IDAF,IODA,A(LDMFA),L2,16,0) GMSDRV * CALL SVDFLD(A(LDMFA),NFRG) GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .TRUE. GMSDRV * CALL DAREAD(IDAF,IODA,A(LDMFB),L2,20,0) GMSDRV * CALL SVDFLD(A(LDMFB),NFRG) GMSDRV * ELSE GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .FALSE. GMSDRV * CALL DAREAD(IDAF,IODA,A(LDMFA),L2,16,0) GMSDRV * CALL SVDFLD(A(LDMFA),NFRG) GMSDRV * END IF GMSDRV * ELSE GMSDRV * IF(UHF) THEN GMSDRV * ALPHA = .TRUE. GMSDRV * BETA = .FALSE. GMSDRV * CALL FEDFLD(A(LDMFA),NFRG) GMSDRV * CALL DAREAD(IDAF,IODA,A(LDMFA+IOFF),L2,16,0) GMSDRV * CALL SVDFLD(A(LDMFA),NFRG) GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .TRUE. GMSDRV * CALL FEDFLD(A(LDMFB),NFRG) GMSDRV * CALL DAREAD(IDAF,IODA,A(LDMFB+IOFF),L2,20,0) GMSDRV * CALL SVDFLD(A(LDMFB),NFRG) GMSDRV * ELSE GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .FALSE. GMSDRV * CALL FEDFLD(A(LDMFA),NFRG) GMSDRV * CALL DAREAD(IDAF,IODA,A(LDMFA+IOFF),L2,16,0) GMSDRV * CALL SVDFLD(A(LDMFA),NFRG) GMSDRV * END IF GMSDRV * END IF GMSDRV C GMSDRV C Modify one-electron integrals and nuclear repulsion: GMSDRV C GMSDRV * CALL INTS1E(DUMMY,DUMMY,0,0) GMSDRV C GMSDRV C Evaluate fragment energy in absence of field using field optimized GMSDRV C density: GMSDRV C GMSDRV * ITS = 1 GMSDRV * CALL GETEHF(EDFF,DUMMY,ITS) GMSDRV C GMSDRV C Rewrite the perturbed density to the dictionary file: GMSDRV C GMSDRV * IF(UHF) THEN GMSDRV * ALPHA = .TRUE. GMSDRV * BETA = .FALSE. GMSDRV * CALL FENEWD(A(LSCR)) GMSDRV * CALL DAWRIT(IDAF,IODA,A(LSCR),L2,16,0) GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .TRUE. GMSDRV * CALL FENEWD(A(LSCR)) GMSDRV * CALL DAWRIT(IDAF,IODA,A(LSCR),L2,20,0) GMSDRV * ELSE GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .FALSE. GMSDRV * CALL FENEWD(A(LSCR)) GMSDRV * CALL DAWRIT(IDAF,IODA,A(LSCR),L2,16,0) GMSDRV * END IF GMSDRV C GMSDRV C Evaluate the energies and densities of the fragment in the absence GMSDRV C of the field: GMSDRV C GMSDRV * CALL GETEHF(EDEF,ECP,MAXIT) GMSDRV * END IF GMSDRV C GMSDRV C Save the optimized, zero field density on the DAF: GMSDRV C GMSDRV * IF(KFRG.EQ.1) THEN GMSDRV * IF(UHF) THEN GMSDRV * ALPHA = .TRUE. GMSDRV * BETA = .FALSE. GMSDRV * CALL DAREAD(IDAF,IODA,A(LDMCA),L2,16,0) GMSDRV * CALL SVDCP(A(LDMCA),NFRG) GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .TRUE. GMSDRV * CALL DAREAD(IDAF,IODA,A(LDMCB),L2,20,0) GMSDRV * CALL SVDCP(A(LDMCB),NFRG) GMSDRV * ELSE GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .FALSE. GMSDRV * CALL DAREAD(IDAF,IODA,A(LDMCA),L2,16,0) GMSDRV * CALL SVDCP(A(LDMCA),NFRG) GMSDRV * END IF GMSDRV * ELSE GMSDRV * IF(UHF) THEN GMSDRV * ALPHA = .TRUE. GMSDRV * BETA = .FALSE. GMSDRV * CALL FEDCP(A(LDMCA),NFRG) GMSDRV * CALL DAREAD(IDAF,IODA,A(LDMCA+IOFF),L2,16,0) GMSDRV * CALL SVDCP(A(LDMCA),NFRG) GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .TRUE. GMSDRV * CALL FEDCP(A(LDMCB),NFRG) GMSDRV * CALL DAREAD(IDAF,IODA,A(LDMCB+IOFF),L2,20,0) GMSDRV * CALL SVDCP(A(LDMCB),NFRG) GMSDRV * ELSE GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .FALSE. GMSDRV * CALL FEDCP(A(LDMCA),NFRG) GMSDRV * CALL DAREAD(IDAF,IODA,A(LDMCA+IOFF),L2,16,0) GMSDRV * CALL SVDCP(A(LDMCA),NFRG) GMSDRV * END IF GMSDRV * END IF GMSDRV C GMSDRV C Compute dipole moments: GMSDRV C GMSDRV * CALL FECOOR(A(LSCR)) GMSDRV * COM(1) = ZERO GMSDRV * COM(2) = ZERO GMSDRV * COM(3) = ZERO GMSDRV * TOTM = ZERO GMSDRV * TOTQ = -FLOAT(NFE(KFRG)) - FLOAT(NFEB(KFRG))*IFUHF GMSDRV * II = LSCR - 4 GMSDRV * DO 10 IAT = 1,NAT GMSDRV * II = II + 3 GMSDRV * IF(ZAN(IAT).NE.ZERO) THEN GMSDRV * TOTM = TOTM + ZMASS(IAT) GMSDRV * TOTQ = TOTQ + ZAN(IAT) GMSDRV * COM(1) = COM(1) + A(II+1)*ZMASS(IAT) GMSDRV * COM(2) = COM(2) + A(II+2)*ZMASS(IAT) GMSDRV * COM(3) = COM(3) + A(II+3)*ZMASS(IAT) GMSDRV * END IF GMSDRV * 10 CONTINUE GMSDRV * COM(1) = TOTQ * COM(1) / TOTM GMSDRV * COM(2) = TOTQ * COM(2) / TOTM GMSDRV * COM(3) = TOTQ * COM(3) / TOTM GMSDRV C GMSDRV C ...electronic contribution: GMSDRV C GMSDRV * DO 20 I = 1,3 GMSDRV * IX = I GMSDRV * CALL FEDXYZ(A(LSCR),IX) GMSDRV * IF(IX.EQ.0) GOTO 60 GMSDRV * CALL PACK(A(LSCR),L1,L1) GMSDRV * DDEF(I) = -TRACEP(A(LDMDA+IOFF),A(LSCR),L1) GMSDRV * DCP(I) = -TRACEP(A(LDMCA+IOFF),A(LSCR),L1) GMSDRV * IF(UHF) DDEF(I) = DDEF(I) - TRACEP(A(LDMDB+IOFF),A(LSCR),L1) GMSDRV * IF(UHF) DCP(I) = DCP(I) - TRACEP(A(LDMCB+IOFF),A(LSCR),L1) GMSDRV * IF(IFFLD.NE.0) THEN GMSDRV * DFLD(I) = -TRACEP(A(LDMFA+IOFF),A(LSCR),L1) GMSDRV * IF(UHF) DFLD(I) = DFLD(I) - TRACEP(A(LDMFB+IOFF),A(LSCR),L1)GMSDRV * END IF GMSDRV * 20 CONTINUE GMSDRV C GMSDRV C ...nuclear contribution: GMSDRV C GMSDRV * CALL FECOOR(A(LSCR)) GMSDRV * II = LSCR - 1 GMSDRV * DO 40 IAT = 1,NAT GMSDRV * DO 30 I = 1,3 GMSDRV * II = II + 1 GMSDRV * DDEF(I) = DDEF(I) + ZAN(IAT) * A(II) GMSDRV * DCP(I) = DCP(I) + ZAN(IAT) * A(II) GMSDRV * IF(IFFLD.NE.0) THEN GMSDRV * DFLD(I) = DFLD(I) + ZAN(IAT) * A(II) GMSDRV * END IF GMSDRV * 30 CONTINUE GMSDRV * 40 CONTINUE GMSDRV C GMSDRV C Convert the dipole to Debye: GMSDRV C GMSDRV * DO 50 I = 1,3 GMSDRV * DDEF(I) = (DDEF(I) - COM(I)) * DEBYE GMSDRV * DFLD(I) = (DFLD(I) - COM(I)) * DEBYE GMSDRV * DCP(I) = (DCP(I) - COM(I)) * DEBYE GMSDRV * DIND(I) = DDEF(I) - DCP(I) GMSDRV * 50 CONTINUE GMSDRV * DDT = SQRT(DDEF(1)*DDEF(1) + DDEF(2)*DDEF(2) + DDEF(3)*DDEF(3)) GMSDRV * DFT = SQRT(DFLD(1)*DFLD(1) + DFLD(2)*DFLD(2) + DFLD(3)*DFLD(3)) GMSDRV * DCT = SQRT( DCP(1)*DCP(1) + DCP(2)*DCP(2) + DCP(3)*DCP(3) ) GMSDRV * DIT = SQRT(DIND(1)*DIND(1) + DIND(2)*DIND(2) + DIND(3)*DIND(3)) GMSDRV * WRITE(LFNPR,940) (DDEF(I),I=1,3),DDT GMSDRV * IF(IFFLD.NE.0) WRITE(LFNPR,950) (DFLD(I),I=1,3),DFT GMSDRV * WRITE(LFNPR,960) (DCP(I), I=1,3),DCT GMSDRV * WRITE(LFNPR,970) (DIND(I),I=1,3),DIT GMSDRV * 60 CONTINUE GMSDRV C GMSDRV C Restore one-electron integrals and nuclear repulsion: GMSDRV C GMSDRV * CALL INTS1E(DUMMY,DUMMY,1,1) GMSDRV * END IF GMSDRV C GMSDRV C Evaluate electrostatic, polarization, exchange, and self-energy GMSDRV C components: GMSDRV C GMSDRV * PLNRG = ZERO GMSDRV * EXNRG = ZERO GMSDRV * DO 70 I = 1,NFRG GMSDRV * SNRG(I) = ZERO GMSDRV * FPNRG(I) = ZERO GMSDRV * FSNRG(I) = ZERO GMSDRV * 70 CONTINUE GMSDRV * IF(NFRG.GT.0.AND.KFRG.EQ.NFRG) THEN GMSDRV * IF(IGROUP.EQ.1.AND.NOPK.EQ.1) THEN GMSDRV * IF(UHF) THEN GMSDRV * ALPHA = .TRUE. GMSDRV * BETA = .FALSE. GMSDRV * CALL FEDCP(A(LDMCA),NFRG) GMSDRV * CALL FEDDEF(A(LDMDA),NFRG) GMSDRV * IF(IFFLD.NE.0) CALL FEDFLD(A(LDMFA),NFRG) GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .TRUE. GMSDRV * CALL FEDCP(A(LDMCB),NFRG) GMSDRV * CALL FEDDEF(A(LDMDB),NFRG) GMSDRV * IF(IFFLD.NE.0) CALL FEDFLD(A(LDMFB),NFRG) GMSDRV * ELSE GMSDRV * ALPHA = .FALSE. GMSDRV * BETA = .FALSE. GMSDRV * CALL FEDCP(A(LDMCA),NFRG) GMSDRV * CALL FEDDEF(A(LDMDA),NFRG) GMSDRV * IF(IFFLD.NE.0) CALL FEDFLD(A(LDMFA),NFRG) GMSDRV * END IF GMSDRV * IF(DIRSCF) THEN GMSDRV * CALL DIRESX(ESNRG,PLNRG,EXNRG,SNRG,A(LDMCA),A(LDMCB), GMSDRV * + A(LDMDA),A(LDMDB),A(LGHOND),MAXG,A(LDDIJ)) GMSDRV * ELSE GMSDRV * CALL DSKESX(ESNRG,PLNRG,EXNRG,SNRG,A(LDMCA),A(LDMCB), GMSDRV * + A(LDMDA),A(LDMDB),A(LBUF),A(LIBUF)) GMSDRV * END IF GMSDRV C GMSDRV C Read the nuclear attraction and field integrals from disk: GMSDRV C GMSDRV * CALL FEVNUC(A(LVNUC),NFRG) GMSDRV * IF(EFLDL) THEN GMSDRV * CALL DAREAD(IDAF,IODA,A(LXINT),L2,95,0) GMSDRV * CALL DAREAD(IDAF,IODA,A(LYINT),L2,96,0) GMSDRV * CALL DAREAD(IDAF,IODA,A(LZINT),L2,97,0) GMSDRV * END IF GMSDRV * IF(ISUM.NE.0) THEN GMSDRV * CALL DAREAD(IDAF,IODA,A(LSCR),L2,89,0) GMSDRV * END IF GMSDRV C GMSDRV C Loop over fragments evaluating nuclear attraction terms: GMSDRV C GMSDRV * DO 90 I = 1,NFRG GMSDRV * IOFF = (I - 1) * L2 GMSDRV * DO 80 J = 1,NFRG GMSDRV * JOFF = (J - 1) * L2 GMSDRV * IF(I.NE.J) THEN GMSDRV * TRACE1 = TRACEP(A(LDMCA+IOFF),A(LVNUC+JOFF),L1) GMSDRV * TRACE2 = TRACEP(A(LDMDA+IOFF),A(LVNUC+JOFF),L1) GMSDRV * ESNRG = ESNRG + TRACE1 GMSDRV * PLNRG = PLNRG - TRACE1 + TRACE2 GMSDRV * SNRG(I) = SNRG(I) + TRACE1 - TRACE2 GMSDRV * IF(UHF) THEN GMSDRV * TRACE1 = TRACEP(A(LDMCB+IOFF),A(LVNUC+JOFF),L1) GMSDRV * TRACE2 = TRACEP(A(LDMDB+IOFF),A(LVNUC+JOFF),L1) GMSDRV * ESNRG = ESNRG + TRACE1 GMSDRV * PLNRG = PLNRG - TRACE1 + TRACE2 GMSDRV * SNRG(I) = SNRG(I) + TRACE1 - TRACE2 GMSDRV * END IF GMSDRV * END IF GMSDRV * 80 CONTINUE GMSDRV C GMSDRV C ...and field terms: GMSDRV C GMSDRV * IF(EFLDL) THEN GMSDRV * IF(EVEC(1).NE.ZERO) THEN GMSDRV * TRACE1 = EVEC(1) * TRACEP(A(LDMCA+IOFF),A(LXINT),L1) GMSDRV * TRACE2 = EVEC(1) * TRACEP(A(LDMDA+IOFF),A(LXINT),L1) GMSDRV * TRACE3 = EVEC(1) * TRACEP(A(LDMFA+IOFF),A(LXINT),L1) GMSDRV * ESNRG = ESNRG + TRACE1 GMSDRV * PLNRG = PLNRG - TRACE1 + TRACE2 GMSDRV * SNRG(I) = SNRG(I) + TRACE1 - TRACE2 GMSDRV * FENRG(I) = FENRG(I) + TRACE1 GMSDRV * FPNRG(I) = FPNRG(I) - TRACE1 + TRACE3 GMSDRV * FSNRG(I) = FSNRG(I) + TRACE1 - TRACE3 GMSDRV * IF(UHF) THEN GMSDRV * TRACE1 = EVEC(1) * TRACEP(A(LDMCB+IOFF),A(LXINT),L1) GMSDRV * TRACE2 = EVEC(1) * TRACEP(A(LDMDB+IOFF),A(LXINT),L1) GMSDRV * TRACE3 = EVEC(1) * TRACEP(A(LDMFB+IOFF),A(LXINT),L1) GMSDRV * ESNRG = ESNRG + TRACE1 GMSDRV * PLNRG = PLNRG - TRACE1 + TRACE2 GMSDRV * SNRG(I) = SNRG(I) + TRACE1 - TRACE2 GMSDRV * FENRG(I) = FENRG(I) + TRACE1 GMSDRV * FPNRG(I) = FPNRG(I) - TRACE1 + TRACE3 GMSDRV * FSNRG(I) = FSNRG(I) + TRACE1 - TRACE3 GMSDRV * END IF GMSDRV * END IF GMSDRV * IF(EVEC(2).NE.ZERO) THEN GMSDRV * TRACE1 = EVEC(2) * TRACEP(A(LDMCA+IOFF),A(LYINT),L1) GMSDRV * TRACE2 = EVEC(2) * TRACEP(A(LDMDA+IOFF),A(LYINT),L1) GMSDRV * TRACE3 = EVEC(2) * TRACEP(A(LDMFA+IOFF),A(LYINT),L1) GMSDRV * ESNRG = ESNRG + TRACE1 GMSDRV * PLNRG = PLNRG - TRACE1 + TRACE2 GMSDRV * SNRG(I) = SNRG(I) + TRACE1 - TRACE2 GMSDRV * FENRG(I) = FENRG(I) + TRACE1 GMSDRV * FPNRG(I) = FPNRG(I) - TRACE1 + TRACE3 GMSDRV * FSNRG(I) = FSNRG(I) + TRACE1 - TRACE3 GMSDRV * IF(UHF) THEN GMSDRV * TRACE1 = EVEC(2) * TRACEP(A(LDMCB+IOFF),A(LYINT),L1) GMSDRV * TRACE2 = EVEC(2) * TRACEP(A(LDMDB+IOFF),A(LYINT),L1) GMSDRV * TRACE3 = EVEC(2) * TRACEP(A(LDMFB+IOFF),A(LYINT),L1) GMSDRV * ESNRG = ESNRG + TRACE1 GMSDRV * PLNRG = PLNRG - TRACE1 + TRACE2 GMSDRV * SNRG(I) = SNRG(I) + TRACE1 - TRACE2 GMSDRV * FENRG(I) = FENRG(I) + TRACE1 GMSDRV * FPNRG(I) = FPNRG(I) - TRACE1 + TRACE3 GMSDRV * FSNRG(I) = FSNRG(I) + TRACE1 - TRACE3 GMSDRV * END IF GMSDRV * END IF GMSDRV * IF(EVEC(3).NE.ZERO) THEN GMSDRV * TRACE1 = EVEC(3) * TRACEP(A(LDMCA+IOFF),A(LZINT),L1) GMSDRV * TRACE2 = EVEC(3) * TRACEP(A(LDMDA+IOFF),A(LZINT),L1) GMSDRV * TRACE3 = EVEC(3) * TRACEP(A(LDMFA+IOFF),A(LZINT),L1) GMSDRV * ESNRG = ESNRG + TRACE1 GMSDRV * PLNRG = PLNRG - TRACE1 + TRACE2 GMSDRV * SNRG(I) = SNRG(I) + TRACE1 - TRACE2 GMSDRV * FENRG(I) = FENRG(I) + TRACE1 GMSDRV * FPNRG(I) = FPNRG(I) - TRACE1 + TRACE3 GMSDRV * FSNRG(I) = FSNRG(I) + TRACE1 - TRACE3 GMSDRV * IF(UHF) THEN GMSDRV * TRACE1 = EVEC(3) * TRACEP(A(LDMCB+IOFF),A(LZINT),L1) GMSDRV * TRACE2 = EVEC(3) * TRACEP(A(LDMDB+IOFF),A(LZINT),L1) GMSDRV * TRACE3 = EVEC(3) * TRACEP(A(LDMFB+IOFF),A(LZINT),L1) GMSDRV * ESNRG = ESNRG + TRACE1 GMSDRV * PLNRG = PLNRG - TRACE1 + TRACE2 GMSDRV * SNRG(I) = SNRG(I) + TRACE1 - TRACE2 GMSDRV * FENRG(I) = FENRG(I) + TRACE1 GMSDRV * FPNRG(I) = FPNRG(I) - TRACE1 + TRACE3 GMSDRV * FSNRG(I) = FSNRG(I) + TRACE1 - TRACE3 GMSDRV * END IF GMSDRV * END IF GMSDRV * END IF GMSDRV * IF(ISUM.NE.0) THEN GMSDRV * TRACE1 = TRACEP(A(LDMCA+IOFF),A(LSCR),L1) GMSDRV * TRACE2 = TRACEP(A(LDMDA+IOFF),A(LSCR),L1) GMSDRV * TRACE3 = TRACEP(A(LDMFA+IOFF),A(LSCR),L1) GMSDRV * ESNRG = ESNRG + TRACE1 GMSDRV * PLNRG = PLNRG - TRACE1 + TRACE2 GMSDRV * SNRG(I) = SNRG(I) + TRACE1 - TRACE2 GMSDRV * FENRG(I) = FENRG(I) + TRACE1 GMSDRV * FPNRG(I) = FPNRG(I) - TRACE1 + TRACE3 GMSDRV * FSNRG(I) = FSNRG(I) + TRACE1 - TRACE3 GMSDRV * IF(UHF) THEN GMSDRV * TRACE1 = TRACEP(A(LDMCB+IOFF),A(LSCR),L1) GMSDRV * TRACE2 = TRACEP(A(LDMDB+IOFF),A(LSCR),L1) GMSDRV * TRACE3 = TRACEP(A(LDMFB+IOFF),A(LSCR),L1) GMSDRV * ESNRG = ESNRG + TRACE1 GMSDRV * PLNRG = PLNRG - TRACE1 + TRACE2 GMSDRV * SNRG(I) = SNRG(I) + TRACE1 - TRACE2 GMSDRV * FENRG(I) = FENRG(I) + TRACE1 GMSDRV * FPNRG(I) = FPNRG(I) - TRACE1 + TRACE3 GMSDRV * FSNRG(I) = FSNRG(I) + TRACE1 - TRACE3 GMSDRV * END IF GMSDRV * END IF GMSDRV * SNRG(I) = SNRG(I) / TWO GMSDRV * FSNRG(I) = FSNRG(I) / TWO GMSDRV * 90 CONTINUE GMSDRV * END IF GMSDRV * END IF GMSDRV C GMSDRV C Save energies on DAF: GMSDRV C GMSDRV * CALL SVE0(ELOC) GMSDRV * IF(NFRG.NE.0) THEN GMSDRV * LEN = 5 + 8 * NFRG GMSDRV * IF(KFRG.EQ.0) THEN GMSDRV * A(LSCR) = ELOC GMSDRV * DO 100 I = 1,LEN-1 GMSDRV * A(LSCR+I) = ZERO GMSDRV * 100 CONTINUE GMSDRV * ELSE GMSDRV * CALL NBREAD(A(LSCR),LEN,70) GMSDRV * A(LSCR+4+ NFRG+KFRG) = EDEF GMSDRV * A(LSCR+4+2*NFRG+KFRG) = ECP GMSDRV * A(LSCR+4+3*NFRG+KFRG) = EFLD GMSDRV * A(LSCR+4+4*NFRG+KFRG) = EDFF GMSDRV * END IF GMSDRV * IF(KFRG.EQ.NFRG) THEN GMSDRV * A(LSCR+1) = ESNRG GMSDRV * A(LSCR+2) = PLNRG GMSDRV * A(LSCR+3) = EXNRG GMSDRV * A(LSCR+4) = ZERO GMSDRV * IF(DFTYPE.NE.ANONE) A(LSCR+4) = ONE GMSDRV * DO 110 I = 1,NFRG GMSDRV * A(LSCR+4 +I) = SNRG(I) GMSDRV * A(LSCR+4+5*NFRG+I) = FENRG(I) GMSDRV * A(LSCR+4+6*NFRG+I) = FPNRG(I) GMSDRV * A(LSCR+4+7*NFRG+I) = FSNRG(I) GMSDRV * 110 CONTINUE GMSDRV * END IF GMSDRV * CALL NBWRIT(A(LSCR),LEN,70) GMSDRV * END IF GMSDRV C GMSDRV C Close the NBO DAF: GMSDRV C GMSDRV * DEL = .FALSE. GMSDRV * CALL NBCLOS(DEL) GMSDRV C GMSDRV C Return memory: GMSDRV C GMSDRV * CALL RETFM(NEED) GMSDRV C GMSDRV C Restore GAMESS print level: GMSDRV C GMSDRV * NPRINT = MPRINT GMSDRV * RETURN GMSDRV C GMSDRV * 900 FORMAT(/1X,'NEDA requires NOPK=1 in $INTGRL and NOSYM=', GMSDRV * + '1 in $CONTRL.') GMSDRV * 910 FORMAT(//1X,'--------------',/1X,' Fragment',I3,':',/1X,'-----', GMSDRV * + '---------',/) GMSDRV * 940 FORMAT(/1X,'Dipole (def): ',F9.4,'(x),',F9.4,'(y),',F9.4,'(z);', GMSDRV * + F9.4,'(tot) Debye') GMSDRV * 950 FORMAT(1X,'Dipole (fld): ',F9.4,'(x),',F9.4,'(y),',F9.4,'(z);', GMSDRV * + F9.4,'(tot) Debye') GMSDRV * 960 FORMAT(1X,'Dipole (cp): ',F9.4,'(x),',F9.4,'(y),',F9.4,'(z);', GMSDRV * + F9.4,'(tot) Debye') GMSDRV * 970 FORMAT(1X,'Dipole (ind): ',F9.4,'(x),',F9.4,'(y),',F9.4,'(z);', GMSDRV * + F9.4,'(tot) Debye') GMSDRV * END GMSDRV C***********************************************************************GMSDRV * SUBROUTINE CHKNBO(T) GMSDRV C***********************************************************************GMSDRV C 11-Feb-93 EDG New subroutine GMSDRV C-----------------------------------------------------------------------GMSDRV * IMPLICIT REAL*8 (A-H,O-Z) GMSDRV C GMSDRV C NBO COMMON Block: GMSDRV C GMSDRV * COMMON/NBINFO/ISPIN,NATOMS,NDIM,NBAS,NNAO,MXBO,MXAO,MXAOLM,MUNIT GMSDRV C GMSDRV C GAMESS COMMON Block: GMSDRV C GMSDRV * COMMON/IOFILE/ IR,IW,IP,IS,IPK,IDAF,NAV,IODA(950) GMSDRV C GMSDRV * DIMENSION T(NDIM,NDIM) GMSDRV C GMSDRV * LEN = NDIM * NDIM GMSDRV * ISPIN = 2 GMSDRV * CALL FECHK(T,ISPIN,IT) GMSDRV * IF(IT.NE.0) CALL DAWRIT(IDAF,IODA,T,LEN,15,0) GMSDRV * ISPIN = -2 GMSDRV * CALL FECHK(T,ISPIN,IT) GMSDRV * IF(IT.NE.0) CALL DAWRIT(IDAF,IODA,T,LEN,19,0) GMSDRV * RETURN GMSDRV * END GMSDRV C***********************************************************************GMSDRV * SUBROUTINE INTS1E(ESNRG,FENRG,IFLG,IFLD) GMSDRV C***********************************************************************GMSDRV C 6-Jun-98 EDG Add electric fields GMSDRV C 13-Jun-97 EDG New subroutine GMSDRV C-----------------------------------------------------------------------GMSDRV * IMPLICIT REAL*8 (A-H,O-Z) GMSDRV * DIMENSION FENRG(*) GMSDRV C GMSDRV C NBO common blocks: GMSDRV C GMSDRV * PARAMETER(MAXATM = 200,MAXBAS = 2000) GMSDRV C GMSDRV * COMMON/NBBAS/LABEL(MAXBAS,6),NBOUNI(MAXBAS),NBOTYP(MAXBAS), GMSDRV * + IATNO(MAXBAS),IBXM(MAXBAS),LOCC(2*MAXBAS),ISCR(2*MAXBAS) GMSDRV * COMMON/NBEDA/KFRG,NFRG,IFRG(MAXBAS),NFE(MAXBAS),NFEB(MAXBAS) GMSDRV * COMMON/NBFLAG/ROHF,UHF,CI,OPEN,COMPLX,ALPHA,BETA,MCSCF,AUHF,ORTHO GMSDRV * LOGICAL ROHF,UHF,CI,OPEN,COMPLX,ALPHA,BETA,MCSCF,AUHF,ORTHO GMSDRV * COMMON/NBINFO/ISPIN,NATOMS,NDIM,NBAS,NNAO,MXBO,MXAO,MXAOLM,MUNIT GMSDRV C GMSDRV C GAMESS common blocks: GMSDRV C GMSDRV * PARAMETER (MXATM=2000, MXFRG=1050, MXDFG=5, MXDPPT=MXFRG*MXDFG*12)GMSDRV C GMSDRV * COMMON /CONV / DENTOL,EN,ETOT,EHF,EHF0,DIFF,ITER,ICALP,ICBET GMSDRV * COMMON /ECPDIM/ NCOEF1,NCOEF2,J1LEN,J2LEN,LLIM,NLIM,NTLIM,J4LEN GMSDRV * COMMON /ECP2 / CLP(400),ZLP(400),NLP(400),KFIRST(MXATM,6), GMSDRV * + KLAST(MXATM,6),LMAX(MXATM),LPSKIP(MXATM), GMSDRV * + IZCORE(MXATM) GMSDRV * COMMON /EFLDC / EVEC(3),EFLDL GMSDRV * COMMON /FMCOM / A(1) GMSDRV * COMMON /FRGINF/ NMPTS(MXFRG),NMTTPT,IEFC,IEFD,IEFQ,IEFO, GMSDRV * * NPPTS(MXFRG),NPTTPT,IEFP, GMSDRV * * NRPTS(MXFRG),NRTTPT,IREP,ICHGP,NFRAG, GMSDRV * * NDPPTS(MXDPPT),NDPTTPT,IEFDP GMSDRV * COMMON /INFOA / NAT,ICH,MUL,NUM,NX,NE,NA,NB, GMSDRV * + ZAN(MXATM),C(3,MXATM),IAN(MXATM) GMSDRV * COMMON /INTOPT/ ISCHWZ,IECP,NECP,IEFLD GMSDRV * COMMON /IOFILE/ IR,IW,IP,IS,IPK,IDAF,NAV,IODA(950) GMSDRV * COMMON /MACHIN/ NWDVAR,MAXFM,MAXSM,LIMFM,LIMSM GMSDRV C GMSDRV * LOGICAL EFLDL,EFLAG GMSDRV C GMSDRV * SAVE NEC,NAC,NBC GMSDRV * SAVE ZERO GMSDRV * DATA ZERO/0.0D0/ GMSDRV C GMSDRV C Recalculate one-electron integrals for NEDA: GMSDRV C GMSDRV * IF(KFRG.NE.0) THEN GMSDRV C GMSDRV C Update number of electrons: GMSDRV C GMSDRV * IF(IFLD.NE.0) THEN GMSDRV * IF(IFLG.EQ.0) THEN GMSDRV * NEC = NE GMSDRV * NAC = NA GMSDRV * NBC = NB GMSDRV * IF(UHF) THEN GMSDRV * NA = NFE(KFRG) GMSDRV * NB = NFEB(KFRG) GMSDRV * ELSE GMSDRV * NA = NFE(KFRG) / 2 GMSDRV * NB = NA GMSDRV * END IF GMSDRV * NE = NA + NB GMSDRV * ELSE GMSDRV * NE = NEC GMSDRV * NA = NAC GMSDRV * NB = NBC GMSDRV * END IF GMSDRV C GMSDRV C Update array ZAN of /INFOA/: GMSDRV C GMSDRV * IF(IFLG.EQ.0) THEN GMSDRV * CALL NBWRIT(ZAN,NAT,75) GMSDRV * DO 10 I = 1,NAT GMSDRV * IF(IATNO(I).LT.0) ZAN(I) = ZERO GMSDRV * 10 CONTINUE GMSDRV * ELSE GMSDRV * CALL NBREAD(ZAN,NAT,75) GMSDRV * END IF GMSDRV C GMSDRV C For ECP runs, update arrays of /ECP2/: GMSDRV C GMSDRV * IF(IECP.NE.0) THEN GMSDRV * IF(IFLG.EQ.0) THEN GMSDRV * CALL NBWRIT(IZCORE,NAT,76) GMSDRV * CALL NBWRIT(LPSKIP,NAT,77) GMSDRV * CALL NBWRIT(LMAX,NAT,78) GMSDRV * CALL BASCHK(NLIM) GMSDRV * NLIM = NLIM + 1 GMSDRV * LLIM = 0 GMSDRV * DO 20 I = 1,NAT GMSDRV * IF(IATNO(I).LT.0) THEN GMSDRV * IZCORE(I) = 0 GMSDRV * LPSKIP(I) = 1 GMSDRV * LMAX(I) = 0 GMSDRV * END IF GMSDRV * NECP = NECP + IZCORE(I) GMSDRV * LLIM = MAX(LLIM,LMAX(I)) GMSDRV * 20 CONTINUE GMSDRV * NLIM = MAX(LLIM,NLIM) GMSDRV * ELSE GMSDRV * CALL NBREAD(IZCORE,NAT,76) GMSDRV * CALL NBREAD(LPSKIP,NAT,77) GMSDRV * CALL NBREAD(LMAX,NAT,78) GMSDRV * NECP = 0 GMSDRV * DO 30 I = 1,NAT GMSDRV * NECP = NECP + IZCORE(I) GMSDRV * 30 CONTINUE GMSDRV * CALL BASCHK(NLIM) GMSDRV * NLIM = NLIM + 1 GMSDRV * LLIM = 0 GMSDRV * DO 40 I = 1,NAT GMSDRV * LLIM = MAX(LLIM,LMAX(I)) GMSDRV * 40 CONTINUE GMSDRV * NLIM = MAX(LLIM,NLIM) GMSDRV * END IF GMSDRV C GMSDRV C Allocate memory: GMSDRV C GMSDRV * NCOEF1 = 8520 GMSDRV * NCOEF2 = 3424 GMSDRV * IF(NLIM.EQ.5) THEN GMSDRV * NCOEF1 = 71660 GMSDRV * NCOEF2 = 10555 GMSDRV * ELSE IF(NLIM.EQ.6) THEN GMSDRV * NCOEF1 = 280000 GMSDRV * NCOEF2 = 28940 GMSDRV * ELSE IF(NLIM.EQ.7) THEN GMSDRV * NCOEF1 = 892584 GMSDRV * NCOEF2 = 60382 GMSDRV * END IF GMSDRV * NTLIM = (NLIM*(NLIM+1)*(NLIM+2))/6 GMSDRV * J1LEN = (NTLIM*NTLIM+NTLIM)/2 + 1 GMSDRV * J2LEN = (LLIM-1)*(LLIM+1)*NTLIM+NTLIM + 1 GMSDRV * J4LEN = NTLIM+NTLIM*(LLIM-1)*(LLIM+1) GMSDRV * NUMDER = 0 GMSDRV * CALL VALFM(LOADFM) GMSDRV * LDCF1 = LOADFM + 1 GMSDRV * LJIN = LDCF1 + NCOEF1 GMSDRV * LLB1 = LJIN + (J1LEN-1)/NWDVAR+1 GMSDRV * LDCF2 = LLB1 + (9*NCOEF1)/NWDVAR GMSDRV * LLJ2 = LDCF2 + NCOEF2 GMSDRV * LLB2 = LLJ2 + (J2LEN-1)/NWDVAR+1 GMSDRV * LFPQR = LLB2 + (6*NCOEF2)/NWDVAR GMSDRV * LZLM = LFPQR + 15625 GMSDRV * LLMF = LZLM + 584 GMSDRV * LLMX = LLMF + 124/NWDVAR GMSDRV * LLMY = LLMX + 584/NWDVAR GMSDRV * LLMZ = LLMY + 584/NWDVAR GMSDRV * LAST = LLMZ + 584/NWDVAR GMSDRV * NEED = LAST - LOADFM GMSDRV * CALL GETFM(NEED) GMSDRV C GMSDRV C Recalculate the ECP integral formulas: GMSDRV C GMSDRV * CALL ECCODR(A(LDCF1),A(LJIN),A(LLB1),A(LDCF2),A(LLJ2), GMSDRV * + A(LLB2),A(LFPQR),A(LZLM),A(LLMF),A(LLMX), GMSDRV * + A(LLMY),A(LLMZ),NUMDER) GMSDRV * LDAF91 = (J1LEN-1)/NWDVAR+1 + (9*NCOEF1-1)/NWDVAR+1 GMSDRV * LDAF93 = (J2LEN-1)/NWDVAR+1 + (6*NCOEF2)/NWDVAR GMSDRV * CALL DAWRIT(IDAF,IODA,A(LDCF1),NCOEF1,90,0) GMSDRV * CALL DAWRIT(IDAF,IODA,A(LJIN) ,LDAF91,91,1) GMSDRV * CALL DAWRIT(IDAF,IODA,A(LDCF2),NCOEF2,92,0) GMSDRV * CALL DAWRIT(IDAF,IODA,A(LLJ2) ,LDAF93,93,1) GMSDRV * CALL RETFM(NEED) GMSDRV * END IF GMSDRV * END IF GMSDRV C GMSDRV C Recalculate the one-electron integrals: GMSDRV C GMSDRV * ISUM = IEFC + IEFD + IEFQ + IEFO + IEFP + NFRAG GMSDRV * EFLAG = EFLDL GMSDRV * IEFCO = IEFC GMSDRV * IEFDO = IEFD GMSDRV * IEFQO = IEFQ GMSDRV * IEFOO = IEFO GMSDRV * IEFPO = IEFP GMSDRV * NFRGO = NFRAG GMSDRV * IF(IFLD.EQ.0) THEN GMSDRV * EFLDL = .FALSE. GMSDRV * IEFC = 0 GMSDRV * IEFD = 0 GMSDRV * IEFQ = 0 GMSDRV * IEFO = 0 GMSDRV * IEFP = 0 GMSDRV * NFRAG = 0 GMSDRV * END IF GMSDRV * CALL ONEEI GMSDRV C GMSDRV C Allocate memory: GMSDRV c GMSDRV * IF(IFLG.EQ.0.AND..NOT.EFLDL) THEN GMSDRV * L2 = NDIM * (NDIM + 1) / 2 GMSDRV * CALL VALFM(LOADFM) GMSDRV * LH = LOADFM + 1 GMSDRV * LT = LH + L2 GMSDRV * LV = LT + L2 GMSDRV * LVNUC = LV + L2 GMSDRV * LAST = LVNUC + NFRG*L2 GMSDRV * NEED = LAST - LOADFM GMSDRV * CALL GETFM(NEED) GMSDRV C GMSDRV C Copy the nuclear attraction integrals to the DAF. The kinetic GMSDRV C energy contributions to obtain V: GMSDRV C GMSDRV * CALL DAREAD(IDAF,IODA,A(LH),L2,11,0) GMSDRV * CALL DAREAD(IDAF,IODA,A(LT),L2,13,0) GMSDRV * CALL VSUB(A(LT),1,A(LH),1,A(LV),1,L2) GMSDRV * IOFF = (KFRG - 1) * L2 GMSDRV * IF(KFRG.GT.1) CALL FEVNUC(A(LVNUC),NFRG) GMSDRV * CALL COPY(A(LV),A(LVNUC+IOFF),L2,L2,1) GMSDRV * CALL SVVNUC(A(LVNUC),NFRG) GMSDRV * CALL RETFM(NEED) GMSDRV * END IF GMSDRV C GMSDRV C Recalculate the nuclear repulsion energy: GMSDRV C GMSDRV * EN = ENUC(NAT,ZAN,C) GMSDRV * EFLDL = EFLAG GMSDRV * IEFC = IEFCO GMSDRV * IEFD = IEFDO GMSDRV * IEFQ = IEFQO GMSDRV * IEFO = IEFOO GMSDRV * IEFP = IEFPO GMSDRV * NFRAG = NFRGO GMSDRV * END IF GMSDRV C GMSDRV C Evaluate the nuclear contribution to the electrostatic energy: GMSDRV C GMSDRV * IF(NFRG.NE.0.AND.IFLG.EQ.0.AND.IFLD.EQ.1) THEN GMSDRV * IF(KFRG.EQ.0) THEN GMSDRV * ESNRG = ENUC(NAT,ZAN,C) GMSDRV * ELSE GMSDRV * EFLAG = EFLDL GMSDRV * IEFCO = IEFC GMSDRV * IEFDO = IEFD GMSDRV * IEFQO = IEFQ GMSDRV * IEFOO = IEFO GMSDRV * IEFPO = IEFP GMSDRV * NFRGO = NFRAG GMSDRV * EFLDL = .FALSE. GMSDRV * IEFC = 0 GMSDRV * IEFD = 0 GMSDRV * IEFQ = 0 GMSDRV * IEFO = 0 GMSDRV * IEFP = 0 GMSDRV * NFRAG = 0 GMSDRV * ETMP = ENUC(NAT,ZAN,C) GMSDRV * ESNRG = ESNRG - ETMP GMSDRV * FENRG(KFRG) = EN - ETMP GMSDRV * EFLDL = EFLAG GMSDRV * IEFC = IEFCO GMSDRV * IEFD = IEFDO GMSDRV * IEFQ = IEFQO GMSDRV * IEFO = IEFOO GMSDRV * IEFP = IEFPO GMSDRV * NFRAG = NFRGO GMSDRV * END IF GMSDRV * END IF GMSDRV * RETURN GMSDRV * END GMSDRV C***********************************************************************GMSDRV * SUBROUTINE GETEHF(EDEF,ECP,ITS) GMSDRV C***********************************************************************GMSDRV C 16-Jun-97 EDG New subroutine GMSDRV C-----------------------------------------------------------------------GMSDRV * IMPLICIT REAL*8 (A-H,O-Z) GMSDRV * DIMENSION IJ(12) GMSDRV C GMSDRV C NBO common blocks: GMSDRV C GMSDRV * COMMON/NBFLAG/ROHF,UHF,CI,OPEN,COMPLX,ALPHA,BETA,MCSCF,AUHF,ORTHO GMSDRV * LOGICAL ROHF,UHF,CI,OPEN,COMPLX,ALPHA,BETA,MCSCF,AUHF,ORTHO GMSDRV C GMSDRV C GAMESS common blocks: GMSDRV C GMSDRV * PARAMETER(MXATM=2000, MXRT=100) GMSDRV C GMSDRV * COMMON /ATHEAT/ ATHEAT GMSDRV * COMMON /CONV / DENTOL,EN,ETOT,EHF,EHF0,DIFF,ITER,ICALP,ICBET GMSDRV * COMMON /DFGRID/ DFTTHR,DFTGTHR,SWOFF,SW0,NDFTFG, GMSDRV * + NRAD,NTHE,NPHI,NRAD0,NTHE0,NPHI0 GMSDRV * COMMON /ENRGYS/ ENUCR,EELCT,ET,SZ,SZZ,ECORE,ESCF,EERD,E1,E2, GMSDRV * + VEN,VEE,EPOT,EKIN,ESTATE(MXRT),STATN,EDFT(2) GMSDRV * COMMON /FUNCT / E,EG(3*MXATM) GMSDRV * COMMON /ENEDA / E0 GMSDRV * COMMON /RUNOPT/ RUNTYP,EXETYP,NEVALS,NGLEVL,NHLEVL GMSDRV * COMMON /SCFOPT/ CONVHF,MAXIT,MCONV,NPUNCH,NPREO(4) GMSDRV * COMMON /WFNOPT/ SCFTYP,CITYP,DFTYPE,CCTYP,MPLEVL,MPCTYP GMSDRV C GMSDRV * SAVE NONE GMSDRV * SAVE ZERO,TOKCAL,TRUDGE GMSDRV * DATA NONE/4HNONE/ GMSDRV * DATA ZERO,TOKCAL/0.0D0,627.52791D0/ GMSDRV * DATA TRUDGE/8HTRUDGE / GMSDRV C GMSDRV C If grid-based DFT, use only fine grid, and do not switch to HF: GMSDRV C GMSDRV * IF(NDFTFG.EQ.1) THEN GMSDRV * NR = NRAD0 GMSDRV * NT = NTHE0 GMSDRV * NP = NPHI0 GMSDRV * NRAD0 = NRAD GMSDRV * NTHE0 = NTHE GMSDRV * NPHI0 = NPHI GMSDRV * XSW = SWOFF GMSDRV * SWOFF = ZERO GMSDRV * END IF GMSDRV C GMSDRV C Evaluate energy: GMSDRV C GMSDRV * MAXIT0 = MAXIT GMSDRV * MAXIT = ITS GMSDRV * NEVALS = 0 GMSDRV * IF(UHF) THEN GMSDRV * CALL UHFOP(SZ,S2) GMSDRV * ELSE GMSDRV * CALL RHFCL GMSDRV * END IF GMSDRV * MAXIT = MAXIT0 GMSDRV * EDEF = E0 GMSDRV * ECP = ETOT GMSDRV C GMSDRV C Reset DFT grid parameters: GMSDRV C GMSDRV * IF(NDFTFG.EQ.1) THEN GMSDRV * NRAD0 = NR GMSDRV * NTHE0 = NT GMSDRV * NPHI0 = NP GMSDRV * SWOFF = XSW GMSDRV * END IF GMSDRV C GMSDRV C Convert energy to kcal/mol for MOPAC wavefunctions: GMSDRV C GMSDRV * IF(MPCTYP.NE.NONE) THEN GMSDRV * E0 = E0 * TOKCAL + ATHEAT GMSDRV * CALL NBREAD(IJ,12,3) GMSDRV * IJ(4) = 2 GMSDRV * CALL NBWRIT(IJ,12,3) GMSDRV * END IF GMSDRV C GMSDRV C Save the energy in common /ENRGYS/ and in the dictionary file GMSDRV C for 'TRUDGE' calculations: GMSDRV C GMSDRV * IF(RUNTYP.EQ.TRUDGE) THEN GMSDRV * EELCT = EHF GMSDRV * ESCF = E0 GMSDRV * ET = E0 GMSDRV * E = E0 GMSDRV *