% The paper is formatted in LaTex Version 2.09.
% Figures are available upon request from the author.
% However, the first four figures are identical with
% the first four figures in ref.[39], Phys. Fluids,
% vol.6, 3063 (1994) and can be obtained there.
% The macro "seceq.sty" is included in the following
% file, which numbers equations sequentially by section.
% It should be included in the directory when formatting
% the LaTex file for the printer. If numbering by section
% is not desired, then the "seceq" may simply be omitted
% from the "documentstyle" command below.
%
%%%%%%%%%%%%%%%%%%%%%%%%%
%% seceq.sty macro %%
%%%%%%%%%%%%%%%%%%%%%%%%%
%
% ripped off by jonathan from euthesis.sty
%
% numbers equations consecutively within sections
% ****************************************
% * MISCELLANEOUS *
% ****************************************
%
% EQUATION and EQNARRAY -- put here because it must follow \chapter definition
%
% \newcounter{equation}
%
\@addtoreset{equation}{section} % Makes \chapter reset 'equation' counter.
\def\theequation{\thesection.\arabic{equation}}
% \jot = 3pt % Extra space added between lines of an eqnarray environment
% The macro \@eqnnum defines how equation numbers are to appear in equations.
%
\def\@eqnnum{(\theequation)}
%\def\@eqnnum{{\rm (\thesection .\theequation)}}
%
BODY
\documentstyle[11pt,seceq]{article}
\newcommand{\be}{\begin{equation}}
\newcommand{\ee}{\end{equation}}
\newcommand{\lb}{\label}
\newcommand{\en}{\epsilon}
\newcommand{\ven}{\varepsilon}
\newtheorem{Th}{Theorem}
\newtheorem{Lem}{Lemma}
\newtheorem{Conj}{Conjecture}
\newtheorem{Prop}{Proposition}
\newtheorem{Ex}{Example}
\newcommand{\bv}{{\bf v}}
\newcommand{\br}{{\bf r}}
\newcommand{\bk}{{\bf k}}
\newcommand{\bl}{{\bf l}}
\newcommand{\bh}{{\bf h}}
\newcommand{\bj}{{\bf j}}
\newcommand{\bm}{{\bf m}}
\newcommand{\bn}{{\bf n}}
\newcommand{\bpd}{{\bf \nabla}}
\newcommand{\bE}{\hat{{\bf e}}}
\newcommand{\bz}{{\bf 0}}
\newcommand{\bR}{{\bf R}}
\newcommand{\bP}{{\bf P}}
\newcommand{\bp}{{\bf p}}
\newcommand{\bD}{{\bf D}}
\newcommand{\bZ}{{\bf Z}}
\newcommand{\bV}{{\bf V}}
\newcommand{\bK}{{\bf K}}
\newcommand{\bT}{{\bf T}}
\newcommand{\bA}{{\bf A}}
\newcommand{\bU}{{\bf U}}
\newcommand{\bG}{{\bf G}}
\newcommand{\ef}{{\rm eff}}
\newcommand{\bC}{{\bf C}}
\newcommand{\bB}{{\bf B}}
\newcommand{\bL}{{\bf L}}
\newcommand{\bI}{{\bf I}}
\newcommand{\med}{{\rm med}}
\newcommand{\bu}{{\bf u}}
\newcommand{\bF}{{\bf f}}
\newcommand{\fs}{{\bf f}^\prime}
\newcommand{\fl}{\overline{{\bf f}}}
\newcommand{\fS}{f^\prime}
\newcommand{\fL}{\overline{f}}
\newcommand{\Fs}{{\bf F}^\prime}
\newcommand{\Fl}{\overline{{\bf F}}}
\newcommand{\BF}{{\bf F}}
\newcommand{\ba}{{\bf a}}
\newcommand{\bs}{{\bf s}}
\newcommand{\hX}{{\hat{X}}}
\newcommand{\hP}{{\hat{P}}}
\newcommand{\hH}{{\hat{H}}}
\newcommand{\hL}{{\hat{L}}}
\newcommand{\hbX}{{\hat{{\bf X}}}}
\newcommand{\hbP}{{\hat{{\bf P}}}}
\newcommand{\bb}{{\bf b}}
\newcommand{\bc}{{\bf c}}
\newcommand{\bx}{{\bf x}}
\newcommand{\bX}{{\bf X}}
\newcommand{\by}{{\bf y}}
\newcommand{\bw}{{\bf w}}
\newcommand{\grad}{{\mbox{\boldmath $\nabla$}}}
\newcommand{\bun}{{\mbox{\boldmath $1$}}}
\newcommand{\st}{$\,_{*}$}
\newcommand{\vl}{\overline{{\bf v}}}
\newcommand{\cvl}{\overline{{\bf V}}}
\newcommand{\taus}{{\tau^r}}
\newcommand{\taul}{{{\tau^s}}}
\newcommand{\vs}{{\bf v}^{\prime}}
\newcommand{\ws}{{\bf w}^{\prime}}
\newcommand{\Vs}{{\bf V}^{\prime}}
\newcommand{\Ws}{{\bf W}^{\prime}}
\newcommand{\ul}{\overline{{\bf u}}}
\newcommand{\Kl}{\overline{{\bf K}}}
\newcommand{\us}{{\bf u}^{\prime}}
\newcommand{\Ks}{{\bf K}^{\prime}}
\newcommand{\pl}{\overline{p}}
\newcommand{\ps}{p^\prime}
\newcommand{\zl}{\overline{z}}
\newcommand{\zs}{z^\prime}
\newcommand{\Zl}{\overline{Z}}
\newcommand{\Zs}{Z^\prime}
\newcommand{\Gl}{\overline{G}}
\newcommand{\Gs}{G^\prime}
\newcommand{\Bl}{\overline{b}}
\newcommand{\Bs}{b^\prime}
\newcommand{\BL}{\overline{B}}
\newcommand{\BS}{B^\prime}
\newcommand{\FL}{\overline{F}}
\newcommand{\FS}{F^\prime}
\newcommand{\Sigmal}{\overline{\Sigma}}
\newcommand{\Sigmas}{\Sigma^\prime}
\newcommand{\phil}{\overline{\phi}}
\newcommand{\Phil}{\overline{\Phi}}
\newcommand{\Phis}{\Phi^\prime}
\newcommand{\el}{\overline{e}}
\newcommand{\btau}{{\mbox{\boldmath $\tau$}}}
\newcommand{\bdot}{{\mbox{\boldmath $\cdot$}}}
\newcommand{\btimes}{{\mbox{\boldmath $\times$}}}
\newcommand{\bleta}{{\mbox{\boldmath $\eta$}}}
\newcommand{\btaul}{{\mbox{\boldmath $\tau$}}^s}
\newcommand{\btaus}{{\mbox{\boldmath $\tau$}}^r}
\newcommand{\bphi}{{\mbox{\boldmath $\phi$}}}
\newcommand{\ben}{{\mbox{\boldmath $\epsilon$}}}
\newcommand{\bsigma}{\overline{\mbox{\boldmath $\sigma$}}}
\textwidth6.25in
\textheight8.5in
\oddsidemargin.25in
\topmargin0in
\renewcommand{\baselinestretch}{1.7}
\begin{document}
\title{Turbulence Noise}
\author{Gregory L. Eyink\\{\em Department of Mathematics, University of Arizona}\\
{\em Building No. 89, Tuscon, AZ, 85721}}
\date{ }
\maketitle
\begin{abstract}
We show that the large-eddy motions in turbulent fluid flow obey a modified hydrodynamic equation with a stochastic
turbulent stress whose distribution is a causal functional of the large-scale velocity field itself. We do so
by means of an exact procedure of ``statistical filtering'' of the Navier-Stokes equations, which formally solves
the closure problem, and we discuss relations of our analysis with the ``decimation theory'' of Kraichnan. We show that
the statistical filtering procedure can be formulated using field-theoretic path-integral methods within the
Martin-Siggia-Rose formalism for classical statistical dynamics. We also establish within the MSR formalism a
``least-effective-action principle'' for mean turbulent velocity profiles, which generalizes Onsager's principle
of least dissipation. This minimum principle is a consequence of a simple realizability inequality and therefore holds
also in any realizable closure. Symanzik's theorem in field-theory---which characterizes the static effective action
as the minimum expected value of the quantum Hamiltonian over all state vectors with prescribed expectations of
fields---is extended to MSR theory with non-Hermitian Hamiltonian. This allows stationary mean velocity profiles and
other turbulence statistics to be calculated variationally by a Rayleigh-Ritz procedure. Finally, we develop
approximations of the exact Langevin equations for large eddies, e.g. a random-coupling DIA model, which yield new
stochastic LES models. These are compared with stochastic subgrid modelling schemes proposed by Rose, Chasnov, Leith,
and others, and various applications are discussed.
\newpage
{\em Key words:} Navier-Stokes, turbulence noise, LES, generalized Langevin equations, variational principles
{\em PACS numbers:} 47.27.Sd, 03.40.Gc, 47.25.Cg, 64.60.Ak
\end{abstract}
\section{Introduction}
The analogy of turbulent motion with the hydrodynamics of molecular fluids has played a central role
in attempts to understand and to model the dynamics of large, turbulent eddies. In ordinary Newtonian hydrodynamics,
at least in the low Mach number, incompressible regime, the local fluid velocity obeys an equation
of the form
\be \partial_t \bv+\grad\cdot(\bv\bv+\btaul)=-\grad p, \lb{1} \ee
in which $p$ is the (kinematic) pressure, required to enforce the incompressibility constraint
$\grad\cdot\bv=0,$ and $\btaul$ is the {\em viscous stress tensor}
\be \btaul= -\nu_0[(\grad\bv)+(\grad\bv)^\top]. \lb{2} \ee
The so-called {\em molecular viscosity} $\nu_0$ represents the residual, dissipative effects of the
``graininess'' of matter at the fluid level. Likewise, in the case of turbulent motion, the
``large-scale velocity'' $\bv_\ell,$ which may be conveniently defined by a smooth filtering
\be \vl_\ell(\br)\equiv \int d^d\br'\,\,G_\ell(\br-\br')\bv(\br'), \lb{3} \ee
obeys an equation of the form
\be \partial_t \vl_\ell+\grad\cdot(\vl_\ell\vl_\ell+\btau_\ell)=-\grad \pl_\ell, \lb{4} \ee
in which $\pl_\ell$ is the filtered pressure field and $\btau_\ell$ is a tensor representing the
{\em turbulent stress} of the eliminated small-scale eddies $<\ell.$ If the Reynolds number of the
fluid is high and $\ell$ is chosen in the long ``inertial range'' of scales, then the molecular
viscosity is believed to play a negligible role in the evolution of the large eddies $>\ell.$
Instead the primary effect is believed to be due to the smaller scale turbulence, which is modelled,
by analogy with the molecular fluids, as
\be \btau_\ell= -\nu_\ell[(\grad\vl_\ell)+(\grad\vl_\ell)^\top], \lb{5} \ee
in which $\nu_\ell$ is now a so-called {\em eddy viscosity}. This representation of
the effect of small-scales as a simple damping is, however, not nearly so accurate as
in the molecular case, where there is generally a large separation of scales between the
atomic and fluid degrees of freedom. In contrast, in the case of turbulence, there is no
such scale-separation and the eddy-viscosity representation is flawed, leading to a variety
of viscoelastic, nonlinear effects \cite{Kr76,Kr87}.
Here we shall be concerned with another over-simplification of the eddy-viscosity representation.
In fact, even in the case of molecular fluids there is an additional influence of the microscopic
degrees of freedom, which is a stochastic effect due to the chaotic molecular motions. Because of
such effects, Landau and Lifshitz \cite{LL} proposed in 1959 that there should be also in the Navier-Stokes fluid equation
Eq.(1) a {\em random stress} $\btaus,$ as
\be \partial_t \bv+\grad\cdot(\bv\bv+\btaul+\btaus)=-\grad p. \lb{6} \ee
Since the molecular degrees of the fluid are, by mixing properties of the dynamics, locally in a
state of thermal equilibrium, the random stress field $\btaus(\br,t)$ should have statistics which
are those of the universal Gibbs states. Motivated by such considerations, Landau and Lifshitz
hypothesized that $\btaus$ for the fluid at temperature $T$ should be a Gaussian random field with zero mean and covariance
\be \langle \tau^r_{ij}(\bx,t)\tau^r_{kl}(\by,s)\rangle=(2k_BT\nu_0/\rho)(\delta_{ik}\delta_{jl}+\delta_{il}\delta_{jk})
\delta^d(\bx-\by)\delta(t-s), \lb{7} \ee
where the delta-functions arise from the fast decay of correlations in the Gibbs states. This ``fluctuation-dissipation
relation'' (FDR) relating the strength of the noisy stresses and the molecular viscosity makes clear that the
two contributions to the momentum flux have the same origin in the microscopic degrees-of-freedom. Generally speaking,
the molecular noise gives a negligible effect on the gross fluid motions, but it does have observable
consequences for fluctuation behavior seen, for example, in light-scattering experiments.
By the same analogy with molecular fluids which motivated the eddy-viscosity model Eq.(\ref{5}), one therefore should expect
that the stress tensor in the equation of motion of large turbulent eddies will consist of both a systematic (or, deterministic)
and a fluctuating part:
\be \partial_t \vl_\ell+\grad\cdot(\vl_\ell\vl_\ell+\btaul_\ell+\btaus_\ell)=-\grad \pl_\ell. \lb{8} \ee
In fact, an exact formula for the turbulent stress can be deduced from the filtering procedure, as
\be \btau_\ell=\overline{(\bv\bv)}_\ell-\vl_\ell\vl_\ell. \lb{9} \ee
This representation of the stress is as a function of the {\em full} velocity field $\bv$, that is, of the small-scale
(or, sub-grid) component $\vs_\ell,$ defined as
\be \vs_\ell(\br)\equiv \bv(\br)-\vl_\ell(\br). \lb{10} \ee
as well as of the large-scale velocity $\vl_\ell$ itself. Therefore, one should expect that, when only the
large-scale field $\vl_\ell$ is specified, the stress $\btau_\ell$ will not be fully determined but that, in fact, a distribution
of values will occur with some given frequency depending upon the precise values of the small-scale velocity $\vs_\ell.$
Unlike the case of laminar fluids where molecular noise represents only a small perturbation to deterministic motion, in turbulent
flow the fluctuations are large and the importance of the random component of the stress, or {\em eddy noise}, is likely to be fully as
significant as the systematic part. On the other hand, if the statistics of the small-scale component $\vs_\ell$ are
{\em universal}, as usually believed, then the distribution law of the fluctuating stress ought to
be fixed by the Kolmogorov ``universal equilibrium state'' when $\ell$ is chosen in the inertial-interval. Nothing so
simple as the FDR in Eq.(\ref{7}) fixing the molecular noise can be expected, again due to lack of scale-separation, but
it ought to be possible to construct a generally applicable model of the statistical distribution.
The recognition of the probable importance of ``eddy noise'' in fully-developed turbulence does not originate with this
work but was pointed out already in the first cited paper of Kraichnan \cite{Kr76} and in the contemporary work of Rose
\cite{Rose}. More recently, a number of stochastic large-eddy models have been proposed to account for such effects
\cite{Lei,Chas,MT,Car}.
Our purpose here is to set up a general formalism for {\em statistical filtering} of the Navier-Stokes turbulence dynamics,
leading to exact ``fluctuating hydrodynamic equations'' for the large eddies of the form proposed. We should point
out that in the present method no finite ``stirring forces'' are added, as in the so-called ``RNG'' method of Yakhot
and Orszag \cite{YO}. Instead, the ``eddy noise'' we study must be dynamically generated by the chaotic fluid motion of the
small-scale eddies. In fact, we believe that a careful analysis of the existence and properties of this turbulent noise
must be made {\em before} any modelling is attempted, for otherwise serious confusions in the physics will result. This work
develops the framework for such an analysis. Our procedure is most similar to that used by Rose \cite{Rose}, but no
iteration schemes
motivated by renormalization group ideas will be employed (since we have not found them to present any significant
advantages). Although our stochastic equations are exact results for the Navier-Stokes dynamics, they are rather formal and
not directly useful for practical modelling. However, subsequent approximations based upon standard ideas (random-coupling
model, multiple-scale expansion and/or Markovianization, functional reversion techniques, etc.) yield simplified forms
which are suitable for numerical solution by computer.
The precise contents of this work are as follows: in the following Section 2 we shall set up the ``statistical filtering''
scheme leading to stochastic equations for the large turbulent eddies and discuss its statistical-dynamical foundations. We
compare there as well the relations of this method with Kraichnan's ``decimation theory'' \cite{Kr85}. His theory leads also
to generalized Langevin equations by a procedure of elimination or ``weeding'' of modes that are then modelled by Langevin
forces subject to realizability constraints. Afterward we develop from first principles the Martin-Siggia-Rose (MSR) field
theory formulation of classical statistical dynamics \cite{MSR}, emphasizing the non-perturbative and exact basis of the
approach. The MSR formalism is shown to yield an elegant and efficient formulation of the statistical filtering
scheme. In the next Section 3 we discuss a rather different---but closely related---subject of variational principles for
mean velocity profiles based upon the MSR effective action. It is shown that a principle of ``least effective-action''
characterizes the mean velocity as a consequence of a simple realizability inequality and thereby generalizes the famous
least-dissipation principle of Onsager \cite{O} for most probable states in nonequilibrium thermodynamics. The
least-action principle is made the basis to calculate mean turbulent profiles and other turbulence statistics
by variational strategies rather than by the solution of effective equations of motion. Finally, in the last
Section 4 we briefy describe approximations of the exact large-eddy equations, more practical for LES,
and we compare these models with those which have been introduced previously \cite{Lei,Chas,MT,Car}. We also
briefly discuss applications of such equations to the problem of atmospheric predictability, where eddy noise can
play a significant role in LES, and to simulation of inhomogeneous turbulence in a complex geometry.
\section{The Statistical Filtering Method and MSR Field Theory}
\noindent{\em (2.1) Formulation of Statistical Filtering}
Let us first recall the usual (deterministic) filtering method, of which an interesting recent account can be
found in the paper of Germano \cite{Ger1} (see also our paper \cite{I-LES}). It was already briefly reviewed in the
Introduction. The basic idea of the method is to convolute the turbulent velocity field with a smooth ``filter function''
$G_\ell(\br)=\ell^{-d}G(\br/\ell),$ as in Eq.(\ref{3}) there, in order to define a ``large-scale'' velocity $\vl_\ell.$
It is then possible to derive an infinite hierarchy of equations for ``generalized central moments,'' of which the
turbulent stress $\btau_\ell$ in Eq.(\ref{9}) is the 2nd order example. The lowest equation in the hierarchy is precisely
the Eq.(\ref{4}) in the Introduction. If the hierarchy is truncated at this order,
then one encounters the ``closure problem'' that $\btau_\ell$ is not a functional
of $\vl_\ell$ alone but also of $\vs_\ell$ and a similar problem occurs at every higher-order
truncation of the hierarchy as well. In the current large-eddy simulation (LES) literature, second-order modelling
is not attempted and only the lowest equation in the hierarchy is retained. Hence a model representation of $\btau_\ell$
must be sought as a (causal) functional of $\vl_\ell$ in order to obtain an autonomous equation for this field.
We shall show that a ``stochastic filtering'' method solves the closure problem in principle (but not in
practice) by yielding an {\em exact} representation of $\btau_\ell$ as
\be \btau_\ell[\vl_\ell]=\btaul_\ell[\vl_\ell]+\btaus_\ell[\vl_\ell], \lb{41} \ee
in which $\btaul_\ell[\vl_\ell]$ is a ``systematic stress,'' some causal functional of $\vl_\ell,$
and $\btaus_\ell[\bv_\ell]$ is a Langevin force, or ``random stress,'' whose distribution given the
past history of $\vl_\ell$ is specified. We shall refer to this exact Langevin dynamical equation
as the ``stochastic large-eddy equation'' (SLE). Although the functional $\btaul_\ell[\vl_\ell]$ and the
distribution of $\btaus_\ell[\vl_\ell]$ are calculable in principle from formulae given below,
in practice the exact expressions cannot be evaluated and approximation schemes must be employed. This
is the subject of the last Section 4.
Just as for the deterministic filtering procedure, a rather wide class of choices for filter function
$G$ may be made. The main requirement is that the Fourier transform $\widehat{G}(\bk)$ must decay
sufficiently rapidly as $k\rightarrow+\infty,$ so that the convolution reallys represents an ``elimination''
of high-wavenumber modes. The most common choices are the {\em Gaussian filter}, usually defined in the LES literature as
\be G(\bx)=\left({{6}\over{\pi}}\right)^{d/2}\exp\left[-{{6|\bx|^2}\over{\pi}}\right], \lb{42} \ee
the {\em tophat filter}
\be G(\bx)= \left\{ \begin{array}{ll}
1 & \mbox{ if $\max_i|x_i|<{{1}\over{2}}$} \cr
0 & \mbox{ otherwise,}
\end{array}
\right. \lb{43} \ee
and the {\em sharp Fourier cutoff filter}, which is most easily defined by its Fourier transform, as:
\be \widehat{G}(\bk)=\left\{ \begin{array}{ll}
1 & \mbox{ if $\max_i|k_i|<\pi$} \cr
0 & \mbox{ otherwise.}
\end{array}
\right. \lb{44} \ee
In physical space the latter gives
\be G(\bx)= 2^d\prod_{i=1}^d {{\sin\pi x_i}\over{x_i}}. \lb{45} \ee
However, as discussed at length in our work \cite{I-LES}, the sharp Fourier filter gives rise
to pathological features in the deterministic method and these can be seen to appear as well in the stochastic
formulation. Therefore, we will always consider a ``graded'' filter obeying the modest condition, Eq.(31) in
\cite{I-LES}, which allows us to avoid the bad features of the sharp Fourier cutoff.
The first steps in the ``statistical filtering'' method are the same as those in the usual deterministic scheme.
The ``large-scale (LS) velocity''
\be \vl_\ell(\br)=(G_\ell*\bv)(\br), \lb{45a} \ee
and ``small-scale (SS) velocity''
\be \vs_\ell(\br)=(H_\ell*\bv)(\br), \lb{45b} \ee
are introduced by convolution with low-pass filter $G_\ell$ and high-pass filter $H_\ell$ satisfying
\be \widehat{G}(\bk)+\widehat{H}(\bk)=1. \lb{45c} \ee
The latter guarantees that $\bv=\vl_\ell+\vs_\ell.$ If the filters are applied to the dynamical equations themselves,
then there result coupled equations:
\be \partial_t \vl_\ell+\grad\cdot(\vl_\ell\vl_\ell+\btau_\ell)=-\grad \pl_\ell+\nu_0\bigtriangleup\vl_\ell, \lb{40} \ee
which we call the ``large-eddy equation'' (LE), and
\be \partial_t \vs_\ell+\grad\cdot(\vl_\ell\vs_\ell+\vs_\ell\vl_\ell+\vs_\ell\vs_\ell-\btau_\ell)=
-\grad \ps_\ell+\nu_0\bigtriangleup\vs_\ell, \lb{40a} \ee
which we call the ``small-eddy equation'' (SE). Here $\btau_\ell$ is the ``turbulent stress-tensor,'' which is given
by an explicit quadratic function of the total velocity, $\btau_\ell=\bT_\ell(\bv,\bv),$
\be \bT_\ell(\bv,\bv)=\overline{(\bv\bv)}_\ell-\vl_\ell\vl_\ell. \lb{45d} \ee
A useful decomposition of this function was made by Germano \cite{Ger2}, according to the following straightforward
procedure. One substitutes $\bv=\vl_\ell+\vs_\ell$ into $\bT_\ell$ to obtain:
\be \bT_\ell(\bv,\bv)=\bL_\ell(\bv,\bv)+\bC_\ell(\bv,\bv)+\bR_\ell(\bv,\bv), \lb{112} \ee
where
\be \bL_\ell(\bv,\bv)\equiv \bT_\ell(\vl_\ell,\vl_\ell), \lb{113} \ee
is the so-called {\em Leonard stress},
\be \bC_\ell(\bv,\bv)\equiv \bT_\ell(\vl_\ell,\vs_\ell)+\bT_\ell(\vs_\ell,\vl_\ell), \lb{114} \ee
is the {\em cross stress}, and
\be \bR_\ell(\bv,\bv)\equiv\bT_\ell(\vs_\ell,\vs_\ell), \lb{115} \ee
is denoted as {\em Reynolds stress}. Note this is a modification of the traditional decomposition of Leonard \cite{Leo}.
It is clear that the LE and SE are coupled equations for the fields $\vl_\ell$ and $\vs_\ell.$
The fundamental new step in ``statistical filtering'' is to add to the SE a weak random force $\fs_\ell$:
\be \partial_t \vs_\ell+\grad\cdot(\vl_\ell\vs_\ell+\vs_\ell\vl_\ell+\vs_\ell\vs_\ell-\btau_\ell)=
-\grad \ps_\ell+\nu_0\bigtriangleup\vs_\ell+\fs_\ell, \lb{45e} \ee
where $\fs_\ell$ has zero mean and spectral support at wavenumbers greater than $\sim 2\pi/\ell.$
For convenience, we choose $\fs_\ell$ Gaussian. A possible choice of the covariance is
\be \langle \fs_\ell(\br t)\fs_\ell(\br' t')\rangle= 2\ven_0 H_\ell(\br-\br')\delta(t-t'). \lb{45f} \ee
The parameter $\ven_0$, when multiplied by the density $\rho$ of the fluid, has units of energy per time.
The role of this force is simply to add some random perturbations to the otherwise deterministic SE
and thus $\rho\ven_0$ should be very small compared to the total energy transfer per time to the small-scale
motions. In fact, $\ven_0$ will momentarily be taken to zero. First, we formally solve the SE part of the dynamics
{\em with the large-scale velocity field considered fixed.} \footnote{ This is similar to Kraichnan's ``clamping method''
in Section 8 of \cite{Kr85}.} The solution may be written as
\be \vs_\ell(\br t)=\Vs_\ell[\br t;\vl_\ell,\fs_\ell], \lb{45g} \ee
in which $\Vs[\br,t; \vl_\ell,\fs_\ell]$ is a causal functional of time-histories of $\vl_\ell$ and $\fs_\ell,$
i.e a functional depending only upon their values for times $From a more physical perspective, one can consider the random force to represent the molecular noise Eq.(\ref{7}).
Just as the molecular viscosity is believed to play no direct role in the dynamics of the large-scale turbulence, so that
the limit $\nu_0\rightarrow 0$ is well-defined, likewise the molecular noise is believed to play no essential role in
that limit except to select the correct measure. This idea was already proposed some time ago by Ruelle \cite{Rue2}. In the
language of field-theory, the bare, molecular noise will be replaced by a ``renormalized'' eddy noise. This is due to the
same physics by means of which molecular viscosity is replaced by an effective or ``renormalized'' eddy viscosity. The
field-theoretic point of view will be exploited later on to come up with approximation schemes and calculational methods
for ``statistical filtering.''
Stochastic descriptions similar to the SLE have been postulated by Hohenberg and Shraiman \cite{HoSh} in the context
of general spatiotemporally chaotic systems and, specifically, for the Kuramoto-Sivashinsky
equation. In the latter example there is an important conjecture of Yakhot \cite{Yak} that the chaotic behavior
generated by small-scale instabilities will produce an ``effective dynamics'' of the large-scales which is just
the noisy Burgers equation. See also \cite{Zal,SKJJB}. Our method generalizes and systematizes this earlier work, which was
based upon weak-coupling perturbation expansions (without a small parameter!). As we shall see shortly, the
``stochastic filtering'' we have formulated non-perturbatively may also be carried out in weak-coupling
expansions and then contains the earlier results. In fact, even for weak-coupling our systematic technique
obtains terms missed in the earlier work.
The ideas proposed here have also some similarities with the ``decimation theory'' of Kraichnan \cite{Kr85}, but
they are really essentially distinct. Kraichnan proposed a general strategy for economical computation of many-mode
systems by a procedure in which modes that are redundant due to underlying statistical symmetries are eliminated
and replaced by suitable Langevin forces. The latter are constructed to enforce certain statistical constraints
imposed by the exact dynamics. In this way efficient numerical approximations might be obtained and, by adding
more and more constraints from the true dynamics, convergence of the approximations to the true values might even
be found. The two methods have in common the elimination of modes with ``generalized Langevin equations'' as the output.
However, the ``statistical filtering'' we have formulated is exact and without any approximation whatsoever.
Of course, the output SLE are, at this stage, completely formal and useless as a tool for numerical computations.
Approximations must be developed for any practical application. Nevertheless, the idea is not at its basis an attempt
to develop {\em statistical approximations}, as is decimation, but is rather an attempt to characterize a physical
{\em statistical law}. This is seen most clearly by comparing the two procedures for ``chaotic'' vs. ``integrable''
systems. Decimation theory might be successfully applied to both cases, yielding generalized Langevin models as
statistical approximations of the true dynamics. In fact, Kraichnan's first example of his decimation method \cite{Kr85}
was for the three-mode dynamical model
\be \partial_t x= A_x yz,\,\,\,\,\,\partial_t y= A_y xz\,\,\,\,\,\partial_t z= A_z xy, \lb{45i} \ee
with $A_x+A_y+A_z=0.$ This dynamics has two quadratic integrals and is thus integrable by quadratures in terms of
elliptic functions, as discussed long ago by Lorenz \cite{Lor}. On the contrary, if the ``stochastic filtering'' were
applied to an integrable system, such as the above three-mode model or Burgers' equation, then the SLE for the explicit
modes would degenerate, or become deterministic, in the vanishing noise limit $\ven_0\rightarrow 0.$ The
chaotic properties of the dynamics are necessary to create a true statistical law as described by the SLE equations.
After these somewhat philosophical remarks, we return to the formal development. It is useful to separate
the turbulent stress $\btau_\ell$ into a ``systematic'' part $\btaul_\ell$ and a ``random'' part
$\btaus_\ell,$
\be \btau_\ell[x;\vl_\ell]=\btaul_\ell[x;\vl_\ell]+\btaus_\ell[x;\vl_\ell]. \lb{45j} \ee
(We represent here the space-time point as $x=(\br,t).$ ) Any such division is essentially arbitrary and
can be made according to various schemes. The one we have found most convenient is as follows: consider
the functional $\Vs_\ell[x;\vl_\ell,\fs_\ell]$ obtained by solving the SE for fixed past histories of $\vl_\ell$
and $\fs_\ell.$ Substituted into the expression for $\bT_\ell$, this can then be averaged over the ensemble of random forces
$\fs_\ell$ {\em treating } $\vl_\ell$ {\em as a deterministic quantity} (or, equivalently, as an independent random
variable). Denoting this average as $\langle\cdot |\vl_\ell\rangle,$ \footnote{Note that for any functional $F$ of
the large-scale velocity alone
\be \langle F[\vl_\ell]G[\vl_\ell,\fs_\ell]|\vl_\ell\rangle= F[\vl_\ell]\langle G[\vl_\ell,\fs_\ell]|\vl_\ell\rangle.
\lb{45m} \ee
Hence, the average $\langle\cdot|\vl_\ell\rangle$ has the appearance of a ``conditional average'' with large-scale
velocity fixed, as its notation also suggests. However, this is not true and rather misleading. The {\em solution}
$\vl_\ell$ {\em of the SLE} for a prescribed past history of random forces $\fs_\ell$ will develop a functional dependence
upon them, so that
\be \langle F[\vl_\ell]G[\vl_\ell,\fs_\ell]\rangle\neq F[\vl_\ell]\langle G[\vl_\ell,\fs_\ell] \rangle, \lb{45n} \ee
when $\langle\cdot\rangle$ denotes the ordinary average over the ensemble of forces $\fs_\ell.$ It would
be possible to consider a true conditional average over $\fs_\ell$ with the {\em solutions} $\vl_\ell$ fixed, but
such a conditioning would change the statistics of the forces $\fs_\ell$ from their {\em a priori} Gaussian
statistics prescribed by Eq.(\ref{45f}). Therefore, the property Eq.(\ref{45m}) is really a consequence just
of the {\em definition} of the linear operation $\langle\cdot |\vl_\ell\rangle$ applied to arbitrary functionals
$F[\vl_\ell,\fs_\ell]$ of the variables $\vl_\ell$ and $\fs_\ell.$ Its usefulness will appear in the context of the MSR
formalism developed later.}
we then set
\be \btaul_\ell[x;\vl_\ell]=\langle \bT_\ell[x;\vl_\ell,\cdot]|\vl_\ell\rangle, \lb{45k} \ee
which no longer depends upon $\fs_\ell,$ and
\be \btaus_\ell[x;\bv_\ell,\fs_\ell]=\bT_\ell[x;\vl_\ell,\fs_\ell]-\btaul[x;\vl_\ell]. \lb{45l} \ee
Note also that the ``conditional mean'' of the SS field,
\be \vs[x;\vl]\equiv \langle \Vs[x;\vl]|\vl\rangle, \lb{72} \ee
will not vanish in general, since the past history of the LS field will set up a SS flow. Therefore, it is
natural to separate this mean contribution from the SS field, by making the definition
\be \Ws[x;\bv]\equiv \Vs[x;\bv]-\vs[x;\vl]. \lb{73} \ee
We can now define a net ``systematic field,''
\be \bv[x;\vl]\equiv \vl(x)+\vs[x;\vl] \lb{74} \ee
as a deterministic functional of the LS field, whereas $\Ws[\vl]$ represents the random part.
As a final remark, let us note that it is possible in principle to ``defilter'' the solutions of the SLE if
\be \fs_\ell=H_\ell*\bF, \lb{45o} \ee
where $\bF$ is a random force whose spectral support is disjoint from the support of $\widehat{G}_\ell(\bk).$
In that case, there is a solution of the full Navier-Stokes equation with the force $\bF$ which, when filtered
by $G_\ell,$ is a solution of the SLE given above. When the filter functions are graded this is not true for the
explicit choice of $\fs_\ell$ suggested in Eq.(\ref{45f}) above. The latter corresponds to filtering with $H_\ell$
the spacetime white-noise force $\bF$ whose covariance is
\be \langle \bF(\br t)\bF(\br' t')\rangle= 2\ven_0 \delta^3(\br-\br')\delta(t-t'). \lb{45p} \ee
Since $\widehat{G}_\ell$ and $\widehat{H}_\ell$ must have overlapping supports for graded filters
and the spatially white-noise force $\bF$ has a flat wavenumber spectrum, a corresponding term $\fl_\ell=G_\ell*\bF$
must then appear in the LE. However, since one does not believe that the results in the limit $\ven_0\rightarrow 0$
will depend too essentially upon the precise form of the force, we expect the same results by using the
``non-defilterable'' force Eq.(\ref{45f}) as a ``defilterable'' one. The possibility of using different types of
forcing in the limit $\ven_0\rightarrow 0$ limit, with hopefully equivalent results, will be exploited further later on.
\noindent{\em (2.2) Martin-Siggia-Rose Formalism for Statistical Dynamics}
We review here the field-theory method of MSR \cite{MSR}. Its fundamental motivation, without technicalities,
is the observation that it is difficult to characterize the probability distribution describing a
stationary ensemble of turbulent flows. In that case, it may be more profitable to work instead with the
ensemble of {\em histories} of the Navier-Stokes dynamics. This can be specified more concretely
and, in many respects, provides an effective substitute for the stationary measure.
To explain the formal principles, we consider first the example of a randomly forced Navier-Stokes fluid. In this case,
the statistical problem to be solved is given by the dynamical equation
\be \partial_t\bv +\bP(\grad)(\bv\cdot\grad)\bv=\nu_0\bigtriangleup \bv+\bF, \lb{11} \ee
where $\bP(\grad)$ is a solenoidal projection required to maintain the incompressibility
condition (replacing a pressure term) and $\bF$ is a random solenoidal force with a Gaussian
distribution and covariance
\be \langle f_i(\bx,t)f_j(\bx',t')\rangle = P_{ij}(\grad_\bx)F(\bx-\bx')\delta(t-t'). \lb{12} \ee
As we shall discuss later, the stochastic nature of the force $\bF$ is not really necessary to our
discussion, or any restriction of the method. The statistical problem posed by Eq.(\ref{11}) may also
be considerably generalized by employing other classes of random forces besides
those Gaussian and white-noise in time. The models with random forcing subsume the problems with random initial data,
since those may be represented by an impulsive random force
\be \bF(\br t)=\bv_0(\br)\delta(t-t_0) \lb{13} \ee
imposed at the initial instant $t_0$ on a quiescent fluid, with the initial datum $\bv_0$ chosen from
some selected random ensemble. Although we do not need to assume that $\langle \bF(\br t)\rangle=0,$
it is more convenient to assume this and to add, if desired, an explicit mean force $\bar{\bF}(\br t)$
to the righthand side of Eq.(\ref{11}). In such cases
\be \bar{\bv}(\br t)\equiv \langle \bv(\br t)\rangle \lb{13a} \ee
will not be zero. If the forcing spectrum $\widehat{F}(\bk)$ is supported only in a small interval near $k_0=1/L,$
then Eqs.(\ref{11}),(\ref{12}) are a relatively realistic model of homogeneous and isotropic turbulence, where the force
is just a convenient way of injecting energy at large scales, i.e. stirring the fluid at length-scale $L.$
Rather than following the somewhat mysterious operator formulation originally developed by MSR \cite{MSR} (see also Phythian
\cite{Phy1,Phy2}), we shall briefly describe the path-integral formulation of Janssen \cite{Jan} and DeDominicis \cite{DeD}.
The simple idea is to write a path-integral representation for the generating functional $Z[\bleta,\hat{\bleta}]$
of correlation and response functions by incorporating the dynamics through a delta functional in its
representation by an integral over an exponential. The objects playing the role of ``momentum'' $p$
in the functional integral analogue of $\delta(x)={{1}\over{2\pi}}\int dp\,e^{ipx},$ are the
fields $\hat{\bv}(\br t),$ whose joint correlations with $\bv$'s turn out to be the response functions.
The expression for the generating functional is just
\begin{eqnarray}
Z[\bleta,\hat{\bleta}]&=&\int{\cal D}\bF \exp\left[-{{1}\over{2}}\langle \bF, F^{-1}\bF\rangle \right] \cr
\,& &\times \int{\cal D}\bv{\cal D}\hat{\bv}\exp\left[ i\int d^d\br\int dt\hat{\bv}(\br t)
\left( (\partial_t-\nu_0\bigtriangleup_\br)\bv(\br t) \right. \right.\cr
\,& &\,\,\,\,\,\,\,\,\,\,\,\left.+\bP(\grad_\br)(\bv(\br t)\cdot \grad_\br)\bv(\br t)-\bar{\bF}(\br t)\right) \cr
\,& &\,\,\,\,\,\,\,\,\,\,\,\left.-i\langle \hat{\bv},\bF\rangle -i\int d^d\br dt \left(\bleta(\br t)\cdot\bv(\br t)
+\hat{\bleta}(\br t)\cdot\hat{\bv}(\br t)\right) \right]J[\bv], \lb{14}
\end{eqnarray}
where $J[\bv]$ is a Jacobian factor which appears in the transformation
\be \delta[\bv-{\bf V}(\bF)]=\delta[(\partial_t-\nu_0\bigtriangleup)\bv-\bar{\bF}+\bP(\grad)(\bv\cdot\grad)\bv-\bF]J[\bv], \lb{15} \ee
and ${\bf V}(\bF)$ is the exact solution of the dynamical equation Eq.(\ref{11}) with specified
force $\bF.$
We note that a discretization of the dynamics in time will always be implicitly assumed in this work to give meaning
to such formulae as those above, even where we use continuum notation as a convenient shorthand. Along with this,
a space-regularization by means of a spatial grid or a Fourier-Galerkin truncation will be employed. Among various
choices of time-discretization, we could use the explicit Euler scheme,
\be (\bv_i-\bv_{i-1})/\tau=\bK(\bv_{i-1})+\bF_i, \lb{15a} \ee
or an implicit numerical scheme, such as:
\be (\bv_i-\bv_{i-1})/\tau=\bK({{\bv_i+\bv_{i-1}}\over{2}})+\bF_i. \lb{15aa} \ee
It is easy to check that the form of the Jacobian $J[\bv]$ and the resulting path-integral expression depend upon the
discretization adopted. In fact, there is a substantial literature on this point \cite{Gr1,Gr2,Wei,Wis}. If the force in
our equation Eq.(\ref{11}) were multiplied by some state-dependent coefficient $\bB(\bv),$ then it is known that the
symmetrical time-splitting as in the implicit rule Eq.(\ref{15aa}) must be used in the noise coefficient in order to yield
the continuum equation in Stratonovich form \cite{Arn}. For a stochastic Langevin equation the Stratonovich interpretation
is usually implicitly assumed since it guarantees validity of the familiar rules of calculus. On the other hand, the
explicit or ``non-anticipating'' scheme leads to the Ito form. We calculate here the Jacobian only for these two
schemes. Beginning with the symmetrical time-splitting, we note if ${\bf V}_i(\bv_{i-1},\bF_i)$ is the
solution of the discretized nonlinear equation Eq.(\ref{15aa}) for $\bv_i$ in terms of $\bv_{i-1}$ and $\bF_i$, it follows
from
\begin{eqnarray}
\,& &\prod_{i}\delta\left((\bv_i-\bv_{i-1})/\tau-\bK({{\bv_i+\bv_{i-1}}\over{2}})-\bF_i\right) \cr
\,& &\,\,\,\,\,\,\,\,\,\,\,\,=\prod_{i}\delta\left(\bv_i-{\bf V}_i(\bv_{i-1},\bF_i)\right)
\left| {{1}\over{\tau}}-{{1}\over{2}}{{\delta\bK}\over{\delta\bv_i}}\left({{\bv_i+\bv_{i-1}}\over{2}}\right) \right|^{-1} \cr
\,& &\,\,\,\,\,\,\,\,\,\,\,\,\approx \prod_{i}\delta\left(\bv_i-{\bf V}_i(\bv_{i-1},\bF_i)\right)
{{1}\over{\tau}}\exp\left[{{1}\over{2}}\tau\cdot {{\delta\bK}\over{\delta\bv_i}}\left({{\bv_i+\bv_{i-1}}\over{2}}\right)\right], \lb{16}
\end{eqnarray}
that the Jacobian factor is
\be J[\bv]=\exp\left[-{{1}\over{2}}\int d^d\br dt\,\, {\rm tr}{{\delta \bK(\br t;\bv)}\over{\delta\bv(\br t)}}\right].
\lb{16a} \ee
On the other hand, it is clear by the same method that $J[\bv]\equiv 1$ for the Ito discretization. We should remark
that for situations where the nonlinear part of the dynamics satisfies a Liouville theorem, the Jacobian in fact
is only a field-independent factor and may always be neglected. It was shown long ago by T. D. Lee \cite{Lee} that the
Fourier-Galerkin truncation of the Navier-Stokes dynamics falls in this category. We shall rederive this result
in Appendix II, where it is shown more generally that the ``effective dynamics'' generated by elimination of the
small-scale modes satisfies such a Liouville property. It thus turns out that in all the cases we consider the nonlinear
dynamics satisfies the Liouville theorem and we are justified in ignoring the Jacobian factor.
Performing finally the Gaussian integration over force $\bF$ by completion of squares, we obtain
\be Z[\bleta,\hat{\bleta}]=\int{\cal D}\bv{\cal D}\hat{\bv}\exp\left(iS[\bv,\hat{\bv}]-i\langle\bv,\bleta\rangle
-i\langle\hat{\bv},\hat{\bleta}\rangle\right), \lb{17} \ee
a path-integral over histories with the ``action''
\begin{eqnarray}
S[\bv,\hat{\bv}] &= & \int d^d\br\int dt \hat{\bv}(\br t)\left( (\partial_t-\nu_0\bigtriangleup_\br)\bv(\br t) \right. \cr
\,& &\,\,\,\,\,\,\,\,\,\,\,\,\left.+\bP(\grad_\br)(\bv(\br t)\cdot \grad_\br)\bv(\br t)-\bar{\bF}(\br t)\right) \cr
\,& &\,\,\,\,\,\,\,\,\,\,\,\,+{{i}\over{2}}\int d^d\br dt\int d^d\br'dt' \hat{v}_i(\br t)F_{ij}(\br t,\br't')\hat{v}_j(\br't'). \lb{18}
\end{eqnarray}
This expression may be used to calculate arbitrary multitime correlations of the velocity field by functional
differentiation with respect to $\bleta(\br t).$ The source $\hat{\bleta}(\br t),$ on the other hand, appears in the
action exactly as an external force in the Eq.(\ref{11}), so that the derivatives ${{\delta^{p+1}Z[\bleta,\hat{\bleta}]}\over
{\delta\eta_i(\br t)\delta\hat{\eta}_{j_1}(\br_1t_1)\cdots\delta\hat{\eta}_{j_p}(\br_pt_p)}}$ are the higher-order average response functions.
It is not difficult to extend the formalism to non-Gaussian or ``colored'' random forces and to forces multiplicative
in the random velocity $\bv.$ Indeed, one can consider any ``generalized Langevin dynamics''
\be \partial_t\bv(x)={\bf K}[x;\bv]+\bF'[x;\bv]. \lb{19} \ee
Here $x=(\br,t)$ is a space-time point and ${\bf K}[x;\bv]$ is an arbitrary functional of $\bv$ which is``non-anticipating''
or ``causal,'' i.e. that depends only upon the ``past'' values of $\bv(x')$ with $t'0$
\be -i\langle\Omega^-|{\rm T}\left[\hX_n(t)\hP_m(0)\right]|\Omega^+\rangle
=\langle {{\partial x_n(t)}\over{\partial x_m(0)}}\rangle. \lb{V56} \ee
Now making a small perturbation $\bx(0)\rightarrow \bx(0)+\ben$ in the initial data is the same as
making a small perturbation $\bF(t)\rightarrow \ben\cdot\delta(t)$ in the external force. By the definition
of functional differentiation it follows that
\be {{\partial x_n(t)}\over{\partial x_m(0)}}=
{{\delta x_n(t)}\over{\delta f_m(0)}}, \lb{V57} \ee
the usual instantaneous response function $\hat{G}_{nm}(t,0)$. The stated identity Eq.(\ref{V53}) immediately follows.
\noindent {\em (3.3) Symanzik-Type Theorems for MSR}
We can now state and prove the version of Symanzik's theorem which holds in MSR field theory. It requires a very
simple modification associated to the non-self-adjoint character of the formal Hamiltonian. More precisely, we have
\begin{Th}
The effective potential
\be V[\bx]=\lim_{T\rightarrow +\infty}{{1}\over{T}}\Gamma[\bx_T], \lb{V58} \ee
for a stationary Markov process is the value at the extremum point of the functional
\be V[\Psi^+,\Psi^-]=-\langle\Psi^-|\hL|\Psi^+\rangle \lb{59a} \ee
varying over all pairs of state vectors $\Psi^+,\Psi^-$ subject to the constraints
\be \langle\Psi^-|\Psi^+\rangle=1 \lb{V60} \ee
and
\be \langle\Psi^-|\hbX|\Psi^+\rangle=\bx. \lb{V61} \ee
\end{Th}
Whereas the original version of the theorem required just one trial state, now there must be {\em two independent
trial states.}
Nevertheless, the proof is similar to the original one of Symanzik \cite{Sym}. Observe first that the
generating functional $W[\bh]$ introduced earlier may be represented in the operator formulation by
\be W[\bh]=\log\langle\Omega^-| {\rm T}\exp\left(\int_0^T dt\,\hL_\bh(t)\right)|\Omega^+\rangle, \lb{V62} \ee
where
\be \hL_\bh(t)=\hL+\bh(t)\bdot\hbX . \lb{V63} \ee
No time-dependence is required for the coordinate operators because the exponential factors
automatically introduce the correct Heisenberg picture operators after differentiating and setting $\bh$ to zero.
We note then that for a {\em static} field $\bh$ in the limit $T\rightarrow +\infty,$
\begin{eqnarray}
\exp(W[\bh_T]) & = & \langle\Omega^-|\exp\left(T\cdot \hL_\bh \right)|\Omega^+\rangle \cr
\,& \approx & \langle\Omega^-|\Omega^+[\bh]\rangle\langle\Omega^-[\bh]|\Omega^+\rangle
\times \exp(T\cdot\lambda[\bh]), \lb{V64}
\end{eqnarray}
where $\lambda[\bh]$ is the eigenvalue of the ``perturbed operator''
\be \hL_\bh=\hL+\bh\bdot\hbX \lb{V65} \ee
with the {\em largest real part} and $\Omega^+[\bh],\Omega^-[\bh]$ are the associated right and left
``ground state'' eigenvectors:
\be \hL_\bh|\Omega^+[\bh]\rangle=\lambda[\bh]|\Omega^+[\bh]\rangle, \lb{V65a} \ee
and
\be \hL_\bh^*|\Omega^-[\bh]\rangle=\lambda^*[\bh]|\Omega^-[\bh]\rangle. \lb{V66} \ee
Furthermore, we can see that
\be {{\partial W[\bh_T]}\over{\partial h_n}}= T\cdot x_n[\bh]+o(T), \lb{V66a} \ee
with
\be x_n[\bh]= \langle\Omega^-[\bh]|\hX_n|\Omega^+[\bh]\rangle. \lb{V67} \ee
This can be obtained from the formula
\begin{eqnarray}
\exp(W[\bh_T]){{\partial W[\bh_T]}\over{\partial h_n}}
& = & \langle\Omega^-|{{\partial}\over{\partial h_n}}\exp\left(T\cdot \hL_\bh \right)|\Omega^+\rangle \cr
\,& = & \langle\Omega^-|\Omega^+[\bh]\rangle\langle\Omega^-[\bh]|\Omega^+\rangle
\langle\Omega^-[\bh]|{{\partial}\over{\partial h_n}}\exp\left(T\cdot \hL_\bh \right)|\Omega^+[\bh]\rangle \cr
\,& & \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,
+O\left(e^{-T\cdot\Delta\lambda}\right), \lb{V67a}
\end{eqnarray}
where $\Delta\lambda$ is the spectral gap between the real parts of the ``ground state'' eigenvalue and
the next highest eigenvalue. We also use the well-known fact that, for any one-parameter family
of operators $\hL(h)$ depending smoothly on a parameter $h,$
\be {{\partial}\over{\partial h}}\exp(\hL(h))=\exp(\hL(h))\varphi(-{\rm Ad}\hL(h))
\left[{{\partial \hL(h)}\over{\partial h}}\right], \lb{V68} \ee
where ${\rm Ad}\hL$ denotes the ``adjoint super-operator'' defined by commutator,
\be ({\rm Ad}\hL)[\hat{O}]=[\hL,\hat{O}], \lb{V69} \ee
and $\varphi(z)$ is the entire function
\be \varphi(z)=(e^z-1)/z=1+{{1}\over{2!}}z+{{1}\over{3!}}z^2\cdots. \lb{V70} \ee
See \cite{HS}. Since
\be \langle\Omega^-[\bh]|[\hL_\bh,\hat{O}]|\Omega^+[\bh]\rangle=0, \lb{V71} \ee
for any operator $\hat{O},$ only the first term survives in the expansion of $\varphi$
when substituted into the first term of formula Eq.(\ref{V67a}). This yields Eq.(\ref{V66a}).
Now let us consider the variational problem. If we incorporate the constraints by suitable
Lagrange multipliers, then the variational equation is just
\be \delta\left[ -\langle\Psi^-|\hL|\Psi^+\rangle-\bh\bdot\langle\Psi^-|\hbX|\Psi^+\rangle
+\lambda\langle\Psi^-|\Psi^+\rangle\right]=0, \lb{V72} \ee
or
\be \langle\delta\Psi^-|\hL_\bh-\lambda|\Psi^+\rangle
+\langle\Psi^-|\hL_\bh-\lambda|\delta\Psi^+\rangle=0. \lb{V72c} \ee
In other words, there are infinitely many stationary points of the functional $V[\Psi^+,\Psi^-]$ subject
to the constraints. They consist precisely of pairs $(\Psi_\alpha^+[\bh],\Psi_\alpha^-[\bh])$ of
eigenvectors of $\hL_\bh$,
\be \hL_\bh|\Psi_\alpha^+[\bh]\rangle=\lambda_\alpha[\bh]|\Psi^+_\alpha[\bh]\rangle, \lb{V72a} \ee
and
\be \hL_\bh^*|\Psi_\alpha^-[\bh]\rangle=\lambda^*_\alpha[\bh]|\Psi_\alpha^-[\bh]\rangle, \lb{V72b} \ee
corresponding to different branches of eigenvalues $\lambda_\alpha[\bh],\alpha=0,1,2,...$ To be precise, we should
consider the stationary point corresponding to the branch with largest real part for each $\bh,$ that is, the pair of
``ground state'' eigenvectors $(\Omega^+[\bh],\Omega^-[\bh])$ introduced above. Applying the left eigenvector to
the eigen-equation of the right vector and using the constraints gives
\be \langle\Omega^-[\bh]|\hL|\Omega^+[\bh]\rangle+ \bh\bdot\bx[\bh]=\lambda[\bh], \lb{V73} \ee
and thus
\begin{eqnarray}
-\langle\Omega^-[\bh]|\hL|\Omega^+[\bh]\rangle & = & \bh\bdot\bx[\bh]-\lambda[\bh] \cr
\, & = & {{1}\over{T}}\left[ \langle\bh_T,{{\delta W}\over{\delta\bh}}[\bh_T]\rangle-W[\bh_T]\right]+o(1), \cr
\, & = & {{1}\over{T}}\Gamma[\bx_T]+o(1). \lb{V74}
\end{eqnarray}
The first quantity is independent of $T,$ so that we see taking the limit $T\rightarrow +\infty$ that
\be -\langle\Omega^-[\bh]|\hL|\Omega^+[\bh]\rangle = V[\bx], \lb{V75} \ee
as was claimed.
We have given only a formal proof of the theorem without a careful statement of the conditions,
which would certainly involve spectral properties of the ``Liouville operator'' $\hL,$ etc. For
deterministic dynamics the existence of a spectral gap in (properly speaking) the Perron-Frobenius operator has
been established only for a few special cases, such as the work of Pollicot and Ruelle on Axiom A systems.
See \cite{Rue3}. The assumption of a spectral gap is probably stronger than required and a fast
polynomial decay, rather than exponential, should suffice. We just make one
remark here on the mathematical aspects, which is that trial states $\Psi^+,\Psi^-$ clearly must be taken
from different spaces. In fact, $\Psi^+$ varies over the space $L^1,$ of integrable functions of coordinates,
while $\Psi^-$ varies over the space $L^\infty$ of bounded functions. Since $L^\infty$ is the Banach space
dual to $L^1,$ the Dirac ``bra-ket'' notation, $\langle\cdot|\cdot\rangle$ above, must be interpreted
as the canonical dual space action of $L^\infty$ vectors on $L^1$ vectors. Later on we will formulate another
more symmetrical version of the result in which both $\Psi^+,\Psi^-$ vary over $L^2.$
The theorem can also be generalized in two important ways. First, there is an analogous variational characterization
of the full, time-dependent effective action $\Gamma[\bx].$ Let us state the result formally as
\begin{Th}
The effective action $\Gamma[\bx]$ for a stationary Markov process is the value at the extremum point of the functional
\be \Gamma[\Psi^+,\Psi^-]= \int^{+\infty}_{-\infty}dt\,\,\langle\Psi^-(t),(\partial_t-\hL)\Psi^+(t)\rangle,
\lb{V25} \ee
when that is independently varied over all pairs of time-dependent state vectors subject to the constraints
for each time $t$:
\be \langle\Psi^-(t),\Psi^+(t)\rangle=1 \lb{V26} \ee
and
\be \langle\Psi^-(t),\hbX\Psi^+(t)\rangle=\bx(t), \lb{V27} \ee
and also to the boundary conditions
\be \lim_{t\rightarrow\mp\infty}|\Psi^\pm(t)\rangle=|\Omega^\pm\rangle. \lb{76} \ee
\end{Th}
We shall not give here the proof of this theorem because it is essentially the same as the proof of a corresponding result
in quantum field theory due to Jackiw and Kerman \cite{JK}. They have shown similarly that the effective
action $\Gamma[\phi]$ of field-theory is the stationary point of the functional
\be \Gamma[\Psi^+,\Psi^-]= \int^{+\infty}_{-\infty}dt\,\,\langle\Psi^-(t),(ih\!\!\!\!^-\partial_t-\hat{H})\Psi^+(t)\rangle,
\lb{V77} \ee
varied over pairs $(\Psi^-(t),\Psi^+(t))$ with constraints
\be \langle\Psi^-(t),\Psi^+(t)\rangle=1 \lb{V78} \ee
and
\be \langle\Psi^-(t),\hat{\Phi}(\br)\Psi^+(t)\rangle=\phi(\br,t). \lb{V79} \ee
Just as the Symanzik theorem is a constrained version of the familiar quantum variational principle for
energy eigenvalues and eigenvectors, the Jackiw-Kerman theorem can be seen as a constrained version of
Dirac's \cite{Dir} variational formulation of the Schr\"{o}dinger equation (a quantum analogue of
Hamilton's principle). According to that principle, the solutions of the Schr\"{o}dinger equation
\be ih\!\!\!\!^-\partial_t\Psi^{\pm}(t)=\hat{H}\Psi^{\pm}(t), \lb{V80} \ee
are obtained as the stationary points of $\Gamma[\Psi^+,\Psi^-]$ subject to independent variations of the
two wavevectors. In addition to providing a basis for time-dependent Rayleigh-Ritz calculations, the Jackiw-Kerman-type
theorem establishes the existence of a Lagrangian functional for the effective action.
A second generalization of Symanzik's theorem has been made by Cornwall, Jackiw, and Tomboulis (CJT)\cite{CJT}
to the static effective action for higher-order statistics. For example, if $\Gamma[\phi,G]$ is the
quantum analogue of $\Gamma[\bv,\bU]$ defined previously, then CJT defined a static version appropriate to discuss
those 2nd-order statistics in the time-invariant ground state. Their definition corresponds to substituting
into $\Gamma[\phi,G]$ time-invariant versions $\phi(\br)$ and $G(\br,t-t';\br',0),$ with this $G$ then eliminated
in terms of the equal-time correlation $G(\br,\br')=G(\br,0;\br',0)$ in a way discussed in detail in \cite{CJT}.
After this, the result was divided by length of the total time interval $T$ and the limit $T\rightarrow +\infty$
taken. Then, CJT showed that the resulting quantity $V[\phi,G]$ is also given by the following
variational prescription:
\be V[\phi,G]=\langle\Psi,\hat{H}\Psi\rangle \lb{V29} \ee
at the minimum point subject to constraints
\be \langle\Psi,\Psi\rangle=1,\langle\Psi,\hat{\Phi}(\br)\Psi\rangle=\phi(\br), \lb{V30} \ee
as before and also
\be \langle\Psi,\hat{\Phi}(\br)\hat{\Phi}(\br')\Psi\rangle=\phi(\br)\phi(\br')+G(\br,\br'). \lb{V31} \ee
This is the straightforward generalization of Symanzik's original result. It can also be generalized to
MSR field theory in a way which should be now obvious. We will refrain here from giving the precise
statement and, instead, turn to the practical implementation of the results.
\newpage
\noindent {\em (3.4) Formulation of Rayleigh-Ritz Variational Methods}
We will just outline here a simple variational method of Rayleigh-Ritz type to approximate the
effective action and, thereby, the ensemble means. To initiate the method a {\em trial weight}
must be chosen,
\be w(\bx)\geq 0,\,\,\,\,\,\,\,\,\,\,\int d\bx\,\,w(\bx)=1, \lb{V81} \ee
as a plausible {\em a priori} guess for the density $\rho(\bx)$ of stationary measure. This
weight will contain a number of parameters which can be optimized by means of the variational
principle in Theorem 1. The most straightforward implementation uses expansions in orthogonal polynomials
with respect to the trial weight,
\be \int P_n(\bx)P_m(\bx) \,w(\bx)\,d\bx=\delta_{n,m}, \lb{V82} \ee
which are assumed to form a complete set in the weighted $L^2$-space associated to $w.$ General convergence
properties of these polynomials are discussed by Kraichnan in \cite{Kr70} and Section 2 of \cite{Kr85}. In the present
application, the pair of state vectors to be varied are expanded to some finite order $N$ as
\be \Psi^+= w\cdot \sum_{n=0}^{N-1} c^+_n P_n, \lb{V83} \ee
and
\be \Psi^-= \sum_{n=0}^{N-1} c^-_n P_n. \lb{V84} \ee
Then, incorporating constraints by Lagrange multipliers, the variational equation becomes
\be \delta\left[\sum_{n,m=0}^{N-1}c_n^- c_m^+(L_{nm}+\bh\bdot\bX_{nm})-\lambda\sum_{n=0}^{N-1}c_n^-c_n^+\right]=0,
\lb{V85} \ee
where for all $n,m=0,1,2,....$
\be L_{nm}=\langle P_n,\hL (wP_m)\rangle,\,\,\,\,\,\,\,\,\,\,\bX_{nm}=\langle P_n,\hbX (wP_m)\rangle. \lb{V86} \ee
The variation may be instituted in two stages. Varying first with respect to the expansion coefficients $\bc^+,\bc^-,$
the $N$-rank eigenvalue equations are obtained
\be \sum_{m=0}^{N-1}(L_{nm}+\bh\bdot\bX_{nm})c_m^+=\lambda c_n^+, \lb{V87} \ee
for $n=0,1,...,N-1,$ and
\be \sum_{n=0}^{N-1}(L_{nm}+\bh\bdot\bX_{nm})c_n^- =\lambda c_m^- \lb{V88} \ee
for $m=0,1,...,N-1.$ Then calculate, for each $\bh,$ the eigenvalue of largest real part and the associated
eigenvectors. With $\bx$ fixed, the value of $\bh$ is determined, yielding values $\lambda_N(\bx),\bc^\pm_N(\bx).$
The intermediate approximation of rank $N$ to the effective potential $V_N(\bx)$ is
\be V_N(\bx)= -\sum_{n,m=0}^{N-1}\,c_{N,n}^-(\bx)c_{N,m}^+(\bx)L_{nm}. \lb{V89} \ee
However, this potential may next be optimized with respect to the remaining parameters in the weight function
itself, giving the final approximation.
A natural way to implement this procedure in calculating mean turbulent velocity profiles would be to
make a Gaussian ansatz for the velocity fluctuations. The mean velocities in the Gaussian trial weight
would be taken as the variational parameters whereas the velocity covariance could be fixed by hypothesis. The covariance
could be determined, for example, from a model energy spectrum which is $k^{-5/3}$ in the inertial range. In this way,
the physics of the K41 theory could be incorporated into the trial state. The orthogonal polynomials with such a Gaussian
trial weight are just appropriate multidimensional Hermite polynomials, or polynomials ``Wick-ordered'' with respect to the
model covariance. The procedure outlined above would produce a set of variational equations for the mean velocity field at
every stage of expansion of state vectors into Hermite polynomials up to the $N$th order. This scheme is essentially a
first-order closure method for mean quantities by postulating 2nd-order statistics but with the advantage that it
may be systematically improved by increasing $N.$ Alternatively, 2nd-order closures could be implemented variationally
as well in which the entire velocity covariance (or some partial 2nd-order statistics such as a local value of the
Kolmogorov constant, or the local turbulence intensity, etc.) would be allowed to range over some trial set and
then varied to optimize the guess. Clearly, similar methods could also be used in time-dependent problems, such as
turbulence growth under a mean constant shear, or freely-decaying turbulence, etc. based upon the variational
characterization of the full time-dependent effective action.
We do not present any such realistic applications at this point. We shall only examine here a couple of very simple, exactly
soluble stochastic models, in order to allow easy comparison of the variational method with reality. First, let us consider
the standard Ornstein-Uhlenbeck process (cf. \cite{Ris}) associated to the linear Langevin equation
\be \partial_t x= -\gamma (x-x_0)+ \sqrt{D}\eta, \lb{V90} \ee
with $\langle\eta(t)\eta(t')\rangle=2\delta(t-t'),$ and Fokker-Planck operator
\be \hL= \gamma{{\partial}\over{\partial x}}\left( (x-x_0)\cdot\right)+ D{{\partial^2}\over{\partial x^2}}. \lb{V91} \ee
As is well-known, the stationary density is just a Gaussian
\be \rho(x)= {{1}\over{\sqrt{2\pi\sigma^2}}}\exp[-(x-x_0)^2/2\sigma^2], \lb{V92} \ee
with $\sigma^2=D/\gamma.$ In particular, the mean position is $\bar{x}=x_0.$ However, suppose we ignore this fact
and instead use a reasonable guess of the 2nd-order statistics to try to obtain the mean value variationally. A natural
choice of trial weight in this case is also Gaussian
\be w(x)= {{1}\over{\sqrt{2\pi\sigma^2}}}\exp[-(x-a)^2/2\sigma^2], \lb{V93} \ee
with arbitrary center $a.$ Furthermore, let us make the simple {\em ansatz} that
\be \Psi^+(x)= w(x)\cdot c_0^+, \lb{V94} \ee
and
\be \Psi^-(x)= 1+ c_1^-(x-a). \lb{V95} \ee
Thus, we use a ``mixed-order'' Hermite expansion in which the series for $\Psi^+$ is taken only to 0th-order but that
for $\Psi^-$ is carried to 1st-order. It is then easy to calculate that
\be \langle\Psi^-,\hL\Psi^+\rangle= -\gamma(a-x_0)\cdot c_1^-c_0^+, \lb{V96} \ee
that
\be \langle\Psi^-,\Psi^+\rangle= c_0^+, \lb{V97} \ee
and that
\be \langle\Psi^-,\hX\Psi^+\rangle= a\cdot c_0^++\sigma^2 \cdot c_1^-c_0^+ . \lb{V98} \ee
Imposing the constraints $\langle\Psi^-,\Psi^+\rangle=1$ and $\langle\Psi^-,\hX\Psi^+\rangle=x,$ it is easy to determine that
\be c_0^+=1 \,\,\,\,\,\,\,\,\,\, c_1^-={{x-a}\over{\sigma^2}}, \lb{V99} \ee
so that
\be \langle\Psi^-,\hL\Psi^+\rangle= -{{\gamma^2}\over{D}}(a-x_0)(x-a). \lb{V100} \ee
Taking $V(x)= -\langle\Psi^-,\hL\Psi^+\rangle$ and varying with respect to $a,$ the value $a={{x+x_0}\over{2}}$ is
obtained. Substituting this stationary point value, the final estimate is obtained:
\be V_{{\rm approx}}(x)= {{\gamma^2}\over{4D}}(x-x_0)^2. \lb{V101} \ee
Clearly, the minimum is just $x_0,$ so that the method has produced the correct mean value. Of course, this
example is very simple and getting the exact mean value at first try is not expected in general. Notice
that the same approximation may be obtained from the ``classical'' action for the process,
\be \Gamma_{OM}[x]={{1}\over{4D}}\int_{-\infty}^{+\infty}dt\,[\dot{x}+\gamma(x-x_0)]^2, \lb{V102} \ee
i.e. the Onsager-Machlup action. For a general diffusion process this is only the exact effective action in the zero-noise
limit $D\rightarrow 0,$ and otherwise it will be just a ``leading'' term in a diagrammatic loop-expansion. It is therefore
similar to the ``Hartree-Fock approximation'' used in \cite{CJT}. However, since
\be {{\Gamma_{OM}[x_T]}\over{T}}={{\gamma^2}\over{4D}}(x-x_0)^2 \lb{V102a} \ee
identically at each finite $T,$ the previous result for $V$ is recovered.
It may be worthwhile giving just one more simple example, to see whether the method can succeed when the
exact statistics are {\em non-Gaussian}. To that end, we take an exactly soluble 1-dimensional gradient dynamics
\be \partial_t x= -\varphi'(x)+ \sqrt{D}\eta, \lb{V103} \ee
with potential
\be \varphi(x)={{\lambda}\over{4!}}(x-x_0)^4. \lb{V104} \ee
As is well-known, the Fokker-Planck operator is now
\be \hL= {{\lambda}\over{6}}{{\partial}\over{\partial x}}\left((x-x_0)^3\cdot\right)
+ D{{\partial^2}\over{\partial x^2}}. \lb{V105} \ee
and the stationary density is just given by
\be \rho(x)\propto \exp\left[-{{\lambda\cdot(x-x_0)^4}\over{4!\cdot D}}\right]. \lb{V106} \ee
Although the fluctuations are non-Gaussian, the mean is still given simply as $\bar{x}=x_0.$
Let us use the same Gaussian trial guess as above for $w$ and the same {\em ansatz} for $\Psi^+,\Psi^-$
as in Eqs.(\ref{V94}),(\ref{V95}). Then Eqs.(\ref{V97}) and (\ref{V98}) are unchanged, while
\begin{eqnarray}
-\langle\Psi^-,\hL\Psi^+\rangle & = & \lambda\left((a-x_0)\cdot{{1}\over{2}}\sigma^2+{{1}\over{6}}(a-x_0)^3\right)
c_1^-c_0^+, \cr
\,& = &{{\lambda}\over{\sigma^2}}
\left((a-x_0)\cdot{{1}\over{2}}\sigma^2+{{1}\over{6}}(a-x_0)^3\right)(x-a), \lb{V107}
\end{eqnarray}
imposing the constraints. The approximate effective potential $V$ can be obtained by substituting the value $a$ at
the stationary point, which is determined by the variational equation $\partial V_a(x)/\partial a=0,$ or
\be \left[ {{1}\over{2}}\sigma^2+{{1}\over{2}}(a-x_0)^2\right](x-a)=
\left[(a-x_0)\cdot{{1}\over{2}}\sigma^2+{{1}\over{6}}(a-x_0)^3\right]. \lb{V108} \ee
It is possible to solve this cubic polynomial equation for $a,$ but it is easier to determine the approximate
mean value without doing so explicitly. In fact, using Eq.(\ref{V108}), it is straightforward to see that
\be V_{{\rm approx}}(x)={{\lambda}\over{\sigma^2}}\left({{1}\over{2}}\sigma^2+{{1}\over{2}}(a-x_0)^2\right)(x-a)^2
\geq 0. \lb{V109} \ee
Thus, the minimum $V_{{\rm approx}}=0$ is achieved if and only if $x=a.$ Substituting back into the
cubic Eq.(\ref{V108}), it is found that
\be {{1}\over{6}}(a-x_0)^3+{{\sigma^2}\over{2}}(a-x_0)=0 \lb{V110} \ee at the minimum. This equation has only one real-valued solution $a=x_0.$ Hence, we conclude that
\be \bar{x}=x_0, \lb{V111} \ee
which is again the exact value. Notice that the variance $\sigma^2$ of the Gaussian trial state never
needed to be specified in this argument. In fact, the Gaussian {\em ansatz} gives a rather poor representation
of the fluctuations for any choice of $\sigma^2,$ but all choices lead here to the same mean value.
Let us make a few remarks on the convergence properties of the systematic expansion procedure
described previously in the limit as $N\rightarrow +\infty.$ Although we give here no rigorous proofs,
it is clear that it should give a convergent scheme if some reasonable properties are satisfied by the
dynamics and the trial weight. A natural condition is that the operator $\hL_w$ defined by the following
formal similarity transformation
\be \hL_w\equiv {{1}\over{\sqrt{w}}}\hL\sqrt{w}, \lb{V112} \ee
be a (generally unbounded) operator on $L^2$ such that
\be \{\sqrt{w}P_n: n=0,1,2,...\}\subset {\rm dom}(\hL_w)\bigcap{\rm dom}(\hL_w^*), \lb{V113} \ee
and that the perturbed operators $\hL_w+\bh\bdot\hbX$ have a complete set of biorthogonal left and
right eigenvectors $\Psi_{w,\alpha}^\pm[\bh],\,\,\,\alpha=0,1,2,...$ in $L^2.$ Notice that if we take
$\Omega_w^\pm[\bh]$ to be the ``ground state'' vectors ($\alpha=0$) with largest real part, then
$\Omega_w^+[\bz]=\rho/\sqrt{w}$ and $\Omega_w^-[\bz]=\sqrt{w}$. \footnote{It follows from this that,
for $w=\rho,$ the two ground states at $\bh=\bz$ of the transformed operator coincide: $\Omega^\pm_\rho=\sqrt{\rho}.$
However, this is probably not helpful, since one must expect the left and right ground states to split
under the perturbation $\bh\bdot\hbX.$} In particular, the first implies the usual condition for orthogonal
expansion methods that
\be \int d\bx\,{{\rho^2(\bx)}\over{w(\bx)}}<\infty. \lb{V114} \ee
The importance of the similarity transformation is that it removes the asymmetry of the problem, replacing the
two different spaces of trial vectors, $L^1$ and $L^\infty,$ by the single Hilbert space $L^2.$ The variational
method may then be reformulated equivalently in terms of the functional
\be V_w[\Psi^+,\Psi^-]= -\langle\Psi^-,\hL_w\Psi^+\rangle, \lb{V115} \ee
varied subject to the same constraints as before but now with both $\Psi^\pm\in L^2.$ It is clear that the
Rayleigh-Ritz method outlined above is equivalent to one for the present principle if the trial states are taken as
\be \Psi^\pm(\bx)=\sqrt{w}\sum_{n=0}^{N-1} c_n^\pm P_n(\bx). \lb{V116} \ee
This reformulation is necessary to justify the procedure, since the trial states in Eqs.(\ref{V83}),(\ref{V84})
do not belong to the proper spaces $L^1,L^\infty.$ It also allows some proofs of convergence of
the $N$th-order approximation $V_N$ pointwise to the true effective potential $V,$ under suitable hypotheses. Note that
we consider here convergence with a fixed choice of trial weight $w,$ although one would expect convergence
to be improved by optimizing the weight at each stage.
Let us just close this section by comparing the variational method we have proposed with some previous ones.
Kraichnan (see Section 4.3 of \cite{Kr58}, also \cite{Kr79}) and Qian \cite{Qi} have devised variational
schemes, which involve satisfying equations of motion in mean-square sense. Instead, we use only the linear ``Liouville
operator'' of the dynamics. More essentially, the scheme advocated by Kraichnan was intended as a scheme to
approximate the entire statistical distribution, whereas the principle discussed here is constructed to
obtain just mean-values or other low-order statistics. The calculations of Qian \cite{Qi} do show that the
$5/3$-spectrum and a reasonable value of the Kolmogorov constant can be obtained by a single-time variational
calculation. There may be a connection of our ideas with those of Castaing \cite{Cas}. The action-principle we
have established here differs from the ``optimum theory'' of Busse \cite{Bus}, in that it characterizes
variationally the true ensemble-averages. Instead, Busse derived variational bounds for turbulent transport
properties and he has proposed that the extremalizing vector fields achieving these bounds will be ``similar''
to the ensemble-averaged turbulent fields. However, in general, they will be distinct. Detailed comparison of the
virtues and failings of these different principles must await future work.
\newpage
\section{Stochastic LES Models and Their Applications}
\noindent {\em (4.1) A Subgrid Random Coupling Model and Comparison With Other Models}
A number of stochastic LES models have already been proposed and implemented for turbulence simulation
\cite{Lei,Chas,MT,Car}. In principle, they might be regarded as simplified approximations of the {\em exact}
SLE derived here. However, there are certain differences.
The model of Chasnov \cite{Chas} was derived by applying the sharp Fourier cutoff filter to the momentum equation
of the ``generalized Langevin equation'' model \footnote{ Properly speaking, these equations differ in type from
what we have called ``generalized Langevin equations'' because they have the property that noise and damping terms
are determined from the past statistics of an infinite ensemble of realizations of the dynamics. A better term to
describe this type of stochastic equation would be ``self-consistent Langevin equation''.} for the EDQNM closure
equations. As we have discussed elsewhere \cite{I-LES}, the use of the sharp Fourier filter is unwise and introduces
non-universal features in the ``subgrid stresses'' which are due to large-scale convection. A more basic
difficulty with constructing LES models in this way was noted a long time ago by Kraichnan and Herring \cite{HKr}
(p.162). The self-consistent Langevin models were introduced to establish realizability for the
closures and, while they are expected to approximate well the low-order {\em statistics} of the turbulence dynamics
(e.g. mean energy spectra and transfer), they cannot be expected to provide faithful approximations for individual
realizations. In particular, Herring and Kraichnan noted that the models ``scramble the dynamics in a way
that makes it implausible that the model could reproduce the build-up of complicated correlations among
large numbers of wavevector modes which occurs in solutions of Navier-Stokes equations.'' This is a particularly
serious concern if the LES model is to be used to study the coherent structures which evolve spontaneously
in the turbulent flow, since production of these structures will be inhibited and their lifetimes reduced.
Another type of model was introduced by Leith \cite{Lei} and applied in a somewhat modified form to
simulation of a turbulent boundary layer by Mason and Thompson \cite{MT}. Those authors modelled the random
accelerations given in our Eq.(\ref{61}) by an application of the same energy-balance ideas used in
deriving the Smagorinsky model and by Kolmogorov-style dimensional reasoning. More precisely, a heuristic dimensional
argument was given that the rate of energy backscatter from turbulent eddies of size $\ell_e$ to the resolved
scale at the filter length $\ell_f$ is of the form $C_B\left({{\ell_e}\over{\ell_f}}\right)^5\ven.$ Here $\ven$ is the
local value of energy dissipation, which was estimated by assuming local energy balance of subgrid flux (from
a Smagorinsky stress model) with dissipation and backscatter:
\be \nu \overline{S}^2=\ven+C_B\left({{\ell_e}\over{\ell_f}}\right)^5\ven. \lb{X1} \ee
The backscatter was then implemented in the scheme by generating a field of independent random vectors
$\bphi=(\phi_x,\phi_y,\phi_z)$ of zero mean at each point of the spacetime grid. This random vector field was used to
generate a ``vector potential'' whose curl gave the random acceleration. The random vector $\bphi$ was filtered on the
scale $\ell_f$ to obtain the components of the potential and rescaled appropriately to give the estimated contribution to
backscatter rate in root mean square sense. Compared with the exact Eq.(\ref{61}), this model acceleration field has one
clearly unnatural feature, that it is given as the curl of a vector $\bA=\grad\btimes\bphi$ rather than the divergence
of a symmetric stress tensor $\bA=-\grad\bdot\btaus$ plus the associated contribution to pressure. Mason and Thompson
argued that only the divergence-free part of the random acceleration needed to be considered. However, the random
backscatter contribution to pressure could make an important contribution, e.g. to isotropization of the small scales,
and their procedure makes the pressure term unrealistically a deterministic quantity (i.e. a function just of resolved
fields). Another difference compared to the exact expression in Eq.(\ref{61}) is that the model acceleration in
this scheme is completely uncorrelated in space-time, whereas the exact expression allows long-range correlations
in space and memory of past history in time. Note that a heuristic ``multifractal model'' suggests the
presence of such correlation effects in energy dissipation \cite{MC}.
A model which contains some of these effects can be constructed by applying Kraichnan's ``random coupling''
method \cite{Kr61} in conjunction with the exact stochastic filtering scheme. \footnote{However, these
models will have the same deficiency as \cite{Chas} in omitting coherent subgrid structures.} We shall not
give full details here but just outline some of the ideas in the simple context of coupled random oscillators. The
dynamical variables are now complex numbers $\zl,\zs$ governed by
\be \dot{\zl}(t)+i\Bl\zl(t)+ic\zs(t)=\fL(t), \lb{X2} \ee
and
\be \dot{\zs}(t)+ic^*\zl(t)+i\Bs\zs(t)=\fS(t), \lb{X3} \ee
Here $\Bs,\Bl$ are real-valued random frequencies with independent distributions of zero means and variances
$\BL=\langle(\Bl)^2\rangle,\BS= \langle(\Bs)^2\rangle,$ also $c$ is a random complex coupling coefficient with mean
$\langle c\rangle$ and variance $C=\langle|c|^2\rangle-|\langle c\rangle|^2,$ and $\fL,\fS$ are external random forces
with means $\langle\fL\rangle,\langle\fS\rangle$ and covariances
\be \langle(\fL(t)-\langle\fL(t)\rangle)(\fL(t')-\langle\fL(t')\rangle)^*\rangle=\FL(t,t')\,\,\,\,\,\,\,\,\,\,
\langle(\fS(t)-\langle\fS(t)\rangle)(\fS(t')-\langle\fS(t')\rangle)^*\rangle=\FS(t,t'). \lb{X4} \ee
This system can be considered a caricature of a passive scalar, in which the $z$ variables mimic
the scalar concentration, $b,c$ mimic the convecting velocity with prescribed statistics, and $f$'s represent scalar sources.
A pair of oscillators are considered to represent the division into ``large-scale'' modes $\zl$ and ``small-scale''
modes $\zs.$ Because the dynamics is linear, it is trivial to implement the stochastic filtering with
\be \Zs[t;\zl,\fS]=\int_{t_0}^t ds\,\,e^{-i\Bs(t-s)}[-ic^*\zl(s)+\fS(s)]. \lb{X5} \ee
Substituted back into the ``LE'', Eq.(\ref{X2}), this yields the exact stochastic equation for $\zl.$ Note that
in this example, which is clearly an ``integrable'' case, the exact SLE dynamics becomes deterministic as
$\BS, C,$ and $\FS$ all vanish. In that case, the filtering just produces the new term
\be \overline{K}[t,\zl]= -|\langle c\rangle|^2\int_{t_0}^t ds\,\,\zl(s)\cdot\zl(t) \lb{X6} \ee
in the effective dynamics.
To make a closure for the general case, we may apply Kraichnan's random coupling scheme. We introduce as in \cite{Kr61}
$N$ independent copies $\zl_{[n]},\zs_{[n]}$ of the above system,$n=0,1,...,N-1,$ and rewrite them coupled together
in ``collective variables'' introduced by a discrete ${\bf Z}(N)$-Fourier transform:
\be \zl_\alpha={{1}\over{\sqrt{N}}}\sum_{n=0}^{N-1}e^{2\pi i\alpha n/N}\zl_{[n]}\,\,\,\,\,\,\,\,\,\,
\zs_\alpha={{1}\over{\sqrt{N}}}\sum_{n=0}^{N-1}e^{2\pi i\alpha n/N}\zs_{[n]}. \lb{X7} \ee
Into the coupled equations are then introduced random interaction phases $\theta_{\alpha,\beta,\gamma},
\theta'_{\alpha,\beta,\gamma}$ {\em but only into the terms involving the SS modes (alone or coupled to the LS modes)}.
The dynamics of the LS modes are left unaltered. For obvious reasons, we refer to this new system as the
``subgrid random coupling model.'' The result upon transferring back to the ``individual variables''
$\zl_{[n]},\zs_{[n]}$ and taking $N\rightarrow +\infty$ is that the individuals become again uncorrelated
in the limit and that each has LS modes now governed by a stochastic effective dynamics of the form
\be \partial_t\zl(t)+i\Bl\zl(t)+\int_{t_0}^t ds\,\,\Sigmal(t,s)\zl(s)=\fL(t)+\phil(t), \lb{X8} \ee
in which $\Sigmal$ is a kernel representing additional damping and $\phil$ is a new random force with mean
$\langle\phil(t)\rangle$ and covariance
\be \langle(\phil(t)-\langle \phil(t)\rangle)(\phil(t')-\langle\phil(t')\rangle)^*\rangle=\Phil(t,t'). \lb{X9} \ee
These additional terms represent the effects of the eliminated SS modes and can be characterized by the statistics
of those modes. Precisely,
\be \Sigmal(t,s)= C\cdot \Gs(t,s) \lb{X10} \ee
and
\be \Phil(t,s)= C\cdot\left(\Zs(t,s)+\langle\zs(t)\rangle\langle\zs(s)\rangle^*\right), \lb{X11} \ee
and
\be \langle\phil(t)\rangle= -i\langle c\rangle\langle\zs(t)\rangle. \lb{X11a} \ee
In addition to the mean $\langle\zs(t)\rangle$ these involve
\be \Gs(t,s)=\langle {{\delta\zs(t)}\over{\delta \fS(s)}}\rangle, \lb{X12} \ee
the average response function of the SS modes, and
\be \Zs(t,s)=\langle (\zs(t)-\langle\zs(t)\rangle)(\zs(s)-\langle\zs(s)\rangle)^*\rangle, \lb{X12a} \ee
the SS correlation function. Note that the new terms in the SLE dynamics Eq.(\ref{X8}) are precisely those
appearing at the 1-loop level in the line-reverted expansion (corresponding to the first and second diagrams
in Section 2.3). All other contributions from the SS modes are eliminated by averaging over the random
phases in the coupling interactions in the large-N limit.
The SS statistics are then determined by a set of self-consistent equations, which we call ``subgrid DIA equations.''
They involve the SS statistics as well as the corresponding statistical quantities $\langle\zl(t)\rangle, \Gl(t,s),
\Zl(t,s)$ of the LS modes. Precisely, the equations are
\be \partial_t\langle\zs(t)\rangle+i\langle c\rangle^*\langle\zl(t)\rangle+
\int_{t_0}^t ds\,\,\Sigmas(t,s)\langle\zs(s)\rangle=\langle \fS(t)\rangle, \lb{X13} \ee
and
\be \partial_t \Gs(t,t')+\int_{t'}^t ds\,\,\Sigmas(t,s)\Gs(s,t')=\delta(t-t'), \lb{X14} \ee
and
\be \partial_t \Zs(t,t')+\int_{t_0}^t ds\,\,\Sigmas(t,s)\Zs(s,t')=
\int_{t_0}^{t'}ds\,\,\left[\FS(t,s)+\Phis(t,s)\right]\Gs(t',s)^*. \lb{X15} \ee
In these equations the new quantities are
\be \Sigmas(t,s)= \BS\cdot\Gs(t,s)+C\cdot\Gl(t,s), \lb{X16} \ee
and
\be \Phis(t,s)= \BS\cdot\left[\Zs(t,s)+\langle\zs(t)\rangle\langle\zs(s)\rangle^*\right]+
C\cdot\left[\Zl(t,s)+\langle\zl(t)\rangle\langle\zl(s)\rangle^*\right]. \lb{X17} \ee
These subgrid DIA equations may be obtained from a set of exact Schwinger-Dyson integral equations
for the SS dynamics ``conditioned'' on a fixed past history of the LS modes. Because of averaging over random
phases and the limit $N\rightarrow +\infty,$ they depend only upon the low-order statistics
of the LS modes and only the 1-loop terms in the ``self-energies'' $\Sigmas$ and $\Phis$ survive. It
may be verified that the combined system of SLE dynamics and subgrid DIA equations satisfy necessary
consistency properties, such as conservation in the mean
\be \partial_t\langle|\zl(t)|^2\rangle+\partial_t\left[\Zs(t,t)+|\langle\zs(t)\rangle|^2\right]=0 \lb{X18} \ee
for $\langle f\rangle,F=0,$ which follow from the existence of a model representation.
This approximation to the exact SLE dynamics clearly does retain some memory effects of the SS modes depending upon
the past history of the mean statistics of the LS modes. A similar ``subgrid random coupling model'' can be
constructed for the nonlinear Navier-Stokes equation, but it is far more complicated and will not be written
here in detail. In this case there is a serious defect of the DIA approximation to the exact SLE equations, which is the
failure of {\em random Galilei invariance} (RGI). This problem is exactly the same as that in ordinary DIA, discussed at
length by Kraichnan in \cite{Kr87,Kr85,Kr77}. The failure of RGI gives rise to qualitatively bad approximations to the
damping, noise, etc. calculated from the DIA equations, or, more generally, from any truncation of the perturbation series
in the ``line-reverted'' form. These defects can be cured by using Lagrangian reformulations of the line-reverted
expansions \cite{Kr77} but only at the price of losing the model representation. Probably the exact ``subgrid DIA''
equations for Navier-Stokes are too complicated anyway to be of practical use in simulations and further
approximations must be made along the lines leading to EDQNM-type closures.
There is one feature of the ``subgrid DIA'' equations for the coupled random oscillator system which should be
mentioned. In the limit as $\BS,C,\FS\rightarrow 0$ these equations degenerate, so that $\Gs(t,t')\rightarrow\theta(t-t'),
\Zs(t,t')\rightarrow 0,$ and the mean obeys the simple equation
\be \partial_t\langle\zs(t)\rangle+i\langle c\rangle^*\langle\zl(t)\rangle=\langle \fS(t)\rangle. \lb{X19} \ee
These results correctly account for the fact that the exact SS dynamics becomes deterministic in that limit.
Likewise, when $C\rightarrow 0,$ the covariance $\Phil\rightarrow 0$ and the force $\phil(t)$ generated by the SS modes
becomes deterministic. Therefore, no spurious stochasticity is generated in the LS dynamics by the DIA approximation.
However, this is possibly an artefact of the linearity of the problem. In the case of nonlinear dynamics, there
is no obvious distinction between the subgrid DIA equations for an ``integrable'' case such as Burgers and a
``chaotic'' case such as Kuramoto-Shivashinsky in the limit $\ven_0\rightarrow 0,$ where external noise is
removed. However, it may be that a qualitative distinction emerges from their detailed solution. The issue may
be formulated as follows: with small randomness $\sim \ven_0$ in the {\em initial data} of the SS modes, the KS
equation should generate an $O(1)$ effective noise for the LS modes in a time proportional to a local ``mixing time''
of the dynamics, with a coefficient which is only weakly dependent on $\ven_0$ (e.g. logarithmically). However, the
Burgers equation will presumably generate an $O(1)$ effective noise for the LS modes from randomness $\sim\ven_0$ in the
initial data of the SS modes in a time which goes rapidly to infinity in the limit as $\ven_0\rightarrow 0.$ It would be
interesting to know whether the DIA approximation to the exact SLE equations is sophisticated enough to make this
distinction between ``chaotic'' and ``integrable'' dynamics.
\noindent {\em (4.2) Atmospheric Predictability and Complex Flows}
Let us just conclude this work by a brief mention of some situations of physical interest where ``turbulent eddy
noise'' may play a significant role. It was already argued by Mason and Thompson \cite{MT} that stochastic backscatter
is necessary to correct errors in the {\em mean} velocity profiles of the turbulent boundary layer produced by
a deterministic LES model of Smagorinsky type. Similarly, the importance of turbulent fluctuations generally in
determining the means lead us to believe that the least-action principle may have some use in calculating those
mean statistics. Mason and Thompson also found in their study significant qualitative differences between individual
realizations in LES with and without stochastic backscatter. For example, with backscatter the velocity fields
were much ``rougher'' or ``irregular.'' Thus, eddy noise may also have important effects on the character of large-scale
structures that evolve in the flow.
Another obvious area of interest is the role of eddy noise in the predictability problem for weather forecasting
and climate change. An early study of this issue within the TFM closure was made by Leith and Kraichnan \cite{LKr}.
It is clear that deterministic LES models of atmospheric flow make a spurious prediction that LS velocities
are completely predictable---in principle---given exact information on just the initial LS fields. A stochastic
LES model with eddy noise from the unresolved SS fields corrects this obvious defect. It is of interest how such
effects of turbulence noise interact with the ``deterministic chaos'' view on the predictability problem \cite{Lor2}.
Since even deterministic LES models show generally ``sensitive dependence'' to initial perturbations, adding
eddy noise may not produce any qualitative differences. However, a quantitative difference, even by a factor of
two or so, would have significant impact on forecasting ability. These issues are currently investigated in the
meteorological literature: see \cite{PNNWZ} for a review. Improved physical understanding of the characteristics
of turbulence noise would be valuable for this inquiry.
\vspace{.3in}
\noindent {\bf Appendix I: The Kraichnan-Wyld Perturbation Theory}
The perturbation expansion for the Navier-Stokes equation Eq.(\ref{11}) can be developed by writing
an exact integral solution, as:
\be \bv(\br t)=\bv^{(0)}(\br t)+\int_{-\infty}^tdt'\int d^d\br'\,G^{(0)}(\br t,\br't')\bP(\grad_{\br'})
(\bv(\br't')\cdot\grad_{\br'})\bv(\br't'), \lb{22} \ee
where $G^{(0)}$ is the bare response function
\be G^{(0)}(\br t,\br't')=(\partial_t-\nu\bigtriangleup)^{-1}(\br t,\br't'), \lb{23} \ee
and $\bv^{(0)}$ is the solution of the linearized problem
\be \bv^{(0)}(\br t)=\int_{-\infty}^tdt'\int d^d\br'\,G^{(0)}(\br t,\br't')\bF(\br't'). \lb{24} \ee
An exact integral solution can also be obtained for the full response tensor,
\be \widehat{G}_{ij}(\br t,\br't')\equiv {{\delta v_i(\br t)}\over{\delta f_j(\br't')}}, \lb{25} \ee
by taking the functional derivative of both sides of Eq.(\ref{22}) with respect to $\bF(r't')$, yielding
\begin{eqnarray}
\,& & \widehat{G}_{ij}(\br t,\br't')=G_{ij}^{(0)}(\br t,\br't') \cr
\,& &\,\,\,\,\,\,\, +\int_{-\infty}^tdt''\int d^d\br''\,G^{(0)}(\br t,\br''t'')P_{il}(\grad_{\br''}) \cr
\,& &\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\times\left[\widehat{G}_{kj}(\br''t'',\br't')(\grad_{\br''})_k v_l(\br''t'')
+(\bv(\br't')\cdot\grad_{\br''})\widehat{G}_{lj}(\br''t'',\br't')\right]. \lb{26}
\end{eqnarray}
Now, iterating the two expressions, Eq.(\ref{22}) and Eq.(\ref{26}), one obtains expansions of
the exact $\bv$ and $\widehat{{\bf G}}$ as power series in $\bv^{(0)}$ and ${\bf G}^{(0)}.$ One may
substitute the series into the ensemble averages for the covariance
\be U_{ij}(\br t,\br't')=\langle v_i(\br t)v_j(\br't')\rangle-\langle v_i(\br t)\rangle \langle v_j(\br't')\rangle, \lb{27} \ee
average response tensor
\be G_{ij}(\br t,\br't')=\langle \widehat{G}_{ij}(\br t,\br't')\rangle, \lb{28} \ee
and mean field $\bar{\bv}(\br t)$ (Eq.(\ref{13a}) above). In the case of Gaussian force $\bF$, the averages over
products of the $\bv^{(0)}$'s may be resolved into products of $\bar{\bv}^{(0)}$'s
and bare correlation functions
\be U^{(0)}_{ij}(\br t,\br't')=\langle v_i^{(0)}(\br t)v_j^{(0)}(\br't')\rangle-\langle v_i^{(0)}(\br t)\rangle \langle v_j^{(0)}(\br't')\rangle, \lb{29} \ee
by the use of Wick's theorem. In this way, series expansions for
$\bar{\bv},{\bf U},$ and ${\bf G}$ may be developed in terms of the corresponding bare
quantities.
A graphical interpretation of this iteration procedure was introduced in \cite{Wyl,Kr61}, in which the
basic elements are the ``propagators,'' which are the bare response and correlation function, along with the mean field
and the bare interaction vertex
\begin{eqnarray}
\gamma_{3;ijk}(\br t,\br't',\br''t'') & = & -{{1}\over{2}}\left(P_{ij}(\grad_\br)(\grad_\br)_k+P_{ik}(\grad_\br)(\grad_\br)_j\right) \cr
\,& &\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\times \delta^d(\br-\br')\delta(t-t')\delta^d(\br-\br'')\delta(t-t''). \lb{30}
\end{eqnarray}
Diagrammatically one represents the bare field as:
\be v^{(0)}_{i}\equiv\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, \lb{31}\ee
the average bare response function as
\be G^{(0)}_{ij}\equiv\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, \lb{32} \ee
and the bare vertex as
\be \gamma_{3;ijk}\equiv\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, \lb{33}.\ee
In terms of these graphical expressions, one generates the series for $\bv$ by first writing down all tree graphs whose ``limbs'' are given by ${\bf G}^{(0)}$'s
and whose ``top branches'' are given by $\bv^{(0)}$'s. To form the series for $\bar{\bv}$ one
eliminates all the $\bv^{(0)}$'s in the single trees for $\bv$, either by replacing them singly by
$\bar{\bv}^{(0)}$'s or joining them in pairs to yield ${\bf U}^{(0)}$'s, which are represented as
\be \bar{v}^{(0)}_{i}\equiv\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\lb{34}\ee
and
\be U^{(0)}_{ij}\equiv\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, \lb{35}\ee
The averaging elimination procedure performed on all pairs of trees (joined by at least one ${\bf U}^{(0)}$)
gives the series for ${\bf U}.$ Finally, the series for ${\bf G}$ is obtained if one
first replaces one ``top branch'' in single trees by a ${\bf G}^{(0)}$ and then performs the
averaging elimination of the $\bv^{(0)}$'s. A convenient compendium of some of the topological
properties of graphs in the resulting diagrammatic perturbation series is contained in \cite{ZL}.
The same perturbation expansion can be obtained as a Feynman diagram expansion for the MSR field-theory.
It is convenient to introduce a doublet field $\Phi_\sigma,\,\sigma=\pm$, following \cite{MSR}, as:
\be \Phi(\br t)=\left( \begin{array}{c}
\Phi_+(\br t) \\ \Phi_-(\br t)
\end{array} \right)
= \left( \begin{array}{c}
\bv(\br t) \\ i\hat{\bv}(\br t)
\end{array} \right) \lb{36}. \ee
Then the path-integral expression for the generating functional $Z[\eta]=Z[-i\eta_+,\eta_-],$ $\eta=(\eta_+,\eta_-)^\top,$
is just
\begin{eqnarray}
Z[\eta]&=&\int{\cal D}\Phi \exp\left[-{{1}\over{2}}\Gamma^{(0)}_2(12) \Phi(1)\Phi(2) \right. \cr
\,& &\,\,\,\,\,\,\,\,\,\,\,\,\left.+\gamma_1(1)\Phi(1)+{{1}\over{3!}}\gamma_3(123)\Phi(1)\Phi(2)\Phi(3)+\eta(1)\Phi(1)\right], \lb{37}
\end{eqnarray}
where we have employed a compact notation $(\sigma_1i_1\br_1t_1)\equiv 1$ along with a summation convention on repeated indices.
Here, $\gamma_1=(0,\bar{\bF})^\top.$ $\Gamma^{(0)}_2$ is the ``wave operator'' given by the $2\times 2$ matrix in doublet space
\be \Gamma^{(0)}_{2;ij}(\br t,\br't')=\left(\begin{array}{ll}
0 & -(\partial_t+\nu_0\bigtriangleup_\br)\delta_{ij}\delta^d(\br-\br')\delta(t-t') \cr
(\partial_t-\nu_0\bigtriangleup_\br)\delta_{ij}\delta^d(\br-\br')\delta(t-t') & -F_{ij}(\br t,\br't')
\end{array} \right), \lb{38} \ee
and $\gamma_3(123)$ is a completely symmetric vertex which is equal to the $\gamma_3$ already defined in Eq.(\ref{30})
when the doublet indices are $(-++),(+-+),$ or $(++-)$ and zero otherwise. Of course, it must be the case---and it is easy
to check directly---that $\Gamma^{(0)}_2=(G^{(0)}_2)^{-1},$ where $G^{(0)}_2$ is the ``bare'' matrix correlation function
$G^{(0)}_{2;\sigma\sigma'}=\langle \Phi_\sigma^{(0)}\Phi_{\sigma'}^{(0)}\rangle-\langle \Phi_\sigma^{(0)}\rangle\langle \Phi_{\sigma'}^{(0)}\rangle$
which has three non-vanishing entries
\be G^{(0)}_2(i\br t,j\br't')=\left(\begin{array}{ll}
U^{(0)}_{ij}(\br t,\br't') & G^{(0)}_{ij}(\br t,\br't') \cr
\bar{G}^{(0)}_{ij}(\br t,\br't') & 0
\end{array} \right), \lb{39} \ee
with $\bar{G}^{(0)}_{ij}(\br t,\br't')=G^{(0)}_{ji}(\br't',\br t)$ a ``(bare) anti-response function'' and with
$\bG^{(0)}$ and $\bU^{(0)}$ already defined in Eqs.(\ref{23}) and (\ref{29}).
\noindent {\bf Appendix II: Proof of the Liouville Property}
It is helpful to begin by giving a proof of Lee's result \cite{Lee}, that the Liouville theorem is valid for
the cutoff Euler dynamics associated to the ``dynamical field''
\be \bK^E(x)=-\bP(\grad_\bx)\grad_\bx\cdot\btau^E(x), \lb{71a} \ee
with
\be \btau^E(x)=\bv(x)\bv(x). \lb{72a} \ee
All of the velocity fields are now defined with a highwavenumber cutoff
\be \bv(\bx,t)=\sum_{\bk\in {\cal K}}\,\widehat{\bv}(\bk,t)e^{i\bk\bdot\bx}, \lb{71aa} \ee
associated to some bounded domain of wavevectors ${\cal K},$ e.g. the sphere ${\cal K}=\{\bk: |\bk|