Re: [AMBER] How to keep Graphene sheet fixed during MD simulation.

From: Jason Swails <jason.swails.gmail.com>
Date: Mon, 7 Nov 2011 08:42:27 -0500

2011/11/7 hai wei <weihy61.gmail.com>

> ”Third, 500 is too strong a force constant for most purposes. Consider
> using
> a value about two orders of magnitude smaller.“
>
> If I want to keep Graphene sheet Constrain in the MD simulaiton, is 500
> still too much? I dont quite understand.
>

Yes, 500 kcal/mol/Ang^2 is a *huge* force constant. It's on the scale of
H-heteroatom bond strengths, which have to be SHAKEn due to their high
frequency motion for a long time step. Thus, with such high force
constants, you'll very likely get integration errors.

You should consider using a force constant between 0 and 10. 10 is very
stiff, and penalizes moving a *single* atom 1 Angstrom away with a 10
kcal/mol penalty! No simulation at a typical temperature can afford that
penalty barring extreme clashes. Even with 1 kcal/mol, you'll get the same
penalty if only 10 atoms move 1 angstrom away from their restrained
position.

Also, are you using ibelly=1, or ntr=1? You definitely shouldn't be using
both.

HTH,
Jason


> THanks
> haiya
>
>
> 在 2011年11月7日 下午2:30,hai wei <weihy61.gmail.com>写道:
>
> > Dear Case:
> >
> > DO I have to set *nscm=0 *in the control. Since I am reading one of the
> > old post, says that "*One thing to try is *
> > *> make sure you set nscm=0 when doing Belly runs*"
> >
> > ***> *
> > *> I thought I would add my 2c to this discussion. Firstly I just want to
> > give *
> > *> a caution with regards to the use of ibelly. That is that you really
> > not use *
> > *> it!!! It simply zeros the forces and is in no way correct or physical.
> > *
> > *> Especially if you start doing things like NPT simulations with it. In
> > this *
> > *> case your virial will be wrong, you pressure will be wrong and thus
> > your *
> > *> density will be wrong. If you must restrain things then you should
> > probably *
> > *> use harmonic restraints. *
> > *> *
> > *> Note, things like seeing a 200K temperature drop on the first step
> > makes *
> > *> sense with using belly for such a huge portion of your syste, If you
> > are *
> > *> starting from a restart file where things were not frozen and then you
> > *
> > *> suddenly freeze a chunk of it you have just magically destroyed a ton
> > of *
> > *> kinetic energy... *
> > *> *
> > *> With regards to the problems you are seeing this could be because of
> > all *
> > *> sorts of issues including trying to do NPT with ibelly. One thing to
> > try is *
> > *> make sure you set nscm=0 when doing Belly runs, Sander may be doing
> > this *
> > *> while pmemd is not. Second make sure you do not try to use langevin
> > dynamics *
> > *> (ntt=3) with ibelly since it will likely not work. *
> > *> *
> > *> All the best *
> > *> Ross *
> >
> > THanks.
> >
> > haiya
> >
> >
> > 2011/11/7 hai wei <weihy61.gmail.com>
> >
> >> Dear Case:
> >>
> >> Thanks for your reply, and apologize for the unclear mail I send, My
> >> system is compsed of Protein (res number 1 582), Graphene (res number
> 583
> >> 605) and counterion and water.
> >>
> >> I am doing simulaiton as following step:
> >>
> >> 1, minimization of the system:
> >> &cntrl
> >> imin = 1, maxcyc = 2500, ncyc = 1000, ntb = 1, ntr =
> >> 1, cut = 10,
> >> /
> >> Hold the protein and graphen sheet fixed
> >> 500.0
> >> RES 1 605
> >> END
> >> END
> >>
> >> 2. Heating the system to 360K:
> >> &cntrl
> >> imin = 0, irest = 0, ntx = 1, ntp = 0, ntb = 1,
> >> cut = 10,
> >> ibelly = 1,
> >> ntc = 2, ntf = 2, tempi = 0.0, temp0 = 360.0,
> >> ntt = 2, (ibelly does not work with ntt=3,)
> >> gamma_ln = 1.0, nstlim = 10000, dt = 0.002, ntpr = 250, ntwx = 250,
> >> ntwr = 10000,
> >> iwrap = 1, ioutfm=1, nscm = 100, /
> >>
> >> Keep protein and Graphen sheet fiexed with ibelly
> >> RES 1 605
> >> END
> >> END
> >>
> >> 3. 1ns MD simulation at 360K ;
> >> &cntrl
> >> imin = 0, irest = 1, ntx = 5, ntb = 2, pres0 = 1.0, ntp = 1, taup
> =
> >> 2.0, cut = 10,
> >> ibelly = 1,
> >> ntc = 2, ntf = 2, tempi = 360.0, temp0 = 360.0,
> >> ntt = 2,
> >> gamma_ln = 1.0, nstlim = 500000, dt = 0.002, ntpr = 250, ntwx = 250,
> >> ntwr = 10000,
> >> nscm = 100, ioutfm = 1, iwrap=1,
> >> /
> >> Keep protein and Graphen sheet fiexed with ibelly
> >> RES 1 605
> >> END
> >> END
> >>
> >> 4. cooling the system to 300k ;
> >> &cntrl
> >> imin = 0, irest = 1, ntx = 5, ntp = 0, ntb = 1,
> >> cut = 10,
> >> ibelly = 1,
> >> ntc = 2, ntf = 2, tempi = 360.0, temp0 = 300.0,
> >> ntt = 2,
> >> gamma_ln = 1.0, nstlim = 100000, dt = 0.002, ntpr = 250, ntwx = 250,
> >> ntwr = 10000,
> >> iwrap = 1, ioutfm=1, nscm = 100,
> >> /
> >> Keep protein and Graphen sheet fiexed with ibelly
> >> RES 1 605
> >> END
> >> END
> >>
> >> 5. MD production run at 300k ;
> >> &cntrl
> >> imin = 0, irest = 1, ntx = 5, ntb = 2, pres0 = 1.0, ntp = 1, taup
> =
> >> 2.0,
> >> cut = 10,
> >> ntr = 0,
> >> ntc = 2, ntf = 2, tempi = 300.0, temp0 = 300.0,
> >> ntt = 2,
> >> gamma_ln = 1.0,
> >> nstlim = 500000, dt = 0.002, ntpr = 250, ntwx = 250, ntwr = 10000,
> >> nscm = 100, ioutfm = 1, iwrap=1,
> >> /
> >> Keep graphene sheet fixed with ibelly
> >> RES 583 605
> >> END
> >> END
> >>
> >> After I read the ambermail list of ibelly question, like, not work with
> >> ntt=3, pememd not work, I think the above control should run normally.
> >>
> >> While I still got a problem here, when I heat the system, How to just
> put
> >> a weak restrain force on protein , while keep graphene sheet ibelly
> >> constrain, Since ntr=1, and ibelly=1 can not be used at the same time.
> >>
> >> THanks.
> >>
> >> haiya
> >>
> >>
> >>
> >>
> >>
> >>
> >>
> >>
> >>
> >>
> >> 2011/11/7 case <case.biomaps.rutgers.edu>
> >>
> >>> On Mon, Nov 07, 2011, hai wei wrote:
> >>> >
> >>> > I am runing a MD simulation with protein on top of Graphene sheet,
> >>> During
> >>> > the simulaiton, I try to fix the graphene sheet while adding the
> >>> following
> >>> > input:
> >>> >
> >>> > hold the protein fixed
> >>> > 500.0
> >>> > RES num1 num2
> >>> > END
> >>> > END
> >>> >
> >>> > While I found that the Graphene sheet does not fixed as it should be,
> >>> There
> >>> > are still some ripples shown.
> >>>
> >>> First, the comment indicates that you are holding the *protein* fixed,
> >>> whereas
> >>> your email text indicates that you wish to keep the graphene fixed.
> >>> Without
> >>> knowing anything about your system (and not knowing what num1 and num2
> >>> are),
> >>> I can't tell what should be happening.
> >>>
> >>> Second, you also need to set ntr=1 in the &cntrl namelist; otherwise
> >>> lines
> >>> like those above will be ignored. The output will summarize the
> >>> restraints,
> >>> and tell you how many atoms were involved. If there is nothing in the
> >>> about
> >>> about restrained atoms, you need to re-examine your inputs.
> >>>
> >>> Third, 500 is too strong a force constant for most purposes. Consider
> >>> using
> >>> a value about two orders of magnitude smaller.
> >>>
> >>> ...good luck....dac
> >>>
> >>>
> >>> _______________________________________________
> >>> AMBER mailing list
> >>> AMBER.ambermd.org
> >>> http://lists.ambermd.org/mailman/listinfo/amber
> >>>
> >>
> >>
> >
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
>



-- 
Jason M. Swails
Quantum Theory Project,
University of Florida
Ph.D. Candidate
352-392-4032
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Mon Nov 07 2011 - 06:00:04 PST
Custom Search