Thanks Justin for your thorough description. Is there any specific tool in gromacs to plot histograms or you mean I follow the customary plotting procedure?<br><br>I have one more question. I am now trying using PR pressure scheme with v-rescale thermostat for NPT dynamics. Before dynamics run, can I apply berendsen barostat combined with v-rescale thermostat to equilibrate or since I am using berendsen barostat, T coupling scheme has to be also berendsen?<br>
<br>Thanks,<br><br><br><div class="gmail_quote">On 16 August 2011 13:17, Justin A. Lemkul <span dir="ltr"><<a href="mailto:jalemkul@vt.edu">jalemkul@vt.edu</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin: 0pt 0pt 0pt 0.8ex; border-left: 1px solid rgb(204, 204, 204); padding-left: 1ex;">
<div class="im"><br>
<br>
Elisabeth wrote:f<br>
<blockquote class="gmail_quote" style="margin: 0pt 0pt 0pt 0.8ex; border-left: 1px solid rgb(204, 204, 204); padding-left: 1ex;">
Hello,<br>
<br>
Thank you. I am looking at potential energies to calculate vaporization heat. I wanted to know how the fact that berendsen does not lead to correct ensemble is affecting total potential energy of the system. I have an unclear image of "correct ensemble". Does this mean whatever output I am getting from my runs are unreliable?<br>
<br>
</blockquote>
<br></div>
The Berendsen coupling algorithms produce very narrow distributions of temperature and pressure. These do not correspond to the correct distribution of a true statistical mechanical ensemble; plot a histogram and you'll find that the results are shockingly different between Berendsen and, say, Nose-Hoover for T. Therefore, the ensemble you're using to measure your properties is not NPT (or NVT, in the case of thermostat only), it's something undefined and likely not real. Applying better thermostats and barostats that have been shown to produce the correct distributions is more rigorously correct.<br>
<br>
The effects are usually not terribly noticeable unless you go looking for them. You'll find that the Berendsen algorithms quite faithfully give you the target T and P that you desire (which does make them very good for initial equilibration) and so you naturally assume that everything is fine. While that's all well and good in one sense, the fluctuations of T (and thus velocities, and thus kinetic energy) are wrong. If KE does not fluctuate properly, then neither does PE (since you're then affecting energy flow between PE and KE, while total energy stays fixed, in theory).<br>
<br>
The larger point is that if you claim to apply an NVT (or NPT) ensemble using Berendsen coupling, you're not correct, strictly speaking.<br>
<br>
-Justin<br>
<br>
<blockquote class="gmail_quote" style="margin: 0pt 0pt 0pt 0.8ex; border-left: 1px solid rgb(204, 204, 204); padding-left: 1ex;"><div class="im">
Appreciate any clarification.<br>
<br>
Thank you,<br>
Best,<br>
<br></div><div><div></div><div class="h5">
On 15 August 2011 19:34, Mark Abraham <<a href="mailto:Mark.Abraham@anu.edu.au" target="_blank">Mark.Abraham@anu.edu.au</a> <mailto:<a href="mailto:Mark.Abraham@anu.edu.au" target="_blank">Mark.Abraham@anu.edu.<u></u>au</a>>> wrote:<br>
<br>
On 16/08/2011 7:05 AM, Elisabeth wrote:<br>
<blockquote class="gmail_quote" style="margin: 0pt 0pt 0pt 0.8ex; border-left: 1px solid rgb(204, 204, 204); padding-left: 1ex;">
Dear all,<br>
<br>
I noticed that applying Parrinello-Rahman (PR) pressure coupling<br>
even after equilibration with berendsen does not lead to target<br>
value for pressure when<br>
<br>
; Bonds<br>
constraints = none<br>
<br>
is used.<br>
</blockquote>
<br>
The use of constraints and the integration step size is linked.<br>
Roughly speaking, no constraints should accompany a 0.5 fs time<br>
step, H-bond constraints with 1fs and all constraints with 2fs.<br>
Haphazard changes to .mdp files have all kinds of these "gotchas".<br>
<br>
<br>
<blockquote class="gmail_quote" style="margin: 0pt 0pt 0pt 0.8ex; border-left: 1px solid rgb(204, 204, 204); padding-left: 1ex;">
I tried Berendsen for to get fixed pressure ( 50 bar) but in the<br>
next run PR even for long time is giving 52 bar. This is the case<br>
for other target pressures too.<br>
</blockquote>
<br>
You need to be sure to collect statistics only after equilibration,<br>
and consider whether the observed variation is consistent with<br>
convergence to a given value.<br>
<br>
<br>
<blockquote class="gmail_quote" style="margin: 0pt 0pt 0pt 0.8ex; border-left: 1px solid rgb(204, 204, 204); padding-left: 1ex;">
<br>
So this made me select Berenden which is giving target pressure<br>
values but my concern is whether my results are reliable because<br>
BR does not give the exact ensemble as PR. I read somewhere on the<br>
list that fluctuation properties can not be calculated when BR is<br>
used. What does "fluctuation property" mean?<br>
</blockquote>
<br>
BR does not produce the correct ensemble. I forget the details about<br>
why, but there are references in the T-coupling section of the<br>
manual you should consider.<br>
<br>
<br>
<blockquote class="gmail_quote" style="margin: 0pt 0pt 0pt 0.8ex; border-left: 1px solid rgb(204, 204, 204); padding-left: 1ex;">
Does this mean that any property calculated form fluctuations of<br>
some other quantity can not be obtained>? like heat capacity which<br>
is defined based on enthalpy fluctuations?<br>
</blockquote>
<br>
IIRC, yes.<br>
<br>
<br>
<blockquote class="gmail_quote" style="margin: 0pt 0pt 0pt 0.8ex; border-left: 1px solid rgb(204, 204, 204); padding-left: 1ex;">
I am interested in potential energy terms (g_energy bonded/non<br>
boned terms) <br>
</blockquote>
<br>
I routinely struggle to see why people think they can learn anything<br>
from these.<br>
<br>
<br>
<blockquote class="gmail_quote" style="margin: 0pt 0pt 0pt 0.8ex; border-left: 1px solid rgb(204, 204, 204); padding-left: 1ex;">
and structural properties like rdf for a number of polymer<br>
molecules, system size around 3000 atoms.<br>
<br>
Thank you for your comments.<br>
Best,<br>
<br>
constraints = none <br>
; Run control integrator = md dt = 0.001 nsteps = 5000000 nstcomm = 100 <br>
<br>
nstenergy = 100 nstxout = 100 nstlist = 10 ns_type = grid <br>
coulombtype = Shift vdw-type = Shift rcoulomb-switch = 0 rvdw-switch = 0.9 ;0 <br>
; Cut-offs<br>
rlist = 1.25 rcoulomb = 1.0 rvdw = 1.0 <br>
Tcoupl = v-rescale tc-grps = System tau_t = 0.1 ref_t = 300 <br>
Pcoupl = berendsen Pcoupltype = isotropic tau_p = 1 compressibility = 3.5e-5 ref_p = 10 <br>
gen_vel = no<br>
gen_temp = 300.0 gen_seed = 173529 <br>
<br>
<br>
<br>
<br>
<br>
<br>
</blockquote>
<br>
<br>
--<br>
gmx-users mailing list <a href="mailto:gmx-users@gromacs.org" target="_blank">gmx-users@gromacs.org</a><br></div></div>
<mailto:<a href="mailto:gmx-users@gromacs.org" target="_blank">gmx-users@gromacs.org</a>><div class="im"><br>
<a href="http://lists.gromacs.org/mailman/listinfo/gmx-users" target="_blank">http://lists.gromacs.org/<u></u>mailman/listinfo/gmx-users</a><br>
Please search the archive at<br>
<a href="http://www.gromacs.org/Support/Mailing_Lists/Search" target="_blank">http://www.gromacs.org/<u></u>Support/Mailing_Lists/Search</a> before posting!<br>
Please don't post (un)subscribe requests to the list. Use the<br>
www interface or send it to <a href="mailto:gmx-users-request@gromacs.org" target="_blank">gmx-users-request@gromacs.org</a><br></div>
<mailto:<a href="mailto:gmx-users-request@gromacs.org" target="_blank">gmx-users-request@<u></u>gromacs.org</a>>.<div class="im"><br>
Can't post? Read <a href="http://www.gromacs.org/Support/Mailing_Lists" target="_blank">http://www.gromacs.org/<u></u>Support/Mailing_Lists</a><br>
<br>
<br>
</div></blockquote>
<br>
-- <br>
==============================<u></u>==========<br>
<br>
Justin A. Lemkul<br>
Ph.D. Candidate<br>
ICTAS Doctoral Scholar<br>
MILES-IGERT Trainee<br>
Department of Biochemistry<br>
Virginia Tech<br>
Blacksburg, VA<br>
jalemkul[at]<a href="http://vt.edu" target="_blank">vt.edu</a> | <a href="tel:%28540%29%20231-9080" value="+15402319080" target="_blank">(540) 231-9080</a><br>
<a href="http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin" target="_blank">http://www.bevanlab.biochem.<u></u>vt.edu/Pages/Personal/justin</a><br>
<br>
==============================<u></u>==========<br><font color="#888888">
-- <br></font><div><div></div><div class="h5">
gmx-users mailing list <a href="mailto:gmx-users@gromacs.org" target="_blank">gmx-users@gromacs.org</a><br>
<a href="http://lists.gromacs.org/mailman/listinfo/gmx-users" target="_blank">http://lists.gromacs.org/<u></u>mailman/listinfo/gmx-users</a><br>
Please search the archive at <a href="http://www.gromacs.org/Support/Mailing_Lists/Search" target="_blank">http://www.gromacs.org/<u></u>Support/Mailing_Lists/Search</a> before posting!<br>
Please don't post (un)subscribe requests to the list. Use the www interface or send it to <a href="mailto:gmx-users-request@gromacs.org" target="_blank">gmx-users-request@gromacs.org</a>.<br>
Can't post? Read <a href="http://www.gromacs.org/Support/Mailing_Lists" target="_blank">http://www.gromacs.org/<u></u>Support/Mailing_Lists</a><br>
</div></div></blockquote></div><br>