Hello Ross. Thank you for the reply. I´ll restart the parameterization
procedure. You suggest I create a new H type. I understand that
setting a new type for this hydrogen would require me to include
vdwradii for this atom. The problem is that I am not sure how I would
do this. I would just have to rename the pdb source file? But then,
how would I set the vdwradii? Can´t I just set the vdwradii for this
particular atom in tleap, after having performed the parameterization?
Wouldn´t the lmodprmtop program be of use in this case?.
Thank you in advance
Fabrício Bracht
2012/6/28 Ross Walker <ross.rosswalker.co.uk>:
> Hi Fabricio,
>
> So you are running into one of the problems with the AMBER force field that
> people have ignored for years and it was never too much of a problem since
> people didn't typically run for long enough to see it. That is that hydroxyl
> hydrogens have a zero VDW radii. In TIP3P water this is not a problem
> because it is rigid and the H sits inside the VDW radii of the oxygen.
> However, when you start moving to highly charged systems, a classic example
> being phosphates you can end up with the H collapsing into one of the other
> negatively charged atoms for which there is not a bond or angle term to
> prevent it and as soon as that distance becomes very small so you divide by
> zero and get NAN. The CPU is a little more tolerant of division by VERY
> small numbers than the GPU and that's why it manifests itself more as a blow
> up (trapped by vlimit) on the CPU and you get NAN on the GPU. Either way the
> simulation is bust because this event is non physical. 1fs can help, since
> the thing won't shoot off as much but it still won't really fix the problem.
>
> I would take a VERY careful look at the parameterization - that is the
> underlying problem. To begin with you probably should not use type HO for
> the hydrogen and instead create a new type. Then make sure this new H type
> has at least a small VDW radii, a very small radii will help (look in the
> frcmod file) but ideally one should probably attempt to refit the radii
> properly for such cases but this is tricky. In your case I would just put a
> small radii on the hydrogen (take a look at what the charmm people do with
> TIP3P - which in charmm should really strictly be called TIP3P-CHARMM since
> it isn't the original TIP3P parameterization due to the addition of small
> radii on the hydroxyl hydrogens.)
>
> What this does highlight though is that it is important that people still
> test things on the CPU and not just rely on the GPU code for everything. If
> you see a problem with a GPU run then switch to CPU and analyse things
> carefully. Don't expect the CPU code to also explode since it is more
> tolerant of marginal stability cases but at least be aware that you might be
> simulating something that is borderline stable at best.
>
> All the best
> Ross
>
>> -----Original Message-----
>> From: Fabrício Bracht [mailto:bracht.iq.ufrj.br]
>> Sent: Thursday, June 28, 2012 12:09 PM
>> To: AMBER Mailing List
>> Subject: Re: [AMBER] NaN with GTX580 pmemd.cuda simulation with and
>> without iwrap = 1
>>
>> Hello. Thank you for the reply. In fact, a time step of 1fs seems to
>> work around the problem. At the same time as i was writing the last
>> email, I had submitted a job with dt = 0.001 and it finished without
>> any problems. I haven´t tested yet with pmemd.cuda, but I´ll do that
>> later on.
>> The two atoms involved in this case are O and H of the hydroxyl
>> residue linked to the zinc atom of the active site. Where do I find
>> the Lennard-Jones rmin values for them? The hole active site was
>> parameterized following the steps on the MTK++ manual.
>> I´ll follow your instructions as to getting a restart file from a
>> previous step and test it out. I´ll inform you as soon as I have the
>> results.
>> Thank you
>> Fabrício Bracht
>>
>> 2012/6/28 David A Case <case.biomaps.rutgers.edu>:
>> > On Thu, Jun 28, 2012, Fabrício Bracht wrote:
>> >
>> >> Hi. I did what you asked. I ran the simulation again, in order to
>> make
>> >> sure that the error happened again. This time, it dies a lot earlier
>> >> with the same warnings as last time (vlimit exceed....etc). I
>> >> continued the simulation using as input the last restrt file written
>> >> using the serial version of pmemd. The simulation crashed after a
>> few
>> >> hours , the last output lines are:
>> >>
>> >>
>> >>
>> >> vlimit exceeded for step 112126; vmax = 91.1788
>> >> vlimit exceeded for step 112176; vmax = 45.6806
>> >>
>> >> Coordinate resetting cannot be accomplished,
>> >> deviation is too large
>> >> iter_cnt, my_bond_idx, i and j are
>> : 2 2690 5315 5316
>> >>
>> >
>> > Thanks for the info. It looks to me like the problem has nothing to
>> do (per
>> > se) with GPU's, since you can reproduce it on a serial CPU using
>> pmemd.
>> > Unfortuantely, problems like this (generally with the force field)
>> can be hard
>> > to track down. Your mdin file looks fine to me, although it is safer
>> to use
>> > dt=0.001 rather than 0.002, especially for constant pressure
>> simulations, and
>> > if you system is not *quite* equilibrated.
>> >
>> > The information above indicates that the problem occurs around atoms
>> 5315
>> > and 5316. You might see what those are, as this might (or might not)
>> give
>> > a clue to the origin of the problem. Look hard (and/or let us know
>> about)
>> > any non-standard residues or molecules that are in your
>> simulation. Do
>> > you have any Lennard-Jones rmin values that are either zero or small
>> (i.e.
>> > less than typical values for that type of atom)? You may have to get
>> a
>> > restart file someplace a little before step 112126; (don't use ig=-1
>> to do
>> > this: copy the actual seed that is printed in your output file, so
>> you will
>> > get an identical trajectory.) Then you can set ntwx=1, and examine
>> the
>> > coordinates at every step as the system crashes....one hopes that
>> this will
>> > identify the source of the problem.
>> >
>> > ...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
>
>
> _______________________________________________
> 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
Received on Thu Jun 28 2012 - 13:30:01 PDT
