Hello, I’m still fairly new to GROMACS but have been managing to make runs of my small guest-host systems (as preparation for RBFE) in water. Apart from some minimisation reluctance, they seemed to work fine.
Now, I’ve tried using a setup with the chloroform (one from Gromacs User Contributions) and for some complexes, it stubbornly claims:
ERROR 1 [file em_rev.mdp]:
The largest distance between non-perturbed excluded atoms is 1.693 nm,
which is larger than the cut-off distance. This will lead to missing
long-range corrections in the forces and energies.
This distance is almost as large as my host structure. It is also overall almost completely rigid, hence I would not expect these ends to interact anyway.
I wouldn’t have thought it is the solvent’s fault, but it is the only difference between the runs.
Otherwise I was using default mdp settings:
; Electrostatics
coulombtype = PME
rcoulomb = 1.2
; van der Waals
vdwtype = cutoff
vdw-modifier = force-switch
rvdw-switch = 1.0
rvdw = 1.2
My guess is that something has gone wrong with periodic boundary conditions and a chloroform molecule is broken over PBC in the structure file used for solvating your system. If that is the case, you need to make it whole again or use another structure file with whole molecules.
Hm, I was under the impression that the chcl3_box indeed had whole molecules and therefore simulations with small molecules (in this case: glycine anhydride, 14 atoms) actually work without this error message. The problem appears for the slightly larger (146 atoms total) system.
In either case, are there any other sources of chloroform solvation models that do not have casually breaking molecules? I was under the impression that the solvation algorithm can only place whole molecules anyway. Is this not so?
I now looked at chcl3_box.gro and there are no PBC issues there. So then I guess the issue is in your host structure. Unfortunately grompp does not print the atom indices in the error message. We should add these.
Have you visualized your host structure in the file you used as input to grompp?
Alright, I had another look at the problem and can still not conclude where the error originates.
Guest+solvent, Host+solvent, and solvent-only simulations work fine; when translating guest relative to host distance changes but there is no possible atom pair that could have such a distance. For the original locations, the largest distance between guest and host atoms is 1.243 nm; 1.719 nm appears but is between two host atoms where it clearly was not a problem.
That said, I also went back to look through the complex input pdb which had a slight discrepancy in that the element type was missing for the guest structure. With that amended PDB file, the calculation was not crashing the “gmx grompp” command to create the .tpr file but instead the execution of it with a more helpful suggestion:
Fatal error:
There are perturbed non-bonded pair interactions beyond the pair-list cutoff
of 1.2 nm, which is not supported. This can happen because the system is
unstable or because intra-molecular interactions at long distances are
excluded. If the latter is the case, you can try to increase nstlist or rlist
to avoid this.The error is likely triggered by the use of couple-intramol=no
and the maximal distance in the decoupled molecule exceeding rlist.
Where only changing the nstlist and rlist made the energy minimization proceed after all.
I was still hesitant to change them, without fully understanding the possible downstream effects of this, or perhaps they’re just some sampling variables?
The system setup, as I previously noticed is also rather sensitive to the LINCS restraints which is why I would run multiple minimisation steps: the first few with a steep integrator, then moving onto a stochastic descent one. With this approach, the error does not show up at all, apparently being fixed somewhere in the steepest descent runs.
Do you perhaps have an idea why that actually worked?
It could be that your initial structure has overlapping atoms which causes excluded atoms to move beyond rlist. Or there could be an issue in the topology.
That is what I concluded, but the recent runs has the original rlist and nstlist value and worked fine also.
Do you know what downstream effect they could have in the simulation and how sensitive it is to changes in them?