[gmx-developers] position restraints, Jarzynski, and virial

[David Mobley]

Dear all,

So, I am trying to come up with an easy way to implement time
dependent position restraints so that I can use the Jarzynski approach
to do some fast-switching calculations. Actually, what I have in mind
is to begin with a protein restrained to a particular structure, and
then turn off the restraints in a time dependent manner and use the
Jarzynski approach to compute the free energy for doing this. But
that's all really beside the point here -- basically, all you need to
know is that I want time dependent position restraints.

It seems to me that the easiest way to do this would be to use the
fact that lambda can already be time-dependent, and just add an input
flag that makes this time-dependent lambda be used as a prefactor for
the position restraint spring constant. Does this sound reasonable, or
is there another way that is obviously better/easier?

>From looking at the code in bondfree.c of version 3.3, it looks to me
like all I should need to change is the call to "harmonic" on line
993:

*dvdlambda += harmonic(pr->posres.fcA[m],pr->posres.fcB[m],
                             0,pr->posres.pos0B[m]-pr->posres.pos0A[m],
                             dx[m],lambda,&v,&fm);

It seems like I should be able to just change the first two arguments
to lambda*posres.fcA[m] and lambda*posres.fcB[m] if my input flag for
doing Jarzynski is turned on. Am I missing something else obvious?

Also, I notice the following comment in the posres section of bondfree.c:

/* Position restraints should not contribute to the virial.
       * We need to subtract this contribution !!!
       * Needs to be implemented, but requires an extra virial matrix.
       */

This seems to be a note for something that needs to be done. Has this
been resolved?

Thanks,
David