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

David van der Spoel spoel at xray.bmc.uu.se
Tue Sep 14 17:40:39 CEST 2010


On 2010-09-14 17.06, Nicolas SAPAY wrote:
>
>>
>> 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?
>
Charmm (and amber) probably use a different definition of the harmonic 
potential. Note that this is not used with constraints, only when you 
have flexible models.

Gromacs uses degrees instead of radians in the user input.

>
>
>>
>> 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
>
>


-- 
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



More information about the gromacs.org_gmx-users mailing list