Re: [AMBER-Developers] Zero VDW radii causing havoc.

From: Gerald Monard <Gerald.Monard.univ-lorraine.fr>
Date: Sun, 18 Aug 2013 10:22:38 +0200

Hi,

It seems to me also that this is more a symptom than the cause of the
problem. This happened to me on a 95k atom system (protein+ligand). At
some time during the MD, one water molecule bumped into the ligand, a
hydrogen from the ligand escaped, SHAKE went mad, and I enjoyed a core
dump. It happened with the GPU and the CPU (but not for the same frame...).
 From what I understand, SHAKE goes mad because some hydrogens move "too
much" between two frames. It cannot put back the hydrogens at the right
distance of their bonded atoms. On my system, I solved the problem by
turning ntf=2 to ntf=1: you add the harmonic potential on the hydrogen
bonds and restrain their movements.
I would suggest a simple workaround in two steps:
1) remove the ntf keyword from the doc and force ntf=1 for all
calculations (time saving is minimal with ntf=2)
2) double (x2 at least) the force constant on all bonds involving
hydrogens with ntc=2. This would prevent the hydrogens to move too much
between t and t+dt before running SHAKE.

Gerald.

On 08/17/2013 06:12 PM, Carlos Simmerling wrote:
> I've never seen this in 20 years, except in free energy calculations where
> the oxygen is being made smaller. In my opinion its probably a symptom, not
> the problem, unless the simulations are really different from anything I've
> ever done.
> On Aug 16, 2013 10:39 PM, "Ross Walker" <ross.rosswalker.co.uk> wrote:
>
>> So that suggests to me that this could indeed be a bug related to having a
>> small water box - has anyone seen this happen with the CPU code? A repro
>> there would be very telling.
>>
>> Also, this is NPT yes? As in it is not some NVT water box at some
>> artificially high compression?
>>
>> All the best
>> Ross
>>
>>
>>
>>
>> On 8/16/13 7:33 PM, "B. Lachele Foley" <lfoley.ccrc.uga.edu> wrote:
>>
>>> I can only add the reminder that our main workaround has been to increase
>>> the size of the water box. That is, it disproportionately seems to
>>> affect smaller systems. The locals who run simulations of really large
>>> systems never seem to see this. Sure, they can't run as long, but there
>>> are many more interactions in each step. I suspect most lipid
>>> simulations are also much larger than some of the systems we run. I
>>> can't promise that this information is useful, but there it is.
>>>
>>> FWIW, we're going to re-run some simulations that failed, but with our HO
>>> vdW terms on the HW, just to see what happens.
>>>
>>> We haven't seen it in CPU. I realize that doesn't mean it couldn't
>>> happen there.
>>>
>>> :-) Lachele
>>>
>>> Dr. B. Lachele Foley
>>> Complex Carbohydrate Research Center
>>> The University of Georgia
>>> Athens, GA USA
>>> lfoley.uga.edu
>>> http://glycam.ccrc.uga.edu
>>>
>>> ________________________________________
>>> From: Ross Walker [ross.rosswalker.co.uk]
>>> Sent: Friday, August 16, 2013 10:02 PM
>>> To: AMBER Developers Mailing List
>>> Subject: Re: [AMBER-Developers] Zero VDW radii causing havoc.
>>>
>>> Scott, I'm still not sure I am convinced on this being the issue. We have
>>> many many simulations done with TIP3P water - probably hundreds of
>>> milliseconds in aggregate, a lot with boosted potentials and other types
>>> of enhanced sampling and haven't seen this issue. Sure we've seen problems
>>> with zero VDW hydroxyl hydrogens before - but those are always when they
>>> are constrained to be close to each other or with massive attraction such
>>> as in some sugars and phosphates. I don't think I've ever seen problems
>>> with TIP3P before since the bonds are rigid and the H's are, I believe,
>>> inside the VDW sphere of the oxygen. Hence I would never expect to see
>>> this for two waters, for a water and a strained hydroxyl perhaps but even
>>> then unlikely.
>>>
>>> Are you certain this is a real problem and not a strange artifact, and or
>>> subtle bug? - Do we have any repro's for the CPU code?
>>>
>>> The only place I'd expect there could be problems with TIP3P is if it is
>>> run without shake. So something must be unique about the examples here
>>> that are showing the repro that doesn't normally occur in protein and/or
>>> dna simulations.
>>>
>>> Ps. we've never seen a problem with the lipid12 force field as far as I
>>> know and we have over 100 microsecs of simulation on these systems.
>>>
>>> Hence I am still a little suspicious...
>>>
>>> All the best
>>> Ross
>>>
>>>
>>>
>>> On 8/16/13 6:20 PM, "Scott Le Grand" <varelse2005.gmail.com> wrote:
>>>
>>>> Yep, that's why I moved the bug over to pmemd and away from pmemd.cuda...
>>>>
>>>>
>>>>
>>>>
>>>> On Fri, Aug 16, 2013 at 6:17 PM, B. Lachele Foley
>>>> <lfoley.ccrc.uga.edu>wrote:
>>>>
>>>>> We only had one, the HO (Ho) term. It's been made non-zero which at
>>>>> least
>>>>> fixed the issue in the highly charged, perfectly aligned situation
>>>>> where we
>>>>> first saw it. BTW, we (specifically Matt) first observed it on CPU.
>>>>>
>>>>> Can changing the LJ A term like that potentially mess up other things
>>>>> like
>>>>> hydrogen bonding?
>>>>>
>>>>> :-) Lachele
>>>>>
>>>>> Dr. B. Lachele Foley
>>>>> Complex Carbohydrate Research Center
>>>>> The University of Georgia
>>>>> Athens, GA USA
>>>>> lfoley.uga.edu
>>>>> http://glycam.ccrc.uga.edu
>>>>>
>>>>> ________________________________________
>>>>> From: Scott Le Grand [varelse2005.gmail.com]
>>>>> Sent: Friday, August 16, 2013 9:11 PM
>>>>> To: AMBER Developers Mailing List
>>>>> Subject: Re: [AMBER-Developers] Zero VDW radii causing havoc.
>>>>>
>>>>> There are a heck of a lot more TIP3P waters running around your test
>>>>> case.
>>>>> Anything with an effective zero VDW radius is a problem.
>>>>>
>>>>> 20 years ago, I used to see these when I did global conformational
>>>>> search
>>>>> of proteins. And when it happened, energies would plummet into a deep
>>>>> abyss and hose the search process. The only thing different with MD is
>>>>> that the finite timestep messes up a water badly enough that the fast
>>>>> Shake
>>>>> recovers it into a bizarre position that now hits a bunch of other
>>>>> atoms
>>>>> and then all you-know-what breaks loose and we're off to NaNdyland...
>>>>>
>>>>> The more things change...
>>>>>
>>>>>
>>>>> On Fri, Aug 16, 2013 at 6:06 PM, B. Lachele Foley <lfoley.ccrc.uga.edu
>>>>>> wrote:
>>>>>
>>>>>> Jason said: "The problem-case hydrogens in the relevant
>>>>> carbohydrates
>>>>>> (and lipid?) force fields should be modified if they represent
>>>>> problems
>>>>> in
>>>>>> normal simulation. "
>>>>>>
>>>>>> FYI: Scott Le Grand, in the bug report Ross mentioned, just said it
>>>>> was
>>>>>> two TIP3P waters doing it.
>>>>>>
>>>>>> :-) Lachele
>>>>>>
>>>>>> Dr. B. Lachele Foley
>>>>>> Complex Carbohydrate Research Center
>>>>>> The University of Georgia
>>>>>> Athens, GA USA
>>>>>> lfoley.uga.edu
>>>>>> http://glycam.ccrc.uga.edu
>>>>>>
>>>>>> ________________________________________
>>>>>> From: Jason Swails [jason.swails.gmail.com]
>>>>>> Sent: Friday, August 16, 2013 3:28 PM
>>>>>> To: AMBER Developers Mailing List
>>>>>> Subject: Re: [AMBER-Developers] Zero VDW radii causing havoc.
>>>>>>
>>>>>> On Aug 16, 2013, at 11:39 AM, David A Case <case.biomaps.rutgers.edu
>>>
>>>>>> wrote:
>>>>>>
>>>>>>> On Fri, Aug 16, 2013, Ross Walker wrote:
>>>>>>>
>>>>>>>> So I am proposing that we have code in rdparm that detects zero
>>>>> VDW
>>>>>> radii,
>>>>>>>> prints a warning message about it and then adjusts these radii in
>>>>> the
>>>>>> code
>>>>>>>> to 0.0000001.
>>>>>>>
>>>>>>> I don't understand how/why this fixes any problem. First, you
>>>>> would
>>>>> have
>>>>>>> to get to distances comparable to the radius for such a term to
>>>>> have
>>>>> any
>>>>>>> effect. You would have already very much messed things up by that
>>>>> point.
>>>>>>> Second, atoms like OH have a zero epsilon as well as a zero radius,
>>>>> so
>>>>>>> just changing the radius would have no effect.
>>>>>>
>>>>>> I think Ross meant changing the A and B coefficients to 1E-6 or so,
>>>>> not
>>>>>> the pre-combined values.
>>>>>>
>>>>>> I'm with Dave on this though. I've run many microseconds of net
>>>>> simulation
>>>>>> and I've never seen a problem with zero vdW terms in HO atom types. I
>>>>>> realize there are some cases where this bug hits you, but as I recall
>>>>> it
>>>>> is
>>>>>> a strange case with carbohydrates. I'm wary about just making these
>>>>> changes
>>>>>> under the assumption it will have only remedial effects. The
>>>>> problem-case
>>>>>> hydrogens in the relevant carbohydrates (and lipid?) force fields
>>>>> should
>>>>> be
>>>>>> modified if they represent problems in normal simulation.
>>>>>>
>>>>>> I can also quite easily write a quick parmed action to implement
>>>>> these
>>>>>> changes without having to hard-code this new convention immediately
>>>>> into
>>>>>> every md engine we have.
>>>>>>
>>>>>> Also, I think that the amber codes should implement the Amber force
>>>>> field
>>>>>> exactly as it is defined in the prmtop. Implicitly changing
>>>>> parameters in
>>>>>> the MD codes is a bad idea, IMO. (How hard would it have been to
>>>>> validate
>>>>>> chamber if CHARMM took that approach?) So if we decide to implement
>>>>> this
>>>>>> change, it should be done in LEaP, not rdparm (but again, I think the
>>>>> LJ
>>>>>> terms should be derivable straight from the combining rules).
>>>>>>
>>>>>>> Note that Istvan's lmodprmtop hack directly changes the LJ A term
>>>>> from
>>>>>>> zero to 1000. This has the advantage of not assuming any
>>>>> particular
>>>>>> mixing
>>>>>>> rules, i.e. it still works when NBFIX like changes have been used.
>>>>>>
>>>>>> As an aside, should lmodprmtop functionality be incorporated into
>>>>> parmed
>>>>>> and/or cpptraj if it does simple things? (I'm not sure what else it
>>>>> does...)
>>>>>>
>>>>>> All the best,
>>>>>> Jason
>>>>>>
>>>>>> _______________________________________________
>>>>>> AMBER-Developers mailing list
>>>>>> AMBER-Developers.ambermd.org
>>>>>> http://lists.ambermd.org/mailman/listinfo/amber-developers
>>>>>>
>>>>>>
>>>>>>
>>>>>> _______________________________________________
>>>>>> AMBER-Developers mailing list
>>>>>> AMBER-Developers.ambermd.org
>>>>>> http://lists.ambermd.org/mailman/listinfo/amber-developers
>>>>>>
>>>>> _______________________________________________
>>>>> AMBER-Developers mailing list
>>>>> AMBER-Developers.ambermd.org
>>>>> http://lists.ambermd.org/mailman/listinfo/amber-developers
>>>>>
>>>>>
>>>>>
>>>>> _______________________________________________
>>>>> AMBER-Developers mailing list
>>>>> AMBER-Developers.ambermd.org
>>>>> http://lists.ambermd.org/mailman/listinfo/amber-developers
>>>>>
>>>> _______________________________________________
>>>> AMBER-Developers mailing list
>>>> AMBER-Developers.ambermd.org
>>>> http://lists.ambermd.org/mailman/listinfo/amber-developers
>>>
>>>
>>>
>>> _______________________________________________
>>> AMBER-Developers mailing list
>>> AMBER-Developers.ambermd.org
>>> http://lists.ambermd.org/mailman/listinfo/amber-developers
>>>
>>>
>>>
>>> _______________________________________________
>>> AMBER-Developers mailing list
>>> AMBER-Developers.ambermd.org
>>> http://lists.ambermd.org/mailman/listinfo/amber-developers
>>
>>
>>
>> _______________________________________________
>> AMBER-Developers mailing list
>> AMBER-Developers.ambermd.org
>> http://lists.ambermd.org/mailman/listinfo/amber-developers
>>
> _______________________________________________
> AMBER-Developers mailing list
> AMBER-Developers.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber-developers
>

-- 
____________________________________________________________________________
  Prof. Gerald MONARD
  Theoretical Chemistry and Biochemistry Group
  SRSMC, Université de Lorraine, CNRS
  Boulevard des Aiguillettes B.P. 70239
  F-54506 Vandoeuvre-les-Nancy, FRANCE
  e-mail : Gerald.Monard.univ-lorraine.fr
  tel.   : +33 (0)383.684.381
  fax    : +33 (0)383.684.371
  web    : http://www.monard.info
____________________________________________________________________________
_______________________________________________
AMBER-Developers mailing list
AMBER-Developers.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber-developers
Received on Sun Aug 18 2013 - 01:30:03 PDT
Custom Search