Diffusion barriers with FHI-aims and aimsChain
Estimated total CPU time: < 5 min when using 70 cores or more
Info
This is a concise tutorial created for the FHI-aims webinar series given in January 2025. It introduces participants to finding minimum energy pathways with aimsChain and covers the most important aspects, but it is far from exhaustive. For learning more about these kinds of calculations, we strongly recommend to also check the detailed description of aimsChain in the FHI-aims manual.
A quick look at diffusion in argyrodite materials
In this tutorial, we will focus on the ordered and most stable Li\(_6\)PS\(_5\)I structure that was described in the previous section. It crystallizes in the cubic F4̅3m spacegroup, but instead of the conventional cell, we will focus on a 2x2x1 supercell of the primitive cell, as done in Ref. 1. The larger supercell is necessary to allow sufficient space for studying the Li diffusion pathway without artificial self-interaction across periodic boundaries.
Here, we will examine one of the possible Li diffusion pathways in Li\(_6\)PS\(_5\)I. The pathway represents inter-cage Li hopping with a jump distance of approximately 3.3 Å between initial and final Li sites. This particular path is also studied in Ref. 1 and was chosen as it represents one of the key diffusion mechanisms identified in previous studies of argyrodite materials.
Warning
In the following exercises, the computational settings, in particular the reciprocal space k-grid and the basis set have been chosen to allow a rapid computation of the exercises. Although the qualitative trends already hold with the present settings, the reciprocal space grid, the basis set, and other settings would all have to be converged with much more care for real production calculations.
Setting up the diffusion pathway
We now want to analyze the chosen diffusion pathway with FHI-aims and aimsChain. Start by creating a new working folder and then navigating to that folder with
In here we need to create files that describe the geometries and control settings of our chosen Li diffusion pathway:
ini.in- this contains the relaxed initial structure of Li\(_6\)PS\(_5\)Ifin.in- this contains the relaxed final structure with the displaced Li ioncontrol.in- this contains the control parameters for all FHI-aims calculations
Warning
The ordering of atoms in the geometry files is crucial. The atoms in the initial and final geometry must establish a one-to-one correspondence, so that the \(i\)th line in the initial and final geometry represents the same atom. Changing the order can change the result from aimsChain dramatically.
Now create both files ini.in and fin.in and copy in their respective blocks below:
Geometry Files
lattice_vector 12.35227200 0.00000000 -0.02798300
lattice_vector 0.00000000 7.13184600 0.00000000
lattice_vector -3.98783100 0.00000000 11.43899700
atom 11.44944670 3.90555555 1.81905657 Li
atom 8.73788177 4.76080482 10.25018107 Li
atom 4.61237925 5.84536513 4.50701941 Li
atom 2.57635956 1.29383338 10.19385703 Li
atom 6.56836403 0.02612659 10.02503012 Li
atom 8.63234725 7.08189692 4.27952525 Li
atom 0.40884066 3.67034667 10.09849200 Li
atom 2.28930378 3.45786468 4.38116762 Li
atom 8.59902355 3.61520438 4.21888385 Li
atom 6.23878689 3.49091584 10.05514666 Li
atom 2.01364876 0.11823780 4.37814653 Li
atom 0.08566444 7.04806818 10.05567843 Li
atom 8.48524078 5.36171711 0.57981653 Li
atom 6.47002911 1.80728649 6.19843469 Li
atom 2.31654889 1.79866782 0.56075555 Li
atom 0.34530081 5.32336220 6.31290477 Li
atom 11.49579618 6.73814848 1.74829175 Li
atom 9.59655528 0.43510331 7.40372738 Li
atom 5.30817298 3.08550990 1.75024947 Li
atom 3.38073205 4.02300416 7.42588009 Li
atom 9.63444526 3.28563296 7.28145970 Li
atom 3.36367896 6.80267265 7.38329597 Li
atom 5.29080224 0.32361158 1.67851901 Li
atom 6.78020686 5.36742252 7.04502873 P
atom 8.68144896 1.74385067 1.37996620 P
atom 0.56244974 1.82926996 7.04679253 P
atom 2.53271148 5.31408815 1.38631375 P
atom 5.78986385 7.05811918 7.68504298 S
atom 7.67539505 0.06615878 2.02088499 S
atom -0.37164048 3.50263955 7.81227211 S
atom 1.54685824 3.63718877 2.04424072 S
atom 7.69090141 3.45477315 1.94974171 S
atom 5.79146903 3.67301660 7.67855663 S
atom 1.57966219 7.03028121 1.99100164 S
atom -0.37133770 0.11187571 7.70987191 S
atom 10.49552565 1.74365289 2.32277168 S
atom 8.66690706 5.40227927 7.87258713 S
atom 4.45530028 5.27863201 2.12413044 S
atom 2.49552903 1.84349870 7.80093778 S
atom 6.92070936 5.35425734 5.00284446 S
atom 4.86500639 1.66905668 10.76736308 S
atom 0.52003346 1.89492018 5.02397032 S
atom -1.35351809 5.26788562 10.77090072 S
atom 8.79613259 1.80887369 5.70252207 S
atom 6.79336844 5.28929156 11.45287351 S
atom 2.62916141 5.39175283 5.67864388 S
atom 0.61580243 1.73929175 11.40658209 S
atom 10.85193142 5.34618083 4.16924709 I
atom 8.78823320 1.88813363 9.84108979 I
atom 4.57121227 1.65198387 4.14629517 I
atom 2.65482040 5.47470060 9.85201522 I
lattice_vector 12.35227200 0.00000000 -0.02798300
lattice_vector 0.00000000 7.13184600 0.00000000
lattice_vector -3.98783100 0.00000000 11.43899700
atom 10.75363197 2.37897437 4.34043094 Li
atom 8.75977785 4.68869299 10.33132916 Li
atom 4.56440813 5.88534366 4.50877747 Li
atom 2.53557745 1.27114863 10.21301831 Li
atom 6.61030961 0.09062558 10.07594822 Li
atom 8.58053095 7.06004835 4.33981703 Li
atom 0.39758715 3.67693025 10.09393514 Li
atom 2.33582676 3.53410477 4.29311712 Li
atom 8.26335126 3.69382454 4.30966063 Li
atom 6.27646702 3.37223316 10.05643481 Li
atom 2.06216647 0.11342163 4.33213151 Li
atom 0.04538877 7.08253418 10.06754877 Li
atom 8.49198240 5.32218886 0.56723452 Li
atom 6.53795499 1.78082201 6.27922545 Li
atom 2.30652166 1.85608297 0.60398362 Li
atom 0.33320669 5.36022145 6.30130757 Li
atom 11.54627829 6.42435822 1.77279368 Li
atom 9.57052172 0.44615834 7.37479799 Li
atom 5.39097198 3.11148591 1.72677215 Li
atom 3.36628757 4.03914027 7.42962913 Li
atom 9.55980302 3.27725497 7.34560571 Li
atom 3.34729642 6.81594915 7.40174153 Li
atom 5.34138287 0.33052365 1.69296170 Li
atom 6.73719257 5.39537905 7.12149385 P
atom 8.77622126 1.70069352 1.36061631 P
atom 0.56696710 1.84696818 7.09495290 P
atom 2.54596907 5.32457745 1.38252956 P
atom 5.75297738 7.07360727 7.77867277 S
atom 7.75383535 0.03955871 2.04574490 S
atom -0.40503344 3.51818050 7.78795665 S
atom 1.45842019 3.72078621 2.00215409 S
atom 7.79025130 3.40908666 1.96086802 S
atom 5.76348242 3.67825786 7.69619500 S
atom 1.58540323 7.04787389 1.97146622 S
atom -0.38869634 0.13503162 7.70826663 S
atom 10.70254817 1.66592107 2.01150709 S
atom 8.63847067 5.41613898 7.91125167 S
atom 4.46784057 5.28109927 2.14237522 S
atom 2.48707371 1.86131413 7.84062310 S
atom 6.83602047 5.46203390 5.06869869 S
atom 4.81094794 1.61705589 10.73757863 S
atom 0.64076547 1.90317274 5.04532224 S
atom -1.33036998 5.28121448 10.76436048 S
atom 8.82279043 1.85768339 5.64643615 S
atom 6.74929161 5.22607908 11.35703279 S
atom 2.60901696 5.41338133 5.68524599 S
atom 0.60252908 1.75454358 11.43085201 S
atom 10.74113587 5.26223681 4.23085701 I
atom 8.82811361 1.85839711 9.81264579 I
atom 4.60875344 1.68351916 4.15700857 I
atom 2.61477744 5.47249390 9.86313130 I
The computational parameters for the FHI-aims calculations need to be written in
a control.in file. Create this new file and copy in the following:
xc pbesol
relativistic atomic_zora scalar
k_grid 1 2 1
compute_forces .true.
charge_mix_param 0.3
#===============================================================================
################################################################################
#
# FHI-aims code project
# VB, Fritz-Haber Institut, 2009
#
# Suggested "light" defaults for Li atom (to be pasted into control.in file)
# Be sure to double-check any results obtained with these settings for post-processing,
# e.g., with the "tight" defaults and larger basis sets.
#
################################################################################
species Li
# global species definitions
nucleus 3
mass 6.941
#
l_hartree 4
#
cut_pot 3.5 1.5 1.0
basis_dep_cutoff 1e-4
#
radial_base 29 5.0
radial_multiplier 1
angular_grids specified
division 0.4484 110
division 0.5659 194
division 0.6315 302
# division 0.6662 434
# division 0.8186 590
# division 0.9037 770
# division 6.2760 974
# outer_grid 974
outer_grid 302
################################################################################
#
# Definition of "minimal" basis
#
################################################################################
# valence basis states
valence 2 s 1.
# ion occupancy
ion_occ 1 s 2.
################################################################################
#
# Suggested additional basis functions. For production calculations,
# uncomment them one after another (the most important basis functions are
# listed first).
#
# Constructed for dimers: 1.80 A, 2.25 A, 2.75 A, 3.50 A, 4.50 A
#
################################################################################
# "First tier" - improvements: -189.23 meV to -6.35 meV
hydro 2 p 1.6
hydro 2 s 2
hydro 3 d 2.6
# "Second tier" - improvements: -4.69 meV to -0.41 meV
# hydro 3 p 4.6
# hydro 2 p 1.8
# hydro 3 s 6.2
# hydro 4 d 4.7
# hydro 4 f 4.1
# "Third tier" - improvements: -0.20 meV to -0.15 meV
# hydro 4 d 0.95
# hydro 3 p 6.2
# hydro 3 s 1.7
# Further functions, listed for completeness
# VB: The following functions are only listed for completeness; test very
# carefully before any kind of production use. From the point of view
# of the basis construction, their main role is to fill up space,
# and they are solely determined by the location of the cutoff potential.
# hydro 3 p 0.1 # -0.15 meV
# hydro 4 d 5 # -0.12 meV
# hydro 4 f 0.1 # -0.14 meV
# hydro 5 g 0.1 # -0.06 meV
#
#
################################################################################
#
# FHI-aims code project
# VB, Fritz-Haber Institut, 2010
#
# Suggested "light" defaults for P atom (to be pasted into control.in file)
# Be sure to double-check any results obtained with these settings for post-processing,
# e.g., with the "tight" defaults and larger basis sets.
#
################################################################################
species P
# global species definitions
nucleus 15
mass 30.973762
#
l_hartree 4
#
cut_pot 3.5 1.5 1.0
basis_dep_cutoff 1e-4
#
radial_base 43 5.0
radial_multiplier 1
angular_grids specified
division 0.5527 50
division 0.8801 110
division 1.1447 194
division 1.3165 302
# division 1.4113 434
# division 1.4781 590
# division 1.5482 770
# division 1.5845 974
# division 2.2606 1202
outer_grid 302
################################################################################
#
# Definition of "minimal" basis
#
################################################################################
# valence basis states
valence 3 s 2.
valence 3 p 3.
# ion occupancy
ion_occ 3 s 1.
ion_occ 3 p 2.
################################################################################
#
# Suggested additional basis functions. For production calculations,
# uncomment them one after another (the most important basis functions are
# listed first).
#
# Constructed for dimers: 1.625 A, 1.875 A, 2.5 A, 3.25 A, 4.0 A
#
################################################################################
# "First tier" - improvements: -726.20 meV to -35.91 meV
ionic 3 d auto
ionic 3 p auto
hydro 4 f 6.2
hydro 5 g 8.6
ionic 3 s auto
# "Second tier" - improvements: -16.02 meV to -1.71 meV
# hydro 4 d 6.2
# hydro 4 p 9.2
# hydro 5 f 9.8
# hydro 1 s 0.7
# hydro 5 g 13.2
# "Third tier" - improvements: -1.82 meV to -0.20 meV
# hydro 3 p 2.5
# hydro 4 d 6.4
# hydro 5 f 11.2
# hydro 2 s 1.5
# "Fourth tier" - improvements: -0.91 meV to -0.17 meV
# hydro 3 d 16.8
# hydro 5 g 18
# hydro 4 p 4.5
# hydro 3 s 2.1
# Further basis functions that fell out of the optimization - < -0.09 meV
# hydro 4 p 10.4
# hydro 4 d 17.6
# hydro 4 s 11.2
################################################################################
#
# FHI-aims code project
# VB, Fritz-Haber Institut, 2009
#
# Suggested "light" defaults for S atom (to be pasted into control.in file)
# Be sure to double-check any results obtained with these settings for post-processing,
# e.g., with the "tight" defaults and larger basis sets.
#
################################################################################
species S
# global species definitions
nucleus 16
mass 32.065
#
l_hartree 4
#
cut_pot 3.5 1.5 1.0
basis_dep_cutoff 1e-4
#
radial_base 44 5.0
radial_multiplier 1
angular_grids specified
division 0.4665 110
division 0.5810 194
division 0.7139 302
# division 0.8274 434
# division 0.9105 590
# division 1.0975 770
# division 1.2028 974
# outer_grid 974
outer_grid 302
################################################################################
#
# Definition of "minimal" basis
#
################################################################################
# valence basis states
valence 3 s 2.
valence 3 p 4.
# ion occupancy
ion_occ 3 s 1.
ion_occ 3 p 3.
################################################################################
#
# Suggested additional basis functions. For production calculations,
# uncomment them one after another (the most important basis functions are
# listed first).
#
# Constructed for dimers: 1.6 A, 1.9 A, 2.5 A, 3.25 A, 4.0 A
#
################################################################################
# "First tier" - improvements: -652.81 meV to -45.53 meV
ionic 3 d auto
hydro 2 p 1.8
hydro 4 f 7
ionic 3 s auto
# "Second tier" - improvements: -30.20 meV to -1.74 meV
# hydro 4 d 6.2
# hydro 5 g 10.8
# hydro 4 p 4.9
# hydro 5 f 10
# hydro 1 s 0.8
# "Third tier" - improvements: -1.04 meV to -0.20 meV
# hydro 3 d 3.9
# hydro 3 d 2.7
# hydro 5 g 12
# hydro 4 p 10.4
# hydro 5 f 12.4
# hydro 2 s 1.9
# "Fourth tier" - improvements: -0.35 meV to -0.06 meV
# hydro 4 d 10.4
# hydro 4 p 7.2
# hydro 4 d 10
# hydro 5 g 19.2
# hydro 4 s 12
#
################################################################################
#
# FHI-aims code project
# VB, Fritz-Haber Institut, 2009
#
# Suggested "light" defaults for I atom (to be pasted into control.in file)
# Be sure to double-check any results obtained with these settings for post-processing,
# e.g., with the "tight" defaults and larger basis sets.
#
################################################################################
species I
# global species definitions
nucleus 53
mass 126.90447
#
l_hartree 4
#
cut_pot 3.5 1.5 1.0
basis_dep_cutoff 1e-4
#
radial_base 64 5.0
radial_multiplier 1
angular_grids specified
division 0.1103 110
division 0.1515 194
division 0.9554 302
# division 1.1196 590
# division 1.1922 770
# division 6.1948 974
# outer_grid 974
outer_grid 302
################################################################################
#
# Definition of "minimal" basis
#
################################################################################
# valence basis states
valence 5 s 2.
valence 5 p 5.
valence 4 d 10.
# ion occupancy
ion_occ 5 s 1.
ion_occ 5 p 4.
ion_occ 4 d 10.
################################################################################
#
# Suggested additional basis functions. For production calculations,
# uncomment them one after another (the most important basis functions are
# listed first).
#
# Constructed for dimers: 2.22, 2.65, 3.25, 4.25 Ang
#
################################################################################
# "First tier" - improvements: -270.74 meV to -21.24 meV
hydro 3 d 4
hydro 4 f 6.4
hydro 2 p 1.6
ionic 5 s auto
# "Second tier" - improvements: -16.39 meV to -0.39 meV
# hydro 5 g 9.4
# hydro 4 f 18.4
# hydro 6 h 12.4
# hydro 4 p 4.5
# hydro 3 d 4.2
# hydro 3 s 3.0
# "Third tier" - max. impr. -0.76 meV, min. impr. -0.06 meV
# hydro 5 f 15.6
# hydro 5 g 12
# hydro 5 d 16
# hydro 4 f 42
# hydro 6 h 15.2
# ionic 5 p auto
# hydro 1 s 6.2
# Further functions that fell out of the optimization: -0.10 meV and below
# hydro 4 f 7
# hydro 6 p 9
# hydro 2 s 6.4
This includes the settings that the FHI-aims calculations will use. Here, we use the PBEsol functional with a 1x2x1 k-grid and light species defaults, that enables us to run the calculations quickly in the scope of this tutorial, but do not represent converged settings. The atomic positions of the initial and final geometry have both been optimized using these settings.
Important
It's important that the control.in file contains the keyword compute_forces .true., since
forces are not computed within FHI-aims by default.
Next we need to write the setup file for aimsChain. Create a new file chain.in
with the following content:
####################################################
## Mandatory settings
####################################################
# Specify the command used for single point calculations with FHI-aims
run_aims mpirun -np 192 aims.x
# Initial geometry
initial_file ini.in
# Final geometry
final_file fin.in
####################################################
## General settings for task
####################################################
# Method to use, neb/string
method string
# Number of images to use
n_images 4
# Convergence criterion for forces, in eV/A
force_thres 0.1
# Optimizer for neb/string method
optimizer dampedBFGS
####################################################
## External starting geometries
####################################################
# List of external geometries
external_geometry extgeo.lst
The first block includes the mandatory options that aimsChain requires, such as
the run command, the initial geometry, and the final geometry. Please adjust the run_aims command,
when running on a system that is different to the AWS cluster we use during the hands-on tutorial
The second block includes the settings related to the minimum energy path algorithm.
Here we use the string method with 4 images between the initial and final point,
and a force convergence criterion of 0.1 eV/Å. For actual production results, one would
use more intermediate images (around 8 or so) plus a tighter convergence threshold of
0.05 eV/Å or lower. Keep in mind that there is no one-stop solution or recipe
for every system, and that these settings need to be tested thoroughly.
For reference, the results in Ref. 1 have been obtained with 8 images and a force criterion of 0.02 eV/Å.
The last block specifies that we are using external geometries as a starting point for our
minimum energy path algorithm, which are specified by an external file extgeo.lst. The number of
external geometries should correspond to the n_images variable. But there is also an option to resample
the external path to obtain more or less images. For more details, please check the description of the resample keyword
for aimsChain in the FHI-aims manual. In general, it is recommended to use (if available) an initial guess for the path
that is not a linear interpolation. This could be a cheaper, less-converged initial aimsChain run or simply physical intuition.
Info
The visualization of the diffusion pathways has been performed with the OVITO package.
A simple linear interpolation between initial and final path looks like the following:
In this tutorial, we will start from a pre-converged path to reduce the amount of computational time required.
Create a new directory interpolation with
and create 4 files in that folder titled image002.in, image003.in, image004.in, and image005.in. Fill in the files with the
following blocks:
Starting point for the diffusion path
lattice_vector 12.3522719999999993 0.0000000000000000 -0.0279830000000000
lattice_vector 0.0000000000000000 7.1318460000000004 0.0000000000000000
lattice_vector -3.9878309999999999 0.0000000000000000 11.4389970000000005
atom 11.4328880792406586 3.8065689889565975 2.1485971368834877 Li
atom 8.7344857797458975 4.7225395827096355 10.2892171022455212 Li
atom 4.6001387319193787 5.8566581495631382 4.5073764056715353 Li
atom 2.5550040690167619 1.2858561626766636 10.1885028273993541 Li
atom 6.5765141649504182 -0.0164434947706597 10.0766553694972920 Li
atom 8.5830373209831485 7.1170753203363946 4.2768415849314074 Li
atom 0.4027567241330697 3.6700058794529533 10.0950502010670800 Li
atom 2.3164271785091044 3.4637518854450162 4.3775746743045669 Li
atom 8.5058733116778171 3.6392011556731312 4.2711765295564321 Li
atom 6.2514033745464710 3.4442611804487839 10.0709397329328727 Li
atom 2.0237209654738906 0.1295468775929066 4.3806012884475818 Li
atom 0.0665933471589836 7.0572648024177305 10.0584663326152448 Li
atom 8.4864956789753236 5.3406557689453225 0.5712538495427752 Li
atom 6.4727071212282912 1.8000424583486279 6.2318507479402632 Li
atom 2.3091141748963127 1.8111352132576395 0.5669680374015594 Li
atom 0.3601133453392645 5.3352471398754506 6.3424211972566384 Li
atom 11.5357483454467520 6.6184142408658930 1.7600020496448139 Li
atom 9.5918286149425462 0.4330320899557707 7.3968955886946954 Li
atom 5.3419318749711602 3.0920834468721359 1.7431752841826833 Li
atom 3.3770443522778928 4.0285793846694533 7.4282438008374578 Li
atom 9.6238463096465861 3.2867362900370605 7.2956960348271860 Li
atom 3.3584241832190629 6.8082202813221295 7.3886575441921094 Li
atom 5.2721359431160506 0.3246861005574877 1.6798370359123214 Li
atom 6.7670132001166614 5.3721294682589500 7.0691743376365812 P
atom 8.6835052787645441 1.6972328768592888 1.3685226141768356 P
atom 0.5600448569594485 1.8410014967442789 7.0448969504428955 P
atom 2.5320365121050079 5.3114699801100631 1.3941856359657354 P
atom 5.7817505361900166 7.0595692341300200 7.7242898501932125 S
atom 7.6490404069049038 0.0384234182081089 2.0101142740310873 S
atom -0.3792897896991582 3.5092215037010739 7.8110302507559135 S
atom 1.5458512546282737 3.6376147978624540 2.0571677975213474 S
atom 7.7376861919345510 3.4272630455852799 1.9611390384927518 S
atom 5.7812250609106561 3.6736243730301377 7.6932040292895003 S
atom 1.5814431911455638 7.0318419377104409 1.9972976746036477 S
atom -0.3748164756660740 0.1256868854167767 7.7083507341058546 S
atom 10.5116035501465888 1.6124827452201114 2.2786979700901795 S
atom 8.6534877331297686 5.4025132070912427 7.8971704174587805 S
atom 4.4539773921466024 5.2743288268472783 2.1317836797018375 S
atom 2.4916202182249672 1.8505912140714362 7.8009745938161856 S
atom 6.8961518185491810 5.3833824090904256 5.0259988690331960 S
atom 4.8472133472969077 1.6343188205672297 10.7550026050724341 S
atom 0.5333917294331738 1.9069756222024628 5.0184890314498070 S
atom -1.3575501023486471 5.2643815834131136 10.7794783158562399 S
atom 8.7898850286375563 1.8218953668910696 5.7081790797193648 S
atom 6.7975296386799329 5.1949637325465972 11.5161153532877023 S
atom 2.6273464526433301 5.3987049928537756 5.6806919133704747 S
atom 0.6080522526425148 1.7459031652326256 11.4123235943683063 S
atom 10.7955951919785047 5.3871419194691574 4.2725459057110342 I
atom 8.7935590394306313 1.8652866623693976 9.8358858274408281 I
atom 4.5880594637160943 1.6583205374556440 4.1514514995754848 I
atom 2.6392339839957812 5.4710243398064016 9.8532402248969877 I
lattice_vector 12.3522719999999993 0.0000000000000000 -0.0279830000000000
lattice_vector 0.0000000000000000 7.1318460000000004 0.0000000000000000
lattice_vector -3.9878309999999999 0.0000000000000000 11.4389970000000005
atom 11.5541796064167244 3.5211265357072548 2.8008821374543702 Li
atom 8.7323168568677740 4.7129234586982438 10.2987882972366940 Li
atom 4.5904390326243050 5.8623431439718647 4.5099484955282430 Li
atom 2.5308926435714998 1.2808256291941778 10.1724926602542354 Li
atom 6.5776402703856416 0.0252122275689895 10.0613435516409240 Li
atom 8.5299783701064555 7.1479899963962890 4.2799390633079950 Li
atom 0.3996207462911313 3.6723920619793780 10.0925807540926211 Li
atom 2.3459313500207410 3.4680645161962289 4.3651794072277239 Li
atom 8.4425363579219059 3.6530593158058258 4.2998152761651074 Li
atom 6.2552554736922641 3.4190970368137124 10.0630020911957381 Li
atom 2.0311398450080564 0.1338823429123098 4.3756147622544850 Li
atom 0.0564329507705924 7.0611938692349163 10.0583670378475816 Li
atom 8.4873523545881771 5.3330815378418359 0.5598232744558344 Li
atom 6.4864873149293460 1.7969496490289079 6.2450278323131077 Li
atom 2.3072107661517323 1.8174785882167270 0.5720540612692624 Li
atom 0.3703180757415313 5.3494096724609852 6.3579804751981719 Li
atom 11.5558522530019179 6.4930371240950278 1.7446058552620960 Li
atom 9.5832455068637739 0.4379396313177699 7.3882392566557114 Li
atom 5.3633112776808431 3.0950723862500591 1.7377313834172159 Li
atom 3.3757677230477090 4.0334633170805523 7.4274639157180395 Li
atom 9.6090246999230544 3.2842227133585875 7.3048478305799494 Li
atom 3.3537769989377280 6.8097125548444888 7.3878914603072268 Li
atom 5.2675849358085109 0.3234620153943970 1.6771115468420117 Li
atom 6.7590017673891065 5.3774859848890966 7.0801565242517448 P
atom 8.6871052645701745 1.6820094019724681 1.3430818733689338 P
atom 0.5511641995035702 1.8565771986993513 7.0382226719372358 P
atom 2.5314595194074503 5.3177812007048395 1.3980973950730256 P
atom 5.7714947885578445 7.0629451700272652 7.7293472261101002 S
atom 7.6443525372602066 0.0332777673975351 1.9994388829855294 S
atom -0.3914497610778513 3.5211933487696547 7.8012724083836211 S
atom 1.5340430627836898 3.6414272302626265 2.0461819190237907 S
atom 7.7724607810029136 3.4181039326884646 1.9681320000255413 S
atom 5.7755171718616944 3.6746674778532871 7.6899757791996066 S
atom 1.5761041507539217 7.0331882286961473 1.9980743863223722 S
atom -0.3805060556734666 0.1416112894530480 7.7000492012566575 S
atom 10.5476819902086092 1.5199207768329948 2.1627092144947859 S
atom 8.6502076694015706 5.4046454805436248 7.8970819587274947 S
atom 4.4508253033266492 5.2708556338102506 2.1386406620291458 S
atom 2.4826478053532774 1.8562884758261848 7.7847647670816471 S
atom 6.8800724559429582 5.4010111236628964 5.0329366025468980 S
atom 4.8122918333356246 1.6335157962837736 10.7284273660948148 S
atom 0.5459587147200216 1.9183261095811055 5.0003412610848459 S
atom -1.3576564312114352 5.2630378867491237 10.7845062770060380 S
atom 8.7913919910320768 1.8359270078714092 5.7025869301680103 S
atom 6.7666296523795202 5.2108479302439550 11.4347748568672785 S
atom 2.6208589843256540 5.4011885028786830 5.6766971417470096 S
atom 0.6020975837130954 1.7498716313880558 11.4135312345924333 S
atom 10.7510298783348404 5.4793007048940510 4.3344402985222228 I
atom 8.7983519317774377 1.8600005552830259 9.8298108347492530 I
atom 4.5966460852178121 1.6610801294787783 4.1501842838791410 I
atom 2.6241150010803818 5.4697764425114554 9.8495128243359655 I
lattice_vector 12.3522719999999993 0.0000000000000000 -0.0279830000000000
lattice_vector 0.0000000000000000 7.1318460000000004 0.0000000000000000
lattice_vector -3.9878309999999999 0.0000000000000000 11.4389970000000005
atom 11.4791139575449765 3.1674707956222390 3.4409096448920535 Li
atom 8.7370839556586635 4.7073441127169584 10.3049763873320348 Li
atom 4.5836522056053601 5.8658120456495064 4.5085907642005303 Li
atom 2.5117556081250285 1.2756942595583112 10.1650074511544393 Li
atom 6.5788684762241179 0.0666596183253879 10.0530737416868359 Li
atom 8.5080775524724803 7.1521009047249766 4.2975735232301826 Li
atom 0.3968071752782256 3.6770502819031425 10.0966939525869392 Li
atom 2.3660125440941551 3.4859080282592734 4.3318465334409195 Li
atom 8.3830432523047236 3.6697655161117613 4.3092255410911795 Li
atom 6.2579116447842722 3.3972043761508250 10.0551179827486479 Li
atom 2.0474111369493952 0.1202214919432862 4.3540649841269463 Li
atom 0.0478010444689524 7.0622533969448202 10.0658182870303978 Li
atom 8.4884267946848606 5.3308904778381914 0.5578986521545202 Li
atom 6.5144706666270951 1.7956702407897533 6.2509519553056023 Li
atom 2.3047057406867970 1.8351432167606303 0.5904453440460866 Li
atom 0.3667701169947584 5.3643467179878295 6.3495783121042848 Li
atom 11.5681151976459926 6.4200910012872328 1.7353222679205904 Li
atom 9.5710601969918141 0.4457529708198725 7.3730799943370133 Li
atom 5.3759352661875504 3.0988966152500934 1.7343148148296481 Li
atom 3.3763932860830028 4.0368865597330617 7.4268054671703752 Li
atom 9.5899878360286941 3.2819445485165986 7.3148959143829213 Li
atom 3.3537349441692128 6.8090831544812467 7.3888747174632234 Li
atom 5.2812941376358369 0.3217052926900954 1.6772468468400055 Li
atom 6.7542497070991372 5.3829253277172606 7.0883500306059561 P
atom 8.7082244908148372 1.6883509585618035 1.3184906310426927 P
atom 0.5450053467414580 1.8623960778435926 7.0565336931777560 P
atom 2.5235657864546956 5.3101787764944639 1.4054911354198341 P
atom 5.7666803757371596 7.0672048492213104 7.7330552126930003 S
atom 7.6771954958748223 0.0420829028843274 2.0046994731572387 S
atom -0.4022907514055987 3.5321445839786030 7.7909786181045018 S
atom 1.4713609009491462 3.6690514049048089 2.0363534447078404 S
atom 7.7926245861063723 3.4157520661053011 1.9668902340535201 S
atom 5.7737778120840586 3.6757752945243318 7.6850068327459260 S
atom 1.5813174327370760 7.0386991843057283 1.9888348398480973 S
atom -0.3905742235112833 0.1492096546938267 7.7008386030909790 S
atom 10.6307424174320193 1.5454397202666856 1.9797125813896603 S
atom 8.6498103327634528 5.4091623112563489 7.8947950976486121 S
atom 4.4500872329047141 5.2664744534795291 2.1440822839599658 S
atom 2.4765470449131581 1.8589436120279585 7.7818131103020143 S
atom 6.8680375009837302 5.4197933215183873 5.0378422334451169 S
atom 4.7777589202679493 1.6361814026285901 10.7004037643332470 S
atom 0.5907026352170606 1.9055464035017555 5.0076152301509076 S
atom -1.3535651350133744 5.2602360130959536 10.7875341121703485 S
atom 8.8173515661936275 1.8573096384349432 5.6754344434454547 S
atom 6.7437270730356396 5.2303807031871914 11.3618709983884045 S
atom 2.6154398693185583 5.4020131763731083 5.6766018109159528 S
atom 0.5989343748767280 1.7523210435455665 11.4159849291704578 S
atom 10.7476111053405408 5.5217094452730633 4.3076667715531789 I
atom 8.8037804394517352 1.8592692853198067 9.8200349081463898 I
atom 4.6037088826168278 1.6646663718033483 4.1512512519383460 I
atom 2.6164593474218423 5.4702051303329453 9.8477698631427675 I
lattice_vector 12.3522719999999993 0.0000000000000000 -0.0279830000000000
lattice_vector 0.0000000000000000 7.1318460000000004 0.0000000000000000
lattice_vector -3.9878309999999999 0.0000000000000000 11.4389970000000005
atom 10.9951366103271813 2.7585388768240011 3.9871501464710994 Li
atom 8.7495495358060875 4.7027471441394173 10.3134174631073439 Li
atom 4.5761218818000104 5.8741579029681033 4.5055549429177058 Li
atom 2.5142316813285479 1.2734513429763117 10.1778984474914456 Li
atom 6.5878811561323047 0.0844468264702590 10.0572521400649482 Li
atom 8.5462153354587151 7.1232370934136506 4.3101680158348463 Li
atom 0.3974883312853609 3.6800830189174567 10.0960563184056902 Li
atom 2.3498609945856366 3.5112252054513160 4.3098369536782020 Li
atom 8.3031201696014083 3.6983518226177230 4.3022018511438036 Li
atom 6.2637992303570194 3.3857203721432776 10.0541257163846449 Li
atom 2.0555791781039616 0.1106609053992260 4.3401205144873511 Li
atom 0.0475289319109858 7.0709432450419936 10.0659516330531815 Li
atom 8.4907896494814725 5.3290370278764030 0.5629503351628176 Li
atom 6.5693180335273276 1.7946029709925908 6.2467383641453464 Li
atom 2.3055804345610045 1.8445336244882575 0.5966507043515331 Li
atom 0.3553031397032321 5.3654945894838821 6.3287740214447314 Li
atom 11.5675283335422066 6.4247402962922795 1.7380795036404739 Li
atom 9.5718848845067352 0.4558711361908368 7.3513368034892324 Li
atom 5.3816001649703793 3.1044891791721971 1.7316017803877954 Li
atom 3.3744631587828327 4.0363499754739269 7.4275870137582762 Li
atom 9.5819531370631932 3.2824616112164646 7.3237224413831976 Li
atom 3.3548174413046015 6.8120267337595157 7.3931180652939172 Li
atom 5.3059735096613290 0.3273713843241502 1.6825622859151959 Li
atom 6.7510561936608768 5.3899767254505457 7.0990954476951282 P
atom 8.7394619629224071 1.6977051767565952 1.3200739950691334 P
atom 0.5561867933830994 1.8536513918787272 7.0736227821948052 P
atom 2.5397866971109675 5.3134303529502880 1.3911064533999791 P
atom 5.7680005717901137 7.0745290087258432 7.7432748010234205 S
atom 7.7173840306309787 0.0469406884682152 2.0189016137120221 S
atom -0.4010336653949791 3.5279428080399065 7.7868227464184079 S
atom 1.4689312390900879 3.6984945770048387 2.0162602831095486 S
atom 7.7941503307922675 3.4138919257708409 1.9581720531785867 S
atom 5.7771290618467734 3.6757820865385318 7.6811508647671358 S
atom 1.5913655612988526 7.0422983989056593 1.9776694983691159 S
atom -0.3924772347306327 0.1430824007084935 7.7051653712988335 S
atom 10.6836026028667703 1.6202275608137389 1.9277397865023926 S
atom 8.6508180199045732 5.4164778609216491 7.8944620459157280 S
atom 4.4653816528338428 5.2744599323483401 2.1407519027690851 S
atom 2.4819707972030920 1.8583821834005061 7.8059026773005513 S
atom 6.8557674181705019 5.4501706358838522 5.0467584380515076 S
atom 4.7783639703490977 1.6315018134273751 10.7001346802145143 S
atom 0.6291007611323032 1.9023419869881701 5.0267007720014414 S
atom -1.3395461866913088 5.2732251354807218 10.7713378293975666 S
atom 8.8956348049328291 1.8894664795208376 5.6054737704225319 S
atom 6.7444987445009161 5.2313314443357442 11.3501233068814553 S
atom 2.6159789948395575 5.4064808324336351 5.6809138570871225 S
atom 0.6003599787743850 1.7546176417262862 11.4204466806324767 S
atom 10.7603982097238742 5.4312244027009777 4.2502932297981166 I
atom 8.8135259167907645 1.8606329110080595 9.8064248139242700 I
atom 4.6118413164369025 1.6739674474587247 4.1545824509195617 I
atom 2.6192664692913055 5.4710996746539404 9.8525031118673621 I
Then create the file extgeo.lst, which specifies that aimsChain uses these geometries:
After creating all these files, your folder should now look like the following:
Li6PS5I_diffusion/
├── chain.in
├── control.in
├── extgeo.lst
├── fin.in
├── ini.in
│
├── interpolation/
│ └── image002.in
│ └── image003.in
│ └── image004.in
│ └── image005.in
With these files in place, we are now ready to run the transition state search using aimsChain. To do so on the provided AWS resources during the hands-on tutorial, create the following submission script:
#!/bin/bash -l
# Standard output and error:
#SBATCH -o ./Output
#SBATCH -e ./Error
# Initial working directory:
#SBATCH -D ./
# Job Name:
#SBATCH -J Jobscript
# Number of nodes and MPI tasks per node:
#SBATCH --nodes=1
#SBATCH --ntasks-per-node=192
#
# Wall clock limit:
#SBATCH --time=00:10:00
module purge
module load fhi-aims
ulimit -s unlimited
export OMP_NUM_THREADS=1
export MKL_NUM_THREADS=1
runchain
Running calculations on other resources
If you are running the calculations on a different machine, make sure to adjust the submission script specific to your system and submission system.
Then submit the job with:
Analyzing the results
When the calculations completed successfully, you should see a folder structure like the following:
Li6PS5I_diffusion/
├── chain.in
├── control.in
├── extgeo.lst
├── fin.in
├── forces.log
├── ini.in
│
├── interpolation/
├── iterations/
├── optimized/
├── path/
│
Inspect the forces.log file which contains the residual forces in eV/Å to make sure that the calculation finished successfully:
Provided analysis scripts
In the following, we will provide python scripts to do the data analysis. We do this to ensure that the tutorial can be completed by any user. Writing scripts (in whatever programming language) is something any user should try to do. There are many good python resources available online, such as the one provided by MolSSI on Python Scripting for Computational Molecular Science.
Next we can plot the optimized energy path. Create a python script plot_energy_path.py with the content from:
Plotting script
"""
Plot energy profiles from NEB (Nudged Elastic Band) or string method calculations from AimsChain.
This script reads energy profiles from one or more calculation directories and creates
a comparative plot. For each directory, it reads:
- ener.lst: Contains energy values and point types (FIXED/NORMAL)
- path.lst: Contains reaction coordinates for each point
- path.xyz: Contains atomic coordinates for calculating reaction coordinates
Usage:
python plot_energy_path.py path/to/calculation1 [path/to/calculation2 ...]
python plot_energy_path.py path/to/calculation1 -o custom_name.png
python plot_energy_path.py calc1 calc2 calc3 -l "Path 1" "Path 2" "Path 3"
Options:
-o, --output Specify custom output filename (default: energy_path.png)
-l, --labels Optional labels for the paths (must match number of folders)
"""
import numpy as np
import matplotlib.pyplot as plt
import argparse
import os
from scipy.interpolate import CubicSpline
def parse_arguments():
parser = argparse.ArgumentParser(description="Plot NEB/String energy profile")
parser.add_argument(
"folders",
type=str,
nargs="+",
help="Paths to folders containing ener.lst and path.lst",
)
parser.add_argument(
"--output",
"-o",
type=str,
default="energy_path.png",
help="Output filename (default: energy_path.png)",
)
parser.add_argument(
"--labels",
"-l",
type=str,
nargs="+",
help="Optional labels for the paths in the legend",
)
args = parser.parse_args()
if args.labels and len(args.labels) != len(args.folders):
parser.error(
f"Number of labels ({len(args.labels)}) must match number of folders ({len(args.folders)})"
)
return args
def read_xyz_trajectory(filename):
"""Read atomic coordinates from an xyz trajectory file."""
trajectory = []
n_atoms = None
current_frame = []
with open(filename, "r") as f:
for line in f:
if n_atoms is None:
n_atoms = int(line.strip())
continue
if "Energy:" in line:
if current_frame:
trajectory.append(np.array(current_frame))
current_frame = []
continue
if len(line.split()) == 4: # Atom line: symbol x y z
x, y, z = map(float, line.split()[1:])
current_frame.append([x, y, z])
if current_frame:
trajectory.append(np.array(current_frame))
return trajectory
def calculate_reaction_coordinate(trajectory):
"""Calculate normalized reaction coordinates based on RMSD between images."""
n_frames = len(trajectory)
distances = np.zeros(n_frames)
# Calculate RMSD between consecutive images
for i in range(1, n_frames):
diff = trajectory[i] - trajectory[i - 1]
rmsd = np.sqrt(np.mean(np.sum(diff**2, axis=1)))
distances[i] = rmsd
# Calculate cumulative distance
cumulative = np.cumsum(distances)
# Normalize to [0,1] range
reaction_coordinates = cumulative / cumulative[-1]
return reaction_coordinates
def read_folder_data(folder_path):
"""Read energy data and fixed points from a folder."""
ener_file = os.path.join(folder_path, "ener.lst")
if not os.path.exists(ener_file):
raise FileNotFoundError(f"Could not find ener.lst in {folder_path}")
# Read the energy data
energies = []
fixed_points = []
with open(ener_file, "r") as f:
next(f) # Skip header
for line in f:
columns = line.split()
energies.append(float(columns[1]))
fixed_points.append(columns[2].strip() == "FIXED")
energies = np.array(energies)
fixed_points = np.array(fixed_points)
# Align energies to minimum
min_energy = np.min(energies)
relative_energies = energies - min_energy
return relative_energies, fixed_points, min_energy
def plot_energy_profiles(folders, output_file, labels=None):
"""Create energy profile plot with smooth interpolation."""
plt.figure(figsize=(10, 6))
colors = plt.cm.tab10(np.linspace(0, 1, len(folders)))
all_coordinates = []
all_energies = []
for i, (folder_path, color) in enumerate(zip(folders, colors)):
try:
# Read energy data
energies, fixed, min_energy = read_folder_data(folder_path)
# Read trajectory and calculate reaction coordinate
xyz_file = os.path.join(folder_path, "path.xyz")
trajectory = read_xyz_trajectory(xyz_file)
reaction_coordinate = calculate_reaction_coordinate(trajectory)
# Create spline interpolation
cs = CubicSpline(reaction_coordinate, energies)
smooth_coords = np.linspace(0, 1, 200)
smooth_energies = cs(smooth_coords)
plt.plot(
smooth_coords,
smooth_energies,
"-",
color=color,
linewidth=2,
label=labels[i] if labels and i < len(labels) else folder_path,
)
plt.plot(
reaction_coordinate,
energies,
"o",
color=color,
markersize=6,
markerfacecolor=color,
)
# Add fixed points with different marker
fixed_coords = reaction_coordinate[fixed]
fixed_energies = energies[fixed]
plt.plot(
fixed_coords,
fixed_energies,
"s",
color=color,
markersize=8,
markeredgecolor="black",
markerfacecolor=color,
)
all_coordinates.extend(reaction_coordinate)
all_energies.extend(energies)
print(f"\nEnergy profile summary for {folder_path}:")
print(f"Minimum energy: {min_energy:.6f} eV")
print(f"Maximum barrier: {np.max(energies):.6f} eV")
except Exception as e:
print(f"Error processing {folder_path}: {e}")
continue
plt.xlabel("Reaction Coordinate", fontsize=15)
plt.ylabel("Relative Energy [eV]", fontsize=15)
plt.grid(True, linestyle="--", alpha=0.7)
plt.legend(loc="upper right", fontsize=15)
# Set axis limits with padding
x_padding = 0.05
y_padding = (max(all_energies) - min(all_energies)) * 0.1
plt.xlim(-x_padding, 1 + x_padding)
plt.ylim(min(all_energies) - y_padding, max(all_energies) + y_padding)
plt.savefig(output_file, dpi=300, bbox_inches="tight")
print(f"\nPlot saved as: {output_file}")
plt.close()
def main():
args = parse_arguments()
try:
plot_energy_profiles(args.folders, args.output, args.labels)
except Exception as e:
print(f"Error: {e}")
return 1
return 0
if __name__ == "__main__":
main()
You can then visualize the results with:
This will take the results from the optimized folder and output the following in your terminal:
> python plot_energy_path.py optimized
Energy profile summary for optimized:
Minimum energy: 0.000000 eV
Maximum barrier: 0.406367 eV
Plot saved as: energy_path.png
and save the following plot in energy_path.png:
We can see that the optimized path has an energy barrier of around 0.4 eV, which agrees well with the converged result from Ref. 1 and also experimental measurements of the room-temperature Li-migration barrier of 0.39 eV 2. Please note that we are not including temperature effects in our simulations here, so we would normally expect a higher energy barrier than room-temperature experiments.
We can also compare the result to a simple linear interpolation between initial and final geometries (lightblue line) which we show here (but did not compute explicitly in the interest of time):
It's clear from this picture that a linear interpolation is not the optimal path, with the energy barrier being more than 1.5 eV. This is also apparent when checking the optimized path (black) compared to the linear path (grey):
Additional steps to properly evaluate the transition state
After obtaining the initial reaction pathway shown in our energy profile, several key steps are omitted in this tutorial, but should follow to properly characterize the transition state (if desired). The climbing image method can be used next, where one image moves uphill along the reaction path to locate an approximate transition state (this can be done with aimsChain as described in the FHI-aims manual). This transition state can then be refined using a dedicated saddle point search algorithm. Finally, vibrational frequency calculations should be performed to verify that the structure is a true transition state, characterized by exactly one imaginary frequency.
Comparison between computational methods
While the previous PBEsol calculation already shows good agreement with experimental data, we want to compare the results with hybrid functionals and different dispersion corrections. The following plot is taken from Fig. 6 of Ref. 1 and compares the same diffusion path we have just calculated for different methods.
Three different approaches are compared:
- PBEsol (yellow line) - the semi-local density functional
- HSE06+TS (pink line) - hybrid functional with the Tkatchenko-Scheffler (pairwise van der Waals) correction
- HSE06+MBD-NL (green line) - hybrid functional with the non-local many-body dispersion method
The inclusion of non-local many-body dispersion effects produces the highest barrier (~0.45 eV), which aligns better with experimental measurements. This higher barrier arises because the MBD-NL method captures important destabilizing effects at the transition state that are missed by simpler methods like semi-local functionals. The pairwise van der Waals correction (TS, dotted pink line) appears to over-stabilize the transition state, leading to artificially lower barriers (~0.25 eV). While the TS method shows significant variation for the contribution from van der Waals interactions along the path, the many-body treatment maintains a more consistent contribution for all nodes.
The higher accuracy of HSE06+MBD-NL stems from its ability to capture the complex interplay between the moving Li ion and the surrounding framework, including effects like electronic polarization that are particularly important at the transition state 1. Similar results were obtained in the other investigated ordered and disordered argyrodite systems.
In case you are interested in learning more about the use of different dispersion corrections in FHI-aims, you can also look at our tutorial on Efficient dispersion-corrected hybrid DFT in FHI-aims.
Solutions
You find the solutions to the above exercise by clicking on the button below.