Nudged Elastic Band (NEB)
*************************

DFTB Catalysis AMS Driver Nudged Elastic Band (NEB) Transition state
Reaction barrier PLAMS

In this tutorial we will use the climbing-image Nudged Elastic Band
method (NEB) to identify the minimum energy path of:

* a chemical reaction in vacuum

* ionic diffusion in a solid


Transition state of a chemical reaction
=======================================

In this first example, we will find the minimum energy reaction path
and transition state of the following reaction:

   [image]

We will be using the DFTB engine with the GFN1-xTB model. This is a
computationally fast method (which is convenient for the purposes of a
tutorial) but it is not very accurate in predicting reaction energies.
For better accuracy, using the DFT engine ADF (or BAND) is
recommended.


Setting up and running the calculation
--------------------------------------

You can set up and run the calculations for this tutorial using either
the Graphical User Interface (GUI), the python library PLAMS or a bash
script.

-[ GUI ]-

   Start **AMSjobs**
   **1.** in the **menu bar** select **SCM → New input**

This will open a new **AMSinput** window.

   In **AMSinput**:
   **1.** Switch to the **DFTB panel**: [image: ADFPanel][image] **→** [image: DFTBPanel][image]
   **2.** Select **Task → NEB**

We will use the GFN1-xTB model, which is the default model as of the
AMS2020 release. Your AMSinput window should look like this:

   [image]

When setting up an NEB calculation we need to specify two systems: the
*initial* and *final* states of the reaction. The NEB algorithm will
then generate a set of *images* by interpolating between the initial
and final systems. This will be the initial approximation of the
reaction path, which will be optimized during the NEB calculation.

   **1.** Download the following xyz files: file "NEB_initial.xyz" and "NEB_final.xyz"
   **2.** Import the coordinates of the initial system: in the **menu bar**, select **File → Import Coordinates (System)...** and select the file **NEB_initial.xyz**
   **3.** Create a new molecule-tab: in the **menu bar**, select **Edit → New Molecule...**
   **4.** Import the coordinates of the final system in the newly created molecule tab: in the **menu bar**, select **File → Import Coordinates (System)...** and select the file **NEB_final.xyz**

You can easily switch between molecule tabs using the arrows at the
bottom of the **molecule drawing area**. Here, you'll also see the
molecule names, such as "C2H4N2O" for the first molecule and
"NEB_final" for the second. Note that the GUI has the following naming
convention:

* **Main System**: The first molecule is considered the main system,
  and is often automatically assigned its chemical formula as the
  name.

* **Additional Molecules**: For any additional molecules, the GUI uses
  by default the name from their original XYZ files.

Important:

  In NEB calculations, the order of the atoms in the initial and final
  system should be the same (if you provide an intermediate system,
  you should use a consistent atom-ordering for that too). The order
  of the atoms should be consistent because the images-interpolation
  algorithm maps the n-th atom of the initial system to the n-th atom
  of the final system.You can see the indices of the atoms by clicking
  in the **menu bar** on **Properties → Atom Info → Name → Show**. It
  is possible to change the order of the atoms in the **Coordinates**
  panel (in the **panel bar**: **Model → Coordinates**) using the
  **Move atom(s)** option.

Now, go to the NEB details panel where we will set up the NEB
calculation:

   **1.** Click on [image: MoreBtn][image] next to **Task → NEB** to go to the NEB details panel
   **2.** From the drop-down menu next to **initial system**, select **C2H4N2O** or **NEB_initial** depending on the name (may differ between AMS versions).
   **3.** From the drop-down menu next to **final system**, select **NEB_final**
   **4.** Check the **Characterize PES point** checkbox

Your AMSinput window should look like this:

   [image]

Tip:

  From most **AMSinput** panels you can jump to the relevant section
  of the user manual by clicking on [image: InfoBtn][image], which is
  located in the top-right corner of the panel.

We are now ready to run the calculation:

   In **AMSinput**:
   **1.** In the **menu bar**, click on **File → Save** and give it the name "NEB_tutorial"
   **2.** In the **menu bar**, click on **File → Run** . This will bring the **AMSjobs** window to the front
   **3.** Wait for the calculation to finish. It should take just a few seconds

In the **logfile** you can monitor the progress of your NEB
calculation:

   In **AMSjobs**:
   **1.** Select the job "NEB_tutorial" and in the **menu bar** click on **SCM → Logfile**. This will open the logfile

