.. index::
   single: Program; SLAPAF
   single: SLAPAF

.. _UG\:sec\:slapaf:

:program:`SlapAf`
=================

.. only:: html

  .. contents::
     :local:
     :backlinks: none

.. xmldoc:: <MODULE NAME="SLAPAF">
            %%Description:
            <HELP>
            This program is a general purpose facility for geometry
            optimization. At present, it is tailored to use analytical
            or numerical gradients produced by ALASKA.
            SLAPAF also computes an approximate Hessian.
            </HELP>

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 MODULE="SLAPAF" NAME="ITERATIONS" APPEAR="Max iterations" KIND="INT" MIN_VALUE="0" DEFAULT_VALUE="500" LEVEL="BASIC">
              %%Keyword: Iterations <basic>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="THRSHLD" APPEAR="Convergence threshold" KIND="REALS" SIZE="2" MIN_VALUE="0.0" DEFAULT_VALUES="0.0,3.0D-4" LEVEL="ADVANCED">
              %%Keyword: Thrshld <advanced>
              <HELP>
              Enter two real numbers
              which specifies the convergence criterion with respect to the
              energy change and the norm of the gradient.
              </HELP>
              The defaults are 0.0 and 3.0D-4 for Gaussian, and 1.0D-6 and 3.0D-4 for Baker.
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="BAKER" APPEAR="Baker style convergence criteria" KIND="SINGLE" LEVEL="ADVANCED">
              %%Keyword: Baker <advanced>
              <HELP>
              Activate convergence criteria according to Baker.
              Default is to use the convergence criteria as in the Gaussian
              program.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="MAXSTEP" APPEAR="Max step" KIND="REAL" MIN_VALUE="0.0" DEFAULT_VALUE="0.3" LEVEL="BASIC">
              %%Keyword: Maxstep <basic>
              <HELP>
              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.
              </HELP>
              The default
              value is 0.3 au or rad.
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="CNWEIGHT" APPEAR="Max constraint weight" KIND="REAL" MIN_VALUE="0.0" DEFAULT_VALUE="1.0" LEVEL="BASIC">
              %%Keyword: CnWeight <basic>
              <HELP>
              Sets the maximum weight assigned to the fulfillment of
              the constraints, relative to the step taken in the
              complementary space for energy minimization.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="TOLE" APPEAR="Constraint tolerance" KIND="REAL" MIN_VALUE="0.0" DEFAULT_VALUE="1.0D10" LEVEL="BASIC">
              %%Keyword: Tolerance <basic>
              <HELP>
              Controls how strictly the constraints must be satisfied at convergence.
              </HELP>
              </KEYWORD>

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:: <GROUP MODULE="SLAPAF" KIND="BOX" NAME="IC" APPEAR="Internal coordinate selection" LEVEL="ADVANCED">

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="CARTESIAN" APPEAR="Cartesian internal coordinates" KIND="SINGLE" EXCLUSIVE="HWRS,NOHWRS,INTERNAL" LEVEL="ADVANCED">
              %%Keyword: Cartesian <advanced>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="INTERNAL" APPEAR="User-defined internal coordinates" KIND="CUSTOM" EXCLUSIVE="HWRS,NOHWRS,CARTESIAN" LEVEL="ADVANCED">
              %%Keyword: Internal <advanced>
              <HELP>
              This marks the start of the definition of the internal
              coordinates.
              </HELP>
              This section is always ended by the keyword "End of Internal".
              Consult the manual for details.
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="HWRS" APPEAR="FWC internal coordinates" KIND="SINGLE" EXCLUSIVE="NOHWRS,CARTESIAN,INTERNAL" LEVEL="ADVANCED">
              %%Keyword: HWRS <basic>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="NOHWRS" APPEAR="Integral coordinates" KIND="SINGLE" EXCLUSIVE="HWRS,CARTESIAN,INTERNAL" LEVEL="ADVANCED">
              %%Keyword: NoHWRS <basic>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="FUZZ" APPEAR="Fuzz" KIND="REAL" MIN_VALUE="0.001" DEFAULT_VALUE="0.5" LEVEL="ADVANCED">
              %%Keyword: Fuzz <advanced>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

  .. xmldoc:: </GROUP>

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:: <GROUP MODULE="SLAPAF" KIND="BOX" NAME="HU" APPEAR="Hessian update options" LEVEL="ADVANCED">

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="HUPDATE" APPEAR="Hessian update method" KIND="CHOICE" LIST="None,BFGS,MSP,EU,TS-BFGS" LEVEL="ADVANCED" DEFAULT_VALUE="BFGS">
              %%Keyword: HUpdate <advanced>
              <HELP>
              Method used for updating the Hessian matrix.
              </HELP>
              It must be one of: None, BFGS, MSP, EU, TS-BFGS.
              </KEYWORD>

