Re: [AMBER] Problem related to the umbrella sampling

From: Sindrila Dutta banik <>
Date: Fri, 12 Jul 2013 02:27:27 +0800 (SGT)

Thanks for your kind reply. 

Yes, these two bonds involve in the proton transfer reaction.  

As you told in your previous mail: "You cannot make or break bonds using a force field" I also agree with this.
But I could not understand how people calculate PMF for proton transfer process using classical mechanics using umbrella sampling. There are many papers where people calculate energy barrier using classical mechanics as well as QM and compare the barrier height. 

I look at atoms 10018 and 10019 in my system. But every thing is fine with this atoms. But the substrate (PPD for which noshake is applied) is deformed during this run. And the above mentioned atoms are not the part of substrate. Is the increase in bond length from 1.1A to 1.7A is too large? Is the problem arises due to large increment in separation? 

With best regards


 From: Jason Swails <>
To: Sindrila Dutta banik <>; AMBER Mailing List <>
Sent: Thursday, 11 July 2013 10:52 PM
Subject: Re: [AMBER] Problem related to the umbrella sampling

On Thu, Jul 11, 2013 at 12:44 AM, Sindrila Dutta banik <> wrote:

> Dear all,
> I am trying to do umbrella sampling considering two reaction coordinates:
> increase one N-H bond length (from 1.1A to 1.8A) and decrease one O-H bond
> (from 1.8 to 1.1A) (which are part of the substrate PPD). To do the same I
> use the following input for equilibration step.

Is this part of a proton transfer in which one bond breaks while another
forms?  If not, what are you trying to learn?

> 10 ps equilibration run
> &cntrl
>  imin  = 0,
>  irest  = 0,
>  ntx    = 1,
>  ntb    = 1,
>  ntc    = 2,
>  noshakemask = ":PPD"
>  ntf    = 2,
>  cut    = 10,
>  tempi  = 300.0,
>  temp0  = 300.0,
>  ntt    = 3,
>  gamma_ln = 1.0,
>  nstlim = 10000, dt = 0.001
>  ntpr = 100, ntwx =100, ntwr =200
>  nmropt = 1, ioutfm=1
>  /
> &wt
>  type='END',
>  &end
> DISANG=comb7
> The comb7 file:
> Two restraints N-H and O-H
>  &rst
>  iat=12230,12245
>  r1=0.0, r2=1.7, r3=1.7, r4=10.0,
>  rk2=1000, rk3=1000
>  &end
>  &rst
>  iat=12219,12245
>  r1=0.0, r2=1.2, r3=1.2, r4=10.0,
>  rk2=1000, rk3=1000,
>  &end
> When I increase N-H bond length from 1.1 upto 1.6, there was no problem.
> When I set it 1.7 as given in the above input file it gives the following
> error.
> Coordinate resetting (SHAKE) cannot be accomplished,
>      deviation is too large
>      NITER, NIT, LL, I and J are :      0      3  4945  10018  10019
>      Note: This is usually a symptom of some deeper
>      problem with the energetics of the system.

You cannot make or break bonds using a force field.  You need to treat a
reacting center with some type of QM method.  You should look at atoms
10018 and 10019 in your system to see if you can figure out why SHAKE may
be failing.  If Amber is trying to SHAKE the bonds you're modifying with
harmonic potentials (which I suspect is happening), this would easily
explain your error.

That said, what your input files are doing doesn't make too much sense to
me (for the reasons mentioned above)...

Also, the usual debugging information here still applies -- try printing
out every step (coordinates, energies, velocities, anything you need) to
isolate the exact step where things break.  Such an approach will often
yield insight into why a problem is occurring (and possibly even how to fix

Good luck,

Jason M. Swails
Quantum Theory Project,
University of Florida
Ph.D. Candidate
AMBER mailing list
AMBER mailing list
Received on Thu Jul 11 2013 - 11:30:03 PDT
Custom Search