A NEB calculation consists of several steps, which are automatically
executed one after the other:

* first, the two end points (the *initial* and *final* molecules) are
  optimized (in the logfile, look for "Optimizing initial state" and
  "Optimizing final state")

* then the NEB reaction path will be optimized (in the logfile, look
  for "NEB Path Optimization"). During the reaction path optimization,
  the highest-energy image on the path will *climb* to the transition
  state

* finally, a single point calculation for the TS is performed (in the
  logfile: "Final calculation for highest-energy image"). If the
  option **Characterize PES point** is on, then the lowest-lying
  normal modes will be calculated in order to validate the TS (the TS
  should have exactly one imaginary frequency). Some information on
  the reaction path is printed at the end of the logfile:

     NEB found a transition state!
     TS barrier height from the left           0.02576078 Hartree
                                              16.165 kcal/mol
                                              67.635 kJ/mol
     TS barrier height from the right          0.08632064 Hartree
                                              54.167 kcal/mol
                                             226.635 kJ/mol

-[ Bash script ]-

The following bash script performs an NEB calculation using the AMS
driver and the DFTB engine. The input options for the AMS driver
described in the AMS driver manual, while the DFTB manual describes
the input options for the **DFTB** "Engine" block.

   #!/bin/sh

   "$AMSBIN/ams" << eor

   Task NEB

   Properties
       PESPointCharacter Yes
   End

   # The initial system (no header)
   System
       Atoms
           N  1.88630912 -0.34204867 -1.59424245
           N  1.14203025 -1.15084766 -1.95458206
           O  0.07639739 -1.27682918 -1.40578886
           C  0.74766633  1.24132374  1.29097263
           C -0.52735637  0.87051548  1.20080269
           H  1.08791192  2.15092371  0.80682956
           H -1.23438622  1.47569873  0.64278673
           H -0.86675778 -0.04061999  1.68265930
           H  1.45534611  0.63457332  1.84647164
       End
   End

   # The final system should be called 'final'
   System final
       Atoms
           N  1.44280525  0.39326405  0.02115802
           N  0.56588608 -0.32983790 -0.53167497
           O -0.68785300 -0.26603404 -0.00725496
           C  0.86550207  1.12517445  1.10513146
           C -0.57889386  0.65549932  1.05974328
           H  0.94852175  2.21696742  0.91684905
           H -1.26384434  1.51200042  0.88137995
           H -0.85746564  0.15513847  2.01191277
           H  1.35845740  0.84559176  2.06071421
       End
   End

   Engine DFTB
       Model GFN1-xTB
   EndEngine
   eor

Important:

  In NEB calculations, the order of the atoms in the initial and final
  system should be the same (if you provide an intermediate system,
  you should use a consistent atom-ordering for that too). The order
  of the atoms should be consistent because the images-interpolation
  algorithm maps the n-th atom of the initial system to the n-th atom
  of the final system.

To run the calculation:

   **1.** Save the above script in a file called "NEB.run"
   **2.** Make the script executable: in a terminal, type "chmod +x NEB.run"
   **3.** Run it interactively and redirect the output to a file: in a terminal, type "./NEB.run > out"

-[ Python ]-

See the NEB scripting example.

To improve the initial approximation of the reaction path in an NEB
calculation, you can (optionally) provide an *intermediate* system.

Another important NEB option is the **number of images**. In case of
problematic NEB path optimization convergence, using more images might
help (note that the computation time increases with the number of
images used).

You can read more about the various NEB options in the AMS manual.


Results of the calculation
--------------------------

-[ GUI ]-

Now, let's visualize the reaction path computed by NEB:

   In **AMSjobs**:
   **1.** Select the job "NEB_tutorial" and, in the **menu bar**, click on **SCM → Movie**. This will open the **AMSmovie** module

   [image]

In AMSmovie, you can click on **play** (or drag the slide-bar) so see
the reaction happening. On right-hand side, the energy and gradients
of the images in the NEB reaction path are plotted.

The transition state is characterized by having one imaginary
frequency. Let's visualize the normal modes of the transition state
geometry with AMSspectra:

   In **AMSmovie**:
   **1.** click on **SCM → Spectra**. This will open the **AMSspectra**

   [image]

Here you will see the computed normal modes for the TS geometry.
Notice that there is one negative frequency (imaginary frequency are
shown as negative numbers).

   In **AMSspectra**:
   **1.** In the table, click on the line with the negative frequency

