[gmx-users] No convergence in Diffusion Coefficient

[Igor Leontyev]

> On 2010-06-29 05.50, Igor Leontyev wrote:
>> See the table bellow, there is no convergence of the Self-Diffusion
>> Coefficient (Dself) over the trajectory length. Dself is obtained for NPT
>> box of 1024 SPC/E water molecules by the Einstein's relation (via RMSD).
>> In
>> Gromacs Manual or Alien&Tildesley's book I didn't find issues related to
>> the
>> problem. Is there any idea why I can not achieve the convergence?
>>
>> Description:
>> To figure out what the trajectory length is needed for an accurate
>> simulation of the self-diffusion coefficient I performed the following
>> test.
>> Split the continuous 10ns long trajectory on parts and calculated Dself
>> for
>> each of the parts, e.g. the splitting on N parts gives N values of Dself
>> obtained on trjlenth=10ns/N trajectory part. For the variety of N values
>> we
>> can calculate the average <Dself> and dispersion Disper. The converged
>> trajectory length is found as the trjlenth value at which the dispersion
>> is
>> sufficiently small and <Dself> equal to the value obtained for the whole
>> 10ns trajectory. But it turned out that Dself does not converge (See the
>> columns 3 and 4 in the Table).
>>
>> Just for the comparison I carried out the same test for the Dielectric
>> Constant Eps (See the columns 5 and 6 in the Table) and the converged
>> trajectory length is about 2.5-5ns which correlates with the length
>> reported
>> in the literature.


> On 2010-06-29 David van der Spoel wrote:
> You don't say how you compute the Dself.

Some details are given on the top of my initial message. The command line
is: "g_msd -trestart 10 -dt 0.5"

The obtained rmds dependence on t is a perfect straight line.

> On 2010-06-29 David van der Spoel wrote:
> Your Dself varies from 2.48 to 2.63,
>  and if you drop the 10 ns and 4.8 ps measurements it varies from
> 2.48 to 2.54. Not too much I would say.

For the fine tuning of water model parameters such a precision is not
enough. If the scattering is between 2.4 and 2.6 it is not clear in which
direction one need to modify Charges or vdW parameters to improve Dself.
Moreover, the approaching of  <Dself> to its infinitely averaged value is
irregular, i.e. as you mentioned the value obtained with the best statistics
(10 ns trajectory) should be dropped while averaging over shorter
trajectories seems to produce better results. I observe no convergence even
for the 100ns trajectory, i.e. improving the statistics to infinity does not
decrease the uncertainty of Dself. By definition it just means there is no
convergence.

The convergence is observed for Dielectric Constant. It regularly approaches
to the converged value 72.

> On 2010-06-29 David van der Spoel wrote: g_msd tries to do something
> semi-intelligent, by fitting the MSD to a
> straight line from 10 to 90% of the length of the graph. You should
> avoid using the very first few ps and the final bit of the graph.

It should not significantly distort results if the trajectory length is in
the range of nanosec. Isn't it?

> Another way of getting statistics is by spitting the system in 2-N
> sub-boxes and computing the Dself for each of these.
>
>>
>> ------------------------------------------------------------
>> trjlenth N <Dself> Disper <Eps> Disper
>> ps 1e-5 cm^2/s
>> ------------------------------------------------------------
>> 10000.0 1 2.6305 0 72.0 0
>> 5000.0 2 2.4591 .0929 71.9 1.1
>> 2500.0 4 2.4554 .0291 71.6 4.7
>> 1250.0 8 2.4816 .0618 71.1 5.2
>> 625.0 16 2.4848 .0786 70.6 7.5
>> 312.5 32 2.4993 .0848 67.7 9.4
>> 156.2 64 2.5128 .1082 63.3 11.6
>> 78.1 128 2.4993 .1119 56.3 14.1
>> 39.0 256 2.5072 .1303 43.6 14.9
>> 19.5 512 2.5133 .1713 31.0 12.4
>> 9.7 1030 2.5449 .2128 19.6 8.0
>> 4.8 2083 2.6318 .2373 12.0 4.9