$DFT group                 (relevant if DFTTYP is chosen)
                         (relevant if SCFTYP=RHF,UHF,ROHF)
 
    Note that if DFTTYP=NONE, an ab initio calculation will
be performed, rather than density functional theory.
 
     This group permits the use of various one electron
(usually empirical) operators instead of the true many
electron Hamiltonian.  Two programs are provided, METHOD=
GRID or GRIDFREE.  The programs have different functionals
available, and so the keyword DFTTYP (which is entered in
$CONTRL) and other associated inputs are documented
separately below.  Every functional that has the same name
in both lists is an identical functional, but each METHOD
has a few functionals that are missing in the other.
 
    The grid free implementation is based on the use of the
resolution of the identity to simplify integrals so that
they may be analytically evaluated, without using grid
quadratures.  The grid free DFT computations in their
present form have various numerical errors, primarily in
the gradient vectors.  Please do not use the grid-free DFT
program without reading the discussion in the 'Further
References' section regarding the gradient accuracy.
 
    The grid based DFT uses a typical grid quadrature to
compute integrals over the rather complicated functionals,
using two possible angular grid types.
 
    Achieving a self-consistent field with DFT is rather
more difficult than for normal HF, so DIIS is the default
converger.
 
    Both DFT programs will run in parallel. See the two
lists below for possible functionals in the two programs.
 
    See also the $TDDFT input group for excited states.
 
METHOD = selects grid based DFT or grid free DFT.
       = GRID     Grid based DFT (default)
       = GRIDFREE Grid free DFT
 
 
 
DFTTYP is given in $CONTRL, not here in $DFT!  Possible
values for the grid-based program are listed first,
 
          ----- options for METHOD=GRID -----
 
DFTTYP = NONE     means ab initio computation (default)
 
Many choices are given below, perhaps the most sensible are
               local DFT: SVWN
            pure DFT GGA: BLYP, PW91, B97-D, PBE/PBEsol
          hybrid DFT GGA: B3LYP, X3LYP, PBE0
       pure DFT meta-GGA: revTPSS
     hybrid DFT meta-GGA: TPSSh, M06
but of course, everyone has their own favorite!
 
 
              pure exchange functionals:
       = SLATER   Slater exchange
       = BECKE    Becke 1988 exchange
       = GILL     Gill 1996 exchange
       = OPTX     Handy-Cohen exchange
       = PW91X    Perdew-Wang 1991 exchange
       = PBEX     Perdew-Burke-Ernzerhof exchange
These will be used with no correlation functional at all.
 
 
              pure correlation functionals:
       = VWN      Vosko-Wilk-Nusair correlation, using
                  their electron gas formula 5 (aka VWN5)
       = VWN3     Vosko-Wilk-Nusair correlation, using
                  their electron gas formula 3
       = VWN1RPA  Vosko-Wilke-Nusair correlation, using
                  their e- gas formula 1, with RPA params.
       = PZ81     Perdew-Zener 1981 correlation
       = P86      Perdew 1986 correlation
       = LYP      Lee-Yang-Parr correlation
       = PW91C    Perdew-Wang 1991 correlation
       = PBEC     Perdew-Burke-Ernzerhof correlation
       = OP       One-parameter Progressive correlation
These will be used with 100% HF exchange, if chosen.
 
 
              combinations (partial list):
       = SVWN     SLATER exchange + VWN5 correlation
                  Called LDA/LSDA in physics for RHF/UHF.
       = SVWN1RPA Slater exchange + VWN1RPA correlation
       = BLYP     BECKE exchange + LYP correlation
       = BOP      BECKE exchange + OP correlation
       = BP86     BECKE exchange + P86 correlation
       = GVWN     GILL exchange + VWN5 correlation
       = GPW91    GILL exchange + PW91 correlation
       = PBEVWN   PBE exchange + VWN5 correlation
       = PBEOP    PBE exchange + OP correlation
       = OLYP     OPTX exchange + LYP correlation
       = PW91     means PW91 exchange + PW91 correlation
       = PBE      means PBE exchange + PBE correlation
There's a nearly infinite set of pairings (well, 6*9), so
we show only enough to give you the idea.  In other words,
pairs are formed by abbreviating the exchange functionals
    SLATER=S, BECKE=B, GILL=G, OPTX=O, PW91X=PW91, PBEX=PBE
and matching them with any correlation functional, of which
only two are abbreviated when used in combinations,
    PW91C==>PW91, PBEC==>PBE
