diffusion coefficient from rmsd or velocity autocorrelation function

[micaela-matta]

I think it would be great to implement the calculation of diffusion coefficients using the Green-Kubo and/or Einstein formulas.

This feature is among gromacs analysis tools and implemented in cpptraj/pytraj [ however, the latter is not well documented].

[mattwthompson]

I can take a stab at doing the Einstein half of this. I don't do much Green-Kubo and I'm not familiar with much more than the basics of it.

Somewhat related, ionic conductivity falls out of MSDs fairly easily either with Nernst-Einstein or Einstein-Helfand. Not sure this would be useful to tons of people but I have been using it lately.

[richardjgowers]

@mattwthompson this would be very welcome! One thing I'm aware of is this: https://github.com/pdebuyl-lab/tidynamics

[kain88-de]

I worked a lot on diffusion in the last three years. Below are some general comments to help implement correct calculation of diffusion coefficients.

---

tidydynamics and how to calculate correlations

The blog post is interesting. It shows that FFT is a good method for calculation of the MSD. However, it means we have to store all coordinates of the trajectory in a numpy array. With the additional requirement of zero-padding for performance, it can be a memory hog. The blocking method is much more memory friendly but it requires to iterate over the trajectory multiple times. The later might not have a huge impact if we calculate the MSD of several molecules at a time.

---

pydiffusion rotational and translational diffusion coefficients. The rotational diffusion module is state of the art but requires preprocessing the trajectories to remove all translational diffusion artifacts. This is true for all rotational diffusion analytic codes. The translational diffusion calculations are similar to tidydynamics. However, I have included a method to determine the translational diffusion in an arbitrary protein fixed coordinate system. This is not interesting for everyone. I only used it to check an implementation of a diffusion algorithm.

A previous PR

We had a PR for the addition of MSD calculations a few years ago. The algorithms looked quite interesting. I never fully understood it but it looked like it can to a blocked correlation analysis with a single pass over the trajectory. It's worth having a look into it.

---

Diffusion coefficients in Periodic Boundary Systems

The Diffusion coefficient is heavily influenced by PBC that we often use in simulations. Any implementation of a diffusion coefficient calculation in MDAnalysis should include PBC corrections. Luckily the corrections only depend on the simulation box size.

membranes translation rotation

For membranes in particular a huge warning sign should be added because for flat boxes the observed diffusion coefficient can grow indefinitely.

---

Some general remarks about the calculation of diffusion coefficients.

Diffusion coefficients are calculated from correlation functions. The MSD is a correlation function. Unfortunately, correlation functions are notoriously difficult to calculate. To reduce sampling noise it should be possible to average the MSD from several simulations or from different proteins of the same simulation.

[kain88-de]

If you have questions about some details you can ping me again on this issue.

[micaela-matta]

thanks a lot. I am only starting to dip my feet in diffusion-related problems.

[ranadevyesh]

I can take a stab at doing the Einstein half of this. I don't do much Green-Kubo and I'm not familiar with much more than the basics of it.

Somewhat related, ionic conductivity falls out of MSDs fairly easily either with Nernst-Einstein or Einstein-Helfand. Not sure this would be useful to tons of people but I have been using it lately.

Yes conductivity would be very useful to me. How have you implemented this calculation?

[orbeckst]

The current development version (and upcoming 2.0.0) has analysis.msd (does not do Green-Kubo, though).

[kain88-de]

Note that there have recently been papers about the systematic errors in self-diffusion coefficient calculation. The authors also provide a python package to remove the errors. The package is GPL3 licensed, but I'm sure if we ask nicely we can include it as GPL2 in MDAnalysis.

[xhgchen]

We now have a working implementation of Green-Kubo self-diffusivity at https://github.com/MDAnalysis/transport-analysis

I will update this comment when the package is officially released, but if anyone needs it before then, I would encourage cloning the repo or using the source code. The package uses standard MDAnalysis units. The docs can be found here and the reproduction blog post can be found here. Please note that the analysis I performed in that blog post did not use standard MDAnalysis units because I was using OpenMM for the simulation and had to enter the velocities manually as a result.

Thanks to everyone at MDAnalysis for working with me and making this possible!