GROMACS

When doing slow-growth FEP (and perhaps other types of FEP), the position restraints are not read from the grompp -rb file if the same filename is passed to both the grompp -r and -rb options.
One might desire the same values for position restraints in the A and B states to keep the position restraint on and invariant of what is being changed in the FEP. I was doing this because I was inserting a ligand to the cognate ligand binding site of a protein and wanted to keep an allosteric binding site from changing conformation.

I used .mdp options like this:

free-energy = yes
init-lambda = 0
couple-lambda1 = vdw
couple-lambda0 = none
couple-intramol = no
couple-moltype = XOL
delta-lambda=0.0002
sc-alpha=0.5
sc-power=1
sc-r-power=6
sc-sigma=0.3

However, the restrained atoms were all pulled to the Cartesian origin as the simulation approached the B-state. This is because the posres_comB array was all zeroes. That's due to the following code in kernel/grompp.c:

static void gen_posres(gmx_mtop_t *mtop, t_molinfo *mi,
char *fnA, char *fnB,
int rc_scaling, int ePBC,
rvec com, rvec comB,
warninp_t wi) {
int i, j;

read_posres  (mtop, mi, FALSE, fnA, rc_scaling, ePBC, com, wi);
    if (strcmp(fnA, fnB) != 0)
{
        read_posres(mtop, mi, TRUE, fnB, rc_scaling, ePBC, comB, wi);
    }
}

I'll also note that if one does not define either -r or -rb, then again the posres_com array is fine (set to the -c positions) but the posres_comB values are all zeroes.

I'm not sure what the original purpose was with the above code, so I won't suggest any particular changes, but I don't think that the current behaviour is optimal, especially since it leads to -rb being ignored if the -r file has the same name as the -rb file, which seems like a common enough desire.

Now that I know what is going on, I am circumventing this by copying a.gro to b.gro and using grompp -r a.gro -rb b.gro

Chris.