.. _li_mgo_vib:

******************************
Vibrational frequency analysis
****************************** 

.. note:: The scripts for this section of the tutorial can be found in ``li_mgo/4_li-mgo_vib``. 
   It is worth highlighting that the transition state structure ``cimg.pun`` was obtained from the ``_dl_find.pun`` file generated in the 
   oxygen 7-5 NEB calculation. 

This section of the tutorial will cover how to perform a vibrational frequency analysis on the equatorial hole migration
transition state. 

First, the transition state is read into ChemShell using the script ``vib_analysis.py``. 
The rest of the script is similar to the geometry optimisation in :ref:`li_doping` section, with the exception of the ``thermal=True`` 
argument in the optimisation method. This specifies that a vibrational frequency analysis should be performed rather than an optimisation. This involves the construction of a finite-difference Hessian from which thermal contributions will be calculated. 

.. caution:: It is necessary to use the same active region as used in :ref:`li_doping` section (regions 1-3) to obtain meaningful results. 
   Therefore, it is highly recommended to use the taskfarmed parallelisation to ensure timely completion of the calculation.

Near the end of the ChemShell output file, DL-FIND will report the symmetrised Hessian matrix, eigenvalues and thermochemical analysis results. 
The key parameters that are of attention are the vibrational frequency modes, zero-point energy and the energetic and entropic vibrational 
contributions. The thermochemical analysis results should look similar to the condensed output illustrated below::

    Thermochemical analysis
    Temperature:     300.00 Kelvin
    Modes assumed to have zero vibrational frequency:  0
     Mode     Eigenvalue Frequency Vib.T.(K)        ZPE (H)   Vib. Ene.(H)      - T*S (H)
        1  -0.0000536602    37.656i
        2   0.0010256750   164.630   236.866   0.0003750549   0.0006238386  -0.0011988288
        3   0.0010404380   165.811   238.564   0.0003777444   0.0006218450  -0.0011923846

    ...
   
    total vibrational energy correction to E_electronic   0.2612139968 H
    total ZPE     864987.58705 J/mol
    total E vib   325126.72265 J/mol
    total S vib     1680.99020 J/mol/K
    Crossover temperature for tunnelling        8.62271 K
     
You should find that the thermochemical analysis results contain one mode with an imaginary frequency. Given that a transition state should only have one 
imaginary frequency, these results are consistent with the climbing image being representative of the equatorial hole migration transition state.

.. hidden:: note for tutorial document builders.
   Defining the active region as regions 1-3 (qm region, qm/mm interface, active mm) is necessary to yield the imaginary frequency, however is significantly
   more computationally expensive. Using solely the qm region (region 1) yields no imaginary frequencies, however is significantly less computationally expensive. 
   For the full sampling of the active region (includes regions 1-3), it takes about 3 hours with taskfarming parallelisation (4 workgroups of 10 processors) and 
   times out (72 hour limit) with non-taskfarming parallelisation (40 processors).