From chemistry-request _-at-_)ccl.net Mon Nov 18 13:43:30 1991 Date: Mon, 18 Nov 91 10:25:01 PST From: ross(+ at +)zeno.mmwb.ucsf.EDU (Bill Ross) To: chemistry "-at-" ccl.net Status: R Here is my promised summary of responses to my question about SHAKE. Another question was raised for implementors of molecular dynamics: does your program have multiple time steps? If not, do you plan to implement them? I'll summarize any answers here. Bill Ross UCSF From: ross ":at:" zeno.mmwb.ucsf.edu (Bill Ross) To: chemistry%!at!%ccl.net Subject: SHAKE failure I'm collecting interesting stories of SHAKE failure in molecular dynamics runs - cases that were never figured out as well as ones that were. References to anything written on this subject would be welcome too. I'll summarize to the reflector. Bill Ross UCSF [The inspiration for this question was occasional SHAKE failure in Amber. Dave Pearlman diagnosed this as stemming from periodic boundary conditions (constant pressure) where ions are treated as part of the solute, all solute-solute interactions are included (no cutoff applied) and so the solute is not imaged with itself: when an ion crosses the edge of the box it is translated to the other side, and if another ion is close by, they suddenly "see" each other and the resulting spike in the virial causes an expansion of the box to equalize the "pressure," which puts some bond lengths beyond the ability of SHAKE to recover. In my opinion, the best way to get around this is to either run with constant volume or treat the ions as part of the solvent so they are imaged. Neither way is completely satisfactory, since in constant volume ions will still suddenly appear in proximity and when ions are treated as part of the solvent cutoffs are applied and long-range electrostatics are lost. I have always run with enough water so that it hasn't happened to me, but this is partially a matter of luck.] >From STOUTEN: at :embl-heidelberg.de Wed Nov 13 01:44:46 1991 Subject: Re: SHAKE failure To: ross.,at,.zeno.mmwb.ucsf.edu Hi Bill, I don't have any interesting stories. I wonder, however, why you do your investigation since SHAKE is being phased out and will be replaced by multiple time step algorithms. Even Wilfred takes this stand. Cheers, Pieter Stouten >From STOUTEN: at :embl-heidelberg.de Thu Nov 14 16:26:56 1991 Subject: Re: SHAKE failure To: ross;at;zeno.mmwb.ucsf.edu Hi Bill, On Wed, 13 Nov 91 09:11:32 PST you wrote: >why even ask? - sociology/folklore of science. > That sounds like a good reason. Still, I don't have exciting stories. When being far from equilibrium I had often problems or when applying heavy torsion constraints (this seems unrelated but is not). Then I just did not use shake. >how soon do you expect multiple time steps to take over? > Hard to say. I am a bit away from the field now. I know that already 3 years ago people were talking about implementing it. As for me, I basically use GROMOS and considering how busy Wilfred c.s. are I don't know when the first official release after GROMOS 87 will see the light. Cheers, Pieter Stouten. #### # # ### # European Molecular Biology Laboratory # ## ## # # # Biocomputing Programme ### # # # ### # Meyerhofstrasse 1, D-6900 Heidelberg, Germany # # # # # # e-mail: stouten%!at!%embl-heidelberg.de #### # # ### #### phone: +49-6221-387 472, fax: 387 517 >From balbes "-at-" osiris.rti.org Wed Nov 13 05:40:27 1991 To: ross;at;zeno.mmwb.ucsf.edu (Bill Ross) Subject: Re: SHAKE failure Okay, this won't be much help, but... I had shake fail after about 14 ps. I had Tom Darden look at it, and he said that because I was saving the steps in binary form, any restart would fail exactly the same way. He said (I think) that he saves things in ascii, so that roundoff on restarting will get around any failures of this type. I have since been using tom's fast amber 3a, and haven't had any more problems of this type. Of course the files have been long since purged. This was quite awhile ago, so details are fuzzy. Lisa %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% standard disclaimer %%%% Lisa M. Balbes, Ph.D. phone: 919-541-6563 Research Triangle Institute, PO Box 12194 vmail: 919-541-6767, xt 6563 Research Triangle Park, NC 27709-2194 email: balbes ":at:" osiris.rti.org - This came directly from a computer and should not be doubted or disbelieved.- >From harris-!at!-athena.mit.edu Wed Nov 13 05:47:02 1991 To: ross /at\zeno.mmwb.ucsf.edu (Bill Ross) Subject: Re: SHAKE failure This example concerns SHAKE's cousin, RATTLE applied to simulation of 50 hexane molecules at liquid density. It is more an example of the effects of the timestep than a complete failure. Using a timestep of .007 ps, if one compares the pressure and temperature computed with the molecular virial with the similar quantities computed using the atomic virial discrepencies are obversed. In a 1.120 nanosecond (160 000 time steps), the average molecular temperature is 303.5 compared to an atomic temperature of 300.5 (there is a thermostat attached to maintain the average temperature at 300K by a weak coupling to an external bath (i.e Berendsen et al. JCP 81, p3684 1984). The pressure from the molecular virial is also about 20 atm higher than that computed from the atomic virial. Cutting the timestep in half reduces the difference between the two methods of computing the temperature and pressure to 1 degree and 5 atm. It appears that the atomic versions of the temperature and pressure are more accurate, but the statistics are too poor to bury the question. Jonathan G. Harris, H. P. Meissner Assistant Professor, Department of Chemical Engineering, MIT Rm 66-450 25 Ames Street, Cambridge, MA 02139 harris -A_T- athena.mit.edu (617)253-5273 Fax 253-9695 ------- From: nobody.,at,.kodak.com To: "amber -8 at 8- cgl.ucsf.edu" -8 at 8- kodak.com Subject: RE: SHAKE nightmares. >From: NAME: Adi M. Treasurywala FUNC: Biophys. & Compu. Chem. TEL: (518)445-7042 To: NAME: Edward P. Jaeger , "amber.,at,.cgl.ucsf.edu".,at,.KODAKR.,at,.MRGATE.,at,.WPC Bill, We have completed a fairly detailed study in AMBER3.0A of this problem. We ran into it because we were trying to do REAL classical dynamics (ie not constant temperature!!!). What happened was that in the initial runs on the molecule that we were interested in we found that the total energy simply did not stabilize but decreased steadily over the run at a fairly precipitous rate. There seemed to be no way to stabilize it. We therefore took a small cyclic peptide WFGLMQ and essentially banged the heck out of it by trying many different conditions to narrow down the reason for this "bug". One of us at least was convinced that it must be something that WE were doing wrong rather than that AMBER could have a bug in it that was so glaring! Anyway, we discovered a LOT of other interesting things along the way such as objective criteria (vs semi religeously based recipies) for knowing when your thermalization was done, and some real reasons to start and finish the collection of MD data (incidentaly it is rather system size dependant...). After much trial and error we found (thanks are due here to George Seibel who REALLY put us on to this) that it was indeed the SHAKE option that had caused us all this grief. If we turned SHAKE OFF on an otherwise identical run that had failed with SHAKE ON, it went just fine!! Everything stabilized and we were able to run real classical dynamics (this is sort of an emotional issue for us here!!). We are in the process of preparing a manuscript covering this work. Let me just ask the gurus on this net... Would that be a good thing to do? Would it help anyone? Would it insult or offend anyone? Incidentaly I learned that the real problem with SHAKE was the integrator (leapfrog just can't handle classical MD with SHAKE turned on I was told by someone.). Hope to hear the other stories soon and to hear any discussions about classical vs constant temp runs! Adi T & Ed Jaeger.