Re: [gmx-users] Inconsistence of SPC/E and TIP4P water in Gromacs

[Justin Lemkul]

On 5/11/17 1:53 AM, Mark Abraham wrote:

Hi,

On Thu, May 11, 2017 at 3:31 AM ZUO Taisen  wrote:



Hi guys:

 I have compared the SPC/E and TIP4P water of the opls force field in
Gromacs.But there is something inconsistent of the potential energy in the
system


The potential energy of SPC/E water is -4.71400e+04/1000=-47.14kJ/mol
which is far from the results from literature (Peter G. Kusalik,
Science,vol,265,26,p1219-1221,1994)
The potential energy TIP4P water is -4.17875e+04/1000=41.787kJ/mol which
is very close to the literature -41.8kJ/mol(Peter G. Kusalik,
Science,vol,265,26,p1219-1221,1994)


I have attached all the files simulating the SPC/E and TIP4P water(.trr
.xtc files are delected) and also the related paper in my last email but
was suspended because of the files were too big.So no files was attached
this time but the .mdp file was posted as following.



We shouldn't have a mailing list where thousands of recipients can receive
large emails :-)



I'm eager to know why! Thank you very much!



Force fields, including water models, are parameterized to reproduce
certain theoretical or experimental data. To the extent that those targets
resemble other observables, you might expect agreement. You should expect
the level of such agreement to differ between models.



Especially when the methods applied are different.  The original Berendsen 1987 
SPC/E paper used very different methods relative to what the OP shows in the 
.mdp file.  In the original paper, the cutoffs were 0.9 nm, no PME, Berendsen's 
weak coupling method for temperature and pressure, temperature at 300 K, SHAKE 
instead of SETTLE (doubt that will make much of a difference, honestly), and the 
simulations were only a whopping 27.5 ps.  So there's no doubt in my mind that 
one will get different results.  Also, every paper that involves SPC/E in some 
form of comparison reports a different diffusion coefficient (though this is not 
unique to SPC/E), and most are not corrected for system size effects. 
Reproducibility is a considerable challenge...


-Justin


Mark



TIP4P simulation results:
Statistics over 3001 steps using 31 frames

   Energies (kJ/mol)
LJ (SR)  Disper. corr.   Coulomb (SR)   Coul. recip.  Potential
  7.71821e+03   -2.05973e+02   -4.93956e+049.58136e+01
 -4.17875e+04
Kinetic En.   Total EnergyTemperature Pres. DC (bar) Pressure (bar)
7.30494e+03   -3.44826e+042.93006e+02   -1.14009e+021.41117e+00


SPC/E simulation results:

Statistics over 3001 steps using 31 frames


   Energies (kJ/mol)
LJ (SR)  Disper. corr.   Coulomb (SR)   Coul. recip.  Potential
9.23352e+03   -2.12206e+02   -5.62576e+049.63633e+01   -4.71400e+04
Kinetic En.   Total EnergyTemperature Pres. DC (bar) Pressure (bar)
7.30494e+03   -3.98350e+042.93005e+02   -1.18040e+021.64654e+00




title= OPLS MD waterethanol siimulation
; Run parameters
integrator= md; leap-frog integrator
nsteps= 3000; 1 * 2000 = 2 ps (20 ns)
dt= 0.001; 1 fs
; Output control
nstxout= 10; save coordinates every 10.0 ps
nstvout= 10; save velocities every 10.0 ps
nstenergy= 10; save energies every 10.0 ps
nstlog= 10; update log file every 10.0 ps
nstxout-compressed  = 10  ; save compressed coordinates every 10.0
ps
define  =-DEFLEXIBLE ; flexible
water
compressed-x-grps   = System; replaces xtc-grps
; Bond parameters
 continuation= yes; Restarting after NPT
 constraint_algorithm= lincs; holonomic constraints
 constraints= none; H bonds (even heavy atom-H bonds)
constrained
 lincs_iter= 1; accuracy of LINCS
 lincs_order= 4; also related to accuracy
; Neighborsearching
cutoff-scheme   = Verlet
ns_type= grid; search neighboring grid cells
nstlist= 10; 20 fs, largely irrelevant with Verlet scheme
rcoulomb= 1.2; short-range electrostatic cutoff (in nm)
rvdw= 1.2; short-range van der Waals cutoff (in nm)
; Electrostatics
coulombtype= PME; Particle Mesh Ewald for long-range electrostatics
pme_order= 4; cubic interpolation
fourierspacing= 0.12; grid spacing for FFT
; Temperature coupling is on
tcoupl= V-rescale; modified Berendsen thermostat
tc-grps= SOL  ; two coupling groups - more accurate
tau_t= 0.1; time constant, in ps
ref_t= 293  ; reference temperature, one for each group, in K
; Pressure coupling is on
pcoupl= Parrinello-Rahman; Pressure coupling on in NPT
pcoupltype= isotropic; uniform scaling of box vectors
tau_p= 2.0; time constant, in ps
ref_p= 1.0; reference pressure, in bar
compressibility = 4.5e-5; isothermal compressibility of
water, bar^-1
; Periodic boundary conditions
pbc= xyz; 3-D PBC
; Dispersion 