The corresponding normal mode will be displayed in the molecule-
visualization area. This normal mode gives you an idea of how the
atoms are moving as they cross the transition state.

-[ Bash script ]-

In the folder where you executed your script you will find a file
"out", which contains the text-output of the calculation, and a folder
called "ams.results" containing binary results of the calculation.

At the end of the output file "out" you will find a section
summarizing the results of the NEB calculation:

   NEB found a transition state!


   ------------------------------------------------------------
   TS barrier height from the left           0.02581511 Hartree
                                            16.199 kcal/mol
                                            67.778 kJ/mol
   TS barrier height from the right          0.08637041 Hartree
                                            54.198 kcal/mol
                                           226.765 kJ/mol
   Reaction energy                          -0.06055530 Hartree
                                           -37.999 kcal/mol
                                          -158.988 kJ/mol
   ------------------------------------------------------------
           Transition state geometry

   --------
   Geometry
   --------

   Atoms
   Index Symbol   x (angstrom)   y (angstrom)   z (angstrom)
      1      N     1.76794468     0.02199401    -0.76530642
      2      N     0.86568320    -0.56021023    -1.14090273
      3      O    -0.28459962    -0.76023024    -0.92163051
      4      C     0.89983351     0.97770099     0.84893165
      5      C    -0.38912967     0.55813520     0.80607959
      6      H     1.18303105     1.94163319     0.44757192
      7      H    -1.15509240     1.14579400     0.31965225
      8      H    -0.73127815    -0.27639812     1.40152984
      9      H     1.61076814     0.51427067     1.51998360

The folder "ams.results" contains:

* a text file called "ams.log" containing a very concise summary of
  the calculation’s progress during the run.

* the main binary results file "ams.rkf",  containing the reaction
  path computed in the NEB calculation.

* the engine results file "NEB_TS_final.rkf" corresponding a single
  point calculation at the transition state geometry. It contains,
  among other properties, the normal modes.

You can explore the content of the "rkf" binary results files using
the **kfbrowser** utility.

   In a terminal, type: "$AMSBIN/kfbrowser ams.results/NEB_TS_final.rkf"

The binary results of the calculation can also be visualized with the
GUI modules:

   In a terminal, type: "$AMSBIN/amsmovie ams.results/ams.rkf"
   In a terminal, type: "$AMSBIN/amsspectra ams.results/NEB_TS_final.rkf"


Ionic diffusion in a solid
==========================

In this second example, we will determine the minimum energy path and
the activation energy of Li diffusion in a typical cathode materials
LiTiS_2.

We will use the ML potential engine with the M3GNet universal machine
learning potential. Although M3GNet was trained to reproduce ground
state structures of the Materials Project database, it can also be
accurate to describe transition states, as we will demonstrate in this
tutorial. This tutorial can be reproduced with the periodic DFT engine
BAND.

Important:

  **First, make sure to install the M3GNet package** inside AMS as
  described in the M3GNet tutorial.

Note:

  When using M3GNet, you can include dispersion via an "engine addon"
  option. In the following, D3 dispersion will be included throughout.


Layer LiTiS_2 cathode material
------------------------------

In some Li-based cathode materials, the diffusion of Li ions is
tightly connected to the content of Li in the materials. We would like
to demonstrate this by computing

* the diffusion barrier of Li in LiTiS_2 with 1 Li vacancy, and

* the diffusion barrier of Li in LiTiS_2 with 2 Li vacancies.

In LiTiS_2, the lithium ions are octahedrally coordinated between the
S−Ti−S sandwich layers.


Li diffusion in LiTiS_2 with 1 vacancy
--------------------------------------

To perform the NEB calculation we will use the graphical user
interface.  Either

* build the systems manually (next step), or

* import them directly (skip to the next section)


Build LiTiS_2 initial and final structures with 1 vacancy
---------------------------------------------------------

We will start with a supercell of crystalline LiTiS_2, constructed
from the optimized lattice parameters, and then introduce a vacancy.


Import the pristine system
~~~~~~~~~~~~~~~~~~~~~~~~~~

   Start AMSinput: **SCM → New Input**

Now, either

* download "LiTiS2_supercell_pristine.in" and import it (**File →
  Import coordinates...**), or

