GROMACS version:2023
GROMACS modification: Yes/No
Here post your question
Hi all,
This is somewhat of a continuous question from Confusing Density from NPT run
Relevant information for this question here
I am running MD on a protein-ligand complex.
My protein is huge with 6 chains and 277.25 kDa.
My forcefield parameterisation was using AmberForce field FF19SB / ambertools. I have decided to use OPC water model to solvate my protein-ligand complex
For any amber users this is how i solvated
solvateoct mol OPCBOX 10.0
addions2 mol Na+ 173
addions2 mol Cl- 76
Now i have two scenarios
Scenario A)
I solvated without crystallography water molecules. My hope was that the water molecules will diffuse into the complex/chains and establish the H-bonding network unfortunately that does seem to be the case. In some ways, the water molecules are compressing the protein imo because they seem to be on the outskirts of the protein and my density appears to be go from ~980 and shoots up to ~1300 during restrained NPT.
Scenario B)
I decided to parameterise with all of the water molecules from the crystallography pdb and tried to energy minimise but got error
Steepest Descents:
Tolerance (Fmax) = 1.00000e+03
Number of steps = 50000
Step= 0, Dmax= 1.0e-01 nm, Epot= 1.54311e+18 Fmax= inf, atom= 222645
Step= 1, Dmax= 1.0e-01 nm, Epot= 1.54311e+18 Fmax= inf, atom= 222645
Step= 18, Dmax= 1.8e-06 nm, Epot= 1.54313e+18 Fmax= inf, atom= 222645
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 19 steps,
but did not reach the requested Fmax < 1000.
Potential Energy = 1.5431128e+18
Maximum force = inf on atom 222645
Norm of force = inf
i have checked the atom number and this pertains to the water molecules on the outskirts of the protein-ligand complex. I was going to attempt removing selected water molecules that may cause this error however, i saw that maybe this is not a good choice Deleting problematic water molecules from .gro files.
EDIT: I forgot to mention that i did check the water molecules and they do not seem to overlapping with other water molecules. I check the output trajectory and i dont see any terrible distortions…
I wonder if this is due to solvating on top of crystallography water molecules and wonder if anyone could suggest a “safe” way to solvate or have i done something wrong in preparing my protein for MD(?).
Could i have only taken particular water molecules and prevent steric clashes? but the issue may still rise where the water molecules do not seem to enter into the protein chains…
Appreciate any advice :)
em.mdp (tried without the define part too)
title = Minimisation
; Parameters describing what to do, when to stop and what to save
define = -DPOSRES
integrator = steep ; Algorithm (steep = steepest descent minimization)
emtol = 1000.0 ; Stop minimisation when the maximum force < 10.0 kJ/mol
emstep = 0.1 ; energy step size
nsteps = 50000 ; Maximum number of (minimization) steps to perform
nstxout = 10
; Parameters describing how to find the neighbors of each atom and how to calculate the interactions
nstlist = 1 ; Frequency to update the neighbour list and long range forces
cutoff-scheme = Verlet
rlist = 1.0 ; Cut-off for making neighbour list (short range forces)
coulombtype = PME ; Treatment of long range electrostatic interactions
rcoulomb = 1.0 ; long range electrostatic cut-off (nm)
vdw-type = Cut-off ; Mdp parameters for non-bonded interaction using the AMBER force field - #4 by dbischoff
vdw-modifier = Potential-shift
rvdw = 1.0 ; long range Van der Waals cut-off (nm)
pbc = xyz ; Periodic Boundary Conditions
DispCorr = EnerPres ; For OPC water model and AmberFF19SB