Questions and problems?

# Energy minimization

Does the parameter NSNB (number of steps between two updates of the nonbonded list) influence the energy gradient? For me it looks as if the energy gradient is reset. If I set NSNB equal to MAXCYC (Maximum number of cycles) the minimization is successful, but if I keep NSNB to the default value of 25 the minimizations ends up in an endless loop with no futher progress in the rms value of the energy gradient.

Every time the nonbonded list is updated, the set of interacting atom pairs is likely to change, which can set the minimization off in a different 'direction'. It is even logically possible that a circular path could occur in this case, such that minimization would never terminate, although it is far more likely that the molecule would wander around indefinitely in a more-or-less infinitesmal space. One approach would be to set a step limit on the minimization - 100-500 is usually sufficient to prepare for running MD. Another solution (especially appropriate for fine-grain minimization for normal mode calculation) would be to set the cutoff to include the entire system and NSNB=9999999.

Bill Ross

My goal is to find a very low gradient accurately on a large system without allowing the system to collapse on itself due to the large number of steps required using conjugate gradient methods. I have run the structure thru 600 steps of steepest descent followed by 1000 steps of conjugate gradient all the while keeping the backbone atoms completely constrained. After releasing the backbones I ran an additional 1000 steps of conjugate gradient. This results in a structure that has not collapsed very much but has a reasonable amount of relaxation from the initial hand built conformation.

• You need to think hard about why you need a "very low gradient", but I'll leave that to you. (Normal mode analysis is about the only thing I can think of that really requires it.)
• What potential model you need to use to find a true minimum that does not have collapsed helices (or whatever) is completely orthogonal to the question of whether to use CG or NR to get there. Current potential functions were designed primarily for MD simulations in water; they typically do not have true minima that are close to X-ray structures (where the latter exist). Example: even for a simple protein like crambin, I have never been able to find a true minimum vaccum structure that is closer than about 0.7 Ang. backbone RMS to the X-ray structure. Yet crystal simulations at finite temperature (including the waters present in the crystal structure) stay within about 0.3-0.4 Ang. of the X-ray result.
• Once you have a potential model, what is the best way to find a true minimum? As you have found out [using nmode], NR is generally not even possible for a large system. CG minimization (probably using the nmode algorithm) is the best I know about. Plan to run 20,000 cycles or more, freezing the non-bonded list at some point, so that non-bonded updates don't provide "noise" to the energy function. The system may well drift further from the starting structure than you would like, but it should be no different with NR (assuming you had enough memory to get the latter to run.) CG methods can get you to RMS gradients of 10**-6 or better, but you do have to be patient.
Dave Case
usig nmode ... minimizing from the initial hand built conformation

You are using a minimizer designed to grind small molecules to the absolute lowest energies for a huge, hand-built model. This is like trying to clean Brooklyn Bridge with a toothbrush. To make such a model reasonable, it is way way better in my experience to run md in the range of 1-100K, possibly using belly or restraints to keep structure in the parts you trust more.

Bill Ross