Testing recent charge-on-spring type polarizable water models. I. Melting temperature and ice properties

We determined the freezing point of eight molecular models of water. All models use the charge-on-spring (COS) method to express polarization. The studied models were the COS/G2, COS/G3 [H. Yu and W. F. van Gunsteren, J. Chem. Phys. 121, 9549 (2004)], the COS/B2 [H. Yu, T. Hansson, and W. F. van Gunsteren, J. Chem. Phys. 118, 221 (2003)], the SWM4-DP [G. Lamoureux, A. D. MacKerell, Jr., and B. Roux, J. Chem. Phys. 119, 5185 (2003)], the SWM4-NDP [G. Lamoureux, E. Harder, I. V. Vorobyov, B. Roux, and A. D. MacKerell, Jr., Chem. Phys. Lett. 418, 245 (2006)], and three versions of our model, the BKd1, BKd2, and BKd3. The BKd1 is the original Gaussian model [P. T. Kiss, M. Darvas, A. Baranyai, and P. Jedlovszky, J. Chem. Phys. 136, 114706 (2012)] with constant polarization and with a simple exponential repulsion. The BKd2 applies field-dependent polarizability [A. Baranyai and P. T. Kiss, J. Chem. Phys. 135, 234110 (2011)], while the BKd3 model has variable size to approximate the temperature-density (T-ρ) curve of water [P. T. Kiss and A. Baranyai, J. Chem. Phys. 137, 084506 (2012)]. We used the thermodynamic integration (TI) and the Gibbs-Helmholtz equation to determine the equality of the free energy for liquid water and hexagonal ice (Ih) at 1 bar. We used the TIP4P and the SPC/E models as reference systems of the TI. The studied models severely underestimate the experimental melting point of ice Ih. The calculated freezing points of the models are the following: COS/G2, 215 K; COS/G3, 149 K; SWM4-DP, 186 K; BKd1, 207 K; BKd2, 213 K; BKd3, 233 K. The freezing temperature of the SWM4-NDP system is certainly below 120 K. It might even be that the water phase is more stable than the ice Ih at 1 bar in the full temperature range. The COS/B2 model melts below 100 K. The best result was obtained for the BKd3 model which indicates that correct description of the (T-ρ) curve improves the estimation of the freezing point. We also determined and compared the densities of high-pressure polymorphs of ice for these models.