Ah, okay, very clear now how this happens. The problem is indeed using
belly dynamics to freeze all the solvent. As I said earlier, for any atom
for which movement is not allowed, the velocities and forces are set to 0,
and you would not expect movement. BUT there is an additional wrinkle with
belly dynamics. That is, that forces and virials are not computed in direct
space for atom pairs that will have frc/vel set to 0 (so basically the
frozen belly atoms interactions are dropped from the pairlist). This screws
up the virial severely if all solvent is marked as nonmoving, so the
pressure always looks low (-). So the pressure adjustment code scales the
coords - and indeed does move everything closer together as the box shrinks.
Well, Dave Case has been advocating using ntr restraints instead of belly
dynamics for a long time... (not that this in any way is a legitimate use of
belly in the first place, but dropping the virial calc in direct space also
has a somewhat unexpected impact). Thanks for the additional info.
- Bob Duke
----- Original Message -----
From: "Bill K" <bill.mercury.chem.pitt.edu>
To: "AMBER Mailing List" <amber.ambermd.org>
Sent: Tuesday, February 10, 2009 2:39 PM
Subject: Re: [AMBER] Stability of a long NPT md trajectory?
Hi,
This is beginning and endpoint of the same trajectory,
which seems very curious to me. I can attach the entire output
(and input files), and you can see the density increasing and volume
decreasing throughout the simulation. I'd think that the ewald
error is a reflection of the use of the belly option for almost
the entire simulation box.
-Bill
> Okay, I am just curious about how the heck this all happened. I am
> presuming that you did some steps between this timestep and the ~350 psec
> timestep under conditions different than the one mdin you have shown us
> (with belly on all the solvent). This intervening period is what is
> interesting - how the heck did the solvent contract to somewhere between a
> quarter and a third of normal volume? - very aphysical. What I have
> objected to is that the whole deal is due to belly - belly should have
> fixed
> the solvent atoms at the current volume, and you would only have the res 1
> atom rattling around. I note that for all the steps you have shown, the
> ewald error estimate is very bad - this in itself indicates something
> anomalous about the setup.
> Regards - Bob
> ----- Original Message -----
> From: "Bill K" <bill.mercury.chem.pitt.edu>
> To: "AMBER Mailing List" <amber.ambermd.org>
> Sent: Tuesday, February 10, 2009 2:37 AM
> Subject: Re: [AMBER] Stability of a long NPT md trajectory?
>
>
> To shed some light on your discussion here, the solute
> was just a single atom, and oddly enough this NPT md
> run started with a typical TIP3P density as you get using
> solvatebox, as I'll show with the initial stats from the output:
>
> NSTEP = 0 TIME(PS) = 40.000 TEMP(K) = 217.61 PRESS = -190.3
> Etot = -72887.4935 EKtot = 0.6487 EPtot = -72888.1421
> BOND = 0.0000 ANGLE = 0.0000 DIHED = 0.0000
> 1-4 NB = 0.0000 1-4 EEL = 0.0000 VDWAALS = -54.6700
> EELEC = -72833.4721 EHBOND = 0.0000 RESTRAINT = 0.0000
> EKCMT = 0.6487 VIRIAL = 189.3043 VOLUME = 45924.9288
> Density = 0.7024
> Ewald error estimate: 0.9839E+00
>
> and terminated, as below, using the input parameters that I provided.
> Certainly, freezing all of the solute was unintentional.
>
>
>> Yes, I think we are in violent agreement here. Something sort of
>> strange,
>> but I suspect the setup as much as how belly is used. The unknown is
>> the
>> size of the solute and whether it is more than 1 residue. But here's
>> the
>> kicker. Only residue 1 moves (makes sense from the low kinetic energy;
>> don't know what I think about no bonded energy terms at all - must be a
>> rather simple solute). But the belly has the solvent completely frozen
>> from
>> the start. So assume the solute is small. Then look at the per atom
>> volume - comes out to 2.77 angstrom**3 per atom. That is consistent
>> with
>> the high density, and way low for water which would more typically be
>> around
>> 10 angstrom**3 per atom at density 1.0. Now, given the belly, the
>> system
>> solvent STARTED out this way... So I would guess a fairly bogus initial
>> setup held that way with a NVT box as well as the belly, and when you
>> switch
>> to NPT the belly still holds you there. How/why the water got stacked
>> this
>> tight in the initial config, I haven't a clue, but I don't see how it
>> could
>> be from solvating with the standard water boxes... The belly
>> effectively
>> takes all solvent-solvent interactions out of the direct force calcs,
>> which
>> is why the pressure/virial is completely bogus (so in addition to the
>> belly
>> zeroing any calc'd forces on the solvent, it also skips their direct
>> force
>> virial contribution).
>> Regards - Bob
>>
>> ----- Original Message -----
>> From: "Carlos Simmerling" <carlos.simmerling.gmail.com>
>> To: "AMBER Mailing List" <amber.ambermd.org>
>> Sent: Monday, February 09, 2009 10:10 PM
>> Subject: Re: [AMBER] Stability of a long NPT md trajectory?
>>
>>
>>> Bob,
>>> I agree that random seed issues can be problematic, but I doubt it
>>> could
>>> ever be anything of the scale seen here. his density is 3 times
>>> normal,
>>> but
>>> for some reason the energy hasn't exploded and the pressure is too low!
>>> This
>>> isn't just correlation issues, or non-ergodicity or something like
>>> that,
>>> it
>>> is simply completely wrong at the most basic level. i bet it's the
>>> belly
>>> option.
>>> carlos
>>>
>>>
>>> On Mon, Feb 9, 2009 at 9:21 PM, Robert Duke <rduke.email.unc.edu>
>>> wrote:
>>>
>>>> I think this may be a really small unit cell (1000 waters - perhaps
>>>> too
>>>> small) with significant periodicity issues that then magnify the
>>>> problems
>>>> associated with usage of the same random seed through multiple
>>>> restarts
>>>> with
>>>> ntt 3, but I am also somewhat skeptical that this is the only problem.
>>>> For
>>>> one thing, I did not get any indication of how many restarts were done
>>>> with
>>>> the same seed - may have been none. I mentioned it purely because it
>>>> is
>>>> something that folks tend to overlook, and it can cause really bad
>>>> problems.
>>>> Regards - Bob
>>>>
>>>> ----- Original Message ----- From: "Carlos Simmerling" <
>>>> carlos.simmerling.gmail.com>
>>>> To: <bill.mercury.chem.pitt.edu>; "AMBER Mailing List"
>>>> <amber.ambermd.org>
>>>> Sent: Monday, February 09, 2009 7:56 PM
>>>>
>>>> Subject: Re: [AMBER] Stability of a long NPT md trajectory?
>>>>
>>>>
>>>> I disagree- I run frequently with ntt=3 and no ig value and have
>>>> never
>>>>> seen the problem you reported.
>>>>>
>>>>> are you sure that just changing ig and nothing else really solved it?
>>>>>
>>>>> beware of problems that go away without explanation- they often come
>>>>> back without explanation too.
>>>>>
>>>>> as I said in my last email, if you use PME, make sure you understand
>>>>> the limits of ibelly as compared to positional restraints, which
>>>>> think
>>>>> are far safer. problems with ibelly and PME, especially during
>>>>> equilibration, are well documented in the archives.
>>>>>
>>>>> On Mon, Feb 9, 2009 at 3:56 PM, Bill K <bill.mercury.chem.pitt.edu>
>>>>> wrote:
>>>>>
>>>>>> Gentlemen,
>>>>>> Thank you for the kind replies to my inquiry. I found your
>>>>>> messages to be confusing initially until I realized my egregious
>>>>>> mistake in applying the ibelly flag. I should have instead put in
>>>>>> my input file (to freeze residue 1):
>>>>>>
>>>>>> Belly Residues (dynamic residues)
>>>>>> RES 2 xxxx
>>>>>> END
>>>>>> END
>>>>>> For edification purposes, though, I found that using ntt=3
>>>>>> and not setting an ig value is what allowed the density to go
>>>>>> such a high value.
>>>>>> Thank you again for your time.
>>>>>> -Bill
>>>>>>
>>>>>> There may be a number of things going on here relating to
>>>>>> nonstandard
>>>>>>> approach to md, but one thing I noticed immediately is the use of
>>>>>>> ntt
>>>>>>> =
>>>>>>> 3.
>>>>>>> This in itself is not a problem, but one must be certain to not use
>>>>>>> the
>>>>>>> same
>>>>>>> random seed at each restart - you must set a new value for "ig",
>>>>>>> and
>>>>>>> I
>>>>>>> didn't notice it in the sample mdin. I don't understand at all
>>>>>>> why
>>>>>>> your
>>>>>>> density would be so high though - there is much more abnormal about
>>>>>>> this
>>>>>>> run
>>>>>>> than reuse of the same random seed in a restart for langevin
>>>>>>> dynamics.
>>>>>>> Regards - Bob Duke
>>>>>>>
>>>>>>> ----- Original Message -----
>>>>>>> From: "Carlos Simmerling" <carlos.simmerling.gmail.com>
>>>>>>> To: <bill.mercury.chem.pitt.edu>; "AMBER Mailing List" <
>>>>>>> amber.ambermd.org>
>>>>>>> Sent: Monday, February 09, 2009 7:08 AM
>>>>>>> Subject: Re: [AMBER] Stability of a long NPT md trajectory?
>>>>>>>
>>>>>>>
>>>>>>> this should not happen.
>>>>>>>> however, I suspect something related to your use of the belly
>>>>>>>> option
>>>>>>>> for nearly everything in the system. sander should have printed a
>>>>>>>> warning about use of ibelly with pme- take a look at that. you may
>>>>>>>> want to use resrtaints instead, but it's hard to give advice when
>>>>>>>> we
>>>>>>>> don't know why you defined only 1 dynamic residue.
>>>>>>>>
>>>>>>>>
>>>>>>>> On Mon, Feb 9, 2009 at 2:07 AM, Bill K
>>>>>>>> <bill.mercury.chem.pitt.edu>
>>>>>>>> wrote:
>>>>>>>>
>>>>>>>>> Hi Everyone,
>>>>>>>>> I'm running a relatively long (~10ns) trajectory of a single
>>>>>>>>> solute in a water box of approx. 1000 tip3p waters, with standard
>>>>>>>>> conditions of 300K and 1Atm. I apologize if this is an
>>>>>>>>> elementary
>>>>>>>>> question, but I've been surfing the mail reflector for some time
>>>>>>>>> now... I've brought the system to temperature
>>>>>>>>> with NVT runs of at least 20ps (I've tried up to 1ns), and then
>>>>>>>>> switched to NPT (using the input at the end of the message), with
>>>>>>>>> ntt = 3. When I do this, I notice that, throughout the run, the
>>>>>>>>> volume continues to decrease, and the pressure remains negative,
>>>>>>>>> even after the system density surpasses 1 by a large margin.
>>>>>>>>> the final output before the calculation crashes is:
>>>>>>>>>
>>>>>>>>> NMR restraints: Bond = 0.000 Angle = 0.000 Torsion =
>>>>>>>>> 0.000
>>>>>>>>>
>>>>>>>>> NSTEP = 313500 TIME(PS) = 353.500 TEMP(K) = 569.51 PRESS
>>>>>>>>> = -3544.4
>>>>>>>>> Etot = -71984.1552 EKtot = 1.6976 EPtot
>>>>>>>>> = -71985.8528
>>>>>>>>> BOND = 0.0000 ANGLE = 0.0000 DIHED =
>>>>>>>>> 0.0000
>>>>>>>>> 1-4 NB = 0.0000 1-4 EEL = 0.0000 VDWAALS
>>>>>>>>> -305.4604
>>>>>>>>> EELEC = -71680.3923 EHBOND = 0.0000 RESTRAINT =
>>>>>>>>> 0.0000
>>>>>>>>> EKCMT = 0.2890 VIRIAL = 636.5491 VOLUME =
>>>>>>>>> 8314.1509
>>>>>>>>> Density =
>>>>>>>>> 3.8796
>>>>>>>>> Ewald error estimate: 0.9533E+00
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> So, instead of trying to do production runs using NTP, should I
>>>>>>>>> stop
>>>>>>>>> the NTP run when my density is in the range of (1 or just over),
>>>>>>>>> and
>>>>>>>>> then
>>>>>>>>> carry on production runs at either NVT, or NVE with ntb=1 and
>>>>>>>>> ntt=0),
>>>>>>>>> or
>>>>>>>>> should this run be self regulating and stable and I've just
>>>>>>>>> completely
>>>>>>>>> missed something? I appreciate your time and consideration.
>>>>>>>>> thank you very much
>>>>>>>>> -Bill
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> # 1ns Dynamic Simulation with Constant Pressure
>>>>>>>>>
>>>>>>>>> # Control section
>>>>>>>>> &cntrl
>>>>>>>>> ntwx = 500, ntpr = 500,
>>>>>>>>> ntt = 3, gamma_ln = 1.0, temp0 = 300.0, tempi = 300.0, tautp =
>>>>>>>>> 1.0,
>>>>>>>>> scnb = 2.0, scee = 1.2, dielc = 1, cut = 8.0,
>>>>>>>>> ntb = 2, ntc = 2, ntf = 2,
>>>>>>>>> nstlim = 999999, dt = 0.0010,
>>>>>>>>> ntp = 1, taup = 5.0, pres0 = 1.0
>>>>>>>>> ibelly = 1, ntr = 0,
>>>>>>>>> imin = 0, irest = 0, ntx = 5, nmropt = 1,
>>>>>>>>> /
>>>>>>>>> &wt
>>>>>>>>> type = 'TEMP0', istep1 = 1, istep2 = 999999, value1 = 300.0,
>>>>>>>>> value2
>>>>>>>>> =
>>>>>>>>> 300.0,
>>>>>>>>> /
>>>>>>>>> &wt
>>>>>>>>> type='END'
>>>>>>>>> /
>>>>>>>>> Belly Residues (dynamic residues)
>>>>>>>>> RES 1
>>>>>>>>> END
>>>>>>>>> END
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> _______________________________________________
>>>>>>>>> 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
>>>>>>
>>>>>>
>>>>> _______________________________________________
>>>>> 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
>>
>>
>
>
>
> _______________________________________________
> 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 Wed Feb 11 2009 - 01:22:02 PST
