I am trying to simulate the adsorption of ammonia onto zeolite. I am using the Vashita pair style for the silicon oxide, and lj/cut/coul/long for the ammonia interactions. I have defined bonds and angles within an ammonium file and have them as harmonic. Essentially I am trying to model the ammonia according to a TraPPE potential which combines bonds angles and lj/coul methods. I eventually want to run a Monte Carlo simulation, but first I want to make sure that I can define a stable environment with the silicon dioxide and the ammonium.
I have run into two main issues:
First that the system contains non-numeric atomic coordinates if I use a delete atoms command with a cutoff equal or greater than 0.8.
If I leave the cutoff great enough, I will run into the “molecule auto special bond generation error”. Reading the documentation, I see that it is uncommon to use to use many-body potentials with explicit bonds, so I set the weighting factors to one. I am not sure if this will be a mute point or not, because I use a fix rigid/nvt command on the ammonia molecules to behave as a rigid body since I am only interested in its interaction with the silicon dioxide.
The silicon dioxide data file contains the initial position and velocity of the silicon dioxide while I just read in an ammonium molecule file.
Any advice or insights would be greatly appreciated.
NH3.data (428 Bytes)
dilatedSiO.data (79.1 KB)
adding_water5.lmp (2.5 KB)
There are a bunch of problems with this just from a cursory look. There are likely more upon closer inspection:
- why use " lj/long/coul/long cut off 2.5"?
- this is the same as lj/cut 2.5
- a cutoff of 2.5 angstrom seems very wrong. it is shorter than your sigma of 3.6 angstrom
- your model of ammonia has charges those are ignored here, but you are using pppm at the same
time. this makes no sense.
- your pair_coeff settings with an epsilon of 1 eV seem very extreme. Where did you get this parameter and the choice of cutoff?
- you are using
pair_coeff * * so you are applying the same LJ parameters to all atom types and for the atoms that use the Vashishta potential this is added, too.
- why special_bonds lj/coul 1.0 1.0 1.0? You do have a molecular system with bond and angle interactions and for those the default setting (0.0 0.0 0.0) is correct. Otherwise you are double counting interactions. The warning can be ignored for as long as you don’t define bonds and angles for the SiO2 part (which you do).
- you have bond and angle interactions added to the Vashishta potential. That is not correct.
- have you checked whether this Vashishta potential is suitable for interfacial systems. Most likely not. The comments in the potential file suggest it has only been tested for bulk systems and that is what makes sense. Interfaces are complex, typically require some kind of charge equilibration and SiO2 will also have exposed oxygen atoms that may have open valences saturated with hydrogen atoms.
With such a bogus input deck all kinds of issues can happen. The problem with delete_atoms is a minor issue in comparison. You need to fix all the other problems first and then I would make some tests where you do not use random placement of ammonia molecules, but a simple fcc or bcc lattice for a sufficient lattice parameter and in a region with sufficient distance to the interface so that there will be no overlaps and then you should see whether you can set up a stable simulation then.