:kword:`UORDer`
  Order the gradients and displacements vectors according to Schlegel prior to
  the update of the Hessian. Default is no reorder.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="UORDER" APPEAR="Ordered update" KIND="SINGLE" EXCLUSIVE="NOUPDATE" LEVEL="ADVANCED">
              %%Keyword: UORDer <basic>
              <HELP>
              Order the gradients and displacements vectors according to Schlegel prior to
              the update of the Hessian. Default is no reorder.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="WINDOW" APPEAR="Update window" KIND="INT" EXCLUSIVE="NOUPDATE" DEFAULT_VALUE="5" LEVEL="ADVANCED">
              %%Keyword: WINDow <basic>
              <HELP>
              Maximum number of previous iterations to include in the Hessian update.
              </HELP>
              Default is 5.
              </KEYWORD>

  .. xmldoc:: </GROUP>

Optional optimization procedure keywords

.. class:: keywordlist

:kword:`NOLIne`
  Disable line search. Default is to use line search for minima.

  .. xmldoc:: <GROUP MODULE="SLAPAF" KIND="BOX" NAME="OP" APPEAR="Optimization options" LEVEL="BASIC">

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="NOLINE" APPEAR="Deactivate line-search" KIND="SINGLE" LEVEL="BASIC">
              %%Keyword: Noline <basic>
              <HELP>
              Disable line search. Default is to use line search for minima.
              </HELP>
              </KEYWORD>

:kword:`RATIonal`
  Activate geometry optimization using the restricted step Rational Functional optimization :cite:`rf,rs-rf`,
  this is the default.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="RATIONAL" APPEAR="RFO option" KIND="SINGLE" EXCLUSIVE="C1-DIIS,C2-DIIS,NEWTON" LEVEL="BASIC">
              %%Keyword: Rational <basic>
              <HELP>
              Activate geometry optimization using the restricted step Rational Functional optimization,
              this is the default.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="C1-DIIS" APPEAR="C1-DIIS option" KIND="SINGLE" EXCLUSIVE="RATIONAL,C2-DIIS,NEWTON" LEVEL="BASIC">
              %%Keyword: C1-diis <advanced>
              <HELP>
              Activate geometry optimization using the C1-GDIIS method.
              The default is to use the Rational Functional approach.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="C2-DIIS" APPEAR="C2-DIIS option" KIND="SINGLE" EXCLUSIVE="RATIONAL,C1-DIIS,NEWTON" LEVEL="BASIC">
              %%Keyword: C2-diis <basic>
              <HELP>
              Activate geometry optimization using the C2-GDIIS method.
              The default is to use the Rational Functional approach.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="DXDX" APPEAR="dxdx DIIS option" KIND="SINGLE" REQUIRE="C1-DIIS.OR.C2-DIIS" EXCLUSIVE="DXG,GDX,GG" LEVEL="ADVANCED">
              %%Keyword: dxdx <basic>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="DXG" APPEAR="dxg DIIS option" KIND="SINGLE" REQUIRE="C1-DIIS.OR.C2-DIIS" EXCLUSIVE="DXDX,GDX,GG" LEVEL="ADVANCED">
              %%Keyword: dxg <basic>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

