Finding the MSD uncertainties from a VASP file#
Using kinisi
to obtain the mean-squared displacement and uncertainty in a VASP Xdatcar type file is straightforward and involves using the DiffusionAnalyzer
class.
[1]:
import numpy as np
from kinisi.analyze import DiffusionAnalyzer
np.random.seed(1)
There the params
dictionary describes details about the simulation, and are fully documented in the parser module. The 'specie'
is the atomic/ionic species that you want to calculate the MSD for, the 'time_step'
is the simulation timestep in your molecular dynamics, and the 'step_skip'
is the frequency which with the data was written in your simulation trajectory.
[2]:
params = {'specie': 'Li',
'time_step': 2.0,
'step_skip': 50,
'progress': False
}
u_params = {'progress': False}
Therefore, for this example, we have a simulation that had a timestep of 2 femtoseconds (note the FAQs about units) but we only stored the trajectory information every 50 steps and we want to investigate only the lithium motion.
[3]:
msd = DiffusionAnalyzer.from_file('example_XDATCAR.gz', parser_params=params, uncertainty_params=u_params)
/home/docs/checkouts/readthedocs.org/user_builds/kinisi/envs/latest/lib/python3.12/site-packages/pymatgen/io/vasp/outputs.py:4496: EncodingWarning: We strongly encourage explicit `encoding`, and we would use UTF-8 by default as per PEP 686
file_len = sum(1 for _ in zopen(filename, mode="rt"))
/home/docs/checkouts/readthedocs.org/user_builds/kinisi/envs/latest/lib/python3.12/site-packages/pymatgen/io/vasp/outputs.py:4499: EncodingWarning: We strongly encourage explicit `encoding`, and we would use UTF-8 by default as per PEP 686
with zopen(filename, mode="rt") as file:
/home/docs/checkouts/readthedocs.org/user_builds/kinisi/envs/latest/lib/python3.12/site-packages/kinisi/parser.py:116: UserWarning: Converting triclinic cell to orthorhombic this may have unexpected results. Triclinic a, b, c are not equilivalent to orthorhombic x, y, z.
warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/kinisi/envs/latest/lib/python3.12/site-packages/kinisi/parser.py:437: UserWarning: UserWarning: Be aware that the expected input unit for 'time_step' is femtoseconds, not picoseconds.
warnings.warn("UserWarning: Be aware that the expected input unit for 'time_step' is femtoseconds, not picoseconds.")
The DiffusionAnalyzer
will perform the bootstrapping process to obtain the displacements and uncertainties (this is detailed in the methodology. To find out how to determine the diffusion coefficient with DiffusionAnalyzer
continue here.
Then the MSD (msd
) as a function of timestep (dt
) can be plotted.
[4]:
import matplotlib.pyplot as plt
[5]:
plt.errorbar(msd.dt, msd.msd, msd.msd_std)
plt.ylabel('MSD/Å$^2$')
plt.xlabel(r'$\Delta t$/ps')
plt.show()