tip4p with rigid epm2

Axel and LC Liu,

Thanks for your kind help,

Actually we are running the same script, pretty much (I am not using shake). May I ask you, why are you reassigning the charges to the atoms? I think that a neighbor 3.0 has solved the problem for me too (I hope).

Otto

The reassigning charge part is because I was testing across water models. You can ignore that if you stick to one model. For a flexible model, it is a different story. You need to keep the time step at least 0.25 fs, or smaller, to catch the vibration modes. Time step is decided by the fastest event in the system. I used neighbor 2.0 before, and there is no problem with the rigid tip4p/2005…so I do not think it would be the root cause. Anyway, there are three parts that may cause the lost atom thing; the input script, initial atom configuration, and LAMMPS itself. I will suggest you to test run my files on your machine. If it works, at least you can have a common ground to talk with others. If, not, even better. The problem lies on the machine side. Otherwise it is truly difficult to debug on the air.

LC Liu

2013/11/6 下午6:52 於 “Moultos, Othon” <othon.moultos@…4184…> 寫道:

Dear all,

Please find attached a density plot for the Tip4p2005 water, simulated with
i) fix rigid/npt
ii) fix shake

This is the only difference in my in files (find them also attached). As you can see in the plots the fix shake reproduces the actual
water density as computed by Vega. The problem is the fix rigid/npt which gives a raised value.

Of course one can (as LC Liu) use the fix shake to run the specific force field but if another want to mix it with another rigid molecule (as me)
you have to use the fix rigid way.

Does anyone have any comments about possible problems I put into the system or any suggestions?

For the simulation of rigid water and co2
Now I am trying to use fix shake only for the water molecules along with fix npt while fix rigid/npt for the co2s.

Thanks a lot,
Otto

density.png

in.rigid (2.69 KB)

in.shake (2.76 KB)

config (868 KB)

PS. The png has a mistake. Where I say lammps it was meant to say shake.

Sorry

Dear all,

Please find attached a density plot for the Tip4p2005 water, simulated with
i) fix rigid/npt
ii) fix shake

This is the only difference in my in files (find them also attached). As you can see in the plots the fix shake reproduces the actual
water density as computed by Vega. The problem is the fix rigid/npt which gives a raised value.

Of course one can (as LC Liu) use the fix shake to run the specific force field but if another want to mix it with another rigid molecule (as me)
you have to use the fix rigid way.

Does anyone have any comments about possible problems I put into the system or any suggestions?

For the simulation of rigid water and co2
Now I am trying to use fix shake only for the water molecules along with fix npt while fix rigid/npt for the co2s.

Thanks a lot,
Otto

density.png

in.rigid (2.69 KB)

in.shake (2.76 KB)

Dear all,

Please find attached a density plot for the Tip4p2005 water, simulated with
i) fix rigid/npt
ii) fix shake

This is the only difference in my in files (find them also attached). As you can see in the plots the fix shake reproduces the actual
water density as computed by Vega. The problem is the fix rigid/npt which gives a raised value.

Of course one can (as LC Liu) use the fix shake to run the specific force field but if another want to mix it with another rigid molecule (as me)
you have to use the fix rigid way.

Does anyone have any comments about possible problems I put into the system or any suggestions?

have you tried a shorter time step? 0.5fs with fix rigid seems a bit
on the large side to me. it may not outright crash, but the
integration of the rotation may not be as accurate as you want (the
major downside of fix rigid vs. fix shake).

axel.

Dear Axel,

Thanks again for the answer. I have tried smaller timesteps (0.25 and 0.1 fs) with similar results. With dt=0.1 fs my system ends up with not such a high density as before but still very high compared to the shake
implementation.

Do you have any other ideas. I really cannot see the way out.

Did anybody succesfully tried a TIP4P2005 with fix rigid/npt before?
Any ideas will be much appreciated.

Otto

You can ask Trung (CCd) about fix rigid/npt in this scenario.

I doubt that it will be any different than fix npt.

Steve

Dear Rajdeep,

Find attached the recent history of our mailing list, concerning my "problems" with tip4p. The input scripts of liu and mine with rigid can for sure predict the density of water tip4p. Try them.

Otto

density.png

in.rigid (2.69 KB)

in.shake (2.76 KB)

config (868 KB)

ATT00001.txt (447 Bytes)

ATT00002.txt (171 Bytes)

Dear Otto,
I am trying to find out the mistake, I have done in my script and examining your script. Thanks a lot for your prompt reply.
I will get back to you I get any difficulties regarding this.

BTW, which was the initial state (i.e. temp, pressure) your data file (config) of water ?

Thanks again.

Rajdeep.

My initial configuration was just a replicated water molecule (with the replication of lammps). So it was a kind of lattice, which very quickly (less than 1ns) becomes fully equilibrated.

Hi Otto,

sorry for my late response. I notice in the data file that there are atoms that are well out of the box boundaries (like atom id’s 37, 38 and 39; 76, 77 and 78; 115, 116 and 117, etc.) without any explicit image flags indicated. These atoms could be related to the deviation in the box volume (and hence density) as you observed, due to the way the virial contribution of fix rigid is determined.

As an evidence for this proposition, I replaced your data file with another file consisting of a single water molecule and replicate it by 16 by 16 by 16 (see the attached file config), the densities from fix shake and fix rigid seem to match (see the log files attached) after equilibration.

Let me know if that resolves the issue.

Best,

-Trung

log.rigid (344 KB)

log.shake (381 KB)

config (441 Bytes)

Dear Trung,

Thanks for your concern. That was one of the mistakes I did there. The other one was a silly typo in the topology (lower HOH angle). My densities are now well computed and I can say for sure that there is no problem with tip4p-like force fields.

Thanks again,

Otto