[gmx-users] EQUILIBRATION

[‪farial tavakoli‬ ‪]

HiI am a new GROMACS user and trying to run md on HDAC2 with PDB ID: 4LY1.  but when I wanted to performance equilibration , i got this ERROR:                      step 1: One or more water molecules can not be settled.
                      Check for bad contacts and/or reduce the timestep if appropriate

Would you please help me to fix this problem?
These are all stages which I did:
I took 4LY1.pdb file from PDB site and removed B and C chains and all ligands and water molecules from it and left just chain A and HETATM ZN A . 
I constructed my protein topology by using gromos 43a1 and this command:pdb2gmx -f 4LY1.pdb -o 4LY1.gro -water spc
then made a complex.gro file which include 4LY1.gro and DRG.gro which i got from PRODRG site. I also copy paste ; Include ligand topology
#include "DRG.itp"

in the topology file and added DRG as a ligand in [ molecules ] directive. 
Then I defined box by using this command :
gmx editconf -f complex.gro -o newbox.gro -bt dodecahedron -d 1.0

then defined solvate by :
gmx solvate -cp newbox.gro -cs spc216.gro -p topol.top -o solv.gro

then used grommp to assemble a .tpr file by using a .mdp file and this command :
gmx grompp -f em.mdp -c solv.gro -p topol.top -o ions.tpr

that was my .mdp afile :
integrator	= steep		; Algorithm (steep = steepest descent minimization)
emtol		= 1000.0  	; Stop minimization when the maximum force < 10.0 kJ/mol
emstep          = 0.01          ; Energy step size
nsteps		= 50000	  	; Maximum number of (minimization) steps to perform
energygrps	= system	; Which energy group(s) to write to disk

; Parameters describing how to find the neighbors of each atom and how to calculate the interactions
nstlist		    = 1		        ; Frequency to update the neighbor list and long range forces
cutoff-scheme       = Verlet
ns_type		    = grid		; Method to determine neighbor list (simple, grid)
rlist		    = 1.0		; Cut-off for making neighbor list (short range forces)
coulombtype	    = PME		; Treatment of long range electrostatic interactions
rcoulomb	    = 1.0		; long range electrostatic cut-off
rvdw		    = 1.0		; long range Van der Waals cut-off
pbc                 = xyz 		; Periodic Boundary Conditions

and added 2 CL ions since my total charge was 2.00 because of zn ion. 

gmx genion -s ions.tpr -o solv_ions.gro -p topol.top -pname NA -nname CL -nn 2

then for minimization, i replaced energygroups from system to protein DRG in em_real.mdp. that was my em_real.mdp file:

integrator	= steep		; Algorithm (steep = steepest descent minimization)
emtol		= 1000.0  	; Stop minimization when the maximum force < 10.0 kJ/mol 
emstep          = 0.01          ; Energy step size
nsteps		= 50000	  	; Maximum number of (minimization) steps to perform
energygrps	= protein DRG   ; Which energy group(s) to write to disk

; Parameters describing how to find the neighbors of each atom and how to calculate the interactions
nstlist		    = 1		        ; Frequency to update the neighbor list and long range forces
cutoff-scheme       = Verlet
ns_type		    = grid		; Method to determine neighbor list (simple, grid)
rlist		    = 1.0		; Cut-off for making neighbor list (short range forces)
coulombtype	    = PME		; Treatment of long range electrostatic interactions
rcoulomb	    = 1.0		; long range electrostatic cut-off
rvdw		    = 1.0		; long range Van der Waals cut-off
pbc                 = xyz 		; Periodic Boundary Conditions

and typed these command: gmx grompp -f em_real.mdp -c solv_ions.gro -p topol.top -o em.tpr    and
gmx mdrun -v -deffnm em

GROMACS replies:
                Steepest Descents converged to Fmax < 1000 in 1833 steps
                Potential Energy  = -5.6356275e+05
                Maximum force     =  9.7815100e+02 on atom 3452
                Norm of force     =  3.1465845e+01

            Simulation ended prematurely, no performance report will be written.
Is there any problem that simulation ended prematurely?

Then for equilibriation i typed gmx genrestr -f jz4.gro -o posre_jz4.itp -fc 1000 1000 1000

and in order to restrain both the protein and the ligand simultaneously, I specify define = -DPOSRES -DPOSRES_LIG in the .mdp file. 
then typed   gmx make_ndx -f em.gro -o index.ndx

to merge the protein and DRG. and selected protein and DRG by :
> 1 | 14
> qthen set tc_grps = Protein_DRG Water_and_ions in NVT.mdp file to achieve my desired "Protein Non-Protein" effect. that was my NVT.mdp file:

title = Protein-ligand complex NVT equilibration define = -DPOSRES -DPOSRES_LIG ; position restrain the protein and ligand
integrator = md ; leap-frog integratornsteps = 50000 ; 2 * 50000 = 100 psdt = 0.002 ; 2 fs; Output controlnstxout = 500 ; save coordinates every 1.0 psnstvout = 500 ; save velocities every 1.0 psnstenergy = 500 ; save energies every 1.0 psnstlog = 500 ; update log file every 1.0 psenergygrps = Protein DRG; Bond parameterscontinuation = no ; first dynamics runconstraint_algorithm = lincs ; holonomic constraints constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrainedlincs_iter = 1 ; accuracy of LINCSlincs_order = 4 ; also related to accuracy; Neighborsearchingcutoff-scheme = Verletns_type = grid ; search neighboring grid cellsnstlist = 10 ; 20 fs, largely irrelevant with Verletrcoulomb = 1.4 ; short-range electrostatic cutoff (in nm)rvdw = 1.4 ; short-range van der Waals cutoff (in nm); Electrostaticscoulombtype = PME ; Particle Mesh Ewald for long-range electrostaticspme_order = 4 ; cubic interpolationfourierspacing = 0.16 ; grid spacing for FFT; Temperature couplingtcoupl = V-rescale ; modified Berendsen thermostattc-grps = Protein_DRG Water_and_ions ; two coupling groups - more accuratetau_t = 0.1 0.1  ; time constant, in psref_t = 300 300 ; reference temperature, one for each group, in K; Pressure couplingpcoupl = no ; no pressure coupling in NVT; Periodic boundary conditionspbc = xyz ; 3-D PBC; Dispersion correctionDispCorr = EnerPres ; account for cut-off vdW scheme; Velocity generationgen_vel = yes ; assign velocities from Maxwell distributiongen_temp = 300 ; temperature for Maxwell distributiongen_seed = -1 ; generate a random seed
then i typed these commands:
gmx grompp -f nvt.mdp -c em.gro -p topol.top -n index.ndx -o nvt.tpr   and 
gmx mdrun -deffnm nvt

but GROMACS replied :
                      step 1: One or more water molecules can not be settled.
                      Check for bad contacts and/or reduce the timestep if appropriate.

                      Back Off! I just backed up step1b.pdb to ./#step1b.pdb.19#

                      Back Off! I just backed up step1c.pdb to ./#step1c.pdb.19#
                      Wrote pdb files with previous and current coordinates
                      Segmentation fault (core dumped)