The pairings shown above only scratch the surface, but
clearly, many possibilities, such as PW91PBE, are nonsense!
 
 
             pure DFT GGA functionals:
       = EDF1     empirical density functional #1, which is
                  a modified BLYP from Adamson/Gill/Pople.
       = PW91     Perdew/Wang 1991
       = PBE      Perdew/Burke/Ernzerhof 1996
       = revPBE   PBE as revised by Zhang/Yang
       = RPBE     PBE as revised by Hammer/Hansen/Norskov
       = PBEsol   PBE as revised by Perdew et al for solids
       = HCTH     Hamprecht/Cohen/Tozer/Handy's 1998 mod
                  to B97, omitting HF exchange (m=4)
       = HCTH93   Hamprecht/Cohen/Tozer/Handy's 1998 mod
                  to B97, omitting HF exchange, fitting to
                  93 atoms and molecules
       = HCTH120  later fit to 120 systems
       = HCTH147  later fit to 147 systems
       = HCTH407  later fit to 407 systems (best)
       = HCTH407P also known as HCTH/407+, designed by
                  Boese/Chandra/Martin/Marx for
                  describing hydrogen bonds
       = HCTH14   also known as HCTH p=1/4, made by
                  Menconi/Wilson/Tozer
       = HCTH76   also known as HCTH p=7/6, made by
                  Menconi/Wilson/Tozer
       = SOGGA    PBE revised by Zhao/Truhlar for solids
       = MOHLYP   metal optimized OPTX, half LYP
       = B97-D    Grimme's modified B97, with dispersion
                  correction (this forces DC=.TRUE.)
       = SOGGA11  optimized with broad applicability for
                  chemistry, by Peverati/Zhao/Truhlar
 
 
             hybrid GGA functionals:
       = BHHLYP   HF and BECKE exchange + LYP correlation
       = B3PW91   Becke's 3 parameter exchange hybrid,
                  with PW91 correlation functional
       = B3LYP    this is a hybrid method combining five
                  functionals: Becke + Slater + HF exchange
                  (B3), with LYP + VWN5 correlation.
                  B3LYPV5 is a synonym for B3LYP.
       = B3LYPV1R use VWN1RPA in place of VWN5, matches the
                  e- gas formula chosen by some programs.
       = B3LYPV3  use VWN3 in place of B3LYP's VWN5
       = B3P86    B3-type exchange, P86 correlation, using
                  VWN3 as the LDA part of the correlation.
                  B3P86V3 is a synonym for B3P86.
       = B3P86V1R use VWN1RPA in place of VWN3
       = B3P86V5  use VWN5 in place of VWN3
       = B97      Becke's 1997 hybrid functional
       = B97-1    Hamprecht/Cohen/Tozer/Handy's 1998
                  reparameterization of B97
       = B97-2    Wilson/Bradley/Tozer's 2001 mod to B97
       = B97-3    Keal/Tozer's 2005 mod to B97
       = B97-K    Boese/Martin's 2004 mod for kinetics
       = B98      Schmider/Becke's 1998 mode to B97,
                  using their best "2c" parameters.
       = PBE0     a hybrid made from PBE
       = X3LYP    HF+Slater+Becke88+PW91 exchange,
                  and LYP+VWN1RPA correlation.
       = SOGGA11X a hybrid based on SOGGA11,
                  with 40.15% HF exchange.
       = APF      mixed functional based on PBE0 and B3PW91
Each includes some Hartree-Fock exchange, and also may use
a linear combination of many DFT parts.
 
 
           range separated functionals:
These are also known as "long-range corrected functionals".
LC-BVWN, LC-BOP, LC-BLYP, or LC-BPBE are available by
selecting BVWN, BOP, BLYP, or BPBE and also setting the
flag LC=.TRUE. (see LC and also MU below).  Others are
selected by their specific name, without using LC:
       = CAMB3LYP coulomb attenuated B3LYP
       = wB97     omega separated form of B97
       = wB97X    wB97 with short-range HF exchange
       = wB97X-D  dispersion corrected wB97X
       = CAMQTP00 Verma/Bartlett's reparametrization of
                  CAMB3LYP for vertical ionization
                  potential
       = CAMQTP01 Jin/Bartlett's reparametrization of
                  CAMB3LYP for vertical ionization
                  potential
