Hi Francois,
I think Jason's reply to this is the direction you need to look. I saw your email and was mostly trying to explain why you could misunderstand and think pmemd ran a 0 cutoff for you when sander wouldn't. The further story here - the code initializes the cut to 0.0 to check and see if you set a value for cut in the namelist, and assumes that after you read the namelist, if the cut value is still 0.0, then the defaults should be set for pme electrostatics cut and vdw cut instead (and it is even more complicated really, in that there are other variables specific for electrostatics vs. vdw cutoffs you could have specified). We could use some weird negative value for cut I suppose, but ultimately someone can surprise you with input, and the only way to be sure that they provided input for a value, at least as far as I know, is to read the namelist twice against different initial values for a variable. Anyway, I defer to Jason regarding the solution to your actual problem here - my intent was merely to try to
be sure you understood what pmemd had done. By the way, there should be clear indications in the mdout output as to the cutoff actually used - in the "ewald parameters" subsection (and of course in an ideal world, you would get the same result using an infinitesimal cut near 0 for the direct space calc because the reciprocal space calc would compensate - but in the real world, if you start messing with the actual pme run parameters, the accuracy of your results varies; messing with the parameters is not turning off pme, so you are always attempting a full pbc system electrostatics evaluation).
Regards - Bob
________________________________________
From: FyD [fyd.q4md-forcefieldtools.org]
Sent: Monday, February 24, 2014 1:31 AM
To: AMBER Mailing List
Subject: Re: [AMBER] PME with cutoff = 0
Dear Robert,
> The cutoff value of "0.d0" in pmemd for a pme run is equivalent to
> requesting that the default values be used - probably 8.0 angstrom,
> if memory serves, assuming igb 0 for a pme simulation.
Yes - this corresponds to what I get...
I compare MD simulations (with PME) with cut = 12 vs cut = 0
the obtained results are very similar.
This means 'cut = O' does not do what I want; I wonder if crashing is
not better than generating data in this case...
> PME does strange things at small-ish cutoff, but does not allow 0
> (in the sense that it says "oh, you wanted 8"), and I presume later
> catches and disallows cut < 0. In my experience, reducing the
> direct space cut below roughly 7.0 angstrom causes unacceptable
> error limits unless you increase grid density for reciprocal space
> (there is basically a balancing act between reciprocal and direct
> space in terms of what is necessary to keep error acceptable). And
> if you think long cutoffs are expensive, play around a bit with
> greatly increasing nfft1,2,3 in a large system :-}. So all I am
> really answering here is why pmemd ran - it basically ignored you...
ok
> I would have to read the paper done with amber 4.1 to have a clue
> what they were really up to,
See
http://pubs.acs.org/doi/full/10.1021/jp9717655 - It is written:
"The heat of vaporization is computed from the average intermolecular
interaction energy Eint via dHvap = -Eint + RT
For molecules that do not have an internal nonbonded interaction
beyond 1?4 (e.g., the chloromethanes), Eint is calculated
straightforwardly as the sum of EELEC and ENONB in the SANDER output,
divided by the number of molecules in the system. For all other
molecules, EELEC and ENONB also contain contributions from
intramolecular nonbonded interactions. To correct for these, we
performed a short simulation with the nonbonded cutoff radius set to
zero. Due to the residue-based cutoff in AMBER, EELEC and ENONB now
accumulate only the intraresidue nonbonded interactions that can be
subtracted from the total to yield the intermolecular portion needed
for the calculation of dHvap."
> but pmemd is designed to always use either
> pme or generalized Born, so resurrecting an old copy of sander may
> be your best bet..
I use PME & would like to continue to use PME.
Resurrecting an old copy of sander? uff...
Does it means this is not possible to carefully study/check solvent
boxes with Amber nowadays?
thank you,
regards, Francois
Dear Ross,
> Both sander and pmemd almost certainly make assumptions about the cutoff.
> Since the days of AMBER 4.1 there have been huge changes in the way the
> pairlist is built etc.
I can understand that ;-)
> For example the code use to just rebuild the
> pairlist at set intervals. These days it is smart and only rebuilds it
> when an atom has moved more than the distance skinnb such that the
> pairlist cutoff is really cut+skinnb. Thus reducing cut to zero probably
> messes up logic in this part of the code. There are likely other issues as
> well, like divisions by zero inside the code that determines the Ewald
> coefficient etc.
My understanding is that division by zero should always be checked.
> The fact one code crashes and the other doesn't is likely
> just dependent on where they hit the first problem. From what you show it
> looks like pmemd just hangs while sander crashes with an error - both
> could be considered a crash.
pmemd does not hang; see above.
> Anyway, to cut a long story short neither code officially supports a value
> of cut=0 for PME runs and it is something that isn't tested. It may be
> possible to run the simulation in GB in sander by setting igb=6 - this
> will run a gas phase simulation. I doubt there is anyway to run a periodic
> simulation with a cutoff of 0. If you really want to then your options are
> likely see if you can find a PRE-PME copy of sander (such as AMBER 4.1) or
> (sander-classic from AMBER 6 perhaps?) or crack open the debugger and see
> if you can trace through the use of the cutoff value and find out where
> the problems are occurring. This may just be a case of going down the
> rabbit hole though.
sander-classic from AMBER 6: do you think this is possible to get
sander-classic from AMBER 6 compiled using 4.4 GNU compiler?
Does it means this is not possible to carefully study/check solvent
boxes with Amber nowadays?
See
http://pubs.acs.org/doi/full/10.1021/jp9717655
> Sorry I can't offer much more help than that. - Note you could also try
> cut=0.001 and see if that works. There is almost certainly a 1.0d0/cut^2
> somewhere in the code and this might work around that and still give you
> what you want from the simulation.
I am testing...
thank you,
regards, Francois
> On 2/23/14, 9:19 AM, "FyD" <fyd.q4md-forcefieldtools.org> wrote:
>
>> Dear All,
>>
>> I ran amber10/sander.MPI with PME + cutoff = 0 -> sander crashes
>> amber12/sander.MPI with PME + cutoff = 0 -> sander crashes
>> (compiled with intel13 + mkl)
>> amber12/pmemd.MPI with PME + cutoff = 0 -> pmemd does not crash...
>>
>> Printed output: printing stops in section 4:
>>
>> --------------------------------------------------------------------------
>> 4. RESULTS
>> --------------------------------------------------------------------------
>>
>> | # of SOLUTE degrees of freedom (RNDFP): 24936.
>> | # of SOLVENT degrees of freedom (RNDFS): 0.
>> | NDFMIN = 24933. NUM_NOSHAKE = 0 CORRECTED RNDFP = 24933.
>> | TOTAL # of degrees of freedom (RNDF) = 24933.
>>
>> -> nothing else is printed...
>>
>> In Fox & Kollman J.Phys.Chem.B 1998, 102, 8070, MD simulation with
>> cutoff = 0 is reported & Sander from Amber 4.1 was used.
>>
>> Why does sander crash, while pmemd does not crash in this case?
>> (with the same input & cut = 12, sander does not crash)
>>
>> thank you, regards,
>> Francois
>>
>>
>>
>> # Cst pressure simulation with cutoff = 0 Angs; 300 K short
>> $PARAL -np $NPROC $EXE -O -i md-cstpres0.in \
>> -o md-cstpres0b.out \
>> -p $mol.parm7 \
>> -c md-cstvol.rst7 \
>> -r md-cstpres0b.rst7 \
>> -x md-cstpres0b.mdcrd
>>
>> my test input below:
>>
>> cst pressure MD, nstlim*dt psec, PME with cutoff
>> &cntrl
>> nmropt = 0, ioutfm = 1, iwrap = 0,
>> ntx = 1, irest = 0, ntrx = 1, ntxo = 1,
>> ntpr = 100, ntwr = 1000, ntwx = 500,
>> ntave = 0, ntwv = 0, ntwe = 0,
>> ntf = 2, ntb = 2,
>> ntc = 2, tol = 0.00001,
>> dielc = 1, ipol = 0,
>> cut = 0, ig = 71277, comp = 52.5,
>> ibelly = 0, ntr = 0, igb = 0,
>> imin = 0, nstlim = 5000,
>> nscm = 1000, t = 0.0, dt = 0.002,
>> temp0 = 300, tempi = 300,
>> ntt = 1, tautp = 0.2, vlimit = 20,
>> ntp = 1, pres0 = 1.0, taup = 0.2,
>> &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
Received on Mon Feb 24 2014 - 10:00:03 PST
