%From schuss@math.tau.ac.il Sun Nov 19 07:55:14 1995
%From: Zeev Schuss <schuss@math.tau.ac.il>
%Date: Sun, 19 Nov 1995 15:57:30 +0200
%To: bob@aix550.phys.rpslmc.edu
%Status: R

\documentstyle[12pt]{article}        
\textwidth=16cm
\textheight=8.5in
\oddsidemargin 0in
\topmargin -0.25in
\begin{document}
\pagestyle{empty}
\begin{center}{\bf
ACTIVATED RATE PROCESSES: A RELATION BETWEEN
HAMILTONIAN AND STOCHASTIC THEORIES}\\

\vspace{5mm}
by\\

\vspace{5mm}
Eli Pollak and Alexander M. Berezhkovskii\footnote{ 
Permanent address: Karpov Institute of Physical Chemistry,
ul. Obukha 10, 103064 Moscow K-64, Russia.}\\
Chemical Physics Department\\
Weizmann Institute of Science, 76100 Rehovot, Israel\\

\vspace{5mm}
Zeev Schuss\\
School of Mathematical Sciences\\
Tel Aviv University, Ramat-Aviv,
69978 Tel Aviv, Israel

\vspace{10mm}
{\bf ABSTRACT}
\end{center}

\vspace{5mm}
Kramers' treatment of activated rate processes is based on the Langevin
equation of motion for the escaping particle. The stochastic dynamics may
be cast equivalently as the dynamics of a particle interacting
bilinearly with a bath of harmonic oscillators.  This paper explores the
connection between the solutions of Kramers' problem (and its generalization
to include memory friction) obtained in the framework of these two approaches.
We demonstrate their equivalence for the specific case of a
parabolic barrier potential.
The Hamiltonian representation is used to
construct (a) a nontrivial eigenfunction of the Fokker-Planck equation which
is generalized to include time dependent friction; (b) the 'Kramers'
stationary flux distribution function; (c) the stochastic separatrix.
\pagestyle{plain}
\newpage
\renewcommand{\theequation}{1.\arabic{equation}}
\setcounter{equation}{0}
\setcounter{section}{1}

\centerline{\bf 1. Introduction}

\vspace{5mm}
The classical representation of activated
rate processes as described by Kramers \cite{kramers}
is in terms of a Fokker-Planck equation for the
probability density function (pdf) $W(p,q,t)$ of the particle
position ($q$) and momentum ($p$),
\begin{eqnarray}
{\partial W(p,q,t)\over \partial t}=LW(p,q,t).\label{fpe}
\end{eqnarray}
The Fokker-Planck operator ($L$) has the form
\begin{eqnarray}
L{\equiv}-{\partial\over\partial q}p+
{\partial\over\partial p}\left[\gamma p+{dV(q)\over \,dq}\right]+
{\gamma\over \beta}{\partial^2\over\partial p^2},\label{fpo}
\end{eqnarray}
where $V(q)$ is the potential acting on the particle, $\gamma$ is the
damping coefficient, and $\beta=1/k_BT$, where $T$ is the temperature
of the bath with which the particle interacts.
The Fokker-Planck equation describes the time ($t$) evolution
of the pdf $W(p,q,t)$ of the particle.  Much effort has
been expended in finding solutions of the partial differential
equation (\ref{fpe}),
as described in standard textbooks \cite{risken}-\cite{schuss}.

It is also well understood that the Fokker-Planck equation
may be derived
from the stochastic Langevin equation of motion for the particle,
\begin{eqnarray}
{\ddot q}+{dV(q)\over \,dq}+\gamma{\dot q}=\xi(t),\label{Langevin}
\end{eqnarray}
where the dot denotes differentiation with respect to time
and $\xi(t)$ is a $\delta$-correlated Gaussian random force
with zero mean, related to the damping by
\begin{eqnarray}
\langle\xi(t)\xi(\tau)\rangle=\frac{2\gamma}{\beta}\delta(t-\tau).
\label{correl}
\end{eqnarray}
For simplicity, we have assumed unit mass for the particle,
or equivalently, the coordinate $q$ and momentum $p={\dot q}$
are mass weighted.

