TIP4P for HOD water

Dear all,

I am working on the dynamics of deuterated water molecule.
For SPC, it is fairly easy to study H2O, D2O and HOD.
But I do not know how to apply TIP4P for HOD, as lj/cut/tipe4p/long only allows the definition of ONE type of hydrogen.
I tried to overcome this by resetting the mass of hydrogen based on atom ID. However, “set mass” only works for particles with per-atom mass attribute.

My question is can LAMMPS do TIP4P for HOD?

Thanks,

Shi

Dear all,

I am working on the dynamics of deuterated water molecule.
For SPC, it is fairly easy to study H2O, D2O and HOD.
But I do not know how to apply TIP4P for HOD, as lj/cut/tipe4p/long only
allows the definition of ONE type of hydrogen.* *
I tried to overcome this by resetting the mass of hydrogen based on atom
ID. However, "set mass" only works for particles with per-atom mass
attribute.

My question is can LAMMPS do TIP4P for HOD?

yes.

there are in fact two ways to do TIP4P water in LAMMPS. one is implementing
a 3-site water and use the tip4p pair styles and fix shake to keep the
water molecules rigid. this is the most efficient approach.
however, you can also define water a 4 site model and then use one of the
rigid fixes. to work around some sanity checks in lammps, you have to
assign the massless site M a tiny mass (say 1.e-10). this is less
efficient, since the rigid integrators typically require a shorter time
step and you have more potentially interacting sites (but there are some
tricks to make that less of a burden).

the first approach would require you to modify the tip4p pair and pppm
styles to handle two hydrogen types. the second can be used right away. an
example input for that can be found here:
http://git.icms.temple.edu/git/?p=lammps-icms.git;a=tree;f=bench-accel/bench_tip4p-rigid;hb=HEAD

axel.

Dear Axel and Ahmed,
Thanks a lot.
But what if I would like to try, q-TIP4P/F, a flexible water potential.
In this case, the position of the massless site M will change during the simulation, and therefore, a rigid molecule can not be used.
Then I guess code modifying is unavoidable for q-TIP4P/F?
Or, we could use “variable” plus “set” command to dynamically change the position of site M?

Hi, Shi:

If you’re not planning to use deuterated hydrogens, then I don’t think you need to modify the code. TIP4P doesn’t require that the bonds be controlled by SHAKE (or at least it shouldn’t). Moreover, since the formula for the “M” site remains unchanged from the original model, there should be no problem with calculating the massless site location.

If you change things to include deuterium, all bets are off.

Hope this helps,

—AEI