( 2 Jun 94)
                      *                               *
                      * Section 2 - Input Description *
                      *                               *
              This section of the manual describes the input to
          GAMESS.  The section is written in a reference, rather
          than tutorial fashion.  However, there are frequent
          reminders that more information can be found on a
          particular input group, or type of calculation, in the
          'Further Information' section of this manual.  There are
          also a number of examples shown in the 'Input Examples'
              The order of this section is chosen to approximate the
          order in which most people prepare their input ($CONTRL,
          $BASIS/$DATA, $GUESS, and so on).  After that comes run
          type related input, then properties input, integral format
          input, and finally CI/MCSCF input.  The next page contains
          a list of all possible input groups, in the order in which
          they can be found in this section.
           name    function                          module:routine
           ----    --------                          --------------
          $CONTRL  chemical control data             INPUTA:START
          $SYSTEM  computer related control data     INPUTA:START
          $BASIS   basis set                         INPUTB:BASISS
          $DATA    molecule, basis set               INPUTB:MOLE
          $ZMAT    coded z-matrix                    ZMATRX:ZMATIN
          $LIBE    linear bend data                  ZMATRX:LIBE
          $SCF     HF-SCF wavefunction control       SCFLIB:SCFIN
          $MP2     2nd order Moller-Plesset          MP2   :MP2INP
          $GUESS   initial orbital selection         GUESS :GUESMO
          $VEC     orbitals              (formatted) GUESS :READMO
          $STATPT  geometry search control           STATPT:SETSIG
          $TRUDGE  nongradient optimization          TRUDGE:TRUINP
          $TRURST  restart data for TRUDGE           TRUDGE:TRUDGX
          $FORCE   hessian, normal coordinates       HESS  :HESSX
          $HESS    force constant matrix (formatted) HESS  :FCMIN
          $GRAD    gradient vector       (formatted) HESS  :EGIN
          $DIPDR   dipole deriv. matrix  (formatted) HESS  :DDMIN
          $VIB     HESSIAN restart data  (formatted) HESS  :HSSNUM
          $MASS    isotope selection                 VIBANL:RAMS
          $IRC     reaction path                     RXNCRD:IRCX
          $FFCALC  finite electric field             FFIELD:FFLDX
          $LOCAL   Boys localization control         LOCAL :LMOINP
          $ELMOM   electrostatic moments             PRPLIB:INPELM
          $ELPOT   electrostatic potential           PRPLIB:INPELP
          $ELDENS  electron density                  PRPLIB:INPELD
          $EFIELD  electric field/gradient           PRPLIB:INPELF
          $POINTS  property calculation points       PRPLIB:INPPGS
          $GRID    property calculation mesh         PRPLIB:INPPGS
          $PDC     MEP fitting mesh                  PRPLIB:INPPDC
          $STONE   distributed multipole analysis    PRPPOP:STNRD
          $SCRF    self consistent reaction field    SCRF  :ZRFINP
          $ECP     effective core pot's              ECPLIB:ECPPAR
          $INTGRL  format for 2e- integrals          INPUTA:START
          $CIINP   control of CI process             GAMESS:MCCI
          $DRT     distinct row table                GUGDRT:ORDORB
          $MCSCF   parameters for MCSCF              MCSCF :MCSCF
          $TRANS   integral transformation           TRFIN :TRANS
          $CISORT  integral sorting                  GUGSRT:GUGSRT
          $GUGEM   Hamiltonian matrix formation      GUGEM :GUGAEM
          $GUGDIA  Hamiltonian eigenvalues/vectors   GUGDGA:GUGADG
          $GUGDM   1e- density matrix                GUGDM :GUGADM
          $LAGRAN  CI lagrangian matrix              LAGRAN:CILGRN
          $GUGDM2  2e- density matrix                GUGDM2:GUG2DM
          $TRFDM2  2e- density back-transformation   TRFDM2:TRF2DM
          $TRANST  transition moments, spin-orbit    TRNSTN:TRNSTX
          * this information is more useful to a programmer than to
          the casual user.
          $CONTRL group          (optional)
          This is a free format group specifying global switches.
          SCFTYP = RHF       Restricted Hartree Fock calculation
                             (default for EVEN number of electrons)
                 = UHF       Unrestricted Hartree Fock calculation
                             (default for ODD number of electrons)
                 = ROHF      Restricted open shell Hartree-Fock.
                             (high spin, see GVB for low spin)
                 = GVB       Generalized valence bond wavefunction
                             or OCBSE type ROHF. (needs $SCF input)
                 = MCSCF     Multiconfigurational SCF wavefunction
                             (this requires $DRT input)
                 = CI        Configuration Interaction calculation
                             (this requires $DRT input)
          RUNTYP = ENERGY    Single point energy. (default)
                 = GRADIENT  Single point energy plus gradient.
                 = OPTIMIZE  Optimize the molecular geometry using
                             analytic energy gradients. See $STATPT.
                 = TRUDGE    Non-gradient total energy minimization.
                             See groups $TRUDGE and $TRURST.
                 = SADPOINT  Locate saddle point (transition state).
                             See the $STATPT group.
                 = HESSIAN   Compute energy second derivatives, and
                             perform harmonic vibrational analysis.
                             See the $FORCE group.
                 = IRC       Follow intrinsic reaction coordinate.
                             See the $IRC group.
                 = PROP      Properties will be calculated.  A $DATA
                             deck and converged $VEC group should be
                             input.  Optionally, Boys localization
                             can be done.  See $ELPOT, etc.
                 = TRANSITN  Find radiative transition moment.
                             See the $TRANST group.
                 = SPINORBT  Find spin-orbit coupling constant.
                             See the $TRANST group.
                 = FFIELD    applies finite electric fields, most
                             commonly to extract polarizabilities.
                             See the $FFCALC group.
          MPLEVL = n         chooses Moller-Plesset perturbation
                             level.  The default is 0 to skip.
                             This is implemented only for n=2,
                             only for RHF, UHF, and ROHF wave
                             functions, and only for ENERGY,
                             TRUDGE, and FFIELD runs.
              * * * * * * * * * * * * * * * * * * * * * * * * *
              Note that RUNTYPs involving the energy gradient
              cannot be used with SCFTYP=CI, or nonzero MPLEVL.
              * * * * * * * * * * * * * * * * * * * * * * * * *
          EXETYP = RUN       Actually do the run. (default)
                 = CHECK     Wavefunction and energy will not be
                             evaluated.  This lets you speedily
                             check input and memory requirements.
                             See the overview section for details.
                 = DEBUG     Massive amounts of output are printed,
                             useful only if you hate trees.
                 = routine   Maximum output is generated by the
                             routine named.  Check the source for
                             the routines this applies to.
          MAXIT  =           Maximum number of SCF iteration cycles.
                             Pertains only to RHF, UHF, ROHF, or
                             GVB runs.  See also MAXIT in $MCSCF.
                             (default = 30)
                           * * * * * * *
          ICHARG =           Molecular charge.  (default=0, neutral)
          MULT   =           Multiplicity of the electronic state
                 = 1         singlet (default)
                 = 2,3,...   doublet, triplet, and so on.
             ICHARG and MULT are used directly for RHF, UHF, ROHF.
             For GVB, these are implicit in the $SCF input, while
             for MCSCF or CI, these are implicit in the $DRT input.
             However, you must still give them correctly here.
                           * * * * * * *
          ECP    =           effective core potential control.
                 = NONE      all electron calculation (default).
                 = READ      read the potentials in $ECP group.
                 = SBK       use Stevens, Basch, Krauss, Jasien,
                             Cundari potentials for all heavy 
                             atoms (Li-Rn are available).
                 = HW        use Hay, Wadt potentials for all the
                             heavy atoms (Na-Xe are available).
                           * * * * * * *
            * * * the next three control molecular geometry * * *
          COORD  =           choice for molecular geometry in $DATA.
                 = UNIQUE    only the symmetry unique atoms will be
                             given, in Cartesian coords (default).
                 = HINT      only the symmetry unique atoms will be
                             given, in Hilderbrandt style internals.
                 = CART      Cartesian coordinates will be input.
                 = ZMT       GAUSSIAN style internals will be input.
                 = ZMTMPC    MOPAC style internals will be input.
                 Note that the final three choices require the input
                 of all atoms in the molecule.  GAMESS will orient
                 the molecule, and determine which atoms are unique.
                 The reorientation is likely to change the order of
                 the atoms from what you input.
                 Note that the final three choices require the use
                 of $BASIS to define the basis set.  The first two
                 choices may or may not use $BASIS, as you wish.
          UNITS  =           distance units, any angles which are
                             entered in $DATA must be in degrees.
                 = ANGS      Angstroms (default)
                 = BOHR      Bohr atomic units
          NZVAR  =           Coordinate switch.
                 = 0         Use Cartesian coordinates (default).
                 = M         If COORD=ZMT or ZMTMPC and a $ZMAT is
                             not given:  the internal coordinates
                             will be those defining the molecule in
                             $DATA.  In this case, $DATA must not
                             contain any dummy atoms.  M is usually
                             3N-6 (or 3N-5 for linear).
                 = M         For other COORD choices, or if $ZMAT is
                             given:  the internal coordinates will
                             be those defined in $ZMAT.  This allows
                             the use of more sophisticated internal
                             coordinate choices.  M is ordinarily
                             3N-6 (3N-5), unless linear bends are
                             used in the $ZMAT.
                 Note that NZVAR refers mainly to the coordinates in
                 which an OPTIMIZE or SADPOINT run is performed, but
                 also to the values printed for any of the other run
                 types.  It is possible to use internals to enter
                 the molecule, but still to use Cartesians during
                 the optimization!
          LOCAL  =           controls orbital localization.
                 = NONE      Skip localization (default).
                 = BOYS      Do Foster-Boys localization.
                 = RUEDNBRG  Do Edmiston-Ruedenberg localization.
                 = POP       Do Pipek-Mezey population localization.
                             See the $LOCAL group.   Localization
                             does not work for SCFTYP's GVB or CI.

                 * * * interfaces to other programs * * *
          MOLPLT = flag that produces an input deck for a molecule
                   drawing program distributed with GAMESS.
                   (default is .FALSE.)
          PLTORB = flag that produces an input deck for an orbital
                   plotting program distributed with GAMESS.
                   (default is .FALSE.)
          AIMPAC = flag to create an input deck for Bader's atoms
                   in molecules properties code. (default=.FALSE.) 
                   For information about this program, contact 
                       Richard F.W. Bader
                       Dept. of Chemistry
                       McMaster University
                       Hamilton, Ontario  L8S-4M1 Canada
          RPAC   = flag to create the input files for Bouman and
                   Hansen's RPAC electronic excitation and NMR
                   shieldings program.  RPAC works only with
                   RHF wavefunctions.  Contact Prof. Aage Hansen
                   in Copenhagen (nahaeh@vm.uni-c.dk) about this 
                   program.  (default is .FALSE.)

          FRIEND = string to prepare input to other quantum
                   programs, choose from
                 = HONDO    for HONDO 8.2
                 = MELDF    for MELDF
                 = GAMESSUK for GAMESS (UK Daresbury version)
                 = GAUSSIAN for Gaussian 9x
                 = ALL      for all of the above

          PLTORB, MOLPLT, and AIMPAC decks are written to file
          PUNCH at the end of the job.  The two binary disk
          files output by RPAC are written at the end of the
          job.  Thus all of these correspond to the final
          geometry encountered during the job.  


          In contrast, selecting FRIEND turns the job into a 
          CHECK run only, no matter how you set EXETYP.  Thus the 
          geometry is that encountered in $DATA.  The input is 
          added to the PUNCH file, and may require some (usually 
          minimal) massaging.
          PLTORB and MOLPLT are written even for EXETYP=CHECK.
          AIMPAC requires at least RUNTYP=PROP.  RPAC requires at
          least RUNTYP=ENERGY, and you must take action to save
          the binary files AOINTS and WORK15.

             The NBO program of Frank Weinhold's group can be
          attached to GAMESS.  The input to control the natural
          bond order analysis is read by the add in code, so is
          not described here.  For information on the NBO program,
          contact Frank Weinhold       WEINHOLD@CHEM.WISC.EDU
                  Theoretical Chemistry Institute
                  Department of Chemistry
                  University of Wisconsin 
                  Madison, WI 53706            (608) 262-0263 

                 * * * computation control switches * * *
             For the most part, the default is the only sensible
          value, and unless you are sure of what you are doing,
          these probably should not be touched.
          NPRINT =           Print/punch control flag
                             See also EXETYP for debug info.
                             (options -7 to 5 are primarily debug)
                 = -7        Extra printing from Boys localization.
                 = -6        debug for geometry searches
                 = -5        minimal output
                 = -4        print 2e-contribution to gradient.
                 = -3        print 1e-contribution to gradient.
                 = -2        normal printing, no punch file
                 =  1        extra printing for basis,symmetry,ZMAT
                 =  2        extra printing for MO guess routines
                 =  3        print out property and 1e- integrals
                 =  4        print out 2e- integrals
                 =  5        print out SCF data for each cycle.
                             (Fock and density matrices, current MOs
                 =  6        same as 7, but narrow terminal output
                             This option isn't perfect.
                 =  7        normal printing and punching (default)
                 =  8        more printout than 7. The extra output
                             is (AO) Mulliken and overlap population
                             analysis, eigenvalues, Lagrangians, ...
                 =  9        everything in 8 plus Lowdin population
                             analysis, final density matrix.

                      * * * restart options * * *
          IREST  =       restart control options
                         (for OPTIMIZE run restarts, see $STATPT)
                         Note that this option is unreliable!
                 = -1    reuse dictionary file from previous run,
                         useful with GEOM=DAF and/or GUESS=MOSAVED.
                         Otherwise, this option is the same as 0.
                 = 0     normal run (default)
                 = 1     2e restart (1-e integrals and MOs saved)
                 = 2     SCF restart (1-,2-e integrls and MOs saved)
                 = 3     1e gradient restart
                 = 4     2e gradient restart
          GEOM   =       select where to obtain molecular geometry
                 = INPUT from $DATA input (default for IREST=0)
                 = DAF   read from DICTNRY file (default otherwise)
              As noted in the first chapter, binary file restart is
          not a well tested option!
          NOSYM  = 0     the symmetry specified in $DATA is used
                         as much as possible in integrals, SCF,
                         gradients, etc.  (this is the default)
                 = 1     the symmetry specified in the $DATA group
                         is used to build the molecule, then
                         symmetry is not used again.   Some GVB
                         or MCSCF runs (those without a totally
                         symmetric charge density) require you
                         request no symmetry.
          INTTYP = POPLE use fast Pople routines for sp integral
                         blocks, and HONDO Rys polynomial code for
                         all other integrals.  (default)
                 = HONDO use HONDO/Rys integrals for all integrals.
                         This option produces slightly more accurate
                         integrals but is also slower.
          NORMF  = 0     normalize the basis functions (default)
                 = 1     no normalization
          NORMP  = 0     input contraction coefficients refer to
                         normalized Gaussian primitives. (default)
                 = 1     the opposite.
          ITOL   =       primitive cutoff factor (default=20)
                 = n     products of primitives whose preexponential
                         factor is less than 10**(-n) are skipped.
          ICUT   = n     integrals less than 10.0**(-n) are not
                         saved on disk. (default = 9)

          $SYSTEM group         (optional)
              This group provides global control information for
          your computer's operation.  This is system related input,
          and will not seem particularly chemical to you!
          TIMLIM =  time limit, in minutes.  Set to about 95 percent
                    of the time limit given to the batch job so that
                    GAMESS can stop itself gently.  (default=600.0)
          MEMORY =  establishes the maximum memory which your job
                    can use.  Some systems allocate just this amount
                    dynamically, others impose a static upper limit.
                    The default causes allocation of a system
                    dependent, moderate amount.  For many systems
                    this amount is 750,000 words.  (default=0)
          KDIAG  =    diagonalization control switch
                 = 0  use a vectorized diagonalization routine
                      if one is available on your machine,
                      else use EVVRSP. (default)
                 = 1  use EVVRSP diagonalization.  This may
                      be more accurate than KDIAG=0.
                 = 2  use GIVEIS diagonalization
                      (not as fast or reliable as EVVRSP)
                 = 3  use JACOBI diagonalization
                      (this is the slowest method)
          COREFL =  a flag to indicate whether or not GAMESS
                    should produce a "core" file for debugging
                    when subroutine ABRT is called to kill
                    a job.  This variable pertains only to
                    UNIX operating systems.  (default=.FALSE.)

            * * * the next four refer to parallel GAMESS * * *

          BALTYP =  parallel load balence scheme for integral
                    sections.  Choose LOOP to pick the inner
                    most loop for parallelization, and NXTVAL
                    to parallelize near the outer loop.  The
                    best strategy for equal speed processors
                    is LOOP, whereas NXTVAL will give better
                    load balance for mixed processors.  The
                    default is NXTVAL, except on a iPSC/860.
          XDR    =  a flag to indicate whether or not messages
                    should be converted into a generic format
                    known as external data representation.
                    If true, messages can exchange between 
                    machines of different vendors, at the cost
                    of performing the data type conversions.
          MEMPAR =  the same meaning as MEMORY, but for the other
                    nodes in a parallel computer.  Most useful for
                    analytic hessian runs, where the master needs
                    much more memory than the other nodes.  The 
                    default is to use the same value as MEMORY.
          PTIME  =  a logical flag to print extra timing info
                    during parallel runs.  This is not currently
          $BASIS group          (optional)
              This group allows certain standard basis sets to be
          easily given.  If this group is omitted, the basis set
          must be given instead in the $DATA group.
          GBASIS =        Name of the Gaussian basis set.
                 = MINI - Huzinaga's 3 gaussian minimal basis set.
                          Available H-Rn.
                 = MIDI - Huzinaga's 21 split valence basis set.
                          Available H-Xe.
                 = STO  - Pople's STO-NG minimal basis set.
                          Available H-Xe, for NGAUSS=2,3,4,5,6.
                 = N21  - Pople's N-21G split valence basis set.
                          Available H-Xe, for NGAUSS=3.
                          Available H-Ar, for NGAUSS=6.
                 = N31  - Pople's N-31G split valence basis set.
                          Available H-Ne,P-Cl for NGAUSS=4.
                          Available H-He,C-F for NGAUSS=5.
                          Available H-Ar, for NGAUSS=6.
                          For Ga-Kr, N31 selects the BC basis.
                 = N311 - Pople's "triple split" N-311G basis set.
                          Available H-Ne, for NGAUSS=6.
                          Selecting N311 implies MC for Na-Ar.
                 = DZV  - "double zeta valence" basis set.
                          a synonym for DH for H,Li,Be-Ne,Al-Cl.
                          a synonym for BC for Ga-Kr.
                 = DH   - Dunning/Hay "double zeta" basis set.
                          (3s)/[2s] for H.
                          (9s,4p)/[3s,2p] for Li.
                          (9s,5p)/[3s,2p] for Be-Ne.
                          (11s,7p)/[6s,4p] for Al-Cl.
                 = BC   - Binning/Curtiss "double zeta" basis set.
                          (14s,11p,5d/[6s,4p,1d] for Ga-Kr.
                 = TZV  - "triple zeta valence" basis set.
                          (5s)/[3s] for H.
                          (10s,3p)/[4s,3p] for Li.
                          (10s,6p)/[5s,3p] for Be-Ne.
                          a synonym for MC for Na-Ar.
                          (14s,9p)/[8s,4p] for K-Ca.
                          (14s,11p,6d)/[10s,8p,3d] for Sc-Zn.
                 = MC   - McLean/Chandler "triple split" basis.
                          (12s,9p)/[6s,5p] for Na-Ar.
                          Selecting MC implies 6-311G for H-Ne.
              additional values for GBASIS are on the next page.
               * * * the next two are ECP bases only * * *
          GBASIS = SBK  - Stevens/Basch/Krauss/Jasien/Cundari
                          valence basis set, for Li-Rn.  This choice
                          implies an unscaled -31G basis for H-He.
                 = HW   - Hay/Wadt valence basis.
                          This is a -21 split, available Na-Xe,
                          except for the transition metals.
                          This implies a 3-21G basis for H-Ne.
               * * * semiempirical basis sets * * *
                   The elements for which these exist can be found
                   in the 'further information' section of this
                   manual.  If you pick one of these, all other data
                   in this group is ignored.
          GBASIS = MNDO - selects MNDO model hamiltonian
                 = AM1  - selects AM1 model hamiltonian
                 = PM3  - selects PM3 model hamiltonian
          NGAUSS = the number of Gaussians (N).   This parameter
                   pertains only to GBASIS=STO, N21, N31, or N311.
          NDFUNC = number of heavy atom polarization functions to
                   be used.  These are usually d functions, except
                   for MINI/MIDI.  The term "heavy" means Na on up
                   when GBASIS=STO, HW, or N21, and from Li on up
                   otherwise.  The value may not exceed 3.  The
                   variable POLAR selects the actual exponents to
                   be used, see also SPLIT2 and SPLIT3. (default=0)
          NPFUNC = number of light atom, p type polarization
                   functions to be used on H-He.  This may not
                   exceed 3, see also POLAR.  (default=0)
          DIFFSP = flag to add diffuse sp (L) shell to heavy atoms.
                   Heavy means Li-F, Na-Cl, Ga-Br, In-I, Tl-At.
                   The default is .FALSE.
          DIFFS  = flag to add diffuse s shell to hydrogens.
                   The default is .FALSE.
          POLAR  = exponent of polarization functions
                 = POPLE     (default for all other cases)
                 = POPN311   (default for GBASIS=N311, MC)
                 = DUNNING   (default for GBASIS=DH, DZV)
                 = HUZINAGA  (default for GBASIS=MINI, MIDI)
                 = HONDO7    (default for GBASIS=TZV)
          SPLIT2 = an array of splitting factors used when NDFUNC
                   or NPFUNC is 2.  Default=2.0,0.5
          SPLIT3 = an array of splitting factors used when NDFUNC
                   or NPFUNC is 3.  Default=4.00,1.00,0.25
          The splitting factors are from the Pople school, and are
          probably too far apart.  See for example the Binning and
          Curtiss paper.  For example, the SPLIT2 value will usually
          cause an INCREASE over the 1d energy at the HF level for
          The actual exponents used for polarization functions, as
          well as for diffuse sp or s shells, are described in the
          'Further References' section of this manual.  This section
          also describes the sp part of the basis set chosen by
          GBASIS fully, with all references cited.
          Note that GAMESS always punches a full $DATA group.  Thus,
          if $BASIS does not quite cover the basis you want, you can
          obtain this full $DATA group from EXETYP=CHECK, and then
          change polarization exponents, add Rydbergs, etc.
          $DATA group          (required)
              This group describes the global molecular data such as
          point group symmetry, nuclear coordinates, and possibly
          the basis set.  It consists of a series of free format
          card images.
          -1-   TITLE     a single descriptive title card.
          -2-   GROUP, NAXIS
          GROUP is the Schoenflies symbol of the symmetry group,
          you may choose from
              C1, CS, CI, CN, S2N, CNH, CNV, DN, DNH, DND,
              T, TH, TD, O, OH.
          NAXIS is the order of the highest rotation axis, and
          must be given when the name of the group contains an N.
          For example, "CNV 2" is C2v.
          For linear molecules, choose either CNV or DNH, and enter
          NAXIS as 4.
              In order to use GAMESS effectively, you must be able
          to recognize the point group name for your molecule.  This
          presupposes a knowledge of group theory at about the level
          of Cotton's "Group Theory", Chapter 3.
              Armed with only the name of the group, GAMESS is able
          to exploit the molecular symmetry throughout almost all of
          the program, and thus save a great deal of computer time.
          GAMESS does not require that you know very much else about
          group theory, although a deeper knowledge (character
          tables, irreducible representations, term symbols, and so
          on) is useful when dealing with the more sophisticated
          Cards -3- and -4- are quite complicated, and are rarely
          given.  A *SINGLE* blank card may replace both cards -3-
          and -4-, to select the 'master frame', which is defined on
          the next page.   If you choose to enter a blank card, skip
          to the bottom of the next page.
          If the point group is C1 (no symmetry), skip over cards
          -3- and -4- (which means no blank card).
          -3-  X1, Y1, Z1, X2, Y2, Z2
          For C1 group, there is no card -3- or -4-.
          For CI group, give one point, the center of inversion.
          For CS group, any two points in the symmetry plane.
          For axial groups, any two points on the principal axis.
          For tetrahedral groups, any two points on a two-fold axis.
          For octahedral groups, any two points on a four-fold axis.
          -4-  X3, Y3, Z3, DIRECT
          third point, and a directional parameter.
          For CS group, one point of the symmetry plane,
                        noncollinear with points 1 and 2.
          For CI group, there is no card -4-.
          For other groups, a generator sigma-v plane (if any) is
          the (x,z) plane of the local frame (CNV point groups).
          A generator sigma-h plane (if any) is the (x,y) plane of
          the local frame (CNH and dihedral groups).
          A generator C2 axis (if any) is the x-axis of the local
          frame (dihedral groups).
          The perpendicular to the principal axis passing through
          the third point defines a direction called D1.  If
          DIRECT='PARALLEL', the x-axis of the local frame coincides
          with the direction D1.  If DIRECT='NORMAL', the x-axis of
          the local frame is the common perpendicular to D1 and the
          principal axis, passing through the intersection point of
          these two lines.  Thus D1 coincides in this case with the
          negative y axis.
              The 'master frame' is just a standard orientation for
          the molecule.  By default, the 'master frame' assumes that
              1.   z is the principal rotation axis (if any),
              2.   x is a perpendicular two-fold axis (if any),
              3.  xz is the sigma-v plane (if any), and
              4.  xy is the sigma-h plane (if any).
          Use the lowest number rule that applies to your molecule.
                  Some examples of these rules:
          Ammonia (C3v): the unique H lies in the XZ plane (R1,R3).
          Ethane (D3d): the unique H lies in the YZ plane (R1,R2).
          Methane (Td): the H lies in the XYZ direction (R2).  Since
                   there is more than one 3-fold, R1 does not apply.
          HP=O (Cs): the mirror plane is the XY plane (R4).
          In general, it is a poor idea to try to reorient the
          molecule.  Certain sections of the program, such as the
          orbital symmetry assignment, do not know how to deal with
          cases where the 'master frame' has been changed.
          Linear molecules (C4v or D4h) must lie along the z axis,
          so do not try to reorient linear molecules.
          You can use EXETYP=CHECK to quickly find what atoms are
          generated, and in what order.  This is typically necessary
          in order to use the general $ZMAT coordinates.
                               * * * *
          Depending on your choice for COORD in $CONTROL,
              if COORD=UNIQUE, follow card sequence U
              if COORD=HINT,   follow card sequence U
              if COORD=CART,   follow card sequence C
              if COORD=ZMT,    follow card sequence G
              if COORD=ZMTMPC, follow card sequence M
          Card sequence U is the only one which allows you to define
          a completely general basis here in $DATA.
          Recall that UNIT in $CONTRL determines the distance units.

          -5U-   Atom input.  Only the symmetry unique atoms are
          input, GAMESS will generate the symmetry equivalent atoms
          according to the point group selected above.
             if COORD=UNIQUE   NAME, ZNUC, X, Y, Z
          NAME  = 10 character atomic name, used only for printout.
                  Thus you can enter H or Hydrogen, or whatever.
          ZNUC  = nuclear charge.  It is the nuclear charge which
                  actually defines the atom's identity.
          X,Y,Z = Cartesian coordinates.
             if COORD=HINT
          NAME = 10 character atomic name (used only for print out).
          ZNUC = nuclear charge.
          CONX = connection type, choose from
            'LC'   linear conn.               'CCPA' central conn.
            'PCC'  planar central conn.              with polar atom
            'NPCC' non-planar central conn.   'TCT'  terminal conn.
            'PTC'  planar terminal conn.             with torsion
          R    = connection distance.
          ALPHA= first connection angle
          BETA = second connection angle
          SIGN = connection sign, '+' or '-'
          POINT1, POINT2, POINT3 =
               connection points, a serial number of a previously
               input atom, or one of 4 standard points: O,I,J,K
               (origin and unit points on axes of master frame).
               defaults:  POINT1='O', POINT2='I', POINT3='J'
          ref- R.L. Hilderbrandt, J.Chem.Phys. 51, 1654 (1969).
          You cannot understand HINT input without reading this.
          Note that if ZNUC is negative, the internally stored
          basis for ABS(ZNUC) is placed on this center, but the
          calculation uses ZNUC=0 after this.  This is useful
          for basis set superposition error (BSSE) calculations.
          * * * If you gave $BASIS, continue entering cards -5U-
                until all the unique atoms have been specified.
                When you are done, enter a " $END " card.
          * * * If you did not, enter cards -6U-, -7U-, -8U-.

          -6U-  GBASIS, NGAUSS, (SCALF(i),i=1,4)
          GBASIS has exactly the same meaning as in $BASIS.  You may
          choose from MINI, MIDI, STO, N21, N31, N311, DZV, DH, BC,
          TZV, MC, SBK, or HW.  In addition, you may choose S, P, D,
          F, G, or L to enter an explicit basis set.  Here, L means
          an s and p shell with a common exponent.
          NGAUSS is the number of Gaussians (N) in the Pople style
          basis, or user input general basis.  It has meaning only
          for GBASIS=STO, N21, N31, or N311, and S,P,D,F,G, or L.
          Up to four scale factors may be entered.  If omitted,
          standard values are used.  They are not documented as
          every GBASIS treats these differently.  Read the source
          code if you need to know more.  They are seldom given.
          * * * If GBASIS is not S,P,D,F,G, or L, either add more
                shells by repeating card -6U-, or go on to -8U-.
          * * * If GBASIS=S,P,D,F,G, or L, enter NGAUSS cards -7U-.
          -7U- IG, ZETA, C1, C2
                IG = a counter, IG takes values 1, 2, ..., NGAUSS.
              ZETA = Gaussian exponent of the IG'th primitive.
                C1 = Contraction coefficient for S,P,D,F,G shells,
                     and for the s function of L shells.
                C2 = Contraction coefficient for the p in L shells.
          * * * For more shells on this atom, go back to card -6U-.
          * * * If there are no more shells, go on to card -8U-.
          -8U-    A blank card ends the basis set for this atom.
          Continue entering atoms with -5U- through -8U- until all
          are given, then terminate the group with a " $END " card.
                 --- this is the end of card sequence U ---

          COORD=CART input:
          -5C- Atom input.
          Cartesian coordinates for all atoms must be entered.  They
          may be arbitrarily rotated or translated, but must possess
          the actual point group symmetry.  GAMESS will reorient the
          molecule into the 'master frame', and determine which
          atoms are the unique ones.  Thus, the final order of the
          atoms may be different from what you enter here.
                NAME, ZNUC, X, Y, Z
          NAME  = 10 character atomic name, used only for printout.
                  Thus you can enter H or Hydrogen, or whatever.
          ZNUC  = nuclear charge.  It is the nuclear charge which
                  actually defines the atom's identity.
          X,Y,Z = Cartesian coordinates.
          Continue entering atoms with card -5C- until all are
          given, and then terminate the group with a " $END " card.
                 --- this is the end of card sequence C ---

          COORD=ZMT input:       (GAUSSIAN style internals)
          -5G-      ATOM
          Only the name of the first atom is required.
          See -8G- for a description of this information.
          -6G-      ATOM  i1 BLENGTH
          Only a name and a bond distance is required for atom 2.
          See -8G- for a description of this information.
          -7G-      ATOM  i1 BLENGTH  i2 ALPHA
          Only a name, distance, and angle are required for atom 3.
          See -8G- for a description of this information.
          -8G-      ATOM  i1 BLENGTH  i2 ALPHA  i3 BETA i4
          ATOM    is the chemical symbol of this atom.  It can be
                  followed by numbers, if desired, for example Si3.
                  The chemical symbol implies the nuclear charge.
          i1      defines the connectivity of the following bond.
          BLENGTH is the bond length "this atom-atom i1".
          i2      defines the connectivity of the following angle.
          ALPHA   is the angle "this atom-atom i1-atom i2".
          i3      defines the connectivity of the following angle.
          BETA    is either the dihedral angle "this atom-atom i1-
                  atom i2-atom i3", or perhaps a second bond
                  angle "this atom-atom i1-atom i3".
          i4      defines the nature of BETA,
                  If BETA is a dihedral angle, i4=0 (default).
                  If BETA is a second bond angle, i4=+/-1.
                  (sign specifies one of two possible directions).
           o  Repeat -8G- for atoms 4, 5, ...
           o  The use of ghost atoms is possible, by using X or BQ
              for the chemical symbol.  Ghost atoms preclude the
              option of an automatic generation of $ZMAT.
           o  The connectivity i1, i2, i3 may be given as integers,
              1, 2, 3, 4, 5,...  or as strings which match one of
              the ATOMs.  In this case, numbers must be added to the
              ATOM strings to ensure uniqueness!

           o  In -6G- to -8G-, symbolic strings may be given in
              place of numeric values for BLENGTH, ALPHA, and BETA.
              The same string may be repeated, which is handy in
              enforcing symmetry.  If the string is preceeded by a
              minus sign, the numeric value which will be used is
              the opposite, of course.  Any mixture of numeric data
              and symbols may be given.  If any strings were given
              in -6G- to -8G-, you must provide cards -9G- and
              -10G-, otherwise you may terminate the group now with
              a " $END " card.
          -9G-   A blank line terminates the Z-matrix section.
          -10G-   STRING VALUE
          STRING is a symbolic string used in the Z-matrix.
          VALUE  is the numeric value to substitute for that string.
          Continue entering -10G- until all STRINGs are defined.
          Note that any blank card encountered while reading -10G-
          will be ignored.  GAMESS regards all STRINGs as variables
          (constraints are sometimes applied in $STATPT).  It is not
          necessary to place constraints to preserve point group
          symmetry, as GAMESS will never lower the symmetry from
          that given at -2-.  When you have given all STRINGs a
          VALUE, terminate the group with a " $END " card.
                 --- this is the end of card sequence G ---
                                * * * *
              The documentation for sequence G above and sequence M
          below presumes you are reasonably familiar with the input
          to GAUSSIAN or MOPAC.  It is probably too terse to be
          understood very well if you are unfamiliar with these.  A
          good tutorial on both styles of Z-matrix input can be
          found in Tim Clark's book "A Handbook of Computational
          Chemistry", published by John Wiley & Sons, 1985.
              Both Z-matrix input styles must generate a molecule
          which possesses the symmetry you requested at -2-.  If
          not, your job will be terminated automatically.

          COORD=ZMTMPC input:       (MOPAC style internals)
          -5M-     ATOM
          Only the name of the first atom is required.
          See -8M- for a description of this information.
          -6M-     ATOM BLENGTH
          Only a name and a bond distance is required for atom 2.
          See -8M- for a description of this information.
          -7M-     ATOM BLENGTH j1 ALPHA j2 
          Only a bond distance from atom 2, and an angle with repect
          to atom 1 is required for atom 3.  If you prefer to hook
          atom 3 to atom 1, you must give connectivity as in -8M-.
          See -8M- for a description of this information.
          -8M-     ATOM BLENGTH j1 ALPHA j2 BETA j3 i1 i2 i3
          ATOM, BLENGTH, ALPHA, BETA, i1, i2 and i3 are as described
          at -8G-.  However, BLENGTH, ALPHA, and BETA must be given
          as numerical values only.  In addition, BETA is always a
          dihedral angle.   i1, i2, i3 must be integers only.
          The j1, j2 and j3 integers, used in MOPAC to signal
          optimization of parameters, must be supplied but are
          ignored here.  You may give them as 0, for example.
          Continue entering atoms 3, 4, 5, ... with -8M- cards until
          all are given, and then terminate the group by giving a
          " $END " card.
                 --- this is the end of card sequence M ---
                         This is the end of $DATA!
          If you have any doubt about what molecule and basis set
          you are defining, or what order the atoms will be
          generated in, simply execute an EXETYP=CHECK job to find

          $ZMAT group      (required if NZVAR is nonzero in $CONTRL)
              This group lets you define the internal coordinates in
          which the gradient geometry search is carried out.  These
          need not be the same as the internal coordinates used in
          $DATA.  See $STATPT to freeze internals.  You must input a
          total of M=3N-6 internal coordinates (M=3N-5 for linear 
          molecules).  NZVAR in $CONTRL can be less than M IF AND 
          ONLY IF you are using linear bends.  
              It is also possible to input more than M coordinates
          if they are used to form exactly M linear combinations 
          for new internals.  These may be symmetry coordinates or 
          natural internal coordinates.  If NZVAR > M, you must 
          input IJS and SIJ below to form M new coordinates.  See 
          DECOMP in $FORCE for the only circumstance in which you 
          may enter a larger NZVAR without giving SIJ and IJS.

          IZMAT is an array of integers defining each coordinate.
          The general form for each internal coordinate is
                code number,I,J,K,L,M,N
          IZMAT =1 followed by two atom numbers. (I-J bond length)
                =2 followed by three numbers. (I-J-K bond angle)
                =3 followed by four numbers. (dihedral angle)
                   Torsion angle between planes I-J-K and J-K-L.
                =4 followed by four atom numbers. (atom-plane)
                   Out-of-plane angle from bond I-J to plane J-K-L.
                =5 followed by three numbers. (I-J-K linear bend)
                   Counts as 2 coordinates for the degenerate bend,
                   normally J is the center atom.  See $LIBE.
                =6 followed by five atom numbers. (dihedral angle)
                   Dihedral angle between planes I-J-K and K-L-M.
                =7 followed by six atom numbers. (ghost torsion)
                   Let A be the midpoint between atoms I and J, and
                   B be the midpoint between atoms M and N.  This
                   coordinate is the dihedral angle A-K-L-B.  The
                   atoms I,J and/or M,N may be the same atom number.
                   (If I=J AND M=N, this is a conventional torsion).
                   Examples: N2H4, or, with one common pair, H2POH.
          Example - a nonlinear triatomic, atom 2 in the middle:
                $ZMAT IZMAT(1)=1,1,2,  2,1,2,3,  1,2,3  $END
          This sets up two bonds and the angle between them.
          The blanks between each coordinate definition are
          not necessary, but improve readability mightily.

                                                         $ZMAT $LIBE

          SIJ is a transformation matrix of dimension NZVAR x M,
              used to transform the NZVAR internal coordinates in
              IZMAT into M new internal coordinates.  SIJ is a
              sparse matrix, so only the non-zero elements are 
              given, by using the IJS array described below.  
              The columns of SIJ will be normalized by GAMESS.
              (Default: SIJ = I, unit matrix)

          IJS is an array of pairs of indices, giving the row and
              column index of the entries in SIJ.

          example - if the above triatomic is water, using
               IJS = 1,1, 3,1,   1,2, 3,2,   2,3
               SIJ = 1.0, 1.0,   1.0,-1.0,   1.0

              gives the matrix S=  1.0   1.0   0.0 
                                   0.0   0.0   1.0
                                   1.0  -1.0   0.0

          which defines the symmetric stretch, asymmetric stretch,
          and bend of water.

          references for natural internal coordinates:
            P.Pulay, G.Fogarasi, F.Pang, J.E.Boggs  
               J.Am.Chem.Soc. 101, 2550-2560(1979)
            G.Fogarasi, X.Zhou, P.W.Taylor, P.Pulay
               J.Am.Chem.Soc. 114, 8191-8201(1992)


          $LIBE group  (required if linear bends are used in $ZMAT)
          A degenerate linear bend occurs in two orthogonal planes,
          which are specified with the help of a point A.  The first
          bend occurs in a plane containing the atoms I,J,K and the
          user input point A.  The second bend is in the plane
          perpendicular to this, and containing I,J,K.  One such
          point must be given for each pair of bends used.
          APTS(1)= x1,y1,z1,x2,y2,z2,...  for linear bends 1,2,...

          Note that each linear bend serves as two coordinates, so
          that if you enter 2 linear bends (HCCH, for example), the
          correct value of NZVAR is M-2, where M=3N-6 or 3N-5, as


          $SCF group         relevant if SCFTYP = RHF, UHF, or ROHF,
                             required if SCFTYP = GVB)
              This group of parameters provides additional control
          over the RHF, UHF, ROHF, or GVB SCF steps.  It must be
          used for GVB open shell or perfect pairing wavefunctions.
          DIRSCF = a flag to activate a direct SCF calculation,
                   which is implemented for all the Hartree-Fock
                   type wavefunctions:  RHF, ROHF, UHF, and GVB.
                   This keyword also selects direct MP2 computation.
                   The default of .FALSE. stores integrals on disk
                   storage for a conventional SCF calculation.
          FDIFF  = a flag to compute only the change in the Fock
                   matrices since the previous iteration, rather
                   than recomputing all two electron contributions.
                   This pertains only to direct SCF, and has a
                   default of .TRUE.  This option is implemented
                   only for the RHF, ROHF, UHF cases.
          ---- The next flags affect convergence rates.
          EXTRAP = controls Pople extrapolation of the Fock matrix.
          DAMP   = controls Davidson damping of the Fock matrix.
          SHIFT  = controls level shifting of the Fock matrix.
          RSTRCT = controls restriction of orbital interchanges.
          DIIS   = controls Pulay's DIIS interpolation.
          DEM    = controls direct energy minimization, which is
                   implemented only for RHF.  (default=.FALSE.)
          defaults for     EXTRAP  DAMP  SHIFT RSTRCT  DIIS
          ab initio:         T      F      F      F      T
          semiempirical:     T      F      F      F      F
               The above parameters are implemented for all SCF
          wavefunction types, except that DIIS will work for GVB
          only for the case NPAIR=1 and NSETO=0.
               Once DIIS is initiated, any other accelerator in
          effect is put in abeyance, but will be turned back on
          again should the DIIS procedure stop.
                    * * * * * * * * * * * * * * * * * * * *
                    For a description of the convergence
                    accelerators, and convergence criteria,
                    see the 'Further Information' section.
                    * * * * * * * * * * * * * * * * * * * *

          ---- These parameters fine tune the various convergers.
          NCONV  = n  SCF density convergence criteria.
                   Convergence is reached when the density change
                   between two consecutive SCF cycles is less than
                   10.0**(-n) in absolute value.  One more cycle
                   is executed after reaching convergence.   Less
                   accuracy in NCONV gives questionable gradients.
                   (default is n=5)
          ETHRSH = energy error threshold for initiating DIIS.  The
                   DIIS error is the largest element of e=FDS-SDF.
                   Increasing ETHRSH forces DIIS on sooner.
                   (default = 0.5 Hartree)
          MAXDII = Maximum size of the DIIS linear equations, so
                   that at most MAXDII-1 Fock matrices are used
                   in the interpolation.  (default=10)
          DEMCUT = Direct energy minimization will not be done
                   once the density matrix change falls below
                   this threshold.  (Default=0.5)
          DMPCUT = Damping factor lower bound cutoff.  The damping
                   damping factor will not be allowed to drop
                   below this value.  This only has effect if
                   DAMP=.T.  (default=0.0)
             note: The damping factor need not be zero to achieve
                   valid convergence (see Hsu, Davidson, and
                   Pitzer, J.Chem.Phys., 65, 609 (1976),
                   especially the section on convergence control),
                   but it should not be astronomical either.
          ---- Miscellaneous options.
          UHFNOS = flag controlling generation of the natural
                   orbitals of a UHF function. (default=.FALSE.)
          MVOQ   = 0  Skip MVO generation (default)
                 = n  Form modified virtual orbitals, using a cation
                      with n electrons removed.   Implemented for
                      RHF, ROHF, and GVB.   If necessary to reach a
                      closed shell cation, the program might remove
                      n+1 electrons.  Typically, n will be about 6.
          NPUNCH = SCF punch option
                 =  0  do not punch out the final orbitals
                 =  1  punch out the occupied orbitals
                 =  2  punch out occupied and virtual orbitals
                       The default is NPUNCH = 2.

              The next parameters define the GVB wavefunction.  Note
          that ALPHA and BETA also have meaning for ROHF.  See also
          MULT in the $CONTRL group.  The GVB wavefunction assumes
          orbitals are in the order core, open, pairs.
          NCO    =   The number of closed shell orbitals.  The
                     default almost certainly should be changed!
          NSETO  =   The number of sets of open shells in the
                     function.  Maximum of 10. (default=0)
          NO     =   An array giving the degeneracy of each open
                     shell set.  Give NSETO values.
          NPAIR  =   The number of geminal pairs in the -GVB-
                     function.  Maximum of 12.  The default
                     corresponds to open shell SCF (default=0).
          CICOEF =   An array of ordered pairs of CI coefficients
                     for the -GVB- pairs.  For example, a two pair
                     case for water, say, might be
                     CICOEF(1)=0.95,-0.05,0.95,-0.05.  If not
                     normalized, as in the default, they will be.
                     This parameter is useful in restarting a GVB
                     run, with the current CI coefficients.
                     (default = 0.90,-0.20,0.90,-0.20,...)
          COUPLE =   A switch controlling the input of F, ALPHA,
                     and BETA.  The default is to use internally
                     stored values for these variables.   Note
                     ALPHA and BETA can be given for -ROHF-, as
                     well as -GVB-.  (Default=.FALSE.)
          F      =   An vector of fractional occupations.
          ALPHA  =   An array of A coupling coefficients given in
                     lower triangular order.
          BETA   =   An array of B coupling coefficients given in
                     lower triangular order.
              Note:  The default for F, ALPHA, and BETA depends on
          the state chosen.  Defaults for the most commonly occuring
          cases are internally stored.

                      * * * * * * * * * * * * * * * * * * *
                      For more discussion of GVB/ROHF input
                      see the 'further information' section
                      * * * * * * * * * * * * * * * * * * *


          VTSCAL =   A flag to request that the virial theorem be
                     satisfied.  An analysis of the total energy
                     as an exact sum of orbital kinetic energies 
                     is printed.  The default is .FALSE.

             This option is implemented for RHF, UHF, and ROHF,

          SCALF  =   initial exponent scale factor when VTSCAL is
                     in use, useful when restarting.  The default
                     is 1.0.

          MAXVT  =   maximum number of iterations (at a single 
                     geometry) to satisfy the energy virial theorem.
                     The default is 20.

          VTCONV =   convergence criterion for the VT, which is 
                     satisfied when 2<T> + <V> + R x dE/dR is less 
                     than VTCONV.  The default is 1.0D-6 Hartree.

          For more information on this option, which is most
          economically employed during a geometry search, see 
          M.Lehd and F.Jensen, J.Comput.Chem. 12, 1089-1096(1991).


          $MP2 group        (not required, relevant if MPLEVL=2)
               Controls 2nd order Moller-Plesset perturbation runs,
          if requested by MPLEVL in $CONTRL.  See also the DIRSCF
          keyword in $SCF to select direct MP2.  MP2 is implemented 
          for RHF, high spin ROHF, or UHF wavefunctions.  Note that 
          the gradient and the properties of the first order wave-
          function cannot be computed, so properties are for the 
          unperturbed wavefunction (the SCF).  The $MP2 group is 
          usually not given.  
          NCORE = n  Omits the first n occupied orbitals from the
                     calculation.  The default for n is the number
                     of chemical core orbitals.

          CUTOFF=    transformed integral retention threshold, the
                     default is 1.0d-9.

          METHOD= n  selects transformation method, 2 being the
                     segmented transformation, and 3 being a more
                     conventional two phase bin sort implementation.
                     3 requires more disk, but less memory.  The
                     default is to attempt method 2 first, and 
                     method 3 second.

          $GUESS group         (optional, relevant for all SCFTYP's)
              This group controls the selection of initial molecular
          GUESS = Selects type of initial orbital guess.
                = HUCKEL   Carry out an extended Huckel calculation
                           using a Huzinaga MINI basis set, and
                           project this onto the current basis.
                           This is implemented for atoms up to Rn,
                           and will work for any all electron or
                           ECP basis set.  (default for most runs)
                = HCORE    Diagonalize the one electron Hamiltonian
                           to obtain the initial guess orbitals.
                           This method is applicable to any basis
                           set, but does not work as well as the
                           HUCKEL guess.
                = MOREAD   Read in formatted vectors punched by an
                           earlier run.  This requires a $VEC group,
                           and you MUST pay attention to NORB below.
                = MOSAVED  (default for restarts)  The initial
                           orbitals are read from the DICTNRY file
                           of the earlier run.
                = SKIP     Bypass initial orbital selection.  The
                           initial orbitals and density matrix are
                           assumed to be in the DICTNRY file.
              All GUESS types except 'SKIP' permit reordering of the
          orbitals, carry out an orthonormalization of the orbitals,
          and generate the correct initial density matrix.  The
          initial density matrix cannot be generated for -CI- and
          -MCSCF-, so property restarts for these wavefunctions will
          require 'SKIP' which is an otherwise seldom used option.
          Note that correct computation of a -GVB- density matrix
          requires CICOEF in $SCF.  Another possible use for 'SKIP'
          is to speed up a EXETYP=CHECK job, or a RUNTYP=HESSIAN job
          where the hessian is supplied.
          PRTMO = a flag to control printing of the initial guess.

                                                         $GUESS $VEC

          NORB   = The number of orbitals to be read in the $VEC
                   group.  This applies only to GUESS=MOREAD.
          For -RHF-, -UHF-, -ROHF-, and -GVB-, NORB defaults to the
          number of occupied orbitals.  NORB must be given for -CI-
          and -MCSCF-.  For -UHF-, if NORB is not given, only the
          occupied alpha and beta orbitals should be given, back to
          back.  Otherwise, both alpha and beta orbitals must
          consist of NORB vectors.
          NORB may be larger than the number of occupied MOs, if you
          wish to read in the virtual orbitals.  If NORB is less
          than the number of atomic orbitals, the remaining orbitals
          are generated as the orthogonal complement to those read.
          NORDER = Orbital reordering switch.
                 = 0  No reordering (default)
                 = 1  Reorder according to IORDER and JORDER.
          IORDER = Reordering instructions.
                   Input to this array gives the new molecular
                   orbital order.  For example, IORDER(3)=4,3 will
                   interchange orbitals 3 and 4, while leaving the
                   other MOs in the original order.  This parameter
                   applies to all orbitals (alpha and beta) except
                   for -UHF-, where it only affects the alpha MOs.
                   (default is IORDER(i)=i )
          JORDER = Reordering instructions.
                   Same as IORDER, but for the beta MOs of -UHF-.
          TOLZ   = level below which MO coefficients will be set
                   to zero.  (default=1.0E-7)
          TOLE   = level at which MO coefficients will be equated.
                   This is a relative level, coefficients are set
                   equal if one agrees in magnitude to TOLE times
                   the other.  (default=5.0E-5)
          MIX    = rotate the alpha and beta HOMO and LUMO orbitals
                   so as to generate inequivalent alpha and beta
                   orbital spaces.  This pertains to UHF singlets
                   only.  (default=.FALSE.)
          $VEC group           (optional, relevant for all SCFTYP's)
                               (required if GUESS=MOREAD)
                This group consists of formatted vectors, as written
          onto file PUNCH in a previous run.  It is considered good
          form to retain the titling comment cards punched before
          the $VEC card, as labeling of what the $VEC is!

          $STATPT group  (optional, for RUNTYP=OPTIMIZE or SADPOINT)
              This group controls the search for stationary points.
          Note that NZVAR in $CONTRL determines if the geometry
          search is conducted in Cartesian or internal coordinates.
          METHOD = optimization selection switch.  You can choose
                   STANDARD or SCHLEGEL (the default is STANDARD).
          OPTTOL = gradient convergence tolerance, in Hartree/Bohr.
                   Convergence of a geometry search requires the
                   largest component of the gradient to be less
                   than OPTTOL, and the root mean square gradient
                   less than 1/3 of OPTTOL.  (default=0.0001)
          NSTEP  = maximum number of steps to take.  Restart data
                   is punched if NSTEP is exceeded.  (default=20)
          DXMAX  = initial trust radius of the step, in Bohr.
                   Steps will be scaled down to this value, if
                   necessary.  (default=0.2, except STANDARD method
                   geometry searches, when it is 0.3)

              the next three apply only to METHOD=STANDARD

          TRUPD  = a flag to allow the trust radius to change as
                   the geometry search proceeds.  (default=.TRUE.)

          TRMAX  = maximum permissible value of the trust radius.
          TRMIN  = minimum permissible value of the trust radius.
              The next two control the hessian matrix quality

          HESS   = selects the initial hessian matrix.
                 = GUESS chooses a positive definite diagonal
                         hessian. (default for RUNTYP=OPTIMIZE)
                 = READ  causes the hessian to be read from a $HESS
                         group. (default for RUNTYP=SADPOINT)
                 = CALC  causes the hessian to be computed, see
                         the $FORCE group.
          IHREP  = the number of steps before the hessian is
                   recomputed.  If given as 0, the hessian will
                   be computed only at the initial geometry if
                   you choose HESS=CALC, and never again.  If
                   nonzero, the hessian is recalculated every
                   IHREP steps, with the update formula used on
                   other steps.  (default=0)

          IFREEZ = array of internal coordinates to freeze.
                   For example, IFREEZ(1)=1,3 would freeze the two
                   bond lengths in the $ZMAT example, while
                   optimizing the angle.  You cannot freeze
                   coordinates using Cartesian coordinates.

          IFOLOW = Mode selection switch, for RUNTYP=SADPOINT.
                   For METHOD=STANDARD, the mode along which the
                   energy is maximized, other modes are minimized.
                   Usually refered to as "eigenvector following".
                   For METHOD=SCHLEGEL, the mode whose eigenvalue
                   is (or will be made) negative.  All other
                   curvatures will be made positive.  
                   (default is 1)

              Let 0 mean the initial geometry, L mean the last
              geometry, and all mean every geometry.
              Let INTR mean the internuclear distance matrix.
              Let HESS mean the approximation to the hessian.
              Note that a directly calculated hessian matrix
              will always be punched, NPUN refers only to the
              updated hessians used by the quasi-Newton step.

          NPRT   =  1  Print INTR at all, orbitals at all
                    0  Print INTR at all, orbitals at 0+L (default)
                   -1  Print INTR at all, orbitals never
                   -2  Print INTR at 0+L, orbitals never

          NPUN   =  3  Punch all orbitals and HESS at all
                    2  Punch all orbitals at all
                    1  Punch occ orbitals after 0, and HESS at all
                    0  Punch occ orbitals after 0, but all 
                       orbitals at geometry 0 (default)
                   -1  Punch occ orbitals at 0+L only
                   -2  Never punch orbitals
          MOVIE  = a flag to create a series of structural data 
                   which can be show as a movie by the MacIntosh
                   program Chem3D.  The data is written to the
                   file IRCDATA.  (default=.FALSE.)

          HSSEND = a flag to control automatic hessian evaluation
                   at the end of a successful geometry search.


           ---- the following parameters are quite specialized ----
          PURIFY = a flag to help eliminate the rotational and 
                   translational degrees of freedom from the 
                   initial hessian (and possibly initial gradient).
                   This is much like the variable of the same name
                   in $FORCE, and will be relevant only if internal
                   coordinates are in use.  (default=.FALSE.)
          ITBMAT = number of micro-iterations used to compute the
                   step in Cartesians which corresponds to the 
                   desired step in internals.  The default is 5.
          UPHESS = SKIP     do not update Hessian (not recommended)
                   BFGS     default for OPTIMIZE using STANDARD
                   POWELL   default for SADPOINT using STANDARD
                   SCHLEGEL only choice for METHOD=SCHLEGEL
          RESTAR = Enables restart of an optimization run.  This
                   can only be used with IREST .ne. 0 in $CONTRL.
                   Use of this variable is discouraged.
           ---- NNEG, RMIN, RMAX, RLIM apply only to SCHLEGEL ----
          NNEG   = The number of negative eigenvalues the force
                   constant matrix should have. If necessary the
                   smallest eigenvalues will be reversed. The
                   default is 0 for RUNTYP=OPTIMIZE, and 1 for

          RMIN   = Minimum distance threshold. Points whose root
                   mean square distance from the current point is
                   less than RMIN are discarded. (default=0.0015)
          RMAX   = Maximum distance threshold. Points whose root
                   mean square distance from the current point is
                   greater than RMAX are discarded. (default=0.1)

          RLIM   = Linear dependence threshold. Vectors from the
                   current point to the previous points must not
                   be colinear.  (default=0.07)
                 * * * * * * * * * * * * * * * * * * * * *
                 See the 'further information' section for
                 some help with OPTIMIZE and SADPOINT runs
                 * * * * * * * * * * * * * * * * * * * * *

          $TRUDGE group (optional, required for RUNTYP=TRUDGE)
              This group defines the parameters for a non-gradient
          optimization of exponents or the geometry.  The TRUDGE
          package is a modified version of the same code from Michel
          Dupuis' HONDO 7.0 system, origially written by H.F.King.
          Presently the program allows for the optimization of 10
              Exponent optimization works only for uncontracted
          primitives, without enforcing any constraints.  Two
          non-symmetry equivalent H atoms would have their p
          function exponents optimized separately, and so would two
          symmetry equivalent atoms!
              Geometry optimization works only in HINT internal
          coordinates (see $CONTRL and $DATA groups).  The total
          energy of all types of SCF wavefunctions can be optimized,
          although this would be extremely stupid as gradient
          methods are far more efficient.  The main utility is for
          MP2 or CI geometry optimizations, which may not be done in
          any other way with GAMESS.
          OPTMIZ = a flag to select optimization of either geometry
                   or exponents of primitive gaussian functions.
                 = GEOMETRY for geometry optimization (default).
                 = BASIS    for basis set optimization.
          CIDRT  = CI selection flag;
                 = NONE means do not perform a CI.   If OPTMIZ is
                        GEOMETRY, GAMESS will optimize the SCFTYP
                        energy, or the MP2 energy if MPLEVL is set.
                 = "STRING" means do perform a CI.  The $STRING
                        group will be read, instead of $DRT.  This
                        permits optimization of a CI function based
                        MCSCF orbitals, which would be defined by
                        the normal $DRT group.  STRING may be any
                        name (e.g. CIDRT), but should not contain
                        a $ sign.  If SCFTYP is RHF, ROHF, or GVB,
                        you may use the normal name (DRT).
              Note that you may optimize the geometry for an excited
              CI state, just specify
                    $GUGDIA   NSTATE=5  $END
                    $GUGDM    IROOT=3   $END
              to find equilibrium geometry of the third state
              in the symmetry implicit in your $STRING group.

          NPAR   = number of parameters to be optimized.
          IEX    = defines the parameters to be optimized.
                   If OPTMIZ=BASIS, IEX declares the serial number
              of the Gaussian primitives for which the exponents
              will be optimized.
                   If OPTMIZ=GEOMETRY, IEX define the pointers to
              the HINT internal coordinates which will be optimized.
              (Note that not all internal coordinates have to be
              optimized.) The pointers to the internal coordinates
              are defined as:  (the number of atom on the input
              list)*10 + (the number of internal coordinate for that
              atom).  For each atom, the HINT internal coordinates
              are numbered as 1, 2, and 3 for BOND, ALPHA, and BETA,
          P  =  Defines the initial values of the parameters to be
                optimized.  You can use this to reset values given
                in $DATA.  If omitted, the $DATA values are used.
                If given here, geometric data must be in Angstroms
                and degrees.
          A complete example is
          CNV 2
          O    8.     LC  0.00  0.0  0.00  -  O  K
          H    1.    PCC  0.94  53.  0.00  +  O  K  I
                   IEX(1)=21,22   P(1)=0.96 $END
           $DRT    GROUP=C2V IEXCIT=2 NFZC=1 NDOC=4 NVAL=14 $END
          which is a CISD/6-31G* geometry optimization of water,
          using RHF orbitals in the CI.  The starting bond length is
          reset to 0.96, while the initial angle will be 106 (twice

          $TRURST group       (optional, relevant for RUNTYP=TRUDGE)
                This  group  specifies restart parameters for TRUDGE
          runs and accuracy thresholds.
          KSTART indicates the conjugate gradient direction in which
          the optimization will proceed. ( default = -1 )
               -1 .... indicates that this is a non-restart run.
                0 .... corresponds to a restart run.
          FNOISE accuracy of function values.
          Variation smaller than FNOISE are not considered to be
          significant (Def. 0.0005)
          TOLF accuracy required of the function (Def. 0.001)
          TOLR accuracy required of conjugate directions (Def. 0.05)
              For geometry optimization, the values which give
          better results (closer to the ones obtained with gradient
          methods) are:  TOLF=0.0001, TOLR=0.001, FNOISE=0.00001

          $FORCE group
          (optional, relevant for RUNTYP=HESSIAN,OPTIMIZE,SADPOINT)
              This group controls the computation of the hessian
          matrix (the energy second derivative tensor, also known
          as the force constant matrix), and an optional harmonic
          vibrational analysis.  This can be a very time consuming
          calculation.  However, given the force constant matrix,
          the vibrational analysis for an isotopically substituted
          molecule is very cheap.  Related input is HESS= in
          $STATPT, and the $MASS, $HESS, $GRAD, $DIPDR, $VIB groups.
          METHOD = chooses the computational method.
                 = ANALYTIC is implemented only for SCFTYPs RHF,
                            ROHF, and GVB (when NPAIR is 0 or 1),
                            if the calculation does not use ECPs.
                            This is the default for these cases.
                 = NUMERIC  is the default for all other cases.
                IR intensities are available only
                for NUMERIC runs at present.
          RDHESS = a flag to read the hessian from a $HESS group,
                   rather than computing it.  This variable pertains
                   only to RUNTYP=HESSIAN.  See also HESS= in the
                   $STATPT group.  (default is .FALSE.)
          PURIFY = controls cleanup
                   Given a $ZMAT, the hessian and dipole derivative
                   tensor can be "purified" by transforming from
                   Cartesians to internals and back to Cartesians.
                   This effectively zeros the frequencies of the
                   translation and rotation "modes", along with
                   their IR intensities.  The purified quantities
                   are punched out.  Purification does change the
                   Hessian slightly, frequencies at a stationary
                   point can change by a wave number or so.  The
                   change is bigger at non-stationary points.
                   (default=.FALSE. if $ZMAT is given)
          PRTIFC = prints the internal coordinate force constants.
                   You MUST have defined a $ZMAT group to use this.

            --- the next four apply only to METHOD=NUMERIC ----
          NVIB   =    Number of displacements in each Cartesian
                      direction for force field computation.
                 = 1  Move one VIBSIZ unit in each positive
                      Cartesian direction.  This requires 3N+1
                      evaluations of the wavefunction, energy, and
                      gradient, where N is the number of SYMMETRY
                      UNIQUE atoms given in $DATA.  (default)
                 = 2  Move one VIBSIZ unit in the positive direction
                      and one VIBSIZ unit in the negative direction.
                      This requires 6N+1 evaluations of the
                      wavefunction and gradient, and gives a small
                      improvement in accuracy.  In particular, the
                      frequencies will change from NVIB=1 results by
                      no more than 10-100 wavenumbers, and usually 
                      much less.  However, the normal modes will be
                      more nearly symmetry adapted, and the residual
                      rotational and translational "frequencies"
                      will be much closer to zero.
          VIBSIZ =    Displacement size (in Bohrs). Default=0.01
                 Let 0 mean the Vib0 geometry, and 
                 D mean all the displaced geometries

          NPRT   = 1  Print orbitals at 0 and D
                 = 0  Print orbitals at 0 only (default)

          NPUN   = 2  Punch all orbitals at 0 and D
                 = 1  Punch all orbitals at 0 and occupied orbs at D
                 = 0  Punch all orbitals at 0 only (default)

            ----- the rest control normal coordinate analysis ----
          VIBANL = flag to activate vibrational analysis.
                   (the default is .TRUE. for RUNTYP=HESSIAN, and
                   otherwise is .FALSE.)
          SCLFAC = scale factor for vibrational frequencies, used
                   in calculating the zero point vibrational energy.
                   Some workers correct for the usual overestimate
                   in SCF frequencies by a factor 0.89.  The output
                   always prints unscaled frequencies, this value
                   is used only in the thermochemical analysis.
                   (Default is 1.0)

          TEMP   = an array of up to ten temperatures at which the
                   thermochemistry should be printed out.  The
                   default is a single temperature, 298.15 K.
          FREQ   = an array of vibrational frequencies.  If the
                   frequencies are given here, the hessian matrix
                   is not computed or read.  You enter any imaginary
                   frequencies as negative numbers, omit the
                   zero frequencies corresponding to translation
                   and rotation, and enter all true vibrational
                   frequencies.  Thermodynamic properties will be
                   printed, nothing else is done by the run.
          PRTSCN = flag to print contribution of each vibrational
                   mode to the entropy.  (Default is .FALSE.)
          DECOMP = activates internal coordinate analysis.
                   Vibrational frequencies will be decomposed into
                   "intrinsic frequencies", by the method of
                   J.A.Boatz and M.S.Gordon, J.Phys.Chem., 93,
                   1819-1826(1989).  If set .TRUE., the $ZMAT group
                   may define more than 3N-6 (3N-5) coordinates.
          PROJCT = controls the projection of the hessian matrix.
                   The projection technique is described by
                   W.H.Miller, N.C.Handy, J.E.Adams in J. Chem.
                   Phys. 1980, 72, 99-112.  At stationary points,
                   the projection simply eliminates rotational and
                   translational contaminants.  At points with
                   non-zero gradients, the projection also ensures
                   that one of the vibrational modes will point
                   along the gradient, so that there are a total of
                   7 zero frequencies.  The other 3N-7 modes are
                   constrained to be orthogonal to the gradient.
                   Because the projection has such a large effect on
                   the hessian, the hessian punched is the one
                   BEFORE projection.  For the same reason, the
                   default is .FALSE. to skip the projection, which
                   is mainly of interest in dynamical calculations.

                                                  $HESS $GRAD $DIPDR

          $HESS group (relevant for RUNTYP=HESSIAN if RDHESS=.TRUE.)
                   (relevant for RUNTYP=IRC if FREQ,CMODE not given)
                (relevant for RUNTYP=OPTIMIZE,SADPOINT if HESS=READ)
              Formatted force constant matrix (FCM), i.e. hessian
          matrix.  This data is punched out by a RUNTYP=HESSIAN job,
          in the correct format for subsequent runs.  The first card
          in the group must be a title card.
              A $HESS group is always punched in Cartesians.  It
          will be transformed into internal coordinate space if a
          geometry search uses internals.  It will be mass weighted
          (according to $MASS) for IRC and frequency runs.
              The initial FCM is updated during the course of a
          geometry optimization or saddle point search, and will be
          punched if a run exhausts its time limit.  This allows
          restarts where the job leaves off.  You may want to read
          this FCM back into the program for your restart, or you
          may prefer to regenerate a new initial hessian.  In any
          case, this updated hessian is absolutely not suitable for
          frequency prediction!
          $GRAD group     (relevant for RUNTYP=OPTIMIZE or SADPOINT)
                    (relevant for RUNTYP=HESSIAN when RDHESS=.TRUE.)
              Formatted gradient vector at the $DATA geometry.  This
          data is read in the same format it was punched out.
              For RUNTYP=HESSIAN, this information is used to
          determine if you are at a stationary point, and possibly
          for projection.  If omitted, the program pretends the
          gradient is zero, and otherwise proceeds normally.
              For geometry searches, this information (if known) can
          be read into the program so that the first step can be
          taken instantly.
          $DIPDR group   (relevant for RUNTYP=HESSIAN if RDHESS=.T.)
          Formatted dipole derivative tensor, punched in a previous
          RUNTYP=HESSIAN job.  If this group is omitted, then a
          vibrational analysis will be unable to predict the IR
          intensities, but the run can otherwise proceed.

                                                          $VIB $MASS
          $VIB group   (relevant for RUNTYP=HESSIAN, METHOD=NUMERIC)
              Formatted card image -restart- data.  This data is
          read in the format it was punched by a previous HESSIAN
          job to the file IRCDATA.  Just add a " $END" card, and if
          the final gradient was punched as zero, delete the last
          set of data.  Normally, IREST in $CONTRL will NOT be used
          in conjunction with a HESSIAN restart.  The mere presence
          of this deck triggers the restart from cards.  This deck
          can also be used to turn a single point differencing run
          into double differencing, as well as recovering from time
          limits, or other bombouts.
          $MASS group           (relevant for RUNTYP=HESSIAN or IRC)
              This group permits isotopic substitution during the
          computation of mass weighted Cartesian coordinates.  Of
          course, the masses affect the frequencies and normal modes
          of vibration.
          AMASS = An array giving the atomic masses, in amu. The
                  default is to use the mass of the most abundant
                  isotope.  Masses through Xenon are stored.
          example - $MASS AMASS(3)=2.0140 $END
          will make the third atom in the molecule a deuterium.

          $IRC group                       (relevant for RUNTYP=IRC)
              This group governs the location of the intrinsic
          reaction coordinate, a steepest descent path in mass
          weighted corrdinates, that connects the saddle point to
          reactants and products.
          ----- there are five integration methods chosen by PACE.
          PACE = GS2    selects the Gonzalez-Schlegel second order
                        method.  This is the default method.  
                        Related input is:
            GCUT   cutoff for the norm of the mass-weighted gradient
                   tangent (the default is chosen in the range from
                   0.00005 to 0.00020, depending on the value for
                   STRIDE chosen below.
            RCUT   cutoff for Cartesian RMS displacement vector.
                   (the default is chosen in the range 0.0005 to
                   0.0020 Bohr, depending on the value for STRIDE)
            ACUT   maximum angle from end points for linear
                   interpolation (default=5 degrees)
            MXOPT  maximum number of contrained optimization steps
                   for each IRC point (default=20)
            IHUPD  is the hessian update formula.  1 means Powell,
                   2 means BFGS (default=2)
            GA     is a gradient from the previous IRC point, and is
                   used when restarting.
            OPTTOL is a gradient cutoff used to determine if the IRC
                   is approaching a minimum.  It has the same meaning
                   as the variable in $STATPT.  (default=0.0001)

          PACE = LINEAR selects linear gradient following (Euler's
                        method).  Related input is:
            STABLZ switches on Ishida/Morokuma/Komornicki reaction
                   path stabilization.  The default is .TRUE.
            DELTA  initial step size along the unit bisector, if
                   STABLZ is on.  Default=0.025 Bohr.
            ELBOW  is the collinearity threshold above which the
                   stabilization is skipped.  If the mass weighted
                   gradients at QB and QC are almost collinear, the 
                   reaction path is deemed to be curving very little, 
                   and stabilization isn't needed.  The default is 
                   175.0 degrees.  To always perform stabilization, 
                   input 180.0.
            READQB,EB,GBNORM,GB are energy and gradient data
                   already known at the current IRC point.  If it
                   happens that a run with STABLZ on decides to skip
                   stabilization because of ELBOW, this data will be
                   punched to speed the restart.

          PACE = QUAD   selects quadratic gradient following.
                        Related input is:
            SAB    distance to previous point on the IRC.
            GA     gradient vector at that historical point.
          PACE = AMPC4  selects the fourth order Adams-Moulton
                        variable step predictor-corrector.
                        Related input is:
            GA0,GA1,GA2 which are gradients at previous points.
          PACE = RK4    selects the 4th order Runge-Kutta variable
                        step method.  There is no related input.

          ----- The next two are used by all PACE choices -----

          STRIDE = Determines how far apart points on the reaction
                   path will be.  STRIDE is used to calculate the
                   step taken, according to the PACE you choose.
                   The default is good for the GS2 method, which is
                   very robust.  Other methods should request much
                   smaller step sizes, such as 0.10 or even 0.05.
                   (default = 0.30 sqrt(amu)-Bohr)
          NPOINT = The number of IRC points to be located in this
                   run. The default is to find only the next point.
                   (default = 1)

          ----- The next two tally the reaction path results.  The
                defaults are appropriate for starting from a saddle
                point, restart values are automatically punched out.
          NEXTPT = The number of the next point to be computed.
          STOTAL = Total distance along the reaction path to next
                   IRC point, in mass weighted Cartesian space.

          ----- The following controls jumping off the saddle point.
                If you give a $HESS group, FREQ and CMODE will be
                generated automatically.
          SADDLE = A logical variable telling if the coordinates
                   given in the $DATA deck are at a saddle point
                   (.TRUE.) or some other point lying on the IRC
                   (.FALSE.).  If SADDLE is true, either a $HESS
                   group or else FREQ and CMODE must be given.
                   (default = .FALSE.)  Related input is:
          TSENGY = A logical variable controlling whether the energy
                   and wavefunction are evaluated at the transition
                   state coordinates given in $DATA.  Since you
                   already know the energy from the transition state
                   search and force field runs, the default is .F.
          FORWRD = A logical variable controlling the direction to
                   proceed away from a saddle point. The forward
                   direction is defined as the direction in which
                   the largest magnitude component of the imaginary
                   normal mode is positive. (default =.TRUE.)
          EVIB   = Desired decrease in energy when following the
                   imaginary normal mode away from a saddle point.
                   (default=0.0005 Hartree)
          FREQ   = The magnitude of the imaginary frequency, given
                   in cm**-1.
          CMODE  = An array of the components of the normal mode
                   whose frequency is imaginary, in Cartesian
                   coordinates.  Be careful with the signs!

             You must give FREQ and CMODE if you don't give a $HESS
             group, when SADDLE=.TRUE.  The option of giving these
             two variables instead of a $HESS does not apply to the
             GS2 method, which must have a hessian input, even for
             restarts.  Note also that EVIB is ignored by GS2 runs.
                      * * * * * * * * * * * * * * * * * *
                      For hints about IRC tracking, see
                      the 'further information' section.
                      * * * * * * * * * * * * * * * * * *

          $FFCALC group                 (relevant for RUNTYP=FFIELD)
              This group permits the study of the influence of an
          applied electric field on the wavefunction.  The most
          common finite field calculation applies a sequence of
          fields to extract the linear polarizability and first and
          second order hyperpolarizability.  The method is general,
          and so works for all ab initio wavefunctions in GAMESS.
          EFIELD      = applied electric field strength
                        (default=0.001 a.u.)
          IAXIS and JAXIS specify the orientation of the applied
                          field.  1,2,3 mean x,y,z respectively.
                          The default is IAXIS=3 and JAXIS=0.
            If IAXIS=i and JAXIS=0, the program computes alpha(ii),
            beta(iii), and gamma(iiii) from the energy changes, and
            a few more components from the dipole changes.  Five
            wavefunction evaluations are performed.
            If IAXIS=i and JAXIS=j, the program computes the cross
            terms beta(ijj), beta(iij), and gamma(iijj) from the
            energy changes, and a few more components from the
            dipole changes.  This requires nine evaluations of the
          AOFF        = a flag to permit evaluation of alpha(ij)
                        when the dipole moment is not available.
                        This is necessary only for MP2, and means
                        the off-axial calculation will do 13, not
                        9 energy evaluations.  Default=.FALSE.
          SYM         = a flag to specify when the fields to be
                        applied along the IAXIS and/or JAXIS (or
                        according to EONE below) do not break the
                        molecular symmetry.  Since most fields do
                        break symmetry, the default is .FALSE.
          ONEFLD      = a flag to specify a single applied field
                        calculation will be performed.  Only the
                        energy and dipole moment under this field
                        are computed.  If this option is selected,
                        only SYM and EONE input is heeded.  The
                        default is .FALSE.
          EONE        = an array of the three x,y,z components of
                        the single applied field.
          There are notes on RUNTYP=FFIELD on the next page.

              Finite field calculations require large basis sets,
          and extraordinary accuracy in the wavefunction.  To
          converge the SCF to many digits is sometimes problematic,
          but we suggest you use the input to increase integral
          accuracy and wavefunction convergence, for example
             $SCF    NCONV=10 FDIFF=.FALSE. $END
              In many cases, the applied fields will destroy the
          molecular symmetry.  This means the integrals are
          calculated once with point group symmetry to do the
          initial field free wavefunction evaluation, and then again
          with point group symmetry turned off.  If the fields
          applied do not destroy symmetry, you can avoid this second
          calculation of the integrals by SYM=.TRUE.  This option
          also permits use of symmetry during the applied field
          wavefunction evaluations.
              Examples of fields that do not break symmetry are a
          Z-axis field for an axial point group which is not
          centrosymmetric (i.e. C2v).  However, a second field in
          the X or Y direction does break the C2v symmetry.
          Application of a Z-axis field for benzene breaks D6h
          symmetry.  However, you could enter the group as C6v in
          $DATA while using D6h coordinates, and regain the prospect
          of using SYM=.TRUE.  If you wanted to go on to apply a
          second field for benzene in the X direction, you might
          want to enter Cs in $DATA, which will necessitate the
          input of two more carbon and hydrogen atom, but recovers
          use of SYM=.TRUE.
          Reference: H.A.Kurtz, J.J.P.Stewart, K.M.Dieter
                     J.Comput.Chem.  11, 82-87 (1990).

          $LOCAL group   (relevant for LOCAL=RUEDNBRG or LOCAL=BOYS)
              This group allows input of additional data to control
          the localization methods.  If no input is provided, the
          valence orbitals will be localized as much as possible,
          while still leaving the wavefunction invariant.
          PRTLOC = a flag to control supplemental printout.  The
                   extra output is the rotation matrix to the
                   localized orbitals, and, for the Boys method,
                   the orbital centroids, for the Ruedenberg
                   method, the coulomb and exchange matrices,
                   for the population method, atomic populations.
          MAXLOC = maximum number of localization cycles.  This
                   applies to BOYS or POP methods only.  If the
                   localization fails to converge, a different
                   order of 2x2 pairwise rotations will be tried.
          CVGLOC = convergence criterion.  The default provides
                   LMO coefficients accurate to 6 figures.
          SYMLOC = a flag to restrict localization so that
                   orbitals of different symmetry types are not
                   mixed.  This option is not supported in all
                   possible point groups.  The purpose of this
                   option is to give a better choice for the
                   starting orbitals for GVB-PP or MCSCF runs,
                   without destroying the orbital's symmetry.
                   This option is compatible with each of the
                   3 methods of selecting the orbitals to be
                   included.  (default=.FALSE.)

              The remaining parameters select the orbitals
              which are to be included in the localization.
              You may select from FCORE, NOUTA/NOUTB, or
              NINA/NINB, but may choose only one of these.
          FCORE  = flag to freeze all the chemical core orbitals
                   present.   All the valence orbitals will be
                   localized.  (default=.TRUE.)
                                 * * *
          NOUTA  = number of alpha orbitals to hold fixed in the
                   localization.  (default=0)
          MOOUTA = an array of NOUTA elements giving the numbers of
                   the orbitals to hold fixed.  For example, the
                   input NOUTA=2 MOOUTA(1)=8,13 will freeze only
                   orbitals 8 and 13.  You must enter all the
                   orbitals you want to freeze, including any cores.
                   This variable has nothing to do with cows.
          NOUTB =  number of beta orbitals to hold fixed in -UHF-
                   localizations.  (default=0)
          MOOUTB = same as MOOUTA, except that it applies to the
                   beta orbitals, in -UHF- wavefunctions only.
                                 * * *
          NINA   = number of alpha orbitals which are to be
                   included in the localization.  (default=0)
          MOINA  = an array of NINA elements giving the numbers of
                   the orbitals to be included in the localization.
                   Any orbitals not mentioned will be frozen.
          NINB   = number of -UHF- beta MOs in the localization.
          MOINB  = same as MOINA, except that it applies to the
                   beta orbitals, in -UHF- wavefunctions only.
          N.B.  Since Boys localization needs the dipole integrals,
                do not turn off dipole moment calculation in $ELMOM.
                      * * * * * * * * * * * * * * * * * *
                      For hints about localizations, see
                      the 'further information' section.
                      * * * * * * * * * * * * * * * * * *

                                                       $ELMOM $ELPOT
          $ELMOM group   (not required)
          This group controls electrostatic moments calculation.
          IEMOM  = 0 - skip this property
                   1 - calculate monopole and dipole (default)
                   2 - also calculate quadrupole moments
                   3 - also calculate octupole moments
          WHERE  = COMASS   - center of mass (default)
                   NUCLEI   - at each nucleus
                   POINTS   - at points given in $POINTS.
          OUTPUT = PUNCH, PAPER, or BOTH (default)
          IEMINT = 0 - skip printing of integrals (default)
                   1 - print dipole integrals
                   2 - also print quadrupole integrals
                   3 - also print octupole integrals
                  -2 - print quadrupole integrals only
                  -3 - print octupole integrals only
              The quadrupole and octupole tensors on the printout
          are formed according to the definition of Buckingham.
          Caution: only the first nonvanishing term in the multi-
          ipole charge expansion is independent of the coordinate
          origin chosen, which is normally the center of mass.
          $ELPOT group   (not required)
          This group controls electrostatic potential calculation.
          IEPOT = 0 skip this property (default)
                  1 calculate electric potential
          WHERE  = COMASS   - center of mass
                   NUCLEI   - at each nucleus (default)
                   POINTS   - at points given in $POINTS
                   GRID     - at grid given in $GRID
                   PDC      - at points controlled by $PDC.
          OUTPUT = PUNCH, PAPER, or BOTH (default)
              This property is the electrostatic potential V(a) felt
          by a test positive charge, due to the molecular charge
          density.  A nucleus at the evaluation point is ignored.
          If this property is evaluated at the nuclei, it obeys the
               sum on nuclei(a)   Z(a)*V(a) = 2*V(nn) + V(ne).
          The electronic portion of this property is called the
          diamagnetic shielding.

                                                     $ELDENS $EFIELD
          $ELDENS group   (not required)
          This group controls electron density calculation.
          IEDEN  = 0 skip this property (default)
                 = 1 compute the electron density.
          MORB   = The molecular orbital whose electron density is
                   to be computed.  If zero, the total density is
                   computed.  (default=0)
          WHERE  = COMASS   - center of mass
                   NUCLEI   - at each nucleus (default)
                   POINTS   - at points given in $POINTS
                   GRID     - at grid given in $GRID
          OUTPUT = PUNCH, PAPER, or BOTH (default)
          IEDINT = 0 - skip printing of integrals (default)
                   1 - print the electron density integrals
          $EFIELD group   (not required)
              This group controls electrostatic field and electric
          field gradient calculation.
          IEFLD  = 0 - skip this property (default)
                   1 - calculate field
                   2 - calculate field and gradient
          WHERE  = COMASS   - center of mass
                   NUCLEI   - at each nucleus (default)
                   POINTS   - at points given in $POINTS
          OUTPUT = PUNCH, PAPER, or BOTH (default)
          IEFINT = 0 - skip printing these integrals (default)
                   1 - print electric field integrals
                   2 - also print field gradient integrals
                  -2 - print field gradient integrals only
          The Hellman-Feynman force on a nucleus is the nuclear
          charge multiplied by the electric field at that nucleus.
          The electric field is the gradient of the electric
          potential, and the field gradient is the hessian of the
          electric potential.  The components of the electric field
          gradient tensor are formed in the conventional way, i.e.
          see D.Neumann and J.W.Moskowitz.

                                                       $POINTS $GRID
          $POINTS group   (not required)
              This group is used to input points at which properties
          will be computed.  This first card in the group must
          contain the string ANGS or BOHR, followed by an integer
          NPOINT, the number of points to be used.  The next NPOINT
          cards are read in free format, containing the X, Y, and Z
          coordinates of each desired point.
          $GRID group     (not required)
              This group is used to input a grid (plane through the
          molecule) on which properties will be calculated.
          ORIGIN(i) = coordinates of the lower left corner of
                      the plot.
          XVEC(i)   = coordinates of the lower right corner of
                      the plot.
          YVEC(i)   = coordinates of the upper left corner of
                      the plot.
          SIZE      = grid increment, default is 0.25.
          UNITS     = units of the above four values, it can be
                      either BOHR or ANGS (the default).
          Note that XVEC and YVEC are not necessarily parallel to
          the X and Y axes, rather they are the axes which you
          desire to see plotted by the MEPMAP contouring program.
                   * * * * * * * * * * * * * * * * * * * *
                   For conversion factors, and references
                   see the 'further information' section.
                   * * * * * * * * * * * * * * * * * * * *

          $PDC group     (relevant if WHERE=PDC in $ELPOT)
               This group determines the points at which to compute
          the electrostatic potential, for the purpose of fitting
          atomic charges to this potential.  Constraints on the fit
          which determines these "potential determined charges" can
          include the conservation of charge, the dipole, and the
          PTSEL  = determines the points to be used, choose from
                   SURFACE to use a set of points on the union of
                         spheres forming a scaled van der Waals
                         surface.  The algorithm is unpublished,
                         from Mark Spackman, with results similar
                         to use of a Connolly surface. (default)
                   CHELPG to use a modified version of the CHELPG
                         algorithm, which produces a symmetric
                         grid of points for a symmetric molecule.
          CONSTR = NONE   - no fit is performed.  The potential at
                            the points is instead output according
                            to OUPUT in $ELPOT.
                   CHARGE - the sum of fitted atomic charges is
                            constrained to reproduce the total
                            molecular charge. (default)
                   DIPOLE - fitted charges are constrained to 
                            exactly reproduce the total charge
                            and dipole.
                   QUPOLE - fitted charges are constrained to 
                            exactly reproduce the charge, dipole,
                            and quadrupole.
              Note: the number of constraints cannot exceed
              the number of parameters, which is the number 
              of nuclei.  Planar molecules afford fewer
              constraint equations, namedly two dipole
              constraints and three quadrupole constraints,
              instead of three and five, repectively.
             * * * the next two pertain to PTSEL=SURFACE * * *
          VDWSCL = scale factor for the atomic VDW radii.  The 
                   default seems to be an empirical best value.
                   Values for VDW radii for most elements up to
                   Z=36 are internally stored, all other elements
                   are assumed to have VDW radii of 1.8 Angstroms.
          NSURF  = number of points on each VDW sphere, choose
                   from 12, 42, 162, or 642.  The default leads
                   to points about 0.1 Angstroms apart on the
                   VDW spheres.  (default=642)

             * * * the next two pertain to PTSEL=CHELPG * * *
          RMAX   = maximum distance from any point to the closest
                   atom.  (default=3.0 Angstroms)
          DELR   = distance between points on the grid.
                   (default=0.8 Angstroms)
          MAXPDC = an estimate of the total number of points whose
                   electrostatic potential will be included in the
                   fit.  For example, an upper bound to the points
                   on the VDW surface is the number of atoms times
                   NSURF.  (default=10000)

          CENTER = an array of coordinates at which the momemts
                   were computed.

          DPOLE  = the molecular dipole.

          QPOLE  = the molecular quadrupole.

          PDUNIT = units for the above values.  ANGS will mean
                   that the coordinates are in Angstroms, the
                   dipole in Debye, and the quadrupole in Buck-
                   inghams.  BOHR implies atomic units for all 3.
            Note: it is easier to compute the moments in the
            current run, by setting IEMOM to at least 2 in
            $ELMOM.  However, you could fit experimental data,
            for example, by reading it in here.

               There is no unique way to define fitted atomic
          charges.  Smaller numbers of points at which the electro-
          static potential is fit, changes in VDW radii, asymmetric
          point location, etc. all affect the results.  A useful
          bibliography is
          U.C.Singh, P.A.Kollman, J.Comput.Chem. 5, 129-145(1984)
          L.E.Chirlain, M.M.Francl, J.Comput.Chem. 8, 894-905(1987)
          R.J.Woods, M.Khalil, W.Pell, S.H.Moffatt, V.H.Smith,
             J.Comput.Chem. 11, 297-310(1990)
          C.M.Breneman, K.B.Wiberg, J.Comput.Chem. 11, 361-373(1990)
          K.M.Merz, J.Comput.Chem. 13, 749(1992)

          $STONE group      (optional)
              This group defines the expansion points for Stone's
          distributed multipole analysis (DMA) of the electrostatic
              The DMA takes the multipolar expansion of each overlap
          charge density defined by two gaussian primitives, and
          translates it from the center of charge of the overlap
          density to the nearest expansion point.  Three references
          for the method are
              Stone, Chem.Phys.Lett. 83, 233 (1981)
              Price and Stone, Chem.Phys.Lett. 98, 419 (1983)
              Buckingham and Fowler, J.Chem.Phys. 79, 6426 (1983)
              The existence of a $STONE group in the input is what
          triggers the analysis.  Enter as many lines as you wish,
          in any order, terminated by a $END record.
          ATOM i name, where
                ATOM     is a keyword indicating that a particular
                         atom is selected as an expansion center.
                i        is the number of the atom
                name     is an optional name for the atom. If not
                         entered the name will be set to the name
                         used in the $DATA input.
          ATOMS          is a keyword selecting all nuclei in the
                         molecule as expansion points.  No other
                         input on the line is necessary.
          BOND i j name, where
                BOND     is a keyword indicating that a bond mid-
                         point is selected as an expansion center.
                i,j      are the indices of the atoms defining the
                         bond, corresponding to two atoms in $DATA.
                name     an optional name for the bond midpoint.
                         If omitted, it is set to 'BOND'.

          CMASS          is a keyword selecting the center of mass
                         as an expansion point.  No other input on
                         the line is necessary.
          POINT x y z name, where
                POINT    is a keyword indicating that an arbitrary
                         point is selected as an expansion point.
                x,y,z    are the coordinates of the point, in Bohr.
                name     is an optional name for the expansion
                         point.  If omitted, it is set to 'POINT'.

          The second and third moments on the printout can be 
          converted to Buckingham's tensors by formula 9 of
            A.D.Buckingham, Quart.Rev. 13, 183-214 (1959)
          These can in turn be converted to spherical tensors
          by the formulae in the appendix of
            S.L.Price, et al.  Mol.Phys. 52, 987-1001 (1984)



          $SCRF group                                     (optional)
              The presence of this group in the input turns on the
          use of the Kirkwood-Onsager spherical cavity model for the
          study of solvent effects.  The method is implemented for
          RHF, UHF, ROHF, and GVB wavefunctions and gradients, and
          so can be used with any RUNTYP involving the gradient. 
          The method is not implemented for MCSCF, MP2, CI, any of
          the semiempirical wavefunction, or with analytic hessians.

          DIELEC = the dielectric constant, 80 is often used for H2O

          RADIUS = the spherical cavity radius, in Angstroms

          G      = the proportionality constant relating the solute
                   molecule's dipole to the strength of the reaction
                   field.  Since G can be calculated from DIELEC and
                   RADIUS, do not give G if they were given.

          Some references on this subject are

          J.G.Kirkwood  J.Chem.Phys. 2, 351 (1934)
          L.Onsager  J.Am.Chem.Soc. 58, 1486 (1936)
          O.Tapia, O.Goscinski  Mol.Phys. 29, 1653 (1975)
          M.M.Karelson, A.R.Katritzky, M.C.Zerner
                Int.J.Quantum Chem., Symp. 20, 521-527 (1986)
          K.V.Mikkelsen, H.Aagren, H.J.Aa.Jensen, T.Helgaker
                J.Chem.Phys. 89, 3086-3095 (1988)
          M.W.Wong, M.J.Frisch, K.B.Wiberg  
                J.Am.Chem.Soc. 113, 4776-4782 (1991)
          M.Szafran, M.M.Karelson, A.R.Katritzky, J.Koput, M.C.Zerner
                J.Comput.Chem. 14, 371-377 (1993)
          M.Karelson, T.Tamm, M.C.Zerner
                J.Phys.Chem. 97, 11901-11907 (1993)

          It is useful to think about Figures 1 and 2 of the Szafran
          reference when you are picking DIELEC and RADIUS values!
          GAMESS implements Zerner's Method A, as described in this
          paper.  The total solute energy includes the Born term, if 
          the solute is an ion.

          It is useful to think about Table VI of the Mikkelsen
          paper, lest you think a solute consists only of a charge
          and a dipole moment.

          $ECP group               (required if ECP=READ in $CONTRL)
              This group lets you read in effective core potentials,
          for some or all of the atoms in the molecule.  You can
          use built in potentials for some of the atoms if you like.
          This is a free format (positional) input group.
          *** Give a card set -1-, -2-, and -3- for each atom ***
          -card 1-    PNAME, PTYPE, IZCORE, LMAX
          PNAME is a 8 character descriptive tag for this potential.
                If it is repeated for a subsequent atom, no other
                information need be given on this card, and cards
                -2- and -3- may also be skipped.  The information
                will be copied from the first atom by this PNAME.
          PTYPE = GEN    a general potential should be read.
                = SBK    look up the Stevens/Basch/Krauss/Jasien/
                         Cundari potential for this type of atom.
                = HW     look up the Hay/Wadt built in potential
                         for this type of atom.
                = NONE   treat all electrons on this atom.
          IZCORE is the number of core electrons to be removed.
          LMAX   is the maximum angular momentum occupied in the
                 core orbitals being removed (usually).  Give
                 IZCORE and LMAX only if PTYPE is GEN.
          *** For the first occurence of PNAME, if PTYPE is GEN, ***
          *** then give cards -2- and -3-.  Otherwise go to -1-. ***
          *** Card sets -2- and -3- are repeated LMAX+1 times    ***
              The potential U(LMAX+1) is given first,
              followed by U(L)-U(LMAX+1), for L=1,LMAX.
          -card 2-    NGPOT
          NGPOT is the number of Gaussians in this part of the
                local effective potential.
          -card 3-    CLP,NLP,ZLP   (repeat this card NGPOT times)
          CLP is the coefficient of this Gaussian in the potential.
          NLP is the power of r for this Gaussian.
          ZLP is the exponent of this Gaussian.
                                * * *
          By far the easiest way to use the SBK potential for all
          atoms in the formic acid molecule is to request ECP=SBK in
          $CONTRL.  But the next page shows two alternatives.

          The first way is to look up the program's internally
          stored SBK potentials one atom at a time:
          C-ECP SBK
          H-ECP NONE
          O-ECP SBK
          H-ECP NONE
          The second oxygen duplicates the first, no core electrons
          are removed for hydrogen.  The order of the atoms must
          follow that generated by $DATA.  Note PTYPE allows you to
          type in one or more atoms explicitly, while using built in
          data for some other atoms.
          The second example reads all SBK potentials explicitly:
          C-ECP GEN 2 1
          1      ----- CARBON U(P) -----
           -0.89371  1  8.56468
          2      ----- CARBON U(S-P) -----
            1.92926  0  2.81497
           14.88199  2  8.11296
          H-ECP NONE
          O-ECP GEN 2 1
          1      ----- OXYGEN U(P) -----
           -0.92550  1 16.11718
          2      ----- OXYGEN U(S-P) -----
            1.96069  0  5.05348
           29.13442  2 15.95333
          H-ECP NONE
          Again, the 2nd oxygen copies from the first.  It is handy
          to use the rest of card -2- as a descriptive comment.
          As a final example, for antimony we have LMAX=3 (there
          are core d's).  One must first enter U(f), followed by
          U(s)-U(f), U(p)-U(f), U(d)-U(f).


          $INTGRL group                                   (optional)
              This group controls AO integral formats.  It should 
          probably never be given, as the program always picks 
          sensible values.

           SCHWRZ = a flag to activate use of the Schwarz inequality
                    to predetermine small integrals.  There is no
                    loss of accuracy when choosing this option, and
                    there are appreciable time savings for bigger
                    molecules.  Default=.TRUE. for over 5 atoms, or
                    for direct SCF, and is .FALSE. otherwise.
           NOPK   = 0 PK integral option on
                      (default for -RHF-, -UHF-, -ROHF-, -GVB-)
                  = 1 PK option off (default for -CI- or -MCSCF-,
                      or for any MP2, analytic hessian, etc. run)
           NORDER = 0 (default)
                  = 1 Sort integrals into canonical order.  There
                      is little point in selecting this option, as
                      no part of GAMESS requires ordered integrals.
                      See also NSQUAR.
           NINTMX =   Maximum no. of integrals in a record block.
                      (default = 2725 for P file, = 1635 for PK)
                      The record length of the file is either
                      (  nwdvar+1) + nintmx/nwdvar + 1  for P files.
                      (2*nwdvar+1) + nintmx/nwdvar + 1 for PK files.
                The following parameters control the integral sort.
                (values given are defaults)

           NSQUAR = 0 Sorted integrals will be in triangular
                      canonical order (default)
                  = 1 instead sort to square canonical order.
           NDAR   = Number of direct access logical records to be
                    used for the integral sort (default=2000)
           LDAR   = Length of direct access records (site dependent)
           NBOXMX =  200   Maximum number of bins.
           NWORD  =    0   Memory to be used (default=all of it).
           NOMEM  =    0   If non-zero, force external sort.
                The following parameters control integral restarts
                (values given are defaults)
           IST=    1      JST=    1    KST=    1    LST=    1
           NREC=   1      INTLOC= 1

               The remaining groups, with the exception of $TRANS,
               apply only to -CI- and -MCSCF- runs.

          $CIINP group                 (optional, relevant for -CI-)
              This group is the control box for Graphical Unitary
          Group Approach (GUGA) -CI- calculations.  Each module
          executed potentially requires a further input group
          described later.  Note that NPFLG (and only NPFLG) will
          also apply to MCSCF runs.
          IREST = n    Restart the -CI- at stage NRNFG(n).
          NRNFG = An array of 10 switches controlling which
                  modules are run during the -CI-.
                  1 means execute the module, 0 means don't.
                  NRNFG(8-10) are not used.
            NRNFG(1) = Generate the distinct row table.  The DRT is
                       the GUGA configuration list. See $DRT.
            NRNFG(2) = Transform the integrals. See $TRANS.
            NRNFG(3) = Sort integrals and calculate the Hamiltonian
                       matrix. See $CISORT and $GUGEM. (default=1)
            NRNFG(4) = Diagonalize the Hamiltonian matrix.
                       See $GUGDIA. (default=1)
            NRNFG(5) = Construct the one electron density matrix,
                       and generate NO's. See $GUGDM. (default=1)
            NRNFG(6) = Construct the two electron density matrix.
                       See $GUGDM2. (default=0)
            NRNFG(7) = Construct the Lagrangian of the CI function.
                       Requires DM2.  See $LAGRAN. (default=0)
          NPFLG = An array of 10 switches to produce debug printout.
                  There is a one to one correspondance to NRNFG, set
                  to 1 for output. (default = 0,0,0,0,0,0,0,0,0,0)
            NPFLG(8) = debug for the MCSCF Newton-Raphson step.
                     * * * * * * * * * * * * * * * * * * * *
                     For hints on how to do -CI- and -MCSCF-
                      see the 'further information' section
                     * * * * * * * * * * * * * * * * * * * *

          $DRT group                  (required for -CI- or -MCSCF-)
              This group describes the -CI- or -MCSCF- wavefunction.
          The distinct row table is the means by which the Graphical
          Unitary Group Approach (GUGA) names the configurations.
              There is no default for GROUP, and you must choose one
          of FORS, FOCI, SOCI, or IEXCIT.
          GROUP = the name of the point group to be used.  This is
                  usually the same as that in $DATA, except for
                  RUNTYP=HESSIAN, when it must be C1.  Choose from
                  the following: C1, C2, CI, CS, C2V, C2H, D2, D2H,
                  C4V, D4, D4H.  If your $DATA group is not listed,
                  choose only C1 here.
          FORS  = flag specifying the Full Optimized Reaction Space
                  set of configuration should be generated.  This
                  is usually set true for MCSCF runs.
          FOCI  = flag specifying first order CI.  In addition to
                  the FORS configurations, all singly excited CSFs
                  from the FORS reference are included.
          SOCI  = flag specifying second order CI.  In addition to
                  the FORS configurations, all singly and doubly
                  excited configurations from the FORS reference
                  are included.  Default=.FALSE.
          IEXCIT= electron excitation level, for example 2 will
                  lead to a singles and doubles CI.  This variable
                  is computed by the program if FORS, FOCI, or
                  SOCI is chosen, otherwise it must be entered.
            Note: there is a known bug (with no known fix) when
                  using IEXCIT with MULT greater than 1.  The
                  program generates extra CSFs, with excitations
                  higher than that given here.  The same thing
                  happens with singly occupied AOS and BOS MOs.
            * * the next variables define the single reference * *
              The single configuration reference is defined by
          filling in the orbitals by each type, in the order shown.
          The default for each type is 0.
                 Core orbitals, which are always doubly occupied:
          NMCC = number of MCSCF core MOs.
          NFZC = number of CI frozen core MOs.

                 Internal orbitals, which are partially occupied:
          NDOC = number of doubly occupied MOs in the reference.
          NAOS = number of alpha occupied MOs in the reference,
                 which are singlet coupled with a corresponding
                 number of NBOS orbitals.
          NBOS = number of beta spin singly occupied MOs.
          NALP = number of alpha spin singly occupied MOs in the
                 reference, which are coupled high spin.
          NVAL = number of empty MOs in the reference.
                 External orbitals, occupied only in FOCI or SOCI:
          NEXT = number of external MOs.  If given as -1, this will
                 be set to all remaining orbitals (apart from any
                 frozen virtual orbitals).
          NFZV = number of frozen virtual MOs, never occupied.
                 * * * the final choices are seldom used * * *
          INTACT= flag to select the interacting space option.
                  The CI will include only those spin couplings
                  which have a nonvanishing matrix element with
                  the reference configuration.
          MXNINT = Buffer size for sorted integrals. (default=10000)
          MXNEME = Buffer size for energy matrix.  (default=1085)
          NPRT   = Configuration printout control switch.
                   This can consume a HUMUNGUS amount of paper!
                   0 = no print (default)
                   1 = print electron occupancies, one per line.
                   2 = print determinants in each CSF.

          $MCSCF group                        (optional for -MCSCF-)
              This group controls the MCSCF Newton-Raphson orbital
          improvement step.
           METHOD = DM2 selects a density driven approach to the
                    construction of the Newton-Raphson matrices.
                  = TEI selects 2e- integral driven NR construction.
                  = FORMULA selects formula driven NR construction.
                    (the FORMULA method is now an inactive option)
                    See the 'further information' section for more
                    details concerning these methods.
           MAXIT  = Maximum number of iterations (default=30)
           MICIT  = Maximum number of microiterations within a
                    single MCSCF iteration. (default=1)
           ACURCY = the major convergence criterion, the maximum
                    permissible asymmetry in the Lagrangian matrix.
           ENGTOL = a secondary convergence criterion, the run is
                    considered converged when the energy change is
                    smaller than this value. (default=1.0E-10)
           DAMP   = damping factor, this is adjusted by the program
                    as necessary.  (default=0.0)
           FORS   = a flag to specify that the MCSCF function is of
                    the Full Optimized Reaction Space type, which is
                    sometimes known as CAS-SCF.  The default is to
                    declare act-act rotations redundant.  This must
                    be set false if FORS was not set in $DRT.
           CANONC = a flag to cause formation of the closed shell
                    Fock operator, and generation of canonical core
                    orbitals.  This will order the MCC core by their
                    orbital energies.  (default=.TRUE.)


           EKT    = a flag to cause generation of extended Koopmans'
                    theorem orbitals and energies.  (Default=.FALSE.)

              For this option, see R.C.Morrison and G.Liu,
              J.Comput.Chem., 13, 1004-1010 (1992).  Note that
              the process generates non-orthogonal orbitals, as
              well as physically unrealistic energies for the
              weakly occupied MCSCF orbitals.  The method is
              meant to produce a good value for the first I.P.
                --- these options are less commonly used ---
                The first 3 pertain only to METHOD=FORMULA.
           FMLFIL = a flag specifying the formula tape has been
                    saved for reuse in this run.  The default is
                    to regenerate this file on the 1st iteration
                    at the 1st geometry run.  (default=.FALSE.)
           LDAR   = length of each direct access record during the
                    sort of the formula tape.  (default is machine
           NDAR   = number of direct access records.  The value used
                    will be the larger of the program's guess of the
                    actual number, and the input value (default=0).
           NWORD  = The maximum memory to be used in the NR step.
                    The default is to use all available memory.
           FCORE  = a flag to freeze optimization of the MCC core
                    orbitals, which is useful in preparation for
                    RUNTYP=TRANSITN jobs.  Setting this flag will
                    automatically force CANONC false.  This option
                    is incompatible with gradients, so can only be
                    used with RUNTYP=ENERGY. (default=.FALSE.)

                    --- the final two are seldom used ---
           NORB   = the number of orbitals to be included in the
                    optimization, the default is to optimize with
                    respect to the entire basis.  This option is
                    incompatible with gradients, so can only be used
                    with RUNTYP=ENERGY.  (default=number of AOs
                    given in $DATA).
           NOROT  = an array of up to 250 pairs of orbital rotations
                    to be omitted from the NR optimization process.
                    The program automatically deletes all core-core
                    rotations, all act-act rotations if FORS=.T.,
                    and all core-act and core-virt rotations if
                    FCORE=.T.  Additional rotations are input as
                    I1,J1,I2,J2... to exclude rotations between
                    orbital I running from 1 to NORB, and J running
                    up to the smaller of I or NVAL in $TRANS.

          $TRANS group                (optional for -CI- or -MCSCF-)
                                     (relevant to analytic hessians)
                                   (relevant to energy localization)
               This group controls the integral tranformation. 
          There is little reason to give any but the first variable.

           DIRTRF = a flag to recompute AO integrals rather than
                    storing them on disk.  The default is .FALSE.
                    for MCSCF and CI runs.  If your job reads $SCF,
                    and you select DIRSCF=.TRUE. in that group, a
                    direct transformation will be done, no matter
                    how DIRTRF is set.

              Note that the transformation may do many passes over 
              the AO integrals for large basis sets, and thus the
              direct recomputation of AO integrals can be very time 

           MPTRAN = method to use for the integral transformation.
                    the default is try 0, then 1, then 2.
                    0 means use the incore method
                    1 means use the segmented method.  This is the
                      only method that works in parallel.
                    2 means use the alternate method, which uses
                      less memory than 2, but requires an extra
                      large disk file.
           NWORD  = Number of words of fast memory to allow.  Zero
                    uses all available memory. (default=0)
           CUTTRF = Threshold cutoff for keeping transformed two
                    electron integrals.  (default= 10**(-9))

                                              $CISORT $GUGEM $GUGDIA

          $CISORT group    (optional, relevant for -CI- and -MCSCF-)
               This group provides further control over the sorting
          of the transformed integrals.
          NDAR   = Number of direct access records.
                   (default = 2000)
          LDAR   = Length of direct access record (site dependent)
          NBOXMX = Maximum number of boxes in the sort.
                   (default = 200)
          NWORD  = Number of words of fast memory to use in this
                   step.  A value of 0 results in automatic use of
                   all available memory.  (default = 0)
          NOMEM  = 0 (set to one to force out of memory algorithm)
          $GUGEM group      (optional, relevant for -CI- or -MCSCF-)
              This group provides further control over the
          calculation of the energy (Hamiltonian) matrix.
          CUTOFF = Cutoff criterion for the energy matrix.
          NWORD  = not used.
          $GUGDIA group     (optional, relevant for -CI- or -MCSCF-)
               This group provides control over the Davidson method
          diagonalization step.
          NSTATE = Number of CI states to be found. (default=1)
                   You can solve for any number of states, but only
                   100 can be saved for subsequent sections, such
                   as state averaging.
          PRTTOL = Printout tolerance for CI coefficients
                   (default = 0.05)
          MXXPAN = Maximum no. of expansion basis vectors used
                   before the expansion basis is truncated.

          ITERMX = Maximum number of iterations (default=50)
          CVGTOL = Convergence criterion for Davidson eigenvector
                   routine.  This value is proportional to the
                   accuracy of the coeficients of the eigenvector(s)
                   found.  The energy accuracy is proportional to
                   its square.  (default = 1.0E-5)
          NWORD  = Number of words of fast memory to use in this
                   step.  A value of zero results in the use of all
                   available memory.  (default = 0)
          MEMMAX = 1000*limgiv + limh where limgiv is the largest
                   matrix diagonalized via Givens-Householder
                   (default=50) and limh is the dimension of the
                   largest Hamiltonian that may be memory resident
          NIMPRV = Maximum no. of eigenvectors to be improved every
                   iteration. (default = nstate)
          NSELCT = Determines initial guess to eigenvectors.
                   = 0 ->  Unit vectors corresponding to the NSTATE
                           lowest diagonal elements and any diagonal
                           elements within SELTHR of them. (default)
                   < 0 ->  First abs(NSELCT) unit vectors.
                   > 0 ->  use NSELCT unit vectors corresponding to
                            the NSELCT lowest diagonal elements.
          SELTHR = Guess selection threshold when NSELCT=0.
          NEXTRA = Number of extra expansion basis vectors to be
                   included on the first iteration.  NEXTRA is
                   decremented by one each iteration.  This may be
                   useful in "capturing" vectors for higher states.
          KPRINT = Print flag bit vector used when
                   NPFLG(4)=1 in the $CIINP group       (default=8)
                   value  1 bit 0 print final eigenvalues
                   value  2 bit 1 print final tolerances
                   value  4 bit 2 print eigenvalues and tolerances
                                  at each truncation
                   value  8 bit 3 print eigenvalues every iteration
                   value 16 bit 4 print tolerances every iteration

                                                      $GUGDM $LAGRAN
          $GUGDM group                 (optional, relevant for -CI-)
               This group provides further control over formation of
          the 1-particle density matrix.
           NFLGDM = Controls each state's density formation.
                    0 -> do not form density for this state.
                    1 -> form density and natural orbitals for this
                         state, print and punch occ.nums. and NOs.
                    2 -> same as 1, plus print density over MOs.
                    (default=1,99*0, meaning G.S. NOs only)
                    See also NSTATE in $GUGDIA.  Note that forming
                    the 1-particle density for a state is negligible
                    against the diagonalization time for that state.
           IROOT  = The -CI- root whose density matrix is saved on
                    the direct access dictionary file for later
                    computation of properties.  (default=1)
           IBLOCK = Density blocking switch. If nonzero, the off
                    diagonal block of the density below row IBLOCK
                    will be set to zero before the (approximate)
                    natural orbitals are found.  One use for this is
                    to keep the internal and external orbitals in a
                    FOCI or SOCI calculation from mixing.
           NWORD  = Number of words of fast memory to use in this
                    step.  A value of zero uses all available memory
          $LAGRAN group                (optional, relevant for -CI-)
               This group provides further control over formation of
          the -CI- Lagrangian.  This is meant, someday, maybe, for
          CI gradient computation, and is seldom used now.
           NDAR  = 4000
           LDAR  = Length of each direct access record
                   (default is site dependent)
           NWORD =   0
           NORB  = number in $TRANS
           NOMEM =   0

                                                     $GUGDM2 $TRFDM2
          $GUGDM2 group     (optional, relevant for -CI- or -MCSCF-)
               This group provides control over formation of the
          2-particle density matrix.
          WSTATE = An array of up to 100 weights to be given to the
                   2 body density of each state in forming the DM2.
                   The default is to optimize a pure ground state.
                   A small amount of the ground state can help the
                   convergence of excited states greatly.
                   Gradient runs are possible only with pure states.
                   Be sure to set NSTATE in $GUGDIA appropriately!
          CUTOFF = Cutoff criterion for the 2nd-order density.
                   (default = 1.0E-9)
          NWORD  = Number of words of fast memory to use in sorting
                   the DM2.  The default uses all available memory.
          NOMEM  = 0 uses in memory sort, if possible.
                 = 1 forces out of memory sort.
          NDAR   = Number of direct access records. (default=4000)
          LDAR   = Length of direct access record (site dependent)
          NBOXMX = Maximum no. of boxes in the sort. (default=200)
          $TRFDM2 group    (optional, relevant for -MCSCF- GRADIENT,
                           IRC, OPTIMIZE, HESSIAN, or SADPOINT runs)
               This group provides control over the transformation
          of the 2-particle density matrix from MO to AO basis.
           NDAR  = 2000
           LDAR  = Length of direct access record (site dependent)
           NBOXMX= 200
           NWORD =   0 Use all available memory or amount specified.
           CUTOFF= 1.0E-9
           NPFLG =   0
           NOMEM =   0 (default - choose method using least memory)
                     1 Transform DM2 in memory, sort out of memory.
                     2 External transformation and sort.

          $TRANST group               (relevant for RUNTYP=TRANSITN)
          (only for SCFTYP=CI)        (relevant for RUNTYP=SPINORBT)
              This group controls the evaluation of the radiative
          transition moment, or spin orbit coupling.  The defaults
          assume that there is one common set of orbitals, all of
          which are occupied.  This would be true for a conventional
          CI-SD run.  The program can also use two separately
          optimized MO sets, provided certain conditions are met.
          NUMVEC = the number of different MO sets. This can be
                   either 1 or 2. (default=1)
          NUMCI  = the number of different CI calculations to do.
                   This should equal NUMVEC for TRANSITN, and is
                   usually 2 for SPINORBT (1 if the multiplicities
                   of the states are the same). (default=1)
          NOCC   = the number of occupied orbitals.  The default is
                   the number of AOs given in $DATA.
          NFZC   = the number of identical core orbitals found in
                   the two MO sets, when NUMVEC=2.  When NUMVEC
                   is 1, it should equal NOCC.  The default is the
                   number of AOs given in $DATA.
          IROOTS = array containing the two CI states for which the
                   transition moments are to be found.  The default
                   is the ground and first excited state, namely
          ZEFF   = an array of effective nuclear charges.  This
                   pertains only to SPINORBT runs.  The defaults
                   are the true nuclear charges.
          PRTCMO = flag to control printout of the corresponding
                   orbitals.  (default is .FALSE.)
          TOLZ   = MO coefficient zero tolerance (as for $GUESS).
          TOLE   = MO coefficient equating tolerance (as for
                   $GUESS).  (default=1.0E-5)
                  * * * * * * * * * * * * * * * * * * * * *
                  For information on TRANSITN and SPINORBT,
                  see the 'further information' section.
                  * * * * * * * * * * * * * * * * * * * * *