fix shake fails when the full-energy option is turned on in fix gcmc

Dear Aidan,

I am currently inserting water molecules into a box with the help of both fix and fix shake. I find that if I turn off the full_energy option in fix , the two fixes work well. However, if I turn on the full_energy option, I always get an error indicating the shake determinant = 0.0, which is really confusing for me.
You can find my test case in the attachment.

Best,
Mingyang

shaketest.in (1.35 KB)

watercell2.molecule (573 Bytes)

Thanks for identifying this problem. I was able to reproduce it. SHAKE is not being implemented correctly for molecules inserted with full_energy. I have not yet figured out why that is.

We have identified the problem. It was a subtle interaction between fix gcmc and fix shake, only when energy_full is used. The solution is:

— MC/fix_gcmc.cpp (old)
+++ MC/fix_gcmc.cpp (new)
@@ -2108,6 +2108,10 @@
int eflag = 1;
int vflag = 0;

  • // force_clear();
  • int nbytes = sizeof(double) * (atom->nlocal + atom->nghost);
  • memset(&atom->f[0][0],0,3*nbytes);