[gmx-users] g_velacc again
Nuno R. L. Ferreira
nunolf at ci.uc.pt
Thu Aug 28 23:17:01 CEST 2003
Dear all gmx users
Few days ago I asked about some hints on how to get diffusion coefficients
from VACF (Cvv(t)).
Read some stuff about the subject (Allen & Tildesley, tkx Erik) and in other
As Anton said :
< the integral of your VACF is exactly the same as the slope of the
positional MSD (or something similar).
If I'm not wrong, "slope of the positional MSD" = Dt (self-diffusion
coefficient) = 1/3 * (integral of VACF) , speaking in 3D.
I've done a 100 ps MD, 2 fs of step and saving velocities every step (tkx
Anton), in a methanol box (400 molecules, tkx Christoph Freudenberger), and
tryed to compute the diffusion coefficient from MSD (g_msd) and also from
VACF (g_velacc). The Cvv(t) vs time graph shows an expected trend for this
kind of correlations (it decreases fast to some negative value, and then a
Cvv tail aproaches to zero)
g_msd gives a good value ( 2.7 +- .2 10-5 cm2 s-1), compared with the
experimental one ( 2.4 10-5 cm2 s-1).
But I "cannot manage" to obtain the value from g_velacc calculation. In
fact, I don't understand very well the output from g_velacc.
I've employed this script with the following flags:
g_velacc -f traj.trr -s traj.tpr -n
being the index file 10 methanol molecules to speed up things.
Here's the output:
COR: Correlation time (plain integral from 0.000 to 99.999 ps) = 0.01533 ps
What is this correlation time?
With -fitfn aexp , I obtain the same information as above plus:
COR: Relaxation times are computed as fit to an exponential:
COR: y = a2 exp(-x/a1)
COR: Fit to correlation function from 0.000 ps to 100.001 ps, results in a
COR: Fit from Integral Tail Value Sum (ps)
a1 (ps) a2 ()
COR: 0.0000e+00 0.0000e+00 6.0026e-02 6.0026e-02 5.1511e-02
I was expecting to get an integral value from here, but it's zero. I'm also
aware that the tail in Cvv(f) function as to be treated with special care
( Solvent diffusion outside macromolecular surfaces, Eric Lindahl & Olle
Edholm, Phys Rev. E ,V57, N1, Y1998, P791).
I used this function (aexp) cause I read somewhere (actually in a NAMD
tutorial, last summer course) that we can fit the Cvv(t) to a first order
exponential, like this :
Cvv(t) = Dt / tau * exp (-t / tau)
This example was applyed to ubiquitin in water, to calculate the Dt.
Can anibody give me some help on this? My aim is to obtain the diffusion
coefficient for a peptide. Since that for calculating MSD's of proteins with
need longer simulations, and the statisticall error with just one molecule
is "big" (MSD vs time does not show a good straight line), I was thinking in
VACF vs time. But till now I did not understood how to calculate Dt from the
information of g_velacc calculation. Where is the integral of VACF? Am I
P.S. I apologise if my explanation is a bit confuse.
Nuno Ricardo Santos Loureiro da Silva Ferreira
Grupo de Química Biológica
Departamento de Química
Faculdade de Ciências e Tecnologia
Universidade de Coimbra
Phone: +351 239 852080
Fax: +351 239 827703
More information about the gromacs.org_gmx-users