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?
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? Thanks, On 16 August 2011 13:17, Justin A. Lemkul <jalem...@vt.edu> wrote: > > > Elisabeth wrote:f > >> Hello, >> >> 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? >> >> > 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. > > 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). > > The larger point is that if you claim to apply an NVT (or NPT) ensemble > using Berendsen coupling, you're not correct, strictly speaking. > > -Justin > > Appreciate any clarification. >> >> Thank you, >> Best, >> >> On 15 August 2011 19:34, Mark Abraham <mark.abra...@anu.edu.au <mailto: >> mark.abra...@anu.edu.**au <mark.abra...@anu.edu.au>>> wrote: >> >> On 16/08/2011 7:05 AM, Elisabeth wrote: >> >>> Dear all, >>> >>> I noticed that applying Parrinello-Rahman (PR) pressure coupling >>> even after equilibration with berendsen does not lead to target >>> value for pressure when >>> >>> ; Bonds >>> constraints = none >>> >>> is used. >>> >> >> The use of constraints and the integration step size is linked. >> Roughly speaking, no constraints should accompany a 0.5 fs time >> step, H-bond constraints with 1fs and all constraints with 2fs. >> Haphazard changes to .mdp files have all kinds of these "gotchas". >> >> >> I tried Berendsen for to get fixed pressure ( 50 bar) but in the >>> next run PR even for long time is giving 52 bar. This is the case >>> for other target pressures too. >>> >> >> You need to be sure to collect statistics only after equilibration, >> and consider whether the observed variation is consistent with >> convergence to a given value. >> >> >> >>> So this made me select Berenden which is giving target pressure >>> values but my concern is whether my results are reliable because >>> BR does not give the exact ensemble as PR. I read somewhere on the >>> list that fluctuation properties can not be calculated when BR is >>> used. What does "fluctuation property" mean? >>> >> >> BR does not produce the correct ensemble. I forget the details about >> why, but there are references in the T-coupling section of the >> manual you should consider. >> >> >> Does this mean that any property calculated form fluctuations of >>> some other quantity can not be obtained>? like heat capacity which >>> is defined based on enthalpy fluctuations? >>> >> >> IIRC, yes. >> >> >> I am interested in potential energy terms (g_energy bonded/non >>> boned terms) >>> >> >> I routinely struggle to see why people think they can learn anything >> from these. >> >> >> and structural properties like rdf for a number of polymer >>> molecules, system size around 3000 atoms. >>> >>> Thank you for your comments. >>> Best, >>> >>> constraints = none >>> ; Run control integrator = md >>> dt = 0.001 nsteps >>> = 5000000 nstcomm = 100 >>> >>> nstenergy = 100 nstxout = >>> 100 nstlist = 10 ns_type >>> = grid >>> coulombtype = Shift vdw-type >>> = Shift rcoulomb-switch = 0 >>> rvdw-switch = 0.9 ;0 >>> ; Cut-offs >>> rlist = 1.25 rcoulomb = >>> 1.0 rvdw = 1.0 >>> Tcoupl = v-rescale tc-grps >>> = System tau_t = 0.1 ref_t >>> = 300 >>> Pcoupl = berendsen Pcoupltype = >>> isotropic tau_p = 1 >>> compressibility = 3.5e-5 ref_p = 10 >>> >>> gen_vel = no >>> gen_temp = 300.0 gen_seed = >>> 173529 >>> >>> >>> >>> >>> >>> >>> >> >> -- >> gmx-users mailing list gmx-users@gromacs.org >> <mailto:gmx-users@gromacs.org> >> >> >> http://lists.gromacs.org/**mailman/listinfo/gmx-users<http://lists.gromacs.org/mailman/listinfo/gmx-users> >> Please search the archive at >> >> http://www.gromacs.org/**Support/Mailing_Lists/Search<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 >> <mailto:gmx-users-request@**gromacs.org<gmx-users-requ...@gromacs.org> >> >. >> >> Can't post? Read >> http://www.gromacs.org/**Support/Mailing_Lists<http://www.gromacs.org/Support/Mailing_Lists> >> >> >> > -- > ==============================**========== > > 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<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<http://lists.gromacs.org/mailman/listinfo/gmx-users> > Please search the archive at http://www.gromacs.org/** > Support/Mailing_Lists/Search<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<http://www.gromacs.org/Support/Mailing_Lists> >
-- 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
