mirror of
https://github.com/nwchemgit/nwchem.git
synced 2026-07-21 14:35:21 -04:00
97f303e6f8
automatically be expanded to include useful information about the checkin (including the file's revision number). With the switch over to SVN this was lost because SVN only does this expansion if you explicitly ask for it (for every single file). I have added a script to the contrib directory that sets the appropriate property to get SVN to do this expansion. This script will make it easy to do this every time new source files are added. It is called svn_expand_Id, the script contains some comments that explain the issue and how it addresses this. This checkin sets this property for a subset of the relevant files (trying to commit all files at once failed with svn crashing). In future the script will only affect those files for which the property was not set before.
1289 lines
51 KiB
TeX
1289 lines
51 KiB
TeX
%
|
|
% $Id$
|
|
%
|
|
\label{sec:scf}
|
|
|
|
The NWChem self-consistent field (SCF) module computes closed-shell
|
|
restricted Hartree-Fock (RHF) wavefunctions, restricted high-spin
|
|
open-shell Hartree-Fock (ROHF) wavefunctions, and spin-unrestricted
|
|
Hartree-Fock (UHF) wavefunctions.
|
|
|
|
The \verb+SCF+ directive provides input to the SCF module and is a
|
|
compound directive that encloses additional directives specific to the
|
|
SCF module:
|
|
\begin{verbatim}
|
|
SCF
|
|
...
|
|
END
|
|
\end{verbatim}
|
|
|
|
\section{Wavefunction type}
|
|
|
|
A spin-restricted, closed shell RHF calculation is performed by
|
|
default. An error results if the number of electrons is inconsistent
|
|
with this assumption. The number of electrons is inferred from the
|
|
total charge on the system and the sum of the effective nuclear
|
|
charges of all centers (atoms and dummy atoms, Section
|
|
\ref{sec:geom}). The total charge on the system is zero by default,
|
|
unless specified at some value by input on the \verb+CHARGE+ directive
|
|
(Section \ref{sec:toplevel}).
|
|
|
|
The options available to define the SCF wavefunction and multiplicity
|
|
are as follows:
|
|
|
|
\begin{verbatim}
|
|
SINGLET
|
|
DOUBLET
|
|
TRIPLET
|
|
QUARTET
|
|
QUINTET
|
|
SEXTET
|
|
SEPTET
|
|
OCTET
|
|
NOPEN <integer nopen default 0>
|
|
RHF
|
|
ROHF
|
|
UHF
|
|
\end{verbatim}
|
|
|
|
The optional keywords \verb+SINGLET+, \verb+DOUBLET+, \ldots,
|
|
\verb+OCTET+ and \verb+NOPEN+ allow the user to specify the number of
|
|
singly occupied orbitals for a particular calculation. \verb+SINGLET+
|
|
is the default, and specifies a closed shell; \verb+DOUBLET+ specifies
|
|
one singly occupied orbital; \verb+TRIPLET+ specifies two singly
|
|
occupied orbitals; and so forth. If there are more than seven singly
|
|
occupied orbitals, the keyword \verb+NOPEN+ must be used, with the
|
|
integer \verb+nopen+ defining the number of singly occupied
|
|
orbitals (sometimes referred to as open shells).
|
|
|
|
If the multiplicity is any value other than \verb+SINGLET+, the
|
|
default calculation will be a spin-restricted, high-spin, open-shell
|
|
SCF calculation (keyword ROHF). The open-shell orbitals must be the
|
|
highest occupied orbitals. If necessary, any starting vectors may be
|
|
rearranged through the use of the \verb+SWAP+ keyword on the
|
|
\verb+VECTORS+ directive (see Section \ref{sec:vectors}) to accomplish
|
|
this.
|
|
|
|
A spin-unrestricted solution can also be performed by specifying the
|
|
keyword \verb+UHF+. In UHF calculations, it is assumed that the
|
|
number of singly occupied orbitals corresponds to the difference
|
|
between the number of alpha-spin and beta-spin orbitals. For example,
|
|
a UHF calculation with 2 more alpha-spin orbitals than beta-spin
|
|
orbitals can be obtained by specifying
|
|
|
|
\begin{verbatim}
|
|
scf
|
|
triplet ; uhf # (Note: two logical lines of input)
|
|
...
|
|
end
|
|
\end{verbatim}
|
|
|
|
The user should be aware that, by default, molecular orbitals are
|
|
symmetry adapted in NWChem. This may not be desirable for fully
|
|
unrestricted wavefunctions. In such cases, the user has the option of
|
|
defeating the defaults by specifying the keywords \verb+ADAPT OFF+
|
|
(see Section \ref{sec:adapt}) and \verb+SYM OFF+ (see Section
|
|
\ref{sec:sym}).
|
|
|
|
The keywords \verb+RHF+ and \verb+ROHF+ are provided in the code for
|
|
completeness. It may be necessary to specify these in order to modify
|
|
the behavior of a previous calculation (see Section \ref{sec:persist}
|
|
for restart behavior).
|
|
|
|
\section{{\tt SYM} --- use of symmetry}
|
|
\label{sec:sym}
|
|
|
|
\begin{verbatim}
|
|
SYM <string (ON||OFF) default ON>
|
|
\end{verbatim}
|
|
|
|
This directive enables/disables the use of symmetry to speed up Fock matrix
|
|
construction (via the petite-list or skeleton algorithm) in the SCF, if
|
|
symmetry was used in the specification of the geometry. Symmetry
|
|
adaptation of the molecular orbitals is not affected by this option.
|
|
The default is to use symmetry if it is specified in the geometry
|
|
directive (Section \ref{sec:geom}).
|
|
|
|
For example, to disable use of symmetry in Fock matrix construction:
|
|
\begin{verbatim}
|
|
sym off
|
|
\end{verbatim}
|
|
|
|
\section{{\tt ADAPT} -- symmetry adaptation of MOs}
|
|
\label{sec:adapt}
|
|
|
|
\begin{verbatim}
|
|
ADAPT <string (ON||OFF) default ON>
|
|
\end{verbatim}
|
|
|
|
The default in the SCF module calculation is to force symmetry
|
|
adaption of the molecular orbitals. This does not affect the speed of
|
|
the calculation, but without explicit adaption the resulting orbitals
|
|
may be symmetry contaminated for some problems. This is especially
|
|
likely if the calculation is started using orbitals from a distorted
|
|
geometry.
|
|
|
|
The underlying assumption in the use of symmetry in Fock matrix
|
|
construction is that the density is totally symmetric. If the orbitals
|
|
are symmetry contaminated, this assumption may not be valid --- which
|
|
could result in incorrect energies and poor convergence of the
|
|
calculation. It is thus advisable when specifying \verb+ADAPT OFF+ to
|
|
also specify \verb+SYM OFF+ (Section \ref{sec:sym}).
|
|
|
|
\section{{\tt TOL2E} --- integral screening threshold}
|
|
\label{sec:tol2e}
|
|
|
|
\begin{verbatim}
|
|
TOL2E <real tol2e default min(10e-7 , 0.01*$thresh$)>
|
|
\end{verbatim}
|
|
|
|
The variable \verb+tol2e+ is used in determining the integral
|
|
screening threshold for the evaluation of the energy and related
|
|
Fock-like matrices. The Schwarz inequality is used to screen the
|
|
product of integrals and density matrices in a manner that results in
|
|
an accuracy in the energy and Fock matrices that approximates the
|
|
value specified for \verb+tol2e+.
|
|
|
|
It is generally not necessary to set this parameter directly. Specify
|
|
instead the required precision in the wavefunction, using the
|
|
\verb+THRESH+ directive (Section \ref{sec:thresh}). The default
|
|
threshold is the minimum of $10^{-7}$ and 0.01 times the requested
|
|
convergence threshold for the SCF calculation (Section
|
|
\ref{sec:thresh}).
|
|
|
|
The input to specify the threshold explicitly within the \verb+SCF+
|
|
directive is, for example:
|
|
|
|
\begin{verbatim}
|
|
tol2e 1e-9
|
|
\end{verbatim}
|
|
|
|
For very diffuse basis sets, or for high-accuracy calculations it
|
|
might be necessary to set this parameter. A value of $10^{-12}$ is
|
|
sufficient for nearly all such purposes.
|
|
|
|
\section{{\tt VECTORS} --- input/output of MO vectors}
|
|
\label{sec:vectors}
|
|
|
|
|
|
\begin{verbatim}
|
|
VECTORS [[input] (<string input_movecs default atomic>) || \
|
|
(project <string basisname> <string filename>) || \
|
|
(fragment <string file1> [<string file2> ...])] \
|
|
[swap [alpha||beta] <integer vec1 vec2> ...] \
|
|
[reorder <integer atom1 atom2> ...] \
|
|
[output <string output_filename default input_movecs>] \
|
|
[lock]
|
|
[rotate <string input_geometry> <string input_movecs>]
|
|
\end{verbatim}
|
|
|
|
The \verb+VECTORS+ directive allows the user to specify the source and
|
|
destination of the molecular orbital vectors. In a startup
|
|
calculation (see Section \ref{sec:start}), the default source for
|
|
guess vectors is a diagonalized Fock matrix constructed from a
|
|
superposition of the atomic density matrices for the particular
|
|
problem. This is usually a very good guess. For a restarted
|
|
calculation, the default is to use the previous MO vectors.
|
|
|
|
The optional keyword \verb+INPUT+ allows the user to specify the
|
|
source of the input molecular orbital vectors as any of the following:
|
|
\begin{itemize}
|
|
\item \verb+ATOMIC+ --- eigenvectors of a Fock-like matrix formed from
|
|
a superposition of the atomic densities (the default guess). See
|
|
Sections \ref{sec:atomscf} and \ref{sec:tolguess}.
|
|
\item \verb+HCORE+ --- eigenvectors of the bare-nucleus Hamiltonian or
|
|
the one-electron Hamiltonian.
|
|
\item \verb+filename+ --- the name of a file containing the MO vectors
|
|
from a previous calculation. Note that unless the path is fully
|
|
qualified, or begins with a dot (``.''), then it is assumed to
|
|
reside in the directory for permanent files (see Section
|
|
\ref{sec:dirs}).
|
|
\item \verb+PROJECT basisname filename+ --- projects the existing MO
|
|
vectors in the file \verb+filename+ from the smaller basis with name
|
|
\verb+basisname+ into the current basis. The definition of the
|
|
basis \verb+basisname+ must be available in the current database,
|
|
and the basis must be smaller than the current basis. In addition,
|
|
the geometry used for the previous calculations must have the atoms
|
|
in the same order and in the same orientation as the current
|
|
geometry.
|
|
\item \verb+FRAGMENT file1 ...+ --- assembles starting MO vectors from
|
|
previously performed calculations on fragments of the system and is
|
|
described in more detail in Section \ref{sec:fragguess}. Even
|
|
though there are some significant restrictions in the use of the
|
|
initial implementation of this method (see Section
|
|
\ref{sec:fragguess}), this is the most powerful initial guess option
|
|
within the code. It is particularly indispensable for open shell
|
|
metallic systems.
|
|
\item \verb+ROTATE input_geometry input_movecs+ --- rotates
|
|
MO vectors generated at a previous geometry
|
|
to the current active geometry.
|
|
\end{itemize}
|
|
|
|
The molecular orbitals are saved every iteration if more than 600
|
|
seconds have elapsed, and also at the end of the calculation. At
|
|
completion (converged or not), the SCF module always canonically
|
|
transforms the molecular orbitals by {\em separately} diagonalizing
|
|
the closed--closed, open--open, and virtual--virtual blocks of the
|
|
Fock matrix.
|
|
|
|
The name of the file used to store the MO vectors is determined as
|
|
follows:
|
|
\begin{itemize}
|
|
\item if the \verb+OUTPUT+ keyword was specified on the \verb+VECTORS+
|
|
directive, then the filename that follows this keyword is used, or
|
|
\item if the input vectors were read from a file, this file is reused
|
|
for the output vectors (overwriting the input vectors); else,
|
|
\item a default file name is generated in the directory for permanent
|
|
files (Section \ref{sec:dirs}) by prepending \verb+".movecs"+ with
|
|
the file prefix, i.e., \verb+"<file_prefix>.movecs"+.
|
|
\end{itemize}
|
|
The name of this file is stored in the database so that a subsequent
|
|
SCF calculation will automatically restart from these MO vectors.
|
|
|
|
Applications of this directive are illustrated in the following
|
|
examples.
|
|
|
|
Example 1:
|
|
\begin{verbatim}
|
|
vectors output h2o.movecs
|
|
\end{verbatim}
|
|
Assuming a start-up calculation, this directive will result in use of
|
|
the default atomic density guess, and will output the vectors to the
|
|
file \verb+h2o.movecs+.
|
|
|
|
Example 2:
|
|
\begin{verbatim}
|
|
vectors input initial.movecs output final.movecs
|
|
\end{verbatim}
|
|
This directive will result in the initial vectors being read from the
|
|
file \verb+"initial.movecs"+. The results will be written to the file
|
|
\verb+final.movecs+. The contents of \verb+"initial.movecs"+ will not
|
|
be changed.
|
|
|
|
Example 3:
|
|
\begin{verbatim}
|
|
vectors input project "small basis" small.movecs
|
|
\end{verbatim}
|
|
This directive will cause the calculation to start from vectors in the
|
|
file \verb+"small.movecs"+ which are in a basis named \verb+"small basis"+.
|
|
The output vectors will be written to the default file
|
|
\verb+"<file_prefix.movecs>"+.
|
|
|
|
Once starting vectors have been obtained using any of the possible
|
|
options, they may be reordered through use of the \verb+SWAP+ keyword.
|
|
This optional keyword requires a list of orbital pairs that will be
|
|
swapped. For UHF calculations, separate \verb+SWAP+ keywords may be
|
|
provided for the alpha and beta orbitals, as necessary.
|
|
|
|
An example of use of the \verb+SWAP+ directive:
|
|
\begin{verbatim}
|
|
vectors input try1.movecs swap 173 175 174 176 output try2.movecs
|
|
\end{verbatim}
|
|
This directive will cause the initial orbitals to be read from the
|
|
file \verb+"try1.movecs"+. The vectors for the orbitals within the
|
|
pairs 173--175 will be swapped with those within 174--176, so the
|
|
resulting order is 175, 176, 173, 174. The final orbitals obtained in
|
|
the calculation will be written to the file \verb+"try2.movecs"+.
|
|
|
|
The swapping of orbitals occurs as a sequential process in the order
|
|
(left to right) input by the user. Thus, regarding each pair as an
|
|
elementary transposition it is possible to construct arbitrary
|
|
permutations of the orbitals. For instance, to apply the permutation
|
|
$(6 7 8 9)$\footnote{The cyclic permutation $(6 7 8 9)$ maps the
|
|
ordered list {\tt 6 7 8 9} into {\tt 9 6 7 8}.} we note that this
|
|
permutation is equal to $(6 7)(7 8)(8 9)$, and thus may be specified
|
|
as
|
|
\begin{verbatim}
|
|
vectors swap 8 9 7 8 6 7
|
|
\end{verbatim}
|
|
|
|
Another example, now illustrating this feature for a UHF calculation,
|
|
is the directive
|
|
\begin{verbatim}
|
|
vectors swap beta 4 5 swap alpha 5 6
|
|
\end{verbatim}
|
|
This input will result in the swapping of the 5--6 alpha orbital pair
|
|
and the 4--5 beta orbital pair. (All other items in the input use the
|
|
default values.)
|
|
|
|
The \verb+LOCK+ keyword allows the user to specify that the ordering
|
|
of orbitals will be locked to that of the initial vectors, insofar as
|
|
possible. The default is to order by ascending orbital energies within
|
|
each orbital space. One application where locking might be desirable
|
|
is a calculation where it is necessary to preserve the ordering of a
|
|
previous geometry, despite flipping of the orbital energies. For such
|
|
a case, the \verb+LOCK+ directive can be used to prevent the SCF
|
|
calculation from changing the ordering, even if the orbital energies
|
|
change.
|
|
|
|
The mapping of the MO's to the nuclei can be changed using the \verb+REORDER+ keyword.
|
|
Once starting vectors have been obtained using any of the possible
|
|
options, the \verb+REORDER+ keyword moves the MO coefficients between atoms
|
|
listed in the integer list. This keyword is particularly useful for calculating localized
|
|
electron and hole states.
|
|
|
|
This optional keyword requires a list containing the new atom ordering.
|
|
It is not necessary to provide separate lists for alpha and
|
|
beta orbitals.
|
|
|
|
An example of use of the \verb+REORDER+ keyword:
|
|
\begin{verbatim}
|
|
vectors input try1.movecs reorder 2 1 output try2.movecs
|
|
\end{verbatim}
|
|
This directive will cause the initial orbitals to be read from the
|
|
file \verb+"try1.movecs"+. The MO coefficients for the basis functions
|
|
on atom 2
|
|
will be swapped with those on atom 1.
|
|
The final orbitals obtained in
|
|
the calculation will be written to the file \verb+"try2.movecs"+.
|
|
|
|
The following example shows how the \verb+ROTATE+ keyword can be used to rotate
|
|
MO vectors calculated at geometry \verb+geom1+ to geometry \verb+geom2+, which has a
|
|
different rotational orientation:
|
|
|
|
\begin{verbatim}
|
|
set geometry geom1
|
|
dft
|
|
vectors input atomic output geom1.mo
|
|
end
|
|
task dft
|
|
|
|
set geometry geom2
|
|
dft
|
|
vectors input rotate geom1 geom1.mo output geom2.mo
|
|
end
|
|
task dft
|
|
|
|
\end{verbatim}
|
|
\subsection{Superposition of fragment molecular orbitals}
|
|
\label{sec:fragguess}
|
|
|
|
The fragment initial guess is particularly useful in the following
|
|
instances:
|
|
\begin{itemize}
|
|
\item The system naturally decomposes into molecules that can be
|
|
treated individually, e.g., a cluster.
|
|
\item One or more fragments are particularly hard to converge and
|
|
therefore much time can be saved by converging them independently.
|
|
\item A fragment (e.g., a metal atom) must be prepared with a specific
|
|
occupation. This can often be readily accomplished with a
|
|
calculation on the fragment using dummy charges to model a ligand
|
|
field.
|
|
\item The molecular occupation predicted by the atomic initial guess
|
|
is often wrong for systems with heavy metals which may have
|
|
partially occupied orbitals with lower energy than some doubly
|
|
occupied orbitals. The fragment initial guess avoids this problem.
|
|
\end{itemize}
|
|
|
|
\begin{verbatim}
|
|
VECTORS [input] fragment <string file1> [<string file2> ...]
|
|
\end{verbatim}
|
|
The molecular orbitals are formed by superimposing the previously
|
|
generated orbitals of fragments of the molecule being studied. These
|
|
fragment molecular orbitals must be in the same basis as the current
|
|
calculation. The input specifies the files containing the fragment
|
|
molecular orbitals. For instance, in a calculation on the water
|
|
dimer, one might specify
|
|
\begin{verbatim}
|
|
vectors fragment h2o1.movecs h2o2.movecs
|
|
\end{verbatim}
|
|
where \verb+h2o1.movecs+ contains the orbitals for the first fragment, and
|
|
\verb+h2o2.movecs+ contains the orbitals for the second fragment.
|
|
|
|
A complete example of the input for a calculation on the water
|
|
dimer using the fragment guess is as follows:
|
|
\begin{verbatim}
|
|
start dimer
|
|
|
|
title "Water dimer SCF using fragment initial guess"
|
|
|
|
geometry dimer
|
|
O -0.595 1.165 -0.048
|
|
H 0.110 1.812 -0.170
|
|
H -1.452 1.598 -0.154
|
|
O 0.724 -1.284 0.034
|
|
H 0.175 -2.013 0.348
|
|
H 0.177 -0.480 0.010
|
|
end
|
|
|
|
geometry h2o1
|
|
O -0.595 1.165 -0.048
|
|
H 0.110 1.812 -0.170
|
|
H -1.452 1.598 -0.154
|
|
end
|
|
|
|
geometry h2o2
|
|
O 0.724 -1.284 0.034
|
|
H 0.175 -2.013 0.348
|
|
H 0.177 -0.480 0.010
|
|
end
|
|
|
|
basis
|
|
o library 3-21g
|
|
h library 3-21g
|
|
end
|
|
|
|
set geometry h2o1
|
|
scf; vectors input atomic output h2o1.movecs; end
|
|
task scf
|
|
|
|
set geometry h2o2
|
|
scf; vectors input atomic output h2o2.movecs; end
|
|
task scf
|
|
|
|
set geometry dimer
|
|
scf
|
|
vectors input fragment h2o1.movecs h2o2.movecs \
|
|
output dimer.movecs
|
|
end
|
|
task scf
|
|
\end{verbatim}
|
|
First, the geometry of the dimer and the two monomers are specified
|
|
and given names. Then, after the basis specification, calculations
|
|
are performed on the fragments by setting the geometry to the
|
|
appropriate fragment (Section \ref{sec:set}) and redirecting the
|
|
output molecular orbitals to an appropriately named file. Note also
|
|
that use of the atomic initial guess is forced, since the default
|
|
initial guess is to use any existing MOs which would not be
|
|
appropriate for the second fragment calculation. Finally, the dimer
|
|
calculation is performed by specifying the dimer geometry, indicating
|
|
use of the fragment guess, and redirecting the output MOs.
|
|
|
|
The following points are important in using the fragment initial guess:
|
|
\begin{enumerate}
|
|
\item The fragment calculations must be in the same basis set as the
|
|
full calculation.
|
|
\item The order of atoms in the fragments and the order in which the
|
|
fragment files are specified must be such that when the fragment
|
|
basis sets are concatenated all the basis functions are in the same
|
|
order as in the full system. This is readily accomplished by first
|
|
generating the full geometry with atoms for each fragment
|
|
contiguous, splitting this into numbered fragments and specifying
|
|
the fragment MO files in the correct order on the \verb+VECTORS+
|
|
directive.
|
|
\item The occupation of orbitals is preserved when they are merged
|
|
from the fragments to the full molecule and the resulting occupation
|
|
must match the requested occupation for the full molecule. E.g., a
|
|
triplet ROHF calculation must be comprised of fragments that have
|
|
a total of exactly two open-shell orbitals.
|
|
\item Because of these restrictions, it is not possible to introduce
|
|
additional atoms (or basis functions) into fragments for the purpose
|
|
of cleanly breaking real bonds. However, it is possible, and highly
|
|
recommended, to introduce additional point charges to simulate the
|
|
presence of other fragments.
|
|
\item MO vectors of partially occupied or strongly polarized systems
|
|
are very sensitive to orientation. While it is possible to specify
|
|
the same fragment MO vector file multiple times in the
|
|
\verb+VECTORS+ directive, it is usually much better to do a separate
|
|
calculation for each fragment.
|
|
\item Linear dependencies which were present in a fragment calculation
|
|
may be magnified in the full calculation. When this occurs,
|
|
some of the fragment's highest virtual orbitals will not be copied to the
|
|
full system, and a warning will be printed.
|
|
|
|
\end{enumerate}
|
|
|
|
A more involved example is now presented. We wish to model the sextet
|
|
state of Fe(III) complexed with water, imidazole and a heme with a net
|
|
unit positive charge. The default atomic guess does not give the
|
|
correct $d^5$ occupation for the metal and also gives an incorrect
|
|
state for the double anion of the heme. The following performs
|
|
calculations on all of the fragments. Things to note are:
|
|
\begin{enumerate}
|
|
\item The use
|
|
of a dummy $+2$ charge in the initial guess on the heme which in part
|
|
simulates the presence of the metal ion, and also automatically forces
|
|
an additional two electrons to be added to the system (the default net
|
|
charge being zero).
|
|
\item The iron fragment calculation (charge +3, $d^5$, sextet) will
|
|
yield the correct open-shell occupation for the full system. If,
|
|
instead, the {\it d}-orbitals were partially occupied (e.g., the doublet
|
|
state) it would be useful to introduce dummy charges around the iron
|
|
to model the ligand field and thereby lift the degeneracy to obtain
|
|
the correct occupation.
|
|
\item $C_s$ symmetry is used for all of the calculations. It is not
|
|
necessary that the same symmetry be used in all of the
|
|
calculations, provided that the order and orientation of the atoms
|
|
is preserved.
|
|
\item The \verb+unset scf:*+ directive is used immediately before
|
|
the calculation on the full system so that the default name for the
|
|
output MO vector file can be used, rather than having to specify it
|
|
explicitly.
|
|
\end{enumerate}
|
|
\begin{verbatim}
|
|
start heme6a1
|
|
title "heme-H2O (6A1) from M.Dupuis"
|
|
|
|
############################################################
|
|
# Define the geometry of the full system and the fragments #
|
|
############################################################
|
|
|
|
geometry full-system
|
|
symmetry cs
|
|
|
|
H 0.438 -0.002 4.549
|
|
C 0.443 -0.001 3.457
|
|
C 0.451 -1.251 2.828
|
|
C 0.452 1.250 2.828
|
|
H 0.455 2.652 4.586
|
|
H 0.461 -2.649 4.586
|
|
N1 0.455 -1.461 1.441
|
|
N1 0.458 1.458 1.443
|
|
C 0.460 2.530 3.505
|
|
C 0.462 -2.530 3.506
|
|
C 0.478 2.844 1.249
|
|
C 0.478 3.510 2.534
|
|
C 0.478 -2.848 1.248
|
|
C 0.480 -3.513 2.536
|
|
C 0.484 3.480 0.000
|
|
C 0.485 -3.484 0.000
|
|
H 0.489 4.590 2.664
|
|
H 0.496 -4.592 2.669
|
|
|
|
H 0.498 4.573 0.000
|
|
H 0.503 -4.577 0.000
|
|
H -4.925 1.235 0.000
|
|
H -4.729 -1.338 0.000
|
|
C -3.987 0.685 0.000
|
|
N -3.930 -0.703 0.000
|
|
C -2.678 1.111 0.000
|
|
C -2.622 -1.076 0.000
|
|
H -2.284 2.126 0.000
|
|
H -2.277 -2.108 0.000
|
|
N -1.838 0.007 0.000
|
|
|
|
Fe 0.307 0.000 0.000
|
|
|
|
O 2.673 -0.009 0.000
|
|
H 3.238 -0.804 0.000
|
|
H 3.254 0.777 0.000
|
|
end
|
|
|
|
geometry ring-only
|
|
symmetry cs
|
|
H 0.438 -0.002 4.549
|
|
C 0.443 -0.001 3.457
|
|
C 0.451 -1.251 2.828
|
|
C 0.452 1.250 2.828
|
|
H 0.455 2.652 4.586
|
|
H 0.461 -2.649 4.586
|
|
N1 0.455 -1.461 1.441
|
|
N1 0.458 1.458 1.443
|
|
C 0.460 2.530 3.505
|
|
C 0.462 -2.530 3.506
|
|
C 0.478 2.844 1.249
|
|
C 0.478 3.510 2.534
|
|
C 0.478 -2.848 1.248
|
|
C 0.480 -3.513 2.536
|
|
C 0.484 3.480 0.000
|
|
C 0.485 -3.484 0.000
|
|
H 0.489 4.590 2.664
|
|
H 0.496 -4.592 2.669
|
|
|
|
Bq 0.307 0.0 0.0 charge 2 # simulate the iron
|
|
end
|
|
|
|
geometry imid-only
|
|
symmetry cs
|
|
H 0.498 4.573 0.000
|
|
H 0.503 -4.577 0.000
|
|
H -4.925 1.235 0.000
|
|
H -4.729 -1.338 0.000
|
|
C -3.987 0.685 0.000
|
|
N -3.930 -0.703 0.000
|
|
C -2.678 1.111 0.000
|
|
C -2.622 -1.076 0.000
|
|
H -2.284 2.126 0.000
|
|
H -2.277 -2.108 0.000
|
|
N -1.838 0.007 0.000
|
|
end
|
|
|
|
geometry fe-only
|
|
symmetry cs
|
|
Fe .307 0.000 0.000
|
|
end
|
|
|
|
geometry water-only
|
|
symmetry cs
|
|
O 2.673 -0.009 0.000
|
|
H 3.238 -0.804 0.000
|
|
H 3.254 0.777 0.000
|
|
end
|
|
|
|
############################
|
|
# Basis set for everything #
|
|
############################
|
|
|
|
basis nosegment
|
|
O library 6-31g*
|
|
N library 6-31g*
|
|
C library 6-31g*
|
|
H library 6-31g*
|
|
Fe library "Ahlrichs pVDZ"
|
|
end
|
|
|
|
##########################################################
|
|
# SCF on the fragments for initial guess for full system #
|
|
##########################################################
|
|
|
|
scf; thresh 1e-2; end
|
|
|
|
set geometry ring-only
|
|
scf; vectors atomic swap 80 81 output ring.mo; end
|
|
task scf
|
|
|
|
set geometry water-only
|
|
scf; vectors atomic output water.mo; end
|
|
task scf
|
|
|
|
set geometry imid-only
|
|
scf; vectors atomic output imid.mo; end
|
|
task scf
|
|
|
|
charge 3
|
|
set geometry fe-only
|
|
scf; sextet; vectors atomic output fe.mo; end
|
|
task scf
|
|
|
|
##########################
|
|
# SCF on the full system #
|
|
##########################
|
|
|
|
unset scf:* # This restores the defaults
|
|
|
|
charge 1
|
|
|
|
set geometry full-system
|
|
|
|
scf
|
|
sextet
|
|
vectors fragment ring.mo imid.mo fe.mo water.mo
|
|
maxiter 50
|
|
end
|
|
|
|
task scf
|
|
\end{verbatim}
|
|
|
|
\subsection{Atomic guess orbitals with charged atoms}
|
|
\label{sec:atomscf}
|
|
|
|
As noted above, the default guess vectors are based on superimposing
|
|
the density matrices of the neutral atoms. If some atoms are
|
|
significantly charged, this default guess may be improved upon by
|
|
modifying the atomic densities. This is done by setting parameters
|
|
that add fractional charges to the occupation of the valence atomic
|
|
orbitals. Since the atomic SCF program does not have its own input
|
|
block, the \verb+SET+ directive (Section \ref{sec:set}) must be used
|
|
to set these parameters.
|
|
|
|
The input specifies a list of tags (i.e., names of atoms in a
|
|
geometry, see Section \ref{sec:geom}) and the charges to be added to
|
|
those centers. Two parameters must be set as follows:
|
|
\begin{verbatim}
|
|
set atomscf:tags_z <string list_of_tags>
|
|
set atomscf:z <real list_of_charges>
|
|
\end{verbatim}
|
|
|
|
\sloppy
|
|
|
|
The array of strings \verb+atomscf:tags_z+ should be set to the list
|
|
of tags, and the array \verb+atomscf:z+ should be set to the list of
|
|
charges which must be real numbers (not integers). All atoms that
|
|
have a tag specified in the list of tags will be assigned the
|
|
corresponding charge from the list of charges.
|
|
|
|
\fussy
|
|
|
|
For example, the following specifies that all oxygen atoms with tag
|
|
\verb+O+ be assigned a charge of \verb+-1+ and all iron atoms with tag
|
|
\verb+Fe+ be assigned a charge of \verb=+2=
|
|
\begin{verbatim}
|
|
set atomscf:z -1 2.0
|
|
set atomscf:tags_z O Fe
|
|
\end{verbatim}
|
|
|
|
There are some limitations to this feature. It is not possible to add
|
|
electrons to closed shell atoms, nor is it possible to remove all
|
|
electrons from a given atom. Attempts to do so will cause the code to
|
|
report an error, and it will not report further errors in the input
|
|
for modifying the charge even when they are detected.
|
|
|
|
Finally, recall that the database is persistent (Section
|
|
\ref{sec:persist}) and that the modified settings will be used in
|
|
subsequent atomic guess calculations unless the data is deleted from
|
|
the database with the \verb+UNSET+ directive (Section
|
|
\ref{sec:unset}).
|
|
|
|
\section{Accuracy of initial guess}
|
|
\label{sec:tolguess}
|
|
|
|
For SCF, the initial Fock-matrix construction from the atomic guess is
|
|
now (staring from version 3.3) performed to a default precision of
|
|
1e-7. However, other wavefunctions, notably DFT, use a lower
|
|
precision. In charged, or diffuse basis sets, this precision may not
|
|
be sufficient and could result in incorrect ordering of the initial
|
|
orbitals. The accuracy may be increased with the following directive
|
|
which should be inserted in the top-level of input (i.e., outside of
|
|
the SCF input block) and before the {\tt TASK} directive.
|
|
\begin{verbatim}
|
|
set tolguess 1e-7
|
|
\end{verbatim}
|
|
|
|
\section{{\tt THRESH} --- convergence threshold}
|
|
\label{sec:thresh}
|
|
|
|
\begin{verbatim}
|
|
THRESH <real thresh default 1.0e-4>
|
|
\end{verbatim}
|
|
|
|
This directive specifies the convergence threshold for the
|
|
calculation. The convergence threshold is the norm of the orbital
|
|
gradient, and has a default value in the code of $10^{-4}$.
|
|
|
|
The norm of the orbital gradient corresponds roughly to the precision
|
|
available in the wavefunction, and the energy should be converged to
|
|
approximately the square of this number. It should be noted, however,
|
|
that the precision in the energy will not exceed that of the integral
|
|
screening tolerance. This tolerance (Section \ref{sec:tol2e}) is
|
|
automatically set from the convergence threshold, so that sufficient
|
|
precision is usually available by default.
|
|
|
|
The default convergence threshold suffices for most SCF energy and
|
|
geometry optimization calculations, providing about 6--8 decimal
|
|
places in the energy, and about four significant figures in the
|
|
density and energy derivative with respect to nuclear coordinates.
|
|
However, greater precision may be required for calculations involving
|
|
weakly interacting systems, floppy molecules, finite-difference of
|
|
gradients to compute the Hessian, and for post-Hartree-Fock
|
|
calculations. A threshold of $10^{-6}$ is adequate for most such
|
|
purposes, and a threshold of $10^{-8}$ might be necessary for very
|
|
high accuracy or very weak interactions. A threshold of $10^{-10}$
|
|
should be regarded as the best that can be attained in most
|
|
circumstances.
|
|
|
|
\section{{\tt MAXITER} --- iteration limit}
|
|
\label{sec:max}
|
|
|
|
\begin{verbatim}
|
|
MAXITER <integer maxiter default 8>
|
|
\end{verbatim}
|
|
|
|
\sloppy
|
|
|
|
The maximum number of iterations for the SCF calculation defaults to
|
|
20 for both ROHF/RHF and UHF calculations. For most molecules, this
|
|
number of iterations is more than sufficient for the quadratically
|
|
convergent SCF algorithm to obtain a solution converged to the default
|
|
threshold (see Section \ref{sec:thresh} above). If the SCF program
|
|
detects that the quad\-ratically con\-ver\-gent algorithm is not
|
|
efficient, then it will resort to a lin\-early con\-ver\-gent
|
|
algorithm and increase the maximum number of iterations by 10.
|
|
|
|
\fussy
|
|
|
|
Convergence may not be reached in the maximum number of iterations for
|
|
many reasons, including input error (e.g., an incorrect geometry or a
|
|
linearly dependent basis), a very low convergence threshold, a poor
|
|
initial guess, or the fact that the system is intrinsically hard to
|
|
converge due to the presence of many states with similar energies.
|
|
|
|
The following sets the maximum number of SCF iterations to 50:
|
|
\begin{verbatim}
|
|
maxiter 50
|
|
\end{verbatim}
|
|
|
|
\section{{\tt PROFILE} --- performance profile}
|
|
|
|
This directive allows the user to obtain timing and parallel
|
|
execution information about the SCF module. It is specified by the
|
|
simple keyword
|
|
|
|
\begin{verbatim}
|
|
PROFILE
|
|
\end{verbatim}
|
|
|
|
This option can be helpful in understanding the computational
|
|
performance of an SCF calculation. However,
|
|
it can introduce a significant overhead
|
|
on machines that have expensive timing routines, such as the SUN.
|
|
|
|
\section{{\tt DIIS} --- DIIS convergence}
|
|
|
|
This directive allows the user to specify DIIS convergence rather than
|
|
second-order convergence for the SCF calculation. The form of the
|
|
directive is as follows:
|
|
|
|
\begin{verbatim}
|
|
DIIS
|
|
\end{verbatim}
|
|
|
|
The implementation of this option is currently fairly rudimentary. It
|
|
does not have level-shifting and damping, and does not support open
|
|
shells or UHF. It is provided on an ``as is'' basis, and should be
|
|
used with caution.
|
|
|
|
When the \verb+DIIS+ directive is specified in the input, the user has
|
|
the additional option of specifying the size of the subspace for the
|
|
DIIS extrapolation. This is accomplished with the \verb+DIISBAS+
|
|
directive, which is of the form:
|
|
\begin{verbatim}
|
|
DIISBAS <integer diisbas default 5>
|
|
\end{verbatim}
|
|
The default of 5 should be adequate for most applications, but may be
|
|
increased if convergence is poor. On large systems, it may be necessary
|
|
to specify a lower value for \verb+diisbas+, to conserve memory.
|
|
|
|
\section{{\tt DIRECT} and {\tt SEMIDIRECT} --- recomputation of integrals}
|
|
\label{sec:semidirect}
|
|
|
|
In the context of SCF calculations direct means that all integrals are
|
|
recomputed as required and none are stored. The other extreme are
|
|
disk- or memory-resident (sometimes termed conventional) calculations
|
|
in which all integrals are computed once and stored. Semi-direct
|
|
calculations are between these two extremes with some integrals being
|
|
precomputed and stored, and all other integrals being recomputed as
|
|
necessary.
|
|
|
|
The default behavior of the SCF module is
|
|
\begin{itemize}
|
|
\item If enough memory is available, the integrals are computed once
|
|
and are cached in memory.
|
|
\item If there is not enough memory to store all the integrals at
|
|
once, then 95\% of the available disk space in the scratch directory
|
|
(see Section \ref{sec:dirs}) is assumed to be available for this
|
|
purpose, and as many integrals as possible are cached on disk (with
|
|
no memory being used for caching). Some attempt is made to store
|
|
the most expensive integrals in the cache.
|
|
\item If there is not enough room in memory or on disk for all the
|
|
integrals, then the ones that are not cached are recomputed in a
|
|
semidirect fashion.
|
|
\end{itemize}
|
|
|
|
The integral file is deleted at the end of a calculation, so it is not
|
|
possible to restart a semidirect calculation when the integrals are
|
|
cached in memory or on disk. Many computer systems (e.g., the EMSL
|
|
IBM SP) clear the fast scratch space at the end of each job, adding a
|
|
further complication to the problem of restarting a {\em parallel}
|
|
semidirect calculation.
|
|
|
|
%Under some situations, it is possible to
|
|
%restart from integrals on disk, but this capability will not be made
|
|
%widely available until a later date.
|
|
|
|
On the IBM SP or any other computer with fast disks local to each
|
|
processor, semidirect calculation offers the best behavior. It can
|
|
result in {\em quadratic speedup} as more processors are added.
|
|
|
|
A fully direct calculation (with recomputation of the integrals at
|
|
each iteration) is forced by specifying the directive
|
|
|
|
\begin{verbatim}
|
|
DIRECT
|
|
\end{verbatim}
|
|
|
|
Alternatively, the \verb+SEMIDIRECT+ directive can be used to control
|
|
the default semidirect calculation by defining the amount of disk
|
|
space and the cache memory size. The form of this directive is as
|
|
follows:
|
|
|
|
\begin{verbatim}
|
|
SEMIDIRECT [filesize <integer filesize default disksize>]
|
|
[memsize <integer memsize default available>]
|
|
[filename <string filename default $file_prefix.aoints$>]
|
|
\end{verbatim}
|
|
|
|
The keyword \verb+FILESIZE+ allows the user to specify the amount of
|
|
disk space to be used per process for storing the integrals in 64-bit
|
|
words. Similarly, the keyword \verb+MEMSIZE+ allows the user to
|
|
specify the number of 64-bit words to be used per process for caching
|
|
integrals in memory. (Note: If the amount of storage space specified
|
|
by the entry for \verb+memsize+ is not available, the code cuts the
|
|
value in half and checks again for available space. This process is
|
|
repeated until the request is satisfied.)
|
|
|
|
By default, the integral files are placed into the scratch directory
|
|
(see Section \ref{sec:dirs}). Specifying the keyword \verb+FILENAME+
|
|
overrides this default. The user-specified name entered in the string
|
|
\verb+filename+ has the process number appended to it, so that each
|
|
process has a distinct file but with a common base-name and directory.
|
|
Therefore, it is not possible to use this keyword to specify different
|
|
disks for different processes. The \verb+SCRATCH_DIR+ directive (see
|
|
Section \ref{sec:dirs}) can be used for this purpose.
|
|
|
|
For example, to force full recomputation of all integrals:
|
|
\begin{verbatim}
|
|
direct
|
|
\end{verbatim}
|
|
|
|
Exactly the same result could be obtained by entering the directive:
|
|
\begin{verbatim}
|
|
semidirect filesize 0 memsize 0
|
|
\end{verbatim}
|
|
|
|
To disable the use of memory for caching integrals and limit disk
|
|
usage by each process to 100 megawords (MW):
|
|
\begin{verbatim}
|
|
semidirect memsize 0 filesize 100000000
|
|
\end{verbatim}
|
|
|
|
The integral records are typically 32769 words long and any non-zero
|
|
value for \verb+filesize+ or \verb+memsize+ should be enough to hold
|
|
at least one record.
|
|
|
|
|
|
\subsection{Integral File Size and Format for the SCF Module}
|
|
|
|
The file format is rather complex, since it accommodates a variety of
|
|
packing and compression options and the distribution of data. This
|
|
section presents some information that may help the user understand
|
|
the output, and illustrates how to use the output information to
|
|
estimate file sizes.
|
|
|
|
If integrals are stored with a threshold of greater than $10^{-10}$,
|
|
then the integrals are stored in a 32-bit fixed-point format (with
|
|
appropriate treatment for large values to retain precision). If
|
|
integrals are stored with a threshold less than $10^{-10}$, however,
|
|
the values are stored in 64-bit floating-point format. If a
|
|
replicated-data calculation is being run, then 8 bits are used for
|
|
each basis function label, unless there are more than 256 functions,
|
|
in which case 16 bits are used. If distributed data is being used,
|
|
then the labels are always packed to 8 bits (the distributed blocks
|
|
always being less than 256; labels are relative to the start of the
|
|
block).
|
|
|
|
Thus, the number ($W$) of 64-bit words required to store $N$
|
|
integrals, may be computed as
|
|
\begin{displaymath}
|
|
W = \left\{ \\
|
|
\begin{array}{c}
|
|
N \mbox{ , 8-bit labels and 32-bit values} \\
|
|
\frac{3}{2}N \mbox{ , 16-bit labels and 32-bit values} \\
|
|
\frac{3}{2}N \mbox{ , 8-bit labels and 64-bit values} \\
|
|
2N \mbox{ , 16-bit labels and 64-bit values}
|
|
\end{array}
|
|
\right.
|
|
\end{displaymath}
|
|
|
|
The actual number of words required can
|
|
exceed this computed value by up to one percent, due to
|
|
bookkeeping overhead, and because the file itself is
|
|
organized into fixed-size records.
|
|
|
|
With at least the default print level, all semidirect (not direct)
|
|
calculations will print out information about the integral file and
|
|
the number of integrals computed. The form of this output is as
|
|
follows:
|
|
|
|
\begin{verbatim}
|
|
Integral file = ./c6h6.aoints.0
|
|
Record size in doubles = 32769 No. of integs per rec = 32768
|
|
Max. records in memory = 3 Max. records in file = 5
|
|
No. of bits per label = 8 No. of bits per value = 32
|
|
|
|
#quartets = 2.0D+04 #integrals = 7.9D+05 direct = 63.6% cached = 36.4%
|
|
\end{verbatim}
|
|
|
|
The file information above relates only to process 0. The line of
|
|
information about the number of quartets, integrals, etc., is a sum
|
|
over all processes.
|
|
|
|
When the integral file is closed, additional information of the following
|
|
form is printed:
|
|
|
|
\begin{verbatim}
|
|
------------------------------------------------------------
|
|
EAF file 0: "./c6h6.aoints.0" size=262152 bytes
|
|
------------------------------------------------------------
|
|
write read awrite aread wait
|
|
----- ---- ------ ----- ----
|
|
calls: 6 12 0 0 0
|
|
data(b): 1.57e+06 3.15e+06 0.00e+00 0.00e+00
|
|
time(s): 1.09e-01 3.12e-02 0.00e+00
|
|
rate(mb/s): 1.44e+01 1.01e+02
|
|
------------------------------------------------------------
|
|
|
|
Parallel integral file used 4 records with 0 large values
|
|
\end{verbatim}
|
|
Again, the detailed file information relates just to process 0, but
|
|
the final line indicates the total number of integral records stored
|
|
by all processes.
|
|
|
|
This information may be used to optimize subsequent calculations, for
|
|
instance by assigning more memory or disk space.
|
|
|
|
\section{SCF Convergence Control Options}
|
|
\label{sec:scfconv}
|
|
|
|
{\em Note to users:} It is desired that the SCF program converge
|
|
reliably with the default options for a wide variety of molecules. In
|
|
addition, it should be guaranteed to converge for any system, with
|
|
sufficient iterations. Please report significant convergence problems
|
|
to \verb+nwchem+-\verb+support@+\-\verb+emsl.pnl.gov+, and include the
|
|
input file.
|
|
|
|
% An understanding of the output of the SCF program and the options
|
|
% controlling convergence requires some knowledge of the convergence
|
|
% scheme.
|
|
|
|
The SCF program uses a preconditioned conjugate gradient (PCG) method
|
|
that is unconditionally convergent. Basically, a search direction is
|
|
generated by multiplying the orbital gradient (the derivative of the
|
|
energy with respect to the orbital rotations) by an approximation to
|
|
the inverse of the level-shifted orbital Hessian. In the initial
|
|
iterations (see Section \ref{sec:nrswitch}), an inexpensive
|
|
one-electron approximation to the inverse orbital Hessian is used.
|
|
Closer to convergence, the full orbital Hessian is used, which should
|
|
provide quadratic convergence. For both the full or one-electron
|
|
orbital Hessians, the inverse-Hessian matrix-vector product is formed
|
|
iteratively. Subsequently, an approximate line search is performed
|
|
along the new search direction. If the exact Hessian is being
|
|
employed, then the line search should require a single step (of
|
|
unity). Preconditioning with approximate Hessians may require
|
|
additional steps, especially in the initial iterations. It is the
|
|
(approximate) line search that provides the convergence guarantee.
|
|
The iterations required to solve the linear equations are referred to
|
|
as micro-iterations. A macro-iteration comprises both the iterative
|
|
solution and a line search.
|
|
|
|
Level-shifting plays the same role in this algorithm as
|
|
it does in the conventional iterative solution of the SCF equations.
|
|
The approximate Hessian used for preconditioning should be positive
|
|
definite. If this is not the case, then level-shifting by a positive
|
|
constant ($\Delta$) serves to make the preconditioning matrix positive
|
|
definite, by adding $\Delta$ to all of its eigenvalues. The
|
|
level-shifts employed for the RHF orbital Hessian should be
|
|
approximately four times (only twice for UHF) the value that one would
|
|
employ in a conventional SCF\footnote{This can be seen by considering
|
|
a one-electron approximation to the closed-shell RHF Hessian in
|
|
canonical orbitals, $A_{ia,jb} = 4 \delta_{ij} \delta_{ab}
|
|
(\epsilon_a - \epsilon_i)$. Similarly, the level shift
|
|
should be twice as large for UHF.}. Level-shifting is automatically enabled
|
|
in the early iterations, and the default options suffice for most test
|
|
cases.
|
|
|
|
So why do things go wrong and what can be done to fix convergence
|
|
problems? Most problems encountered so far arise either poor initial
|
|
guesses or from small or negative eigenvalues of the orbital Hessian.
|
|
The atomic orbital guess is usually very good. However, in
|
|
calculations on charged systems, especially with open shells,
|
|
incorrect initial occupations may result. The SCF might then converge
|
|
very slowly since very large orbital rotations might be required to
|
|
achieve the correct occupation or move charge large distances in the
|
|
molecule. Possible actions are
|
|
\begin{itemize}
|
|
\item Modify the atomic guess by assigning charges to the atoms
|
|
known to carry substantial charges (Section \ref{sec:atomscf})
|
|
\item Examining an analysis of the initial orbitals (Section
|
|
\ref{sec:scfprint}) and then swapping them to attain the desired
|
|
occupation (Section \ref{sec:vectors}).
|
|
\item Converging the calculation in a minimal basis set, which is
|
|
usually easier, and then projecting into a larger basis set (Section
|
|
\ref{sec:vectors}).
|
|
\item Using the fragment orbital initial guess (Section
|
|
\ref{sec:fragguess}).
|
|
\end{itemize}
|
|
|
|
Small or negative Hessian eigenvalues can occur even though the
|
|
calculation seem to be close to convergence (as measured by the
|
|
gradient norm, or the off-diagonal Fock matrix elements). Small
|
|
eigenvalues will cause the iterative linear equation solver to
|
|
converge slowly, resulting in an excessive number of micro-iterations.
|
|
This makes the SCF expensive in terms of computation time, and it is
|
|
possible to exceed the maximum number of iterations without achieving
|
|
the accuracy required for quadratic convergence --- which causes more
|
|
macro-iterations to be performed.
|
|
|
|
Two main options are available when a problem will not converge:
|
|
Newton-Raphson can be disabled temporarily or permanently (see Section
|
|
\ref{sec:nrswitch}), and level-shifting can be applied to the matrix
|
|
(see Section \ref{sec:level}). In some cases, both options may be
|
|
necessary to achieve final convergence.
|
|
|
|
If there is reason to suspect a negative eigenvalue, the first course
|
|
is to disable the Newton-Raphson iteration until the solution is
|
|
closer to convergence. It may be necessary to disable it completely.
|
|
At some point close to convergence, the Hessian will be positive
|
|
definite, so disabling Newton-Raphson should yield a solution with
|
|
approximately the same convergence rate as DIIS.
|
|
|
|
If temporarily disabling Newton-Raphson is not sufficient to achieve
|
|
convergence, it may be necessary to disable it entirely and apply a
|
|
small level-shift to the approximate Hessian. This should improve the
|
|
convergence rate of the micro-iterations and stabilize the
|
|
macro-iterations. The level-shifting will destroy exact quadratic
|
|
convergence, but the optimization process is automatically adjusted to
|
|
reflect this by enforcing conjugacy and reducing the accuracy to which
|
|
the linear equations are solved. The net result of this is that the
|
|
solution will do more macro-iterations, but each one should take less
|
|
time than it would with the unshifted Hessian.
|
|
|
|
The following sections describe the directives needed to disable the
|
|
Newton-Raphson iteration and specify level-shifting.
|
|
|
|
\section{{\tt NR} --- controlling the Newton-Raphson}
|
|
\label{sec:nrswitch}
|
|
|
|
\begin{verbatim}
|
|
NR <real nr_switch default 0.1>
|
|
\end{verbatim}
|
|
|
|
The exact orbital Hessian is adopted as the preconditioner when the
|
|
maximum element of the orbital gradient is below the value specified
|
|
for \verb+nr_switch+. The default value is 0.1, which means that
|
|
Newton-Raphson will be disabled until the maximum value of the orbital
|
|
gradient (twice the largest off-diagonal Fock matrix element) is less
|
|
than 0.1. To disable the Newton-Raphson entirely, the
|
|
value of \verb+nr_switch+ must be set to zero. The directive to accomplish
|
|
this is as follows:
|
|
\begin{verbatim}
|
|
nr 0
|
|
\end{verbatim}
|
|
|
|
\section{{\tt LEVEL} --- level-shifting the orbital Hessian}
|
|
\label{sec:level}
|
|
|
|
This directive allows the user to specify level-shifting to obtain a
|
|
positive-definite preconditioning matrix for the SCF solution
|
|
procedure. Separate level shifts can be set for the first-order
|
|
convergent one-electron approximation to the Hessian used with the
|
|
preconditioned conjugate gradient (PCG) method, and for the full
|
|
Hessian used with the Newton-Raphson (NR) approach. It is also
|
|
possible to change the level-shift automatically as the solution
|
|
attains some specified accuracy. The form of the directive is as
|
|
follows:
|
|
\begin{verbatim}
|
|
LEVEL [pcg <real initial default 20.0> \
|
|
[<real tol default 0.5> <real final default 0.0>]] \
|
|
[nr <real initial default 0.0> \
|
|
[<real tol default 0.0> <real final default 0.0>]]
|
|
\end{verbatim}
|
|
|
|
This directive contains only two keywords: one for the PCG method and
|
|
the other for the exact Hessian (Newton Raphson, or NR). Use of PCG
|
|
or NR is determined by the input specified for \verb+nr_switch+ on the
|
|
\verb+NR+ directive, Section \ref{sec:nrswitch} above.
|
|
|
|
Specifying the keyword \verb+pcg+ on the \verb+LEVEL+ directive allows
|
|
the user to define the level shifting for the approximate (i.e., PCG)
|
|
method. Specifying the keyword \verb+nr+ allows the user to define
|
|
the level shifting for the exact Hessians. In both options, the
|
|
initial level shift is defined by the value specified for the variable
|
|
\verb+initial+. Optionally, \verb+tol+ can be specified independently
|
|
with each keyword to define the level of accuracy that must be
|
|
attained in the solution before the level shifting is changed to the
|
|
value specified by input in the real variable \verb+final+. Level
|
|
shifts and gradient thresholds are specified in atomic units.
|
|
|
|
For the PCG method (as specified using the keyword \verb+pcg+), the
|
|
defaults for this input are 20.0 for \verb+initial+, 0.5 for
|
|
\verb+tol+, and 0.0 for \verb+final+. This means that the
|
|
approximate Hessian will be shifted by 20.0 until the maximum element
|
|
of the gradient falls below 0.5, at which point the shift will be set
|
|
to zero.
|
|
|
|
For the exact Hessian (as specified using the keyword \verb+nr+), the
|
|
defaults are all zero. The exact Hessian is usually not shifted since
|
|
this destroys quadratic convergence. An example of an input directive
|
|
that applies a shift of 0.2 to the exact Hessian is as follows:
|
|
\begin{verbatim}
|
|
level nr 0.2
|
|
\end{verbatim}
|
|
|
|
To apply this shift to the exact Hessian only until the maximum
|
|
element of the gradient falls below 0.005, the required input
|
|
directive is as follows:
|
|
\begin{verbatim}
|
|
level nr 0.2 0.005 0
|
|
\end{verbatim}
|
|
|
|
Note that in both of these examples, the parameters for the PCG method
|
|
are at the default values. To obtain values different from the
|
|
defaults, the keyword \verb+pcg+ must also be specified. For example,
|
|
to specify the level shifting in the above example for the exact
|
|
Hessian {\em and} non-default shifting for the PCG method, the
|
|
directive would be something like the following:
|
|
\begin{verbatim}
|
|
level pcg 20 0.3 0.0 nr 0.2 0.005 0.0
|
|
\end{verbatim}
|
|
|
|
This input will cause the PCG method to be level-shifted by 20.0 until
|
|
the maximum element of the gradient falls below 0.3, then the shift
|
|
will be zero. For the exact Hessian, the level shifting is initially
|
|
0.2, until the maximum element falls below 0.005, after which the
|
|
shift is zero.
|
|
|
|
The default options correspond to
|
|
\begin{verbatim}
|
|
level pcg 20 0.5 0 nr 0 0 0
|
|
\end{verbatim}
|
|
|
|
\section{Orbital Localization}
|
|
\label{orbloc}
|
|
The SCF module includes an {\em experimental} implementation of
|
|
orbital localization, including Foster-Boys and Pipek-Mezey which only
|
|
works for closed-shell (RHF) wavefunctions. There is currently no
|
|
input in the SCF block to control this so the \verb+SET+ directive
|
|
(Section \ref{sec:set}) must be used.
|
|
|
|
The directive
|
|
\begin{verbatim}
|
|
set scf:localize t
|
|
\end{verbatim}
|
|
will separately localize the core, valence, and virtual orbital spaces
|
|
using the Pipek-Mezey algorithm. If the additional directive
|
|
\begin{verbatim}
|
|
set scf:loctype FB
|
|
\end{verbatim}
|
|
is included, then the Foster-boys algorithm is used. The partitioning
|
|
of core-orbitals is performed using the atomic information described
|
|
in Section \ref{mp2:core}.
|
|
|
|
In the next release, this functionality will be extended to included all
|
|
wavefunctions using molecular orbitals.
|
|
|
|
|
|
\newpage
|
|
\section{Printing Information from the SCF Module}
|
|
\label{sec:scfprint}
|
|
|
|
All output from the SCF module is controlled using the \verb+PRINT+
|
|
directive described in Section \ref{sec:printcontrol}. The following
|
|
list describes the items from SCF that are currently under direct
|
|
print control, along with the print level for each one.
|
|
|
|
\begin{table}[htbp]
|
|
\begin{center}
|
|
\begin{tabular}{lcc}
|
|
{\bf Name} & {\bf Print Level} & {\bf Description} \\
|
|
``atomic guess density'' & debug & guess density matrix \\
|
|
``atomic scf'' & debug & details of atomic SCF \\
|
|
``mo guess'' & default & brief info from mo guess \\
|
|
``information'' & low & results \\
|
|
``initial vectors'' & debug & \\
|
|
``intermediate vectors'' & debug & \\
|
|
``final vectors'' & debug & \\
|
|
``final vectors analysis'' & default & \\
|
|
``initial vectors analysis'' & never & \\
|
|
``intermediate evals'' & debug & \\
|
|
``final evals'' & default & \\
|
|
``schwarz'' & high & integral screening info \& stats at completion\\
|
|
``screening statistics'' & debug & display stats after every Fock build \\
|
|
``geometry'' & high & \\
|
|
``symmetry'' & debug & detailed symmetry info \\
|
|
``basis'' & high & \\
|
|
``geombas'' & debug & detailed basis map info \\
|
|
``vectors i/o'' & default & report vectors I/O \\
|
|
``parameters'' & default & convergence parameters \\
|
|
``convergence'' & default & info each iteration \\
|
|
``mulliken ao'' & never & Mulliken population of basis functions
|
|
\end{tabular}
|
|
\end{center}
|
|
\caption{SCF Print Control Specifications}
|
|
\end{table}
|
|
|
|
\newpage
|
|
\section{Hartree-Fock or SCF, MCSCF and MP2 Gradients}
|
|
|
|
\input{scfgrad}
|
|
|
|
% LocalWords: MO's
|