lammps got stuck when keyword shake used in fix gcmc combined with fix nvt

Hi, I’m using fix gcmc to fill an empty box with spc water.

Here are some lines that use fix shake and shake keyword in fix gcmc , combined with fix nvt to simulate the grand canonical ensemble.

fix wshake water shake 0.0001 50 0 b 1 a 1 mol h2o

variable temp equal 300.0

fix 1 water gcmc 1000 1000 1000 0 2949499 ${temp} -0.48 0.5 mol h2o shake wshake maxangle 180 #with keyword shake wshake

fix 2 water nvt temp {temp} {temp} 100

fix_modify 2 dynamic/dof yes

thermos 100

However, Lammps did not output the results for a long time, it got stuck:

Setting up Verlet run …

Unit style : real

Current step : 0

Time step : 1

Per MPI rank memory allocation (min/avg/max) = -5302 | -5302 | -5302 Mbytes

Step Atoms Temp Press Density E_vdwl E_coul c_T4

0 0 0 0 0 0 0 0

—No output in a time of minutes!

Instead, when the shake keyword is excluded from the fix gcmc:

fix 1 water gcmc 1000 1000 1000 0 2949499 ${temp} -0.48 0.5 mol h2o maxangle 180 # shake wshake

the Lammps run normally and immediately output the results!

I’m confused about the reason, does the keyword shake in fix gcmc cannot be used with fix nvt?

Reply to Andrew Jewett:

  1. I have simplified and attached my LAMMPS files. During the checking process, I still cannot solve this problem, i.e. the shake keyword in fix gcmc combined with fix nvt causes LAMMPS to stuck to output results.

  2. My initial intention to use fix gcmc to fill an empty box of water is to learn and test gcmc through the use of water chemical. I expected to fill the empty box with water of about 0.9 g/ml. My final simulation is to get the adsorption of water in clay minerals.

Also, I have checked the software tools suggested by Andrew to create water files and the moltemplate example files, these are useful for me.

in.water-gcmc.txt (1.62 KB)

mol.water.txt (634 Bytes)

Hi, I’m using fix gcmc to fill an empty box with spc water.

​there can be all kinds of strange problems when starting from an empty
box. things will be much easier to manage for LAMMPS, if you start with a
system, that contains at least one water molecule.

axel.

I posted a similar question on fix gcmc as I’m too trying to adsorb water into a porous substrate and having issues with it. I’m using 31 March 2017 version of lammps. The current version has examples on gcmc but that’s for a single system of CO2 where the bond types and angle types are similar in the host system and the new molecules. For such a system I have no trouble running gcmc and it works fine for just water molecules in my case (as Axel mentioned with preexisting water in the box). However, when I try to adsorb water in a different system (for example Ca-silicate) then lammps gives ERROR: Molecule auto special bond generation overflow. For the host system I have atom types 1-4 with no bonds or angles. For the incoming SPC water molecules I have bonds and angles with types 5 and 6 for O and H of water respectively. I kept the provision for these extra types in the substrate data file. But I’m still seeing this error. Is there something else that need to be mentioned here? I am copying my inputfile with the water molecule and the data file of the host system that I’m reading with read_data. Please let me know what I’m doing wrong about this. I have very little experience with fix gcmc. I tried looking up for similar questions in the lammps user archive but couldn’t find a solution for my problem. So now I definitely need help on this from the lammps users who have experience with fix gcmc. Any example for adsorbing water or some other molecule into a different substrate would be also fine. Thanks in advance.

LAMMPS INPUT DECK

units metal
atom_style full
bond_style harmonic
dihedral_style none
improper_style none
angle_style harmonic

pair_style lj/cut 12.00

read_data data.anhydrous

molecule watermol mol.water

group water type 5 6

mass 1 28.065 #Si

mass 2 16.00 #O
mass 3 40.078 #Ca1
mass 4 40.078 #Ca2
mass 5 16.00 #Ow
mass 6 1.00 #Hw

set type 1 charge 1.72
set type 2 charge -1.1608219
set type 3 charge 1.43
set type 4 charge 1.70
set type 5 charge -0.82
set type 6 charge 0.41

pair_coeff * * 0.000000 0.100000 0.0000000
pair_coeff 1 2 0.2433e-4 3.6716
pair_coeff 2 2 0.5392e-1 3.0687
pair_coeff 2 3 0.3773e-4 4.898
pair_coeff 2 4 0.6331e-4 4.898
pair_coeff 1 5 0.2299e-4 3.6298
pair_coeff 3 5 0.3802e-4 4.898
pair_coeff 4 5 0.2617e-4 5.0168
pair_coeff 5 5 6.74E-03 3.5532 #from clayff
pair_coeff 2 5 0.2278e-3 4.7557