:kword:`GDX`
  See above.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="GDX" APPEAR="gdx DIIS option" KIND="SINGLE" REQUIRE="C1-DIIS.OR.C2-DIIS" EXCLUSIVE="DXDX,DXG,GG" LEVEL="ADVANCED">
              %%Keyword: gdx <basic>
              <HELP>
              See the dxg keyword.
              </HELP>
              </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 MODULE="SLAPAF" NAME="GG" APPEAR="gg DIIS option" KIND="SINGLE" REQUIRE="C1-DIIS.OR.C2-DIIS" EXCLUSIVE="DXDX,DXG,GDX" LEVEL="ADVANCED">
              %%Keyword: gg <basic>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

:kword:`NEWTon`
  Activate geometry optimization using the standard quasi-Newton approach.
  The default is to use the Rational Functional approach.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="NEWTON" APPEAR="Newton-Raphson optimization" KIND="SINGLE" EXCLUSIVE="RATIONAL,C1-DIIS,C2-DIIS" LEVEL="BASIC">
              %%Keyword: Newton <basic>
              <HELP>
              Activate geometry optimization using the standard quasi-Newton approach.
              The default is to use the Rational Functional approach.
              </HELP>
              </KEYWORD>

:kword:`RS-P-rfo`
  Activate RS-P-RFO :cite:`rs-rf` as default for TS-search. Default is RS-I-RFO.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="RS-P-RFO" APPEAR="RS-P-RFO option" KIND="SINGLE" REQUIRE="TS.OR.FINDTS" LEVEL="BASIC">
              %%Keyword: RS-P-RFO <basic>
              <HELP>
              Activate RS-P-RFO as default for TS-search. Default is RS-I-RFO.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="TS" APPEAR="Brute force TS optimization" KIND="SINGLE" EXCLUSIVE="C1-DIIS,C2-DIIS,NEWTON,NOUPDATE,MEYER,BPUPDATE,BFGS" LEVEL="BASIC">
              %%Keyword: TS <basic>
              <HELP>
              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).
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="MODE" APPEAR="Mode specification" KIND="INT" REQUIRE="TS" LEVEL="BASIC">
              %%Keyword: Mode <advanced>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="FINDTS" APPEAR="Find TS algorithm" KIND="SINGLE" EXCLUSIVE="TS" LEVEL="BASIC">
              %%Keyword: FindTS <basic>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="TSCONSTRAINTS" APPEAR="TS constraints" KIND="CUSTOM" REQUIRE="FINDTS" LEVEL="BASIC">
              %%Keyword: TSConstraints <basic>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="GNRM" APPEAR="Gradient norm threshold" KIND="REAL" MIN_VALUE="0.0" DEFAULT_VALUE="0.2" REQUIRE="FINDTS" LEVEL="ADVANCED">
              %%Keyword: GNRM <basic>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

  .. xmldoc:: </GROUP>

:kword:`MEP-search` or :kword:`MEP`
  Enable a minimum energy path (MEP) search.

  .. xmldoc:: <GROUP MODULE="SLAPAF" KIND="BOX" NAME="ADVANCED" APPEAR="Advanced PES exploration options" LEVEL="BASIC">

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="MEP-SEARCH" APPEAR="MEP-search" KIND="SINGLE" LEVEL="BASIC" ALSO="MEP">
              %%Keyword: MEP-search <basic>
              <HELP>
              Enable a minimum energy path (MEP) search.
              </HELP>
              MEP is a valid synonym.
              </KEYWORD>

  .. xmldoc:: %%Keyword: MEP <basic>
              Enable a minimum energy path (MEP) search.
              Synonym of MEP-search.

