CHARMM Development Project Forums User Discussion & Questions General Chemistry Discussions TIP3 Bulk water hydrogen bonding

Hello,

I am trying to test a short analysis script I wrote to calculate the number of Hydrogen bonds per molecule using bulk TIP3 water. I found a few references about hbond criteria, but I am not finding much about what the correct answer should be or how accurately I should be able to replicate that with a simple D--H...A angle and DA length definition.

I would appreciate any references you know of and any advice on how accurate!

Thanks,

Jason

r(H..A)<2.4A with no angle restriction is a reasonable criterion (De Loof, H., Nilsson, L., and Rigler, R. (1992). JACS 114, 4028-4035); this is also the default for COOR HBOND (corman.doc). Note that in contrast to the usual crystallographic situation when hydrogen positions are not known it is easier to define hydrogen bonds in simulations using all-atom representations.

In the next CHARMM release (c33) COOR HBONd can handle periodic systems. Without this, surface effects can influence your result. A workaround in order to get the correct average number of hbonds/water is to use TRAJ READ and COOR HBOND in a CHARMM loop. For the first atom selection to COOR HBOND use a selection that leaves out the outermost layer of water, and for the second selection all waters should be included.

Thanks lennart!

I appreicate all the help you provide.

Lennart,

What about TIP3P self diffusion coefficient? I have found references that say 5.06 * 10^-5 cm2/s, but also ~2.8?

Thanks,
Jason

The calculated self-diffusion depends on the non-bond treatment, esp. the electrostatics; Ewald and spherical cutoff methods give fairly different results. TIP3P diffusion is, in general, too fast; TIP4P is better in this regard, but adding the additional 'atom' decreases the simulation rate (fewer ps/day). If you search the forums for TIP4P, you'll find a ref to a paper describing a TIP4Pew model, modified for Ewald summation use by a group at IBM Almaden and one of California schools (UCLA?).

TIP3P self diffusion is around 5.7 * 10^-5 cm2/s with Ewald and with (sensible) spherical truncation. See
J. Phys. Chem. A, 2001. 105(43): p. 9954-9960
J. Comp. Chem, 2002. 23(13): p. 1211-1219.

We observed about 5.1 for Ewald, and roughly 3.9 for a typical 12 A cutoff in

J. Phys. Chem. 1996, 100, 17011-17020

where the expt value is about 3x10^5 cm^2/sec. The results are method and cutoff dependent, and are always too fast for TIP3P. The shear viscosity showed similar and consistent trends, i.e. it is too small (not viscous enough), and further from expt for Ewald.

I think the bottom line is this-- TIP3P has reasonable solvation energetics with CHARMM parameter sets (by design) and the density isn't too bad, but transport properties are incorrect (too fast). This isn't necessarily surprising, as the water model was originally developed via Monte Carlo methods (not dynamics). If you want to seriously study water properties, you should look into other models, esp. polarizable ones. Given the development of the CHARMM force field wrt. TIP3P and the highly similar nature of TIP4P parameters, it seems reasonable to use TIP4P related models, esp. if other CHARMM parameters will be used.

Again Thanks! I will look at using TIP4P.