[gmx-developers] multicomponent lambda free energy calculations.

[Sander Pronk]

Hi Michael,

This sounds really useful, and I personally think something like this should be included in the main distribution: it is the natural generalization of the current code to something more usable.

One question:
Does this mean that the couple-lambda0 and couple-lambda1 keywords are essentially obsolete?

I have two suggestions:
1) there should probably be a pull code lambda

2) A small suggestion: the naming of these mdp variables is slightly confusing. I'd suggest naming the lambda sets 'lambdas': that way the idea that this is a set is reinforced. Also, the main setting (fep-lambda) can be named 'lambdas'. The lambda point identifier 'init_fep_state' should probably be something like 'init_lambda_state'. 

init_lambda_state = 
lambdas = ...
coul-lambdas = ...
vdw-lambdas = 

Sander


On 15 Jan 2010, at 05:46 , Michael Shirts wrote:

> I've put together a reworking of the foreign lambda calculation in the
> CVS to make it easier to use best practices for free energy
> calculations.  This code is caught up to CVS master (as of today).
> This can be downloaded for testing with the command:
> 
> git clone git://git.gromacs.org/~mrshirts/gromacs_mrs.git
> 
> What is described below sounds userful, You should 1) try to see if
> you can rewrite your topologies in a way that can use fewer separate
> topologies with this new code and 2) see if you can set up tests that
> -verify- that what is being done is actually correct, comparing
> between the old and new codes.  I've checked for several test cases,
> but there may always be some new test case that I missed.
> 
> Key points:
> 
> 1.  Previous functionality is preserved (including couple keywords),
> with the exception that instead of using the key word
> 'foreign_lambda', the key word 'fep-lambda' is used instead.  A bug in
> the free energy calculations code 1,4 interactions without LJ terms
> was found -- this has also been fixed in the main branch as well.
> 
> 2. Lambda can now be expressed as a vector of the different
> interaction terms; currently 'coul-lambda', 'vdw-lambda',
> 'bonded-lambda', 'restraint-lambda' and 'mass-lambda' are possible,
> but finer/different divisions are possible in the future as needed.
> 
> 
<snip>

> Examples:
> 
> Example A:
> fep-lambda      = 0.0 0.0 0.0
> coul-lambda     = 0.0 0.8 1.0
> vdw-lambda      = 0.0 0.2 1.0
> 
> Coul and vdw terms change, all other terms remain in A state; the same
> as if fep-lambda was omitted.
> 
> Example B:
> fep-lambda        = 0.0 0.0 0.0
> coul-lambda       = 0.0 0.8 1.0
> vdw-lambda        = 0.0 0.2 0.4
> Restraint-lambda  = 0.0 0.4 1.0
> Bonded-lambda     = 0.0 0.1 1.0
> Mass-lambda       = 0.0 0.0 1.0
> 
> Fep-lambda does nothing, since all defined components are specified.
> 
> Example C:
> fep-lambda        = 0.0 0.2 0.8
> coul-lambda       = 0.0 0.8 1.0
> 
> Vdw, restraint, bonded, mass all have 0.0 in the first state, 0.2 in
> the second state, and 0.8 in the last state, copied from fep-lambda.
> 
> This scheme should be sufficient to handle almost all possible
> transformations, but can be extended fairly easily as needed.