I am checking the fix gcmc feature on a very simple problem, before I try it in a more complicated setup.
My simulation is rather simple: a cubic box of lj atoms, with pair coefficients eps=1.0 and sigma=1.0 (all in lj units), and cut-off 2.5. If I run an NPT simulation of that system at pressure 1.0 and temperature 1.0 (lj units), the resulting equilibrium density is 0.8.
Now, if I run an NVT simulation of the same system at T=1.0, starting at a lower density (let's say 0.5) and coupling my system with a reservoir at T=1.0 and P=1.0 through the command:
fix 3 all gcmc 100 100 500 1 4839843 1.0 1.0 0.01 pressure 1.0
The density of my system can go up to 0.84 (and even higher). From the top of my head, I can't remember what is the jamming density for a LJ fluid, but I must be getting closer. How is that possible? Am I missing something? Is there a bug in the insertion/deletion algorithm? Have I used the wrong P in the GCMC reservoir?
Any help will be highly appreciated.