[gmx-users] Significant differences in Ewald Vs PME for clay sheet
Gaurav Goel
gauravgoeluta at gmail.com
Mon Feb 2 14:43:08 CET 2015
Hi Mark,
We have done repeat calculations with version 5.0.4. To briefly give you
the main results--
1. Both versions, 4.5.1 and 5.0.4, with either Ewald and PME give the same
values for potential energy and interparticle forces.
2. Ewald (gromacs version 5.0.4), PME (5.0.4) and PME (4.5.1) give the
same result for the virials also. Only Ewald (4.5.1) is different. However.
only Ewald (4.5.1) matches the sheet size and Young's modulus as obtained
by force field developers.
3. Use of 3d or 3dc does not change the above determination.
Details:
Fourierspacing =0.08; rc=1.2 nm; ewald_rtol=1e-5; geometry=3d in all cases.
Increasing accuracy of PME does not affect results.
Case 1. Ewald / 4.5.1
Case 2: PME / 4.5.1
case 3: Ewald / 5.0.4
case 4; PME / 5.0.4
Ew/4.5.1 PME/4.5.1 Ew/5.0.4 PME/5.0.4
Total PE: -404119 -404118 -404112 -404113
Note 1: Individual contributions to the total PE in all 4 cases are
similarly matching.
Note 2: Force versus atom plots as obtained from .trr file is again
matching in the four cases.
However, the various virial contributions differ significantly between case
1 and rest:
Ew/4.5.1 PME/4.5.1 Ew/5.0.4 PME/5.0.4
Vir-xx: 7444.71 -43309.3 -43317.7 -43315.8
Vir-yy: -6327.84 -56862.3 -56869.5 -56869.5
Above numbers clearly indicate a large compressive pressure in cases 2-4
(PME), that is commensurate with observed increase in sheet size along x-
and y- directions in (semi-isotropic) NPT ensemble at zero pressure.
So in absence of any other data we can say something funky was happening in
virial calculation with Ewald/4.5.1 which was rectified in version 5.0.4.
However, Ewald/4.5.1 is the only case that matches with simulations of
Heinz et al. done using the molecular studio package. Any suggestions?
Thanks,
G
On Sun, Jan 25, 2015 at 6:25 PM, Mark Abraham <mark.j.abraham at gmail.com>
wrote:
> On Sat, Jan 24, 2015 at 6:48 AM, Gaurav Goel <gauravgoeluta at gmail.com>
> wrote:
>
> > Hi mark,
> >
> > We have used version 4.5.1.
> >
> > Ewald_geometry = 3d gives exactly same results as 3dc. E.g., viral_xx is
> > -43303.7 in former and -43303.8 in latter.
> >
>
> OK, thanks. In practice, probably only Berk has any clue what this code
> does and how, and if there's a problem then we'd need to know whether it
> still exists in the latest versions (e.g. 5.0.4). Can you try that on your
> inputs, please?
>
> Mark
>
>
> > Regards,
> > G
> > On 23-Jan-2015 10:28 pm, "Mark Abraham" <mark.j.abraham at gmail.com>
> wrote:
> >
> > > Hi,
> > >
> > > What do you see with ewald_geometry=3d? In what GROMACS version?
> > >
> > > Mark
> > >
> > > On Fri, Jan 23, 2015 at 5:47 PM, Gaurav Goel <gauravgoeluta at gmail.com>
> > > wrote:
> > >
> > > > A (brief) background:
> > > >
> > > > We have prepared a periodic Na-MMT clay sheet using force-field
> > > parameters
> > > > as developed by Heinz and co-workers. The box (and clay sheet) are
> > > 2.6-2.7
> > > > nm in x- and y- directions. Sheet thickness along z is 0.7nm and we
> > add a
> > > > vacuum layer to get a box_z=4nm.
> > > >
> > > > On using Ewald summation for electrostatic interactions, we get an
> > exact
> > > > correspondence with results of Heinz et al. The equilibrium box-size
> in
> > > > 2.6nm in x and y, as reported by Heinz using Ewald summations.
> Further
> > > the
> > > > Young's modulus values are in close agreement with other simulation
> > > studies
> > > > as well as experimental data.
> > > >
> > > > However when we use PME electrostatics, we see the equilibrium box
> size
> > > > increases to 2.7nm. Also, the Young's modulus values are doubled
> > compared
> > > > to Ewald summations. This prompted us to dig more into the energy and
> > > force
> > > > values.
> > > >
> > > > *** ***
> > > > As test case, we took a single snapshot of Ewald summation
> equilibrated
> > > > clay sheet configuration (box dimensions: 2.6 x 2.6 x 4 nm) and using
> > > > 'mdrun -rerun' determined the potential energy, forces and virial
> > > > components. Details are below:
> > > >
> > > > Case 1. Ewald (fourierspacing=0.12 ~ 40 wave vectors; ewald_rtol:
> > 1e-5;
> > > > rc=1.2nm)
> > > > Case 2: PME (fourierspacin=0.08; ewald_rtol: 1e-7; pme_order 8;
> > > ewald_3dc;
> > > > rc=1.2nm)
> > > > We have tried a few other combinations, such as standard Ewald-PME
> > > instead
> > > > of ewald_3dc, but all the numbers are same as case 2.
> > > > Case 1(Ewald) case 2 (PME)
> > > >
> > > > Total PE: -404119 -404118
> > > > Note 1: Individual contributions to the total PE in both cases are
> > > > similarly matching.
> > > > Note 2: Force versus atom plots as obtained from .trr file is again
> > > exactly
> > > > matching in the two cases.
> > > >
> > > > However, the various Virial contributions differ significantly:
> > > > Vir-xx: 7444.71 -43303.8
> > > > Vir-yy: -6327.84 -56858.7
> > > >
> > > > Above numbers clearly indicate a large compressive pressure in case 2
> > > > (PME), that is commensurate with observed increase in sheet size
> along
> > x-
> > > > and y- directions in (semi-isotropic) NPT ensemble at zero pressure.
> > > >
> > > > ** ** **
> > > > Since the virial components are determined from interparticle forces,
> > why
> > > > we are seeing this difference in terms of 'force versus atom' data
> > > matching
> > > > in Ewald versus PME, but the virials showing a significant
> difference.
> > > >
> > > > We will appreciate any help/suggestions on above issue. Thank you for
> > > your
> > > > time.
> > > >
> > > > Thanks,
> > > > G
> > > >
> > > > --
> > > > Gaurav Goel, PhD
> > > > Assistant Professor
> > > > Department of Chemical Engineering
> > > > Indian Institute of Technology, Delhi
> > > > Hauz Khas, New Delhi 110016
> > > > India
> > > > --
> > > > Gromacs Users mailing list
> > > >
> > > > * Please search the archive at
> > > > http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> > > > posting!
> > > >
> > > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> > > >
> > > > * For (un)subscribe requests visit
> > > > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users
> or
> > > > send a mail to gmx-users-request at gromacs.org.
> > > >
> > > --
> > > Gromacs Users mailing list
> > >
> > > * Please search the archive at
> > > http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> > > posting!
> > >
> > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> > >
> > > * For (un)subscribe requests visit
> > > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> > > send a mail to gmx-users-request at gromacs.org.
> > >
> > --
> > Gromacs Users mailing list
> >
> > * Please search the archive at
> > http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> > posting!
> >
> > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> >
> > * For (un)subscribe requests visit
> > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> > send a mail to gmx-users-request at gromacs.org.
> >
> --
> Gromacs Users mailing list
>
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.
>
--
Gaurav Goel, PhD
Assistant Professor
Department of Chemical Engineering
Indian Institute of Technology, Delhi
Hauz Khas, New Delhi 110016
India
More information about the gromacs.org_gmx-users
mailing list