$FMOPRP group (optional for FMO runs)
Options setting up SCF convergers, parallelization and
properties are given here.
I. Parameters for SCF convergers and initial guess
MAXIT = the maximum number of monomer SCF iterations.
(default 30)
CONV = monomer SCF energy convergence criterion.
It is considered necessary to set CONV in $SCF to
a value less or equal to the CONV in $FMO.
Usually 1e-7 works well, but for poorly converging
monomer SCF (frequently seen with DFT) one order,
smaller value for CONV in $SCF is recommended,
(1e-7 in $FMO and 1e-8 in $SCF) (default: 1e-7).
NGUESS = controls initial guess (cumulative options, add
all options desired) (default=2):
1 run free monomer SCF
2 if set, dimer density/orbitals are constructed
from the "sum" of monomer quantities, otherwise
Huckel guess will be used for dimers,
and the same applies to trimers.
4 insert HMO projection operator in Huckel guess
8 apply dimer HO projection to dimer initial guess
16 do RHF for each dimer and trimer, then run DFT.
128 do not use orbitals from the previous geometry
during geometry optimization. This is mostly
useful for multilayer optimizations, when this
choice must always be set if basis sets differ.
256 if SCF does not converge, try the alternative
converger (flip between SOSCF and DIIS).
If 2048 is not added, the alternative converger
will start up with the final orbitals of the
unconverged SCF.
512 reorder initial orbitals manually using $GUESS
options (IORDER), applies to MCSCF layers only.
2048 a modifier of 256 - when both 256 and 2048 are
set, the alternative converger will use the
same initial set of initial orbitals as the
unconverged SCF.
IJVEC = Index array enabling reading $VEC inputs defining
initial orbitals for individual n-mer runs.
This consists of quintuplets:
ifg1,jfg1,kfg1,ilay1,norb1, ifg2,jfg2, ...
The first pair indexes $VEC1 with ifg1,jfg1,kfg1
for layer ilay1 expecting norb1 MOs.
The second quintuplet handles $VEC2 etc.
ifg,0,0,ilay1,norb1 is used for monomer IFG,
ifg,jfg,0,ilay1,norb1 is used for dimer IFG,JFG.
ifg equal to 0 in a quintuplet ends the list.
$VEC groups must be used consequently from $VEC1.
(default: all 0s; at most 100 can be given)
MODORB = controls whether orbitals and energies are
exchanged between fragments (additive options).
1 exchange orbitals if set, otherwise densities
2 exchange energies
DFT, ROHF, UHF and MCSCF, SCZV and PIEDA require
MODORB=3. MODORB=3 is in general more robust,
because it provides a better initial guess.
(Default: 0 for RHF, 3 for DFT/ROHF/UHF/MCSCF.)
MCONV = an array specifying SCF convergers for each FMO
step. Individually (MCONV(2) is for monomers,
MCONV(4) for dimers, MCONV(9) for trimers). Each
array element is set to A1+A2+A3, where A1
determines SCF and A2 MCSCF convergers, and A3 is
the direct/conventional bit common for all SCF
methods. MCONV is an additive option:
A1(SCF): A2(MCSCF): A3(direct)
1 EXTRAP 1024 FOCAS 256 FDIFF
2 DAMPH 2048 SOSCF 512 DIRSCF
4 VSHIFT 4096 DROPC
8 RSTRCT 8192 CANONC
16 DIIS 16384 FCORE
32 DEM 32768 FORS
64 SOSCF 65536 n/a
128 LOCOPT 131072 EKT
262144 LINSER
524288 JACOBI
1048576 QUD
There are some limitations on joint usage for each
that can be understood from $SCF or $MCSCF.
If set to -1, the defaults given in $SCF or $MCSCF
are used. See MCONFG. (default: all -1's).
MCONFG = an array specifying SCF convergers for each
fragment during the monomer SCF runs. The value -1
means use the default (defined by MCONV).
The priority in which convergers are chosen is:
MCONFG (highest), if not defined MCONV,
if not defined, $SCF (lowest).
This option is useful in case of poor convergence
caused by charge fluctuations and SCF converger
problems in particular, SOSCF instability for poor
initial guess. Default: all -1.
ESPSCA = scale factors for up to nine initial monomer SCF
iterations. ESPs will be multiplied by these
factors, to soften the effect of environment and
help convergence. At most nine factors can be
defined. (default: all 1.0's)
CNVDMP = damping of SCF convergence, that is, loosen
convergence during the initial monomer SCF
iterations to gain speed. CONV in $SCF and ITOL
and ICUT in $CONTRL are modified.
CONV is set roughly to min(DE/CNVDMP,1e-4), where
DE is the convergence in energy at the given
monomer SCF iteration. It is guaranteed that
CONV,ITOL and ICUT at the end will be set to the
values given in $SCF. Damping is disabled if
CNVDMP is 0. Reasonable values are 10-100.
Care should be taken for restart jobs: since
restart jobs do not know how well FMO converged,
restart jobs start out at the same rough values as
nonrestart jobs, if CNVDMP is used. Therefore for
restart jobs either set CNVDMP appropriately for
the restart (i.e., normally 10-100 times larger
than for the original run) or turn this option
off, otherwise regressive convergence can incur
additional iterations (default: 0).
COROFF = parameter turning off DFT in initial monomer SCF,
similar to SWOFF. COROFF is used during monomer
SCF, and it turns off DFT until monomer energies
converge to this threshold. If COROFF is nonzero,
SWOFF is ignored during monomer SCF, but is used
for dimers and trimer iterations. Setting
COROFF=1e-3 and SWOFF=0 usually produces good DFT
convergence. COROFF may be thought as a macro-
analogue of SWOFF. If monomer SCF converges poorly
(>25 iterations), it is also recommended to raise
CONV in $SCF to 1e-8 (if CONV in $FMO is 1e-7).
Default:1.0E-3 (0.0 skips this option).
NPCMIT = the maximum number of FMO/PCM[m] iterations,
applicable to m>1 only (for m=1, $FMOPRP MAXIT is
used). NPCMIT=2 can be thought as having special
meaning: it is used to define FMO/PCM[l(m)] runs
by forcing the FMO/PCM loop run only twice, which
corresponds to determining PCM charges during the
first iteration (and the m-body level) and then
using them during the second iteration (l-body).
For FMO/PCM[l(m)] only l=1 is implemented and "m"
is given in $PCM IFMO. Default: 30.
CNVPCM = convergence threshold for FMO/PCM[m] iterations,
applicable to m>1 only (for m=1, $FMOPRP CONV is
used). CNVPCM is applied to the total FMO energy
Default: 1.0D-07 Hartree.
PCMOFF = parameter turning PCM off in initial monomer SCF
iterations, analogous to COROFF. PCM is turned
off, until convergence reaches PCMOFF. PCMOFF=0
disables this feature. Default: 0.0
NCVSCF = an array of 2 elements to alter SCF convergers.
After NCVSCF(1) monomer SCF iterations the SCF
converger will switch between SOSCF <-> FULLNR.
This option is useful in converging difficult
cases in the following way:
$SCF diis=.t. soscf=.f. $end
$FMOPRP NCVSCF(1)=2 mconv(4)=65 $end
This results in the initial 2 monomer SCF
iterations being done with DIIS, then a switch to
SOSCF occurs. mconv(4)=65 switches to SOSCF for
dimers.
Note that NCVSCF(1) will only overwrite MCONV, but
not MCONFG. The SCF converger in MCONV(2) will be
enforced after NCVSCF(2) monomer SCF iterations,
overwriting MCONFG as well. This is useful for
the most obnoxiously converging cases. See other
FMO documentation.
Default: 9999,9999 (which means do not use).
NAODIR = a parameter to decide whether to enforce DIRSCF.
Useful for incore integral runs in parallel.
NAODIR is the number of AO orbitals that is
expected to produce 100,000,000 non-zero
integrals. Using this and assuming NAO**3.5
dependence, the program will then guess how many
integrals will each n-mer have and whether they
will fit into the available memory. If they are
determined not to fit, DIRSCF will be set true.
This option overwrites MCONV but not MCONFG.
If set to 0, then the default in-core integral
strategy is used. (default=0)
VDWRAD = array of van der Waals radii in Angstrom, one for
each atom in the periodic table. Reasonable values
are set only for a few light atoms and otherwise a
value of 2.5 is used. VDWRAD values are used only
to compute distance between fragments and thus
somewhat affect all distance-based approximations.
II. Parameters defining parallel execution
MODPAR = parallel options (additive options)
(default: 13, which is 1+4+8)
1 turns on/off heavy job first strategy (reduces
waiting on remaining jobs at barrier points)
(see also 8)
4 broadcast all fragments done by a group at once
rather than fragment by fragment.
8 alters the behavior of fragment initialixation:
if set, fragments are always done in the reverse
order (nfg, nfg-1, ...1) because distance
calculation costs decrease in the same order and
they usually prevail over making Huckel orbitals
or running free monomer SCF. Note that during
SCC (monomer SCF) iterations the order in which
monomers are done is determined by MODPAR=1.
16 if set, hybrid orbital projectors will not be
parallelized (may be useful on slow networks)
32 reserved
64 Broadcast F40 for FMO restarts. F40 should only
be precopied to the grand master scratch
directory and it should NOT exist on all slaves.
256 Replace I/O to fragment density file by
parallel broadcasts from group masters
512 Use DDI memory to store fragment densities
during the monomer step, using supervector
(smallest memory).
1024 Use DDI memory to store fragment densities
during the monomer step, using supervector
(smallest communications).
Only one option out of 512 and 1024 may be used.
NGRFMO = an array that sets the number of GDDI compute
groups during various stages of the calculation.
The first ten elements are used for layer 1, the
next 10 for layer 2, etc.
ngrfmo(1) monomer SCF
ngrfmo(2) dimers
ngrfmo(3) trimers
ngrfmo(4) correlated monomers
ngrfmo(5) separated dimers
ngrfmo(6) SCF monomers in FMO-MCSCF (MCSCF
monomer will be done with ngrfmo(1) groups)
ngrfmo(7) SCF dimers in FMO-MCSCF (MCSCF dimer
be done with ngrfmo(2) groups)
ngrfmo(8-10) reserved
If any of them is zero, the corresponding stage
runs with the previously defined number of groups.
If NGRFMO option is used, it is recommended to set
NGROUP in $GDDI to the total number of nodes.
(default: 0,0,0,0).
MANNOD = manually define node division into groups.
Contrary to MANNOD in $GDDI and here it is defined
for each FMO stage (see NGRFMO) in each layer.
If MANNOD values are set at all, it is required
that they be given corresponding to the first
nonzero NGRFMO value. The MANNOD values should be
given for each nonzero NGRFMO.
E.g. ngrfmo(1)=6,3,0,0,0, 0,0,0,0,0, 4,3
mannod(1)=4,2,2,2,2,2, 5,5,4, 4,4,3,3, 6,6,2
where 6 groups are defined for monomers in layer
1, then 3 for dimers in layer 1, and 4 and 3
groups for monomers and dimers in layer 2.
(default: all -1 which means do not use).
Note that nodes with very large counts may be too large for
good scaling with certain kinds of FMO runs. Any such fat
nodes can be divided into "logical nodes" by using the
kickoff option :cpus= for TCP/IP, or DDI_LOGICAL_NODE_SIZE
for MPI runs. See the DDI instructions.
LOADBF = an array for semi-dynamic load balancing,
specifying the basis set sizes. If it is exceeded,
static load balancing is used. LOADBF has the same
structure as NGRFMO.
(Default: all 0's (disabling this feature).
LOADGR = an array for semi-dynamic load balancing,
specifying the group sizes to be used with LOADBF.
LOADGR has the same structure as NGRFMO.
LOADGR is normally used with MANNOD, because it is
only useful with uneven group sizes. The main
purpose is to handle cases when only a few very
large fragments are mixed with many small ones.
An example: LOADBF(1)=200,400 LOADGR(1)=1,4
NGRFMO(1)=3,5 MANNOD(1)=11,1,1, 3,3,3,3,1
Monomers whose size exceeds 200 basis functions
are executed on 1 group consisting of 11 nodes
with static load balancing. Other monomers are
computed on 2 groups consisting of 1 node with
dynamic load balancing. Dimers whose size exceeds
400 basis functions are executed on 4 groups
consisting of 3 nodes with static load balancing.
Other monomers are computed on 1 group consisting
of 1 node with dynamic load balancing.
Note that after finishing static work load, large
groups will join the dynamic load balancing pool.
(Default: all 0's (disabling this feature).
III. Orbital conversion
File F40 that contains orbital density can be manipulated
in some way to change the information stored in it without
running any FMO calculations. Such conversion requires
irest=2 and the basis sets in the input should define the
old (before conversion) format. The results will be stored
in F30. You should then rename it to F40 and use in a
consequent run (with irest>=2).
Two basic conversion types are supported: A) changing RHF
into MCSCF and B) changing basis sets for RHF. RHF and
MCSCF use different stucture of the restart file (F40) and
therefore conversion is necessary.
For type A the following orbital reordering manipulation
before storing the results can be done, for example
$guess guess=modaf norder=1 iorder(28)=34,28
Type B is typically used for preparing good initial
orbitals for hard to converge cases. E.g., you can use
something like 6-21G to converge the orbitals and then
convert F40 to be used with 6-311G*. At present there is a
limitation that only density based (MODORB=0) files may be
converged, i.e. you cannot do it for DFT and MCSCF.
MAXAOC = The new (i.e., after conversion) maximum number of
AOs per fragment. If you don't know what it should
be you can run a CHECK job with the new basis set
and find the number in "Max AOs per frg:".
If this number is equal to the old value, then
type A is chosen.
IBFCON = the array giving pairs of the old and new numbers
of AOs for each atom in $DATA (type B only).
MAPCON = maps determining how to copy old orbitals into new
(type B only). See the example.
Example: $DATA contains only H and O (in this order), F40
was computed with 6-31G and you want to convert to 6-31G**.
One water per fragment.
MAXAOC=25 25=5*2+15=new basis size for 6-31G**
IBFCON(1)=2,5, 9,15
2 and 5 for H (6-31 and 6-31G**), 9 and 15 for O
MAPCON(1)=1,2,0,0,0,
1,2,3,4,5,6,7,8,9,0,0,0,0,0,0
Here we copy the two s functions of each H, and add p
polarization p to each H (3 0's), and similarly we copy
nine s,p functions for O, and add d polarization (6 0's)
In order to construct MAPCON, you should know in what order
Gaussian primitives are stored. The easiest way to learn
this is to run a simple calculation and check the output
(SHELL information).
IV. Printing, properties, restart, and dimensions.
NPRINT = controls print-out (bit additive)
bits 1-2
0 normal output
1 reduced output (recommended for single points)
2 small output (recommended for MD)
3 minimal output (for large scale production MD)
4 print interfragment distances. Note: any of
RESPAP, RESPPC, or RESDIM must be non-zero or
otherwise nothing will be printed. If you only
want the distances but no approximations, set
the thresholds to huge values, e.g. resdim=1000.
8 print Mulliken charges.
Note: RESPPC must be set (non-zero), see above.
16 special test run to check for errors in $FMOBND.
32 increase the print-out/punch-out level a) on the
last SCC iteration for monomers and b) for all
dimers/trimers. It can be employed to plot
fragment MOs while using $CONTRL NPRINT=-5 to
reduce monomer output.
64 print atomic coordinates for each fragment.
128 skip printing ES dimer energies during their
calculation (but a summary will list them).
Note that choosing 128+2 will reduce the
NFG^2 memory requirement for DFTB in half.
PRTDST = array of four print-out thresholds:
1. print all pairs of fragments separated by less
than PRTDST(1).
2. print a warning if two fragments are closer
than PRTDST(2), intended mostly to monitor
suspicious geometries during optimization.
3. print a warning if two fragments are closer
than PRTDST(3) and have no detached bond between
them, intended to check input.
PRTDST(3) values should slightly exceed the
longest detached bond in the system.
4. PRTDST(4) has a completely different meaning.
In the summary of pair interactions, only those
values will be printed, which are larger than
PRTDST(4), in the units of kcal/mol.
Using zero for PRTDST(1) and PRTDST(2) turns them
off. Similarly, use PRTDST(3)=-1 to turn it off.
PRTDST has no units, as it applies to unitless FMO
distances (e.g., 0.5 means half the sum of van der
Waals radii for the closest pair of atoms).
(default: 0.0, 0.5, 0.6, 0.0)
IREST = restart level (to use it, you must copy .F40.000
to the scratch disk of either the first node only
(if MODPAR=64) or all nodes (otherwise).
Multilayer FMO can be restarted if either all
layers use the same basis set or if you save
the .F40.000 file for layer 1 and use it for the
same layer.
0 no restart
2 restart monomer SCF (SCC).
4 restart dimers. Requires monomer energies be
given in $FMOENM. Some or no dimer energies
may also be given in $FMOEND, in which case
those dimers with energies will not be run.
Usually the only property that can be obtained
with IREST=4 is the energy. The only exception
is: a) 1024 was added to IREST when monomer SCF
was run and b) property restart files (*.F38*)
from each node were saved and copied to the
scratch directory for the IREST=1028 job. If
these two conditions are met, gradient and ES
moments can be restarted with IREST=1028.
1024 write property restart files during monomer SCF
and/or use them to restart gradient and/or ES
moments. No other property may be restarted.
Note that monomer restarts (irest=2) do not need adding
1024, as the properties are recomputed. 1024 should only
be used for IREST=4 (or for IREST=0 to save restart data).
MODPRP = some extra FMO properties (bit additive)
Default: 0.
1 total electron density (AO-basis matrix, written
to F10: useful to create initial orbitals for ab
initio).
2 reserved.
4 electron density on a grid.
8 if set, the grid output is a "sparse cube
file", otherwise a "Gaussian cube file".
16 automatically generate grid for modprp = 4 or 8.
32 molecular electrostatic potential on a grid.
128 spin density on a grid
Only one bit out of 4 and 8 may be set.
Default: 0.
NGRID = three integers, giving the number of 3D grid
points for monomers with NOPFRG=4 in x,y and z
directions (default 0,0,0).
GRDPAD = Grid padding. Contributions to density on grid
will be restricted to the box surrounding an n-mer
with each atom represented by a sphere of GRDPAD
vdW radii. In general the finer effects one is
interested in, the larger GRDPAD should be. For
example, if one plots not density, but density
differences and a very small cutoff is used, then
a larger value of GRDPAD (2.5 or 3.0) may be
preferred.
Default: 2.0.
IMECT = The partitioning method for the interfragment
charge transfer (computed from Mulliken charges).
IMECT pertains only to those dimers between which
a bond is detached.
IMECT=0,1,2,3,4 are supported (see source code).
(default: 4)
V. Interaction analysis (PIEDA)
IPIEDA = 0 skip the analysis (default)
1 perform brief PL-state analysis (FMO pair
interactions)
2 perform full PL-state analysis with the PL0-
state data.
PCM may be used with IPIEDA=0 or 1.
N0BDA = gives the number of detached bonds. This
parameter should be set to a nonzero value only in
runs that produce BDA pair energies. (default: 0)
R0BDA = array of the detached bond lengths, whose number
is N0BDA. R0BDA must be given if E0BDA is used.
E0BDA = the array of BDA pair energies, whose number is
N0BDA*4.
EFMO0 = the array of the free state fragment energies,
first NFRAG correlated, then NFRAG uncorrelated
values.
EPL0DS = monomer polarization energies, first NFRAG values
of PL0d, then NFRAG values of PL0s, then NFRAG
values of PL0DI.
EINT0 = the total components for the PL0 state:
ES0, EX0, CT+mix0, DI0.
None of the PIEDA input values (except IPIEDA) are to be
manually prepared, all should come from the punch file of
preceeding calculations.
The brief order of IPIEDA=2 execution is:
1. run FMO0.
2. compute BDA energies (if detached bonds are present),
using sample files in tools/fmo/pieda. To do this, one
needs only R0BDA for a given system. R0BDA is punched by
any FMO run at the very beginning, so NBODY=0 type of run
might be used to generate it.
3. The results of (1) are EFMO0; the results of (2) are
E0BDA; use them to run PL0, whose results will be EPL0DS
and EINT0.
4. Run PL with the results of (1),(2) and (3).
The alternative is to run IPIEDA=1, which requires none of
the above data, but it will use E0BDA is available.
==========================================================
==========================================================
529 lines are written.
Edited by Shiro KOSEKI on Mon Feb 13 10:50:16 2017.