C
C $Id: setup.F,v 1.4 1998/07/16 16:40:44 jjv5 Exp $
C
C------------------------------------------------------------------------
      SUBROUTINE SETUP(IERROR)
C
C     DOES ALL THE PRELIMINARY SETUP TO PREPARE THE PROGRAM FOR THE
C     FIRST ENERGY EVALUATION.
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
#include "divcon.dim"
#include "divcon.h"
C
      IERROR = 0

C
C     READ AND PROCESS ALL INFORMATION SUPPLIED IN FILES.
C
      CALL RDALL(IERROR)
      IF(IERROR.NE.0) RETURN
C
C     CHECK INTER-ATOMIC DISTANCES FOR EACH RESIDUE IF REQUESTED.
C
      IF(INDEX(KEYWRD,'CHKRES').NE.0.AND.NRES.GT.0) CALL CHKRES
C
C     DETERMINE TOTAL NUMBER OF ELECTRONS AND ORBITALS IN SYSTEM.
C
!     IICHG = INDEX(KEYWRD,'CHARGE=')
!     IF(IICHG.EQ.0)THEN
!       CHG = 0.0D0
!     ELSE
!       CALL RDNUM(KEYWRD,IICHG,CHG,IERR)
!     ENDIF
      NORBS = 0
      NETCHG = NETCHARGE
      NZCHG = 0
      DO 20 I=1,NATOMS
        IAI = IATNUM(I)
        IF(IAI.EQ.0) GO TO 20
        NZCHG = NZCHG + ZCHG(IAI)
        NORBS = NORBS + NATORB(IAI)
 20   CONTINUE
      NELEC = NZCHG - NETCHG
      if (setch) then
C-RDC          write(iout,'(/)')
         do icr=1,ncores
            nzchg1 = 0
            i0 = icorel1(icr)-1
            nresn = icorel1(icr+1)-icorel1(icr)
            DO IR=1,NRESN
               ires = icorel(ir+i0)
               i1 = irpnt(ires)
               i2 = irpnt(ires+1)-1
               do i=i1,i2
                  iai = iatnum(i)
                  if (iai.ne.0) then
                     nzchg1 = nzchg1 + zchg(iai)
                  endif
               enddo
            enddo
            nelecef(icr) = nzchg1 - nelecef(icr)
C-RDC             write(iout,'(" NUMBER OF ELECTRONS IN CLUSTER GROUP ",
C-RDC      &           I4," IS ",i5)') icr, nelecef(icr)
         enddo
      else
         nelecef(1) = nelec
      endif

C
C-RDC       WRITE(IOUT, '(//" NUMBER OF ATOMS                 = ",I5,
C-RDC      .     /" TOTAL NUMBER OF BASIS FUNCTIONS = ",I5,
C-RDC      .     /" TOTAL NUMBER OF ELECTRONS       = ",I5)') NATOMS,
C-RDC      .     NORBS,
C-RDC      .     NELEC
      IF(MOD(NELEC,2).NE.0)THEN
         IERROR = 1
         WRITE(IOUT,'(/" ODD NUMBER OF ELECTRONS IS NOT ALLOWED FOR",
     .        " A CLOSED-SHELL CALCULATION")')
      ENDIF
      IF(NELEC.GT.2*NORBS)THEN
         IERROR = 1
         WRITE(IOUT,'(/" NUMBER OF ELECTRONS IS MORE THAN TWICE THE",
     .        " NUMBER OF ATOMIC ORBITALS")')
      ENDIF
C
C     CHECK STORAGE FOR PAIRWISE QUANTITIES IF 'DIRECT' KEYWORD HAS
C     NOT BEEN SPECIFIED.
C
      IF(.not.DIRECT)THEN
         NPAIRS = (NATOMS*(NATOMS-1))/2
         IF(NPAIRS.GT.MXPAIR)THEN
            IERROR = 1
            WRITE(IOUT,'(/" STORAGE FOR ALL POSSIBLE ATOM PAIRS ",
     .           "EXCEEDED",
     .           /" -- INCREASE MXPAIR PARAMETER IN divcon.dim",
     .           " TO AT LEAST ",I7,
     .           /" OR USE ''DIRECT'' KEYWORD")') NPAIRS
         ENDIF
      ENDIF
  
      IF(IERROR.NE.0) RETURN
C
C     ASSIGN DEFAULT PROGRAM SETTINGS FOR ANYTHING THE USER HAS NOT
C     SPECIFIED.
C
      CALL SETDEF
C
C     READIN EXTERNAL PARAMETERS
      if (extrn) then
         inf=iunit(26)
         open(unit=inf,file=fname(26),status='old',err=1000)
         call intpr(inf,ierror)
         close(inf)
     
         if (ierror.ne.0) return
      endif
C
C     ASSIGN SEMIEMPIRICAL PARAMETERS BASED ON THE CHOSEN HAMILTONIAN.
C     these are the intramolecular contributions, needed to build
C     up the 1 electron part of the fock Matrix, so getpar needs
C     to be called in any case 
C
      CALL GETPAR
      if (prtpar) call printpar
C
C     check if some parameters are missing (based on Uss)
C
      do i = 1, natoms
        iat = iatnum(i)
        if (UCORE(0,IAT).eq.0.0d0) then
           write(iout,'("")')
           write(iout,'(" ERROR: NO PARAMETER FOR ATOM ",i6,
     .               " (ATOMIC NUMBER =",i3,")")') i, iat
           call mexit(iout,1)
        end if
      end do
C
C     IDENTIFY PEPTIDE LINKAGES FOR MOLECULAR MECHANICS CORRECTION IF
C     USER HAS REQUESTED IT.
C
      N2PEP = 0
!      IF(INDEX(KEYWRD,'ADDMM').NE.0) CALL GETPEP(IERROR)
!      CALL GETPEP(IERROR)
      RETURN
!1000 continue
 1000 write(iout,'(/" ERROR: COULD NOT OPEN ",A20," FOR READING ",
     &     "OF EXTERNAL PARAMETERS",
     &     /"        (NOTE THAT THE NAME OF THE EXTERNAL FILE ",
     &     "SHOULD BE CAPITALIZED)")') fname(26)
      ierror = 1
      return
      END