Polarizable Water Models (Technical)

One frustrating thing is that I often have ideas about stuff that I can't write up as a paper because it isn't substantial enough. Because theirs no real new results or theories I have its just opinion and speculation. So I figure I could put it here, so its recorded.

So one thing that puzzles me a lot about Molecular Dynamics water models is that as far as I'm aware no one has built one that incorporates induced multipole moments higher than the dipole moment.

This doesn't make sense to me, it has been shown that the electric field varies substantially over the range of a single water molecule.   That means it will have significant derivatives and hence significant induced quadrupole and octupole moments. There are models that use Drude oscillators but this is not the same thing as I'm pretty sure these will be dominated by the dipole moment as the point charges are so close together. 

I think this causes the common problem with water models. That the polarisability has to be artificially reduced from 1.44 cubic angstroms to about 1 as otherwise the dielectric constant of water comes out too large at around 100. (paper here) This paper argues that on average the electric field is lower around a water molecule than it is at the centre, and that that is why you should lower the polarisability. But strictly speaking the correct way to account for that is to use higher order induced multipole moments as you are taking a taylor expansion about the centre of the molecule. Additional evidence for this comes from a relatively newer water model which has Drude oscillators on the hydrogens as well as the oxygens. This model seems to able to reproduce the dielectric constant of water with the vacuum polarisability value of 1.44 angstroms. That implies that it is the variation in the electric field that explains this effect, and the correct way to account for this is with a multipole expansion. I guess it's too computationally expensive though, although it looks like a quantum drude oscillator model might work.

They also apply their model to ion solvation energies and get nice results for the salt pairs. But they don't calculate the surface potential of their model. So there is no way to know what they're getting for the real or intrinsic values. Their single ion values in periodic boundary conditions are really close to Tissandier's values, which is surprising but I guess you can't read too much into that as the quadrupole trace of their water molecule could be anything. Surely it's not too hard to calculate this, would be really useful for working out what the real and intrinsic values are, which is super important.