* copy-paste the below coordinates into AMSinput.

  -[ Reveal copy-pastable coordinates ]-

     System
         Atoms
             Li 1.71915 2.9776551458 6.174905 
             Ti 1.71915 2.9776551458 3.0874525 
             S 1.71915 4.9627585764 4.5513282748 
             S 3.4383 3.9702068611 7.7984817252 
             Li 3.4383 5.9553102917 0.0 
             Ti 3.4383 5.9553102917 9.2623575 
             S 3.4383 7.9404137222 10.7262332748 
             S 5.15745 6.947862007 1.6235767252 
             Li 3.4383 5.9553102917 6.174905 
             Ti 3.4383 5.9553102917 3.0874525 
             S 3.4383 7.9404137222 4.5513282748 
             S 5.15745 6.947862007 7.7984817252 
             Li 6.8766 0.0 0.0 
             Ti 6.8766 0.0 9.2623575 
             S 6.8766 1.9851034305 10.7262332748 
             S 8.595750000000001 0.9925517153 1.6235767252 
             Li 6.8766 0.0 6.174905 
             Ti 6.8766 0.0 3.0874525 
             S 6.8766 1.9851034305 4.5513282748 
             S 8.595750000000001 0.9925517153 7.7984817252 
             Li 5.15745 2.9776551458 0.0 
             Ti 5.15745 2.9776551458 9.2623575 
             S 5.15745 4.9627585764 10.7262332748 
             S 6.8766 3.9702068611 1.6235767252 
             Li 5.15745 2.9776551458 6.174905 
             Ti 5.15745 2.9776551458 3.0874525 
             S 5.15745 4.9627585764 4.5513282748 
             S 6.8766 3.9702068611 7.7984817252 
             Li -3.4383 5.9553102917 0.0 
             Ti -3.4383 5.9553102917 9.2623575 
             S -3.4383 7.9404137222 10.7262332748 
             S -1.71915 6.947862007 1.6235767252 
             Li -3.4383 5.9553102917 6.174905 
             Ti -3.4383 5.9553102917 3.0874525 
             S -3.4383 7.9404137222 4.5513282748 
             S -1.71915 6.947862007 7.7984817252 
             Li 0.0 0.0 0.0 
             Ti 0.0 0.0 9.2623575 
             S 0.0 1.9851034305 10.7262332748 
             S 1.71915 0.9925517153 1.6235767252 
             Li 0.0 0.0 6.174905 
             Ti 0.0 0.0 3.0874525 
             S 0.0 1.9851034305 4.5513282748 
             S 1.71915 0.9925517153 7.7984817252 
             Li -1.71915 2.9776551458 0.0 
             Ti -1.71915 2.9776551458 9.2623575 
             S -1.71915 4.9627585764 10.7262332748 
             S 0.0 3.9702068611 1.6235767252 
             Li -1.71915 2.9776551458 6.174905 
             Ti -1.71915 2.9776551458 3.0874525 
             S -1.71915 4.9627585764 4.5513282748 
             S 0.0 3.9702068611 7.7984817252 
             Li 0.0 5.9553102917 0.0 
             Ti 0.0 5.9553102917 9.2623575 
             S 0.0 7.9404137222 10.7262332748 
             S 1.71915 6.947862007 1.6235767252 
             Li 0.0 5.9553102917 6.174905 
             Ti 0.0 5.9553102917 3.0874525 
             S 0.0 7.9404137222 4.5513282748 
             S 1.71915 6.947862007 7.7984817252 
             Li 3.4383 0.0 0.0 
             Ti 3.4383 0.0 9.2623575 
             S 3.4383 1.9851034305 10.7262332748 
             S 5.15745 0.9925517153 1.6235767252 
             Li 3.4383 0.0 6.174905 
             Ti 3.4383 0.0 3.0874525 
             S 3.4383 1.9851034305 4.5513282748 
             S 5.15745 0.9925517153 7.7984817252 
             Li 1.71915 2.9776551458 0.0 
             Ti 1.71915 2.9776551458 9.2623575 
             S 1.71915 4.9627585764 10.7262332748 
             S 3.4383 3.9702068611 1.6235767252 
         End
         Lattice
             10.3149 0.0 0.0
             -5.15745 8.9329654375 0.0
             0.0 0.0 12.34981
         End
     End



Create the "initial" system for NEB with a Li vacancy
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

