Andrew DeYoung wrote:
Hi,
I have run a simulation of 254 SPC/E water molecules in a cubic box of edge
length 19.660 Angstrom. I have used the gromacs utility g_velacc to
calculate the velocity autocorrelation function for the motion of only the
oxygen atoms. When I plot the velocity autocorrelation function versus
time, my result is a curve with more than one local minimum.
My question is, should I expect having more than one local minimum? I am
not sure that my result is reasonable.
For example, the Democritus molecular dynamics website gives one example of
a velocity autocorrelation function, although they do not specifically say
what system (water or something else?) or what state (liquid or gas?) was
used: http://www.cse.scitech.ac.uk/ccg/software/Democritus/Theory/vaf.html .
Their plot shows one "strong" local minimum and one "weaker" one; both of
these local minima are in the negative regime. I am new to this field, but
I believe that the "strong" local minimum corresponds conceptually to the
recoil of the atoms when they collide.
If you have time, could you please look at my result at
http://www.andrew.cmu.edu/user/adeyoung/velacc_may23.pdf and see if you
think it is a qualitatively reasonable result? My main concern is the first
local minimum, which I have labeled with an arrow on the plot. This minimum
is in the _positive_ regime, which I believe (but I am not sure) cannot
correspond to the recoil of atoms when they collide.
In my pdf file I have displayed the parameters I have used for equilibration
and dynamics. To summarize, I have used the leap frog algorithm with 2 ns
equilibration in 1 fs steps, and 5 ns dynamics in 1 fs steps. For dynamics,
I have saved configurations every 2 fs. Although my result from g_velacc
has velocity autocorrelation frunction data from 0 to 1000 ps, the function
decays rapidly to zero (as it should, I think, because the velocity rapidly
decorrelates since the atoms/molecules "forget" their initial velocities),
so I have plotted only from 0 to 1 ps. I am running Gromacs version 4.5.4.
I don't know if this explains your result or not, but you do realize that if you
have "gen_vel = yes" at the beginning of your simulation for data collection, it
negates any value in previous equilibration, right? You've destroyed the
ensemble you created before. You should pass the checkpoint file from the
equilibration run to grompp -t and set "gen_vel = no" for your simulation. I'd
start by doing that and seeing if it affects the result, because at this point,
it may just be that you haven't properly re-equilibrated, or at the very least,
the original frames are messing with the VACF calculation.
-Justin
--
========================================
Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
========================================
--
gmx-users mailing list gmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
