[gmx-users] Stuck --> unstuck ...

[David van der Spoel]

Stéphane Téletchéa wrote:
> I've found out where the problem lie (at least one of ...).
> 
> Small reminder: i've tried to setup a system consisting of HIV protease 
> + one ligand on a water box. I used the excellent protocol from John E. 
> Kerrigan (Drug-enzyme tutorial) in order to set up my simulation and 
> encountered some problems with it (basically: Box was shifted at least 
> 10 times. Please see log-file.).
> 
> I add some constructive comments on ways to explore where the problem 
> was (changing the protocol, including ions into the water to get a 
> better coupling bath, etc.). Thanks to all of them, even if sometimes 
> humour is hard to catch by email!
> 
> I'm not interested in the system per se (i've done other dynamics 
> working fine with this protocol) but was rather interested in providing 
> an updated way of doing a drug-enzyme tutorial with latest force fields 
> (g45a3, g53a6) and using the latest prodrg (since it allow the output of 
>  g43a1-type atoms, mostly transferable to the newer ones).
> I'll put this on the wiki soon.
> 
> I've checked many parameters except one: all-bons vs hbonds.
> 
> I do not recall a 'rationale' reason for using the all-bonds parameter 
> (instead of hbonds where it worked for other dynamics simulations) but 
> the fact is that if i run all-bonds simulations (whatever the starting 
> system) with NVT+NPT, i see a crash with the error message like this:
> 
> Step 1  Warning: pressure scaling more than 1%, mu: -1.59647e+20 
> -1.59647e+20 -1.59647e+20
> Correcting invalid box:
> old box (3x3):
>    old box[    0]={-9.82630e+20,  0.00000e+00, -0.00000e+00}
>    old box[    1]={ 0.00000e+00, -8.90833e+20, -0.00000e+00}
>    old box[    2]={ 0.00000e+00,  0.00000e+00, -1.12456e+21}
> new box (3x3):
>    new box[    0]={-9.82630e+20,  0.00000e+00, -0.00000e+00}
>    new box[    1]={ 0.00000e+00, -8.90833e+20, -0.00000e+00}
>    new box[    2]={ 0.00000e+00,  8.90833e+20, -1.12456e+21}
> 
> using the *same* protocol (same starting structures, etc) and only 
> changing all-bonds to h-bonds leads to a 'normal' trajectory (no error, 
> no warning, no crash).
> 
> I remember gromacs 3.3.0 had a problem where pme=4 (that lead to the 
> 3.3.1 rapid release). There are a lot of reports on this mailing-list 
> with crash/segfaults when using all-bonds, and it seems more when 
> all-bonds are used instead of h-bonds.
> 
> Of course it may come from not completely converged dynamics (my setup 
> is very rapid, the 4 PR steps last only 100 ps), from not enough 
> minimsed structures, etc.
> 
> Could it be also possible there is a problem while using all-bonds (even 
> with a 2fs time step) instead of in gromacs, and pme_order was just 
> hiding it ?

I don't think so. All-bonds with 2 fs should work reliably with both 
lincs and shake, provided your system is equilibrated to some extent.



-- 
David.
________________________________________________________________________
David van der Spoel, PhD, Assoc. Prof., Molecular Biophysics group,
Dept. of Cell and Molecular Biology, Uppsala University.
Husargatan 3, Box 596,  	75124 Uppsala, Sweden
phone:	46 18 471 4205		fax: 46 18 511 755
spoel at xray.bmc.uu.se	spoel at gromacs.org   http://folding.bmc.uu.se
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++