GROMACS version: 2022
GROMACS modification: No
Hello everyone, I am trying to make free energy calculation for the first time and have a problem. I am mutating aminoacids in protein and using fast-growth method to calculate free energy of stability. Using pmx web server, I’ve created hybrid topology(which contains dummy atoms), then placed protein in a box, neutralized it with ions, and conducted energy minimization procedure.
Then after command
gmx grompp -f nvt.mdp -c em.gro -r em.gro -p hybrid.top -o nvt.tpr -maxwarn 2
I’ve got this warning:
The bond in molecule-type Protein_chain_A between atoms 5324 DOH and 5325
DHH has an estimated oscillational period of 6.5e-03 ps, which is less
than 5 times the time step of 2.0e-03 ps.
Maybe you forgot to change the constraints mdp option.
And then, after gmx mdrun -deffnm nvt
I got this error:
There are perturbed non-bonded pair interactions beyond the pair-list cutoff
of 1.111 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.
I have no error, if I just try to simulate this molecule without hybrid topology, so it seems that dummy particles is a source of this problem. Could please somebody help me to solve this problem. Thanks.
Here are my mdp files
minimization.mdp
;Run parameters
define = -DFLEXIBLE
nsteps = 1500
integrator = steep
emtol = 50.0
emstep = 0.01
;Neighbour searching
cutoff-scheme = Verlet
ns-type = grid
rcoulomb = 1.0
vdwtype = cut-off
rvdw = 1.0
nstlist = 10
rlist = 1
;Eletrostatics
coulombtype = PME
;Periodic boundary conditions
;pbc = xyz
;ld-seed = 11
; FREE ENERGY
;----------------------------------------------------
free-energy = yes
init-lambda = 0 ; start from state A
delta-lambda = 0 ; stay in state A
;
nstcomm = 100
constraints = none
coulomb-modifier = Potential-shift-Verlet
pme-order = 4
fourierspacing = 0.12
ewald-rtol = 1e-5
vdw-modifier = Potential-shift-Verlet
verlet-buffer-tolerance = 0.005
DispCorr = EnerPres
Tcoupl = no
Pcoupl = no
gen_vel = no
nvt.mdp
;This mdp file has been created by the pymdsetup.gromacs_wrapper.grompp.create_mdp()
;Type of MDP: nvt
;Position restrain
Define = -DPOSRES
;Run parameters
nsteps = 50000
;integrator = md
integrator = sd
dt = 0.002
comm-mode = Linear ; remove center of mass translation
;Output control
nstxout = 500
nstvout = 500
nstenergy = 500
nstlog = 500
nstcalcenergy = 100
nstcomm = 100
nstxout-compressed = 1000
compressed-x-precision = 1000
compressed-x-grps = System
;Bond parameters
constraint-algorithm = lincs
constraints = h-bonds
lincs-iter = 1
lincs-order = 4
continuation = no
;Neighbour searching
cutoff-scheme = Verlet
ns-type = grid
rcoulomb = 1.0
vdwtype = cut-off
rvdw = 1.0
nstlist = 10
rlist = 1
;Eletrostatics
coulombtype = PME
pme-order = 4
fourierspacing = 0.12
fourier-nx = 0
fourier-ny = 0
fourier-nz = 0
ewald-rtol = 1e-5
;Temperature coupling
tcoupl = V-rescale
tc-grps = Protein Non-Protein
tau-t = 0.1 0.1
ref-t = 300 300
;Pressure coupling
pcoupl = no
;Dispersion correction
DispCorr = EnerPres
;Velocity generation
gen-vel = yes
gen-temp = 300
gen-seed = -1
;Periodic boundary conditions
pbc = xyz
ld-seed = 17
; FREE ENERGY
;----------------------------------------------------
free-energy = yes
init-lambda = 0 ; start from state A
delta-lambda = 0 ; stay in state A
;
vdw-modifier = Potential-shift-Verlet