M11 is also range-separated, but is listed below with the
other meta-GGAs.
 
 
              "double hybrid" GGA:
       = B2PLYP   mixes BLYP, HF exchange, and MP2!
                  See related inputs CHF and CMP2 below.
          "double hybrid" and "range separated":
       = wB97X-2  intended for use with GBASIS=CCT,CCQ,CC5
       = wB97X-2L intended for use with GBASIS=N311
                  NGAUSS=6 NDFUNC=3 NFFUNC=1 NPFUNC=3
                  DIFFSP=.T. DIFFS=.T.
Note: there are no analytic gradients for "double hybrids".
Note: the B2PLYP family uses the conventional MP2 energy
and may be used for closed shell or spin-unrestricted open
shell cases.  The wB97X-2 family uses the SCS-MP2 energy,
and thus is limited to closed shell cases at present.
 
 
              meta-GGA functionals:
These are not hybridized with HF exchange, unless that is
explicitly stated below.
       = VS98     Voorhis/Scuseria, 1998
       = PKZB     Perdew/Kurth/Zupan/Blaha, 1999
       = tHCTH    Boese/Handy's 2002 metaGGA akin to HCTH
       = tHCTHhyb tHCTH's hybrid with 15% HF exchange
       = BMK      Boese/Martin's 2004 parameterization of
                  tHCTHhyb for kinetics
       = TPSS     Tao/Perdew/Staroverov/Scuseria, 2003
       = TPSSh    TPSS hybrid with 10% HF exchange
       = TPSSm    TPSS with modified parameter, 2007
       = revTPSS  revised TPSS, 2009
       = dlDF     a reparameterized M05-2X, reproducing
                  interaction energies which have had all
                  dispersion removed.  This MUST be used
                  with a special -D correction to recover
                  dispersion.  See 'Further References'.
       = M05      Minnesota exchange-correlation, 2005
                  a hybrid with 28% HF exchange.
       = M05-2X   M05, with doubled HF exchange, to 56%
       = M06      Minnesota exchange-correlation, 2006
                  a hybrid with 27% HF exchange.
       = M06-L    M06, with 0% HF exchange (L=local)
       = M06-2X   M06, with doubled HF exchange, to 54%
       = M06-HF   M06 correlation, using 100% HF exchange
       = M08-HX   M08 with 'high HF exchange'
       = M08-SO   M08 with parameters that enforce the
                  correct second order gradient expansion.
       = M11      M11 range-separated hybrid
       = M11-L    M11 local (0% HF exchange) with
                  dual-range exchange
       = MN12-L   uses a nonseparable functional form aiming to provide
                  balanced peformance for both chemistry and solid-state
                  physics applications
       = MN12-SX  screened-exchanged (SX) hybrid functional with 20% HF
                  exchange for short-range and 0% HF exchange for
                  long-range
       = MN15     uses a nonseparable functional form with 44% HF
                  exchange. This fuctional is a global hybride (no
                  range-separation).
       = MN15-L   a local functional with 0% HF exchange.
       = REVM06   revised hybrid M06
       = REVM06-L revised local M06-L
 
When the M06 family was created, Truhlar recommended M06
for the general situation, but see his "concluding remarks"
in the M06 reference about which functional is best for
what kind of test data set.  The most recent M11 family is
probably a better choice, and two functionals fit all the
needs of the older M05/M06/M08 families.
 
https://en.wikipedia.org/wiki/Minnesota_functionals
 
An extensive bibliography for all functionals can be
found in the 'Further References' section of this manual.
 
Note that only a subset of these functionals can be used
for TD-DFT energy or gradients.  These subsets are listed
in the $TDDFT input group.
 
* * * dispersion corrections * * *
 
Many exchange-correlation functionals fail to compute
intra- and inter-molecular dispersion interactions
accurately.  Two possible correction schemes are provided
below.  The first uses empirically chosen C6 and C8
coefficients, while the latter obtains these from the
molecular DFT densities.  At most, only one of the LRDFLG
or DC options below may be chosen.
 
DC     = a flag to turn on Grimme's empirical dispersion
         correction, involving scaled R**(-6) terms.
         N.B.  This empiricism may also be added to plain
         Hartree-Fock, by choosing DFTTYP=NONE with DC=.T.
         Three different versions exist, see IDCVER.
         (default=.FALSE., except if DFTTYP=B97-D, wB97X-D)
 
