I have been trying to minimize the energy of a coarse-grained MARTINI ion channel- membrane system that I constructed with the CHARMM membrane GUI builder. I keep running into the following error even though I have increased rvdw/rcoulomb from 1.1 to 2 and increased nstlist from 20 to 40. I even tried increasing rvdw/rcoulomb from 1.1 to 5 and nstlist to 50, and I ended up with the same error message even though the system energy minimized to step 1103 with a potential energy of -5.48e6. I have read the other issues on this forum with error messages similar to this one, but it wasn’t clear to me what the solution is or what I should do to move past this error.
Program: gmx mdrun, version 2021.3
Source file: src/gromacs/gmxlib/nonbonded/nb_free_energy.cpp (line 854)
Fatal error:
There are 1 perturbed non-bonded pair interactions beyond the pair-list cutoff
of 2.1 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.
Does the simulation run when you’re not performing any kind of perturbation? It seems perhaps your system simply isn’t well equilibrated enough to start making transformations. Minimize and equilibrate using normal MD methods, then apply the free energy features.
;soft-core-minimization so that single precision GROMACS works here
; Free energy parameters
free-energy = yes
init-lambda = 0.01
sc-alpha = 4
sc-power = 2
sc-coul = yes
nstdhdl = 0
couple-moltype = system
; we are changing both the vdw and the charge. In the initial state, both are on
couple-lambda0 = vdw-q
; in the final state, both are off.
couple-lambda1 = none
couple-intramol = yes
Also because the system is around 30 by 30 by 20 nm, does it still make sense to continue increase the rvdw and rcoulomb at all beyond 10 nm? I am not sure exactly what the limit of the cutoffs should be.
Well, you are by virtue of setting a non-zero value of λ, which I presume is being used to soften interactions and help resolve clashes. This is a strategy I would generally only recommend for problematic systems, but in your case it seems to be causing problems. I would try a normal minimization and equilibration protocol before trying to do anything unorthodox.
You should not mess with cutoffs; they are (generally considered) part of the force field.
There are other options that govern DD cell size and such that could be tweaked in special cases, but again I wouldn’t be messing with free energy options here.
Thank you. I was using the default CHARMM-GUI file, and then the only things I initially tried to vary were rcoulomb and rvdw since it said in the README file “step6.0 - soft-core minimization. If you encountered ‘There are 1 perturbed non-bonded pair interaction …’ error message, please modify rvdw and rcoulomb values from 1.1 to 2.0 in the step6.0_minimization.mdp file.” I have modified the mdp file to remove the free energy section, but now I am running into a different error regarding infinite forces.
Output:
Steepest Descents:
Tolerance (Fmax) = 1.00000e+01
Number of steps = 5000
Energy minimization has stopped because the force on at least one atom is not
finite. This usually means atoms are overlapping. Modify the input
coordinates to remove atom overlap or use soft-core potentials with the free
energy code to avoid infinite forces.
You could also be lucky that switching to double precision is sufficient to
obtain finite forces.
writing lowest energy coordinates.
Steepest Descents converged to machine precision in 15 steps,
but did not reach the requested Fmax < 10.
Potential Energy = 8.7619801e+23
Maximum force = inf on atom 15123
Norm of force = inf
It seems your system has some very nasty contacts. Try minimizing with the free energy options enabled, then minimize those coordinates again without them to see if you need to continue using the softened potential for any future steps.
Thank you! I will try that. This particular ion channel is known to induce a membrane indentation/dome around it at 1 bar pressure, but since I can only initialize the lipids in a planar configuration, I am expecting that through minimization the lipids will start to fill in around the blades of the ion channel to form the dome. However, I am starting to think that I should instead take the ion channel pdb file and relax the protein structure in a vacuum until it is in its planar configuration (by adding negative pressure in the lateral dimension to stretch it out) before building the membrane around it. If I energy minimize and equilibrate a protein on its own in a vacuum, is is possible to add negative pressure in the lateral dimension without completely deforming the channel?
Hi @kyrstyno and @jalemkul, I’m working on Piezo1 ion channel embedded in membrane system and encountered the same problem as @kyrstyno.
“Fatal error:
There are 1 perturbed non-bonded pair interactions beyond the pair-list cutoff
of 2.1 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.”
My .mdp file is similar to @kyrstyno as well because I also used CHARMM-GUI for building the system.
As suggested, I tried couple-intramol=yes, increase rlist and nstlist, but still get the same error.
When I tried minimizing without free energy options, or minimizing with the free energy options enabled then minimize those coordinates again without them, I also got the same error as @kyrstyno regarding infinite forces.
I’m not sure what’s happening here and what should I do to get this fix and get my minimization running.
I really appreciate any advice/suggestions.
@kyrstyno: Have you by any chance fixed this error?
I’m sorry for any confusion this may cause. I came across this thread from @kyrstyno with the same system with the exact same errors, so I just wanted to ask if he/she has figured out the solution, which would be applicable for my system as well.
