[gmx-users] Re: v-rescale
Peter C. Lai wrote Generally it's probably not a good idea to rely on tutorials designed around 3.3 when a google search for gromacs tutorial shows a series of 4.5.x tutorials written by Justin himself, with explanations of why certain steps are conducted. (also when certain features may be implemented differently, such as the introduction of the v-rescale thermostat). Well, it's not like I went looking for obsolete information... I started with GROMACS 3.3, and I wasn't keeping abreast of the changes as the software was upgraded. -- View this message in context: http://gromacs.5086.n6.nabble.com/v-rescale-tp5001066p5001187.html Sent from the GROMACS Users Forum mailing list archive at Nabble.com. -- gmx-users mailing listgmx-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
[gmx-users] Re: v-rescale
Hi Peter, Thanks for your response. Rather than dragging this thread too far off-topic, I'll direct you back to my thread, where I have just posted additional details. I took a warning message from GROMACS a bit too literally and it caused me to use conditions that blew up my simulations. I am interested in your protocol for the typical equilibration. If this is in fact standardized, do you have a reference? It doesn't match up with anything in the tutorial files I have been using to run my own simulations. Admittedly, those files are from GROMACS 3.3, and the procedures may be a bit out of date. -- View this message in context: http://gromacs.5086.n6.nabble.com/v-rescale-tp5001066p5001136.html Sent from the GROMACS Users Forum mailing list archive at Nabble.com. -- gmx-users mailing listgmx-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
Re: [gmx-users] Re: v-rescale
On 2012-09-20 12:18:02AM -0700, Ladasky wrote: Hi Peter, Thanks for your response. Rather than dragging this thread too far off-topic, I'll direct you back to my thread, where I have just posted additional details. I took a warning message from GROMACS a bit too literally and it caused me to use conditions that blew up my simulations. I am interested in your protocol for the typical equilibration. If this is in fact standardized, do you have a reference? It doesn't match up with anything in the tutorial files I have been using to run my own simulations. Admittedly, those files are from GROMACS 3.3, and the procedures may be a bit out of date. Generally it's probably not a good idea to rely on tutorials designed around 3.3 when a google search for gromacs tutorial shows a series of 4.5.x tutorials written by Justin himself, with explanations of why certain steps are conducted. (also when certain features may be implemented differently, such as the introduction of the v-rescale thermostat). -- == Peter C. Lai| University of Alabama-Birmingham Programmer/Analyst | KAUL 752A Genetics, Div. of Research | 705 South 20th Street p...@uab.edu| Birmingham AL 35294-4461 (205) 690-0808 | == -- gmx-users mailing listgmx-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
Re: [gmx-users] Re: v-rescale
On 20/09/12 01:35, Peter C. Lai wrote: then switching to nose-hoover for production runs (as nose-hoover chains result in the correct canonical distribution)? I was under the impression that v-rescale resulted in the correct canonical distribution as well. Is this incorrect? -- Massimo Sandal, Ph.D. http://devicerandom.org -- gmx-users mailing listgmx-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
Re: [gmx-users] Re: v-rescale
On 20/09/2012 9:35 AM, Peter C. Lai wrote: I am not sure where the idea of using berendsen barostat with the v-rescale thermostat for equilibration came from, however. Doesn't the typical equilibration begin with v-rescale for temperature equilibration then adding parinello-rahman barostat then switching to nose-hoover for production runs (as nose-hoover chains result in the correct canonical distribution)? N-H does have known issues, see http://link.aip.org/link/doi/10.1063/1.2989802 and http://link.aip.org/link/doi/10.1063/1.2408420. I am not aware of any shortcomings of the Bussi v-rescale thermostat in GROMACS. Mark On 2012-09-19 04:24:27PM -0700, Ladasky wrote: Dear Sara, I just had a problem with my simulations that I traced to the use of the V-rescale temperature algorithm. Here is my recent post: http://gromacs.5086.n6.nabble.com/Re-Water-molecules-cannot-be-settled-why-td4999302.html;cid=1348087067061-71#a5001121 V-rescale may be appropriate in certain simulations, but it is apparently NOT appropriate when used with Berendsen pressure coupling during the initial equilibration. I don't know if that is related to your problem, but it's something that I just discovered the hard way. -- View this message in context: http://gromacs.5086.n6.nabble.com/v-rescale-tp5001066p5001122.html Sent from the GROMACS Users Forum mailing list archive at Nabble.com. -- gmx-users mailing listgmx-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 -- gmx-users mailing listgmx-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
Re: [gmx-users] Re: v-rescale
I've done some extensive testing (paper on testing method in the works) and vrescale gives a very accurate ensemble very well for NVT. Parrinello-Rahman and MTTK are the only algorithms that are correct for NPT. Berendsen barostat is not. Note that there is a bug with vrescale + md-vv + that is fixed in 4.6 and (hopefully) for the next patch of 4.5. For equilibrating a system, Berendsen for both temperature and pressure is the best bet. They artificially minimize fluctuations, which is great for equilibration, bad for data collection. On Thu, Sep 20, 2012 at 6:05 PM, Mark Abraham mark.abra...@anu.edu.au wrote: On 20/09/2012 9:35 AM, Peter C. Lai wrote: I am not sure where the idea of using berendsen barostat with the v-rescale thermostat for equilibration came from, however. Doesn't the typical equilibration begin with v-rescale for temperature equilibration then adding parinello-rahman barostat then switching to nose-hoover for production runs (as nose-hoover chains result in the correct canonical distribution)? N-H does have known issues, see http://link.aip.org/link/doi/10.1063/1.2989802 and http://link.aip.org/link/doi/10.1063/1.2408420. I am not aware of any shortcomings of the Bussi v-rescale thermostat in GROMACS. Mark On 2012-09-19 04:24:27PM -0700, Ladasky wrote: Dear Sara, I just had a problem with my simulations that I traced to the use of the V-rescale temperature algorithm. Here is my recent post: http://gromacs.5086.n6.nabble.com/Re-Water-molecules-cannot-be-settled-why-td4999302.html;cid=1348087067061-71#a5001121 V-rescale may be appropriate in certain simulations, but it is apparently NOT appropriate when used with Berendsen pressure coupling during the initial equilibration. I don't know if that is related to your problem, but it's something that I just discovered the hard way. -- View this message in context: http://gromacs.5086.n6.nabble.com/v-rescale-tp5001066p5001122.html Sent from the GROMACS Users Forum mailing list archive at Nabble.com. -- gmx-users mailing listgmx-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 -- gmx-users mailing listgmx-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 -- gmx-users mailing listgmx-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
[gmx-users] Re: v-rescale
Dear Sara, I just had a problem with my simulations that I traced to the use of the V-rescale temperature algorithm. Here is my recent post: http://gromacs.5086.n6.nabble.com/Re-Water-molecules-cannot-be-settled-why-td4999302.html;cid=1348087067061-71#a5001121 V-rescale may be appropriate in certain simulations, but it is apparently NOT appropriate when used with Berendsen pressure coupling during the initial equilibration. I don't know if that is related to your problem, but it's something that I just discovered the hard way. -- View this message in context: http://gromacs.5086.n6.nabble.com/v-rescale-tp5001066p5001122.html Sent from the GROMACS Users Forum mailing list archive at Nabble.com. -- gmx-users mailing listgmx-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
Re: [gmx-users] Re: v-rescale
I am not sure where the idea of using berendsen barostat with the v-rescale thermostat for equilibration came from, however. Doesn't the typical equilibration begin with v-rescale for temperature equilibration then adding parinello-rahman barostat then switching to nose-hoover for production runs (as nose-hoover chains result in the correct canonical distribution)? On 2012-09-19 04:24:27PM -0700, Ladasky wrote: Dear Sara, I just had a problem with my simulations that I traced to the use of the V-rescale temperature algorithm. Here is my recent post: http://gromacs.5086.n6.nabble.com/Re-Water-molecules-cannot-be-settled-why-td4999302.html;cid=1348087067061-71#a5001121 V-rescale may be appropriate in certain simulations, but it is apparently NOT appropriate when used with Berendsen pressure coupling during the initial equilibration. I don't know if that is related to your problem, but it's something that I just discovered the hard way. -- View this message in context: http://gromacs.5086.n6.nabble.com/v-rescale-tp5001066p5001122.html Sent from the GROMACS Users Forum mailing list archive at Nabble.com. -- gmx-users mailing listgmx-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 -- == Peter C. Lai| University of Alabama-Birmingham Programmer/Analyst | KAUL 752A Genetics, Div. of Research | 705 South 20th Street p...@uab.edu| Birmingham AL 35294-4461 (205) 690-0808 | == -- gmx-users mailing listgmx-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
[gmx-users] Re: v-rescale - harmonic oscillator
Dear Servaas,
the problem that you found is not related to the stochastic rescaling
itself, but to the way the effective energy is calculated in gromacs.
In particular, the increment in the integral part is not consistent
with the applied scaling. You can correct for this by changing
src/mdlib/coupling.c, line 462:
// NEW:
if (Ek + dEk = 0) {
ekind-tcstat[i].lambda = 0.0;
therm_integral[i] -= -Ek;
} else {
ekind-tcstat[i].lambda = sqrt((Ek + dEk)/Ek);
ekind-tcstat[i].lambda = max(min(1.25,ekind-tcstat[i].lambda),0.8);
therm_integral[i] -=
Ek*(ekind-tcstat[i].lambda*ekind-tcstat[i].lambda-1.0);
}
// OLD:
// if (Ek + dEk = 0) {
// ekind-tcstat[i].lambda = 0.0;
// } else {
// ekind-tcstat[i].lambda = sqrt((Ek + dEk)/Ek);
// ekind-tcstat[i].lambda = max(min(1.25,ekind-tcstat[i].lambda),0.8);
// }
// therm_integral[i] -= dEk;
Note that this is a problem only when lambda is very different from
one, and this happens in very small systems. The effect for large
systems should be negligible. Furthermore, you should combine this
with the fix that Berk posted a few days ago, which is in file
src/mdlid/update.c, line 165:
// OLD:
// vv = lg*(vn + f[n][d]*w_dt);
// NEW:
vv = lg*vn + f[n][d]*w_dt;
These fixes will give you a much better energy conservation.
Consider also the following comments:
* For very small systems, the special cases for lambda (max(min...))
introduce small artifacts in the ensemble which are out of control.
Thus, I don't think that you can really use this scheme to perform
rigorous hybrid Monte Carlo. A possible solution is to use a better
integrator for the thermostat (the one discussed in the appendix of
Bussi, Donadio and Parrinello, JCP 2007), which works for every choice
of taut and any number of degrees of freedom, thus it does not need
any special case for lambda. I can send you an implementation of it
for gromacs. Berk is also working on that, and will probabily
introduce it soon in the official gromacs (4.1 or the cvs).
* Even with all these fixes (including the better integrator), there
still a (small) drift in the effective energy for the harmonic
oscillator. I did not test, but I guess that it originates from the
fact that gromacs uses leapfrog and not velocity-Verlet, and leapfrog
is not easily combined with our stochastic scheme which is ultimately
based on Trotter decomposition. The drift should *completely*
disappear for a harmonic system integrated with velocity-Verlet. If
you are interested, we have proven this for a different thermostat
(Bussi Parrinello, PRE 2007, for Langevin dynamics). The demonstration
can be straightforwardly ported to the velocity-rescaling thermostat.
* A part from the two previous implementation issues, I do not expect
any problem with the stochastic rescaling applied to the harmonic
oscillator. In particular, lack of ergodicity in harmonic systems
(which is the big pain in Nose-Hoover) should not be a problem at all.
Furthermore, the algorithm should work and provide the proper ensemble
also on very small systems.
Best regards,
Giovanni Bussi
___
gmx-users mailing listgmx-users@gromacs.org
http://www.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] Re: v-rescale - harmonic oscillator
Dear Dr. Bussi,
Thank you very much for the explanation of this problem.
So I understand that a problem remains with the integrator used in
GROMACS for small systems and that Berk is fixing this (and also a small
error is comming from using leap frog instead of velocity verlet).
I don't use this thermostat in my hybrid monte carlo algoritm I just
compared the results from a hybrid monte carlo run the results with
v-rescale. So for me there is no hurry the implement the different
integrator, since Berk is already working on it I think it is not very
useful if I also start doing this. I am prepared to help with it or test
it if this is required.
Kind regards,
Servaas
-
Message: 1
Date: Tue, 12 May 2009 09:06:52 +0200
From: Giovanni Bussi giovanni.bu...@gmail.com
Subject: [gmx-users] Re: v-rescale - harmonic oscillator
To: gmx-users@gromacs.org
Message-ID:
3297c3ca0905120006l2f8e630dv34c2ed4b6cd02...@mail.gmail.com
Content-Type: text/plain; charset=ISO-8859-1
Dear Servaas,
the problem that you found is not related to the stochastic rescaling
itself, but to the way the effective energy is calculated in gromacs.
In particular, the increment in the integral part is not consistent
with the applied scaling. You can correct for this by changing
src/mdlib/coupling.c, line 462:
// NEW:
if (Ek + dEk = 0) {
ekind-tcstat[i].lambda = 0.0;
therm_integral[i] -= -Ek;
} else {
ekind-tcstat[i].lambda = sqrt((Ek + dEk)/Ek);
ekind-tcstat[i].lambda = max(min(1.25,ekind-tcstat[i].lambda),0.8);
therm_integral[i] -=
Ek*(ekind-tcstat[i].lambda*ekind-tcstat[i].lambda-1.0);
}
// OLD:
// if (Ek + dEk = 0) {
// ekind-tcstat[i].lambda = 0.0;
// } else {
// ekind-tcstat[i].lambda = sqrt((Ek + dEk)/Ek);
// ekind-tcstat[i].lambda = max(min(1.25,ekind-tcstat[i].lambda),0.8);
// }
// therm_integral[i] -= dEk;
Note that this is a problem only when lambda is very different from
one, and this happens in very small systems. The effect for large
systems should be negligible. Furthermore, you should combine this
with the fix that Berk posted a few days ago, which is in file
src/mdlid/update.c, line 165:
// OLD:
// vv = lg*(vn + f[n][d]*w_dt);
// NEW:
vv = lg*vn + f[n][d]*w_dt;
These fixes will give you a much better energy conservation.
Consider also the following comments:
* For very small systems, the special cases for lambda (max(min...))
introduce small artifacts in the ensemble which are out of control.
Thus, I don't think that you can really use this scheme to perform
rigorous hybrid Monte Carlo. A possible solution is to use a better
integrator for the thermostat (the one discussed in the appendix of
Bussi, Donadio and Parrinello, JCP 2007), which works for every choice
of taut and any number of degrees of freedom, thus it does not need
any special case for lambda. I can send you an implementation of it
for gromacs. Berk is also working on that, and will probabily
introduce it soon in the official gromacs (4.1 or the cvs).
* Even with all these fixes (including the better integrator), there
still a (small) drift in the effective energy for the harmonic
oscillator. I did not test, but I guess that it originates from the
fact that gromacs uses leapfrog and not velocity-Verlet, and leapfrog
is not easily combined with our stochastic scheme which is ultimately
based on Trotter decomposition. The drift should *completely*
disappear for a harmonic system integrated with velocity-Verlet. If
you are interested, we have proven this for a different thermostat
(Bussi Parrinello, PRE 2007, for Langevin dynamics). The demonstration
can be straightforwardly ported to the velocity-rescaling thermostat.
* A part from the two previous implementation issues, I do not expect
any problem with the stochastic rescaling applied to the harmonic
oscillator. In particular, lack of ergodicity in harmonic systems
(which is the big pain in Nose-Hoover) should not be a problem at all.
Furthermore, the algorithm should work and provide the proper ensemble
also on very small systems.
Best regards,
Giovanni Bussi
--
___
gmx-users mailing listgmx-users@gromacs.org
http://www.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
RE: [gmx-users] Re: v-rescale - harmonic oscillator
Hi,
I have an optimal leap-frog v-rescale thermostat implemented for CVS head.
This will be in Gromacs 4.1.
Gromacs 4.1 will probably also have a velocity verlet integrator,
for which we can also remove the last very small drift in the v-rescale
implementation.
Berk
From: servaas.michielss...@student.kuleuven.be
To: gmx-users@gromacs.org
Date: Tue, 12 May 2009 13:25:02 +0200
Subject: [gmx-users] Re: v-rescale - harmonic oscillator
Dear Dr. Bussi,
Thank you very much for the explanation of this problem.
So I understand that a problem remains with the integrator used in
GROMACS for small systems and that Berk is fixing this (and also a small
error is comming from using leap frog instead of velocity verlet).
I don't use this thermostat in my hybrid monte carlo algoritm I just
compared the results from a hybrid monte carlo run the results with
v-rescale. So for me there is no hurry the implement the different
integrator, since Berk is already working on it I think it is not very
useful if I also start doing this. I am prepared to help with it or test
it if this is required.
Kind regards,
Servaas
-
Message: 1
Date: Tue, 12 May 2009 09:06:52 +0200
From: Giovanni Bussi giovanni.bu...@gmail.com
Subject: [gmx-users] Re: v-rescale - harmonic oscillator
To: gmx-users@gromacs.org
Message-ID:
3297c3ca0905120006l2f8e630dv34c2ed4b6cd02...@mail.gmail.com
Content-Type: text/plain; charset=ISO-8859-1
Dear Servaas,
the problem that you found is not related to the stochastic rescaling
itself, but to the way the effective energy is calculated in gromacs.
In particular, the increment in the integral part is not consistent
with the applied scaling. You can correct for this by changing
src/mdlib/coupling.c, line 462:
// NEW:
if (Ek + dEk = 0) {
ekind-tcstat[i].lambda = 0.0;
therm_integral[i] -= -Ek;
} else {
ekind-tcstat[i].lambda = sqrt((Ek + dEk)/Ek);
ekind-tcstat[i].lambda =
max(min(1.25,ekind-tcstat[i].lambda),0.8);
therm_integral[i] -=
Ek*(ekind-tcstat[i].lambda*ekind-tcstat[i].lambda-1.0);
}
// OLD:
// if (Ek + dEk = 0) {
// ekind-tcstat[i].lambda = 0.0;
// } else {
// ekind-tcstat[i].lambda = sqrt((Ek + dEk)/Ek);
// ekind-tcstat[i].lambda =
max(min(1.25,ekind-tcstat[i].lambda),0.8);
// }
// therm_integral[i] -= dEk;
Note that this is a problem only when lambda is very different from
one, and this happens in very small systems. The effect for large
systems should be negligible. Furthermore, you should combine this
with the fix that Berk posted a few days ago, which is in file
src/mdlid/update.c, line 165:
// OLD:
// vv = lg*(vn + f[n][d]*w_dt);
// NEW:
vv = lg*vn + f[n][d]*w_dt;
These fixes will give you a much better energy conservation.
Consider also the following comments:
* For very small systems, the special cases for lambda (max(min...))
introduce small artifacts in the ensemble which are out of control.
Thus, I don't think that you can really use this scheme to perform
rigorous hybrid Monte Carlo. A possible solution is to use a better
integrator for the thermostat (the one discussed in the appendix of
Bussi, Donadio and Parrinello, JCP 2007), which works for every choice
of taut and any number of degrees of freedom, thus it does not need
any special case for lambda. I can send you an implementation of it
for gromacs. Berk is also working on that, and will probabily
introduce it soon in the official gromacs (4.1 or the cvs).
* Even with all these fixes (including the better integrator), there
still a (small) drift in the effective energy for the harmonic
oscillator. I did not test, but I guess that it originates from the
fact that gromacs uses leapfrog and not velocity-Verlet, and leapfrog
is not easily combined with our stochastic scheme which is ultimately
based on Trotter decomposition. The drift should *completely*
disappear for a harmonic system integrated with velocity-Verlet. If
you are interested, we have proven this for a different thermostat
(Bussi Parrinello, PRE 2007, for Langevin dynamics). The demonstration
can be straightforwardly ported to the velocity-rescaling thermostat.
* A part from the two previous implementation issues, I do not expect
any problem with the stochastic rescaling applied to the harmonic
oscillator. In particular, lack of ergodicity in harmonic systems
(which is the big pain in Nose-Hoover) should not be a problem at all.
Furthermore, the algorithm should work and provide the proper ensemble
also on very small systems.
Best regards,
Giovanni Bussi
--
___
gmx-users mailing listgmx-users@gromacs.org
http
RE: [gmx-users] Re: v-rescale - harmonic oscillator
Hi, I have come to the conclusion that the v-rescale thermostat is only correct in the limit of large numbers of degrees of freedom, which is equivalent to a small scaling limit. The conservation of the conserved energy quantity improves as the system gets bigger. But I noticed that in 4.0 the energy conservation is not good. It can be improved significantly by changing line 164 or src/mdlib/update.c from vv = lg*(vn + f[n][d]*w_dt); to vv = lg*vn + f[n][d]*w_dt; There are some other similar lines in case you are PR p-coupling or cosine acceleration. I will probably not fix this for 4.0.5, but only for 4.1, so we do not change the reproducibility in minor revisions. Berk From: servaas.michielss...@student.kuleuven.be To: gmx-users@gromacs.org Date: Thu, 7 May 2009 15:27:37 +0200 Subject: [gmx-users] Re: v-rescale - harmonic oscillator Hi Berk, Thank you for your concern on this issue. It did not expect it to be a problem from the paper of Bussi et al., but it is of course possible. I am very curious to hear what they have to say about this. Are you contacting him or will I do it? For the rotational degrees of freedom would removing them with comm_mode=angular be a problem then because you are using restraints? (you use orire than?) Thanks again, Servaas Hi, With a diatomic molecule there is the issue that the rotation of the molecule has two degrees of freedom which are completely uncoupled from the rest of the system. Anyhow, I made a 1D harmonic oscillator using a position restraint and even there the energy does not seem to be conserved. I am not sure if a harmonic oscillator is a pathetic case, like with Nose-Hoover, or if this is an issue with the integrator implementation in Gromacs. I might need to ask Bussi about this. Berk From: servaas.michielss...@student.kuleuven.be To: gmx-users@gromacs.org Date: Thu, 7 May 2009 14:03:30 +0200 Subject: [gmx-users] v-rescale - harmonic oscillator Hi, It is just a diatomic molecule without other interactions (so only bonded interaction, LJ=0 and charges=0 on the atoms). What other details would you like to know? Should I send you the tpr file? Servaas Hi, No, only for small values is should be off, maybe tau_t 10 or 100 * delta_t. How did you set up your harmonic oscillator? Berk From: servaas.michielss...@student.kuleuven.be To: gmx-users@gromacs.org Date: Thu, 7 May 2009 12:54:54 +0200 Subject: [gmx-users] Re: v-rescale - harmonic oscillator Hi, Thank you for your reply. Yes by effective energy I mean the Gromacs conserved energy term. So you would say the cause of my problem here is that I need a small tau_t to thermostat this system but for small tau_t values the algorithm is not 100% ok yet? And what are very small tau_t values? I tested with quite a big range here. The problem occurred to me when I was experimenting with a hybrid monte carlo algorithm. I first tested it on a noble gas and compared to MD simulation with nose-hoover and v-rescale. Those results were 100% identical. Then I tried it on a dipeptide in vacuum, results were a little off here (I would expect this for nose-hoover thermostat, but not for v-rescale). So I decided to go the the extreme simple case of a harmonic oscillator. Those results were seriously off, and then I checked the conserved energy term, which seems impossible to keep it conserved for this system. (of course there can also be an error in my hybrid monte carlo code, but the fact that the conserved energy is not conserved is disturbing me here) Servaas Hi, What do you mean with effective energy? The Gromacs conserved energy term? For very small tau_t the current implementation does not work well. Bussi mailed me a proper implementation that I will put in when I have time. Berk From: servaas.michielssens at student.kuleuven.be To: gmx-users at gromacs.org Date: Thu, 7 May 2009 09:02:52 +0200 Subject: [gmx-users] v-rescale - harmonic oscillator Hi, I did some experiments with a harmonic oscillator (diatomic molecule without charge en LJ parameters) using the v-rescale thermostat. First I ran a simulation in the NVE ensemble with a time-step of 0.0001ps, total energy was constant here. Then I tried simulating the system with the v-rescale thermostat, same time-step of 0.0001ps and tau_t=100,10,1,0.1,0.01,0.001. For none of those simulations the effective energy was conserved. (effective energy should be a conserved quantity for NVT simulations
[gmx-users] Re: v-rescale - harmonic oscillator
Hi, Thanks for your confirmation of this problem. One thing is not clear to me from your answer. Has v-rescale fundamental problems with a Harmonic oscillator or is it a problem with the GROMACS implementation? Kind regards, Servaas Message: 3 Date: Fri, 8 May 2009 10:12:45 +0200 From: Berk Hess g...@hotmail.com Subject: RE: [gmx-users] Re: v-rescale - harmonic oscillator To: Discussion list for GROMACS users gmx-users@gromacs.org Message-ID: col113-w10f68330ea113929cbc3518e...@phx.gbl Content-Type: text/plain; charset=iso-8859-1 Hi, I have come to the conclusion that the v-rescale thermostat is only correct in the limit of large numbers of degrees of freedom, which is equivalent to a small scaling limit. The conservation of the conserved energy quantity improves as the system gets bigger. But I noticed that in 4.0 the energy conservation is not good. It can be improved significantly by changing line 164 or src/mdlib/update.c from vv = lg*(vn + f[n][d]*w_dt); to vv = lg*vn + f[n][d]*w_dt; There are some other similar lines in case you are PR p-coupling or cosine acceleration. I will probably not fix this for 4.0.5, but only for 4.1, so we do not change the reproducibility in minor revisions. Berk From: servaas.michielss...@student.kuleuven.be To: gmx-users@gromacs.org Date: Thu, 7 May 2009 15:27:37 +0200 Subject: [gmx-users] Re: v-rescale - harmonic oscillator Hi Berk, Thank you for your concern on this issue. It did not expect it to be a problem from the paper of Bussi et al., but it is of course possible. I am very curious to hear what they have to say about this. Are you contacting him or will I do it? For the rotational degrees of freedom would removing them with comm_mode=angular be a problem then because you are using restraints? (you use orire than?) Thanks again, Servaas Hi, With a diatomic molecule there is the issue that the rotation of the molecule has two degrees of freedom which are completely uncoupled from the rest of the system. Anyhow, I made a 1D harmonic oscillator using a position restraint and even there the energy does not seem to be conserved. I am not sure if a harmonic oscillator is a pathetic case, like with Nose-Hoover, or if this is an issue with the integrator implementation in Gromacs. I might need to ask Bussi about this. Berk From: servaas.michielss...@student.kuleuven.be To: gmx-users@gromacs.org Date: Thu, 7 May 2009 14:03:30 +0200 Subject: [gmx-users] v-rescale - harmonic oscillator Hi, It is just a diatomic molecule without other interactions (so only bonded interaction, LJ=0 and charges=0 on the atoms). What other details would you like to know? Should I send you the tpr file? Servaas Hi, No, only for small values is should be off, maybe tau_t 10 or 100 * delta_t. How did you set up your harmonic oscillator? Berk From: servaas.michielss...@student.kuleuven.be To: gmx-users@gromacs.org Date: Thu, 7 May 2009 12:54:54 +0200 Subject: [gmx-users] Re: v-rescale - harmonic oscillator Hi, Thank you for your reply. Yes by effective energy I mean the Gromacs conserved energy term. So you would say the cause of my problem here is that I need a small tau_t to thermostat this system but for small tau_t values the algorithm is not 100% ok yet? And what are very small tau_t values? I tested with quite a big range here. The problem occurred to me when I was experimenting with a hybrid monte carlo algorithm. I first tested it on a noble gas and compared to MD simulation with nose-hoover and v-rescale. Those results were 100% identical. Then I tried it on a dipeptide in vacuum, results were a little off here (I would expect this for nose-hoover thermostat, but not for v-rescale). So I decided to go the the extreme simple case of a harmonic oscillator. Those results were seriously off, and then I checked the conserved energy term, which seems impossible to keep it conserved for this system. (of course there can also be an error in my hybrid monte carlo code, but the fact that the conserved energy is not conserved is disturbing me here) Servaas Hi, What do you mean with effective energy? The Gromacs conserved energy term? For very small tau_t the current implementation does not work well. Bussi mailed me a proper implementation that I will put in when I have time. Berk From: servaas.michielssens at student.kuleuven.be To: gmx-users at gromacs.org
RE: [gmx-users] Re: v-rescale - harmonic oscillator
Hi, That is a fundamental problem of most global thermostats, just like Nose-Hoover. Nose-Hoover chains solve this. The only thermostats that guarantee proper ensembles for any system are local ones such as Langevin dynamics. Berk From: servaas.michielss...@student.kuleuven.be To: gmx-users@gromacs.org Date: Fri, 8 May 2009 11:17:33 +0200 Subject: [gmx-users] Re: v-rescale - harmonic oscillator Hi, Thanks for your confirmation of this problem. One thing is not clear to me from your answer. Has v-rescale fundamental problems with a Harmonic oscillator or is it a problem with the GROMACS implementation? Kind regards, Servaas Message: 3 Date: Fri, 8 May 2009 10:12:45 +0200 From: Berk Hess g...@hotmail.com Subject: RE: [gmx-users] Re: v-rescale - harmonic oscillator To: Discussion list for GROMACS users gmx-users@gromacs.org Message-ID: col113-w10f68330ea113929cbc3518e...@phx.gbl Content-Type: text/plain; charset=iso-8859-1 Hi, I have come to the conclusion that the v-rescale thermostat is only correct in the limit of large numbers of degrees of freedom, which is equivalent to a small scaling limit. The conservation of the conserved energy quantity improves as the system gets bigger. But I noticed that in 4.0 the energy conservation is not good. It can be improved significantly by changing line 164 or src/mdlib/update.c from vv = lg*(vn + f[n][d]*w_dt); to vv = lg*vn + f[n][d]*w_dt; There are some other similar lines in case you are PR p-coupling or cosine acceleration. I will probably not fix this for 4.0.5, but only for 4.1, so we do not change the reproducibility in minor revisions. Berk From: servaas.michielss...@student.kuleuven.be To: gmx-users@gromacs.org Date: Thu, 7 May 2009 15:27:37 +0200 Subject: [gmx-users] Re: v-rescale - harmonic oscillator Hi Berk, Thank you for your concern on this issue. It did not expect it to be a problem from the paper of Bussi et al., but it is of course possible. I am very curious to hear what they have to say about this. Are you contacting him or will I do it? For the rotational degrees of freedom would removing them with comm_mode=angular be a problem then because you are using restraints? (you use orire than?) Thanks again, Servaas Hi, With a diatomic molecule there is the issue that the rotation of the molecule has two degrees of freedom which are completely uncoupled from the rest of the system. Anyhow, I made a 1D harmonic oscillator using a position restraint and even there the energy does not seem to be conserved. I am not sure if a harmonic oscillator is a pathetic case, like with Nose-Hoover, or if this is an issue with the integrator implementation in Gromacs. I might need to ask Bussi about this. Berk From: servaas.michielss...@student.kuleuven.be To: gmx-users@gromacs.org Date: Thu, 7 May 2009 14:03:30 +0200 Subject: [gmx-users] v-rescale - harmonic oscillator Hi, It is just a diatomic molecule without other interactions (so only bonded interaction, LJ=0 and charges=0 on the atoms). What other details would you like to know? Should I send you the tpr file? Servaas Hi, No, only for small values is should be off, maybe tau_t 10 or 100 * delta_t. How did you set up your harmonic oscillator? Berk From: servaas.michielss...@student.kuleuven.be To: gmx-users@gromacs.org Date: Thu, 7 May 2009 12:54:54 +0200 Subject: [gmx-users] Re: v-rescale - harmonic oscillator Hi, Thank you for your reply. Yes by effective energy I mean the Gromacs conserved energy term. So you would say the cause of my problem here is that I need a small tau_t to thermostat this system but for small tau_t values the algorithm is not 100% ok yet? And what are very small tau_t values? I tested with quite a big range here. The problem occurred to me when I was experimenting with a hybrid monte carlo algorithm. I first tested it on a noble gas and compared to MD simulation with nose-hoover and v-rescale. Those results were 100% identical. Then I tried it on a dipeptide in vacuum, results were a little off here (I would expect this for nose-hoover thermostat, but not for v-rescale). So I decided to go the the extreme simple case of a harmonic oscillator. Those results were seriously off, and then I checked the conserved energy term, which seems impossible to keep it conserved for this system. (of course there can also
RE: [gmx-users] Re: v-rescale - harmonic oscillator
Hi, No, only for small values is should be off, maybe tau_t 10 or 100 * delta_t. How did you set up your harmonic oscillator? Berk From: servaas.michielss...@student.kuleuven.be To: gmx-users@gromacs.org Date: Thu, 7 May 2009 12:54:54 +0200 Subject: [gmx-users] Re: v-rescale - harmonic oscillator Hi, Thank you for your reply. Yes by effective energy I mean the Gromacs conserved energy term. So you would say the cause of my problem here is that I need a small tau_t to thermostat this system but for small tau_t values the algorithm is not 100% ok yet? And what are very small tau_t values? I tested with quite a big range here. The problem occurred to me when I was experimenting with a hybrid monte carlo algorithm. I first tested it on a noble gas and compared to MD simulation with nose-hoover and v-rescale. Those results were 100% identical. Then I tried it on a dipeptide in vacuum, results were a little off here (I would expect this for nose-hoover thermostat, but not for v-rescale). So I decided to go the the extreme simple case of a harmonic oscillator. Those results were seriously off, and then I checked the conserved energy term, which seems impossible to keep it conserved for this system. (of course there can also be an error in my hybrid monte carlo code, but the fact that the conserved energy is not conserved is disturbing me here) Servaas Hi, What do you mean with effective energy? The Gromacs conserved energy term? For very small tau_t the current implementation does not work well. Bussi mailed me a proper implementation that I will put in when I have time. Berk From: servaas.michielssens at student.kuleuven.be To: gmx-users at gromacs.org Date: Thu, 7 May 2009 09:02:52 +0200 Subject: [gmx-users] v-rescale - harmonic oscillator Hi, I did some experiments with a harmonic oscillator (diatomic molecule without charge en LJ parameters) using the v-rescale thermostat. First I ran a simulation in the NVE ensemble with a time-step of 0.0001ps, total energy was constant here. Then I tried simulating the system with the v-rescale thermostat, same time-step of 0.0001ps and tau_t=100,10,1,0.1,0.01,0.001. For none of those simulations the effective energy was conserved. (effective energy should be a conserved quantity for NVT simulations with v-rescale thermostat, JOURNAL OF CHEMICAL PHYSICSVolume: 126 Issue: 1 Article Number: 014101Published: JAN 7 2007 ) So what is going wrong here? Would one expect this thermostat to fail for such system (e.g. a simple Nose-Hoover is known to fail for the harmonic oscillator)? Thanks in advance, Servaas ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php _ Express yourself instantly with MSN Messenger! Download today it's FREE! http://messenger.msn.click-url.com/go/onm00200471ave/direct/01/___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php
[gmx-users] Re: v-rescale - harmonic oscillator
Hi Berk, Thank you for your concern on this issue. It did not expect it to be a problem from the paper of Bussi et al., but it is of course possible. I am very curious to hear what they have to say about this. Are you contacting him or will I do it? For the rotational degrees of freedom would removing them with comm_mode=angular be a problem then because you are using restraints? (you use orire than?) Thanks again, Servaas Hi, With a diatomic molecule there is the issue that the rotation of the molecule has two degrees of freedom which are completely uncoupled from the rest of the system. Anyhow, I made a 1D harmonic oscillator using a position restraint and even there the energy does not seem to be conserved. I am not sure if a harmonic oscillator is a pathetic case, like with Nose-Hoover, or if this is an issue with the integrator implementation in Gromacs. I might need to ask Bussi about this. Berk From: servaas.michielss...@student.kuleuven.be To: gmx-users@gromacs.org Date: Thu, 7 May 2009 14:03:30 +0200 Subject: [gmx-users] v-rescale - harmonic oscillator Hi, It is just a diatomic molecule without other interactions (so only bonded interaction, LJ=0 and charges=0 on the atoms). What other details would you like to know? Should I send you the tpr file? Servaas Hi, No, only for small values is should be off, maybe tau_t 10 or 100 * delta_t. How did you set up your harmonic oscillator? Berk From: servaas.michielss...@student.kuleuven.be To: gmx-users@gromacs.org Date: Thu, 7 May 2009 12:54:54 +0200 Subject: [gmx-users] Re: v-rescale - harmonic oscillator Hi, Thank you for your reply. Yes by effective energy I mean the Gromacs conserved energy term. So you would say the cause of my problem here is that I need a small tau_t to thermostat this system but for small tau_t values the algorithm is not 100% ok yet? And what are very small tau_t values? I tested with quite a big range here. The problem occurred to me when I was experimenting with a hybrid monte carlo algorithm. I first tested it on a noble gas and compared to MD simulation with nose-hoover and v-rescale. Those results were 100% identical. Then I tried it on a dipeptide in vacuum, results were a little off here (I would expect this for nose-hoover thermostat, but not for v-rescale). So I decided to go the the extreme simple case of a harmonic oscillator. Those results were seriously off, and then I checked the conserved energy term, which seems impossible to keep it conserved for this system. (of course there can also be an error in my hybrid monte carlo code, but the fact that the conserved energy is not conserved is disturbing me here) Servaas Hi, What do you mean with effective energy? The Gromacs conserved energy term? For very small tau_t the current implementation does not work well. Bussi mailed me a proper implementation that I will put in when I have time. Berk From: servaas.michielssens at student.kuleuven.be To: gmx-users at gromacs.org Date: Thu, 7 May 2009 09:02:52 +0200 Subject: [gmx-users] v-rescale - harmonic oscillator Hi, I did some experiments with a harmonic oscillator (diatomic molecule without charge en LJ parameters) using the v-rescale thermostat. First I ran a simulation in the NVE ensemble with a time-step of 0.0001ps, total energy was constant here. Then I tried simulating the system with the v-rescale thermostat, same time-step of 0.0001ps and tau_t=100,10,1,0.1,0.01,0.001. For none of those simulations the effective energy was conserved. (effective energy should be a conserved quantity for NVT simulations with v-rescale thermostat, JOURNAL OF CHEMICAL PHYSICSVolume: 126 Issue: 1 Article Number: 014101Published: JAN 7 2007 ) So what is going wrong here? Would one expect this thermostat to fail for such system (e.g. a simple Nose-Hoover is known to fail for the harmonic oscillator)? Thanks in advance, Servaas ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/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/mailing_lists/users.php ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search