Re: [gmx-users] Inconsistence of SPC/E and TIP4P water in Gromacs

[Mark Abraham]

Hi,

On Thu, May 11, 2017 at 3:31 AM ZUO Taisen  wrote:

>
> Hi guys:
>
>  I have compared the SPC/E and TIP4P water of the opls force field in
> Gromacs.But there is something inconsistent of the potential energy in the
> system
>
>
> The potential energy of SPC/E water is -4.71400e+04/1000=-47.14kJ/mol
> which is far from the results from literature (Peter G. Kusalik,
> Science,vol,265,26,p1219-1221,1994)
> The potential energy TIP4P water is -4.17875e+04/1000=41.787kJ/mol which
> is very close to the literature -41.8kJ/mol(Peter G. Kusalik,
> Science,vol,265,26,p1219-1221,1994)
>
>
> I have attached all the files simulating the SPC/E and TIP4P water(.trr
> .xtc files are delected) and also the related paper in my last email but
> was suspended because of the files were too big.So no files was attached
> this time but the .mdp file was posted as following.
>

We shouldn't have a mailing list where thousands of recipients can receive
large emails :-)


> I'm eager to know why! Thank you very much!
>

Force fields, including water models, are parameterized to reproduce
certain theoretical or experimental data. To the extent that those targets
resemble other observables, you might expect agreement. You should expect
the level of such agreement to differ between models.

Mark


> TIP4P simulation results:
> Statistics over 3001 steps using 31 frames
>
>Energies (kJ/mol)
> LJ (SR)  Disper. corr.   Coulomb (SR)   Coul. recip.  Potential
>   7.71821e+03   -2.05973e+02   -4.93956e+049.58136e+01
>  -4.17875e+04
> Kinetic En.   Total EnergyTemperature Pres. DC (bar) Pressure (bar)
> 7.30494e+03   -3.44826e+042.93006e+02   -1.14009e+021.41117e+00
>
>
> SPC/E simulation results:
>
> Statistics over 3001 steps using 31 frames
>
>
>Energies (kJ/mol)
> LJ (SR)  Disper. corr.   Coulomb (SR)   Coul. recip.  Potential
> 9.23352e+03   -2.12206e+02   -5.62576e+049.63633e+01   -4.71400e+04
> Kinetic En.   Total EnergyTemperature Pres. DC (bar) Pressure (bar)
> 7.30494e+03   -3.98350e+042.93005e+02   -1.18040e+021.64654e+00
>
>
>
>
> title= OPLS MD waterethanol siimulation
> ; Run parameters
> integrator= md; leap-frog integrator
> nsteps= 3000; 1 * 2000 = 2 ps (20 ns)
> dt= 0.001; 1 fs
> ; Output control
> nstxout= 10; save coordinates every 10.0 ps
> nstvout= 10; save velocities every 10.0 ps
> nstenergy= 10; save energies every 10.0 ps
> nstlog= 10; update log file every 10.0 ps
> nstxout-compressed  = 10  ; save compressed coordinates every 10.0
> ps
> define  =-DEFLEXIBLE ; flexible
> water
> compressed-x-grps   = System; replaces xtc-grps
> ; Bond parameters
>  continuation= yes; Restarting after NPT
>  constraint_algorithm= lincs; holonomic constraints
>  constraints= none; H bonds (even heavy atom-H bonds)
> constrained
>  lincs_iter= 1; accuracy of LINCS
>  lincs_order= 4; also related to accuracy
> ; Neighborsearching
> cutoff-scheme   = Verlet
> ns_type= grid; search neighboring grid cells
> nstlist= 10; 20 fs, largely irrelevant with Verlet scheme
> rcoulomb= 1.2; short-range electrostatic cutoff (in nm)
> rvdw= 1.2; short-range van der Waals cutoff (in nm)
> ; Electrostatics
> coulombtype= PME; Particle Mesh Ewald for long-range electrostatics
> pme_order= 4; cubic interpolation
> fourierspacing= 0.12; grid spacing for FFT
> ; Temperature coupling is on
> tcoupl= V-rescale; modified Berendsen thermostat
> tc-grps= SOL  ; two coupling groups - more accurate
> tau_t= 0.1; time constant, in ps
> ref_t= 293  ; reference temperature, one for each group, in K
> ; Pressure coupling is on
> pcoupl= Parrinello-Rahman; Pressure coupling on in NPT
> pcoupltype= isotropic; uniform scaling of box vectors
> tau_p= 2.0; time constant, in ps
> ref_p= 1.0; reference pressure, in bar
> compressibility = 4.5e-5; isothermal compressibility of
> water, bar^-1
> ; Periodic boundary conditions
> pbc= xyz; 3-D PBC
> ; Dispersion correction
> DispCorr= EnerPres; account for cut-off vdW scheme
> ; Velocity generation
> gen_vel= no; Velocity generation is off
>
>
>
> --
>
>
> Taisen Zuo
>
> China Spallation Neutron Source，Institute of High Energy Physics， Chinese
> Academy of Science
> A1-510, Zhongziyuan road NO.1, Dongguan, Guangdong, PR China. 523770
> Tel: 86-0769-89156495
> Cell: 13650469795
>
> --
> 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
> 