Orient the crystal along the x-axis (Ctrl+1), then rotate the view
slightly the molecule for a better view of the middle Li atom in the
top layer (index 25).

   **Select → Select By Index → 25**
   In the bottom right corner, you see the coordinates for this atom. Write them down because you will need them later: "5.157 2.978 6.175"
   **Atoms → Delete Atom(s)** (or press the 'Delete' key on your keyboard)
   In the bottom right corner, verify that the stoichiometry is Li17 S36 Ti18  (one fewer Li than Ti).
 [image]

The Li atom that will diffuse has index 1. Let's highlight it with a
region for better visibility. Regions can be used to affect the
calculation, but here we just use it to better see the diffusing atom.

   Click outside any atoms so that no atom is selected
   **Select → Select by Index → 1**
   **Atoms → New Region From Selected Atoms**

This will highlight atom 1 with a colored orb around it.

[image]

Now export this structure to a ".in" file, so that we do not
accidentally lose it:

   File → Export coordinates → .in
   Save it with the filename "LiTiS2_vacancy_initial.in"


Create the "final" system for NEB with a Li vacancy
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

NEB calculations have both a starting point (initial system) and an
end point (final system).

We will create the final system by changing the coordinates of atom
"1" to the coordinates of the atom that was removed.

You can either do this

* manually by first selecting atom 1 and then right-click-dragging the
  atom into the approximate position of the previous vacancy, or

* just set the coordinates directly:

     In the panel bar, go to **Model → Coordinates**
     Change the coordinates of the first atom to "5.157 2.978 6.175" (the values for the atom you removed earlier)
     **File → Export coordinates → .in**
     Use the file name "LiTiS2_vacancy_final.in"

[image]


Set up the NEB calculation
--------------------------


Import the initial and final systems
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

Note:

  * By default, the NEB task will optimize both the initial and final
    points before interpolating between them.

  * Make sure that the initial and final structures are close to their
    respective local minima!

  * Best practice is to explicitly optimize them and verify that the
    structures look good before running the NEB.

  * When you perform NEB in a crystal it is recommended to keep the
    cell dimension between the initial and final states.

Create a new clean input:

   **SCM → New Input**

Now, either

* **File → Import Coordinates →** "LiTiS2_vacancy_initial.in"
  (download or use your own), or

