[AMBER] How to define reaction coordinate in Umbrella sampeling

From: Sindrila Dutta banik <sindrila.duttabanik.yahoo.com>
Date: Mon, 1 Jul 2013 19:42:59 +0800 (SGT)

Thanks a lot for explaining all my questions so nicely. 

I have another question. As I described in my previous mail I want to study a reaction which involve a shifting of a H-atom from a N-atom to O-atom. For the same first I equilibrate the structure using SHAKE algorithm for the bond involving H-atom. Now for the umbrella sampling what shall I do? If I use SHAKE for the umbrella sampling then the N-H bond will be constrain. 

With best regards


 From: Jason Swails <jason.swails.gmail.com>
To: Sindrila Dutta banik <sindrila.duttabanik.yahoo.com>
Cc: AMBER Mailing List <amber.ambermd.org>
Sent: Saturday, 29 June 2013 6:43 PM
Subject: Re: [AMBER] How to define reaction coordinate in Umbrella sampeling

On Sat, Jun 29, 2013 at 1:52 AM, Sindrila Dutta banik <sindrila.duttabanik.yahoo.com> wrote:

Thanks again for your reply. 
>In the previous mail you suggest to use rstwt = -1.0, 1.0;  can I use another value for the same. The -1 and 1 indicates that simple difference of two separation, if I use 2.0 then it will be difference of double separation. Am I right?

-2.0, 2.0 is functionally equivalent to -1.0,1.0, so there's not much point in 'doubling' IMO.  -2.0,1.0 is a different reaction coordinate (and may be preferable if the first reaction coordinate is expected to change only half as much as the second reaction coordinate, for instance).

If you use 2.0 for both, then that generalized distance would double as well, yes. (I encourage you not to think of the generalized distance as a 'separation', since it's a linear combination of two separations).

I think by using (x_2 - x_1) you indicate the separation between atom1 and atom2. But in the input I have to specify r1, r2, r3, r4; which is the difference of separation, but what will be the value of individual separation? 
r1, r2, r3, and r4 are points in the flat-well potential.  They correspond not to separations, but to the value of the generalized distance restraint.  The values I gave in my example

LCOD of 2 distances
   iat=1,2,3,4, rk2=5.0, rk3=5.0, r1=-10, r2=0, r3=0, r4=10, rstwt=-1.0,1.0,

correspond to aharmonic potential with a force constant of 5.0 kcal/mol (in the energy equation U=k(x-x_0)^2 -- there is no 1/2 factor) that is harmonic between -10 and 10 with the center of the potential at 0.  The generalized coordinate is zero when the 2--1 separation is equal to the 3--4 separation.

Note that a generalized distance restraint imposes a single restraint on two degrees of freedom.  In general, a single value of the generalized restraint cannot be broken down into specific individual separations.  There are an infinite number of individual separations that can yield a specific value for their linear combination.  For instance, if 1--2 were separated by 2 Å and so was 3--4, the generalized coordinate would have a value of 0.  Likewise if both were separated by 1 Å.  If you want to control both coordinates independently, you will need to specify two restraints and generate a free energy surface using 2D-WHAM.


>Actually, in my case one H-atom is bonded with one N-atom and during the reaction the H-atom will shift towards another O-atom. I want to increase the N-H bond length and decrease H-O separation and set the reaction coordinate as the difference of two separation. So my doubt is that if I specify the r1, r2... (difference of separation) then two separation with change simultaneously, or it will change only one.

They will change simultaneously in the way that minimizes the energy (in principle).
>It is not clear to me what is the difference in sander and pmemd? 

This is described in the Amber manual.  sander and pmemd give the same answers for the same input, in general.  But pmemd does not support all of the features that sander supports.  Therefore, since pmemd fails when rstwt is added to the &rst namelist, that means that pmemd does not support LCOD restraints, and you have to use sander.

In general, pmemd is ~2x faster than sander on a single core (in explicit solvent), and the parallel scaling is much better.  As a result, when a feature you want is implemented in pmemd, you should use that program.  (But when it's not, you must stick to sander)

Jason M. Swails
Quantum Theory Project,
University of Florida
Ph.D. Candidate
AMBER mailing list
Received on Mon Jul 01 2013 - 05:00:02 PDT
Custom Search