Issue with fix qeq/point for SPC/E water simulation

Hello everyone,
I am simulating SPC/E water model. I assigned charges with fix qeq method instead of set command. I observed several issues I am not understanding.

  1. I used fix qeq/point to assign charge to 1 molecule of water without invoking any pair potential for run 0. The charge on O and H were -0.2 and 0.1, respectively for this simulation. And, after NVE dynamics run with shake, the potential energy was -13 Kcal/mol. The parameter files are from Goddard’s (1991) http://pubs.acs.org/doi/pdf/10.1021/j100161a070 paper,

param file:

1 8.741 13.364 0 0 0

2 4.528 13.89 0 0 0

  1. I used fix qeq/point to assign charge on O and H. But, this time I used chi and eta parameters from ffield.reax.CHO file. I multiplied eta by 2. The charge on O and H were -0.125 and 0.0625, respectively, and the potential energy was -5 Kcal/mol.

param file (ffield.reax.CHO):
1 8.5 17.9978 0 0 0

2 5.32 14.87 0 0 0

  1. I used set command to assign charges O=-0.8476 and H=0.4238. After NVE dynamics run with shake, the potential energy was -211 Kcal/mol.

According to Goddard’s paper, charge on H atom is 0.35, even then it’s well below SPC/E model’s charge. Am I applying fix qeq incorrectly? Potential energy is very low in method (1) and (2) compared to method 3. If we take, SPC/E model’s charge as reference point, the first two methods are clearing giving me wrong answers.

Here’s my input file:

You can’t have a cutoff that is only 3 angstroms; at least 10 angstroms or the same as the Coulomb cutoff is reasonable.

Ray

Hello Ray,

Thank you for replying. I did run the simulation with 10 Angstrom and with 3 times the box length. Still the results are the same.

Read the Rappe & Goddard paper and you will see their charge model (charge-charge interaction) is not a point charge model. Find the corresponding charge model.

Ray