* -[ Reveal copy-pasteable coordinates (LiTiS2_vacancy_initial.in) ]-

     System
         Atoms
             Li 1.71915 2.9776551458 6.174905 region=Region_1
             Ti 1.71915 2.9776551458 3.0874525 
             S 1.71915 4.9627585764 4.5513282748 
             S 3.4383 3.9702068611 7.7984817252 
             Li 3.4383 5.9553102917 0.0 
             Ti 3.4383 5.9553102917 9.2623575 
             S 3.4383 7.9404137222 10.7262332748 
             S 5.15745 6.947862007 1.6235767252 
             Li 3.4383 5.9553102917 6.174905 
             Ti 3.4383 5.9553102917 3.0874525 
             S 3.4383 7.9404137222 4.5513282748 
             S 5.15745 6.947862007 7.7984817252 
             Li 6.8766 0.0 0.0 
             Ti 6.8766 0.0 9.2623575 
             S 6.8766 1.9851034305 10.7262332748 
             S 8.595750000000001 0.9925517153 1.6235767252 
             Li 6.8766 0.0 6.174905 
             Ti 6.8766 0.0 3.0874525 
             S 6.8766 1.9851034305 4.5513282748 
             S 8.595750000000001 0.9925517153 7.7984817252 
             Li 5.15745 2.9776551458 0.0 
             Ti 5.15745 2.9776551458 9.2623575 
             S 5.15745 4.9627585764 10.7262332748 
             S 6.8766 3.9702068611 1.6235767252 
             Ti 5.15745 2.9776551458 3.0874525 
             S 5.15745 4.9627585764 4.5513282748 
             S 6.8766 3.9702068611 7.7984817252 
             Li -3.4383 5.9553102917 0.0 
             Ti -3.4383 5.9553102917 9.2623575 
             S -3.4383 7.9404137222 10.7262332748 
             S -1.71915 6.947862007 1.6235767252 
             Li -3.4383 5.9553102917 6.174905 
             Ti -3.4383 5.9553102917 3.0874525 
             S -3.4383 7.9404137222 4.5513282748 
             S -1.71915 6.947862007 7.7984817252 
             Li 0.0 0.0 0.0 
             Ti 0.0 0.0 9.2623575 
             S 0.0 1.9851034305 10.7262332748 
             S 1.71915 0.9925517153 1.6235767252 
             Li 0.0 0.0 6.174905 
             Ti 0.0 0.0 3.0874525 
             S 0.0 1.9851034305 4.5513282748 
             S 1.71915 0.9925517153 7.7984817252 
             Li -1.71915 2.9776551458 0.0 
             Ti -1.71915 2.9776551458 9.2623575 
             S -1.71915 4.9627585764 10.7262332748 
             S 0.0 3.9702068611 1.6235767252 
             Li -1.71915 2.9776551458 6.174905 
             Ti -1.71915 2.9776551458 3.0874525 
             S -1.71915 4.9627585764 4.5513282748 
             S 0.0 3.9702068611 7.7984817252 
             Li 0.0 5.9553102917 0.0 
             Ti 0.0 5.9553102917 9.2623575 
             S 0.0 7.9404137222 10.7262332748 
             S 1.71915 6.947862007 1.6235767252 
             Li 0.0 5.9553102917 6.174905 
             Ti 0.0 5.9553102917 3.0874525 
             S 0.0 7.9404137222 4.5513282748 
             S 1.71915 6.947862007 7.7984817252 
             Li 3.4383 0.0 0.0 
             Ti 3.4383 0.0 9.2623575 
             S 3.4383 1.9851034305 10.7262332748 
             S 5.15745 0.9925517153 1.6235767252 
             Li 3.4383 0.0 6.174905 
             Ti 3.4383 0.0 3.0874525 
             S 3.4383 1.9851034305 4.5513282748 
             S 5.15745 0.9925517153 7.7984817252 
             Li 1.71915 2.9776551458 0.0 
             Ti 1.71915 2.9776551458 9.2623575 
             S 1.71915 4.9627585764 10.7262332748 
             S 3.4383 3.9702068611 1.6235767252 
         End
         Lattice
             10.3149 0.0 0.0
             -5.15745 8.9329654375 0.0
             0.0 0.0 12.34981
         End
     End


Create a new "molecule" in AMSinput for the final state:

   **Edit → New Molecule**

Now, either

* **File → Import Coordinates →** "LiTiS2_vacancy_final.in" (download
  or use your own), or

