Molecular Properties and Conversion Factors
 
These two papers are of general interest:
 A.D.Buckingham, J.Chem.Phys. 30, 1580-1585(1959).
 D.Neumann, J.W.Moskowitz J.Chem.Phys. 49, 2056-2070(1968).
The first deals with multipoles, and the second with other 
properties such as electrostatic potentials.
 
All units are derived from the atomic units for distance 
and the monopole electric charge, as given below.
 
distance                 1 au = 5.291771E-09 cm
 
monopole                 1 au = 4.803242E-10 esu
                        1 esu = sqrt(g-cm**3)/sec
 
dipole                   1 au = 2.541766E-18 esu-cm
                      1 Debye = 1.0E-18 esu-cm
 
quadrupole               1 au = 1.345044E-26 esu-cm**2
                 1 Buckingham = 1.0E-26 esu-cm**2
 
octopole                 1 au = 7.117668E-35 esu-cm**3
 
electric potential       1 au = 9.076814E-02 esu/cm
 
electric field           1 au = 1.715270E+07 esu/cm**2
                  1 esu/cm**2 = 1 dyne/esu
 
electric field gradient  1 au = 3.241390E+15 esu/cm**3
 
The atomic unit for the total electron density is 
electron/bohr**3, but 1/bohr**3 for an orbital density.
 
The atomic unit for spin density is excess alpha spins per 
unit volume, h/4*pi*bohr**3.  Only the expectation value is 
computed, with no constants premultiplying it.
 
IR intensities are printed in Debye**2/amu-Angstrom**2. 
These can be converted into intensities as defined by 
Wilson, Decius, and Cross's equation 7.9.25, in km/mole, by 
multiplying by 42.255.  If you prefer 1/atm-cm**2, use a 
conversion factor of 171.65 instead.  A good reference for 
deciphering these units is A.Komornicki, R.L.Jaffe 
J.Chem.Phys. 1979, 71, 2150-2155.  A reference showing how 
IR intensities change with basis improvements at the HF 
level is Y.Yamaguchi, M.Frisch, J.Gaw, H.F.Schaefer, 
J.S.Binkley, J.Chem.Phys. 1986, 84, 2262-2278.  

Raman activities in A**4/amu multiply by 6.0220E-09 for
units of cm**4/g.  One of the many sources explaining how 
activity relates to intensity is D.Michalska, R.Wysokinski 
Chem.Phys.Lett. 403, 211-217(2005)


Polarizabilities

Static polarizabilities are named alpha, beta, and gamma; 
and are called the polarizability, hyperpolarizability, and 
second hyperpolarizability.  They are the 2nd, 3rd, and 4th 
derivatives of the energy with respect to applied uniform 
electric fields, with the 1st derivative being the dipole 
moment!  It might be worth mentioning that a uniform field 
can be applied using $EFIELD, if you wish to develop custom 
usages, but $EFIELD is not to be given for any kind of run 
discussed below.

The general approach to computing static polarizabilities 
is numerical differentiation, namely RUNTYP=FFIELD, which 
should work for any energy method provided by GAMESS.  A 
sequence of computations with fields applied in the x, y, 
and/or z directions will generate the tensors.  See $FFCALC 
for details.  Analytic computation of all three tensors is 
available for SCFTYP=RHF only, see RUNTYP=TDHF and $TDHF.  
If you need to know just the alpha polarizability, see 
POLAR in $CPHF during any analytic hessian job.

A break down of the static alpha polarizability in terms of 
contributions from individual localized orbitals can be 
obtained by setting POLDCM=.TRUE. in $LOCAL.  Calculation 
will be by analytic means, unless POLNUM in that group is 
selected.  This option is available only for SCFTYP=RHF.  
The keyword LOCHYP in $FFCALC gives a similar analysis for 
all three static polarizabilities, determined by numerical 
differentiation.

Polarizabilities in a static electric field differ from 
those in an oscillating field such as a laser produces.  
For RHF cases only, a variety of frequency dependent 
polarizabilities alpha, beta, and gamma can be generated, 
depending on the experiment.  A particularly easy one to 
describe is 'second harmonic generation', governed by a 
beta tensor describing the absorption of two photons and 
the emission of one at doubled frequency.  See RUNTYP=TDHF, 
and papers listed under $TDHF, for many other non-linear 
optical experiments.

Nuclear derivatives of the dipole moment and the various 
polarizabilities are also of interest.  For example, 
knowledge of the derivative of the dipole with respect to 
nuclear coordinates yields the IR intensity.  Similarly, 
the nuclear derivative of the static alpha polarizability 
gives Raman intensities: see RUNTYP=RAMAN.  Analytic of 1st 
or 2nd nuclear derivatives of static or frequency dependent 
polarizabilities are available for SCFTYP=RHF only, see 
RUNTYP=TDHFX and $TDHFX, giving rise to experimental 
observations such as resonance Raman and hyper-Raman.

