NEMD water shear viscosity, TIP4P/2005: possible? (LAMMPS version 16 Mar 2018)

I am applying NEMD (SLLOD) to obtain shear viscosity of water at 293K (20°C) for water models TIP3P, SPC/E & TIP4P/2005.

The setup has worked fine for TIP3P and SPC/E: after applying unit conversions, viscosity values are in good range of the experimental value/literature (slightly more accurate than reported simulations, which are mostly at 300K).

However, for TIP4P/2005 I just can’t seem to make it work. The system setups are the same, the only differences are parameters specific to the water models.

‘’ERROR: Out of range atoms – cannot compute PPPM’’

comes up every time after PPPM Initialization; the simulation never actually gets to carrying out the equilibration run. The error refers to the kspace_style command line. I used pair_style ‘lj/cut/tip4p/long’ with kspace_style ‘pppm/tip4p’.

However, when running on a single processor (instead of 16), the simulation does proceed, with sensible output and viscosity values (accuracy > TIP3P & SPC/E). Dangerous builds: 1 in the equilibration run (10k fs) and 0 for the data gathering run (1ns).

In various attempts to get through this error message I used slightly different potentials with ‘pppm/disp/tip4p’ which gave the error of ‘KSpace style does not yet support triclinic geometries’. Could this also be the case for the ‘pppm/tip4p’ style?

Out of curiosity, I also tried removing the kspace_style line, changing pair_style to ‘lj/cut/tip4p/cut’ - the simulation does proceed, but viscosities come out ~ 80x larger than expected (just to mention in case this helps).

I have run simulations (orthogonal box) for all three water models with no problems before attempting viscosity calculations. I’d ideally like to use NEMD, rather than shear wall or rNEMD, if possible. So far I haven’t come across a direct literature comparison in terms of whether this method with TIP4P/2005 can work or not.

I’m thinking this is some sort of parallel computing problem, however I’m not sure what sort of things I could change/investigate to try and resolve it.
I would be really grateful for some help with this, or any wisdom/opinions/pointers on what might be happening.

Olivera

I am applying NEMD (SLLOD) to obtain shear viscosity of water at 293K (20°C) for water models TIP3P, SPC/E & TIP4P/2005.

The setup has worked fine for TIP3P and SPC/E: after applying unit conversions, viscosity values are in good range of the experimental value/literature (slightly more accurate than reported simulations, which are mostly at 300K).

However, for TIP4P/2005 I just can’t seem to make it work. The system setups are the same, the only differences are parameters specific to the water models.

‘’ERROR: Out of range atoms – cannot compute PPPM’’

comes up every time after PPPM Initialization; the simulation never actually gets to carrying out the equilibration run. The error refers to the kspace_style command line. I used pair_style ‘lj/cut/tip4p/long’ with kspace_style ‘pppm/tip4p’.

this is a known bug, please see: https://github.com/lammps/lammps/issues/869

axel.

Hi Axel,

Thank you for responding, so quickly too J It’s comforting to finally know what was causing the problem.

I’ll have to do some replanning - what month might the Spring release likely be, assuming the bug can be fixed?

Olivera​

Hi Axel,

Thank you for responding, so quickly too J It’s comforting to finally know what was causing the problem.

that is why we have moved LAMMPS development to GitHub, so the majority of the development process is recorded in public, and that includes plans for changes and documenting known bugs. We encourage all LAMMPS users to report all known issues in the issue tracker, so people don’t have to wonder about known issues. You can also find resolved (and invalid) issues in the list of closed issues.

I’ll have to do some replanning - what month might the Spring release likely be, assuming the bug can be fixed?

there is no guarantee whether this will happen, as the majority of LAMMPS development is done in people’s “free” time. while this specific issue is know, we don’t know the exact origin and without that, there is no path to resolve it. as you may see from the history of the issue, it has been moved from milestone to milestone. that assignment is primarily a record of the urgency.

Please note, that lj/cut/tip4p/long and pppm/tip4p are not the only way to simulate the TIP4P model. You can also model it as an explicit 4-site rigid model with fix rigid/small or one of its variants instead of fix shake plus time integration. This is less efficient, because you have more interacting pairs and it will require a reduction in the time step to avoid divergence of the time integration, but it can be used right away, and works with lj/cut/coul/long and pppm. Since LAMMPS checks for zero masses (as they would cause a division by zero for regular time integration), you have to use a tiny mass (say 1.e-10) for the massless site or comment out the check, since for rigid objects only the total mass matters.

axel.

Thank you, that’s good to know.

I had been considering an explicit 4-site TIP4P/2005. Admittedly I’m not sure if I can totally figure it out in practice, being a newbie with a lot to get done in very little time.

Also not sure how recommendable it would be down the line when it comes to viscosities of complex aqueous solutions; I’m thinking it may create complications/sub-optimal results. I might just leave my TIP4P/2005 efforts on the backburner for now, maybe contemplate SPC/E or possibly a different viscosity method.

Thanks again for your helpful responses (and Happy New Year :slight_smile: )

Olivera

Dear Olivera,

Thank you, that’s good to know.

I had been considering an explicit 4-site TIP4P/2005. Admittedly I’m not sure if I can totally figure it out in practice, being a newbie with a lot to get done in very little time.

Also not sure how recommendable it would be down the line when it comes to viscosities of complex aqueous solutions; I’m thinking it may create complications/sub-optimal results. I might just leave my TIP4P/2005 efforts on the backburner for now, maybe contemplate SPC/E or possibly a different viscosity method.

Seeing that you qualify yourself as a newbie, I would indeed recommend switching to viscosity methods that have been used in published work on the viscosity of TIP4P/2005 water and don’t require triclinic boxes. Here are a few examples:

https://aip.scitation.org/doi/abs/10.1063/1.3697977

http://iopscience.iop.org/article/10.1088/0953-8984/24/28/284117/meta
https://aip.scitation.org/doi/abs/10.1063/1.3330544
https://aip.scitation.org/doi/10.1063/1.5042209
https://aip.scitation.org/doi/10.1063/1.3515262
http://advances.sciencemag.org/content/3/8/e1700399
https://pubs.rsc.org/en/Content/ArticleLanding/2017/CP/C6CP07863J#!divAbstract

Best regards,

Laurent

Hello Laurent,
It’s very kind of you to have put together this reading list - I have come across a few of these but will look into the others - I’ve been reading a great deal since starting a couple of months ago.

I’ll be looking into non-Newtonian fluids down the line, hence the non-equilibrium approach.

The main query I had has been clarified by Axel, but thanks again for the recommendations.

All the best,

Olivera