.. _geom_opt:

*********************
Geometry Optimisation
*********************

.. note:: The scripts for this step of the tutorial can be found in 
   ``mgo_surface/3_mgo_opt``.

In this section we demonstrate an energy minimisation of the MgO cluster.
The input script ``mgo_opt_energy.py`` is very similar in form to the single
point calculation and you can refer back to the previous step for 
explanations of the identical parts.

To run a geometry optimisation we need to know which atoms are active, that 
is can be relaxed (corresponding to regions 1-3 of the ionic embedded 
cluster model). We can use the ``getRegion()`` method to get a list of 
active atoms from the input fragment in the same way as for the QM region.

.. literalinclude:: ../../samples/mgo_surface/3_mgo_opt/mgo_opt_energy.py
    :lines: 17-19

The QM/MM settings are the same as for the single point calculation with 
one exception: ``shl_maxcycles`` is reduced to 5. This makes sure that the
geometry optimisation does not get stuck in a long sequence of shell 
relaxations at any particular step, while it should be sufficient to 
ensure that the shells are fully relaxed by the end of the geometry 
optimisation.

The energy minimisation is run by calling DL-FIND, the standard geometry 
optimisation module for ChemShell. An instance of the ``Opt()`` class 
is first set up containing DL-FIND options, and the optimisation 
run using the ``run()`` method.

.. literalinclude:: ../../samples/mgo_surface/3_mgo_opt/mgo_opt_energy.py
    :lines: 39-43

Here the geometry optimisation has several variables:

 * ``algorithm`` is used to specify the optimisation algorithm used by DL-FIND. 
   Here we use the default L-BFGS algorithm which is suitable for both 
   small and large systems.
 * ``maxcycle`` determines the maximum number of optimisation cycles. 
   The optimiser will cease after convergence or at the limit set here.
 * ``maxene`` determines the maximum number of energy and gradient calculations
   during the optimisation. The optimiser will cease if this is reached before
   convergence.
 * ``active`` specifies the atoms that will be relaxed during the optimisation.
   Here, we pass in the ``active_region`` list of region 1-3 atoms.

With NWChem running on 4 cores, the optimisation will run in approximately
45 minutes.

After the optimisation completes, the optimised structures are saved
for later use in punch file format, and for visualisation in an XYZ file.

.. literalinclude:: ../../samples/mgo_surface/3_mgo_opt/mgo_opt_energy.py
    :lines: 47-49

Finally, the result is validated against the pre-calculated optimised energy.
