[gmx-users] Regular vs. CHARMM TIP3P water model

Nicolas SAPAY nicolas.sapay at cermav.cnrs.fr
Tue Sep 14 17:06:50 CEST 2010 

>
> Good.
>
> I just calculated tip3p energies for a Charmm tip3p trajectory.
> Standard tip3p gives 1.1 kJ/mol per water molecule higher energies.
> Thus the LJ on H's gives significant extra attraction which increases the
> density.

Great, I obtain a similar result: -40.0 kJ/mol for regular tip3p vs -41.2
kJ/mol for Charmm tip3p (so a difference of 1.2 kJ/mol in favor of Charmm
tip3p).

By the way, I have noted that the SHAKE force constants are different
between Gromacs and Charmm. In the Gromacs implementation, the constant is
376560 kJ/mol/A^2 for the O-H bond and 460,24 kJ/mol/rad^2 for the H-O-H
angle. In the Charmm parameter file, the force constant is 188280
kJ/mol/A^2 for the bond (2 times lower) and 23012 kJ/mol/rad^2 for the
angle (50 times higher). Is there a justification for that?



>
> Berk
>
>> Date: Tue, 14 Sep 2010 15:49:03 +0200
>> Subject: RE: [gmx-users] Regular vs. CHARMM TIP3P water model
>> From: nicolas.sapay at cermav.cnrs.fr
>> To: gmx-users at gromacs.org
>>
>>
>> >
>> > Hi,
>> >
>> > We did some checking.
>> > The value of the density for tip3p reported in the Gromacs Charmm ff
>> > implementation of 1001.7 is incorrect.
>> > This should have been 985.7. The number of 1014.7 for Charmm tip3p is
>> > correct.
>>
>> I have just done a quick test using dispersion correction
>> (DispCorr=EnerPres) and the density value extracted from the edr file by
>> g_energy. I obtain fairly similar results: 984.8 +/-4.5 for TIP3P and
>> 1015.0 +/- 3.6 for Charmm TIP3P.
>>
>> >
>> > I would expect that the difference with your number is mainly due to
>> the
>> > shorter LJ cut-off you are using.
>> >
>> > Berk
>> >
>> > From: gmx3 at hotmail.com
>> > To: gmx-users at gromacs.org
>> > Subject: RE: [gmx-users] Regular vs. CHARMM TIP3P water model
>> > Date: Tue, 14 Sep 2010 11:33:16 +0200
>> >
>> >
>> >
>> >
>> >
>> >
>> >
>> >
>> > Hi,
>> >
>> > I don't understand what exactly you want to reproduce.
>> > Standard tip3p and Charmm tip3p are different models, so the density
>> does
>> > not have to be identical.
>> > The Gromacs Charmm FF implementation paper:
>> > http://pubs.acs.org/doi/full/10.1021/ct900549r
>> > gives 1002 for tip3p and 1015 for charmm tip3p (this is with LJ
>> switched
>> > from 1 to 1.2 nm).
>> > I have in an old paper 986 for tip3p, with 0.85 nm cut-off plus
>> dispersion
>> > correction.
>> > These numbers are quite different. I'll do some checking.
>> >
>> > Berk
>> >
>> >> Date: Tue, 14 Sep 2010 11:00:55 +0200
>> >> From: spoel at xray.bmc.uu.se
>> >> To: gmx-users at gromacs.org
>> >> Subject: Re: [gmx-users] Regular vs. CHARMM TIP3P water model
>> >>
>> >> On 2010-09-14 10.23, Nicolas SAPAY wrote:
>> >> > Hello everybody,
>> >> >
>> >> > I have many difficulties to reproduce TIP3P simulation results with
>> >> CHARMM
>> >> > TIP3P. Regular TIP3P gives systematically a lower density than its
>> >> CHARMM
>> >> > counterpart, independantly from the cutoff for non-bonded
>> >> interactions,
>> >> > the version of GROMACS (4.5.1 or 4.0.7) or the double precision.
>> >> > Regular 961,067 +/-0,756 g/L
>> >> > CHARMM  980,860 +/-0,492 g/L
>> >> >
>> >> > The Enthalpy of vaporization follows a similar scheme:
>> >> > regular -39,992 +/-0,021 kJ/mol
>> >> > CHARMM  -40,665	+/-0,009 kJ/mol
>> >>
>> >> How about the dispersion correction?
>> >> If that is not turned on densities will be too low (in both cases).
>> >>
>> >> >
>> >> > In fact, CHARMM gives results closer to what I should obtain at
>> 300K
>> >> and 1
>> >> > bar for TIP3P. I suspect an issue with the bond constraints, but I
>> >> can't
>> >> > locate precisely where it is. Settle parameters are exactly the
>> same
>> >> (as
>> >> > well as constraints for flexible water). Did someone ever face a
>> >> similar
>> >> > problem?
>> >> >
>> >> > Rgards,
>> >> > Nicolas
>> >> >
>> >> > P.S. I My system is just a box of 1728 water molecules
>> >> pre-equilibrated
>> >> > for 500 ps at 300 K and 1 bar.
>> >> >
>> >> > P.S. II I'm using the following simulation parameters:
>> >> > constraints              = hbonds
>> >> > constraint_algorithm     = LINCS
>> >> > continuation             = no
>> >> > lincs-order              = 4
>> >> > lincs-iter               = 1
>> >> > lincs-warnangle          = 30
>> >> > morse                    = no
>> >> > (LINCS or SHAKE should not make any difference since TIP3P is
>> normally
>> >> > treated by SETTLE).
>> >> >
>> >> > Other parameters are:
>> >> > ; COUPLING
>> >> > ; Temperature coupling
>> >> > tcoupl                   = Berendsen
>> >> > nsttcouple               = -1
>> >> > nh-chain-length          = 10
>> >> > tc_grps                  = System
>> >> > tau_t                    = 0.1
>> >> > ref_t                    = 300
>> >> > ; Pressure coupling
>> >> > pcoupl                   = Berendsen
>> >> > pcoupltype               = isotropic
>> >> > nstpcouple               = -1
>> >> > tau_p                    = 1
>> >> > compressibility          = 4.5e-5
>> >> > ref_p                    = 1
>> >> >
>> >> > ; NEIGHBORSEARCHING PARAMETERS
>> >> > ; nblist update frequency
>> >> > nstlist                  = 10
>> >> > ns_type                  = grid
>> >> > pbc                      = xyz
>> >> > rlist                    = 1.15
>> >> >
>> >> > ; OPTIONS FOR ELECTROSTATICS AND VDW
>> >> > ; Method for doing electrostatics
>> >> > coulombtype              = PME
>> >> > rcoulomb-switch          = 0
>> >> > rcoulomb                 = 1.15
>> >> > vdwtype                  = Cutoff
>> >> > rvdw                     = 1.0
>> >> > fourierspacing           = 0.10
>> >> > pme_order                = 6
>> >> > ewald_rtol               = 1.0e-6
>> >> > ewald_geometry           = 3d
>> >> > epsilon_surface          = 0
>> >> > optimize_fft             = yes
>> >> >
>> >> > ; RUN CONTROL PARAMETERS
>> >> > integrator               = md
>> >> > dt                       = 0.002
>> >> > nsteps                   = 200000
>> >> > comm_mode                = Linear
>> >> > nstcomm                  = 1
>> >> > comm_grps                = System
>> >> >
>> >> >
>> >>
>> >>
>> >> --
>> >> David van der Spoel, Ph.D., Professor of Biology
>> >> Dept. of Cell & Molec. Biol., Uppsala University.
>> >> Box 596, 75124 Uppsala, Sweden. Phone:	+46184714205.
>> >> spoel at xray.bmc.uu.se    http://folding.bmc.uu.se
>> >> --
>> >> gmx-users mailing list    gmx-users at gromacs.org
>> >> http://lists.gromacs.org/mailman/listinfo/gmx-users
>> >> Please search the archive at
>> >> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
>> >> Please don't post (un)subscribe requests to the list. Use the
>> >> www interface or send it to gmx-users-request at gromacs.org.
>> >> Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>> >
>> >
>> > --
>> > gmx-users mailing list    gmx-users at gromacs.org
>> > http://lists.gromacs.org/mailman/listinfo/gmx-users
>> > Please search the archive at
>> > http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
>> > Please don't post (un)subscribe requests to the list. Use the
>> > www interface or send it to gmx-users-request at gromacs.org.
>> > Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>> > --
>> > gmx-users mailing list    gmx-users at gromacs.org
>> > http://lists.gromacs.org/mailman/listinfo/gmx-users
>> > Please search the archive at
>> > http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
>> > Please don't post (un)subscribe requests to the list. Use the
>> > www interface or send it to gmx-users-request at gromacs.org.
>> > Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>>
>>
>> --
>> [ Nicolas Sapay - Post-Doctoral Fellow ]
>> CERMAV - www.cermav.cnrs.fr
>> BP53, 38041 Grenoble cedex 9, France
>> Phone: +33 (0)4 76 03 76 44/53
>>
>> --
>> gmx-users mailing list    gmx-users at gromacs.org
>> http://lists.gromacs.org/mailman/listinfo/gmx-users
>> Please search the archive at
>> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
>> Please don't post (un)subscribe requests to the list. Use the
>> www interface or send it to gmx-users-request at gromacs.org.
>> Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>  		 	   		  --
> gmx-users mailing list    gmx-users at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-users
> Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
> Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-users-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


-- 
[ Nicolas Sapay - Post-Doctoral Fellow ]
CERMAV - www.cermav.cnrs.fr
BP53, 38041 Grenoble cedex 9, France
Phone: +33 (0)4 76 03 76 44/53

More information about the gromacs.org_gmx-users mailing list