High pressure in simulation of carbon nano-particle inside water.

hi everyone
I am trying to simulate 2 layer graphene sheet(15nmX15nm) inside water box(30nmX30nm) (using lammps version released on 17 Aug 2017). Modelling is done using vmd and merging between water and graphene sheet following the axels tutorial(https://sites.google.com/site/akohlmey/software/topotools/topotools-tutorial—various-tips-tricks).
The interaction potential parameters are taken from references as shown in system.in.settings file. I want to simulate the system at lab conditions of (p = 1 atm and t = 298K), so I use npt after equilibrating the system for 5000 time steps of size 2 femtoseconds. At the end of equilibration pressure is of 5 orders and when npt is used in next phase the size of the box start to increase. I tried following methods:

  1. Graphene sheet is fixed using rigid single, and langevin thermostat along with nve/limit 0.1 on water for equilibration and then for next phase npt with dilate all option.
  2. Graphene sheet is fixed using rigid/nve single, and langevin thermostat along with nve/limit 0.1 on water for equilibration and then for next phase npt with dilate all option.
  3. Graphene sheet is fixed using rigid/nvt single, and langevin thermostat along with nve/limit 0.1 on water for equilibration and then for next phase npt with dilate all option.
  4. Then in second set, above three cases along with rigid I’ve also used neigh_modify to exclude the interaction between C atoms of graphene sheet and in third set, apart from exclude pair interactions delete_bond multi command for C atoms.
    Box expansion was less in cases when only rigid was used and no pair exclusion and delete bond instruction.

I also simulated water box without graphene sheet and found that pressure is too high for the equilibration phase involving langevin and nve/limit, but comes down to 1 atm.
Please go through my input files and suggest where is the glitch.

I came across similar querry on list, where reply was that the reason could be ill geometry or overlapping. if so, how could I figure out this? and solve it.

contents of system.in.init:

units real
dimension 3
boundary p p p
atom_style full

(Hybrid force fields were not necessary but are used for portability.

pair_style hybrid/overlay lj/cut/coul/long 10.0 lj/cut 10.0 # hybrid lj/charmm/coul/long 9.0 10.0 10.0 lj/cut/coul/long 10.0 #hybrid lj/charmm/coul/long 9.0 10.0 10.0
pair_modify tail yes
bond_style hybrid harmonic
angle_style hybrid harmonic
#dihedral_style charmm
kspace_style pppm 0.00001


lj parameters are taken from

Werder, T., Walther, J.H., Jaffe, R.L., Halicioglu, T. and Koumoutsakos, P., 2003. On the water− carbon interaction for use in molecular dynamics simulations of graphite and carbon nanotubes. The Journal of Physical Chemistry B, 107(6), pp.1345-1352.

We can not go through every step line by line. Generally, if you are trying to adjust volume of cell filled with rigid water and stiff graphene, it is a good idea to not rescale the positons of the graphene atoms. The reason is that the rescaling can result in large changes in bond length, and large forces on the atoms. To achieve this, use “dilate water”. But there can be issues with this too. In order to make progress, I suggest simplifying your script as much as possible until you get things working. Also, start out by doing pure water. If you can’t do that, why bother trying to add graphene. Conversely, if the graphene case is working well, it should behave similarly to the pure water case.


Thanks Thompson

I was successful in simulating water box using fix npt temp 298.0 298.0 100 iso 1.0 1.0 1000 and confirmed the density of water with literature though there were pressure fluctuations, unavoidable for liquids as per previous threads. For combined system I used dilate water option but to no avail. Finally, I even tried simulating graphene sheet, without water, with fix rigid/nvt and again the box size increased to a large value with large fluctuations in temperature. Apart from this I have one doubt that when using rigid on graphene and fix npt on water does the virial part of pressure takes into account the pair interactions between the graphene sheet atoms? If neede I can upload the data file.