One of the most important quantities in thermodynamics is the reversible
work needed to change the thermodynamic state of a system. Under
certain conditions, this quantity will be independent of the path taken
between the initial and final states and will,
therefore, be related to the free energy difference between these two states.
As a result, considerable effort has been invested in the development of
methods to compute such free energy differences (see, e.g., Ref. [1]
for a review and comparison).
Of particular importance is the free energy profile along a reaction
path characterized by a reaction coordinate q. Not only does the
free energy profile provide a thermodynamic picture along the reaction
path, but it also permits determination of activation energies and
associated rate constants via classical or quantum transition state theory.
Statistical mechanics provides a means whereby such free energy profiles
can be determined directly from ensemble averages. In the case of a
classical N-particle system at temperature T,
the reaction coordinate q is
expressible as a function 
of the N Cartesian position
vectors 
of the N particles.
The free energy profile, F(q'), is then defined by
where
is the probability density for the reaction coordinate to take on the
value q' and 
.
In principle, the probability distribution in Eq. (2) can be computed directly from a molecular dynamics (MD) or Monte Carlo simulation. However, if the reaction path described by q corresponds to a rare event with a high activation energy, then the use of direct simulation techniques is infeasible due to the presence of very low probability regions in the configuration space. Consequently special techniques such as umbrella sampling [2, 3, 4] and the so called ``bluemoon ensemble'' approach [5, 6] have been developed. In these approaches, the reaction coordinate is restrained to lie within certain windows or rigidly fixed at certain values and individual simulations in each window or at each fixed value are performed. The free energy profile is then reconstructed a posteriori via weighted histogram (see, e.g. Ref. [7]) or thermodynamic integration techniques. Despite their success, the necessity of performing numerous individual simulations is nevertheless computationally intensive. In addition, the need to change the reaction coordinate by hand in these methods may involve many nontrivial changes in the other degrees of freedom in a system which require long equilibration times to relax or may, itself, be nontrivial depending on the complexity of the reaction coordinate. Finally, both the umbrella and bluemoon approaches require nontrivial postprocessing of the output data.
In this paper, an alternative approach to the calculation of free energy
profiles along reaction paths is presented. The new method is based on
the creation of a dynamical adiabatic separation between the
reaction coordinate and remaining degrees of freedom. In particular,
a dynamical scheme is constructed in which the reaction coordinate
evolves slowly relative to the other degrees of freedom and is simultaneously
maintained at a high temperature. The latter condition, which has also
been employed on entire molecules by other authors to enhance
sampling of configuration space [8, 9] and which has been used to
perform approximate quantum dynamical simulations [10, 11],
ensures that all activation barriers along the reaction path can be easily crossed
and can be enforced by coupling the reaction coordinate to its
own heat bath or thermostat. The former condition permits the remaining
degrees of freedom to fully relax in response to the motion
of the reaction coordinate and, thereby, sample a large portion of
their configuration space as the reaction coordinate slowly evolves.
A careful analysis of the resulting dynamics reveals that the free energy
profile will then be given by Eq. (1) with 
replaced
by 
, where 
is the temperature of the reaction
coordinate (see Sec. 2).
It is then clear that the use of such an approach eliminates all postprocessing
of the simulation data. It is also found that the new
adiabatic free energy dynamics (AFED)
method allows the free energy profile to be determined with greater
efficiency than constrained/restrained methods such as the bluemoon
and umbrella sampling schemes.
This paper is organized as follows. In Sec. 2, the adiabatic
dynamics method is analyzed in detail and shown to yield the
free energy profile directly from Eq. (1) with
replaced by 
. The analysis is based on the
use of the Liouville operator technique and a formally exact multiple
time scale breakup of the adiabatic classical propagator. In Sec. 3,
several model problems are examined using the new approach.
In these examples, direct comparison to analytical solutions
(when available) and the bluemoon ensemble method are made.
Conclusions and prospects for extending this work to more
complex situations, other types of free energy calculations, and
quantum free energy profiles are discussed in Sec. 4.
