How to Obtain Partial Correlation Function from XDATCAR

[rmz4ed]

Hi,

I am trying to obtain the partial correlation function after running a MD simulation in VASP. In an older thread, it was suggested to use the XDATCAR file to obtain the partial function. However, as the XDATCAR VASP wiki page is very sparse, I am wondering if anyone has any advice on how to use the XDATCAR file to calculate the partial PCFs, or if there are any existing scripts for doing so.

Thanks!

[marie-therese.huebsch]

Hi,
so roughly speaking the pair-correlation function is the probability of finding a particle at a given distance from the center of another particle. Formally, this can be written as

Code: Select all

g(r) = (V/N^2) < \sum_i \sum_{j\neq i} \delta( \vec{r} - \vec{R}_{ij} ) >,

where V is the real space volume, N is the number of particles, \vec{R}_{ij}= \vec{R}_i - \vec{R}_{j} is the distance of particles i and j, \delta is the delta function and < > takes the ensemble average. (Sorry there is no math mode) What you need from your MD run is thus the structure data and the ionic position at all time steps.

At the beginning of XDATCAR, you can find the structure parameters in the same format as in the POSCAR file. Then, for each MD step the ionic positions are listed

Code: Select all

Direct configuration=     <no of MD step>

"Direct" means that fractional coordinates and not cartesian coordinates are written. The order in which the atomic species are listed is again the same as in the POSCAR file. With the definition above and the data in XDATCAR you can write a script in your preferred language. I am sorry I don't have one at hand right now.

In case you are running VASP version 6.3.0 with HDF5 support and py4vasp version 0.4 is installed, you can also use the following python code:

Code: Select all

from py4vasp import Calculation
import mdtraj
calc = Calculation.from_path(pathname)
traj = calc.structure[:].to_mdtraj()
mdtraj.compute_rdf(traj, other_settings)

However, currently this approach has some issues: mdtraj assumes a certain standard format. The first lattice vector must be along x, the second can have x and y components and the third is arbitrary. For orthorhombic cells that should work, but please be aware that for others it might fail. Also mind that mdtraj follows a minimum image convention. So all periodic images of an ion are ignored.

You can also use py4vasp to read the data and then implement your own code to compute the partial paircorrelation function. That way you are fully aware of all conventions that have been used. I hope this was helpful. Let me know if there are further questions and please consider sharing your script if you write one.

Best regards,
Marie-Therese