:kword:`rMEP-search` 
  Enable a reverse minimum energy path (MEP) search.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="RMEP-SEARCH" APPEAR="Reverse MEP-search" KIND="SINGLE" LEVEL="BASIC">
              %%Keyword: RMEP-search <basic>
              <HELP>
              Enable a reverse minimum energy path (MEP) search.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="IRC" APPEAR="IRC analysis" KIND="SINGLE" LEVEL="BASIC">
              %%Keyword: IRC <basic>
              <HELP>
              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!
              </HELP>
              </KEYWORD>

:kword:`NMEP` or :kword:`NIRC`
  Maximum number of points to find in a minimum energy path search or intrinsic reaction coordinate analysis.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="NMEP" APPEAR="Max points on a MEP" KIND="INT" MIN_VALUE="1" REQUIRE="MEP-SEARCH.OR.RMEP-SEARCH.OR.IRC" LEVEL="BASIC" ALSO="NIRC">
              %%Keyword: NMEP <basic>
              <HELP>
              Maximum number of points to find in a minimum energy path search or intrinsic reaction coordinate analysis.
              </HELP>
              NIRC is a valid synonym.
              </KEYWORD>

              %%Keyword: NIRC <basic>
              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 MODULE="SLAPAF" NAME="MEPSTEP" APPEAR="MEP Step" KIND="STRING" REQUIRE="MEP-SEARCH.OR.RMEP-SEARCH.OR.IRC" LEVEL="BASIC" ALSO="IRCSTEP">
              %%Keyword: MEPStep <basic>
              <HELP>
              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).
              </HELP>
              IRCStep is a valid synonym.
              </KEYWORD>

              %%Keyword: IRCStep <basic>
              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 MODULE="SLAPAF" NAME="MEPTYPE" APPEAR="MEP Type" KIND="STRING" REQUIRE="MEP-SEARCH.OR.RMEP-SEARCH.OR.IRC" LEVEL="ADVANCED" ALSO="IRCTYPE">
              %%Keyword: MEPType <advanced>
              <HELP>
              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.
              </HELP>
              IRCType is a valid synonym.
              </KEYWORD>

              %%Keyword: IRCType <advanced>
              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 MODULE="SLAPAF" NAME="MEPALGORITHM" APPEAR="MEP Algorithm" KIND="STRING" REQUIRE="MEP-SEARCH.OR.RMEP-SEARCH.OR.IRC" LEVEL="ADVANCED" ALSO="IRCALGORITHM">
              %%Keyword: MEPAlgorithm <advanced>
              <HELP>
              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.
              </HELP>
              IRCAlgorithm is a valid synonym.
              </KEYWORD>

              %%Keyword: IRCAlgorithm <advanced>
              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 MODULE="SLAPAF" NAME="MEPCONVERGENCE" APPEAR="MEP gradient convergence" KIND="REAL" REQUIRE="MEP-SEARCH.OR.RMEP-SEARCH.OR.IRC" LEVEL="ADVANCED" ALSO="IRCCONVERGENCE">
              %%Keyword: MEPConvergence <advanced>
              <HELP>
              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.
              </HELP>
              IRCConvergence is a valid synonym.
              </KEYWORD>

              %%Keyword: IRCConvergence <advanced>
              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 MODULE="SLAPAF" NAME="REFERENCE" APPEAR="MEP reference structure" KIND="REALS_LOOKUP" SIZE="DEG_FREEDOM" LEVEL="BASIC">
              %%Keyword: REFErence <basic>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="GRAD" APPEAR="Gradient of reference" KIND="REALS_LOOKUP" SIZE="DEG_FREEDOM" LEVEL="BASIC">
              %%Keyword: GRADient of reference <basic>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="REAC" APPEAR="Reaction vector" KIND="REALS_LOOKUP" SIZE="DEG_FREEDOM" LEVEL="BASIC">
              %%Keyword: REACtion vector <basic>
              <HELP>
              The keyword is followed by the reaction vector specified as the Cartesian vector components
              on each of the symmetry unique atoms.
              </HELP>
              </KEYWORD>

  .. xmldoc:: </GROUP>

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:: <GROUP MODULE="SLAPAF" KIND="BOX" NAME="FCO" APPEAR="Force constant options" LEVEL="ADVANCED">

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="OLDFORCE" APPEAR="External Hessian in internals" KIND="SINGLE" LEVEL="ADVANCED">
              %%Keyword: Oldforce <advanced>
              <HELP>
              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!
              </HELP>
              </KEYWORD>

