.. index::
single: Program; SLAPAF
single: SLAPAF
.. _UG\:sec\:slapaf:
:program:`SlapAf`
=================
.. only:: html
.. contents::
:local:
:backlinks: none
.. xmldoc::
%%Description:
This program is a general purpose facility for geometry
optimization. At present, it is tailored to use analytical gradients
produced by ALASKA.
SLAPAF also computes an approximate Hessian.
Provided with the first order derivative with respect to nuclear displacements
the program is capable to optimize molecular structures with or
without constraints for minima or
transition states. This will be achieved with a quasi-Newton approach
in combination with 2nd ranks updates of the approximate Hessian or
with the use of a computed (analytic or numerical) Hessian.
Note that *if* a computed Hessian is available on the
:file:`RUNFILE` then it will be used rather than the approximate Hessian generated by :program:`Slapaf`.
On completion of an optimization :program:`SlapAf` will automatically execute a single energy evaluation,
this can be disabled with the :kword:`NOLAst` keyword.
.. _UG\:sec\:slapaf_description:
Description
-----------
:program:`SlapAf` has three different ways in selecting the
basis for the displacements during the optimization.
The first format require user input (not recommended), whereas the two other options are totally black-boxed.
The formats are:
#. the old format as in |molcasiii|, which is user specified.
The internal coordinates
are here represented as linear combination of internal coordinates
(such as bonds, angles, torsions, out of plane angles, Cartesian coordinates
and separation of centers of mass) and the linear combinations are totally defined
by user input.
This format does also require the user to specify the
Hessian (default a diagonal matrix).
This option *allows* for frozen internal coordinates.
#. the second format is an automatic
option which employs the Cartesian eigenvectors of the approximative Hessian (generated by the
Hessian model functional :cite:`HMF`).
#. the third format (this is the recommend and default) is an automatic option which utilizes linear combinations
of some curvilinear coordinates (stretches, bends, and torsions).
This implementation
has two variations. The first can be viewed as the conventional use of
non-redundant internal coordinates :cite:`nric1,nric2,nric3`.
The second variation is a force constant weighted (FCW)
redundant space (the HWRS option) version of the former
implementation :cite:`Lindh:97`.
All three formats of internal coordinates can be used in combinations with
constraints on the molecular parameters or other type of constraints as for
example energy differences.
The displacements are symmetry adapted
and any rotation and translation if present is deleted.
The relaxation is symmetry preserving.
.. _UG\:sec\:slapaf_dependencies:
Dependencies
------------
:program:`SlapAf` depends on the results of :program:`ALASKA` and also possibly
on :program:`MCKINLEY` and :program:`MCLR`.
.. _UG\:sec\:slapaf_files:
Files
-----
Input files
...........
Apart from the standard input file :program:`SlapAf` will use the following input
files.
.. class:: filelist
:file:`RUNFILE`
File for communication of auxiliary information. If a computed Hessian is available on this file it will be used rather than
the approximate Hessian generated by :program:`Slapaf`.
:file:`RUNFILE2`
File for communication of auxiliary information of the "ground state" in case of minimum energy cross point optimizations.
:file:`RUNOLD`
File for communication of auxiliary information for reading an old Hessian matrix from a previous geometry optimization.
Output files
............
In addition to the standard output file :program:`SlapAf` will use the following output
files.
.. class:: filelist
:file:`RUNFILE`
File for communication of auxiliary information.
:file:`RUNFILE2`
File for communication of auxiliary information of the "ground state" in case of minimum energy cross point optimizations.
:file:`MD_GEO`
Molden input file for geometry optimization analysis.
:file:`MD_MEP`
Molden input file for minimum energy path (MEP).
:file:`MD_SADDLE`
Molden input file for energy path (MEP) of a Saddle TS optimization.
:file:`MD_IRC`
Molden input file for intrinsic reaction coordinate analysis of a TS.
:file:`MD_FREQ`
Molden input file for harmonic frequency analysis.
:file:`UNSYM`
ASCII file where all essential information, like geometry, Hessian normal modes and dipole
derivatives are stored.
:file:`STRUCTURE`
Output file with a statistics of geometry optimization convergence.
.. _UG\:sec\:slapaf_input:
Input
-----
:program:`SlapAf` will as standard
provided with an energy and a corresponding gradient
update the geometry (optimize).
Possible update methods include different quasi-Newton methods.
The program will also provide for updates of the Hessian.
The program has a number of different variable metric methods available for
the Hessian update.
This section describes the input to the :program:`SlapAf` program.
This section describes the input to the
:program:`SLAPAF` program in the |molcas| program system. The input starts
with the program name ::
&SLAPAF
There are no compulsory keywords
Optional convergence control keywords
.. class:: keywordlist
:kword:`ITERations`
Maximum number of iterations which
will be allowed in the relaxation procedure. Default is 500
iterations, however, if environment variable :variable:`MOLCAS_MAXITER` has been exported by the user
this is the assumed default value.
.. xmldoc::
%%Keyword: Iterations
Specify the max number of iterations which
will be allowed in the relaxation procedure. Default is 500
iterations however, if MOLCAS_MAXITER has been exported by the user
this is the assumed default value.
:kword:`THRShld`
Enter two real numbers which specifies the convergence criterion with respect to the
energy change and the norm of the gradient. The defaults are
0.0 and 3.0D-4 au for Gaussian convergence criteria
(which normally do not consider the energy change), and
1.0D-6 and 3.0D-4 for Baker criteria (see the :kword:`BAKER` keyword).
.. xmldoc::
%%Keyword: Thrshld
Enter two real numbers
which specifies the convergence criterion with respect to the
energy change and the norm of the gradient.
The defaults are 0.0 and 3.0D-4 for Gaussian, and 1.0D-6 and 3.0D-4 for Baker.
:kword:`BAKEr`
Activate convergence criteria according to Baker :cite:`Baker`.
Default is to use the convergence criteria as in the Gaussian
program :cite:`GAUSSIAN94`.
.. xmldoc::
%%Keyword: Baker
Activate convergence criteria according to Baker.
Default is to use the convergence criteria as in the Gaussian
program.
:kword:`MAXStep`
This keyword is followed by the value which defines the seed of largest
change of the internal coordinates which will be accepted. A
change which is larger is reduced to the max value. The value is dynamically modified each iterations.
The default value is 0.3 au or rad.
.. xmldoc::
%%Keyword: Maxstep
Enter the value which defines the seed of largest
change of the internal coordinates which will be accepted. A
change which is larger is reduced to the max value. The value is dynamically modified each iterations.
The default
value is 0.3 au or rad.
:kword:`CNWEight`
Sets the maximum weight assigned to the fulfillment of the constraints, relative to the step taken in the
complementary space for energy minimization. The step in the constraint space is truncated to be at most as
large as the step in the minimization space, or half the maximum total step, whichever is larger, multiplied
by this value. Default is 1.0.
.. xmldoc::
%%Keyword: CnWeight
Sets the maximum weight assigned to the fulfillment of
the constraints, relative to the step taken in the
complementary space for energy minimization.
:kword:`TOLErance`
Controls how strictly the constraints (if any) must be satisfied at convergence. The default value
is very large, such that this criterion is always met, and only the gradient and maximum step (or
energy difference) control convergence. If you set this keyword to some value, a constrained optimization
will only converge if the maximum error in any constraint is lower than this number (in atomic units,
and radians).
.. xmldoc::
%%Keyword: Tolerance
Controls how strictly the constraints must be satisfied at convergence.
Optional coordinate selection keywords
.. class:: keywordlist
:kword:`CARTesian`
Activate :program:`SlapAf` to use the eigenvectors
of the approximative Hessian expressed in Cartesian as the
definition of the internal coordinates. The default is to
use the FCW non-redundant internal coordinates.
The Hessian will be modeled by the Hessian Model Functional.
.. xmldoc::
.. xmldoc::
%%Keyword: Cartesian
Activate SlapAf to use the eigenvectors
of the approximative Hessian expressed in Cartesian as the
definition of the internal coordinates. The default is to
use the FCW non-redundant internal coordinates.
The Hessian will be modeled by the Hessian Model Functional.
:kword:`INTErnal`
This marks the start of the definition of the internal
coordinates. This section is always ended by the keyword
:kword:`End of Internal`.
For a complete description of this
keyword see
:numref:`UG:sec:definition_of_internal_coordinates`.
This option will also use a diagonal matrix as default for
the Hessian matrix.
The default is to
use the FCW non-redundant internal coordinates.
.. xmldoc::
%%Keyword: Internal
This marks the start of the definition of the internal
coordinates.
This section is always ended by the keyword "End of Internal".
Consult the manual for details.
:kword:`HWRS`
Use the force constant weighted (FCW) redundant space version of the
nonredundant internal coordinates. This is the default.
The Hessian will be modeled by the Hessian Model Functional.
.. xmldoc::
%%Keyword: HWRS
Use the force constant weighted (FCW) redundant space version of the
nonredundant internal coordinates.
The Hessian will be modeled by the Hessian Model Functional.
This is the default.
:kword:`NOHWrs`
Disable the use of the force constant weighted redundant space version of the
nonredundant internal coordinates. The default is to use the HWRS option.
The Hessian will be modeled by the Hessian Model Functional.
.. xmldoc::
%%Keyword: NoHWRS
Disable the use of the force constant weighted redundant space version of the
nonredundant internal coordinates. The default is to use the HWRS option.
The Hessian will be modeled by the Hessian Model Functional.
:kword:`FUZZ`
When automatically generating the primitive internal coordinates, the system may
end up in disconnected fragments, in which case additional bonds are defined
between the fragments.
This keyword controls how many inter-fragment bonds are added. Bonds are generated
between the closest atoms of two fragments, and all pairs of atoms in separate
fragments at a distance up to the specified value longer.
The value can be followed with the unit BOHR or ANGSTROM. The default is 0.5 a.u.
.. xmldoc::
%%Keyword: Fuzz
When automatically generating the primitive internal coordinates, the system may
end up in disconnected fragments, in which case additional bonds are defined
between the fragments.
This keyword controls how many inter-fragment bonds are added. Bonds are generated
between the closest atoms of two fragments, and all pairs of atoms in separate
fragments at a distance up to the specified value longer.
The value can be followed with the unit BOHR or ANGSTROM. The default is 0.5 a.u.
.. xmldoc::
Optional Hessian update keywords
.. class:: keywordlist
:kword:`HUPDate`
Method used for updating the Hessian matrix. It must be one of:
* ``None`` --- No update is applied.
* ``BFGS`` --- Activate update according to Broyden--Fletcher--Goldfarb--Shanno.
This is the default.
* ``MSP`` --- Activate the Murtagh--Sargent--Powell update according to Bofill :cite:`MSP`.
This update is preferred for the location of transition states.
* ``EU`` --- Activate the EU update according to Bofill :cite:`EU`.
This update can be used for the location of transition states.
* ``TS-BFGS`` --- Activate the TS-BFGS update according to Bofill :cite:`EU`.
This update can be used for the location of minima or transition states.
.. xmldoc::
.. xmldoc::
%%Keyword: HUpdate
Method used for updating the Hessian matrix.
It must be one of: None, BFGS, MSP, EU, TS-BFGS.
:kword:`UORDer`
Order the gradients and displacements vectors according to Schlegel prior to
the update of the Hessian. Default is no reorder.
.. xmldoc::
%%Keyword: UORDer
Order the gradients and displacements vectors according to Schlegel prior to
the update of the Hessian. Default is no reorder.
:kword:`WINDow`
Maximum number of previous iterations to include in the Hessian update.
When using RVO (see :kword:`KRIGing` keyword), the maximum number of sample points used is twice this value.
Default is 5.
.. xmldoc::
%%Keyword: WINDow
Maximum number of previous iterations to include in the Hessian update.
Default is 5.
.. xmldoc::
Optional optimization procedure keywords
.. class:: keywordlist
:kword:`NOLIne`
Disable line search. Default is to use line search for minima.
.. xmldoc::
.. xmldoc::
%%Keyword: Noline
Disable line search. Default is to use line search for minima.
:kword:`RATIonal`
Activate geometry optimization using the restricted step Rational Functional optimization :cite:`rf,rs-rf`,
this is the default.
.. xmldoc::
%%Keyword: Rational
Activate geometry optimization using the restricted step Rational Functional optimization,
this is the default.
:kword:`C1-Diis`
Activate geometry optimization using the C1-GDIIS method :cite:`gdiis,diis1,diis2`.
The default is to use the Rational Functional approach.
.. xmldoc::
%%Keyword: C1-diis
Activate geometry optimization using the C1-GDIIS method.
The default is to use the Rational Functional approach.
:kword:`C2-Diis`
Activate geometry optimization using the C2-GDIIS method :cite:`c2-diis`.
The default is to use the Rational Functional approach.
.. xmldoc::
%%Keyword: C2-diis
Activate geometry optimization using the C2-GDIIS method.
The default is to use the Rational Functional approach.
:kword:`DXDX`
This option is associated to the use of the C1- and C2-GDIIS
procedures. This option will activate the computation of the
so-called error matrix elements as :math:`e=\delta x^{\text{T}}\delta x`,
where :math:`\delta x` is the displacement vector.
.. xmldoc::
%%Keyword: dxdx
This option is associated to the use of the C1- and C2-GDIIS
procedures. This option will activate the computation of the
so-called error matrix elements as e=dx(T)dx,
where dx is the displacement vector.
:kword:`DXG`
This option is associated to the use of the C1- and C2-GDIIS
procedures. This option will activate the computation of the
so-called error matrix elements as :math:`e=\delta x^{\text{T}}g`,
where :math:`\delta x` is the displacement vector and :math:`g` is the
gradient vector.
.. xmldoc::
%%Keyword: dxg
This option is associated to the use of the C1- and C2-GDIIS
procedures. This option will activate the computation of the
so-called error matrix elements as e=dx(T)g,
where dx is the displacement vector and g is the
gradient vector.
:kword:`GDX`
See above.
.. xmldoc::
%%Keyword: gdx
See the dxg keyword.
:kword:`GG`
This option is associated to the use of the C1- and C2-GDIIS
procedures. This option will activate the computation of the
so-called error matrix elements as :math:`e=g^{\text{T}}g`,
where :math:`g` is the gradient vector. This is the default.
.. xmldoc::
%%Keyword: gg
This option is associated to the use of the C1- and C2-GDIIS
procedures. This option will activate the computation of the
so-called error matrix elements as e=g(T)g,
where g is the gradient vector. This is the default.
:kword:`NEWTon`
Activate geometry optimization using the standard quasi-Newton approach.
The default is to use the Rational Functional approach.
.. xmldoc::
%%Keyword: Newton
Activate geometry optimization using the standard quasi-Newton approach.
The default is to use the Rational Functional approach.
:kword:`RS-P-rfo`
Activate RS-P-RFO :cite:`rs-rf` as default for TS-search. Default is RS-I-RFO.
.. xmldoc::
%%Keyword: RS-P-RFO
Activate RS-P-RFO as default for TS-search. Default is RS-I-RFO.
:kword:`TS`
Keyword for optimization of transition states. This flag will activate
the use of the mode following rational functional approach :cite:`mfrf`.
The mode to follow can either be the one with the lowest eigenvalue (if positive
it will be changed to a negative value) or by the eigenvector which index
is specified by the :kword:`MODE` keyword (see below). The keyword will also
activate the Murtagh--Sargent--Powell update of the Hessian and inactivate
line search. This keyword will also enforce that the Hessian has the
right index (i.e. one negative eigenvalue).
.. xmldoc::
%%Keyword: TS
Keyword for optimization of transition states. This flag will activate
the use of the mode following rational functional approach.
The mode to follow can either be the one with the lowest eigenvalue (if positive
it will be changed to a negative value) or by the eigenvector which index
is specified by the MODE keyword. The keyword will also
activate the Murtagh-Sargent-Powell update of the Hessian and inactivate
line search. This keyword will also enforce that the Hessian has the
right index (i.e. one negative eigenvalue).
:kword:`MODE`
Specification of the Hessian eigenvector index, this mode will be followed
by the mode following RF method for optimization of transition states.
The keyword card is followed by a single card specifying the eigenvector index.
.. xmldoc::
%%Keyword: Mode
Specification of the Hessian eigenvector index, this mode will be followed
by the mode following RF method for optimization of transition states.
The keyword card is followed by a single card specifying the eigenvector index.
:kword:`FINDTS`
Enable a constrained optimization to release the constraints and locate
a transition state if negative curvature is encountered and the
gradient norm is below a specific threshold (see the :kword:`GNRM` option).
Keyword :kword:`TSCOnstraints` should be used in combination with :kword:`FINDTS`.
.. xmldoc::
%%Keyword: FindTS
Enable a constrained optimization to release the constraints and locate
a transition state if negative curvature is encountered and the
gradient norm is below a specific threshold (see the GNRM option).
Keyword TSCOnstraints should be used in combination with FINDTS.
:kword:`TSCOnstraints`
Specify constraints that will be active during the initial stage of an
optimization with :kword:`FINDTS`. When negative curvature and low
gradient are encountered, these constraints will be released and
other constraints will remain active. If this block is not given in
the input, all constraints will be released. The syntax of this
keyword is exactly like normal constraints, and it must be ended with
:kword:`End of TSConstraints`
(see :numref:`UG:sec:definition_of_internal_coordinates` below).
.. xmldoc::
%%Keyword: TSConstraints
Specify constraints that will be active during the initial stage of an
optimization with FINDTS. When a transition state region is reached
these constraints will be released. If this keyword is not used,
all constraints will be released.
:kword:`GNRM`
Modify the gradient norm threshold associated with the :kword:`FINDTS` option.
The actual threshold is specified on the subsequent line. The default
value is 0.2.
.. xmldoc::
%%Keyword: GNRM
Modify the gradient norm threshold associated with the FINDTS option.
The actual threshold is specified on the subsequent line. The default
value is 0.2.
.. xmldoc::
:kword:`MEP-search` or :kword:`MEP`
Enable a minimum energy path (MEP) search.
.. xmldoc::
.. xmldoc::
%%Keyword: MEP-search
Enable a minimum energy path (MEP) search.
MEP is a valid synonym.
.. xmldoc:: %%Keyword: MEP
Enable a minimum energy path (MEP) search.
Synonym of MEP-search.
:kword:`IRC`
The keyword is used to perform an intrinsic reaction coordinate (IRC) analysis of a
transition state structure. The analysis will follow the reaction path forward and
backward until the energy increases. The keyword requires that the starting structure be
that of a transition state and that the reaction vector be specified explicitly
(check the keyword :kword:`REACtion vector`) or implicitly if it can be found on :file:`RUNOLD`.
Note that the user should not specify any explicit constraints!
.. xmldoc::
%%Keyword: IRC
The keyword is used to perform an intrinsic reaction coordinate (IRC) analysis of a
transition state structure. The analysis will follow the reaction path forward and
backward until the energy increase. The keyword require that the starting structure is
that of a transition state and that the reaction vector is specified explicitly
(check the keyword "REACtion vector") or implicitly can be found on RUNOLD.
Note that the user should not specify any explicit constraints!
:kword:`NMEP` or :kword:`NIRC`
Maximum number of points to find in a minimum energy path search or intrinsic reaction coordinate analysis.
.. xmldoc::
%%Keyword: NMEP
Maximum number of points to find in a minimum energy path search or intrinsic reaction coordinate analysis.
NIRC is a valid synonym.
%%Keyword: NIRC
Maximum number of points to find in an intrinsic reaction coordinate analysis or minimum energy path search.
Synonym of NMEP.
:kword:`MEPStep` or :kword:`IRCStep`
The keyword is used to specify the step length done in the MEP search or IRC analysis.
The step length can be followed with the unit BOHR or ANGSTROM. The default is 0.1 a.u.
(in normalized mass-weighted coordinates).
.. xmldoc::
%%Keyword: MEPStep
The keyword is used to specify the step length done in the MEP search or IRC analysis.
The step length can be followed with the unit BOHR or ANGSTROM. The default is 0.1 a.u.
(in normalized mass-weighted coordinates).
IRCStep is a valid synonym.
%%Keyword: IRCStep
The keyword is used to specify the step length done in the IRC analysis or MEP search.
The step length can be followed with the unit BOHR or ANGSTROM. The default is 0.1 a.u.
(in normalized mass-weighted coordinates).
Synonym of MEPStep.
:kword:`MEPType` or :kword:`IRCType`
Specifies what kind of constraint will be used for optimizing the points during the MEP search or IRC analysis.
The possibilities are SPHERE, the default, which uses the Sphere constraint (each structure is at a given distance in coordinate space from the reference),
or PLANE which uses the Transverse constraint (each structure is at a given distance from the hyperplane defined by the reference and the path direction).
The reference structure changes at each step, according to the :kword:`MEPAlgorithm` keyword.
.. xmldoc::
%%Keyword: MEPType
Specifies what kind of constraint will be used for optimizing the points during the MEP search or IRC analysis.
The possibilities are SPHERE, the default, which uses the Sphere constraint (each structure is at a given distance in coordinate space from the reference),
or PLANE which uses the Transverse constraint (each structure is at a given distance from the hyperplane defined by the reference and the path direction).
The reference structure changes at each step, according to the MEPAlgorithm keyword.
IRCType is a valid synonym.
%%Keyword: IRCType
Specifies what kind of constraint will be used for optimizing the points during the IRC analysis or MEP search.
The possibilities are SPHERE, the default, which uses the Sphere constraint (each structure is at a given distance in coordinate space from the reference),
or PLANE which uses the Transverse constraint (each structure is at a given distance from the hyperplane defined by the reference and the path direction).
The reference structure changes at each step, according to the IRCAlgorithm keyword.
Synonym of MEPType.
:kword:`MEPAlgorithm` or :kword:`IRCAlgorithm`
Selects the algorithm for a MEP search or IRC analysis.
The possibilities are GS for the González--Schlegel algorithm, the default, or MB for the Müller--Brown algorithm.
.. xmldoc::
%%Keyword: MEPAlgorithm
Selects the algorithm for a MEP search or IRC analysis.
The possibilities are GS for the Gonzalez-Schlegel algorithm, the default, or MB for the Mueller-Brown algorithm.
IRCAlgorithm is a valid synonym.
%%Keyword: IRCAlgorithm
Selects the algorithm for a MEP search or IRC analysis.
The possibilities are GS for the Gonzalez-Schlegel algorithm, the default, or MB for the Mueller-Brown algorithm.
Synonym of MEPAlgorithm.
:kword:`MEPConvergence` or :kword:`IRCConvergence`
Sets the gradient convergence for a MEP search or IRC analysis.
The path will be terminated when the gradient norm at an optimized point is below this threshold.
By default it is the same as the gradient threshold for the normal iterations, specified with :kword:`THRShld`,
it may be necessary to reduce it to follow a path on a very flat surface.
.. xmldoc::
%%Keyword: MEPConvergence
Sets the gradient convergence for a MEP search or IRC analysis.
The path will be terminated when the gradient norm at an optimized point is below this threshold.
By default it is the same as the gradient threshold for the normal iterations, specified with THRShld,
it may be necessary to reduce it to follow a path on a very flat surface.
IRCConvergence is a valid synonym.
%%Keyword: IRCConvergence
Sets the gradient convergence for a MEP search or IRC analysis.
The path will be terminated when the gradient norm at an optimized point is below this threshold.
By default it is the same as the gradient threshold for the normal iterations, specified with THRShld,
Synonym of MEPConvergence.
:kword:`REFErence`
The keyword is followed by a list of the symmetry unique coordinates (in au)
of the origin of the hyper sphere. The default origin is the structure
of the first iteration.
.. xmldoc::
%%Keyword: REFErence
The keyword is followed by a list of the symmetry unique coordinates (in au)
of the origin of the hyper sphere. The default origin is the structure
of the first iteration.
:kword:`GRADient of reference`
The keyword is followed by a list of the gradient vector components. This keyword is
compulsory when using the Transverse kind of constraint. The optimization is performed in
a space orthogonal to the given vector.
.. xmldoc::
%%Keyword: GRADient of reference
The keyword is followed by a list of the gradient vector components. This keyword is
compulsory when using the Transverse kind of constraint. The optimization is performed in
a space orthogonal to the given vector.
:kword:`REACtion vector`
The keyword is followed by the reaction vector specified as the Cartesian vector components
on each of the symmetry unique atoms.
.. xmldoc::
%%Keyword: REACtion vector
The keyword is followed by the reaction vector specified as the Cartesian vector components
on each of the symmetry unique atoms.
.. xmldoc::
Optional force constant keywords
.. class:: keywordlist
:kword:`OLDForce`
The Hessian matrix is read from the file :file:`RUNOLD`.
This Hessian is either
an analytic or approximative Hessian updated by Slapaf.
Note that for this option to work properly the type of
internal coordinates must be the same!
.. xmldoc::
.. xmldoc::
%%Keyword: Oldforce
The Hessian matrix is read from the file RUNOLD.
This Hessian is either
an analytic or approximative Hessian updated by Slapaf.
Note that for this option to work properly the type of
internal coordinates must be the same!
:kword:`FCONstant`
Input of Hessian in internal coordinates.
There are two different syntaxes.
#. The keyword is followed by an entry with
the number of elements which will be set (observe that the
update will preserve that the elements :math:`H_{ij}` and :math:`H_{ji}` are
equal). The next entries will contain the value and the indices of
the elements to be replaced.
#. The keyword if followed by the label :kword:`Square` or
:kword:`Triangular`. The subsequent line specifies the rank of the
Hessian. This is then followed by entries specifying the Hessian
in square or lower triangular order.
.. xmldoc:: %%Keyword: Fconstant
Input of Hessian in internal coordinates.
Note this is
There are two different syntaxes.
1) The keyword is followed by an entry with
the number of elements which will be set (observe that the
update will preserve that the elements Hij and Hji are
equal). The next lines will contain the value and the indices of
the elements to be replaced.
2) The keyword if followed by the label "Square" or
"Triangular". The subsequent entry specifies the rank of the
Hessian. This is then followed by entries specifying the Hessian
in square or lower triangular order.
:kword:`XFCOnstant`
Input of an external Hessian matrix in cartesian coordinates. The
syntax is the same as for the :kword:`FCONSTANT` keyword.
.. xmldoc::
%%Keyword: XFConstant
Input of an external Hessian matrix in cartesian coordinates. The
syntax is the same as for the FCONSTANT keyword.
:kword:`NUMErical`
This invokes as calculation of the force constant matrix by a
two-point finite difference formula. The resulting force
constant matrix is used for an analysis of the harmonic
frequencies. **Observe** that in case of the use of internal
coordinates defined as Cartesian coordinates that these has to be
linear combinations which are free from translational and
rotational components for the harmonic frequency analysis to be
valid. **Alternative:** see keyword :kword:`RowH` in the section
about Internal coordinates.
.. xmldoc::
%%Keyword: Numerical
This invokes as calculation of the force constant matrix by a
two-point finite difference formula. The resulting force
constant matrix is used for an analysis of the harmonic
frequencies. Observe that in case of the use of internal
coordinates defined as Cartesian coordinates that these has to be
linear combinations which are free from translational and
rotational components for the harmonic frequency analysis to be
valid.
:kword:`CUBIc`
This invokes a calculation of the 2nd and the 3rd order
force constant matrix by finite difference formula.
.. xmldoc::
%%Keyword: Cubic
This invokes a calculation of the 2nd and the 3rd order
force constant matrix by finite difference formula.
:kword:`DELTa`
This keyword is followed by a real number which defines the
step length used in the finite differentiation. Default: 1.0D-2.
.. xmldoc::
%%Keyword: Delta
This keyword is followed by a real number which defines the
step length used in the finite differentiation. Default: 1.0D-2.
:kword:`PRFC`
The eigenvalues and eigenvectors of the Hessian matrix
are printed. The internal coordinates definitions are also printed.
.. xmldoc::
%%Keyword: PrFC
The eigenvalues and eigenvectors of the Hessian matrix
are printed. The internal coordinates definitions is also printed.
:kword:`RHIDden`
Define the hidden atoms selection radius in order to improve a QM/MM Hessian. It can be followed by :kword:`Angstrom`.
.. xmldoc::
%%Keyword: rHid
Define the hidden atoms selection radius in order to improve a QM/MM Hessian.
.. xmldoc::
Optional miscellaneous keywords
.. class:: keywordlist
:kword:`CTOF`
Coordinates TO Follow defines an internal coordinate whose values
will be printed in the output during the optimization. Both
the original and the new values will be printed.
The keyword must be followed by the definition on the primitive
coordinate.
.. xmldoc::
%%Keyword: CTOF
Coordinates TO Follow defines an internal coordinate whose values
will be printed in the output during the optimization. Both
the original and the new values will be printed.
The keyword must be followed by the definition on the primitive
coordinate.
:kword:`RTRN`
Maximum number of atoms for which bond lengths, angles and dihedral
angles are listed, and
the radius defining the maximum length of a bond follows.
The latter is used as a threshold when printing out
angles and dihedral angles. The length can be followed by
:kword:`Bohr` or
:kword:`Angstrom` which indicates the unit in which the length
was specified, the default is
:kword:`Bohr`.
The default values are 15 and 3.0 au.
.. xmldoc::
%%Keyword: RTRN
Maximum number of atoms for which bond lengths, angles and dihedral
angles are listed, and
the radius defining the maximum length of a bond follows on
the next line. The latter is used as a threshold when printing out
angles and dihedral angles. The length can be followed by
"Bohr" or "Angstrom" which indicates the unit in which the length
was specified, the default is "Bohr".
:kword:`THERmochemistry`
Request frequencies to be computed followed by an user specified thermochemical analysis.
The keyword must be followed by different entries containing the Rotational Symmetry Number,
the Pressure (in atm), and one entry per Temperature (in K)
for which the thermochemistry will be calculated.
The section is ended by the keyword :kword:`End of PT`.
.. xmldoc::
%%Keyword: THER
Request frequencies to be computed followed by an user specified thermochemical analysis.
The keyword must be followed by different entries containing the Rotational Symmetry Number,
the Pressure (in atm), and one entry per Temperature (in K)
for which the thermochemistry will be calculated.
The section is ended by the keyword "End of PT".
:kword:`DISOtope`
Calculates frequencies modified for double isotopic substitution.
.. xmldoc::
%%Keyword: DISOtope
Calculates frequencies modified for double isotopic substitution.
:kword:`TRACk`
Tries to follow electronic states during an optimization, by computing state overlaps with :program:`RASSI`
at each step. Root numbers selected with :kword:`RlxRoot` in :program:`RASSCF` or with the "EDiff" constraint
are only fixed in the first iteration, then the best-matching states are chosen.
.. xmldoc::
%%Keyword: Track
Tries to follow electronic states during an optimization, by computing state overlaps with RASSI.
:kword:`LASTenergy`
Specifies the quantum chemical method requested for the Last_Energy module (e.g., SCF, CASSCF, CASPT2, etc.)
The keyword must be followed by the name of the module. Moreover, the EMIL command COPY needs to be used
in the global input to provide a file named LASTEN, containing the input for the specified module.
.. xmldoc::
%%Keyword: LAST
Specifies the quantum chemical method requested for the Last_Energy module (e.g., SCF, CASSCF, CASPT2, etc.)
The keyword must be followed by the name of the module. Moreover, the EMIL command COPY needs to be used
in the global input to provide a file named LASTEN, containing the input for the specified module.
:kword:`NOLAst energy`
Disables the call to the :program:`Last_Energy` module when convergence is achieved.
.. xmldoc::
%%Keyword: NoLastEnergy
Disables the call to the Last_Energy module when convergence is achieved.
Optional restricted variance optimization (RVO) :cite:`Raggi2020` keywords
.. class:: keywordlist
:kword:`KRIGing`
Activate RVO using gradient-enhanced Kriging (GEK) to describe the surrogate model.
The maximum number of sample points (energies and gradients) is twice the value indicated by the :kword:`WINDow` keyword (i.e. 10 by default).
.. xmldoc::
%%Keyword: Kriging
Activate restricted variance optimization (RVO) using gradient-enhanced Kriging to describe the surrogate model.
:kword:`TFOFfset`
Trend function or baseline offset for the GEK surrogate model.
The surrogate model will tend to the maximum energy among the sample points plus this value (in au).
The default value is 10.0 au.
.. xmldoc::
%%Keyword: TFOFfset
Trend function or baseline offset for the GEK surrogate model.
The surrogate model will tend to the maximum energy among the sample points plus this value.
Default: 10.0 au.
:kword:`MAXDisp`
Maximum energy dispersion allowed during each macro iteration of the RVO procedure.
A real value is read from the input, the maximum dispersion is this value times the maximum Cartesian gradient.
The default value is 0.3 au.
.. xmldoc::
%%Keyword: MAXDISP
Maximum energy dispersion allowed during each macro iteration of the RVO procedure.
A factor, multiplied by the maximum Cartesian Gradient.
Default: 0.3 au.
:kword:`MXMI`
Maximum number of micro iterations in each macro iteration of the RVO procedure.
The default value is 50.
.. xmldoc::
%%Keyword: MXMI
Maximum number of micro iterations in each macro iteration of the RVO procedure.
Default: 50.
Example: A complete set of input decks for a CASSCF geometry
optimization. These are the input decks for the optimization
of the enediyne molecule.
.. extractfile:: ug/SLAPAF.input
&GATEWAY
Title= Enediyne
Coord= $MOLCAS/Coord/enediyne.xyz
Basis= ANO-L-VQZP
Group= x z
> DoWhile
&SEWARD
&SCF
ITERATIONS= 30; Occupied= 9 8 2 1; Thresholds= 1.0d-8 1.0d-3 1.5d-3 0.2d-3; IVO
&RASSCF
Symmetry= 1; Spin= 1
NactEl= 12 0 0; Inactive= 7 7 0 0; Ras2= 3 3 3 3
Iterations= 50 50; CiRoot= 1 1; 1; Thrs= 1.0e-08 1.0e-05 1.0e-05
Lumorb
&SLAPAF; Iterations= 20
> EndDo
Example: Thermochemistry for an asymmetric top (Rotational Symmetry Number
= 1), at 1.0 atm and 273.15, 298.15, 398.15 and 498.15 K. ::
&SLAPAF; THERmochemistry= 1; 1.0; 273.15; 298.15; 398.15; 498.15; End of PT
End of input
.. _UG\:sec\:definition_of_internal_coordinates:
Definition of internal coordinates or constraints
.................................................
The input section defining the internal coordinates always start with the
keyword :kword:`Internal coordinates` and the definition of the constraints
starts with the keyword :kword:`Constraints`. Note that the latter
is an input section for the :program:`GATEWAY` module.
The input is always sectioned into two
parts where the first section defines a set of primitive internal
coordinates
and the second part defines the actual internal coordinates as
any arbitrary linear combination of the primitive internal coordinates
that was defined in the first section.
In case of constraints the second part does also assign values to the
constraints.
In the first section we will refer to the atoms by their atom label
(:program:`SEWARD` will make sure that there is no redundancy). In case of
symmetry one will have to augment the atom label with a symmetry operation
in parenthesis in order to specify a symmetry related center.
Note that the user only
have to specify distinct internal coordinates (:program:`ALASKA` will make the
symmetry adaptation).
In the specification below *rLabel* is a user defined label with no more
than 8 (eight) characters. The specifications atom1, atom2, atom3, and atom4
are the unique atom labels as specified in the input to :program:`SEWARD`.
**The primitive internal coordinates** are defined as
.. class:: primlist
*rLabel* = Bond atom1 atom2
a primitive internal coordinate *rLabel* is defined as the bond
between center atom1 and atom2.
*rLabel* = Angle atom1 atom2 atom3
a primitive internal coordinate *rLabel* is defined as the angle
between the bonds formed from connecting atom1 to atom2 and
connecting atom2 to atom3.
*rLabel* = LAngle(1) atom1 atom2 atom3
a primitive internal coordinate *rLabel* is defined as the linear angle
between the bonds formed from connecting atom1 to atom2 and
connecting atom2 to atom3. To define the direction of the angle the following
procedure is followed.
#. --- *the three centers are linear*,
#. form a reference axis, :math:`R_1`, connecting atom1 and atom3,
#. compute the number of zero elements, *nR*, in the reference vector,
#. --- *nR=0*,
a first perpendicular direction to the reference axis is formed by
.. compound::
.. math:: R=(R_{1x},R_{1y},-R_{1z})
followed by the projection
.. math:: R_2=R-\frac{R \cdot R_1}{R_1 \cdot R_1} R_1.
The second perpendicular direction completes the right-handed system.
#. --- *nR=1*,
a first perpendicular direction to the reference axis is defined by setting the element in :math:`R_2`
corresponding to the zero entry in :math:`R_1` to unity.
The second perpendicular direction completes the right-handed system.
#. --- *nR=2*,
a first perpendicular direction to the reference axis is defined by setting the element
corresponding to the first zero entry in :math:`R_1` to unity.
The second perpendicular direction completes the right-handed system.
#. --- *the three centers are nonlinear*,
the first perpendicular direction is the one which is in the plane formed by atoms atom1, atom2, and atom3.
The second perpendicular direction is taken as the direction perpendicular to the same plane.
The direction of the bend for **LAngle(1)** is taken in the direction of the first perpendicular direction, etc.
*rLabel* = LAngle(2) atom1 atom2 atom3
a primitive internal coordinate *rLabel* is defined as the linear angle
between the bonds formed from connecting atom1 to atom2 and
connecting atom2 to atom3. The definition of the perpendicular directions
is as described above. The direction of the bend for **LAngle(2)** is taken in the direction of
the second perpendicular direction.
*rLabel* = Dihedral atom1 atom2 atom3 atom4
a primitive internal coordinate *rLabel* is defined as the angle
between the planes formed of atom1, atom2 and atom3, and atom2,
atom3 and atom4, respectively.
*rLabel* = OutOfP atom1 atom2 atom3 atom4
a primitive internal coordinate *rLabel* is defined as the angle
between the plane formed by atom2, atom3, and atom4 and the
bond formed by connecting atom1 and atom4.
*rLabel* = Dissoc (n1+n2) atom1 atom2 atom3 ... atomN
a primitive internal coordinate *rLabel* is defined as the distance
between the center of masses of two sets of centers. The first
center has n1 members and the second has n2.
The input contains the labels of the atoms of the first group followed
immediately by the labels of the second group.
This option is not available for constraints.
*rLabel* = Cartesian i atom1
a primitive internal coordinate *rLabel* is defined as the pure
Cartesian displacement of the center labeled atom1. The label
i is selected to x, y, or z to give the appropriate component.
*rLabel* = Ediff [i j]
the energy difference between states i and j (if provided, the brackets indicate they
are optional, do not include the brackets).
If i and j are not provided, the difference is between the "current" state and
the state provided on :file:`RUNFILE2`.
This is only used in constrained optimization in which crossings or conical intersections
are located. If this constraint is used, the average energy of the two states will
be optimized, subject to the constraint. If the value is 0.0 and the spin and spatial
symmetry of both states is the same, a conical intersection will be searched.
In this case, :program:`SLAPAF` will request an analytical calculation of the nonadiabatic
coupling vector, if available. If it is not available, or if :file:`RUNFILE2` is being used
(i and j not provided), the branching space update method of Maeda et al. will be used :cite:`Maeda2010`.
*rLabel* = Sphere
the radius of the hypersphere defined by two different molecular structures
(the origin is the first structure) in relative mass-weighted coordinates.
This is only used in constrained optimization in which minimum reaction paths (MEP) or intrinsic reaction
coordinate (IRC) paths are followed. The units of the radius is in mass-weighted coordinates
divided with the square root of the total mass of the molecule.
*rLabel* = Transverse
a level of "orthogonality". This is used to perform an optimization in a space
orthogonal to a given vector. Recommended value 0.0. Requires usage of GRAD keyword.
*rLabel* = Fragment atom1 atom2 atom3 ... atomN
a dummy internal coordinate *rLabel* is defined. This translates to
that a set of internal coordinates are generated automatically according
to a standard Z-matrix format to define all degrees of freedom
of the fragment defined by the list of atoms on the same line. These
internal coordinates will be automatically fixed in the geometry optimizations to
the values of starting structure. Note, the values of these do not need to
be explicitly defined and set in the :kword:`Values` section. Note, too, that
the generation of the internal coordinates is done according to the order
in which atom1, atom2, etc. are given; for some systems, especially with
linear angles, it may be preferable to define the coordinates manually.
The second section starts with the label :kword:`Vary` or in the case of constraints
with the label :kword:`Values`.
.. compound::
In case of a definition of **internal coordinates** in this section the user
specifies all symmetric internal coordinates excluding translation and rotation
using a list of expressions like
*label* = f1 *rLabel1* + f2 *rLabel2* + ...
which defines an internal coordinate *label* as the linear combination of the
primitive internal coordinates *rLabel1*, *rLabel2*, ... with the coefficients
f1, f2, ..., respectively. If the internal coordinate just corresponds to
the primitive internal coordinate, the same label can be used
*label*
If some internal coordinates are chosen to be fixed they should be defined after
the label :kword:`Fix`. The fixed internal coordinate are defined with
expressions as in the section :kword:`Vary`. Observe: using expression can
introduce linear dependence and/or undefined nuclear coordinates, so use with care.
For the internal coordinates defined after :kword:`Vary` (and :kword:`Fix`, if present)
a numerical estimation of rows and columns of the Hessian matrix can be performed. The
*label* of internal coordinates (max 10) must be specified after keyword :kword:`RowH`.
Keywords :kword:`NUMErical` and :kword:`RowH` are mutually exclusive.
.. compound::
In case of a definition of **constraints** the sections contains either a
direct reference to a *rLabel* as in
*rLabel* = *rValue* [Angstrom,Degrees] [Soft,Hard] [Phantom]
or one can also use expressions like
f1 *rLabel1* |+-| f2 *rLabel2* |+-| ... = *Value* [Angstrom,Degrees] [Soft,Hard] [Phantom]
where *rValue* is the desired value of the constraint in au or rad, or in
angstrom or degrees if the corresponding keyword is added. The "Hard" and "Soft"
keywords are only meaningful for numerical differentiation: the coordinates corresponding
to soft constraints are differentiated, those of hard constraints are not :cite:`Stenrup2015`.
By default almost all constraints are hard, only constraints of the type "Sphere", "Transverse"
and "Ediff" default to soft. The "Hard" and "Soft" keywords override the default.
When using constraints in combination with the :kword:`FINDTS` keyword, one should use
soft constraints, at least for the constraint most similar to the expected reaction vector.
Constraints defined in :kword:`TSCOnstraints` (recommended) are automatically considered
soft.
The "Phantom" modifier can be used to ignore a constraint in the optimization. A phantom
constraint will only be considered for numerical differentiation. Phantom constraints are
useful in combination with the :kword:`KEEPOldGradient` keyword of :program:`ALASKA`.
Using :kword:`NGEXclude` in :program:`GATEWAY` is equivalent to phantom constraints,
and it is the preferred way to set up composite gradients :cite:`Stenrup2015`.
Alternatively, if the current value of an internal coordinate is to be used, i.e.
no change is to be allowed (frozen), this is expressed as
*rLabel* = Fix [Soft,Hard] [Phantom]
Note that a coordinate of type "Fragment" does not need to appear in the :kword:`Values`
section, but if it does it must be assigned the value "Fix".
Example: A definition of user specified internal coordinates of benzene. The molecule is
in :math:`D_{6h}` and since |molcas| only uses up to :math:`D_{2h}` the
:kword:`Fix` option is used to
constrain the relaxation to the higher point group. **Observe** that this will
only restrict the nuclear coordinates to :math:`D_{6h}`. The electronic wavefunction,
however, can have lower symmetry. ::
Internal coordinates
r1 = Bond C1 C2
r2 = Bond C1 H1
r3 = Bond C2 H2
r4 = Bond C2 C2(x)
f1 = Angle H1 C1 C2
f2 = Angle H2 C2 C1
Vary
a = 1.0 r1 + 1.0 r4
b = 1.0 r2 + 1.0 r3
c = 1.0 f1 + 1.0 f2
Fix
a = 1.0 r1 + -1.0 r4
b = 1.0 r2 + -1.0 r3
c = 1.0 f1 + -1.0 f2
End of Internal
Example: A input for the optimization of water constraining the structure to be linear
at convergence.
.. extractfile:: ug/SLAPAF.constrains.input
&GATEWAY
Title= H2O geom optim, using the ANO-S basis set.
Coord=$MOLCAS/Coord/Water.xyz
Basis=ANO-S-VDZ
Group= c1
Constraints
a1 = langle(1) H2 O1 H3
Values
a1 = 179.99 degrees
End of Constraints
>>> DO WHILE <<<
&SEWARD; &SCF
&SLAPAF
>>> END DO <<<
Example: A complete set of input decks for a UHF transition
structure geometry optimization of an identity hydrogen
transfer reaction (:math:`\ce{HO + H_2O -> H_2O + OH}`).
.. extractfile:: ug/SLAPAF.Zmat.input
&GATEWAY
ZMAT
O.STO-3G....
H.STO-3G....
H1
Z2 1 1.0
O3 1 1.15 2 92.
O4 1 1.15 2 92. 3 180.
H5 3 0.98 4 105.4 2 120.
H6 4 0.98 3 105.4 2 120.
>>> DO WHILE <<<
&SEWARD;
&SCF; UHF
&SLAPAF; TS; PRFC
Internal
bOO4 = Bond O3 O4
bOH5 = Bond H5 O3
bOH6 = Bond H6 O4
bOH1 = Bond O3 H1
aOOH5 = Angle O4 O3 H5
aOOH6 = Angle O3 O4 H6
aHOH1 = Angle H5 O3 H1
dH6 = Dihedral H6 O4 O3 H5
dH1 = Dihedral O4 H5 O3 H1
Vary; bOO4; bOH5; bOH6; bOH1; aOOH5; aOOH6; aHOH1; dH6; dH1
RowH; bOH1
End of Internal
>>> ENDDO <<<
Example: Optimization of a minimum energy conical intersection point,
using automatic calculation of analytical gradients and nonadiabatic coupling.
.. extractfile:: ug/SLAPAF.CI.input
&GATEWAY
Coord = acrolein.xyz
Basis = cc-pVDZ
Group = NoSymm
Constraints
a = Ediff 1 2
Values
a = 0.0
End of constraints
>>> DoWhile
&SEWARD
>>> If (iter = 1)
&SCF
&MBPT2
PrPt
>>> EndIf
&RASSCF
FileOrb = $Project.MP2Orb
Charge = 0
NActEl = 6 0 0
RAS2 = 5
CIRoot = 4 4 1
&SLAPAF
>>> EndDo
.. xmldoc::
.. xmldoc::
.. xmldoc::
.. xmldoc::
.. xmldoc::