This notebook demonstrates how to perform diffusivity and ionic conductivity analyses starting from a series of VASP AIMD simulations using Python Materials Genomics (pymatgen) and its add-on package pymatgen-diffusion. These notebooks are described in detail in
Deng, Z.; Zhu, Z.; Chu, I.-H.; Ong, S. P. Data-Driven First-Principles Methods for the Study and Design of
Alkali Superionic Conductors. Chem. Mater. 2017, 29 (1), 281–288 DOI: 10.1021/acs.chemmater.6b02648.
If you find these notebooks useful and use the functionality demonstrated, please consider citing the above work.
Let's start by importing some modules and classes that we will be using.
In [1]:
from IPython.display import Image
%matplotlib inline
import matplotlib.pyplot as plt
import json
import collections
from pymatgen import Structure
from pymatgen.analysis.diffusion_analyzer import DiffusionAnalyzer, \
get_arrhenius_plot, get_extrapolated_conductivity
from pymatgen_diffusion.aimd.pathway import ProbabilityDensityAnalysis
from pymatgen_diffusion.aimd.van_hove import VanHoveAnalysis
In [2]:
# files = ["run1/vasprun.xml", "run2/vasprun.xml", "run3/vasprun.xml"]
# analyzer = DiffusionAnalyzer.from_files(files, specie="Li", smoothed=False)
In this work, all trajectories are stored in an efficient document-based MongoDB database. The format of the documents in the database is a binary JSON format. Here, we will instead instantiate the DiffusionAnalyzer from a pre-serialized DiffusionAnalyzer for each temperature.
In [2]:
temperatures = [600, 800, 1000, 1200]
analyzers = collections.OrderedDict()
for temp in temperatures:
with open("aimd_data/%d.json" % temp) as f:
d = json.load(f)
analyzers[temp] = DiffusionAnalyzer.from_dict(d)
In [3]:
plt = analyzers[1000].get_msd_plot()
title = plt.title("1000K", fontsize=24)
In [4]:
diffusivities = [d.diffusivity for d in analyzers.values()]
plt = get_arrhenius_plot(temperatures, diffusivities)
From the temperatures and diffusivities, one may obtained the extrapolated room-temperature conductivity as follows.
In [5]:
rts = get_extrapolated_conductivity(temperatures, diffusivities,
new_temp=300, structure=analyzers[800].structure,
species="Li")
print("The Li ionic conductivity for Li6PS5Cl at 300 K is %.4f mS/cm" % rts)
We can compute the probability density function from the AIMD trajectories using the ProbabilityDensityAnalysis class implemented in the pymatgen-diffusion add-on. We will use the calculation at 800K as an example. The probability density function can then be output to a CHGCAR-like file for visualization in VESTA.
In [6]:
structure = analyzers[800].structure
trajectories = [s.frac_coords for s in analyzers[800].get_drift_corrected_structures()]
pda = ProbabilityDensityAnalysis(structure, trajectories, species="Li")
pda.to_chgcar("aimd_data/CHGCAR.vasp") # Output to a CHGCAR-like file for visualization in VESTA.
The VESTA visualization software can be used to visualize isosurfaces in the probability density. The 800K probability density function at an isosurface of 0.002 is shown below.
In [8]:
Image(filename='Isosurface_800K_0.png')
Out[8]:
In [7]:
vha = VanHoveAnalysis(analyzers[800])
We can then plot the self ($G_s$) and distinct ($G_d$) parts of the van Hove correlation function as follows.
In [8]:
vha.get_3d_plot(type="self")
vha.get_3d_plot(type="distinct")
Out[8]:
In [ ]: