[gmx-users] Using g_dipoles to calculate dipole auto-correlation function

[James Zhu]

Hello!

I am a user from china. I am using gromacs-4.6.5. My problem is I can't get
the dipole auto-correlation function I want using g_dipoles.

I did a simulation of a system consisting of 6912 SPC water molecules. The
box size is 15.05872 * 3.76468 * 3.76468 nm^3. The simulation is all right.
Then I used g_dipoles to calculate the auto-correlation function averaged
over all molecules, which is defined as:
  C_{mol}(t)=\frac{1}{N}\sum_{i=1}^N<\mu_i(0) \cdot \mu_i(t)>
where \mu_i is the dipole of the i-th water molecule, < \cdot > means the
average over the simulation.

The command is:
g_dipoles -corr mol -f traj.allatom.TS7.528.newboxsize.trr -nonormalize
As the gromacs-4.6.5 manual said, "mol" means the correlation functions can
be averaged over all molecules.

The square mod of the dipole of a water molecule is 0.00229 calculate from
the configuration file below:
Pure Water
20736
    1SOL     OW    1  14.447   0.308   2.226  0.0391 -0.5260  0.6260
    1SOL    HW1    2  14.368   0.263   2.269  0.6952 -0.3402  2.0665
    1SOL    HW2    3  14.419   0.344   2.137 -2.2182  1.6884  2.0829
    2SOL     OW    4  14.245   1.162   0.750 -0.6480  0.1202  0.3786
    2SOL    HW1    5  14.326   1.103   0.744 -1.0285 -0.3643 -0.0118
    2SOL    HW2    6  14.218   1.190   0.658 -0.8392  0.4103  0.5220
    3SOL     OW    7   0.207   3.480   2.694 -0.4147 -0.2898  0.1138
    3SOL    HW1    8   0.180   3.432   2.778  0.8112  0.0259  0.7189
    3SOL    HW2    9   0.305   3.498   2.697  0.5814 -3.7135 -3.0765
......
Since the shape of each water molecule is almost the same, C(0) should be
very close to 0.00229, i.e. <\mu^2>. But the result of g_dipoles says C(0)
is 5.17, which is far away. Below is the output in the file dipcorr.xvg:
# This file was created Tue Jul 29 13:21:24 2014
# by the following command:
# g_dipoles -corr mol -f traj.allatom.TS7.528.newboxsize.trr -nonormalize
#
# g_dipoles is part of G R O M A C S:
#
# Gyas ROwers Mature At Cryogenic Speed
#
@    title "Dipole Autocorrelation Function"
@    xaxis  label "Time (ps)"
@    yaxis  label "C(t)"
@TYPE xy
     0.000     5.17053
     0.020     4.82886
     0.040     4.70246
     0.060     4.69387
     0.080     4.66471
......

So I wonder if I am not using the proper command or there is something
wrong with gromacs.

Thanks a lot for your help! If these information are not enough, please
contact me.

Best regards,
Jinglong Zhu



Jinglong Zhu
PhD Candidate in Computational Mathematics
School of Mathematical Sciences
Peking University
Beijing,100871
P.R.China