<html>
<head>
<style><!--
.hmmessage P
{
margin:0px;
padding:0px
}
body.hmmessage
{
font-size: 10pt;
font-family:Tahoma
}
--></style>
</head>
<body class='hmmessage'>
Hi,<br><br>The choice of constraint algorithm is irrelevant for the results.<br><br>As I understood Charmm tip3p has less structure than standard tip3p<br>and it would not surprise me if the density is different.<br>But I have never closely looked at Charmm tip3p properties myself.<br><br>Berk<br><br>&gt; Date: Tue, 14 Sep 2010 13:11:47 +0200<br>&gt; Subject: RE: [gmx-users] Regular vs. CHARMM TIP3P water model<br>&gt; From: nicolas.sapay@cermav.cnrs.fr<br>&gt; To: gmx-users@gromacs.org<br>&gt; <br>&gt; &gt;<br>&gt; &gt; Hi,<br>&gt; &gt;<br>&gt; &gt; I don't understand what exactly you want to reproduce.<br>&gt; &gt; Standard tip3p and Charmm tip3p are different models, so the density does<br>&gt; &gt; not have to be identical.<br>&gt; &gt; The Gromacs Charmm FF implementation paper:<br>&gt; &gt; http://pubs.acs.org/doi/full/10.1021/ct900549r<br>&gt; &gt; gives 1002 for tip3p and 1015 for charmm tip3p (this is with LJ switched<br>&gt; &gt; from 1 to 1.2 nm).<br>&gt; &gt; I have in an old paper 986 for tip3p, with 0.85 nm cut-off plus dispersion<br>&gt; &gt; correction.<br>&gt; &gt; These numbers are quite different. I'll do some checking.<br>&gt; &gt;<br>&gt; &gt; Berk<br>&gt; <br>&gt; Actually, TIP3P and CHARMM TIP3P are supposed to give (almost) identical<br>&gt; results as they are not exactly different models but different<br>&gt; implementations of a same model. That is why I'm surprised of the rather<br>&gt; large difference in the density values observed in my simulations or in<br>&gt; the JCTC paper.<br>&gt; <br>&gt; In CHARMM, the additionnal Lennard-Jones parameters for the TIP3P hydrogen<br>&gt; are not here to better reproduce water properties but for practical<br>&gt; purpose - normally, the hydrogens are almost entirely included in the<br>&gt; oxygen VDW radius. As far as I understand, those parameters are here to<br>&gt; improve the stability of simulations by limiting artificial close contacts<br>&gt; between a large negative charge and the point charge on the hydrogen<br>&gt; atoms.<br>&gt; <br>&gt; In fact, I was wondering if the Settle algorithm is appropriate to CHARMM<br>&gt; TIP3P. In the original CHARMM RTF file, there is a SHAKE constraint on the<br>&gt; H-H distance which is not present in the GROMACS tips3p.itp file. So, a<br>&gt; better use of the CHARMM TIP3P model in GROMACS would be to add this<br>&gt; constraint in tips3p.itp and use Shake instead of Settle, isn't it?<br>&gt; <br>&gt; &gt;<br>&gt; &gt;&gt; Date: Tue, 14 Sep 2010 11:00:55 +0200<br>&gt; &gt;&gt; From: spoel@xray.bmc.uu.se<br>&gt; &gt;&gt; To: gmx-users@gromacs.org<br>&gt; &gt;&gt; Subject: Re: [gmx-users] Regular vs. CHARMM TIP3P water model<br>&gt; &gt;&gt;<br>&gt; &gt;&gt; On 2010-09-14 10.23, Nicolas SAPAY wrote:<br>&gt; &gt;&gt; &gt; Hello everybody,<br>&gt; &gt;&gt; &gt;<br>&gt; &gt;&gt; &gt; I have many difficulties to reproduce TIP3P simulation results with<br>&gt; &gt;&gt; CHARMM<br>&gt; &gt;&gt; &gt; TIP3P. Regular TIP3P gives systematically a lower density than its<br>&gt; &gt;&gt; CHARMM<br>&gt; &gt;&gt; &gt; counterpart, independantly from the cutoff for non-bonded<br>&gt; &gt;&gt; interactions,<br>&gt; &gt;&gt; &gt; the version of GROMACS (4.5.1 or 4.0.7) or the double precision.<br>&gt; &gt;&gt; &gt; Regular 961,067 +/-0,756 g/L<br>&gt; &gt;&gt; &gt; CHARMM  980,860 +/-0,492 g/L<br>&gt; &gt;&gt; &gt;<br>&gt; &gt;&gt; &gt; The Enthalpy of vaporization follows a similar scheme:<br>&gt; &gt;&gt; &gt; regular -39,992 +/-0,021 kJ/mol<br>&gt; &gt;&gt; &gt; CHARMM  -40,665        +/-0,009 kJ/mol<br>&gt; &gt;&gt;<br>&gt; &gt;&gt; How about the dispersion correction?<br>&gt; &gt;&gt; If that is not turned on densities will be too low (in both cases).<br>&gt; &gt;&gt;<br>&gt; &gt;&gt; &gt;<br>&gt; &gt;&gt; &gt; In fact, CHARMM gives results closer to what I should obtain at 300K<br>&gt; &gt;&gt; and 1<br>&gt; &gt;&gt; &gt; bar for TIP3P. I suspect an issue with the bond constraints, but I<br>&gt; &gt;&gt; can't<br>&gt; &gt;&gt; &gt; locate precisely where it is. Settle parameters are exactly the same<br>&gt; &gt;&gt; (as<br>&gt; &gt;&gt; &gt; well as constraints for flexible water). Did someone ever face a<br>&gt; &gt;&gt; similar<br>&gt; &gt;&gt; &gt; problem?<br>&gt; &gt;&gt; &gt;<br>&gt; &gt;&gt; &gt; Rgards,<br>&gt; &gt;&gt; &gt; Nicolas<br>&gt; &gt;&gt; &gt;<br>&gt; &gt;&gt; &gt; P.S. I My system is just a box of 1728 water molecules<br>&gt; &gt;&gt; pre-equilibrated<br>&gt; &gt;&gt; &gt; for 500 ps at 300 K and 1 bar.<br>&gt; &gt;&gt; &gt;<br>&gt; &gt;&gt; &gt; P.S. II I'm using the following simulation parameters:<br>&gt; &gt;&gt; &gt; constraints              = hbonds<br>&gt; &gt;&gt; &gt; constraint_algorithm     = LINCS<br>&gt; &gt;&gt; &gt; continuation             = no<br>&gt; &gt;&gt; &gt; lincs-order              = 4<br>&gt; &gt;&gt; &gt; lincs-iter               = 1<br>&gt; &gt;&gt; &gt; lincs-warnangle          = 30<br>&gt; &gt;&gt; &gt; morse                    = no<br>&gt; &gt;&gt; &gt; (LINCS or SHAKE should not make any difference since TIP3P is normally<br>&gt; &gt;&gt; &gt; treated by SETTLE).<br>&gt; &gt;&gt; &gt;<br>&gt; &gt;&gt; &gt; Other parameters are:<br>&gt; &gt;&gt; &gt; ; COUPLING<br>&gt; &gt;&gt; &gt; ; Temperature coupling<br>&gt; &gt;&gt; &gt; tcoupl                   = Berendsen<br>&gt; &gt;&gt; &gt; nsttcouple               = -1<br>&gt; &gt;&gt; &gt; nh-chain-length          = 10<br>&gt; &gt;&gt; &gt; tc_grps                  = System<br>&gt; &gt;&gt; &gt; tau_t                    = 0.1<br>&gt; &gt;&gt; &gt; ref_t                    = 300<br>&gt; &gt;&gt; &gt; ; Pressure coupling<br>&gt; &gt;&gt; &gt; pcoupl                   = Berendsen<br>&gt; &gt;&gt; &gt; pcoupltype               = isotropic<br>&gt; &gt;&gt; &gt; nstpcouple               = -1<br>&gt; &gt;&gt; &gt; tau_p                    = 1<br>&gt; &gt;&gt; &gt; compressibility          = 4.5e-5<br>&gt; &gt;&gt; &gt; ref_p                    = 1<br>&gt; &gt;&gt; &gt;<br>&gt; &gt;&gt; &gt; ; NEIGHBORSEARCHING PARAMETERS<br>&gt; &gt;&gt; &gt; ; nblist update frequency<br>&gt; &gt;&gt; &gt; nstlist                  = 10<br>&gt; &gt;&gt; &gt; ns_type                  = grid<br>&gt; &gt;&gt; &gt; pbc                      = xyz<br>&gt; &gt;&gt; &gt; rlist                    = 1.15<br>&gt; &gt;&gt; &gt;<br>&gt; &gt;&gt; &gt; ; OPTIONS FOR ELECTROSTATICS AND VDW<br>&gt; &gt;&gt; &gt; ; Method for doing electrostatics<br>&gt; &gt;&gt; &gt; coulombtype              = PME<br>&gt; &gt;&gt; &gt; rcoulomb-switch          = 0<br>&gt; &gt;&gt; &gt; rcoulomb                 = 1.15<br>&gt; &gt;&gt; &gt; vdwtype                  = Cutoff<br>&gt; &gt;&gt; &gt; rvdw                     = 1.0<br>&gt; &gt;&gt; &gt; fourierspacing           = 0.10<br>&gt; &gt;&gt; &gt; pme_order                = 6<br>&gt; &gt;&gt; &gt; ewald_rtol               = 1.0e-6<br>&gt; &gt;&gt; &gt; ewald_geometry           = 3d<br>&gt; &gt;&gt; &gt; epsilon_surface          = 0<br>&gt; &gt;&gt; &gt; optimize_fft             = yes<br>&gt; &gt;&gt; &gt;<br>&gt; &gt;&gt; &gt; ; RUN CONTROL PARAMETERS<br>&gt; &gt;&gt; &gt; integrator               = md<br>&gt; &gt;&gt; &gt; dt                       = 0.002<br>&gt; &gt;&gt; &gt; nsteps                   = 200000<br>&gt; &gt;&gt; &gt; comm_mode                = Linear<br>&gt; &gt;&gt; &gt; nstcomm                  = 1<br>&gt; &gt;&gt; &gt; comm_grps                = System<br>&gt; &gt;&gt; &gt;<br>&gt; &gt;&gt; &gt;<br>&gt; &gt;&gt;<br>&gt; &gt;&gt;<br>&gt; &gt;&gt; --<br>&gt; &gt;&gt; David van der Spoel, Ph.D., Professor of Biology<br>&gt; &gt;&gt; Dept. of Cell &amp; Molec. Biol., Uppsala University.<br>&gt; &gt;&gt; Box 596, 75124 Uppsala, Sweden. Phone:        +46184714205.<br>&gt; &gt;&gt; spoel@xray.bmc.uu.se    http://folding.bmc.uu.se<br>&gt; &gt;&gt; --<br>&gt; &gt;&gt; gmx-users mailing list    gmx-users@gromacs.org<br>&gt; &gt;&gt; http://lists.gromacs.org/mailman/listinfo/gmx-users<br>&gt; &gt;&gt; Please search the archive at<br>&gt; &gt;&gt; http://www.gromacs.org/Support/Mailing_Lists/Search before posting!<br>&gt; &gt;&gt; Please don't post (un)subscribe requests to the list. Use the<br>&gt; &gt;&gt; www interface or send it to gmx-users-request@gromacs.org.<br>&gt; &gt;&gt; Can't post? Read http://www.gromacs.org/Support/Mailing_Lists<br>&gt; &gt;                                                --<br>&gt; &gt; gmx-users mailing list    gmx-users@gromacs.org<br>&gt; &gt; http://lists.gromacs.org/mailman/listinfo/gmx-users<br>&gt; &gt; Please search the archive at<br>&gt; &gt; http://www.gromacs.org/Support/Mailing_Lists/Search before posting!<br>&gt; &gt; Please don't post (un)subscribe requests to the list. Use the<br>&gt; &gt; www interface or send it to gmx-users-request@gromacs.org.<br>&gt; &gt; Can't post? Read http://www.gromacs.org/Support/Mailing_Lists<br>&gt; <br>&gt; <br>&gt; -- <br>&gt; [ Nicolas Sapay - Post-Doctoral Fellow ]<br>&gt; CERMAV - www.cermav.cnrs.fr<br>&gt; BP53, 38041 Grenoble cedex 9, France<br>&gt; Phone: +33 (0)4 76 03 76 44/53<br>&gt; <br>&gt; -- <br>&gt; gmx-users mailing list    gmx-users@gromacs.org<br>&gt; http://lists.gromacs.org/mailman/listinfo/gmx-users<br>&gt; Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting!<br>&gt; Please don't post (un)subscribe requests to the list. Use the <br>&gt; www interface or send it to gmx-users-request@gromacs.org.<br>&gt; Can't post? Read http://www.gromacs.org/Support/Mailing_Lists<br>                                               </body>
</html>