Finally, instead of considering polarizabilities as a 
function of real frequencies, they can be considered to be 
dependent on the imaginary frequency.  The imaginary 
frequency dependent alpha polarizability can be computed 
analytically for SCFTYP=RHF only, using POLDYN=.TRUE. in 
$LOCAL.  Integration of this quantity over the imaginary 
frequency domain can be used to extract C6 dispersion 
constants.

Polarizabilities are tensor quantities.  There are a number 
of different ways to define them, and various formulae to 
extract "scalar" and "vector" quantites from the tensors.  
A good reference for learning how to compare the output of 
a theoretical program to experiment is
    A.Willetts, J.E.Rice, D.M.Burland, D.P.Shelton
    J.Chem.Phys. 97, 7590-7599(1992)




Localized Molecular Orbitals
 
    Three different orbital localization methods are
implemented in GAMESS.  The energy and dipole based
methods normally produce similar results, but see
M.W.Schmidt, S.Yabushita, M.S.Gordon in J.Chem.Phys.,
1984, 88, 382-389 for an interesting exception.  You can
find references to the three methods at the beginning of
this chapter.
 
    The method due to Edmiston and Ruedenberg works by
maximizing the sum of the orbitals' two electron self
repulsion integrals.  Most people who think about the
different localization criteria end up concluding that
this one seems superior.  The method requires the two
electron integrals, transformed into the molecular orbital
basis.  Because only the integrals involving the orbitals
to be localized are needed, the integral transformation is
actually not very time consuming.
 
    The Boys method maximizes the sum of the distances
between the orbital centroids, that is the difference in
the orbital dipole moments.
 
    The population method due to Pipek and Mezey maximizes
a certain sum of gross atomic Mulliken populations.  This
procedure will not mix sigma and pi bonds, so you will not
get localized banana bonds.  Hence it is rather easy to
find cases where this method give different results than
the Ruedenberg or Boys approach.
 
    GAMESS will localize orbitals for any kind of RHF, UHF,
ROHF, or MCSCF wavefunctions.  The localizations will
automatically restrict any rotation that would cause the
energy of the wavefunction to be changed (the total
wavefunction is left invariant).  As discussed below,
localizations for GVB or CI functions are not permitted.
 
    The default is to freeze core orbitals.  The localized
valence orbitals are scarcely changed if the core orbitals
are included, and it is usually convenient to leave them
out.  Therefore, the default localizations are:  RHF
functions localize all doubly occupied valence orbitals.
UHF functions localize all valence alpha, and then all
valence beta orbitals.  ROHF functions localize all valence
doubly occupied orbitals, and all singly occupied orbitals, 
but do not mix these two orbital spaces.  MCSCF functions
localize all valence MCC type orbitals, and localize all
active orbitals, but do not mix these two orbital spaces.
To recover the invariant MCSCF function, you must be using
a FORS=.TRUE. wavefunction, and you must set GROUP=C1 in
$DRT, since the localized orbitals possess no symmetry.
 
    In general, GVB functions are invariant only to
localizations of the NCO doubly occupied orbitals.  Any
pairs must be written in natural form, so pair orbitals
cannot be localized.  The open shells may be degenerate, so
in general these should not be mixed.  If for some reason
you feel you must localize the doubly occupied space, do a
RUNTYP=PROP job.  Feed in the GVB orbitals, but tell the
program it is SCFTYP=RHF, and enter a negative ICHARG so
that GAMESS thinks all orbitals occupied in the GVB are
occupied in this fictitous RHF.  Use NINA or NOUTA to
localize the desired doubly occupied orbitals.  Orbital
localization is not permitted for CI, because we cannot
imagine why you would want to do that anyway.
 
    Boys localization of the core orbitals in molecules
having elements from the third or higher row almost never
succeeds.  Boys localization including the core for second
row atoms will often work, since there is only one inner
shell on these.  The Ruedenberg method should work for any
element, although including core orbitals in the integral
transformation is more expensive.
 
    The easiest way to do localization is in the run which
generates the wavefunction, by selecting LOCAL=xxx in the
$CONTRL group.  However, localization may be conveniently
done at any time after determination of the wavefunction,
by executing a RUNTYP=PROP job.  This will require only
$CONTRL, $BASIS/$DATA, $GUESS (pick MOREAD), the converged
$VEC, possibly $SCF or $DRT to define your wavefunction,
and optionally some $LOCAL input.
 
    There is an option to restrict all rotations that would
