Hydrogen Bond Autocorrelation Functions

User discussions
1

GROMACS version: 2024
GROMACS modification: No

Hello everyone!

I recently found out that GROMACS can compute the forward and backward constants of the theory of Luzar and Chandler. I have a simulation box containing 1414 molecules of water (TIP4P) and I did gmx hbond-legacy -f file.xtc -s file.tpr -ac ac.xvg. The aucotorrelation function is ok, I compared this with a Python code that I have for this. However, the constant values are weird:

image

This is for 298 K. According to Luzar and Chandler, k = 0.6 and k’ = 0.9. Why don’t I have these values? The simulation length is of 50 ps, maybe it should be longer?

Also, where can I find the code for the calculation of these constants?

Edit: I am also calculating it for simulations of 500 ps, but still get weird results (negative values for example).

Thank you very much!
Best,
Maria

2

There are a few factors that could be relevant.

  1. Do you know the expected values for the water model you are using?
  2. How often are you writing the simulation trajectory output? In Resolving the hydrogen bond dynamics conundrum | The Journal of Chemical Physics | AIP Publishing Luzar writes “This is why it is important to use sampling intervals well below 10 fs.” I haven’t checked if his analyses are the same, but it’s quite possible that you should have a very high trajectory output frequency. He also writes “Therefore, to obtain good statistics for long time tails, one must collect data over hundreds of picoseconds (hundreds of nanoseconds at low T).” So, 500 ps sounds like a reasonable lower limit.
3

My .mdp file for production is:

; Run control
integrator = md
tinit = 0
dt = 0.002
nsteps = 250000 ; 500 ps
nstcomm = 100

; Output control
nstxout = 0
nstvout = 50
nstfout = 0
nstlog = 50
nstenergy = 50
nstxout-compressed = 50

I will try to run a 5 ns simulation and improve the trajectory output frequency.

4

You are writing coordinates and velocities every 100 fs. Luzar said that sampling intervals well below 10 fs should be used.

I’m usually suggesting to write coordinates (and velocities, if any) at low frequencies (every 10-200 ps depending on use case). But in this case, it seems like it’s important to write frequently.

Hopefully it will help.