[gmx-users] Re: restraint-lambdas for position restraints in hamiltonian exchange

Dejun Lin dejun.lin at gmail.com
Mon Aug 12 23:57:17 CEST 2013


Hi Michael,

After looking at the codes, it seems to me that the exchange scheme I
posted previously incurs violation of detailed balance, which seems to be a
bug. I'll really appreciate it if you can clarify this.

Let me re-state what I was trying to do. I have the following 7-component
restraint-lambdas:

0 0.1 0.2 0.3 0.09 0.19 0.29

The restraints are position restraints defined in top file of molecule I:
#define posA    0       0
0
#define posB    0       0
3
[ position_restraints ]
;       ai      funct   parA
parB
        1       1       posA
posB
        2       1       posA
posB
        3       1       posA
posB
        4       1       posA    posB


I initialized replicas 0 to 3 in structure C and 4 to 6 in structure D and
hoped that replica 0 would sample structure C and D according to the
correct thermodynamics.

Now suppose an exchange is attempt between replica i and j, which are
running in Hamiltonian Ui and Uj respectively.

Then dE should be Ui(xj) - Uj(xj) + Uj(xi) - Ui(xi), where Ui(xj) =
(1-lambda_i)*POSRES_A_i(xj) + lambda_i*POSRES_B_i(xj) =
lambda_i*POSRES_B_i(xj) (last equality holds since the force constants are
all 0 in state A as shown above) for any pair of i and j, where POSRES_B_i
is the state-B restraint potential energy function with reference position
defined in replica i.
So dE = lambda_i*POSRES_B_i(xj) - lambda_j*POSRES_B_j(xj) +
lambda_j*POSRES_B_j(xi) - lambda_i*POSRES_B_i(xi)

But the respective terms Ui(xj) - Uj(xj) and Uj(xi) - Ui(xi) are actually
calculated in the 2 simulations indepedently (function
test_for_replica_exchange() in kernel/repl_ex.c around line 954) and
POSRES_B_i(xj) is actually POSRES_B_j(xj) in simulation j while
POSRES_B_j(xi) is actually POSRES_B_i(xi) in simulation i because
simulation i and j don't communicate for reference position in this case as
far as I can tell from the function real posres() as defined in
include/bondf.h. The codes was able to run but I think this will give wrong
detail balance and it would be nice for the program to issue a
warning/error in this case.

Again, thanks a lot for your help!

-Dejun


2013/8/7 Dejun Lin <dejun.lin at gmail.com>

