[gmx-users] "B" state for free energy corrupts minimization?

[Berk Hess]

>From: mernst at tricity.wsu.edu
>Reply-To: Discussion list for GROMACS users <gmx-users at gromacs.org>
>To: gmx-users at gromacs.org
>Subject: [gmx-users] "B" state for free energy corrupts minimization?
>Date: Tue, 31 Jan 2006 00:16:11 -0800 (PST)
>
>I am attempting to perform a free energy simulation that involves the
>transformation of the lesioned DNA base 8-oxoguanine to guanine. Before I
>set up the B state in my .itp file, I was able to simulate 8-oxoguanine
>with MD and its behavior looked to be in reasonable agreement with
>simulations I have performed in another package.
>
>Now that I have defined B-state parameters, MD and even energy
>minimization give strange and undesirable results, even before I begin the
>free energy calculations or set lambda to a non-zero value.
>
>I've narrowed it down to a simple test case where I supply a B state for
>just one atom.
>
>If I use this line in my .itp file, which changes nothing between A and B,
>the energy minimization gives a reasonable-looking structure:
>
>46 amber99_35 2 8oG  N7   46   -0.5129    14.01   amber99_35 -0.5129 14.01
>
>This is what the output looks like:
>
>http://www.sciencemadness.org/scipics/goodtrimer.png
>
>If I use this line, using a different atom type in the B state,
>
>46 amber99_35 2 8oG  N7   46   -0.5129    14.01   amber99_36 -0.5129 14.01
>
>the minimization gives me a hydrogen that sticks out at a strange angle:
>
>http://www.sciencemadness.org/scipics/badtrimer.png
>
>The amber99_35 atom is AMBER type NA, corresponding to sp2 N in 5
>memb.ring w/H atom (HIS). The amber99_36 atom is AMBER type NB,
>corresponding to sp2 N in 5 memb.ring w/LP (HIS,ADE,GUA). In the A state
>there is a hydrogen that becomes a dummy in the B state, so I would like
>to go from NA to NB. The AMBER force field port to Gromacs is from the
>Pande group.
>
>I cannot figure out why this effect shows up when lambda=0, and how I
>might set up the free energy simulation so that it is not horribly
>distorted and prone to crashing in the equilibrium period.

This sounds very strange.
B-state parameters should have no influence on the forces at lambda=0.
I have done many free-energy simulations with the Amber force-field
without problems.

Three questions:

Have you tried the free-energy topology with free-energy in the mdp
file turned on and off, and is there a difference?

Could you compare the normal and free-energy potential energies
for the first step of the minimization?
This should tell you in which term the difference is located.

You can compare tpr files with gmxcheck -s1 ... -s2 ...
You should do this for your normal and free-energy tpr.
You might get a lot of output if the 36 atom type is not at all
present in the A-state, but we are probably looking for
differences in the bonded interactions.

Berk.