TIP4P/Ice (2005) force field for simulation of ice 1h

LAMMPS LAMMPS Mailing List Mirror
1

Dear lammps Users,

I am trying to simulate ice 1h structure which is hexagonal in nature using TIP4P/Ice (2005) force field. My objective is to obtain phonon density of states of water at T = 100 K (though temperature does not effect DOS, but still I am doing at low temperature to avoid diffusion effects) using a force field which properly describes this phase. And compare it with the behaviour of water at similar T & P using a reaxFF implementation for a system, which also has water in it. This way I can compare how good the reaxFF implementation is.

Since ice 1h is the most stable phase at my required condition of temperature and pressure, hence I am trying to obtain its phonon DOS.

I have constructed the input script based on the information related to charges, and LJ parameters on the lammps web site (attached with this email).

But after a trying a lot of things, I am still unable to obtain a stable structure at even 10K.

So here I list all the problems I encountered and solutions I tried:

  1. I encountered the first problem in kspace_style, since TIP4P/Ice require kspace_style pppm/tip4p for coloumb interaction calculations, and this style only works with orthogonal strcutures.

So I changed my hexagonal structure to a orthogonal one, by some basic transformations.

  1. Then I tried to simulate this orthogonal structure and minimize it using style CG/hftn but both did not work, and the system blew up in less than 50 MD steps with a time step of 1 fs.

  2. Then I tried other ways to minimize it, using suggestions on lammps page, i.e. putting viscous damping using fix viscous, and fix nve/limit. I tried both of them alone and combined initially, to obtain a stable system. But even that did not work, as system continuously increases in temperature and results in blow up.

  3. I have also tried making initial time step extremely small, so as to make integration very fine, but even that does not control system temperature.

I really do not know now what to do now.

I think there can be some problem with my data file but I do not know, what specific problem. As structure is okay on visualization. There are no free H or O atoms in the system. As everything is bonded. There are 1152 water molecules in my system, and the total number of atoms is 3456. Also I have maintained the order of O atom id < H atom id and there arrangement in the data file.

So I do not know if there is anything wrong with it, as I have tried my best to avoid any mistakes in the data file which I am sending to everybody’s scrutiny, to save your time and effort. Though any suggestion related to possible problems in it will be really helpful.

Can someone recommend me what to do based on my objective and problems that I encountered.

I am attaching my input and data file for clarity.

I will be really grateful for suggestions.

Thanks,
Ankit

in.H2O (2.72 KB)

data.water (203 KB)