|Subject:||Re: [ESPResSo-devel] Somewhat strange trajectory|
|Date:||Tue, 8 Apr 2014 17:36:34 +0200|
there is one question I have:
Do you observe a diffusive (a) or a directed (b) motion along 1,1,1? From the trajectory it looks like a directed motion but that would lead to a quadratic MSD. To quantify (a), it should be possible to derive a drift velocity, and for (b) you could calculated the off-diagonal elements of the diffusion tensor for <(x(t)-x(t+tau))(y(t)-y(t+tau))>. This would mean implementing the correlation operation "dyadic_product" which I want to do now for some time :-).
Regarding the reason: Variant a would be related to a correlation of subsequent numbers, and Variant b would be a nonzero mean. In principle, both should be easy to distinguish and in both cases, the same should apply for a single isolated particle. Maybe this would be an even easier test case. Use cellsystem nsquare.
Maybe that helps….
> Pokud je tento e-mail součástí obchodního jednání, Přírodovědecká fakulta Univerzity Karlovy v Praze:
On Apr 8, 2014, at 12:32 PM, Peter Košovan <address@hidden> wrote:
> Hi Espressies,
> Simulating a single polymer chain of the Kremer-Grest model, we observed that the centre of mass "diffuses" in narrow a slab, with a normal vector approx. (1,1,1). In order to see that, one needs to run about 10^8 LJ time units (~10^10 time steps, ~1cpu-day), plot the centre of mass trajectory and tilt it in an appropriate angle. With shorter trajectories it's not obvious if the bias is in the trajectory or in the observer's eye.
> I could partly fix the problem replacing the call to l_random() (ran1 RNG) in d_random() of random.hpp with the built-in rand() function. The trajectory might still be biased, but at least it's not a slab.
> Attached files: sample tcl script, gnuplot script to plot what I have described and its source data: centre-of-mass and monomer msds (average over 20 runs), com trajectory with the l_random() and with rand() RNG.
> I think that we should replace the l_random() by some other RNG. Surely, rand() is not a good choice either, it was just the first hack I tried. Any suggestions / opinions are welcome. I intend to investigate it further: test other RNGs and check if/how they affect dynamical and thermodynamic properties.
> Some more detail of what I learned:
> I asked a student to simulate the Kremer-Grest polymer in Langevin thermostat and compute msd(t), from 20 independent starting configurations, to get 20 independent msds's and estimate the error. Plotting msd(t)/t, it first converges to approx. 2/N, as it should for Rouse-like dynamics, but at t~3*10^4, it drops to ~3/4 of this value, while the error is much smaller. We suspected a mistaken input or a subtle correlator bug, so I repeated the simulations with my own script, various chain lengths, box sizes, gammas, timesteps, correlator parameters, initial conditions, and the result was robust. Finally, I printed the centre of mass coordinates and asked a colleague to compute the msd(t) using his own independent tool. Plotting the output trajectory revealed the problem. The tilt of the plane appears robust with respect to initial conditions, random seed and chain length. I have not yet checked other parameters, since other runs did not write the trajectory. I also checked the original Grest+Kremer paper (Phys.Rev.A,33,3628,1986), and their msds also change slope on long times (g3 in Fig.3), but they have not specified their RNG.
> Looking forward to your comments.
> Dr. Peter Košovan
> Departmtent of Physical and Macromolecular Chemistry
> Faculty of Science, Charles University in Prague, Czech Republic
> Katedra fyzikální a makromolekulární chemie
> Přírodovědecká fakulta Univerzity Karlovy v Praze
> Tel. +420221951290
> Fax +420224919752
> a) si vyhrazuje právo jednání kdykoliv ukončit a to i bez uvedení důvodu,
> b) stanovuje, že smlouva musí mít písemnou formu,
> c) vylučuje přijetí nabídky s dodatkem či odchylkou,
> d) stanovuje, že smlouva je uzavřena teprve výslovným dosažením shody na všech náležitostech smlouvy.
|[Prev in Thread]||Current Thread||[Next in Thread]|