[lammps-users] Multiple fix gcmc for a binary mixture

Dear All,

I’m trying to run a GCMC simulation of a binary mixture. I got very strange results even for a binary lennard-jones mixture. I then tested a minimum system for “fix gcmc”: a binary LJ mixture but with the two components having the same interaction, so effectively a pure fluid. This too gave strange results. In all cases I’m using LAMMPS v. (3 Mar 2020).

The state point is from the LAMMPS lj gcmc example “in.gcmc.lj”, and is T=2.0 and mu=-1.25. The variables in “fix gcmc” are mu=mu1=mu2=-1.25, Ttarget=2.0 and disp=1.0.
After 10^6 MC cycles the averages are:

Pure liquid (one gcmc fix):
fix gcmc1 1 gcmc 1 100 100 1 1818 ((v_Ttarget)) ((v_mu1)) $((v_disp) group 1

T = 2.0038342
P = 1.4356461
density = 0.53376733

Pure liquid but binary (two gcmc fixes):
fix gcmc1 1 gcmc 1 100 100 1 1818 ((v_Ttarget)) ((v_mu1)) $((v_disp)) group 1

fix gcmc2 2 gcmc 1 100 100 2 1818 ((v_Ttarget)) ((v_mu2)) $((v_disp)) group 2
T = 2.0036288
P = 1.6200632
density = 0.55261573
mole fraction = 0.86737019

For reference, a quick NPT MD simulation at T=2.0 and P=1.45 gives:
T = 2.0000051
P = 1.4500375
density = 0.53636195

It seems as though using two “fix gcmc” fixes gives a pressure and density that are slightly off, and a mole fraction very different from x1=x2=0.5.

What am I doing wrong? Is it possible to do a GCMC simulation of a binary mixture using LAMMPS? I’ve attached the input and log files for the two GCMC runs.

Any help would be much appreciated,

Kind Regards,

James

in.gcmc.ljbinary (2.93 KB)

ljbinary.log (191 KB)

in.gcmc.ljpure (2.93 KB)

ljpure.log (191 KB)

Each of the two gcmc fixes does not know what the other is doing.
Since they have the same settings and random number seeds, they will be executed on the steps and they are active one after the other.
So you are breaking the “balance” with the second always executing when the first may have already pushed the system to a high(er) potential energy.

Please also note the comments in https://matsci.org/t/how-to-perform-grand-canonical-monte-carlo-simulation-for-binary-mixture-by-lammps/28922/2
where something similar is attempted by one of the authors/maintainers of fix gcmc but there are significant differences in how the gcmc fixes are issued.

Thank you Axel for the quick response and helpful insight.

Shifting the chemical potential by kbTln(2), giving the two gcmc fixes different seeds and executing them on different timesteps gives the averages expected from the pure fluid.

GCMC settings

variable mushift equal -v_Ttarget*ln(2.0)
variable mu1 equal -1.25+v_mushift
variable mu2 equal -1.25+v_mushift
variable disp equal 1.0

fix gcmc1 1 gcmc 2 100 100 1 1818 ((v_Ttarget)) ((v_mu1)) ((v_disp)) group 1 run 1 fix gcmc2 2 gcmc 2 100 100 2 111 ((v_Ttarget)) ((v_mu2)) ((v_disp)) group 2

averages

T = 2.0048075
P = 1.4404807
density = 0.53450412
mole fraction = 0.50093096

I will now move onto a more complicated LJ system.

Kind regards,

James