IDCVER = 1 means 1st 2004 implementation.
       = 2 means 2nd 2006 implementation DFT-D2,
           default for B97-D, wB97X-D.
       = 3 means 3rd 2010 implementation DFT-D3.
           Default if DC is chosen and IDCVER isn't given.
       = 4 means modified 3rd implementation DFT-D3(BJ).
           (-4 is used for DFT-D3(BJ) for HF-3c).
         Setting IDCVER will force DC=.TRUE.
 
 
GCP     = a flag for the geometric counterpoise scheme
          correction in HF-3c.
 
SRB     = a flag for short-range basis set incompleteness
          (SRB) correction in HF-3c.
 
DCCHG  = a flag to use Chai-Head-Gordon damping function
         instead of Grimme's 2006 function. Pertinent only
         for the DFT-D2 method.  Forces DC=.TRUE.
         (default=.FALSE. except for wB97X-D)
 
DCABC  = a flag to turn on the computation of the E(3) non-
         additive energy term. Pertinent only for DFT-D3,
         it forces DC=.TRUE.  (default=.FALSE.)
 
    The following parameters govern Grimme's semiempirical
dispersion term. They are basis set and functional
dependent, so they exist for only a few DFTTYP. Default
values are automatically selected and printed out in the
output file for many common density functionals.
    The following keywords are for entering non-standard
values. For DFT-D2 values, see also:
      R.Peverati and K.K.Baldridge
      J.Chem.Theory Comput. 4, 2030-2048 (2008).
For DFT-D3 values, and a detailed explanation of each
parameter, see:
      S. Grimme, J. Antony, S. Ehrlich and H. Krieg,
      J.Chem.Phys. 132, 154104/1-19(2010)
and for DFT-D3(BJ):
      S. Grimme, S. Ehrlich and L. Goerigk,
      J.Comput.Chem. 32, 1456-1465 (2011)
 