In an activated rate process, the potential $V$ is assumed to
have a well at some location $q=q_0$ and a barrier at a point
$q=q^{\ddag}$ that separates the well from a different well
or a continuum.  The theory of activated
rate processes describes the effect of damping and temperature on the
rate of escape of the particle.
If $V(q)\rightarrow\infty$ sufficiently fast as
$|q|\rightarrow\infty$, the Fokker-Planck operator has well
defined eigenfunctions, $W_n$, and eigenvalues, $-\mu_n$. The
time dependence is then exponential, with modes
$W_n(p,q,t)=e^{-\mu_nt}W_n(p,q)$.
The eigenfunction $W_0$, associated with the zero eigenvalue
$\mu_0=0$, is the thermal equilibrium
density $W_0=\exp\left[-\beta({p^2/2}+V(q)\right]$.
The standard methodology \cite{risken,gardiner,schuss} for obtaining
the rate is based on various approximations for nontrivial
density functions, use of variational principles such as
generalized Ritz principles for
the Fokker-Planck operator \cite{talkner,dekker},
and more recently, estimates of the mean first passage time to the stochastic
separatrix \cite{ryter1}-\cite{klosek}.
All these approaches were based on studying the properties of
the Fokker-Planck equation or the Langevin equation.

An alternate representation of the stochastic differential equation is
obtained by describing the particle dynamics in terms of a Hamiltonian
in which the particle is bilinearly coupled to a harmonic
bath \cite{zwanzig}. The equivalent Hamiltonian has the form
\begin{eqnarray}
H={p^2\over 2}+V(q)+\sum_j \left[{p_{x_j}^2\over 2}+{1\over 2}
(\omega_jx_j-{c_j\over \omega_j}q)^2\right],\label{Hamiltonian}
\end{eqnarray}
where the bath (mass weighted)
modes have coordinates $x_j$, momenta $p_j$, and
frequencies $\omega_j$.  Each bath oscillator is bilinearly coupled to
the system through the coupling coefficient $c_j$.
This Hamiltonian is related to the stochastic Langevin equation
by noting that the equation of motion of each bath oscillator
is just that of a forced oscillator in which the forcing
function is dependent only on the system coordinate $q$.
This allows for an explicit equation of motion for
the system coordinate which has the same form as the (generalized)
Langevin equation,
\begin{eqnarray}
{\ddot q}+{dV(q)\over \,dq}+\int^td\tau\gamma(t-\tau){\dot q}
(\tau)=\xi(t), \label{GLE}
\end{eqnarray}
in which the time dependent friction is
\begin{eqnarray}
\gamma(t)=\sum_j{c_j^2\over \omega_j^2}\cos(\omega_jt).
\label{gammat}
\end{eqnarray}
The random force, $\xi(t)$, depends on the initial positions
and momenta of all the bath variables.
If this initial distribution is chosen from the
canonical ensemble, $\exp(-\beta H)$, then the random force is
Gaussian and its autocorrelation function is related to the
time dependent friction by the fluctuation-dissipation theorem.
For the Langevin equation, the time dependent friction is Ohmic,
$\gamma(t)=2\gamma\delta(t)$.
Obtaining the Langevin dynamics from the Hamiltonian implies a
representation of the friction function as a periodic function,
whose cosine Fourier series is well defined, for example
$\gamma_{\tau}(t)=2\gamma\sum_{n=-\infty}^{n=\infty}\delta(t+n\tau)$,
and then taking an ever increasing period $\tau$.  For a given
period, $\tau$, the Hamiltonian representation
of the dynamics is accurate for times $t$ less than $\tau$.

Although this representation would seem at first glance to
cause enormous difficulties, it has been shown in recent
work \cite{pollak1,pollak2} that direct analysis of the Hamiltonian
dynamics can lead to vast simplifications. Specifically,
one can use classical variational transition state theory to
estimate the rate for finite barrier heights\cite{pollak1,pollak2},
cusped potentials
\cite{pollak2}, space dependent friction problems \cite{voth1},
and anisotropic friction in the multidimensional case \cite{bpz}. The
Hamiltonian dynamics has been used to derive a Kramers
turnover theory for the rate \cite{pgh},
valid for all values of the damping constant $\gamma$.
It has also been used as a starting point for a perturbation theory
expansion to derive exact finite barrier corrections
to the rate of reaction \cite{pollak3}.

What has not been done to date is to establish the relation
between the detailed dynamics of the Hamiltonian
representation on the one hand and the
dynamics as determined from the Fokker-Planck equation or
the stochastic differential equation on the other.
While it has been shown that the rates, when known,
are identical, it has not been shown for example how one
may derive the eigenfunctions of the Fokker-Planck operator
directly from the Hamiltonian representation.  Similarly,
one would think that there should be a way to relate
the stochastic separatrix and the Hamiltonian dynamics.

In this paper we start exploring the connection between the
two approaches. We deal
exclusively with the dynamics of a parabolic barrier potential, for
which everything is presumably known.  The Fokker-Planck equation
can be solved analytically \cite{risken} and
the Hamiltonian dynamics is separable \cite{pollak4}.
>From the point of view of each separate approach,
there is nothing new, however, the relationship between the \
two very different solutions is new and we believe interesting.
For example, the Hamiltonian approach can deal with memory
friction as easily as with Ohmic friction.
The Fokker-Planck equation is well defined only
for Ohmic friction.  However, from the Hamiltonian approach
we can derive a distribution function in the $q,p$ phase
space which decays exponentially for {\it any} time dependent
friction.  On a conceptual basis, the analysis presented in
this paper also indicates why the
variational transition state theory for the rate is so much more
powerful and useful than other methods.

In Section 2, we review briefly the Hamiltonian dynamics for
the parabolic barrier potential.  We then show in Section
3 how these dynamics may
be exploited for a direct construction of a nontrivial eigenfunction.
In Section 4, the ``Kramers" stationary flux distribution
function is derived, and in Section 5, we construct
the stochastic separatrix. We end with a discussion, indicating how the
present results may be used for further development of the theory of
activated rate processes.\\

\renewcommand{\theequation}{2.\arabic{equation}}
\setcounter{equation}{0}
\setcounter{section}{2}

\centerline{\bf 2. Hamiltonian dynamics for a parabolic barrier}

\vspace{5mm}
When the potential function is that of a parabolic barrier
\begin{eqnarray}
V(q)=-{1\over 2}\omega^{\ddag^2}q^2,\label{V}
\end{eqnarray}
where $\omega^{\ddag}$ is the barrier frequency, the Hamiltonian
(\ref{Hamiltonian})
is a quadratic form and thus may be cast in terms of
normal modes
\begin{eqnarray}
H={1\over 2}p_\rho^2-{1\over 2}\lambda^{\ddag^2}\rho^2+
\sum_j{1\over 2}\left[p_{y_j}^2+\lambda_j^2y_j^2\right].\label{Hnm}
\end{eqnarray}
Here, $\rho$ denotes the unstable normal mode, which is the analog of
the system coordinate $q$, and reduces to it in the zero damping limit.
Similarly, $y_j$ denotes the j-th stable normal mode.
The normal modes
are related to the original system and bath coordinates,
$q,{\underline x}$, via an orthogonal transformation matrix
${\bf U}$, such that
\begin{eqnarray}
q=u_{00}\rho+\sum_ju_{j0}y_j,\ \ \ \ \ \ \
p=u_{00}p_{\rho}+\sum_ju_{j0}p_{y_j}.\label{qu00}
\end{eqnarray}
The details of this transformation have been worked out in various
places \cite{pollak4,levine1,levine2}, so that
here we summarize the main results which will be of interest in the
next sections.

Denoting the Laplace transform by a ``hat", one may easily
see that the Laplace transform of the time dependent friction,
as represented in Eq.(\ref{gammat}), is
\begin{eqnarray}
{\hat \gamma}(s)=s\sum_j{c_j^2\over \omega_j^2}{1\over
\omega_j^2+s^2}.\label{gammal}
\end{eqnarray}
One finds that the unstable mode parabolic barrier frequency,
$\lambda^{\ddag}$, is just the usual Kramers-Grote-Hynes
frequency \cite{kramers,grote}
\begin{eqnarray}
{\lambda^{\ddag}\over \omega^{\ddag}}=
\left[1+{{\hat \gamma}(\lambda^{\ddag})\over
\lambda^{\ddag}}\right]^{-{1\over 2}}.\label{lambdag}
\end{eqnarray}
The transformation matrix element that gives the projection
of the system coordinate on the unstable mode may also
be expressed in terms of the Laplace transform of the
time dependent friction as \cite{pgh}
\begin{eqnarray}
u_{00}^2=\left[1+{1\over 2}({{\hat \gamma}(\lambda^{\ddag})\over
\lambda^{\ddag}}+{\partial{\hat \gamma}(s)\over {\partial}s}
\vert_{s=\lambda^{\ddag}})\right]^{-1}.  \label{u00}
\end{eqnarray}
One may define an effective bath frequency, $\Omega$,
as \cite{pollak5}
\begin{eqnarray}
\Omega^{-2}\equiv\sum_j{u_{j0}^2\over \lambda_j^2}=
{u_{00}^2\over \lambda^{\ddag^2}}-{1\over \omega^{\ddag^2}},
\label{Omega}
\end{eqnarray}
where the last equality on the right hand side of the equation may
be proved from properties of the normal mode
transformation \cite{levine2}.
Since the three quantities, $\lambda^{\ddag}$, $u_{00}$, and
$\Omega$, are all expressed in terms of Laplace transforms of
the time dependent friction, they are well defined in the
continuum limit. We also note that, as shown in \cite{pollak4},
the product over all normal mode frequencies is identical to
the product over all original frequencies,
\begin{eqnarray}
\lambda^{\ddag}\prod_j\lambda_j=\omega^{\ddag}\prod_j\omega_j.
\label{lambdap}
\end{eqnarray}

For the specific case of Ohmic friction,
\begin{eqnarray}
\gamma(t)=2\gamma\delta(t), \label{gammad}
\end{eqnarray}
which is the main case of interest in this paper, the explicit solution
for all the three quantities $\lambda^{\ddag}$, $u_{00}$, and
$\Omega$, is simple, specifically,
\begin{eqnarray}
u_{00}^{-2}=1+{\gamma\over 2\lambda^{\ddag}}\label{u00g}
\end{eqnarray}
and
\begin{eqnarray}
u_1^2\Omega^2=\omega^{\ddag^2}, \label{u1}
\end{eqnarray}
where we have used the notation
\begin{eqnarray}
u_1^2{\equiv}1-u_{00}^2.\label{u12}
\end{eqnarray}

Since the Hamiltonian in the normal mode form is separable,
the solution for the time dependence of each of the modes
in terms of the initial conditions (denote by a 0 subscript)
is simple.  For the unstable mode one has
\begin{eqnarray}
\rho(t)=\rho_0\cosh(\lambda^{\ddag}t)+{p_{\rho_0}\over \lambda^{\ddag}}
\sinh(\lambda^{\ddag}t)\label{rho}
\end{eqnarray}
and for the j-th stable mode
\begin{eqnarray}
y_j(t)=y_{j_0}\cos(\lambda_jt)+{p_{y_{j_0}}\over
\lambda_j}\sin(\lambda_jt) \label{yj}
\end{eqnarray}

In the normal mode form, the question of crossing the barrier
becomes trivial. If the energy in the unstable mode is greater
than zero and the momentum is positive, the trajectory proceeds
from reactants (negative $\rho$) to products (positive $\rho$).
However, the system is observed in the physical phase space $q,p$
and thus the question remains of how does the motion of
the normal modes reflect in the physical space. This is
dealt with in detail in the next section.\\

\renewcommand{\theequation}{3.\arabic{equation}}
\setcounter{equation}{0}
\setcounter{section}{3}

\centerline{\bf 3. An Eigenfunction}

\vspace{5mm}
As noted in the Introduction, the time dependence of the eigenfunctions
of the Fokker-Planck operator is that of a single exponent.  For the
parabolic barrier, the lowest nonzero eigenvalue is just the normal
mode barrier frequency $\lambda^{\ddag}$.  It is shown below that
this is not a mere coincidence. In fact, when considering the
normal mode representation, one understands that if
the initial conditions of the unstable mode are chosen such that
\begin{eqnarray}
\lambda^{\ddag}\rho_0-p_{\rho_0}=0 \label{lambdar}
\end{eqnarray}
then for any time t,
\begin{eqnarray}
\rho(t)=\rho_0e^{+\lambda^{\ddag}t}, \ \ \ \ \ \ \ \ p_{\rho}(t)=
p_{\rho_0}e^{+\lambda^{\ddag}t}.\label{rhot}
\end{eqnarray}
The time dependence of the physical coordinate and momentum
is on the average the same exponential time dependence as
implied by eqs. (\ref{qu00}). Thus one should expect
the eigenfunction to be represented by the projection onto
the physical phase space of all points which obey the fundamental
relation as in Eq.(\ref{lambdar}). In this section, we
quantify these ideas and show that they are correct.

A function $P(q,p,t)$ is defined as
\begin{eqnarray}
P({\tilde q},{\tilde p},t)=
\mbox{Tr}\left\{\delta(q(t)-{\tilde q})\delta(p(t)-{\tilde p})
\delta(\lambda^{\ddag}\rho_0-p_{\rho_0})e^{-\beta H}
\over \mbox{Tr}\left\{
 \delta(q-{\tilde q})\delta(p-{\tilde p})e^{-\beta H}
\right\}\right\}.  \label{P}
\end{eqnarray}
The time dependence of $q(t),p(t)$ is known exactly for the Hamiltonian
dynamics in terms of the initial conditions
and is given by a combination of eqs.(\ref{qu00}), (\ref{rho}), and
(\ref{yj}).  The trace implies an integration over all initial
coordinates and momenta and may be performed either in
the $q,{\underline x}$ representation, using the
form of the Hamiltonian given in Eq.(\ref{Hamiltonian})
or in the $\rho,{\underline y}$
representation using the normal mode Hamiltonian (\ref{Hnm}).

Instead of performing the trace operation over the initial
conditions, one may use the invariance of the volume element
in phase space to perform the integration over the volume
element at time $t$.
One must then express the initial conditions for the unstable mode
$\rho_0,p_{\rho_0}$
in terms of the unstable coordinate and momentum at time t.
By virtue of eqs.(\ref{lambdar}) and (\ref{rhot}) and properties
of the $\delta$ function, one finds that
\begin{eqnarray}
\delta(\lambda^{\ddag}\rho_0-p_{\rho_0})=e^{-\lambda^{\ddag}t}
\delta\left[\lambda^{\ddag}\rho(t)-p_{\rho}(t)\right].  \label{d}
\end{eqnarray}
It immediately follows that the time dependence of the distribution
is as designed,
\begin{eqnarray}
P(q,p,t)=e^{-\lambda^{\ddag}t}P(q,p,0).\label{Pp}
\end{eqnarray}
It is thus reasonable to expect that $P(q,p,0)$
is related to the eigenfunction associated with the eigenvalue
$-\lambda^{\ddag}$ of
the Fokker-Planck equation in the case of Ohmic friction.
It is just a matter of some algebra to demonstrate the explicit
relation.

The denominator in Eq.(\ref{P}) is easily evaluated in
the $q,{\underline x}$ representation, since the integration
over all the bath modes is just a simple Gaussian integration.
One finds that
\begin{eqnarray}
\mbox{Tr}\left\{ \delta(q-{\tilde q})\delta(p-{\tilde p})e^{-\beta H}
\right\}=
\prod_j\left({2\pi\over \beta\omega_j}\right)e^{-\beta {1\over 2}
({\tilde p}^2-\omega^{\ddag^2}{\tilde q}^2)}.\label{Trace}
\end{eqnarray}
For the numerator, it is convenient to use the normal
mode representation, this implies expressing the
system coordinate $q$ and momentum $p$ in terms of the
normal modes, as given in eqs.(\ref{qu00}).
The integrals may be reduced to simple Gaussian integration by using
the Fourier expansion of the $\delta$ function.  For example,
\begin{eqnarray}
\delta(q-{\tilde q})={1\over 2\pi}\int_{-\infty}^{\infty}d\kappa
e^{i\kappa(u_{00}\rho+\sum_ju_{j0}y_j-{\tilde q})}.
\label{dq}
\end{eqnarray}
Integration over all bath coordinates and momenta leaves one with
the interim result
\begin{eqnarray}
&&P({\tilde q},{\tilde p},0)=
{\beta\Omega\over {2\pi u_1}}\prod_j\left({\omega_j\over \lambda_j}
\right)\times\label{Pt}\\&&
 \int_{-\infty}^{\infty}\,d\rho \,dp_{\rho}\,
\delta(\lambda^{\ddag}\rho-p_{\rho})
\exp\left(-{\beta\over 2}\left[\Omega^2
{\lambda^{\ddag^2}\over \omega^{\ddag^2}}
(\rho-u_{00}{\tilde q}{\omega^{\ddag^2}\over
\lambda^{\ddag^2}})^2+{1\over u_1^2}(p_{\rho}-
u_{00}{\tilde p})^2\right]\right) .\nonumber
\end{eqnarray}
The remaining two integrations are now straightforward.
Using the identities (\ref{lambdag}, (\ref{Omega}),
(\ref{lambdap}), (\ref{u00g})-(\ref{u12}), and rearranging,
one finds
\begin{eqnarray}
P({\tilde q},{\tilde p},0)=({\beta\over 2\pi})^{1\over 2}
\left[{\lambda^{\ddag}\over u_{00}^2{\hat
\gamma}(\lambda^{\ddag})}\right]^{1\over 2}
\exp\left[-{\beta\over 2}{\lambda^{\ddag}\over
{\hat \gamma}(\lambda^{\ddag})}
({\tilde p}-{\omega^{\ddag^2}\over
\lambda^{\ddag}}{\tilde q})^2\right].\label{Ptq}
\end{eqnarray}

The eigenfunction of the Fokker-Planck operator is
the product of $P$ with the equilibrium distribution (see
\cite{ef})
\begin{eqnarray}
W(q,p)=P(q,p,0)e^{-{\beta\over 2}(p^2-\omega^{\ddag^2}q^2)}
\label{W1}
\end{eqnarray}
It is a matter of algebra to show that this function is indeed the
eigenfunction of a Fokker-Planck like operator,
\begin{eqnarray}
L=-{\partial\over\partial q}p+{\partial\over\partial p}
\left[{\hat \gamma}(\lambda^{\ddag}) p+{dV(q)\over \,dq}\right]
+{{\hat \gamma}
(\lambda^{\ddag})\over \beta}{\partial^2\over\partial p^2}
\label{L}
\end{eqnarray}

For Ohmic friction ${\hat \gamma}(\lambda^{\ddag})=\gamma$ and so
Eq.(\ref{W1}) is just the standard Fokker-Planck operator.
However, by utilizing the microscopic dynamics which is
valid also for memory friction, we find that we may define
a new operator, whose eigenfunction is
$W(q,p)$ with eigenvalue $-\lambda^{\ddag}$ which in
turn obeys the Kramers-Grote-Hynes equation (\ref{lambdag}).\\


\renewcommand{\theequation}{4.\arabic{equation}}
\setcounter{equation}{0}
\setcounter{section}{4}

\centerline{\bf 4.Kramers' stationary flux distribution function}

\vspace{5mm}
In his original paper, Kramers determined the escape rate
of the activated particle in the spatial diffusion limit by
finding a distribution function
$F(q,p)$ which has the property that the flux associated with
it is stationary and that it obeys the boundary condition
which is that deep in the reactants region it is just the
equilibrium distribution function. For the Fokker-Planck
operator,  Eq.(\ref{fpo}), the current operator has components
\begin{eqnarray}
J_q=p,\ \ \ \ J_p=-(\gamma p+{dV(q)\over \,dq})-{\gamma\over \beta}
{\partial \over \partial p}\label{Jq}
\end{eqnarray}
The stationarity condition takes the form
\begin{eqnarray}
{\partial F\over \partial t}=0=-\nabla\cdot{\underline J}F.
\label{F}
\end{eqnarray}
As shown by Kramers \cite{kramers}, for the parabolic barrier
\begin{eqnarray}
F(q,p)=e^{-{\beta\over 2}(p^2-\omega^{\ddag^2}q^2)}
{1\over \surd\pi}
\int_{-\infty}^{({\beta\lambda^{\ddag}\over 2\gamma})^{1\over 2}
(p-{\omega^{\ddag^2}\over \lambda^{\ddag}}q)}e^{-y^2}\,dy.
\label{Fqp}
\end{eqnarray}

When considering the Hamiltonian dynamics, the flux operator in
the full phase
space, in the normal mode representation, has the components
\begin{eqnarray}
J_{\rho}={\dot \rho},\ \ \ \ \,J_{p_{\rho}}={\dot p_{\rho}},
\label{Jr}
\end{eqnarray}
and similarly for each of the stable modes,
\begin{eqnarray}
J_{y_j}={\dot y_j},\ \ \ \ \ J_{p_{y_j}}={\dot p}_{y_j}.
\label{Jyj}
\end{eqnarray}
Here the dot denotes time differentiation, the time
dependence is determined from Hamilton's equations of motion.
The stationarity condition takes the form
\begin{eqnarray}
{\nabla}\cdot{\underline J}F=0,
\label{nabla}
\end{eqnarray}
where F is a distribution function and the gradient operator
is in the full phase space.

Of course, all this is just notation for the fact
that any function in phase space which is constant along a classical
trajectory will be stationary.  Thus for example
the Hamiltonian itself is stationary.  However, Kramers was interested
not only in a stationary flux, but he imposed the boundary
condition which is that the flux is associated with particles
which were initiated in the infinite past in the reactant region.
Following Pechukas \cite{pechukas}, we define the characteristic
function of points in phase space $\chi_r$ which is unity
on all phase space points of a
trajectory which was initiated in the infinite past at reactants
and is zero otherwise.  By definition, $\chi_r$
is stationary.  The distribution function associated with the
characteristic function $\chi_r$, projected onto the physical phase
space is then
\begin{eqnarray}
F_r({\tilde q},{\tilde p})=N\mbox{Tr}\left[\delta(q-{\tilde q})\delta(p-{\tilde
 p})
e^{-\beta H}\chi_r\right], \label{Fr}
\end{eqnarray}
where $N$ is a normalization constant and the Tr operation
is over the whole phase space.

For the purely parabolic barrier the dynamics
is trivial and so we know the explicit form of the
characteristic function
\begin{eqnarray}
\chi_r=\theta(p_{\rho}-\lambda^{\ddag}\rho).\label{chi}
\end{eqnarray}
where $\theta(x)$ is the unit step function. It is now a
matter of repeating the same algebra as in the previous
section to find that the desired function is simply
\begin{eqnarray}
F_r({\tilde q},{\tilde p})=e^{-{\beta\over 2}({\tilde p}^2
-\omega^{\ddag^2} {\tilde q}^2)}{1\over \surd\pi}
\int_{-\infty}^{\left[{\beta\over 2}
{\lambda^{\ddag}\over {\hat \gamma}(\lambda^{\ddag})}\right]^{1\over 2}
({\tilde p}-{\omega^{\ddag^2}\over
\lambda^{\ddag}}{\tilde q})}e^{-y^2}\,dy.\label{Frt}
\end{eqnarray}
It is easy to see that for Ohmic friction this result reduces
to Kramers function. Moreover, using the same generalization
of the Fokker-Planck operator to include memory friction as
in the previous section one finds that
\begin{eqnarray}
LF(q,p)=0,\label{LF}
\end{eqnarray}
where the operator $L$ is defined in Eq.(3.11).
It is the same Fokker-Planck operator as in Eq.(\ref{fpo}),
but generalized to include the memory friction case.

For the purely parabolic barrier, the eigenvalue $\mu=0$
is degenerate, as both the equilibrium and the stationary
flux distributions are eigenfunctions.
Of course the equilibrium distribution is a nonnormalizable
function.  In addition,
the equilibrium net flux
through any surface vanishes while the net flux associated with the
Kramers distribution is constant but not zero. As is well known, in the
presence of nonlinearities, the degeneracy is removed.\\

\renewcommand{\theequation}{5.\arabic{equation}}
\setcounter{equation}{0}
\setcounter{section}{5}

\centerline{\bf 5. The stochastic separatrix}

\vspace{5mm}
The stochastic separatrix is defined as a curve in the
particle phase space $(q,p)$ from which the particle
reaches products with a probability of
$1\over 2$ \cite{ryter1}-\cite{klosek}.
To construct the stochastic separatrix we, construct a
function $T({\tilde q},{\tilde p})$ which is defined as
the probability that a particle initiated at phase space point
${\tilde q},{\tilde p}$ reaches products in the limit that
the time goes to $\infty$.  The stochastic separatrix
is a curve in phase space defined by the condition
\begin{eqnarray}
T({\tilde q},{\tilde p})={1\over 2}.
\label{T}
\end{eqnarray}

To obtain the ``transmission probability function'' $T$,
one only needs to consider the dynamics of the unstable mode.
Positive values of $\rho$ are defined as
the products region, negative values are defined as reactants.
Denote the energy in the unstable mode as $E_{\rho}$,
\begin{eqnarray}
2E_{\rho}\equiv p_{\rho}^2-\lambda^{\ddag^2}\rho^2.
\label{2E}
\end{eqnarray}
Any trajectory initiated with $E_{\rho}>0, p_{\rho}>0$ necessarily
makes it to products as time tends to $\infty$. Similarly,
any trajectory initiated with $E_{\rho}<0, \rho>0$, makes it
to products at large times.
Barring the set of measure zero, $E_{\rho}=0$, no other
initial condition will lead to products as the time goes to
$\infty$.  Therefore the transmission
probability function is the projection of these two conditions
onto the physical phase space.  Explicitly,
\begin{eqnarray}
T({\tilde q},{\tilde p})=\mbox{Tr}\left\{\delta(q(t)-{\tilde q})
\delta(p(t)-{\tilde p})\left[\theta(E_{\rho})\theta(p_{\rho})+
\theta(-E_{\rho})\theta(\rho)\right]e^{-\beta H}
\over \mbox{Tr} \delta(q-{\tilde q})\delta(p-{\tilde p})e^{-\beta
H}\right\}. \label{Tt}
\end{eqnarray}
As in the
previous section, it is convenient to carry out the
integrations in the normal mode space.
The integration over all the bath variables is identical to that in the
previous section so that one can immediately write down the
interim result as in Eq.(\ref{Pt}), except that the $\delta$
function is replaced by the sum of the
products of step functions appearing in Eq.(\ref{Tt}). The
result is
\begin{eqnarray}
T({\tilde q},{\tilde p})=({\beta\over 2\pi})^{1\over 2}
({\lambda^{\ddag}\over {\hat \gamma}(\lambda^{\ddag})})^{1\over 2}
\int_{-({\tilde p}+{\omega^{\ddag^2}\over
\lambda^{\ddag}}{\tilde q})}^{\infty}
dve^{-{\beta\over 2}{\lambda^{\ddag}\over {\hat
\gamma}(\lambda^{\ddag})}v^2}.\label{Ttq}
\end{eqnarray}

>From the definition of the stochastic separatrix we find that
$T={1\over 2}$ on the line
\begin{eqnarray}
{\tilde p}_{ss}+{\omega^{\ddag^2}\over \lambda^{\ddag}}{\tilde
q}_{ss}=0.\label{T1}
\end{eqnarray}
This result is correct for arbitrary memory friction, since the
normal mode barrier frequency ($\lambda^{\ddag}$) is a solution
of the Kramers-Grote-Hynes equation (2.5).
In the Ohmic limit, this is identical to the results derived
in  \cite{klosek}.
Note again the ease with which we have derived the stochastic
separatrix for arbitrary memory friction.\\

\renewcommand{\theequation}{6.\arabic{equation}}
\setcounter{equation}{0}
\setcounter{section}{6}

\centerline{\bf 6. Discussion}

\vspace{5mm}
The Hamiltonian representation of the dynamics of a particle moving
on a parabolic barrier in the presence of time dependent
friction and an external Gaussian random force has been used to derive
stationary distribution functions and separatrices in the
physical ($q,p$) phase space. In the limit of Ohmic friction,
these become the known separatrices and eigenfunctions of
the Fokker-Planck operator. We have seen that knowledge of
the deterministic separatrix of the unstable normal
mode ($\rho$) is really all that counts.  The rest is
just projection of the normal mode dynamics onto the physical
phase space, taking into account
the correct thermal averaging over initial conditions.

Interestingly, the Hamiltonian approach for the parabolic barrier is
more efficient than the direct approach based on the
Fokker-Planck operator since it allows immediate extension
of all results to include
memory friction. This is not a mere artifact of the parabolic barrier
case.  An important aspect of the parabolic barrier is that it
provides a useful point of departure for a perturbation theory which
takes into account the nonparabolic part of the potential.  For example,
one may use the perturbation theory presented in \cite{pollak3} together with
the approach presented here, to obtain perturbation expansions of Kramers'
function and the stochastic separatrix even in the presence of memory and
space dependent friction.

The Hamiltonian approach for the parabolic
barrier, seems to suggest a useful generalization of the
Fokker-Planck operator to the memory friction case (cf Eq. 3.11).
This operator does not in general describe the full dynamics for the
memory friction case.  However, it does suggest that one can usefully
define Fokker-Planck like operators which provide a good representation
for specific properties.

Perhaps the most interesting result of the present paper
is the insight it gives into the dynamics underlying the
stationary distribution functions. We have seen that it
is necessary to understand the dynamics in all of phase
space, in order to construct the eigenfunctions.  Equivalently,
the stochastic separatrix may be found only if one knows
the fate of each trajectory.  This is a much more
severe restriction than the one used in bounding the rate
with the use of variational transition state theory (VTST),
In VTST, all that is needed is the instantaneous velocity of
each trajectory as it crosses a dividing surface.  This is
a much weaker demand leading to a substantially
more powerful theory of reaction rates.  Of course,
the weaker demand has its price, the VTST approach only
gives an upper bound to the rate.
However, experience has shown that often this upper bound is
very good.\\

\noindent
{\bf Acknowledgments:}\
We are grateful to Professor P. Talkner for many
stimulating discussions.  This work has been supported by a grant from
the Minerva foundation and the Einstein Center of the
Weizmann Institute of Science and by a grant from the
Foundation for Basic Reserch, administered by the Israel
Academy of Science.\\

\begin{thebibliography}{99}
\bibitem{kramers} H.A. Kramers, Physica {\bf 7}, 284 (1940).
\bibitem{risken} H.Risken, {\it The Fokker-Planck Equation} (Springer,
Belin, New York, 1989), 2nd ed., Springer Series in Synergetics, Vol. 18.
\bibitem{gardiner} C.W. Gardiner, {\it Handbook of Stochastic Methods}
(Springer, New York, 1983).
\bibitem{schuss}Z. Schuss, {\it Theory and Applications of
Stochastic Differential Equations}, Wiley, N.Y. 1980.
\bibitem{talkner} P. Talkner, Ber. Bunsenges. Phys. Chem. {\bf 95}, 327
(1991).
\bibitem{dekker} H. Dekker, Phys. Lett. {\bf 112A}, 197 (1985).
\bibitem{ryter1} D. Ryter, Physica A {\bf 142} 103 (1987).
\bibitem{ryter2} D. Ryter, J. Stat. Phys. {\bf 49}, 751 (1987).
\bibitem{klosek} M.M. K\l osek, B.J. Matkowsky and Z. Schuss, Ber. Bunsenges.
Phys. Chem. {\bf 95}, 331 (1991).
\bibitem{zwanzig} R. Zwanzig, J. Stat. Phys. {\bf 9}, 215 (1973).
\bibitem{pollak1} E. Pollak, S.C. Tucker and B.J. Berne, Phys. Rev. Lett.
{\bf 65}, 1399 (1990).
\bibitem{pollak2} E. Pollak, J. Chem. Phys. {\bf 93}, 1116 (1990).
\bibitem{voth1} G.R. Haynes, G.A. Voth and E. Pollak, Chem. Phys. Lett.,
in press.
\bibitem{bpz} A.M. Berezhkovskii, E. Pollak and V. Yu. Zitserman, J. Chem.
Phys. {\bf 97}, 2422 (1992).
\bibitem{pgh} E. Pollak, H.Grabert and P.H{\" a}nggi, J. Chem. Phys.
{\bf 91}, 4073 (1989).
\bibitem{pollak3} E. Pollak and P. Talkner, Phys. Rev. E {\bf 47}, 922 (1993).
\bibitem{pollak4} E. Pollak, J. Chem. Phys. {\bf 85}, 865 (1986).
\bibitem{levine1} A.M. Levine, M. Shapiro and E. Pollak, J. Chem. Phys.
{\bf 88}, 1959 (1988).
\bibitem{levine2} A.M. Levine, W. Hontscha and E. Pollak, Phys. Rev. B
{\bf 40}, 2138 (1989).
\bibitem{grote} R.F. Grote and J.T. Hynes, J. Chem. Phys. {\bf 73},
2715 (1980).
\bibitem{pollak5} E. Pollak, J. Chem. Phys. {\bf 95}, 533 (1991).
\bibitem{ef}Z. Schuss, (preprint).
\bibitem{pechukas} P. Pechukas, in {\it Dynamics of Molecular Collisions},
Part B,
edited by W.H. Miller, (Plenum Press, N.Y., 1976), p. 269.
\end{thebibliography}
\end{document}


