Re: LJ Parameters on hydrogen for TIP4P water models

From: Shayon Bhattacharya (shayon.bhattacharya_at_gmail.com)
Date: Wed Sep 30 2015 - 16:39:45 CDT

Thanks Brian and Peter for your valuable inputs. Finally got rid of the hydrogen vdw sites (epsilon 0, Rmin/2 0). For the force constant, I received a warning and then an error while using a zero Kb value. Hence, using a non-zero (450) value.

regards,
Shayon
 
On 30-Sep-2015, at 2:14 PM, Peter Freddolino <petefred_at_umich.edu> wrote:

> On the tail coefficient point, I wonder if the most proper thing would be to exclude water hydrogens altogether from the averaging process, since they're basically not supposed to have vdw sites anyway...
>
> On the nonzero force constants, NAMD does not impose constraints during minimization, so having nonzero force constants helps maintain a reasonable-ish geometry during minimization. As soon as dynamics start, the constraints snap into place and the force constants don't matter. If I recall correctly, NAMD also uses the lengths from the parameter file to figure out the geometry to constrain the water geometry for some water models.
>
> Best,
> Peter
>
>
> On Tue, Sep 29, 2015 at 3:52 PM, Brian Radak <brian.radak.accts_at_gmail.com> wrote:
> I believe Lennart Nilsson (sp?) published a paper comparing the thermodynamics and dynamics of TIP3P and the modified CHARMM TIP3P. It is an old paper, but I believe the conclusion was that the bulk properties are statistically unchanged.
>
> I would think the main difference with the added hydrogen parameters would be in water-solute interactions. On my end that is just speculation, but the CHARMM force field community would seem to believe it enough that they continue to use the modified parameter set. The added parameters might also impact any tail corrections by altering the average LJ coefficients.
>
> The non-zero force constants might have a different meaning in CHARMM vs NAMD. I believe NAMD rather explicitly throws out length parameters for constrained bonds. For lone pairs I would hope they are never used, as that would be a bit of waste of operations. CHARMM also uses length parameters when building structures from the IC table (clearly Psfgen does not), so sometimes the lengths are put in for that purpose.
>
> Regards,
> Brian
>
>
> On Tue, Sep 29, 2015 at 2:04 PM, Shayon Bhattacharya <shayon.bhattacharya_at_gmail.com> wrote:
> Hi Peter,
>
> Very reassuring to have that response from you, and the link to the paper - I was thinking along the same lines. Thank you! Since I have already run a few simulations using your tip4 parameter file and Jordi Faraudo's tip4p/2005 par file, I don't really want to start all over again by revamping the parameters. Though you have pointed out that hydrogen vdw parameters are the least of my problems, do you think having a small value on hydrogen epsilon or Rmin/2 will significantly affect simulation results as opposed to keeping them as zeroes?
>
> Also the MacKerell group has a Kb value of 0 for bond OT LP, whereas you used 450 for the bond stretching force constant for the dummy atom with oxygen - same value as for HT OT force constant. How does it matter?
>
> Thank you again in advance.
>
> cheers,
> Shayon
>
> On 29-Sep-2015, at 7:25 PM, Peter Freddolino <petefred_at_umich.edu> wrote:
>
>> Hi Shayon,
>> From what I've always seen, the effects of the hydrogen vdw site is basically non-existent because of how small it is relative to the oxygen, as you note. It was also shown in http://pubs.acs.org/doi/abs/10.1021/ct900549r that the hydrogen vdw parameter makes almost no difference. In all of my TIP4 files I used the original TIP4 definition of no hydrogen vdW site (which is what you get with epsilon=0); the comments in that tip4p.par file that you link to are somewhat misleading in that respect (oops). I don't know what Jordi Faraudo's logic is for having hydrogen vdw sites there; people put the hydrogen vdw site on charmm style tip3 water (TIPS3P), but I am not aware of anyone advocating that for any of the TIP4 models.
>>
>> From my understanding the safest recommendation is to use the TIPS3P model with charmm force fields (that is what they were parameterized with) even though TIP3P is virtually equivalent (and might be preferred for performance reasons). Things get dicey very fast if you try using ANY TIP4 models with charmm forcefields, because of how deeply the water model is baked into the charmm parameterization procedure. People do it, sometimes it seems to work well, but there's no reason to expect that it actually SHOULD give reasonable results. So I think the best advice there is (1) think three times about whether that's really a good idea, and (2) if you do, the water hydrogen vdw parameters are probably the least of your problems, but if MacKerell's group has them turned off in his TIP4 files, that's probably the place to be.
>>
>> Best,
>> Peter
>>
>> On Tue, Sep 29, 2015 at 11:20 AM, Shayon Bhattacharya <shayon.bhattacharya_at_gmail.com> wrote:
>> Dear NAMD list,
>>
>> I am running simulations with various TIP4P water models (TIP4P, TIP4P-Ew, TIP4P/2005, etc) for my protein to test a certain feature and choose the right model for my simulation - all using CHARMM36 FF for proteins. I have used the tip4p kit from Peter - http://www.ks.uiuc.edu/Research/namd/mailing_list/namd-l.2011-2012/0322.html with the pre-equilibrated water box, and also using this for the other water models (generating psfs of boxes with different topologies).
>>
>> I see that the Mackerell website have these tip4p water model parameters (in toppar stream files) with no values (they are zeroes) on the Hydrogen Lennard-Jonnes parameters - epsilon and Rmin/2 as used in charmm convention, whereas the parameter file from Peter has a Rmin/2 value of 0.0250, while the epsilon is 0. Same applies for tip4p/2005, where the parameter file from Jordi Faraudo has hydrogen LJ values of -0.046 and 0.2245 respectively for epsilon and Rmin/2, whereas the one from Mackerell website had zeroes in place.
>>
>> Knowing that the oxygen LJ parameters are what matters, and the vdw radii for oxygen is bigger than hydrogen, my question is do the LJ parameters for hydrogen really matter? If so, why such discrepancies in values for NAMD compatible Charmm conventions?
>>
>> thank you,
>> Shayon
>> --
>> Shayon Bhattacharya
>> PhD student
>> Department of Physics and Energy
>> University of Limerick, Ireland
>>
>
>
>
>
> --
> Brian Radak
> Postdoctoral Scholar
> University of Chicago
> Department of Biochemistry & Molecular Biology
> Gordon Center for Integrative Science, W323A
> 929 E. 57th St.
> Chicago, IL 60637-1454
> Tel: 773/834-2812
> e-mail: radak_at_uchicago.edu
>

This archive was generated by hypermail 2.1.6 : Tue Dec 27 2016 - 23:21:21 CST