Incorrect handling of “Distance between pull groups” with PBC
While doing some AWH scans of the ligand complexes with proteins and DNA, I sometimes get the errors at some unpredictable points in the simulation
Fatal error: Distance between pull groups 3 and 10 (5.204208 nm) is larger than 0.49 times the box size (5.203962).
In previous cases I was able to work around it by including an additional restraining potential that supposedly prevented the ligand from going too far. However, with the recent system (DNA 28-mer + organic ligand in water box) I was not so lucky.
Moreover, I cannot understand why the error occurs when it does. I also wonder what the correct options should be and/or if the logic in pull code might actually be incorrect.
So, the box is actually quite large (10.6 nm including 1 nm margin on all sides) and the flat-bottom restraint is set at 5 nm. However, the system never reaches that distance. In one run, it started at 1.9 nm and crashed soon after 3.35 nm. I guess that, due to some random drift of molecules, the pull code suddenly finds another periodic image of one or both molecules more “attractive”, leading to this discontinuity.
I tried playing with
pull-groupN-pbcatom, and COM motion removal options (even tried setting the ligand
pbcatom to another molecule’s reference) but to no avail. Increasing the box size obviously would not help either, as the desired distance is not reached anyway (but the simulation would take longer to crash)))).
However, while trying larger box (12.6 nm), I have luckily stumbled upon a seed combination that reproduces this crash quickly (see attachment pulldist.7z). Unfortunately, VMD crashes on the incomplete trajectory.
Step Time 841 1.68200 Program: gmx mdrun, version 2020.2 Source file: src/gromacs/pulling/pull.cpp (line 497) Fatal error: Distance between pull groups 3 and 10 (6.173351 nm) is larger than 0.49 times the box size (6.172868).
The starting value for this distance (t=0) is 1.951 and the last recorded
pullx.xvg state (t=1.5420) is 1.93284, so obviously it could not possibly have reached 6.17 at t=1.682 in normal way.