Does anyone know how I can use TIP5P water model on CHARMM?
I just want to test the performance of different water model in molecular dynamics.
I guess, we can run with TIP4P/TIP3P water models.
I'm not sure that anyway has actually confirmed that using TIP4P
in CHARMM gives to correct properties, though people are working
on this. If I get additional information I'll post it. Concerning TIP5P,
it has also not yet been implemented, although it should function.
Again, I'll post info as it comes in.
I had a bash on both SPCE TIP 3 n 4 n 5.
The first 3 ive got to work fine in charmm and reproduce literature data....TIP5P does funny stuff.
All sims were 1296 waters in a 33.9 A box (correct physical density). For the implimentation of TIP3,4 I did not use the literature values for the forcefield. In the literature only O has LJ params. The FF i used for TIP3, 4 had LJ for the H as well (mckerrel standard params). TIP5, irrespective of if H has LJ params has v. short HO hbonds (down to about 1.3A), and as a result large voids occur in the water box (clearly a V. unhappy simulation!). All the water mol. geometries are rigid and correct. PSF readout from charmm says all the correct charges are on all the correct atoms. Im outta ideas. Ive been over the LJ params for TIP5 (there is only one!) and it is correct and v. similar to those used for TIP3. I suspect this may be a prob. with the Lonepair command in charmm.
Phil Mason, Brady group, Cornell
as a matter of spam evasion i never post my email anywhere public. I do check the spam account.. for an immediate response check the cornell web directory.
We've done some TIP4P testing and TIP3P re-evaluation ourselves, and can offer the following:
One can get improved density for water (closer to expt and Jorgensen's original result) for TIP3P by removing the VDW radii from the H atoms; they were apparently added at Harvard eons ago to account for problems with model building. There's a 2-3% increase in density as a consequence of the H atom VDW radii. Other properties such as rotational correlation time and the bulk dielectric do not seem to be affected.
There's been a change in the LONEPAIR code for how the forces are calculated, which should be available in c31a2 (and perhaps c30b2 as well) when they are released in February. With the change, we get much better agreement for TIP4P density compared to a group at IBM Almaden using independent code (not CHARMM).
This is v e r y interesting, and what worries me a bit is:
Given the CHARMM27 parameterization against the modified TIP3P, what would the effects be if one changes (or returns to) solvent params vs protein params?.
Perhaps not a big issue, but I've tested a protein (CHARMM27) in SPC and SPC/E. The results points to some unphysically short h-bonds between solvent and carboxylate/carbonyl oxygens when using these two models..as jugded from 8x2 ns independent simulations. TIP3P (modified) did not exhibit such anomalies.Thus, it would be nice to establish one physically consistent-at-any-level-of-theory good water-model and then keep this model for a further force-field optimization of parameters.
Another problem with TIP3P has been the viscosity; compared to expt, it's off by a factor of 3 or so with Ewald methods, which is related to why translational diffusion is wrong by a similar amount when calculated from <r**2> vs. time, either for water itself, or for small molecules in explicit TIP3P water. If you're trying to compute NMR time decay properties or the rate of a conformation change that may have a dependence on solvent viscosity, you have to take into account that the TIP3P model gives a viscosity close to 0.3 cp, rather than 1.0. For spherical truncation methods, it's a bit better, more like 0.5 or 0.6, but Ewald methods are the preferred electrostatic treatment. See the following ref:
Effect of Electrostatic Force Truncation on Interfacial and Transport Properties of Water
Scott E. Feller, Richard W. Pastor, Atipat Rojnuckarin,
Stephen Bogusz, and Bernard R. Brooks
J. Phys. Chem. 100, 17011-17020 (1996)
The viscosity problem was one issue that led to our closer look at TIP4P, especially a new variant adjusted for use with Ewald methods. I believe the O-H bond length, H-O-H angle, and O atom VDW radius are all the same as for TIP3P, so it's a minimal perturbation (just the dipole moment). The down side is the cost-- for solvated systems, simulations take more computer time because of the extra particle on each water molecule.
This is in a way "good news", because I noticed the same trend in using spherical vs PBC/Ewald summation. When it comes to NMR exhange correlation (order parameters), it turned out to yield rather similar patterns/residue on the protein in (TIP3P/SPC or SPCE/E). But as you stated, smaller systems are probably more sensitive to such issues.
Since we're now in a wealth of computer power, I feel it is worth the extra CPU-cost to get answers closer to reality..
Please keep me informed about the TIP4P progress.
Thanks for the ref btw.
We have also looked at the difference between the original TIP3P and the CHARMM implementation, and found the differences to be very small:
Mark P, Nilsson L. 2001. Structure and Dynamics of the TIP3P, SPC and SPC/E Water Models at 298K. J. Phys. Chem. A 105(43):9954-9960.
Mark P, Nilsson L. 2002. Structure and dynamics of liquid water with different long-range interaction truncation and temperature control methods in molecular dynamics simulations. J. Comp. Chem 23(13):1211-1219.