#H2O SPC/Fw
bond_coeff 1 45.93 1.012
angle_coeff 1 3.29136 113.24

kspace_style none

neighbor 1.0 bin
neigh_modify delay 4 every 2 check yes

timestep 0.001
dump id all atom 1000 dump.${name}
thermo 10
thermo_style custom step pe ke temp press density vol atoms

variable tfac equal 5.0/3.0

fix g1 all gcmc 10 100 100 0 27868 300 0.1 0.01 mol watermol tfac_insert ${tfac} group water #rigid m1 # maxangle 180 full_energy

fix r3 all nvt temp 300.0 300.0 1.0
run 1000

DATA FILE (data.anhydrous)

LAMMPS data file written using python script

70 atoms
6 atom types
0 bonds
1 bond types
0 angles
1 angle types

-0.849377 28.998874 xlo xhi
-0.018042 29.630416 ylo yhi
-1.064126 29.490964 zlo zhi

Masses

1 26.980000
2 16.000000
3 40.078000
4 40.078000
5 16.000000
6 1.000000

Atoms

1 0 1 0.000000 11.944669 5.947721 18.383305
2 0 2 0.000000 13.540550 6.356191 19.344793
3 0 2 0.000000 12.390029 5.412430 16.601149
4 0 2 0.000000 10.828668 7.498327 18.311030
5 0 2 0.000000 11.015977 4.531274 19.260118
6 0 1 0.000000 26.455698 5.114050 9.789356
7 0 2 0.000000 24.884309 4.843987 8.742483
8 0 2 0.000000 26.518044 3.819125 11.196568
9 0 2 0.000000 26.420799 6.866727 10.552225
10 0 2 0.000000 27.998874 4.927879 8.684124
11 0 1 0.000000 10.155061 15.784010 6.449462
12 0 2 0.000000 9.578487 17.072908 7.731821
13 0 2 0.000000 8.630115 14.982098 5.617330
14 0 2 0.000000 11.177258 14.426079 7.324791

MOLECULE FILE

A single water molecule file

3 atoms
2 bonds
1 angles

Masses

1 1.000000
2 1.000000
3 16.00000

Coords

1 0.000000 0.000000 0.000000
2 0 0.000000 1.633045 0.000000
3 0.577415 0.816522 0.000000

Types

1 6
2 6
3 5

Charges

1 0.41
2 0.41
3 -0.82

Bonds

1 1 1 3
2 1 2 3

Angles

1 1 1 3 2

The ERROR: Molecule auto special bond generation overflow described by Mohammad, may be related with the read_data command and the data file.

The LAMMPS manual gives a detail description of read_data command and the format of data file.

In this situation, there should be some “extra bond per atom” etc in the data file as I quoted from the LAMMPS Users Manual 13 Apr 2017 version:

· extra bond per atom = leave space for this many new bonds per atom

· extra angle per atom = leave space for this many new angles per atom

· extra dihedral per atom = leave space for this many new dihedrals per atom

· extra improper per atom = leave space for this many new impropers per atom

· extra special per atom = leave space for this many new special bonds per atom

I have not tested Mohammad’s system, if it does not work, please inform me.

Hi Youzhi,

Thanks for your reply. I tried using extra/bond/types and extra/angle/types arg with read_data but then I get “ERROR: Molecule topology/atom exceeds system topology atom”. As I’m reading the data and not using create_box its not possible for me to use extra bond per atom or extra angle per atom as shown in the example for CO2. I need to use read_data for the substrate as the system under study is amorphous and cant be made using lattice command. In the substrate I don’t have any bonds unlike the adsorbing SPC water molecules where I have both bonds and angles. I have attached my input script along with the molecule file and the substrate data file in case someone wants to test it. Thanks again for the help.

Best regards,

data.anhydrous (3.54 KB)

mol.water (338 Bytes)

in.gcmc (1.62 KB)

Hi Mohammad,

The ERROR: “Molecule topology/atom exceeds system topology atom", you may need to add this:

extra special per atom = leave space for this many new special bonds per atom

in your data file. Please try this and inform me if this error occurs.

Hi Mohammad,

The extra section should be added in the header section of the data file when using read_data command:

· extra bond per atom = leave space for this many new bonds per atom

· extra angle per atom = leave space for this many new angles per atom

· extra dihedral per atom = leave space for this many new dihedrals per atom