> Hi Micheal,
>
> Sorry for keep bugging you about this but I want to make sure I'm doing
> what I think I'm doing :)
> I did some test on the following setup.
>
> restraint-lambdas = 0.3 0.299999 (essentially no difference) on position
> restraints.
> replica 0 starts in structure A and replica 1 in structure B. The 2
> position restraints have exactly the same set of force constants, i.e., the
> 2 replica used the same top/itp file.
> But each of this 2 replicas should have one distinct set of restraints
> because the difference in initial structure (gro file) are huge.
>
> Replica 0:  structure *A*, restraints A (ref. position A, force constant
> C)
> Replica 1:  structure *B*, restraints B (ref. position B, force constant
> C)
>
> Swapped between the 2 replica were observed and I looked at the trajectory
> from *replica 0* in VMD and noticed that there were *structure B* in it.
>
> Can assume that swap resulted in the following states ?
>
> Replica 0:  structure *B*, restraints A (ref. position A, force constant
> C)
> Replica 1:  structure *A*, restraints B (ref. position B, force constant
> C)
>
> i.e. It was the coordinates instead of lambda values were swapped between
> replicas, which is essentially the same as if the reference positions were
> swapped.
>
> Thanks again for your help!
>
> -Dejun
>
>
> 2013/8/3 Michael Shirts-2 [via GROMACS] <
> ml-node+s5086n5010326h5 at n6.nabble.com>
>
>  That is correct. Such a functionality wouldn't be that hard to
>> implement - but there are a long list of easy functionalities to be
>> implemented. You can submit a request to redmine.gromacs.org so that
>> the request is archived.
>>
>> On Sat, Aug 3, 2013 at 6:12 PM, Dejun Lin <[hidden email]<http://user/SendEmail.jtp?type=node&node=5010326&i=0>>
>> wrote:
>>
>> > So I take it that in the position restraint case (not COM-pulling),
>> where
>> > the reference positions are determined by the starting structure
>> instead of
>> > a B-state topology, the reference positions won't be swapped ?
>> >
>> >
>> > 2013/8/3 Michael Shirts-2 [via GROMACS] <
>> > [hidden email] <http://user/SendEmail.jtp?type=node&node=5010326&i=1>>
>> >
>> >> Short answer is anything that has a B state parameter can be included
>> >> in in Hamiltonian exchange.
>> >>
>> >> If it's pull code or explicit restraints, it's controlled by restraint
>> >> lambda.
>> >>
>> >> > I went through the manual and couldn't find any definite answers to
>> the
>> >> > following questions.
>> >> >
>> >> > First, I wonder if the reference positions of position restraints,
>> not
>> >> just
>> >> > the force constants, of different replicas are exchanged in
>> hamiltonian
>> >> > exchange based on restraint-lambdas? For example, if I have 1
>> molecule
>> >> that
>> >> > has 2 structures, say, A and B and the following 2-component
>> >> > restraint-lambdas:  1.0 1.0 with replica 0 starting in A and replica
>> 1
>> >> in B,
>> >> > would the exchange between replica 0 and 1 yield anything meaningful
>> ?
>> >>
>> >> Right now, the pull code doesn't have B state positions, just force
>> >> constants. It CAN be done using distance restraints, but then you have
>> >> to make the atoms involved part of the same molecule.  Adding B state
>> >> distances to pull code was recently requested, and will probably be
>> >> straightforward to put in soon.
>> >>
>> >> > In
>> >> > other word, would the relaxation towards the other structure occur
>> in
>> >> both
>> >> > states once an exchange was accepted?
>> >> >
>> >> > Another question is what types of restraints can restraint-lambdas
>> act
>> >> on?
>> >> > For example, bond/angle/distance/position restraints, COM-pulling
>> >> potential
>> >> > ?
>> >>
>> >> Anything with a B state in gromacs topologies.
>> >> --
>> >> gmx-users mailing list    [hidden email]<
>> http://user/SendEmail.jtp?type=node&node=5010324&i=0>
>> >> 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 [hidden email]<
>> http://user/SendEmail.jtp?type=node&node=5010324&i=1>.
>> >>
>> >> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>> >>
>> >>
>> >> ------------------------------
>> >>  If you reply to this email, your message will be added to the
>> discussion
>> >> below:
>> >>
>> >>
>> http://gromacs.5086.x6.nabble.com/restraint-lambdas-for-position-restraints-in-hamiltonian-exchange-tp5010322p5010324.html
>> >>  To unsubscribe from restraint-lambdas for position restraints in
>> >> hamiltonian exchange, click here<
>>
>> >> .
>> >> NAML<
>> http://gromacs.5086.x6.nabble.com/template/NamlServlet.jtp?macro=macro_viewer&id=instant_html%21nabble%3Aemail.naml&base=nabble.naml.namespaces.BasicNamespace-nabble.view.web.template.NabbleNamespace-nabble.view.web.template.NodeNamespace&breadcrumbs=notify_subscribers%21nabble%3Aemail.naml-instant_emails%21nabble%3Aemail.naml-send_instant_email%21nabble%3Aemail.naml>
>>
>> >>
>> >
>> >
>> >
>> >
>> > --
>> > View this message in context:
>> http://gromacs.5086.x6.nabble.com/restraint-lambdas-for-position-restraints-in-hamiltonian-exchange-tp5010322p5010325.html
>> > Sent from the GROMACS Users Forum mailing list archive at Nabble.com.
>> > --
>> > gmx-users mailing list    [hidden email]<http://user/SendEmail.jtp?type=node&node=5010326&i=2>
>> > 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 [hidden email]<http://user/SendEmail.jtp?type=node&node=5010326&i=3>.
>>
>> > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>> --
>> gmx-users mailing list    [hidden email]<http://user/SendEmail.jtp?type=node&node=5010326&i=4>
>> 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 [hidden email]<http://user/SendEmail.jtp?type=node&node=5010326&i=5>.
>>
>> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>>
>>
>> ------------------------------
>>  If you reply to this email, your message will be added to the
>> discussion below:
>>
>> http://gromacs.5086.x6.nabble.com/restraint-lambdas-for-position-restraints-in-hamiltonian-exchange-tp5010322p5010326.html
>>  To unsubscribe from restraint-lambdas for position restraints in
>> hamiltonian exchange, click here<http://gromacs.5086.x6.nabble.com/template/NamlServlet.jtp?macro=unsubscribe_by_code&node=5010322&code=ZGVqdW4ubGluQGdtYWlsLmNvbXw1MDEwMzIyfDE0MzAxNDIxNA==>
>> .
>> NAML<http://gromacs.5086.x6.nabble.com/template/NamlServlet.jtp?macro=macro_viewer&id=instant_html%21nabble%3Aemail.naml&base=nabble.naml.namespaces.BasicNamespace-nabble.view.web.template.NabbleNamespace-nabble.view.web.template.NodeNamespace&breadcrumbs=notify_subscribers%21nabble%3Aemail.naml-instant_emails%21nabble%3Aemail.naml-send_instant_email%21nabble%3Aemail.naml>
>>
>
>




--
View this message in context: http://gromacs.5086.x6.nabble.com/restraint-lambdas-for-position-restraints-in-hamiltonian-exchange-tp5010322p5010480.html
Sent from the GROMACS Users Forum mailing list archive at Nabble.com.



More information about the gromacs.org_gmx-users mailing list