The bar module of GROMACS makes the calculation of ΔG_AB very simple. Simply collect all of the md*.xvg files that were produced by mdrun (one for each value of λ) in the working directory and run gmx bar:

gmx bar -f Lambda_*/Production_MD/md*.xvg -o -oi

The program will print lots of useful information to the screen, in addition to producing three output files: bar.xvg, barint.xvg, and histogram.xvg. The screen output from gmx bar should look something like:

Temperature: 298 K

Detailed results in kT (see help for explanation):

 lam_A  lam_B      DG   +/-     s_A   +/-     s_B   +/-   stdev   +/- 
     0      1    2.53  0.01    0.19  0.01    0.19  0.01    0.66  0.01
     1      2    2.12  0.00    0.21  0.02    0.21  0.02    0.65  0.00
     2      3    1.73  0.02    0.19  0.01    0.18  0.01    0.61  0.01
     3      4    1.36  0.01    0.18  0.01    0.18  0.01    0.59  0.00
     4      5    1.03  0.01    0.15  0.01    0.14  0.01    0.54  0.00
     5      6    0.76  0.01    0.13  0.01    0.12  0.01    0.49  0.00
     6      7    0.54  0.00    0.10  0.01    0.10  0.01    0.45  0.00
     7      8    0.36  0.00    0.09  0.00    0.08  0.00    0.41  0.00
     8      9    0.20  0.00    0.08  0.00    0.07  0.00    0.38  0.00
     9     10    0.06  0.00    0.06  0.00    0.06  0.00    0.37  0.00
    10     11    0.27  0.00    0.04  0.00    0.05  0.00    0.30  0.00
    11     12    0.26  0.00    0.05  0.00    0.06  0.00    0.31  0.00
    12     13    0.25  0.00    0.04  0.00    0.05  0.00    0.32  0.00
    13     14    0.22  0.00    0.06  0.00    0.07  0.00    0.33  0.00
    14     15    0.20  0.00    0.06  0.00    0.07  0.00    0.35  0.00
    15     16    0.17  0.00    0.06  0.00    0.07  0.00    0.37  0.00
    16     17    0.14  0.00    0.07  0.01    0.09  0.01    0.39  0.01
    17     18    0.10  0.00    0.07  0.01    0.09  0.01    0.41  0.00
    18     19    0.06  0.01    0.09  0.01    0.11  0.01    0.44  0.00
    19     20   -0.01  0.00    0.11  0.01    0.14  0.01    0.48  0.00
    20     21   -0.11  0.01    0.14  0.01    0.18  0.01    0.54  0.00
    21     22   -0.22  0.01    0.15  0.01    0.19  0.01    0.59  0.00
    22     23   -0.36  0.01    0.20  0.00    0.25  0.00    0.65  0.01
    23     24   -0.56  0.03    0.25  0.02    0.32  0.02    0.73  0.02
    24     25   -0.78  0.02    0.27  0.04    0.33  0.04    0.80  0.01
    25     26   -0.95  0.03    0.27  0.02    0.29  0.02    0.75  0.01
    26     27   -0.91  0.02    0.13  0.01    0.13  0.01    0.55  0.00
    27     28   -0.69  0.01    0.08  0.01    0.07  0.01    0.39  0.01
    28     29   -0.43  0.00    0.04  0.00    0.04  0.00    0.27  0.00
    29     30   -0.16  0.00    0.01  0.00    0.01  0.00    0.19  0.00


Final results in kJ/mol:

point      0 -      1,   DG  6.26 +/-  0.02
point      1 -      2,   DG  5.26 +/-  0.01
point      2 -      3,   DG  4.28 +/-  0.04
point      3 -      4,   DG  3.37 +/-  0.04
point      4 -      5,   DG  2.56 +/-  0.02
point      5 -      6,   DG  1.89 +/-  0.02
point      6 -      7,   DG  1.34 +/-  0.01
point      7 -      8,   DG  0.89 +/-  0.01
point      8 -      9,   DG  0.49 +/-  0.01
point      9 -     10,   DG  0.16 +/-  0.01
point     10 -     11,   DG  0.68 +/-  0.00
point     11 -     12,   DG  0.64 +/-  0.01
point     12 -     13,   DG  0.61 +/-  0.01
point     13 -     14,   DG  0.56 +/-  0.01
point     14 -     15,   DG  0.48 +/-  0.01
point     15 -     16,   DG  0.43 +/-  0.01
point     16 -     17,   DG  0.35 +/-  0.01
point     17 -     18,   DG  0.26 +/-  0.01
point     18 -     19,   DG  0.15 +/-  0.01
point     19 -     20,   DG -0.02 +/-  0.01
point     20 -     21,   DG -0.27 +/-  0.02
point     21 -     22,   DG -0.54 +/-  0.02
point     22 -     23,   DG -0.90 +/-  0.03
point     23 -     24,   DG -1.40 +/-  0.06
point     24 -     25,   DG -1.94 +/-  0.05
point     25 -     26,   DG -2.36 +/-  0.07
point     26 -     27,   DG -2.25 +/-  0.05
point     27 -     28,   DG -1.70 +/-  0.02
point     28 -     29,   DG -1.06 +/-  0.01
point     29 -     30,   DG -0.41 +/-  0.01

total      0 -     30,   DG 17.80 +/-  0.20

The value of ΔG_AB I obtained is 17.80 ± 0.20 kJ mol^-1, or 4.25 ± 0.05 kcal mol^-1. Since the process I conducted for this demonstration was the transformation of ethanol into a dummy molecule, the reverse process (the introduction of uncharged ethanol into water, thus the actual hydration free energy of the process) corresponds to a ΔG_AB of -4.25 ± 0.05 kcal mol^-1 (assuming reversibility). This value differs slightly from the experimental value of -5.01 kcal mol^-1 and the value of -5.34 ± 0.12 kcal mol^-1 previously reported with CHARMM22. Deviations are likely due to the fact that the topology used here came from the CGenFF parameter set, rather than the ethanol model compound optimized for the CHARMM protein force field; there are minor differences in partial charge assignments and Lennard-Jones parameters of the constitutent atom types.

As in the example of methane in water, one can plot the output from the BAR program to see how the free energy differences between neighboring λ windows resulted in the final output. We will not plot them in this tutorial.