: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 <advanced>
              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 MODULE="SLAPAF" NAME="XFCONSTANT" APPEAR="External Hessian in Cartesians" KIND="SINGLE" LEVEL="ADVANCED">
              %%Keyword: XFConstant <basic>
              <HELP>
              Input of an external Hessian matrix in cartesian coordinates. The
              syntax is the same as for the FCONSTANT keyword.
              </HELP>
              </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 MODULE="SLAPAF" NAME="NUMERICAL" APPEAR="Numerical Hessian" KIND="SINGLE" LEVEL="ADVANCED">
              %%Keyword: Numerical <basic>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

:kword:`CUBIc`
  This invokes a calculation of the 2nd and the 3rd order
  force constant matrix by finite difference formula.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="CUBIC" APPEAR="Numerical anharmonic force constants" KIND="SINGLE" LEVEL="ADVANCED">
              %%Keyword: Cubic <basic>
              <HELP>
              This invokes a calculation of the 2nd and the 3rd order
              force constant matrix by finite difference formula.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="DELTA" APPEAR="Numerical displacement value" KIND="REAL" REQUIRE="NUMERICAL.OR.CUBIC" LEVEL="ADVANCED">
              %%Keyword: Delta <basic>
              <HELP>
              This keyword is followed by a real number which defines the
              step length used in the finite differentiation. Default: 1.0D-2.
              </HELP>
              </KEYWORD>

:kword:`PRFC`
  The eigenvalues and eigenvectors of the Hessian matrix
  are printed. The internal coordinates definitions are also printed.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="PRFC" APPEAR="Print eigen vectors and values of H" KIND="SINGLE" LEVEL="ADVANCED">
              %%Keyword: PrFC <basic>
              <HELP>
              The eigenvalues and eigenvectors of the Hessian matrix
              are printed. The internal coordinates definitions is also printed.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="RHID" APPEAR="Hidden atoms selection radius" KIND="REAL" LEVEL="ADVANCED">
              %%Keyword: rHid <advanced>
              <HELP>
              Define the hidden atoms selection radius in order to improve a QM/MM Hessian.
              </HELP>
              </KEYWORD>

  .. xmldoc:: </GROUP>

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 MODULE="SLAPAF" NAME="CTOF" KIND="STRINGS" SIZE="2" LEVEL="ADVANCED">
              %%Keyword: CTOF <advanced>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="RTRN" KIND="CUSTOM" LEVEL="ADVANCED">
              %%Keyword: RTRN <advanced>
              <HELP>
              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".
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="THERMO" APPEAR="Thermochemical analysis" KIND="CUSTOM" LEVEL="ADVANCED">
              %%Keyword: THER <advanced>
              <HELP>
              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".
              </HELP>
              </KEYWORD>

:kword:`DISOtope`
  Calculates frequencies modified for double isotopic substitution.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="DISO" APPEAR="Double isotopic substitutions" KIND="SINGLE" LEVEL="ADVANCED">
              %%Keyword: DISOtope <advanced>
              <HELP>
              Calculates frequencies modified for double isotopic substitution.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="TRACK" APPEAR="Track states" KIND="SINGLE" LEVEL="BASIC">
              %%Keyword: Track <basic>
              <HELP>
              Tries to follow electronic states during an optimization, by computing state overlaps with RASSI.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="LASTENERGY" APPEAR="Last Energy method" KIND="STRING" LEVEL="ADVANCED">
              %%Keyword: LAST <advanced>
              <HELP>
              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.
              </HELP>
              </KEYWORD>