· extra improper per atom = leave space for this many new impropers per atom

· extra special per atom = leave space for this many new special bonds per atom

Hi Youzhi,

I tried putting the following in the headers of my data file:

LAMMPS data written by python

70 atoms
6 atom types
0 bonds
1 bond types
0 angles
1 angle types

1 extra bond per atom
1 extra angle per atom
1 extra special per atom

It still gives me "ERROR: Molecule topology/atom exceeds system topology atom". I tried increasing the no. of extra bonds and angles but the result is still the same. However if I try increasing the extra special to too high then the simulation doesn’t proceed and seems to be stuck without any progress. The water molecules have just one bond types which is O-H and one angle type H-O-H and they are already specified in the molecule file. The error is generated right after reading the fix gcmc command (I used a simple print command to confirm that). So it has to be something with the fix gcmc and the way its trying to add the molecules into the substrate. I did add the special bond counts and special bonds section in my molecule file as shown below. Please let me know if I’m doing something wrong. Thanks again for the help.

A single water molecule file

3 atoms
2 bonds
1 angles

Masses

1 1.000000
2 1.000000
3 16.00000

Coords

1 0.000000 0.000000 0.000000
2 0 0.000000 1.633045 0.000000
3 0.577415 0.816522 0.000000

Types

1 6
2 6
3 5

Charges

1 0.41
2 0.41
3 -0.82

Bonds

1 1 1 3
2 1 2 3

Angles

1 1 1 3 2

Special Bond Counts

1 1 1 0
2 1 1 0
3 2 0 0

Special Bonds

1 3 2
2 3 1
3 1 2

Hi Mohammad,

The "ERROR: Molecule topology/atom exceeds system topology atom" still relates with the “extra” in the data file.

In your system, for water molecule, there are:

2 extra bond per atom

1 extra angle per atom

2 extra special per atom

I tested this correction and found it work on my machine (Attached as anhydrous+water.tar.gz ).

anhydrous+water.tar.gz (51.8 KB)

(1) Thanks to Axel Kohlmeyer’s comment, I tried to add some water molecules in the intial box,

then LAMMPS works properly with shake keyword used in fix gcmc.

So it was a bad idea to start gcmc using shake in an empty box.

(2) Above situation also suits the simulation box filled with clay pores, which have no water, then after add some waters to the system,
the fix gcmc with shake starts works properly.
So, it seems to me that the fix gcmc with shake keyword need at least one shake molecule in the box to work.

(3) The rigid keyword in fix gcmc always works (even in empty box), and I’m confused about the rigid and the shake used with water.

If I want to simulate water adsorption on clay minerals, should I use fix shake of the water molecules, or just use the rigid water molecule ?

Hi Youzhi,

Thanks for sending me the output. Very strange, after implementing the changes it still does’t work for my version of lammps. Which version are you using? I’m using the pre-built Ubuntu Linux executable lammps-daily and I have updated it as well. At the beginning of the run it says “19 May 2017”. Could it be something due to the pre-built lammps version? It just gets stuck after saying WARNING: Molecule attributes do not match system attributes.

UPDATE:

Youzhi, I ran your script and they seem to work fine. So there must be something wrong in my script (still cant find the reason). Anyways, this is great but I tried to extend my run with preexisting water molecules in the substrate data file and it gives me segmentation fault. I did add the extra bonds/angles/special same as before. Please let me know if you have any insight on this. I have attached the data file with pre-existing water molecules. Could it be related to some bug in fix gcmc? Please let me know your suggestion. I’m still using 19 May 2017 version. Thanks!

data.silicate_5-water (2.9 KB)

Hi Mohammad,

(1) I’m using lammps version: lammps-stable-31Mar2017.

(2) The “WARNING: Molecule attributes do not match system attributes” is suspected to be related with the duplicate definitions of atom properties both in the input file and the water molecule file (in this case, the atom mass). By removing the “Masses” in the water molecule file, this warning vanished.

(3) In your data.silicate_5-water file, one of the water molecule ID (line 63-65) is the same with that of silicate which causes molecule ID conflict.
By modifying this water molecule ID from “1” (same with the silicate molecule ID “1” ) to a different value, such as “6”, which is different from all the other molecule IDs, this “segmentation fault error” is avoided.

Hi Youzhi,

Its working perfect now. Thanks a lot for staying with me this long and help me solve all the problems.

Best regards,

Mohammad

Hi Mohammad,

Thanks is mutual between us! I have learnt a lot from discussions with you, and especially, the “set style ID charge command” to assign atom type charges through the input file conveniently instead of setting each atom in the data file.