mix orbitals of different symmetries.  SYMLOC=.TRUE. yields
only partially localized orbitals, but these still possess
symmetry.  They are therefore very useful as starting
orbitals for MCSCF or GVB-PP calculations.  Because they
still have symmetry, these partially localized orbitals run
as efficiently as the canonical orbitals.  Because it is
much easier for a user to pick out the bonds which are to
be correlated, a significant number of iterations can be
saved, and convergence to false solutions is less likely.

                          * * *

    The most important reason for localizing orbitals is to 
analyze the wavefunction.  A simple example is to look at 
shapes of the orbitals with the MacMolPlt program.  Or, you 
might read the localized orbitals in during a RUNTYP=PROP 
job to examine their Mulliken populations.

    Localized orbitals are a particularly interesting way 
to analyze MCSCF computations.  The localized orbitals may 
be oriented on each atom (see option ORIENT in $LOCAL) to 
direct the orbitals on each atom towards their neighbors 
for maximal bonding, and then print a bond order analysis. 
The orientation procedure is newly programmed by J.Ivanic 
and K.Ruedenberg, to deal with the situation of more than 
one localized orbital occuring on any given atom.  Some 
examples of this type of analysis are
    D.F.Feller, M.W.Schmidt, K.Ruedenberg
       J.Am.Chem.Soc.  104, 960-967 (1982)
    T.R.Cundari, M.S.Gordon
       J.Am.Chem.Soc.  113, 5231-5243 (1991)
    N.Matsunaga, T.R.Cundari, M.W.Schmidt, M.S.Gordon
       Theoret.Chim.Acta  83, 57-68 (1992).

    In addition, the energy of your molecule can be 
partitioned over the localized orbitals so that you may
be able to understand the origin of barriers, etc.  This
analysis can be made for the SCF energy, and also the MP2
correction to the SCF energy, which requires two separate
runs.  An explanation of the method, and application to
hydrogen bonding may be found in J.H.Jensen, M.S.Gordon,
J.Phys.Chem. 99, 8091-8107(1995).

    Analysis of the SCF energy is based on the localized
charge distribution (LCD) model: W.England and M.S.Gordon,
J.Am.Chem.Soc. 93, 4649-4657 (1971).  This is implemented
for RHF and ROHF wavefunctions, and it requires use of 
the Ruedenberg localization method, since it needs the
two electron integrals to correctly compute energy sums.
All orbitals must be included in the localization, even
the cores, so that the total energy is reproduced.

    The LCD requires both electronic and nuclear charges
to be partitioned.  The orbital localization automatically
accomplishes the former, but division of the nuclear
charge may require some assistance from you.  The program
attempts to correctly partition the nuclear charge, if you
select the MOIDON option, according to the following: a 
Mulliken type analysis of the localized orbitals is made.
This determines if an orbital is a core, lone pair, or 
bonding MO.  Two protons are assigned to the nucleus to 
which any core or lone pair belongs.  One proton is 
assigned to each of the two nuclei in a bond.  When all 
localized orbitals have been assigned in this manner, the 
total number of protons which have been assigned to each 
nucleus should equal the true nuclear charge.

    Many interesting systems (three center bonds, back-
bonding, aromatic delocalization, and all charged species)
may require you to assist the automatic assignment of
nuclear charge.  First, note that MOIDON reorders the
localized orbitals into a consistent order: first comes
any core and lone pair orbitals on the 1st atom, then
any bonds from atom 1 to atoms 2, 3, ..., then any core
and lone pairs on atom 2, then any bonds from atom 2 to
3, 4, ..., and so on.  Let us consider a simple case 
where MOIDON fails, the ion NH4+.  Assuming the nitrogen
is the 1st atom, MOIDON generates
     NNUCMO=1,2,2,2,2
       MOIJ=1,1,1,1,1
              2,3,4,5
        ZIJ=2.0,1.0,1.0,1.0,1.0,
                1.0,1.0,1.0,1.0
The columns (which are LMOs) are allowed to span up to 5 
rows (the nuclei), in situations with multicenter bonds.
MOIJ shows the Mulliken analysis thinks there are four
NH bonds following the nitrogen core.  ZIJ shows that
since each such bond assigns one proton to nitrogen, the
total charge of N is +6.  This is incorrect of course,
as indeed will always happen to some nucleus in a charged
molecule.  In order for the energy analysis to correctly
reproduce the total energy, we must ensure that the 
charge of nitrogen is +7.  The least arbitrary way to
do this is to increase the nitrogen charge assigned to
each NH bond by 1/4.  Since in our case NNUCMO and MOIJ
and much of ZIJ are correct, we need only override a
small part of them with $LOCAL input:
       IJMO(1)=1,2,  1,3,  1,4,  1,5
       ZIJ(1)=1.25, 1.25, 1.25, 1.25
which changes the charge of the first atom of orbitals 
2 through 5 to 5/4, changing ZIJ to
        ZIJ=2.0,1.25,1.25,1.25,1.25,
                1.0, 1.0, 1.0, 1.0
