Nudged Elastic Band (NEB) Examples¶

Here are a few examples showing how the NEB method can be used to obtain the the path and transition state of a reaction.

HCN isomerization reaction with NEB¶

Download NEB_HCN.run

#!/bin/sh

# This example demonstrates the use of the Nudged Elastic Band method in AMS for
# finding a transition state of the HCN isomerisation reaction. A shell script
# used to run the example calculation is shown below.

$AMSBIN/ams <<eor Task NEB System Atoms C 0.000000 0.000000 0.000000 N 1.180000 0.000000 0.000000 H 2.196000 0.000000 0.000000 End End System final Atoms C 0.000000 0.000000 0.000000 N 1.163000 0.000000 0.000000 H -1.078000 0.000000 0.000000 End End NEB Images 9 Iterations 100 End Engine DFTB Model DFTB3 ResourcesDir DFTB.org/3ob-3-1 DispersionCorrection D3-BJ EndEngine eor  H2 dissociation on graphene¶ Download NEB_H2_on_graphene.run #!/bin/sh for map1 in no yes do for map2 in no yes do AMS_JOBNAME=test.map1=$map1.map2=$map2$AMSBIN/ams << eor

System
Atoms
C -1.23079526  1.45426666 -0.60065817
C  1.25628125  1.47378921 -0.00438520
C -2.44885253 -2.15934252 -0.51047194
C -0.01262403 -2.15933012 -0.51045137
C -1.23076689 -0.05455308 -0.60068696
C  1.25630128 -0.07404160 -0.00441877
C -2.44887495 -0.74178904 -0.51048970
C -0.01264534 -0.74177652 -0.51049375
H  1.25629472 -0.34738650  1.07475852
H  1.25630533  1.74708724  1.07481086
End
Lattice
4.97414302 0.0 0.0
0.0 4.30083513 0.0
End
MapAtomsToUnitCell $map1 End System final Atoms C -1.24330284 1.42711444 -0.25240444 C 1.23494804 1.42739453 -0.25493930 C -2.48233069 -2.12769513 -0.25358631 C -0.00426996 -2.12770161 -0.25349945 C -1.24330266 0.00893844 -0.25240978 C 1.23494831 0.00864800 -0.25496555 C -2.48233386 -0.70957077 -0.25358531 C -0.00426696 -0.70956458 -0.25349840 H 1.27176634 0.72925445 2.26688554 H 1.28137254 0.73053214 3.01215519 End Lattice 4.95651766 0.0 0.0 0.0 4.27331845 0.0 End MapAtomsToUnitCell$map2
End

Properties
NormalModes Yes
End

GeometryOptimization
HessianFree
End
End

NEB
ClimbingThreshold 0.01
OptimizeLattice Yes
End

Engine DFTB
Model DFTB3
ResourcesDir DFTB.org/3ob-3-1
DispersionCorrection D3-BJ
KSpace Quality=Basic
EndEngine

eor
done
done


Running multiple NEB calculations using PLAMS¶

Download NEB.plms

This example should be executed using PLAMS.

PLAMS documentation and tutorial

import os

# This PLAMS script locates TS for several reactions using the Nudged Elastic Band (NEB) method.
# For each reaction two xyz files are required: one corresponding to the reactant state, called
# '{}_initial.xyz', and one corresponding to the product state, called '{}_final.xyz' (e.g.
# 'MyMolecule_initial.xyz' and 'MyMolecule_final.xyz')

# The folder containing the xyz files:
xyz_folder = os.path.join(os.environ['AMSHOME'],'examples','AMS','NEB_PLAMS','xyz')

names = [name.rsplit('_initial.xyz')[0] for name in os.listdir(xyz_folder) if name.endswith('_initial.xyz')]

# Settings for the AMS driver
neb_sett = Settings()
neb_sett.input.ams.Properties.NormalModes = 'Yes'
neb_sett.input.ams.GeometryOptimization.Convergence.Step = 1.0e-3

# Settings for the engine (here we use the DFTB engine with semiempirical GFN1-xTB method)
engine_sett = Settings()
engine_sett.input.DFTB.Model='GFN1-xTB'

for name in sorted(names):

# For NEB we need two system blocks in the AMS input (the initial and final molecule).
# In PLAMS you can have multiple system blocks by passing the the AMSJob a dictionary of molecules.
# The 'keys' of the dictionary will be used as the headers of the System block
mols = {}
# The initial molecule should be in the main 'System' block (the main system block has no header).
# The key of the mols dictionary should therefore be an empty string
mols[''] = Molecule(os.path.join(xyz_folder,name+'_initial.xyz'))
# The final molecule should be specified in a system block with the header 'final'
mols['final'] = Molecule(os.path.join(xyz_folder,name+'_final.xyz'))

# Create and run the job:
job = AMSJob(molecule=mols, name=name, settings=neb_sett+engine_sett)
job.run()

print('')
print("System name: {}".format(name))
if job.ok():
pes_point_character = job.results.readrkf('AMSResults', 'PESPointCharacter', file='NEB_TS_final')
print("NEB calculation converged. PES Point character: {}".format(pes_point_character))
print("Left  TS barrier: {0:.6f} [Hartree]".format(job.results.readrkf('NEB', 'LeftBarrier')))
print("Right TS barrier: {0:.6f} [Hartree]".format(job.results.readrkf('NEB', 'RightBarrier')))
else:
print("Unsuccesful NEB calculation")
print('')