[gmx-users] Crashes when doing free energy due to pressure and LINCS warnings

[Matteo Aldeghi]

Dear users/developers,

I am trying to carry out some alchemical free energy calculations for the
binding of a small molecule. I run energy minimisation (10000 steps), NVT
(10 ps) and NPT (100 ps) equilibration steps before the NVT production run
(1 ns). The problem is that I get crashes at the beginning of the NPT, after
receiving a series of pressure scaling and LINCS warnings. I am trying to
modify the restraints only.

I am having troubles figuring out what I'm doing wrong, since the crashes do
not happen at low lambda values (e.g. lambda <= 0.1) or if I run the same
simulation with full restraints but without using the free energy code
(free-energy = no). I ran the simulations more times obtaining the same
results. I apply the distance restraint using the pull code, whereas I
specify 2 angle and 3 dihedral restraints in the topology (as shown below).

Also, I tried running simulations for the decoupling of the ligand charges
(defining a B state with zero charges) using:
fep-lambdas              = 0.0 0.25 0.5 0.75 1.0
restraint-lambdas	 = 1.0 1.0 1.0 1.0 1.0
Now, I ran a simulation at lambda = 0.5 where I commented out all the
restraints from the topology and mdp files, and it crashed. On the other
hand, I ran the same simulation (restraints commented out) twice, but
filling the restraint-lambdas vector with 0.0, and they both ran fine until
completion.

I noticed that the crashed simulations have large negative pressures (-10^5
bar) and large positive potentials.

I'm using: 
- gromacs 4.6.3 and run mdrun_mpi with particle decomposition. 
- Amber99SB-ILDN and GAFF with AM1-BCC charges for the ligand.

Any suggestion would be much appreciated,

Matteo



===============================================================================
                                                                       
Additional notes 
===============================================================================
* I tried playing around with the tau_p value with no success.
* I ran a couple of simulations with 4.6.1 and they crashed too at the same
point.
* the B-state constant for the pull code is not defined as it should, as I
mean to have k1=0 and kB1=4184 (rather than only k1=4184). This is a
mistake, however I believe this shouldn't affect the crashes. 
* I couldn't find information about whether there should be a B-state
specified for the angle and dihedral restraints, am I right to assume the
lambda value scales the force constant directly?

================================================================================
                                                                     
Example of warnings
================================================================================

Step 17, time 0.034 (ps)  LINCS WARNING
relative constraint deviation after LINCS:
rms 349.859101, max 17941.974609 (between atoms 2297 and 2298)
bonds that rotated more than 30 degrees:
 atom 1 atom 2  angle  previous, current, constraint length
   2084   2086   89.6    0.1330   0.1240      0.1335
   2086   2087   89.8    0.1006   0.1440      0.1010
# many more bonds follow

Step 18  Warning: pressure scaling more than 1%, mu: 156.752 156.752 156.752

# a mix of these 2 warnings continues until the simulation crashes with a
LINCS fatal error due to too many warnings.

================================================================================
                                                                Restraints
in the topology file
================================================================================

[ angle_restraints ]
; ai     aj    ak    al    type     th0    fc           multiplicity
 1393   1391   2615  1391  1        88.8   4.1840e+01   1
 1391   2615   2614  2615  1        32.9   4.1840e+01   1

[ dihedral_restraints ]
; ai    aj    ak    al       type   phi     dphi   kfac
  1410  1393  1391  2615     1      -159.7  0.0    4.1840e+01
  1393  1391  2615  2614     1      122.6   0.0    4.1840e+01
  1391  2615  2614  2610     1      12.8    0.0    4.1840e+01

mdp_files.tar.gz
<http://gromacs.5086.x6.nabble.com/file/n5013820/mdp_files.tar.gz>  

--
View this message in context: http://gromacs.5086.x6.nabble.com/Crashes-when-doing-free-energy-due-to-pressure-and-LINCS-warnings-tp5013820.html
Sent from the GROMACS Users Forum mailing list archive at Nabble.com.