DCALP  = alpha parameter in the DFT-D damping function
         (same as alpha6 in Grimme's DFT-D3 notation).
         Note also that alpha8 and alpha10 in DFT-D3 have
         constrained values of:
         alpha8 = alpha6 + 2, alpha10 = alpha8 + 2.
         Default=14.0 for DFT-D3
                =20.0 for DFT-D2
                =23.0 for DFT-D1
                =6.00 for DCCHG=.TRUE.
 
DCSR   = sR exponential parameter to scale the van der
         Waals radii (same as sR,6 in Grimme's DFT-D3
         notation). Note also that sR,8 in DFT-D3 have a
         fixed value of 1.0.
         Optimized values are automatically selected for
         some of the more common functionals, otherwise,
         the default is 1.00 for DFT-D3, 1.10 for DFT-D2,
         and 1.22 for DFT-D1.
 
DCS6   = s6 linear parameter for scaling the C6 term.
         Optimized values are automatically selected for
         some of the more common functionals, otherwise,
         the default is 1.00.
 
DCS8   = s8 linear parameter for scaling the C8 term of
         DFT-D3. Pertinent only for DFT-D3.
         Optimized values are automatically selected for
         some of the more common functionals, otherwise,
         the default is 1.00.
 
DCA1   = a1 parameter appearing in the -D3(BJ) dispersion
         model. Optimized values are automatically
         selected for a set of known functionals,
         otherwise the default is 0.50.
 
DCA2   = a2 parameter appearing in the -D3(BJ) dispersion
         model. Optimized values are automatically
         selected for a set of known functionals,
         otherwise the default is 4.00.
 
The old keywords DCPAR and DCEXP were replaced by DCS6 and
DCSR in 2010.  Similarly, DCOLD has morphed into IDCVER.
 
                          - - -
 
The Local Response Dispersion (LRD) correction includes
atomic pair-wise -C6/R**6, -C8/R**8, and -C10/R**10 terms,
whose coefficients are computed from the molecular system's
electron density and its nuclear gradient.  The nuclear
gradient assumes the dispersion coefficients do not vary
with geometry, which causes only a very small error in the
gradient.  Optionally, 3 and 4 center terms may be added,
at the 1/R**6 level; in this case, nuclear gradients may
not be computed at all.
 
Since the three numerical parameters are presently known
only for the long-range exchange corrected BOP functional,
calculations may specify simply DFTTYP=LCBOPLRD.  The
"LCBOPLRD" functional will automatically select the
following:
     DFTTYP=BOP   LC=.TRUE. MU=0.47
     LRDFLG=.TRUE. LAMBDA=0.232 KAPPA=0.600 RZERO=3.22
leaving only the choice for MLTINT up to you.
 
References for LRD are
   T.Sato, H.Nakai J.Chem.Phys. 131, 224104/1-12(2009)
   T.Sato, H.Nakai J.Chem.Phys. 133, 194101/1-9(2010)
 
LRDFLG = flag choosing the Local Response Dispersion (LRD)
         C6, C8, and C10 corrections.  Default=.FALSE.
 
MLTINT = flag to add the 3 and 4 center 6th order terms,
         the default=.FALSE.  Note that nuclear gradients
         are not available if these multi-center terms
         are requested.
 
Three numerical parameters may be input.  The defaults
shown are optimized for the BOP functional with the LC
correction for long-range exchange.
 
LAMBDA = parameter adjusting the density gradient
         correction for the atomic and atomic pair
         polarizabilities.  (default=0.232)
KAPPA  = parameter in the damping function (default=0.600)
RZERO  = parameter in the damping function (default=3.22)
 
It may be interesting to see a breakdown of the total
dispersion correction, using these keywords:
 
PRPOL  = print out atomic effective polarizabilities
         (default=.FALSE.)
PRCOEF = N  (default N=0)
         print out dispersion coefficient to N-th order.
PRPAIR = print out atomic pair dispersion energies
         (default=.FALSE.)
 
 
              * * * range separation * * *
 
LC     = flag to turn on the long range correction (LC),
         which smoothly replaces the DFT exchange by the
         HF exchange at long inter-electron distances.
         (default=.FALSE.)
         This option can be used only with the Becke
         exchange functional (Becke) and a few correlation
         functionals: DFTTYP=BVWN, BOP, BLYP, BPBE only.
         For example, B3LYP has a fixed admixture of HF
         exchange, so it cannot work with the LC option.
         See H.Iikura, T.Tsuneda, T.Yanai, and K.Hirao,
         J.Chem.Phys. 115, 3540 (2001).
 
MU     = A parameter for the long range correction scheme.
         Increasing MU increases the HF exchange used,
         very small MU produces the DFT limit.
         (default=0.33)
 
Other range-separated options exist, invoked by naming the
functional, such as DFTTYP=CAMB3LYP (see the DFTTYP keyword
for a full list).
 
 
    * * * B2x-PLYP double hybrid functionals * * *
 
B2xPLYP Double Hybrid functionals have the general formula:
    Exc = (1-cHF) * ExGGA + cHF * ExHF
       + (1-cMP2) * EcGGA + cMP2 * E(2)
 
The next keywords allow the choice of cHF and cMP2. Both
values must be between 0 and 1 (0-100%).
 
CHF    = amount of HF exchange. (default=0.53)
 
CMP2   = amount of MP2. (default=0.27)
 
Some other common double hybrid functionals are available
simply by choosing DFTTYP=B2PLYP, and changing the CHF and
CMP2 parameters. Popular parametrizations are:
                             CHF       CMP2
     ------------------------------------------
     B2-PLYP (default)   |   0.53  |   0.27   |
     ------------------------------------------
     B2K-PLYP            |   0.72  |   0.42   |
     ------------------------------------------
     B2T-PLYP            |   0.60  |   0.31   |
     ------------------------------------------
     B2GP-PLYP           |   0.65  |   0.36   |
     ------------------------------------------
 
 
                 * * * Grid Input * * *
 
Only one of the three grid types may be chosen for the run.
The default (if no selection is made) is the Lebedev grid.
In order to duplicate results obtained prior to April 2008,
select the polar coordinate grid NRAD=96 NTHE=12 NPHI=24.
Energies can be compared if and only if the identical grid
type and density is used, analogous to needing to compare
with the identical basis set expansions.  See REFS.DOC for
more information on grids.  See similar inputs in $TDDFT.
 
RADTYP = type of radial quadrature
       = MHL  : Murray, Handy and Laming radial grid (default)
       = MK3  : Mura and Knowles Log3 radial grid
       = TA   : Treutler and Ahlrichs radial grid
 
Lebedev grid:
 
NRAD   = number of radial points in the Euler-MacLaurin
         quadrature. (default=96)
 
NLEB   = number of angular points in the Lebedev grids.
         (default=302).  Possible values are 86, 110, 146,
         170, 194, 302, 350, 434, 590, 770, 974, 1202,
         1454, 1730, 2030...
         The same switch can be used to select spherical grids
         developed by A.S.Popov:
         246, 264, 342, 432 : rotaional octahedral symmetry grids
         350 and 398 : full octahedral symmetry grids (like Lebedev grids)
         212 : rotational icosahedral grid
         rotational icosahedral grids with 192 and 242 points
         are currently disabled in the source code
 
 
Meta-GGA functionals require a tighter grid to achieve the
same accuracy. For this reason a tighter default grid of
NRAD=99 and NLEB=590 is chosen by default with all meta-GGA
functionals.
 
The default for NLEB means that nuclear gradients will be
accurate to about the default OPTTOL=0.00010 (see $STATPT),
590 approaches OPTTOL=0.00001, and 1202 is "army grade".
 
The next two specify radial/angular in a single keyword:
 
SG1    = a flag to select the "standard grid 1", which has
         24 radial points, and various pruned Lebedev
         grids, from 194 down to 6.  (default=.FALSE.
         This grid is very fast, but produces gradients
         whose accuracy reaches only OPTTOL=0.00050.
         This grid should be VERY USEFUL for the early
         steps of a geometry optimization.
 
JANS   = two unpublished grids due to Curtis Janssen,
         implemented here differently than in MPQC:
       = 1 uses 95 radial points for all atoms, and prunes
           from a Lebedev grid whose largest size is 434,
           thus using about 15,000 grid points/atom.
       = 2 uses 155 radial points for all atoms, and prunes
           from a Lebedev grid whose largest size is 974,
           thus using about 71,000 grid points/atom.
           This is a very accurate grid, e.g. "army grade".
       The information for pruning exists only for H-Ar,
       so heavier elements will use the large radial/
       Lebedev grid without any pruning.
 
polar coordinate grid:
 
NRAD   = number of radial points in the Euler-MacLaurin
         quadrature. (96 is reasonable)
 
NTHE   = number of angle theta grids in Gauss-Legendre
         quadrature (polar coordinates). (12 is reasonable)
 
NPHI   = number of angle phi grids in Gauss-Legendre
         quadrature.  NPHI should be double NTHE so points
         are spherically distributed. (24 is reasonable)
 
The number of angular points will be NTHE*NPHI.  The values
shown give a gradient accuracy near the default OPTTOL of
0.00010, while NTHE=24 NPHI=48 approaches OPTTOL=0.00001,
and "army grade" is NTHE=36 NPHI=72.
 
 
              * * * Grid Switching * * *
 
At the first geometry of the run, pure HF iterations will
be performed, since convergence of DFT is greatly improved
by starting with the HF density matrix.  After DFT engages,
most runs (at all geometries, except for PCM or numerical
Hessians) will use a coarser grid during the early DFT
iterations, before reaching some initial convergence.
After that, the full grid will be used.  Together, these
switchings can save considerable CPU time.
 
SWOFF =  turn off DFT, to perform pure SCF iterations,
         until the density matrix convergence falls below
         this threshold.  This option is independent of
         SWITCH and can be used with or without it. It is
         reasonable to pick SWOFF > SWITCH > CONV in $SCF.
         SWOFF pertains only to the first geometry that the
         run computes, and is automatically disabled if you
         choose GUESS=MOREAD to provide initial orbitals.
         The default is 5.0E-3.
 
SWITCH = when the change in the density matrix between
         iterations falls below this threshhold, switch
         to the desired full grid (default=3.0E-4)
         This keyword is ignored if the SG1 grid is used.
         Setting to zero disables DFT grid switching.
 
NRAD0  = same as NRAD, but defines initial coarse grid.
         default = smaller of 24 and NRAD/4
 
NLEB0  = same as NLEB, but defines initial coarse grid.
         default = 110
 
NTHE0  = same as NTHE, but defines initial coarse grid.
         default = smaller of 8, NTHE/3
 
NPHI0  = same as NPHI, but defines initial coarse grid.
         default = smaller of 16, NPHI/3
 
 
molecular grid construction parameters:
 
BFCTYP = specify algorithm of molecular grid construction
         from atomic grids (Becke's fuzzy cell algorithm)
       = BECKE  : use Becke's fuzzy cell method from the original
                  paper.
       = SSF    : use algorithm by Stratmann, Scuseria and Frisch.
                  It usually screens out more points than Becke's
                  algorithm
       = LEGACY : use legacy DFT code (default)
       BECKE and SSF switches allow multithreaded DFT run, but
       does not support DC, LRD and analytical Hessian calculation
       types, as well as polar spherical grid with non-C1 symmetry.
 
PARTFN = select grid partitioning function
         relevant only if BFCTYP=SSF or BECKE
       = BECKE   : use polynomial proposed in Becke's original work.
                   Default when BFCTYP=BECKE
       = SSF     : use polynomial proposed by Stratmann, Scuseria
                   and Frisch. Default when BFCTYP=SSF
       = ERF     : use erf-based partitioning function based on
                   unpublished ERF1 algorithm in NWChem
       = SMSTEP2 : use smoothstep-2 polynomial (very effective,
                   but aggressive partitioning)
       = SMSTEP3 : use smoothstep-3 polynomial
       = SMSTEP4 : use smoothstep-4 polynomial
       = SMSTEP5 : use smoothstep-5 polynomial
 
 
technical parameters:
 
THRESH = threshold for ignoring small contributions to the
         Fock matrix.  The default is designed to produce
         no significant energy loss, even when the grid is
         as good as "army grade".  If for some reason you
         want to turn all threshhold tests off, of course
         requiring more CPU, enter 1.0e-15.
         default: 1.0e-4/Natoms/NRAD/NTHE/NPHI
 
GTHRE  = threshold applied to gradients, similar to THRESH.
         < 1 assign this value to all thresholds
         = 1 use the default thresholds (default).
         > 1 divide default thresholds by this value.
         If you wish to increase accuracy, set GTHRE=10.
         The default introduces an error of roughly 1e-7
         (a.u./bohr) in the gradient.
 
WTDER  = switch on/off grid weigths derivative contribution to
         the nuclear gradient.
       = .TRUE.  : compute grid weights derivatives
       = .FALSE. : do not compute grid weights derivatives
         Default is .TRUE. if BFCTYP=LEGACY, and .FALSE.
         otherwise.
         Note: BFCTYP=BECKE and BFCTYP=SSF currently do not
         support grid weight derivatives.
         Grid weights are only important for low-quality grids.
         The default grid settings are good enough to disable
         grid weight derivatives calculation. This switch can
         be turned off to improve convergence of geometry
         optimization, if the optimized structure is not the
         one with lowest energy.
 
 
 
The keyword $DFTTYP is given in $CONTRL, and may have these
values if the grid-free program is chosen:
 
        ----- options for METHOD=GRIDFREE -----
 
DFTTYP = NONE     means ab initio computation (default)
                     exchange functionals:
       = XALPHA   X-Alpha exchange (alpha=0.7)
       = SLATER   Slater exchange (alpha=2/3)
       = BECKE    Becke's 1988 exchange
       = DEPRISTO Depristo/Kress exchange
       = CAMA     Handy et al's mods to Becke exchange
       = HALF     50-50 mix of Becke and HF exchange
                     correlation functionals:
       = VWN      Vosko/Wilke/Nusair correlation, formula 5
       = PWLOC    Perdew/Wang local correlation
       = LYP      Lee/Yang/Parr correlation
                     exchange/correlation functionals:
       = BVWN     Becke exchange + VWN5 correlation
       = BLYP     Becke exchange + LYP correlation
       = BPWLOC   Becke exchange + Perdew/Wang correlation
       = B3LYP    hybrid HF/Becke/LYP using VWN formula 5
       = CAMB     CAMA exchange + Cambridge correlation
       = XVWN     Xalpha exchange + VWN5 correlation
       = XPWLOC   Xalpha exchange + Perdew/Wang correlation
       = SVWN     Slater exchange + VWN5 correlation
       = SPWLOC   Slater exchange + PWLOC correlation
       = WIGNER   Wigner exchange + correlation
       = WS       Wigner scaled exchange + correlation
       = WIGEXP   Wigner exponential exchange + correlation
 
AUXFUN = AUX0  uses no auxiliary basis set for resolution
               of the identity, limiting accuracy.
       = AUX3  uses the 3rd generation of RI basis sets,
               These are available for the elements H to
               Ar, but have been carefully considered for
               H-Ne only.  (DEFAULT)
 
THREE  = a flag to use a resolution of the identity to
         turn four center overlap integrals into three
         center integrals.  This can be used only if
         no auxiliary basis is employed. (default=.FALSE.)
==========================================================
 
 
 
==========================================================
715 lines are written.
Edited by Shiro KOSEKI on Thu Mar 5 10:25:38 2020.