[gmx-developers] doing energy evaluations at alternate lambda values (i.e. for bennett acceptance ratio)

[Berk Hess]

David Mobley wrote:

> Dear all,
>
> We're attempting to modify mdrun to be able to write out the system
> total potential energy evaluated at a number of alternate lambda
> values during a free energy calculation, when energies are written to
> output. (The goal is to use this output for doing free energy
> calculations with Bennett Acceptance Ratio or other methods which
> require potential energies evaluated at a variety of alternate lambda
> values).
>
> I'm assuming this would be a straightforward (perhaps 20 minutes?)
> task for anyone who already knows the inner workings of the code, but
> since I'm not that familiar with it, it seems a bit more difficult for
> me. Is there anyone who wants to volunteer?
>
> If not, can anyone help by telling me the complete set of function
> calls that need to be done in the inner loop of do_md(), in order to
> safely:
> 1) Change the current lambda value to some other lambda'
> 2) Compute the new total energy
> 3) Restore lambda to its original value so the simulation can continue.
>
> Thanks very much in advance.
>
The main problem is that you want the forces from the correct lambda value.
You probably also want efficiency, otherwise you it would be more efficient
to do separate runs.
I think both things can be combined by:
first calling do_force as normal,
then calling do_force without neighbor searching and without 
calculcating forces,
this can be accomplished with parameters passed to do_force.
 From your mail it seems as if your bonded interactions do not influence 
dvdl,
then you can also not calculate bonded interactions, this also makes 
sure the bonded
forces are not double counted.

The only thing I am not sure of is if the virial is not overwritten.

For efficiency you also want to call only the free-energy energy loops 
in nonbonded.c
for the other lambda values.

Berk.