* -[ Reveal copy-pasteable coordinates (LiTiS2_vacancy_final.in) ]-

     System
         Atoms
             Li 5.157 2.9776551458 6.174905 region=Region_1
             Ti 1.71915 2.9776551458 3.0874525 
             S 1.71915 4.9627585764 4.5513282748 
             S 3.4383 3.9702068611 7.7984817252 
             Li 3.4383 5.9553102917 0.0 
             Ti 3.4383 5.9553102917 9.2623575 
             S 3.4383 7.9404137222 10.7262332748 
             S 5.15745 6.947862007 1.6235767252 
             Li 3.4383 5.9553102917 6.174905 
             Ti 3.4383 5.9553102917 3.0874525 
             S 3.4383 7.9404137222 4.5513282748 
             S 5.15745 6.947862007 7.7984817252 
             Li 6.8766 0.0 0.0 
             Ti 6.8766 0.0 9.2623575 
             S 6.8766 1.9851034305 10.7262332748 
             S 8.595750000000001 0.9925517149999999 1.6235767252 
             Li 6.8766 0.0 6.174905 
             Ti 6.8766 0.0 3.0874525 
             S 6.8766 1.9851034305 4.5513282748 
             S 8.595750000000001 0.9925517149999999 7.7984817252 
             Li 5.15745 2.9776551458 0.0 
             Ti 5.15745 2.9776551458 9.2623575 
             S 5.15745 4.9627585764 10.7262332748 
             S 6.8766 3.9702068611 1.6235767252 
             Ti 5.15745 2.9776551458 3.0874525 
             S 5.15745 4.9627585764 4.5513282748 
             S 6.8766 3.9702068611 7.7984817252 
             Li -3.4383 5.9553102917 0.0 
             Ti -3.4383 5.9553102917 9.2623575 
             S -3.4383 7.9404137222 10.7262332748 
             S -1.71915 6.947862007 1.6235767252 
             Li -3.4383 5.9553102917 6.174905 
             Ti -3.4383 5.9553102917 3.0874525 
             S -3.4383 7.9404137222 4.5513282748 
             S -1.71915 6.947862007 7.7984817252 
             Li 0.0 0.0 0.0 
             Ti 0.0 0.0 9.2623575 
             S 0.0 1.9851034305 10.7262332748 
             S 1.71915 0.9925517149999999 1.6235767252 
             Li 0.0 0.0 6.174905 
             Ti 0.0 0.0 3.0874525 
             S 0.0 1.9851034305 4.5513282748 
             S 1.71915 0.9925517149999999 7.7984817252 
             Li -1.71915 2.9776551458 0.0 
             Ti -1.71915 2.9776551458 9.2623575 
             S -1.71915 4.9627585764 10.7262332748 
             S 0.0 3.9702068611 1.6235767252 
             Li -1.71915 2.9776551458 6.174905 
             Ti -1.71915 2.9776551458 3.0874525 
             S -1.71915 4.9627585764 4.5513282748 
             S 0.0 3.9702068611 7.7984817252 
             Li 0.0 5.9553102917 0.0 
             Ti 0.0 5.9553102917 9.2623575 
             S 0.0 7.9404137222 10.7262332748 
             S 1.71915 6.947862007 1.6235767252 
             Li 0.0 5.9553102917 6.174905 
             Ti 0.0 5.9553102917 3.0874525 
             S 0.0 7.9404137222 4.5513282748 
             S 1.71915 6.947862007 7.7984817252 
             Li 3.4383 0.0 0.0 
             Ti 3.4383 0.0 9.2623575 
             S 3.4383 1.9851034305 10.7262332748 
             S 5.15745 0.9925517149999999 1.6235767252 
             Li 3.4383 0.0 6.174905 
             Ti 3.4383 0.0 3.0874525 
             S 3.4383 1.9851034305 4.5513282748 
             S 5.15745 0.9925517149999999 7.7984817252 
             Li 1.71915 2.9776551458 0.0 
             Ti 1.71915 2.9776551458 9.2623575 
             S 1.71915 4.9627585764 10.7262332748 
             S 3.4383 3.9702068611 1.6235767252 
         End
         Lattice
             10.3149 0.0 0.0
             -5.15745 8.9329654375 0.0
             0.0 0.0 12.34981
         End
     End


You can toggle between the initial and final structures using the
arrows below the 3D area.


Set up the engine M3GNet settings
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

   In the yellow dropdown, switch to [image: MLPotentialPanel][image]
   Set **Model → M3GNet → M3GNet-UP-2022**
   On the panel bar, go to **Details → Engine Addons**
   Tick the **D3 Dispersion → Yes**


Set up the NEB settings
~~~~~~~~~~~~~~~~~~~~~~~

   **1.** From the **Main** tab, select the **Task → NEB**
   **2.** Click the [image: MoreBtn][image] next to **Task → NEB**
   **3.** Set the **Maximum number of iterations** to 100
   **4.** Set the **Number of images** to 19
   **5.** Set **Initial system** to **LiTiS2_vacancy_initial** (or it may be called Mol-1)
   **6.** Set **Final system** to  **LiTiS2_vacancy_final** (or it may be called Mol-2)

The NEB is now ready. Save the calculation. This will automatically
rename the molecules as **initial** and **final**.

   [image]

Run the calculation. This might take a few more minutes as AMS will
perform a simultaneous optimization of all the 17 images by minimizing
the forces parallel to the reaction path.

Once the calculation is finished, select **LiTiS2_NEB** in **AMSjobs**
and click **SCM → Movie**. You can re-orient the crystal and use the
slidebar under the molecule area to appreciate the diffusion along the
minimum energy path. If you double click the y-axis label of the
graph, you can switch to your prefered units. We found the activation
energy E_A = 0.59 eV.

   [image]


NEB with Python
---------------

See the Li diffusion NEB scripting example.


Li diffusion in Li_xTiS_2
-------------------------

To explore the effect of lithium concentration on diffusion, begin
with the structures from step 3. Remove the Li atom with index 9 for
both the initial and final structures and at the end increase the NEB
iterations to 500.

This change reduce the activation energy by approximately a factor of
2.

* In the first case, the lithium migration trajectory was roughly
  linear.

* The presence of the additional Li vacancy causes the path to curve,
  passing through a tetrahedral site and creating a shallow energy
  minimum.

   [image]
