Dear Gromacs users
I am doing a test calculation in MD simulation between ligand and a protein residue for bond angle and dihedral restraints. I have set them in my topology file as below:-
; distance restraints [ bonds ] ; i j type r0A r1A r2A fcA r0B r1B r2B fcB 71 29 10 0.403 0.403 10.0 0.0 0.403 0.403 10.0 4184.000
[ angle_restraints ] ; ai aj ak al type thA fcA multA thB fcB multB 29 71 69 71 1 104.92 0.0 1 104.92 41.840 1 71 29 31 29 1 147.02 0.0 1 147.02 41.840 1
[ dihedral_restraints ] ; ai aj ak al type phiA dphiA fcA phiB dphiB fcB 69 71 29 31 1 146.72 0.0 0.0 146.72 0.0 41.840 72 71 29 31 1 17.55 0.0 0.0 17.55 0.0 41.840 72 71 29 25 1 173.10 0.0 0.0 173.10 0.0 41.840
However, after energy minimization, I observed that my restraints are not working and I find the following change:-
bonds:-
0.403 changed to 0.413
angles:-
104.92 changed to 109.24
147.02 changed to 147.65
dihedrals:-
146.72 changed to 133.39
17.55 changed to 12.94
173.10 changed to -166.16
How should I work with these so that they may not change much? Dihedrals seem to be changing more than bonds and angles. I have to apply these in the free energy calculation which already crashed because of restraint failure.
I would really be thankful for your kind suggestions, please.
Thanks.
Sadaf
Sadaf Rani It may be that the position/dihedral restraints just aren't being turned on. Does your *.log file report an energy value for the position restraints? If not, that may be an indicator that your position restraints, though defined in your topology, are not being read during the simulation.
Anther way to check: Does your *.mdp file define the position restraints to tell GROMACS that the restraints are on?
For instance, anytime I have position restraints, my *.mdp file will start with a line like:
define = -DPOSRES_ON
Dear JD Tamucci and Anwesh Pandey I have set following in my mdp file:-
Run control
integrator = sd ; stochastic leap-frog integrator
define = -DPOSRES
I have edited my topology file as well as mentioned above and defined state B of my system as well.
I have added position restraint for some of residue atoms and for ligand I kept distance, angle and dihedral angles only
I am getting following in log file:-
Step Time
0 0.00000
Energies (kJ/mol)
Bond Restraint Pot. Angle Proper Dih. Improper Dih.
4.27459e+01 8.76032e-03 1.87138e+02 1.46976e+02 1.62133e+00
LJ-14 Coulomb-14 LJ (SR) Disper. corr. Coulomb (SR)
6.29960e+01 6.96927e+02 8.20600e+03 -2.80515e+02 -6.44091e+04
Coul. recip. Position Rest. Angle Rest. Dih. Rest. Potential
3.01125e+02 5.95019e-02 2.17501e+00 7.48658e+00 -5.50344e+04
Kinetic En. Total Energy Temperature Pres. DC (bar) Pressure (bar)
9.98715e+03 -4.50472e+04 2.93974e+02 -1.14243e+02 -8.16560e+00
dVcoul/dl dVvdw/dl dVrestraint/dl Constr. rmsd
1.31273e+03 -4.29347e+02 9.67035e+00 2.79463e-06
How can I know that my restraints are being read?
Thanks in advance.
Sadaf
- Badges
- Science topic