[gmx-users] Re: problem in g_density

[cdalgicdir]

Hi Mark,

The set of frames is arbitrary, overall I do not expect any sudden 
significant change in density.

I calculated the density profiles by coinciding end of intervals, i.e 
instead of 5-10, 10-15, 15-20 and 5-20 ns, I calculated 4999-10000, 
9999-15000, 14999-end and 4999 - end ps. (simulation lasts for 20 ns)
where I expect the data sets do not skip any intermediate frames. 
However the results don't change.
I have also tried using traj.trr instead of traj.xtc, results don't 
change again.

What I don't understand is considering that the data sets represent the 
intervals,
x1: 5-10 ns
x2: 10-15 ns
x3: 15-20 ns

and X represents 5-20 ns.

X = x1+x2+x3
shouldn't the average of subset averages be equal to the average of the 
main set?
mean(X) = mean(mean(x1)+mean(x2)+mean(x3))


Thanks,
Cahit

On 06/08/2011 05:40 PM, Mark Abraham wrote:

Date: Wed, 08 Jun 2011 22:45:39 +1000
From: Mark Abraham<Mark.Abraham at anu.edu.au>
Subject: Re: [gmx-users] problem in g_density
To: Discussion list for GROMACS users<gmx-users at gromacs.org>
Message-ID:<4DEF6EF3.2060000 at anu.edu.au>
Content-Type: text/plain; charset=ISO-8859-9; format=flowed

On 8/06/2011 9:41 PM, cdalgicdir wrote:

> >  Hi,
> >
> >  The simulation box consists of layers of cyclohexane and vacuum.
> >  Ensemble is NVT.
> >
> >  The density profiles of cyclohexane for 5-10 ns, 10-15 ns and 15-20 ns
> >  coincide. However the density profile for 5-20 ns is significantly
> >  higher than the ones calculated.
> >  g_density -n -b 5000 -e 20000
> >
> >  If I calculate the 5-20 ns interval using every 5th frame the density
> >  profile resembles that of the 5-10, 10-15, 15-20 ns's.
> >  g_density -n -b 5000 -e 20000 -dt 5
> >
> >  Attached is the graph of these profiles.
> >
> >  I have tried this in 4.07, 4.5.3 and 4.5.4. The results are similar
> >  for all these versions.
> >
> >  I am aware that such inconsistencies may occur with g_density in NPT
> >  ensembles. However I am using a NVT ensemble. Is there something am I
> >  missing? And more importantly,  how should I decide which one is the
> >  correct density?
I'd guess this is an artefact of the fact that the times in the
trajectory file are not exactly multiples of the time step. That is just
life - see
http://www.gromacs.org/Documentation/Floating_Point_Arithmetic. If so,
one or more of these ranges is probably not computing over the set of
frames you think it is. Consider using suitable non-integers for -b and
-e, and avoiding -dt unless you know it will work on the times that
exist. trjcat -settime can help, also.

Mark