Hi LAMMPS list,
I was trying this “fix gcmc” with the “fix nvt” to sample the system with the MD integrator, here is what it looks like in the simplest implementation:
compute_modify thermo_temp dynamic yes
timestep 2.5
run_style verlet
dump 1 all xyz 10000 test-gcmc-*.xyz
fix 1 all nvt temp 303.0 303.0 10.0
fix 2 all gcmc 2000 1 0 1 29494 303.0 -1.5 0.1 molecule yes
run 1000000
So it works to insert and delete molecules, and I have a version that uses a molecule ID “mol” pointing to type 1 instead of “all” in the “fix gcmc”, but the system starts as one molecule and it seems to work the same either way.
I had played around with the Tdamp value but it seems like no matter what I do the system seems to lose kinetic energy proportionally to the number of molecules added to the system so maybe the NVT fix is not being applied to the new molecules. There is a bit of code in the fix_gcmc.cpp file that looks like this:
int nfix = modify->nfix;
Fix **fix = modify->fix;
for (int j = 0; j < nfix; j++) {
if (strcmp(modify->fix[j]->style,“shake”) == 0) {
fix[j]->update_arrays(m,atom_offset);
} else if (fix[j]->create_attribute) fix[j]->set_arrays(m);
}
I was wondering if there was possibly a way to set it to update the “fix nvt” in that section and if that might be the problem here? It seems to proceed fine with a syntax like this in the input file, as a test of the possible source of the problem:
label top
fix 1 all nvt temp 303.0 303.0 10.0
fix 2 all gcmc 2000 1 0 1 29494 303.0 -1.5 0.1 molecule yes
run 10000
jump SELF top
But it seems that this is just a stopgap measure and I am just comparing the source here to another GCMC implementation.
I’m not sure what’s happening in the GCMC algorithm. I see the temperature is initialized to the target 303K for new molecules in the fix_gcmc.c code, but I’m guessing the “fix nvt” is not being applied to the new molecules that are inserted and they slowly lose velocity because of no fix (so they are poorly equilibrated).
So I was next wondering if hybrid MD/MC can be applied to “reservoir” molecules in the GCMC so they are properly equilibrated when a molecule insertion is attempted? I have been working on a set of MC moves for angle change and dihedral change of molecules (it’s not necessary for gases and water molecules but it is needed for more complex molecules and cyclic molecules have their own MC moves too) and was wondering if there are any of these moves in LAMMPS? Using the hybrid MD seems like a good alternative to MC moves but these molecules are non-interacting (i.e., they only interact intramolecularly and wouldn’t be in the system yet).
So is there a way to reference a reservoir of non-interacting molecules for applying either hybrid MD/MC moves or regular MC moves in the GCMC algorithm to equilibrate molecules for insertion into the system? I was not able to find any examples in the LAMMPS source but I have made some source code examples that should fit into the GCMC paradigm.
Thanks,
Luke