:kword:`NOLAst energy`
  Disables the call to the :program:`Last_Energy` module when convergence is achieved.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="NOLASTENERGY" APPEAR="No Last Energy" KIND="SINGLE" LEVEL="BASIC">
              %%Keyword: NoLastEnergy <basic>
              <HELP>
              Disables the call to the Last_Energy module when convergence is achieved.
              </HELP>
              </KEYWORD>

Optional restricted variance optimization (RVO) :cite:`Raggi2020,FdezGalvan2021,FdezGalvan2023` 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 MODULE="SLAPAF" NAME="KRIGING" APPEAR="Restricted variance optimization" KIND="SINGLE" LEVEL="ADVANCED">
              %%Keyword: Kriging <advanced>
              <HELP>
              Activate restricted variance optimization (RVO) using gradient-enhanced Kriging to describe the surrogate model.
              </HELP>
              </KEYWORD>

: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 MODULE="SLAPAF" NAME="TFOFFSET" APPEAR="Trend function offset" KIND="REAL" DEFAULT_VALUE="10.0" REQUIRE="KRIGING" LEVEL="ADVANCED">
              %%Keyword: TFOFfset <advanced>
              <HELP>
              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.
              </HELP>
              Default: 10.0 au.
              </KEYWORD>

: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.
  During the constrained phase of an optimization with :kword:`FindTS`, the default is 0.1 au.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="MAXDISP" APPEAR="Maximum dispersion factor" KIND="REAL" MIN_VALUE="0.0" DEFAULT_VALUE="0.3" REQUIRE="KRIGING" LEVEL="ADVANCED">
              %%Keyword: MAXDISP <advanced>
              <HELP>
              Maximum energy dispersion allowed during each macro iteration of the RVO procedure.
              A factor, multiplied by the maximum Cartesian Gradient.
              </HELP>
              Default: 0.3 au.
              </KEYWORD>

:kword:`MXMI`
  Maximum number of micro iterations in each macro iteration of the RVO procedure.
  If you set this to a small value, you may want to set :kword:`NOFAllback` too.
  The default value is 50.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="MXMI" APPEAR="Micro iterations" KIND="INT" MIN_VALUE="1" DEFAULT_VALUE="50" REQUIRE="KRIGING" LEVEL="ADVANCED">
              %%Keyword: MXMI <advanced>
              <HELP>
              Maximum number of micro iterations in each macro iteration of the RVO procedure.
              </HELP>
              Default: 50.
              </KEYWORD>

:kword:`NOFAllback`
  Disable fallback to conventional optimization if RVO microiterations do not converge.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="NOFALLBACK" APPEAR="No fallback to conventional" KIND="SINGLE" REQUIRE="KRIGING" LEVEL="ADVANCED">
              %%Keyword: NOFALLBACK <advanced>
              <HELP>
              Disable fallback to conventional optimization if RVO microiterations do not converge.
              </HELP>
              </KEYWORD>

:kword:`NDELta`
  Activate partial gradient enhanced Kriging, PGEK. This integer number determines for how many fewer iterations the gradients will
  be included in the PGEK procedure.
  The default value is 0, that is standard GEK.

  .. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="NDELTA" APPEAR="Samples without gradient" KIND="INT" MIN_VALUE="0" DEFAULT_VALUE="0" REQUIRE="KRIGING" LEVEL="ADVANCED">
              %%Keyword: NDELta <advanced>
              <HELP>
              Activate partial gradient enhanced Kriging, PGEK. This integer number determines for how many fewer iterations the gradients will
              be included in the PGEK procedure.
              </HELP>
              Default: 0.
              </KEYWORD>

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:: <KEYWORD MODULE="SLAPAF" NAME="REDUNDANT" KIND="SINGLE" LEVEL="UNDOCUMENTED" />

.. xmldoc:: <KEYWORD MODULE="SLAPAF" NAME="NOEMEP" KIND="SINGLE" LEVEL="UNDOCUMENTED" />

.. xmldoc:: </MODULE>

.. xmldoc:: <MODULE NAME="LAST_ENERGY">
            </MODULE>