I did check the forces (as well as the relative EEL energy for different snapshots), but they don't match.
Here is the forces (first line) for cut=8:
-7.515 -1.297 -6.394 20.476 -1.337 4.897 -5.195 -2.233 3.519 -6.623
And for cut=10:
-7.619 -0.880 -6.132 21.938 -2.793 2.594 -4.984 -2.038 4.672 -6.990
Again, I don't think the difference is neglectable.
The EEL energy for the first two snapshot under cut=8 are: -8608.5697 and -8703.3207
those under cut=10 are: -8714.0624 and -8907.5023.
Besides, I did the calculation for the same system with the same force field and water model in OpenMM, and the results for
nonbondedCutoff=8*angstrom generally agree well with nonbondedCutoff=10*angstrom.
Here is the 8 A result:
MM energy is:
-7455.212402905567
MM force is:
[[-8.87632470e+00 -2.60438588e-03 -6.16627323e+00]
[ 2.23216476e+01 -1.52421599e+00 2.67599205e-01]
[-5.95707899e+00 -2.09336558e+00 4.51639169e+00]
...
[-2.25740521e+01 -1.54278686e+00 -7.58745082e+00]
[ 9.38813232e+00 1.24806778e+00 2.62724138e+00]
[ 1.60604423e+01 -5.04632512e-01 6.92133250e+00]]
And here is the 10 A result:
MM energy is:
-7456.282483798095
MM force is:
[[-8.87784094e+00 -1.70950930e-03 -6.16671149e+00]
[ 2.23083739e+01 -1.53055808e+00 2.54194937e-01]
[-5.95628932e+00 -2.09862229e+00 4.51622240e+00]
...
[-2.25861562e+01 -1.54085941e+00 -7.58845137e+00]
[ 9.40012347e+00 1.26370078e+00 2.62059393e+00]
[ 1.60632034e+01 -4.96608480e-01 6.92517787e+00]]
It's also weird that they don't match with any of the amber results (probably because OpenMM also applies PME for LJ by default while Amber does not).
If I turned off the PBC, they could perfectly match.
OpenMM looks a lot more than I would have expected Amber to do as well, with more or less roundoff differences.It looks like either something is wrong, or at least different in what the cutoffs mean in Amber's PME implementation.
Best regards,
Ge
________________________________
From: Steinbrecher, Thomas <thomas.steinbrecher.roche.com>
Sent: Tuesday, August 20, 2024 11:05 AM
To: Ge Song <g.song.duke.edu>; AMBER Mailing List <amber.ambermd.org>
Subject: Re: [Sender Not Verified] [AMBER] Cutoff value significantly impacts the electrostatic energy in PME
Hi Ge,
Im not a PME expert, but maybe what you see is more like a shifting of the (arbitrary) EEL zero point with a different cutoff. I think it would be more meaningful to compare the forces you get with the different cutoffs, which is what would lead to different dynamics and hence different results of your MD.
Kind Regards,
Thomas
On Tue, Aug 20, 2024 at 3:44 PM Ge Song via AMBER <amber.ambermd.org<mailto:amber.ambermd.org>> wrote:
**Warning** The sender address (Ge Song via AMBER ) can not be verified, sender email address could be spoofed. Please take care to proceed.
Dear community,
In sander, it performs PME for calculation of electrostatic interaction in periodic system, and the 'cut' keyword is used only for the cutoff of direct space summation. Therefore, in principle, setting it as 8 Angstrom or 10 Angstrom won't significantly change the EEL energy in the output, this is also in a previous discussion at here RE: Re: [AMBER] cutoff from Ross Walker on 2010-06-05 (Amber Archive Jun 2010) (ambermd.org<https://urldefense.com/v3/__http://ambermd.org__;!!OToaGQ!scrgBH-NJzRrvrn4sh-YXYux6I-ZF9uIQb72MzSYwJeB5xcRT8BcBe4crFca46KIgTg3rqU3b0917uR1Wd16qo6HGw$>)<http://archive.ambermd.org/201006/0138.html<https://urldefense.com/v3/__http://archive.ambermd.org/201006/0138.html__;!!OToaGQ!scrgBH-NJzRrvrn4sh-YXYux6I-ZF9uIQb72MzSYwJeB5xcRT8BcBe4crFca46KIgTg3rqU3b0917uR1Wd3zgjQbng$>>. On the other hand, I found it does for my alanine-in-12 A water box system when I rerun the trajectory.
Here is my sander input:
# md_rerun.in<
https://urldefense.com/v3/__http://md_rerun.in__;!!OToaGQ!scrgBH-NJzRrvrn4sh-YXYux6I-ZF9uIQb72MzSYwJeB5xcRT8BcBe4crFca46KIgTg3rqU3b0917uR1Wd0MFIiQbw$>
&cntrl
imin=5,
ntx=5,
maxcyc=1,
irest=0,
temp0=300.0,
ioutfm=0,
ntwf=1,
ntwx=1,
cut=10.0,
ntb=1,
/
and in the output for cut=10:
NSTEP ENERGY RMS GMAX NAME NUMBER
1 -7.5647E+03 1.4040E+01 1.0199E+02 O 2168
BOND = 1.6571 ANGLE = 8.7221 DIHED = 12.0941
VDWAALS = 1079.3543 EEL = -8714.0624 HBOND = 0.0000
1-4 VDW = 3.0208 1-4 EEL = 44.5127 RESTRAINT = 0.0000
minimization completed, ENE= -.75647012E+04 RMS= 0.140397E+02
changing the cut=10 to cut=8 gives:
NSTEP ENERGY RMS GMAX NAME NUMBER
1 -7.4416E+03 1.3970E+01 9.7058E+01 O 2168
BOND = 1.6571 ANGLE = 8.7221 DIHED = 12.0941
VDWAALS = 1097.0112 EEL = -8608.5697 HBOND = 0.0000
1-4 VDW = 3.0208 1-4 EEL = 44.5127 RESTRAINT = 0.0000
minimization completed, ENE= -.74415517E+04 RMS= 0.139702E+02
As you can see, the EEL part changes from -8714.0624 to -8608.5697 kcal/mol, which is definitely not neglectable.
I wonder if any experts know what could be the probable reason.
Best regards,
Ge
_______________________________________________
AMBER mailing list
AMBER.ambermd.org<mailto:
AMBER.ambermd.org>
http://lists.ambermd.org/mailman/listinfo/amber<
https://urldefense.com/v3/__http://lists.ambermd.org/mailman/listinfo/amber__;!!OToaGQ!scrgBH-NJzRrvrn4sh-YXYux6I-ZF9uIQb72MzSYwJeB5xcRT8BcBe4crFca46KIgTg3rqU3b0917uR1Wd2d3yjmOA$>
--
Dr. Thomas Steinbrecher
Principal Scientist CADD
Roche Pharma Research and Early Development
Roche Innovation Center Basel
F. Hoffmann-La Roche Ltd
Bldg. 06/OG07
Grenzacherstrasse 124
4070 Basel
Switzerland
Phone +41 61 682 1319
mailto: thomas.steinbrecher.roche.com<mailto:thomas.steinbrecher.roche.com>
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Tue Aug 20 2024 - 11:30:02 PDT