[gmx-users] Inconsistence of SPC/E and TIP4P water in Gromacs

[ZUO Taisen]

Hi guys:
   
 I have compared the SPC/E and TIP4P water of the opls force field in 
Gromacs.But there is something inconsistent of the potential energy in the 
system 


The potential energy of SPC/E water is -4.71400e+04/1000=-47.14kJ/mol which is 
far from the results from literature (Peter G. Kusalik, 
Science,vol,265,26,p1219-1221,1994) 
The potential energy TIP4P water is -4.17875e+04/1000=41.787kJ/mol which is 
very close to the literature -41.8kJ/mol(Peter G. Kusalik, 
Science,vol,265,26,p1219-1221,1994)


I have attached all the files simulating the SPC/E and TIP4P water(.trr .xtc 
files are delected) and also the related paper in my last email but was 
suspended because of the files were too big.So no files was attached this time 
but the .mdp file was posted as following.


I'm eager to know why! Thank you very much!


TIP4P simulation results:
Statistics over 3001 steps using 31 frames

   Energies (kJ/mol)
LJ (SR)  Disper. corr.   Coulomb (SR)   Coul. recip.  Potential
  7.71821e+03   -2.05973e+02   -4.93956e+049.58136e+01   
-4.17875e+04
Kinetic En.   Total EnergyTemperature Pres. DC (bar) Pressure (bar)
7.30494e+03   -3.44826e+042.93006e+02   -1.14009e+021.41117e+00


SPC/E simulation results:

Statistics over 3001 steps using 31 frames


   Energies (kJ/mol)
LJ (SR)  Disper. corr.   Coulomb (SR)   Coul. recip.  Potential
9.23352e+03   -2.12206e+02   -5.62576e+049.63633e+01   -4.71400e+04
Kinetic En.   Total EnergyTemperature Pres. DC (bar) Pressure (bar)
7.30494e+03   -3.98350e+042.93005e+02   -1.18040e+021.64654e+00




title= OPLS MD waterethanol siimulation 
; Run parameters
integrator= md; leap-frog integrator
nsteps= 3000; 1 * 2000 = 2 ps (20 ns)
dt= 0.001; 1 fs
; Output control
nstxout= 10; save coordinates every 10.0 ps
nstvout= 10; save velocities every 10.0 ps
nstenergy= 10; save energies every 10.0 ps
nstlog= 10; update log file every 10.0 ps
nstxout-compressed  = 10  ; save compressed coordinates every 10.0 ps
define  =-DEFLEXIBLE ; flexible water
compressed-x-grps   = System; replaces xtc-grps
; Bond parameters
 continuation= yes; Restarting after NPT 
 constraint_algorithm= lincs; holonomic constraints 
 constraints= none; H bonds (even heavy atom-H bonds) constrained
 lincs_iter= 1; accuracy of LINCS
 lincs_order= 4; also related to accuracy
; Neighborsearching
cutoff-scheme   = Verlet
ns_type= grid; search neighboring grid cells
nstlist= 10; 20 fs, largely irrelevant with Verlet scheme
rcoulomb= 1.2; short-range electrostatic cutoff (in nm)
rvdw= 1.2; short-range van der Waals cutoff (in nm)
; Electrostatics
coulombtype= PME; Particle Mesh Ewald for long-range electrostatics
pme_order= 4; cubic interpolation
fourierspacing= 0.12; grid spacing for FFT
; Temperature coupling is on
tcoupl= V-rescale; modified Berendsen thermostat
tc-grps= SOL  ; two coupling groups - more accurate
tau_t= 0.1; time constant, in ps
ref_t= 293  ; reference temperature, one for each group, in K
; Pressure coupling is on
pcoupl= Parrinello-Rahman; Pressure coupling on in NPT
pcoupltype= isotropic; uniform scaling of box vectors
tau_p= 2.0; time constant, in ps
ref_p= 1.0; reference pressure, in bar
compressibility = 4.5e-5; isothermal compressibility of water, 
bar^-1
; Periodic boundary conditions
pbc= xyz; 3-D PBC
; Dispersion correction
DispCorr= EnerPres; account for cut-off vdW scheme
; Velocity generation
gen_vel= no; Velocity generation is off 



--


Taisen Zuo

China Spallation Neutron Source，Institute of High Energy Physics， Chinese 
Academy of Science
A1-510, Zhongziyuan road NO.1, Dongguan, Guangdong, PR China. 523770
Tel: 86-0769-89156495
Cell: 13650469795

-- 
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-requ...@gromacs.org.