<div>Hi,</div><div>what I am trying to do is implementing the REST2 tecnique [1]</div><div>The goal is enanching the conforlational sampling splitting the system in two parts: solvent and solute.</div><div>The two groups have a different APPARENT temperature because the bonded interaction within the solute atoms are rescaled by a factor Bi/B ( Bi=1/Ti*Kb, Ti = apparent temperature for the i-th replica, B=Boltzmann constant) and the protein-water interactions are rescaled by sqrt(Bi/B). So the system temperature, regulated by the thermostat, is equal along the different replicas but the dynamic of the different replicas</div>
<div>changes becaus of the different Boltzmann factor.</div><div>To implement the tecnique, the naive approach (adopted in [2]) is:</div><div>1) Generating N-different topology (N= number of replicas).</div><div>2) Removing same check at the beginning of the Replica Exchange code, in order to run N replicas with equal run parameters (temperature, pression, ecc. ecc.) </div>
<div>3) Changing the acceptance ratio formula.</div><div><br></div><div>This approach has been developped for gmx4.0.3 but I need to use gmx4.5.5 because it has a native implementation of the CHARMM force-field. Since the code between the two versions is different I am trying to figure out how to reinvent the wheel.</div>
<div><br></div><div>As shown in [3] and adopted in [2] the same result can be obtained using the lamda dynamics in gromacs.In this case, only one topology has to be generated. The state A correspond to the replica at the lowest apparent temperature, while the state B represents the highest temperature replica. Intermediate replicas are generated changing the init_lambda value. </div>
<div><br></div><div>What has been found, this approach results in slowest run.</div><div><br></div><div><br></div><div>So, right now I got stuck at the point 3. In fact to find:</div><div><br></div><div>delta = B( Va(Xb) + Vb(Xa) - Va(Xa) - Vb(Xb)) </div>
<div><br></div><div>I need to calculate Va(Xb) and Vb(Xa)</div><div><br></div><div>[1] <a href="http://pubs.acs.org/doi/abs/10.1021/jp204407d">http://pubs.acs.org/doi/abs/10.1021/jp204407d</a></div><div>[2]<a href="http://pubs.acs.org/doi/abs/10.1021/ct100493v">http://pubs.acs.org/doi/abs/10.1021/ct100493v</a> </div>
<div>[3] <a href="http://onlinelibrary.wiley.com/doi/10.1002/jcc.21703/abstract">http://onlinelibrary.wiley.com/doi/10.1002/jcc.21703/abstract</a></div><div><br></div><div><br></div><div><br></div><div><br></div><br><div class="gmail_quote">
2012/7/16 Berk Hess <span dir="ltr"><<a href="mailto:hess@kth.se" target="_blank">hess@kth.se</a>></span><br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div class="HOEnZb">
<div class="h5">On 7/16/12 14:18 , Shirts, Michael (mrs5pt) wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
Actually I did it but I the performance decrease a lot (2x).<br>
In particular I run the same run with free_energy = yes and free_energy =<br>
no and in the second case the run is 2times faster.<br>
I found in the mailing list that this performance drop is known.<br>
</blockquote>
What are the particular parameters you are trying to change? I wrote the<br>
new Hamiltonian Exchange code, so I'm interested in making sure it's as<br>
useful as possible. There might be some cases where this performance drop<br>
can be avoided, and It's likely the PME parameters, since doing repeated<br>
that's the only thing that really makes free energy slow.<br>
<br>
UNLESS of course you are changing ALL of the atoms, which which case all of<br>
them will go through the free energy code, and it certainly will be slower.<br>
There are some possible ways to fix this, but that will not be until 5.0.<br>
Best,<br>
~~~~~~~~~~~~<br>
Michael Shirts<br>
Assistant Professor<br>
Department of Chemical Engineering<br>
University of Virginia<br>
<a href="mailto:michael.shirts@virginia.edu" target="_blank">michael.shirts@virginia.edu</a><br>
(434)-243-1821<br>
<br>
<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
From: francesco oteri <<a href="mailto:francesco.oteri@gmail.com" target="_blank">francesco.oteri@gmail.com</a>><br>
Reply-To: Discussion list for GROMACS development <<a href="mailto:gmx-developers@gromacs.org" target="_blank">gmx-developers@gromacs.org</a>><br>
Date: Mon, 16 Jul 2012 14:08:49 +0200<br>
To: Discussion list for GROMACS development <<a href="mailto:gmx-developers@gromacs.org" target="_blank">gmx-developers@gromacs.org</a>><br>
Subject: Re: [gmx-developers] calculating poteintial energy into do_md()<br>
<br>
Dear Berk<br>
<br>
2012/7/16 Berk Hess <<a href="mailto:hess@kth.se" target="_blank">hess@kth.se</a>><br>
<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
On 7/16/12 13:05 , francesco oteri wrote:<br>
<br>
Dear gromacs developers,<br>
I am trying to implement in gromacs 4.5.5 a particular Hemiltonian Replica<br>
Exchange tecnique.<br>
Right now, at every exchange attempt in do_md, I figured out how to<br>
access at the potential<br>
energy of replica A (=configuration A at temperature A) and B<br>
(=configuration B at temperature B) and so on.<br>
<br>
In my case, each replica has the same temperature, but there is a<br>
different Hemiltonian equation for every replica.<br>
The different Hemiltonian are obtained simply changing the force field<br>
parameters in the input topology so I dont<br>
need to modify anything in gromacs.<br>
<br>
But at every exchange attempt I have to test if the configuration B can<br>
exist in the state A so I need to calculate<br>
its potential energy using the force field data of replica A.<br>
<br>
I found that function sum_epot calculates the potential energy but I<br>
suspect that it uses values calculated in<br>
do_force since sum_epot is called by do_force_lowlevel in turn called by<br>
do_force.<br>
<br>
<br>
<br>
So my question is, should I call do_force in replica A with coordinates<br>
from replica B reach my goal?<br>
<br>
Yes.<br>
Are you changing bonded and/or non-bonded parameters?<br>
Some non-bonded parameters might be preprocessed, so you might need<br>
reprocess those<br>
and them reprocesses again to get back to the original state.<br>
<br>
<br>
</blockquote>
Which parameters need such a preprocessing?<br>
<br>
<br>
<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
Note that 4.6 will have proper Hamiltonian replica exchange, but that will<br>
use the lambda coupling<br>
parameter approach. If you need to do something similar, it might be much<br>
simpler to use this code.<br>
<br>
<br>
</blockquote>
Actually I did it but I the performance decrease a lot (2x).<br>
In particular I run the same run with free_energy = yes and free_energy =<br>
no and in the second case the run is 2times faster.<br>
I found in the mailing list that this performance drop is known.<br>
</blockquote></blockquote></div></div>
I think the performance drop is due to the, temporarily, missing sse non-bonded kernels.<br>
They should be back somewhere this week.<br>
<br>
Cheers,<br>
<br>
Berk<div class="HOEnZb"><div class="h5"><br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
Since I dont need to change lambda along the simulation, I prefer using a<br>
modified input topology and running a normal REMD changing something in<br>
repl_ex.c and md.c.<br>
<br>
Of course, any suggestions to improve performance are wellcome :)<br>
<br>
<br>
<br>
<br>
Cheers,<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
Berk<br>
<br>
<br>
Thanks in advance,<br>
<br>
Francesco<br>
<br>
<br>
<br>
<br>
<br>
<br>
--<br>
gmx-developers mailing list<br>
<a href="mailto:gmx-developers@gromacs.org" target="_blank">gmx-developers@gromacs.org</a><br>
<a href="http://lists.gromacs.org/mailman/listinfo/gmx-developers" target="_blank">http://lists.gromacs.org/<u></u>mailman/listinfo/gmx-<u></u>developers</a><br>
Please don't post (un)subscribe requests to the list. Use the<br>
www interface or send it to <a href="mailto:gmx-developers-request@gromacs.org" target="_blank">gmx-developers-request@<u></u>gromacs.org</a>.<br>
<br>
</blockquote>
<br>
<br>
-- <br>
Cordiali saluti, Dr.Oteri Francesco<br>
-- <br>
gmx-developers mailing list<br>
<a href="mailto:gmx-developers@gromacs.org" target="_blank">gmx-developers@gromacs.org</a><br>
<a href="http://lists.gromacs.org/mailman/listinfo/gmx-developers" target="_blank">http://lists.gromacs.org/<u></u>mailman/listinfo/gmx-<u></u>developers</a><br>
Please don't post (un)subscribe requests to the list. Use the<br>
www interface or send it to <a href="mailto:gmx-developers-request@gromacs.org" target="_blank">gmx-developers-request@<u></u>gromacs.org</a>.<br>
</blockquote></blockquote>
<br>
<br>
-- <br>
gmx-developers mailing list<br>
<a href="mailto:gmx-developers@gromacs.org" target="_blank">gmx-developers@gromacs.org</a><br>
<a href="http://lists.gromacs.org/mailman/listinfo/gmx-developers" target="_blank">http://lists.gromacs.org/<u></u>mailman/listinfo/gmx-<u></u>developers</a><br>
Please don't post (un)subscribe requests to the list. Use the www interface or send it to <a href="mailto:gmx-developers-request@gromacs.org" target="_blank">gmx-developers-request@<u></u>gromacs.org</a>.<br>
</div></div></blockquote></div><br><br clear="all"><div><br></div>-- <br>Cordiali saluti, Dr.Oteri Francesco<br>