The purpose of the IJMO sparse matrix pointer is to let
you give only the changed parts of ZIJ and/or MOIJ.

    Another way to resolve the problem with NH4+ is to
change one of the 4 equivalent bond pairs into a "proton".  
A "proton" orbital AH treats the LMO as if it were a
lone pair on A, and so assigns +2 to nucleus A.  Use of
a "proton" also generates an imaginary orbital, with
zero electron occupancy.  For example, if we make atom
2 in NH4+ a "proton", by
     IPROT(1)=2
     NNUCMO(2)=1
     IJMO(1)=1,2,2,2   MOIJ(1)=1,0   ZIJ(1)=2.0,0.0
the automatic decomposition of the nuclear charges will be
     NNUCMO=1,1,2,2,2,1
       MOIJ=1,1,1,1,1,2
                3,4,5
        ZIJ=2.0,2.0,1.0,1.0,1.0,1.0
                    1.0,1.0,1.0
The 6th column is just a proton, and the decomposition
will not give any electronic energy associated with
this "orbital", since it is vacant.  Note that the two ways
we have disected the nuclear charges for NH4+ will both
yield the correct total energy, but will give very
different individual orbital components.  Most people
will feel that the first assignment is the least arbitrary,
since it treats all four NH bonds equivalently.
 
    However you assign the nuclear charges, you must
ensure that the sum of all nuclear charges is correct.
This is most easily verified by checking that the energy
sum equals the total SCF energy of your system.  
 
    As another example, H3PO is studied in EXAM26.INP.
Here the MOIDON analysis decides the three equivalent
orbitals on oxygen are O lone pairs, assigning +2 to
the oxygen nucleus for each orbital.  This gives Z(O)=9,
and Z(P)=14.  The least arbitrary way to reduce Z(O)
and increase Z(P) is to recognize that there is some
backbonding in these "lone pairs" to P, and instead 
assign the nuclear charge of these three orbitals by
1/3 to P, 5/3 to O.

    Because you may have to make several runs, looking
carefully at the localized orbital output before the
correct nuclear assignments are made, there is an 
option to skip directly to the decomposition when the
orbital localization has already been done.  Use
   $CONTRL RUNTYP=PROP
   $GUESS  GUESS=MOREAD  NORB=
   $VEC containing the localized orbitals!
   $TWOEI
The latter group contains the necessary Coulomb and
exchange integrals, which are punched by the first
localization, and permits the decomposition to begin
immediately.

    SCF level dipoles can also be analyzed using the
DIPDCM flag in $LOCAL.  The theory of the dipole 
analysis is given in the third paper of the LCD
sequence.  The following list includes application of
the LCD analysis to many problems of chemical interest:

W.England, M.S.Gordon  J.Am.Chem.Soc. 93, 4649-4657 (1971)
W.England, M.S.Gordon  J.Am.Chem.Soc. 94, 4818-4823 (1972)
M.S.Gordon, W.England  J.Am.Chem.Soc. 94, 5168-5178 (1972)
M.S.Gordon, W.England  Chem.Phys.Lett. 15, 59-64 (1972)
M.S.Gordon, W.England  J.Am.Chem.Soc. 95, 1753-1760 (1973)
M.S.Gordon             J.Mol.Struct. 23, 399 (1974)
W.England, M.S.Gordon, K.Ruedenberg,
                       Theoret.Chim.Acta 37, 177-216 (1975)
J.H.Jensen, M.S.Gordon, J.Phys.Chem. 99, 8091-8107(1995)
J.H.Jensen, M.S.Gordon, J.Am.Chem.Soc. 117, 8159-8170(1995)
M.S.Gordon, J.H.Jensen, Acc.Chem.Res. 29, 536-543(1996)

                       * * *

    It is also possible to analyze the MP2 correlation
correction in terms of localized orbitals, for the RHF
case.  The method is that of G.Peterssen and M.L.Al-Laham,
J.Chem.Phys., 94, 6081-6090 (1991).  Any type of localized
orbital may be used, and because the MP2 calculation
typically omits cores, the $LOCAL group will normally 
include only valence orbitals in the localization.  As
mentioned already, the analysis of the MP2 correction must 
be done in a separate run from the SCF analysis, which must
include cores in order to sum up to the total SCF energy.

                       * * *

    Typically, the results are most easily interpreted
by looking at "the bigger picture" than at "the details".
Plots of kinetic and potential energy, normally as a
function of some coordinate such as distance along an
IRC, are the most revealing.  Once you determine, for
example, that the most significant contribution to the
total energy is the kinetic energy, you may wish to look
further into the minutia, such as the kinetic energies
of individual localized orbitals, or groups of LMOs
corresponding to an entire functional group.