\documentclass[reqno]{amsart}
\usepackage{hyperref}


\AtBeginDocument{{\noindent\small
{\em Electronic Journal of Differential Equations},
Vol. 2005(2005), No. 136, pp. 1--14.\newline
ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu
\newline ftp ejde.math.txstate.edu  (login: ftp)}
\thanks{\copyright 2005 Texas State University - San Marcos.}
\vspace{9mm}}

\begin{document}
\title[\hfilneg EJDE-2005/136\hfil Residual models for nonlinear PDEs]
{Residual models for nonlinear partial differential equations}
\author[G. Pantelis\hfil EJDE-2005/136\hfilneg]
{Garry Pantelis}

\address{Garry Pantelis \hfill\break
Department of Computer Science, Intercollege,
P.O. Box 42572, 6500 Larnaca, Cyprus}
\email{pantelis.g@adm.lar.intercol.edu}


\date{}
\thanks{Submitted November 14, 2005. Published November 30, 2005.}
\subjclass[2000]{35G20, 35G25, 53D10, 93A30}
\keywords{Partial differential equations; nonlinear evolution;
\hfill\break\indent  contact manifolds; mathematical modelling}


\begin{abstract}
 Residual terms that appear in nonlinear PDEs that are constructed
 to generate filtered representations of the variables of the
 fully resolved system are examined by way of a consistency condition.
 It is shown that certain commonly used empirical gradient models
 for the residuals fail the test of consistency and therefore cannot
 be validated as approximations in any reliable sense.
 An alternate method is presented for computing the residuals.
 These residual models are independent of free or artificial parameters
 and there direct link with the functional form of the system of PDEs
 which describe the fully resolved system are established.
\end{abstract}

\maketitle
\numberwithin{equation}{section}
\newtheorem{theorem}{Theorem}[section]

\section{Introduction}

Nonlinear systems of PDEs can generate solutions that exhibit
multiscale fluctuations. When numerical methods are used as the
solution method the discretization of the solution domain will
result in a loss of resolution of the finer scale fluctuations.
The macroscopic equations that are constructed to generate some
kind of filtered representation of the fully resolved variables
contain residual terms that attempt to model mechanisms that are
manifestations of some microscopic process. Here we shall examine
a number of residual models for nonlinear PDE systems in common
use today, which we shall refer to as empirical gradient models.
It will be shown on the basis of a consistency condition alone
that these models cannot be regarded as residual approximations in
any reasonable sense.

Some of the results presented here appear in \cite{Pan04} but have
been recast in a setting of contact manifolds where the base
manifold is space/time/scale. This is done because in such a
setting analysis of nonlinear PDEs appears more natural and a new
level of insight into the derivations is gained from a geometric
perspective. We start with a few preliminary details of contact
manifolds that will be adequate for our purposes. For a more
rigorous development of integral manifolds of distributions on
contact manifolds see for instance \cite{War01}.

Let $x = ( x^1 , \dots , x^n )$ represent the cartesian
coordinates  of the $n$-dimensional space domain $\mathbb{R}^n$
($n=1,2,3$), $t \in (0,t_0)$ the time, where $t_0 >0$, and $\eta
\in (0,\eta_0)$ the scale parameter, where $0 < \eta_0 \ll 1$. Set
$M \subset \mathbb{R}^m$, $m=n+2$, such that $M= \mathbb{R}^n
\times (0,t_0) \times (0,\eta_0)$ with coordinates on $M$
represented by $(x,t,\eta)=( x^1 , \dots , x^n ,t,\eta)$, where
$x^k$ ($k=1, \dots ,n$) are the spatial coordinates, $x^{n+1} = t$
is the time coordinate and $x^{n+2}=\eta$ is the scale coordinate.
It is important to keep in mind that we use the symbol $x$ to
refer only to the spatial component of $(x^i) = (x,t,\eta)=( x^1 ,
\dots , x^n ,t,\eta)$.

Let $K \subseteq \mathbb{R}^{m+N+Nm+Nm^2}$ denote the second order
contact  manifold whose coordinates are given by
$(x^i,u^\alpha,u_i^\alpha,u_{ij}^\alpha)$, $1 \leq \alpha \leq N$,
$1 \leq i,j \leq m=n+2$. We associate $u^\alpha$ as placeholders
for the solutions of PDEs and $u_i^\alpha$ and $u_{ij}^\alpha$ as
the placeholders for the first and second partial derivatives of
$u^\alpha$. As components of the coordinates of $K$ they are
independent variables.

The natural basis for the tangent space $T(M)$ of the base
manifold $M$ is $\{ \partial_i \}_{1 \leq i \leq m}$ and the natural
basis of the tangent space $T(K)$ of $K$ is $\{ \partial_i ,
\dot{\partial}_\alpha , \dot{\partial}_\alpha^i , \dot{\partial}_\alpha^{ij} \}_{1
\leq i,j \leq m; 1 \leq \alpha \leq N}$, where we use the notation
\begin{equation}
\partial_i = \frac{\partial}{\partial x^i} , \quad \dot{\partial}_\alpha = \frac{\partial~}{\partial
u^\alpha} , \quad \dot{\partial}_\alpha^i = \frac{\partial~}{\partial u_i^\alpha} ,
\quad \dot{\partial}_\alpha^{ij} = \frac{\partial~}{\partial u_{ij}^\alpha}
\label{cm110}
\end{equation}
Let $\mathcal{D}$ be an $m$-dimensional distribution on
$K$ (i.e. an $m$-dimensional subset of  $T(K)$). Let the
distribution admit as a basis  $\{ \mathcal{V}_i \}_{1 \leq i \leq m}$
where the vector fields $\mathcal{V}_i$ take the form
\begin{equation}
\mathcal{V}_i = \partial_i + u_i^\alpha \dot{\partial}_\alpha + u_{ij}^\alpha
\dot{\partial}_\alpha^j +
A_{ijk}^\alpha \dot{\partial}_\alpha^{jk} \, , \quad 1 \leq i \leq m
\label{cm150}
\end{equation}
for some $A_{ijk}^\alpha \in \mathcal{F} (K)$ and
where the convention of summation over repeated upper and lower
indices is used throughout. The distribution $\mathcal{D}$ is
involutive, or completely integrable, if $[\mathcal{D},\mathcal{D}] \in
\mathcal{D}$, where $[ \cdot , \cdot ]$ denotes the Lie brackets. With
respect to the basis $\{ \mathcal{V}_i \}_{1 \leq i \leq m}$ the
involutive condition for the distribution $\mathcal{D}$ is satisfied
only if
\begin{equation}
u_{ij}^\alpha = u_{ji}^\alpha \, , \quad
A_{ijk}^\alpha = A_{jik}^\alpha = A_{jki}^\alpha\, , \quad \mathcal{V}_i
\langle A_{jkl}^\alpha \rangle = \mathcal{V}_j \langle A_{ikl}^\alpha
\rangle \label{cm170}
\end{equation}
which leads to $[ \mathcal{D},\mathcal{D} ] =0$. By
the Frobenius theorem we know that if the distribution $\mathcal{D}$ is
involutive then through each point $p \in K$ there passes an
integral manifold of the same dimension as the distribution
$\mathcal{D}$.

Let $\phi : M \to K$ be a smooth ($C^\infty$) map from $M$ to $K$.
The smooth map $\phi : M \to K$ induces the pullback map
$\phi^* : \mathcal{F}(K) \to \mathcal{F}(M)$ which maps smooth functions on
$K$ to smooth functions on $M$.
We will write $\mathcal{F}(A)$ to mean the set of smooth functions on
the manifold $A$ and wherever the map $\phi : M \to K$ has been
specified and there is no confusion as to which map we a dealing
with we will sometimes use the shorthand notation
$\,^*F = \phi^* F = F \circ \phi$ for $F \in \mathcal{F} (K)$.

The smooth map $\phi$ also induces a smooth differential map
$\phi_*$ such that at each point $p \in M$ the differential map
$\phi_*$ maps members of $T_p(M)$ (the tangent space of $M$ at $p
\in M$) to members of $T_{\phi(p)} (K)$ (the tangent space of $K$
at $\phi(p) \in K$). Suppose that $(M,\phi)$ is an integral
manifold of the involutive distribution $\mathcal{D}$ and that at each
point $p \in M$ the differential map $\phi_*$ maps each member of
the natural basis of $T_p(M)$ to each member of the basis of
$\mathcal{D}$ at $\phi(p) \in (M,\phi)$, i.e. $\phi_* ( \partial_i |_{p}) =
\mathcal{V}_i |_{\phi(p)}$. Consider the function $F \in \mathcal{F}(K)$ which
upon the action of the pullback of $\phi$ is mapped to $\,^*F \in
\mathcal{F} (M)$. By a straightforward calculation
\begin{equation}
\partial_i (^*F ) =
\,^*(\partial_i F) + \frac{\partial \,^*u^\alpha}{\partial x^i} \,^*(\partial_\alpha F) +
\frac{\partial \,^*u_j^\alpha}{\partial x^i} \,^*(\partial_\alpha^j F) + \frac{\partial
\,^*u_{jk}^\alpha}{\partial x^i} \,^*(\partial_\alpha^{jk} F) \label{cm120}
\end{equation}
Since we are using $u_i^\alpha$ and $u_{ij}^\alpha$ as place
holders for the first and second partial derivatives of
$\,^*u^\alpha$ we immediately write
\begin{equation}
\,^*u_i^\alpha = \frac{\partial
\,^*u^\alpha}{\partial x^i} ~, \quad \,^*u_{ij}^\alpha = \frac{\partial
\,^*u_j^\alpha}{\partial x^i} = \frac{\partial^2 \,^*u^\alpha}{\partial x^i \partial x^j}
\label{cm130}
\end{equation}
Because we have no place holders for third
partial derivatives and we are confined to $K$ we require that
under the action of the pullback of $\phi$,
\begin{equation}
\,^*A_{ijk}^\alpha
= \frac{\partial \,^*u_{jk}^\alpha}{\partial x^i} = \frac{\partial^3 \,^*u^\alpha}{\partial x^i
\partial x^j \partial x^k} \label{cm140}
\end{equation}
In shorthand notation
(\ref{cm120}) can be written \begin{equation} \partial_i (^*F ) = \,^*(\mathcal{V}_i F)
\label{cm160} \end{equation} We shall call a smooth map $\phi : M \to K$ a
\emph{regular map} if $(M,\phi)$ is an $m$-dimensional integral
manifold of the involutive distribution $\mathcal{D}$. The distribution
$\mathcal{D}$ will admit as a basis $\{ \mathcal{V}_i \}_{1 \leq i \leq m}$,
where the vector fields $\mathcal{V}_i$ are given in (\ref{cm150}).
Since $\mathcal{D}$ is involutive we can assume the symmetry relations
(\ref{cm170}). Under the action of the pullback of the map $\phi$
we have the identities (\ref{cm130})  and (\ref{cm140}) which
allow us also to apply (\ref{cm160}). We shall say that the map
$\phi :M \to K$ is a \emph{bounded regular map} if it is a regular
map and $u^\alpha,u_i^\alpha,u_{ij}^\alpha$ are bounded on the
image $\phi(M)$ contained in $K$.

Consider some $N$ functions $P^\alpha \in \mathcal{F}(K)$ which have the
representation given by $P^\alpha (x^i,u^\beta,u_i^\beta,u_{ij}^\beta)$.
Suppose that the map $\phi : M \to K$ is associated with an integral
 manifold of dimension $m$ of the involutive distribution $\mathcal{D}$ such
that $\phi$ annihilates $P^\alpha$, i.e. $\phi^* P^\alpha =0$.
Then by the pullback identities (\ref{cm130}) $\phi^* P^\alpha =0$
has the representation
\begin{equation}
P^\alpha \left ( x^i, \,^*u^\beta,
\frac{\partial \,^*u^\beta}{\partial x^i} , \frac{\partial^2 \,^*u^\beta}{\partial x^j \partial x^i}
\right ) =0 \label{cm180}
\end{equation}
which can be regarded as a system of
$N$ PDEs on the base manifold $M$ for the $N$ unknown dependent
variables $\,^*u^\alpha \in \mathcal{F}(M)$. On the other hand consider
the point set $\Xi \in K$ defined by $\Xi = \{ p \in K :
P^\alpha = 0 \}$. The image $\phi (M)$ will be contained in $\Xi$
and existence theorems for solutions can then be transferred to
the study of the properties of the constraints $P^\alpha =0$ on
the point set $\Xi$. Crucial to this endeavour is the
establishment of the symmetry conditions (\ref{cm170}) on the
image $\phi(M)$ contained in the point set $\Xi$.

The question of existence of solutions to PDEs is often attacked by way
of the so called fundamental ideal which contains the contact ideal as
a subideal (see for instance \cite{Ede90}, \cite{Ede01}).
It turns out that $\{ \mathcal{V}_i \}_{1 \leq i \leq m}$ forms a basis for
the $m$-dimensional module of Cartan annihilators of the closed contact
ideal.
In \cite{Ede01} yet another ideal is constructed called the horizontal
ideal which contains the contact ideal as a subideal.
The set of all solution maps of the horizontal ideal contains all the
solution maps of the closed contact ideal.
It is shown that the module of Cauchy characteristics of the horizontal
ideal admits as a basis $\{ \mathcal{V}_i \}_{1 \leq i \leq m}$ and that the
horizontal ideal is stable under the Lie transport with respect to the
vector fields $\mathcal{V}_i$.
It follows that the closure of the horizontal ideal
(integrability condition) is gauranteed by the symmetry relations
(\ref{cm170}).

It is worth noting that the analysis that follows can be recast on a 
setting of the first order contact manifold.
To do this one needs to define solution maps as annihilators of the so 
called balance ideal whose generators consist of the generators of the 
horizontal ideal and the balance $m$-forms associated with the PDEs under 
consideration \cite{Ede01}.
Our objectives here are not concerned with establishing existence
theorems for PDEs but rather to derive certain results for filtered
variables of nonlinear PDEs using space/time/scale as a base manifold.
For this purpose the concept of integral manifolds in relation to
involutive distributions will be adequate.
In the most part we will examine the properties of solution maps that
we know exist and in particular solution maps associated with the heat
equation.
We do this because the heat equation recast on space/scale rather than
space/time can be directly linked to spatial filters of the fully resolved
variables.

\section{Macroscopic Equations}

We set $N=2 N'$ and decompose each of $u^\alpha$,
$u_i^\alpha$ and $u_{ij}^\alpha$ into two parts
\begin{equation}
\begin{alignedat}{3}
 u^A &= v^A ,  &\quad u^{N' + A} &= r^A , &\quad 1 &\leq A \leq N'  \\
u_i^A &= v_i^A , & u_i^{N' + A} &= r_i^A ,&    1 &\leq A \leq N',\;
  1 \leq i \leq n+2  \\
u_{ij}^A &= v_{ij}^A , & u_{ij}^{N' + A} &= r_{ij}^A , &
 1 &\leq A \leq N' , \;  1 \leq i,j \leq n+2
\end{alignedat} \label{me100}
\end{equation}
where the $v^A$ will be associated with place holders for the filtered
 variables and the $r^A$ will be associated with place holders for the
residuals which will be defined below.
In expanded form the coordinates of $K$ become
\begin{equation}
(x^i,u_i^\alpha,u_{ij}^\alpha)=(x^i,v^A,r^A,v_i^A,r_i^A,v_{ij}^A,r_{ij}^A)
\label{me101}
\end{equation}
Throughout, unless otherwise stated, we use the Latin indices for the
range $1 \leq i,j,k \leq m=n+2$, the Greek indices for the range
$1 \leq \alpha ,\beta, \gamma \leq N=2N'$ and the capital Latin indices
for the range $1 \leq A,B,C \leq N'$.
For easier identification of certain important operators we use the
indices $t$ instead of $n+1$ and $\eta$ instead of $n+2$, e.g.
we write $\mathcal{V}_t$ instead of $\mathcal{V}_{n+1}$ and $\mathcal{V}_\eta$
instead of $\mathcal{V}_{n+2}$.

Let $P^A$, ($1 \leq A \leq N'$), have the representation $P^A =
P^A (x^i , v^B , v_i^B , v_{ij}^B )$ indicating that $P^A$ will be
independent of the residuals and their partial derivatives.
Suppose that under the pullback of some regular map $\phi : M \to
K$, there exists a limiting solution $\phi^*  v^A |_{\eta=0^+} =
\widetilde{v}^A$, $\widetilde{v}^A \in \mathcal{F} (\mathbb{R}^n \times (0,t_0) )$,
such that 
\begin{equation} 
P^A \left ( x^i , \widetilde{v}^B , \frac{\partial \widetilde{v}^B}{\partial x^i}
, \frac{\partial^2 \widetilde{v}^B}{\partial x^i \partial x^j} \right ) = 0 
\label{me105} \end{equation}
We define (\ref{me105}) as the system of PDEs which define the
fully resolved system and can be represented as the limiting
system $\phi^* P^A |_{\eta=0^+} = 0$. Note that given our agreed
range of the Latin indices that (\ref{me105}) will contain terms
involving $\partial/\partial \eta$. In most cases of interest such terms will
be redundant. Since these terms will not effect the calculations
in the derivations we leave the general representation of $P^A$ as
it is to avoid introducing more notation.

Let $\triangle$ denote the spatial Laplacian operator acting on
members of $\mathcal{F} (M)$
\begin{equation}
\triangle = \sum_{b=1}^n \frac{\partial^2}{\partial x^b \partial x^b}
\label{me140}
\end{equation}
and the associated operator acting on members of $\mathcal{F} (K)$
\begin{equation}
\mathcal{L} = \sum_{b=1}^n \mathcal{V}_b \mathcal{V}_b
\label{me150}
\end{equation}
We shall make repeated use of the vector field
\begin{equation}
\mathcal{W} = \partial_\eta + (\mathcal{L} u^\alpha ) \dot{\partial}_\alpha
+ (\mathcal{L} u_i^\alpha) \dot{\partial}_\alpha^i
+ (\mathcal{L} u_{ij}^\alpha) \dot{\partial}_\alpha^{ij}
\label{me160}
\end{equation}
where the convention of summation over repeated upper and lower
indices of $\alpha$ and $i,j$ is maintained.

Henceforth, all variables are to be regarded as nondimensional through
some appropriate scaling and we shall assume, without further mention,
that for any bounded regular map $\phi : M \to K$, boundedness includes
 the auxiliary condition such that on the image $\phi(M)$ in $K$
 we have $\lim_{|x| \to \infty} u_i^\alpha,u_{ij}^\alpha, A_{ijk}^\alpha =0$.
Although in some cases we need only $\lim_{|x| \to \infty} u_b^\alpha =0$
($b=1, \dots ,n$), we do this because we wish, on occasion, to make use
of certain fundamental solutions of the heat equation which, under this
auxiliary condition in combination with suitable Cauchy data,
will gaurantee that they are uniquely defined.

Consider $\phi : M \to K$ to be a bounded regular map that
annihilates the $N'$ functions $F^A \in \mathcal{F} (K)$ given by
\begin{equation}
F^A = v_\eta^A - \sum_{b=1}^n v_{bb}^A \label{me310}
\end{equation}
Under the action of the pullback of $\phi$ we have
\begin{equation}
\frac{\partial \,^*v^A}{\partial
\eta} - \triangle \,^*v^A = 0 \quad \text{on }M \label{me320}
\end{equation}
where $\,^*v^A = \phi^* v^A$. Suppose also that
\begin{equation}
 \,^*v^A |_{\eta=0^+} = \widetilde{v}^A \quad \text{on }\mathbb{R}^n
 \times (0,t_0) \label{me330}
\end{equation}
for some bounded $\widetilde{v}^A \in \mathcal{F}
(\mathbb{R}^n \times (0,t_0))$. The system (\ref{me320}) and
(\ref{me330}) defines a Cauchy problem for the heat equation on
space/scale $\mathbb{R}^n \times (0,\eta_0)$ whose coordinates are
$(x,\eta)$ and where the time $t \in (0,t_0)$ enters the problem
only as a parameter. Along with the requirement of the boundedness
of $\phi$ and its auxiliary condition, for each $t \in (0,t_0)$
the unique solution $\,^*v^A$ of (\ref{me320}) and (\ref{me330})
can be written
\begin{equation}
\,^*v^A (x,t,\eta) = \int_{\mathbb{R}^n}
G(x-x',\eta) \widetilde{v}^A (x',t) \, d^n x' \label{me340}
\end{equation}
where
\begin{equation}
G(x,\eta) = (4 \pi \eta)^{-n/2} \exp [- |x|^2 / (4 \eta )]
\label{me350}
\end{equation}
is given in normalized form such that
 \begin{equation}
\int_{\mathbb{R}^n} G(x,\eta) \, d^n x = 1~. \label{me360}
\end{equation}
We use the the notation $d^n x = dx^1 \dots dx^n$. The solution map
of (\ref{me310}) along with the Cauchy data (\ref{me330}) and the
boundedness of $\phi$ defines each $\,^*v^A$ as the spatial
Gaussian filter of each $\widetilde{v}^A$. The scale parameter $\eta$ can
be expressed as
\begin{equation} \eta= \beta \delta^2 \label{me370}
\end{equation}
where $\delta$ is a (nondimensional) characteristic space scale
associated with the resolution and $\beta$ is a parameter that
controls the rate of damping of the filter.

Suppose that for the filter map $\phi$ defined above, $( \widetilde{v}^A
)$ turns out to be a limiting solution of the system $\phi^* P^A
|_{\eta=0^+} = 0$, which has the representation given by
(\ref{me105}). The $N'$ equations $\phi^* P^A |_{\eta=0^+} = 0$
define the fully resolved system of PDEs which is of interest to
us if it generates solutions $( \widetilde{v}^A )$ that are highly
fluctuating on $\mathbb{R}^n \times (0,t_0)$. We cannot in general
expect that $\,^*P^A = P^A (x^i , \,^*v^B , \,^*v_i^B ,
\,^*v_{ij}^B )$ will vanish everywhere on $M$. We therefore
introduce the residuals $\,^* r^A$ that quantify the deviation
from zero of $\phi^* P^A$ when $\eta > 0$. The following defines
the consistency condition from which immediately follow the exact
macroscopic equations that are satisfied by the filtered variables
and the exact residuals.

\begin{theorem}[Consistency] \label{thme1}
Let $P^A \in \mathcal{F} (K)$ have the representation
$P^A = P^A (x^i , v^B , v_i^B , v_{ij}^B )$ and suppose that there
 exist $\widetilde{v}^A \in \mathcal{F} (\mathbb{R}^n \times (0,t_0) )$,
  bounded on
$\mathbb{R}^n \times (0,t_0)$, such that
\begin{equation}
P^A \left ( x^i ,
\widetilde{v}^B , \frac{\partial \widetilde{v}^B}{\partial x^i} , 
\frac{\partial^2 \widetilde{v}^B}{\partial x^i \partial
x^j} \right ) = 0 \label{me410}
\end{equation}
Let $F^A \in \mathcal{F} (K)$ be
given by $F^A = v_\eta^A - \sum_{b=1}^n v_{bb}^A$ and let $\phi: M
\to K$ be a bounded regular map such that
\begin{equation}
\,^*F^A = 0 \quad
\text{on }M~, \quad \,^*v^A |_{\eta=0^+} = \widetilde{v}^A \quad \text{on
}\mathbb{R}^n \times (0,t_0) ~. \label{me420}
 \end{equation}
Let $E^A \in \mathcal{F} (K)$ be defined by
\begin{equation} E^A = r_\eta^A - \sum_{b=1}^n
r_{bb}^A - S^A \label{me421}
\end{equation}
where $S^A \in \mathcal{F} (K)$ is given by
\begin{equation}
S^A = (\mathcal{L} - \mathcal{W} ) P^A (x^i , v^B , v_i^B ,
v_{ij}^B ) ~. \label{me460}
 \end{equation}
and suppose that in addition
$\,^*r^A|_{\eta=0^+} =0$. Then on $M$
\begin{equation}
P^A \left ( x^i ,
\,^*v^B ,\frac{\partial \,^*v^B}{\partial x^i} ,
\frac{\partial^2 \,^*v^B}{\partial x^i \partial x^j}
\right ) + \,^*r^A = \int_0^\eta \int_{\mathbb{R}^n} G(x-x',\eta -
\eta') \,^*E^A (x',t,\eta') \, d^n x' \, d \eta' \label{me461}
\end{equation}
and hence
\begin{equation}
\Big| P^A \Big( x^i , \,^*v^B ,\frac{\partial \,^*v^B}{\partial
x^i} ,\frac{\partial^2 \,^*v^B}{\partial x^i \partial x^j} \Big) + \,^*r^A \Big| \leq
\eta \sup_{ x \in \mathbb{R}^n, \, \eta \in (0,\eta_0)} |^*E^A |
\label{me430}
\end{equation}
\end{theorem}


\begin{proof}
The system (\ref{me420}) defines the Cauchy problem for the heat
equation on space/scale $\mathbb{R}^n \times (0,\eta_0)$ whose
coordinates are $(x,\eta)$ and where the time $t \in (0,t_0)$
appears only as a parameter. Given the existence of $\widetilde{v}^A \in
\mathcal{F} (\mathbb{R}^n \times (0,t_0))$, bounded on $\mathbb{R}^n
\times (0,t_0)$, we are guaranteed that a bounded regular solution
map $\phi : M \to K$ associated with (\ref{me420}) exists and that
$(\,^*v^A)$ is unique. We also know that each $\,^*v^A$ is
explicitly given by (\ref{me340}) and hence defines $\,^*v^A$ as
the spatial Gaussian filter of $\widetilde{v}^A$. We note that the
regular map $\phi : M \to K$ itself is not unique because it does
not as yet place any constraints on the coordinate components $r^A
,r_i^A, r_{ij}^A$ of $K$ other than through the pullback
identities (\ref{cm130}), (\ref{cm140}) and the associated
involutive conditions on the distribution $\mathcal{D}$ (\ref{cm170}).

Consider the point set $\Xi = \{ p \in K :  F^A =0 \}$.
The map $\phi : M \to K$ will map $M$ to $(n+2)$-dimensional integral 
manifolds of the distribution $\mathcal{D}$ on $K$ contained in the point 
set $\Xi$.
We have on $K$
\begin{equation}
\begin{aligned}
E^A & =  r_\eta^A - \sum_{b=1}^n r_{bb}^A - S^A  \\
    & =  (\mathcal{V}_\eta - \mathcal{L} ) r^A - (\mathcal{L} - \mathcal{W}) P^A  \\
    & =  (\mathcal{V}_\eta - \mathcal{L} ) (P^A+r^A) 
    - (\mathcal{V}_\eta - \mathcal{W}) P^A
\end{aligned}\label{me500}
\end{equation} On $\Xi$ we can set $\mathcal{V}_\eta = \mathcal{W}$ and hence
(\ref{me500}) reduces to
\begin{equation}
(\mathcal{V}_\eta - \mathcal{L} ) (P^A+r^A) = E^A
\quad \text{on } \Xi \label{me501}
\end{equation}
Introducing $Z^A \in \mathcal{F}
(K)$ such that $Z^A = P^A + r^A$ it follows that under the action
of the pullback of $\phi$, (\ref{me501}) becomes
\begin{equation}
\frac{\partial \,^*Z^A}{\partial \eta} - \triangle \,^*Z^A = \,^*E^A \quad \text{on }M
\label{me570}
\end{equation}
Since $\phi^* P^A |_{\eta=0^+} = 0$ and $\,^*r^A
|_{\eta=0^+} = 0$ we also have
\begin{equation}
\,^*Z^A |_{\eta=0^+} = 0 \quad
\text{on }\mathbb{R}^n \times (0,t_0) \label{me580}
\end{equation}
The system (\ref{me570}), (\ref{me580}) defines a Cauchy problem for the
nonhomogeneous heat equation on space/scale $\mathbb{R}^n \times
(0,\eta_0)$ whose coordinates are $(x,\eta)$ and where the time $t
\in (0,t_0)$ enters the problem only as a parameter. Given that
$\phi$ is a bounded regular map we can write the solution of
(\ref{me570}) and (\ref{me580}) as \cite{Sob64}
 \begin{equation}
 \,^*Z^A (x,t,\eta) = \int_0^\eta \int_{\mathbb{R}^n} G(x-x',\eta - \eta')
\,^*E^A (x',t,\eta') \, d^n x' \, d \eta' \label{me590}
\end{equation}
where $G(x,\eta)$ is defined by (\ref{me350}). This is the identity
(\ref{me461}). It follows that
\begin{equation}
\begin{aligned}
|^*Z^A | & =  | \int_0^\eta \int_{\mathbb{R}^n} G(x-x',\eta -\eta')
\,^*E^A (x',t,\eta') \, d^n x' \, d \eta' | \\
 & \leq  \int_0^\eta \int_{\mathbb{R}^n} G(x-x',\eta -\eta') |^*E^A
 (x',t,\eta')| \, d^n x' \, d \eta'  \\
 & \leq\Big( \sup_{x \in \mathbb{R}^n ,\, \eta \in (0,\eta_0)}|^*E^A |\Big)
\int_0^\eta \int_{\mathbb{R}^n} G(x-x',\eta -\eta') \, d^n x' \, d \eta'
\end{aligned}\label{me600}
\end{equation}
and the inequality (\ref{me430}) follows from (\ref{me360})
\end{proof}

The choice of $r^A$ that renders $E^A =0$ on the point set $\Xi$ leads
us to the exact macroscopic equations that are satisfied by the
filtered variables $\,^*v^A$ and the exact residuals $\,^*r^A$:
\begin{gather}
P^A \left ( x^i ,^*v^B ,\frac{\partial \,^*v^B}{\partial x^i} ,
\frac{\partial^2 \,^*v^B}{\partial x^i \partial x^j} \right ) + \,^*r^A  
=  0 \quad \text{on }M
\label{me700} \\
\frac{\partial \,^*r^A}{\partial \eta} - \triangle \,^*r^A  =  \,^*S^A \quad \text{on }
M \label{me710} \\
\,^*r^A |_{\eta=0+}  =  0 \quad \text{on } \mathbb{R}^n \times
(0,t_0) \label{me720}
\end{gather}
where $S^A \in \mathcal{F} (K)$ is given by
\begin{equation}
S^A = (\mathcal{L} - \mathcal{W} ) P^A (x^i , v^B , v_i^B , v_{ij}^B ) ~.
\label{me730}
\end{equation}
The fully resolved variables are not explicitly contained in the macroscopic
equations but are implied as a limiting solution of (\ref{me700}) as
$\eta \to 0^+$.
As a result the macroscopic system (\ref{me700}), (\ref{me710}) and
(\ref{me720}), along with the identity (\ref{me730}), have certain
useful properties for the purposes of application.
In practical application we wish to generate solutions of the filtered
variables only on a single scale slice $M_{\eta=\text{const}}$ without
the need to access the fully resolved variables.
The presence of the term $\partial_\eta \,^*r^A$ in (\ref{me710}) is the only
remaining obstacle to reach this objective.
To overcome this obstacle we resort to approximation methods.

\section{Approximation of the Residual}

As before consider the point set $\Xi = \{ p \in K :  F^A =0 \}$ and let
$\phi : M \to K$ be the bounded regular map defined in the consistency
theorem.
Since $\phi$ is regular it will map $M$ to $(n+2)$-dimensional integral
manifolds of the distribution $\mathcal{D}$ on $K$ contained in the point 
set $\Xi$.
As mentioned in the proof of the the consistency theorem, $\phi$ is not
unique because it does not as yet place any constraints on the coordinate
components $r^A ,r_i^A, r_{ij}^A$ of $K$ other than through the pullback
identities (\ref{cm130}), (\ref{cm140}) and the associated involutive
conditions on the distribution $\mathcal{D}$ (\ref{cm170}).
Consider the point set $\Xi' \subset \Xi$ such that
$\Xi' = \{ p \in \Xi : r^A - \eta \mathcal{L} r^A - \eta S^A =0 \}$.
Under the action of the pullback of $\phi$ to the constraints
$r^A - \eta \mathcal{L} r^A - \eta S^A =0$ on $\Xi'$ we have
\begin{equation}
\triangle \,^*r^A - \frac{\,^*r^A}{\eta} + \,^*S^A = 0 \quad \text{on }M
\label{a010}
\end{equation}
Under the definition of $E^A \in \mathcal{F} (K)$ given by (\ref{me421})
we obtain
\begin{equation}
E^A = r_\eta^A - \frac{r^A}{\eta}  \quad \text{on }\Xi'
\label{a050}
\end{equation}
Noting that with the additional constraint on $\phi$ which requires
 that $\,^*r^A |_{\eta=0^+} = 0$, under the action of the pullback of
$\phi$ (\ref{a050}) becomes
\begin{equation}
\,^*E^A = \frac{\partial \,^*r^A}{\partial \eta} - \frac{\,^*r^A}{\eta}
= \frac{\partial \,^*r^A}{\partial \eta} - \frac{\,^*r^A 
- \,^*r^A |_{\eta=0^+} }{\eta}
= \frac{\eta}{2} \left [ \frac{\partial^2 \,^*r^A}{\partial \eta^2} \right]_{\eta
\in (0, \eta_0 )}
\label{a160}
\end{equation}
Hence
\begin{equation}
|^*E^A | \leq \frac{\eta}{2} \sup_{\eta \in (0,\eta_0)} | \,^*r_{\eta \eta}^A |
\label{a170}
\end{equation}
Since $\phi$ is a bounded regular map, as defined in the Section 1,
it follows from the consistency theorem that there is a constant
$C>0$ such that
\begin{equation}
| P^A \left ( x^i , \,^*v^B ,\frac{\partial \,^*v^B}{\partial x^i} ,
\frac{\partial^2 \,^*v^B}{\partial x^i \partial x^j} \right ) + \,^*r^A | \leq C \eta^2
\label{a070}
\end{equation}
This demonstrates that for the residual approximation based on (\ref{a010})
the exact filtered variables will satisfy the macroscopic system of
 PDEs (\ref{me700}) to a consistency error $O(\eta^2)$.
On the other hand if we force (\ref{me700}) using a residual that is
not exact (i.e. not satisfying (\ref{me710}) and (\ref{me720})) we
cannot expect that $\phi$ will be the filter map associated with
(\ref{me420}).

Let $\widehat{\phi} : M \to K$ be a bounded regular map such that 
\begin{equation}
P^A \left ( x^i , \,^*\widehat{v}^B ,\frac{\partial \,^*\widehat{v}^B}{\partial x^i}
,\frac{\partial^2 \,^*\widehat{v}^B}{\partial x^i \partial x^j} \right ) 
+ \,^*\widehat{r}^A = 0
\quad \text{on }M ~, \label{a100} 
\end{equation} 
where $\,^*\widehat{v}^A =
\widehat{\phi}^* v^A$, $\,^*\widehat{r}^A = \widehat{\phi}^* r^A$ and
$\,^*\widehat{r}^A$ is an approximation of the exact residual. Suppose
that in addition $\,^*\widehat{r}^A |_{\eta=0^+} = 0$ on $\mathbb{R}^n
\times (0,t_0)$. The system (\ref{a100}) can be thought of as the
$N'$ PDEs which generate the $N'$ dependent variables
$\,^*\widehat{v}^A$. We note that the $\,^*\widehat{v}^A$ cannot be the
exact filters of the of the
 fully resolved variables $\widetilde{v}^A$, i.e. the map $\widehat{\phi}$ cannot
annihilate (\ref{me310}), because we have forced (\ref{a100}) while using
residuals that are not exact.
To generate the $N'$ approximations for the residuals for the above
approximation
\begin{equation}
\triangle \,^*\widehat{r}^A - \frac{\,^*\widehat{r}^A}{\eta} + \,^*\widehat{S}^A =0
 \quad \text{on } M
\label{a110}
\end{equation}
where $\,^*\widehat{S}^A = \widehat{\phi}^* S^A$ and
\begin{equation}
S^A = (\mathcal{L} - \mathcal{W}) P^A
\label{a130}
\end{equation}
We have seen above that with a residual equation error of $O(\eta)$ we
obtain a consistency error of $O(\eta^2)$.
While $\,^*\widehat{v}^A$ cannot be the exact filters of the of the fully
resolved variables $\widetilde{v}^A$ we expect that with a consistency error
of $O(\eta^2)$ $\,^*\widehat{v}^A$ will be reasonably good approximations of
the exact filtered variables $\,^*v^A$.
To obtain an estimate for the error $\,^*\widehat{v}^A - \,^*v^A$ is of course
a stronger validation of the residual approximation and is very much related,
from a geometrical point of view, to the magnitude of the vector field
$\mathcal{V}_\eta - \mathcal{W}$ on the image $\widehat{\phi} (M)$ contained 
in $K$.
Our main objective here is to demonstrate the usefulness of the consistency
error alone to examine certain empirically based residual models.
It will be seen that consistency will be adequate for this purpose and the
 stronger validation by way of the error just mentioned will be presented
 elsewhere.

\section{Application}

We consider the equations that describe the motion of an
incompressible and inviscid fluid. Because we also wish to examine
certain empirical gradient models in the next section, we will
augment these equations with an equation that describes the
transport of a single conservative solute within the fluid medium.
Before proceeding a few points need to be made: Any consideration
of the Euler or Navier-Stokes equations in the present context may
appear problematic because of the current open question on the
existence of regular solutions, particularly in three spatial
dimensions. Because of their wide use in modelling of complex
flows it seems appropriate that they receive some attention
although the motivation for the ideas presented here is much
wider. It is important to keep in mind throughout that the
restriction $\widetilde{v}^A \in \mathcal{F} (\mathbb{R}^n \times (0,t_0))$ is
made for brevity rather than by necessity. If working in the class
of smooth solutions is too restrictive the consistency theorem
could be modified as follows: We replace $\widetilde{v}^A \in \mathcal{F}
(\mathbb{R}^n \times (0,t_0))$ with $\widetilde{v}^A \in L_1^{loc}
(\mathbb{R}^n)$ for each $t \in (0,t_0)$ and require that
(\ref{me410}) be satisfied in a generalized sense (see for
instance \cite{Lad85}). The solution of (\ref{me420}) will still
be given by (\ref{me340}) but the Cauchy data condition $\,^*v^A
|_{\eta=0^+} = \widetilde{v}^A$ will be satisfied a.e. on $\mathbb{R}^n$
in the limit as $\eta \to 0^+$ for each $t \in (0,t_0)$. The
macroscopic system of equations (\ref{me700})-(\ref{me720}), along
with the identity (\ref{me730}) for the source term of the
residual equation, will still hold.

Let $v^b$ ($b=1, \dots ,n$) be the placeholders on $K$ for the $n$
filtered velocity components and $v^{n+1}=p$ be the placeholder on $K$
for the filtered pressure.
We set $v^{n+2} = \omega$ to denote the mass fraction of some conservative
 solute in the fluid medium.
Here $N' = n+2$. Define
\begin{gather}
P^a  =  v^a_t + \sum_{b=1}^n v^b v_b^a + p_a \quad 1 \leq a \leq n
\label{gm100} \\
P^{n+1}  =  \sum_{b=1}^n v_b^b \label{gm110} \\
P^{n+2}  =  \omega_t + \sum_{b=1}^n v^b \omega_b \label{gm111}
\end{gather}
where we shall use the notation $p_i=v_i^{n+1}$, $p_{ij}=v_{ij}^{n+1}$,
$\omega_i=v_i^{n+2}$, and $\omega_{ij}=v_{ij}^{n+2}$.
A calculation based on (\ref{me460}) and (\ref{gm100})-(\ref{gm111}) yields
\begin{gather}
S^a  =  2 \sum_{b,c=1}^n v_c^b v_{bc}^a  \quad 1 \leq a \leq n \label{gm120} \\
S^{n+1}  =  0 \label{gm116} \\
S^{n+2}  =  2 \sum_{b,c=1}^n v_c^b \omega_{bc} \label{gm121}
\end{gather}
We should note the following: For a viscous fluid we would introduce
the term $- \sum_{b=1}^n v_{bb}^a/Re$ on the right hand side of (\ref{gm100}),
where $Re$ is the Reynolds number.
Similarly we may also include a molecular diffusion term
$- \kappa \sum_{b=1}^n \omega_{bb}$ on the right hand side of (\ref{gm111}),
where $\kappa$ is the molecular diffusion coefficient.
Both these terms have no effect on the residual equation source terms
and (\ref{gm120})-(\ref{gm121}) will remain unchanged.

Since the source term $S^{n+1}=0$ on $K$ the residual for the continuity
 equation will vanish.
Under the action of the pullback of $\widehat{\phi}$, the system (\ref{a100})
and (\ref{a110}) for this application
becomes
\begin{gather}
\frac{\partial \,^*\widehat{v}^a}{\partial t} + \sum_{b=1}^n \,^*\widehat{v}^b
 \frac{\partial \,^*\widehat{v}^a}{\partial x^b} 
 + \frac{\partial \,^*\widehat{p}}{\partial x^a}
 + \,^*\widehat{r}^a =  0 \quad 1 \leq a \leq n \label{a620} \\
\frac{\partial \,^*\widehat{v}^b}{\partial x^b}  =  0 \label{a621} \\
\frac{\partial \,^*\widehat{\omega}}{\partial t} + \sum_{b=1}^n \,^*\widehat{v}^b
\frac{\partial \,^*\widehat{\omega}}{\partial x^b} + \,^*\widehat{r}^{n+2} 
 =  0 \label{a622} \\
\triangle \,^*\widehat{r}^a - \frac{\,^*\widehat{r}^a}{\eta} + \,^*\widehat{S}^a
 =  0  \quad 1 \leq a \leq n \label{a623} \\
\triangle \,^*\widehat{r}^{n+2} - \frac{\,^*\widehat{r}^{n+2}}{\eta}
+ \,^*\widehat{S}^{n+2}  =  0 \label{a624}
\end{gather}
where
\begin{gather}
\,^*\widehat{S}^a  =  2 \sum_{b,c=1}^n \frac{\partial ~}{\partial x^b}
\left ( \frac{\partial \,^*\widehat{v}^b}{\partial x^c}
\frac{\partial \,^*\widehat{v}^a}{\partial x^c} \right )
 \quad 1 \leq a \leq n \label{a640} \\
\,^*\widehat{S}^{n+2}  =  2 \sum_{b,c=1}^n \frac{\partial ~}{\partial x^b}
 \left ( \frac{\partial \,^*\widehat{v}^b}{\partial x^c}
 \frac{\partial \,^*\widehat{\omega}}{\partial x^c} \right )
  \label{a641}
\end{gather}
and use has been made of the form invariance of  the continuity
equation. Since the system (\ref{a620})-(\ref{a641}) contains no
terms involving $\partial_\eta$ and no reference to the fully resolved
variables, we can seek a solution on any desired scale slice
$M_{\eta=\text{const}}$. The choice of the value of the scale
parameter $\eta$ will be dictated by the level of refinement in
the spatial discretization used in the numerical solution scheme.
The numerical scheme involves also a temporal discretization of
the evolution equations (\ref{a620}) and (\ref{a622}) from which
one obtains an update of the velocities and the solute mass
fraction at each timestep. For incompressible flows a staggered
grid is often used and the velocities updated along with the
pressure using a projection method (see for instance
\cite{Pey83}). Within each timestep a finite difference
approximation of elliptic equations (\ref{a623}) and (\ref{a624})
are solved iteratively to obtain an update of the residuals. The
procedure is repeated until a desired time is reached.

Numerical experiments have been conducted for applications  in
fluid mechanics by solving (\ref{a100}) and (\ref{a110}) on single
scale slices using finite difference methods \cite{Pan04}. The
computations using this approach are found to be stable despite
the complex flow patterns that emerge during the breakdown of
hydrodynamic stability. For compressible flows, energy balance
studies on numerical solutions based on the system (\ref{a100})
and (\ref{a110}) indicate the presence of some intrinsic property
that conserves the total energy of the fluid system. This is
particularly interesting given the evidence that the model
accommodates the flow of energy both to and from the smaller
scales, i.e. internal energy can flow up from the unresolvable
scale into the macroscopic scale where it appears as kinetic
energy. It is also observed that numerical solutions generated
independently on increasing scale slices result in a corresponding
increase in the smoothing of finer scale fluctuations in the
complex flow patterns, indicating that filtering is occurring.

\section{Empirical Gradient Models}

We shall investigate a class of residual models in common  use
today in the field of fluid mechanics. We consider again the flow
of an incompressible and inviscid fluid with a single solute. As
such we use the prescription given by (\ref{gm100})-(\ref{gm111}).
The identities (\ref{gm120})-(\ref{gm121}) still hold under any
residual approximation used.

Empirical gradient models take the form
\begin{gather}
r^a  =  - \eta \nu \mathcal{L} v^a \quad 1 \leq a \leq n \label{gm200} \\
r^{n+1}  =  0 \label{gm201} \\
r^{n+2}  =  - \eta \kappa \mathcal{L} \omega \label{gm202}
\end{gather}
 where $\nu, \kappa \in \mathcal{F} (K)$. The coefficients $\nu$ and $\kappa$
are empirically based and are dependent on constant parameters
(assumed to be measurable). In application these residual models
are explicitly defined and hence can be inserted directly into the
system (\ref{a620}) and (\ref{a622}). The equations (\ref{a623})
and (\ref{a624}) are of course not applicable here.

The philosophy behind these models is that large scale
fluctuations interact in a diffusion like fashion analagous to
those at the molecular level. The residuals (\ref{gm200}) are
meant to capture  the turbulence stresses in the fluid and the
residual (\ref{gm202}) is meant to capture the solute dispersion
in the turbulent fluid medium. The form invariance of the
continuity equation is assumed (note that the form invariance of
the continuity equation of the previous section is not assumed but
follows immediately from the vanishing of the source term
$S^{n+1}$). A useful coverage on current practices and
applications of these type of residual models can be found in
\cite{Liu01}.

In early applications of these models the coefficients $\mu$ and
$\kappa$ were assumed to be constants. Due to their failure to
predict observations later variants of these models were proposed
such that $\mu$ and $\kappa$ are some functions of the dependent
variables. We need not investigate any particular case here
because it can be shown by consistency alone that empirical
gradient type models are flawed in the general case given above
and cannot be regarded as approximations of the residuals in any
reasonable sense.

Consider the point set $\Xi = \{ p \in K :  F^A =0 \}$. The map
$\phi : M \to K$ of the consistency theorem  will map $M$ to
$(n+2)$-dimensional integral manifolds of the distribution $\mathcal{D}$
on $K$ contained in the point set $\Xi$. On the basis of
(\ref{gm200})-(\ref{gm202}) we have on $K$
\begin{gather}
E^a  =  - \nu \mathcal{L} v^a - S^a - \eta (\mathcal{V}_\eta 
- \mathcal{L}) (\nu \mathcal{L} v^a) \label{gm210} \\
E^{n+1}  =  0  \label{gm211} \\
E^{n+2}  =  - \kappa \mathcal{L} \omega - S^{n+2}
- \eta (\mathcal{V}_\eta - \mathcal{L}) (\kappa \mathcal{L} \omega) \label{gm212}
\end{gather}
On the point set $\Xi$, (\ref{gm210}) and (\ref{gm212}) can be written
\begin{gather}
E^a  =  - \nu \mathcal{L} v^a - S^a - \eta (\mathcal{W} - \mathcal{L}) 
(\nu \mathcal{L} v^a) \label{gm220} \\
E^{n+2}  =  - \kappa \mathcal{L} \omega - S^{n+2} - \eta (\mathcal{W} 
- \mathcal{L}) (\kappa \mathcal{L} \omega) \label{gm222}
\end{gather}
where we use the fact that  $\mathcal{V}_\eta =\mathcal{W}$ on the point
set $\Xi$. Note that in the case that the coefficients $\nu$ and
$\kappa$ are constants the $O(\eta)$ terms vanish under the action of 
the pullback of $\phi$. However,
whether they are assumed as constants or not, the troublesome
zeroth order terms remain. Noting the identities (\ref{gm120}) and
(\ref{gm121}) we require that
\begin{gather}
\nu \sum_{b=1}^n v_{bb}^a  \sim  - 2 \sum_{b,c=1}^n v_c^b v_{bc}^a
 \quad 1 \leq a \leq n \label{gm300} \\
\kappa \sum_{b=1}^n \omega_{bb} \sim  - 2 \sum_{b,c=1}^n v_c^b
\omega_{bc} \label{gm310}
\end{gather}
where we use the symbol $\sim$ to denote  equality to within an
error $O(\eta)$. We see that there do not exist choices for the
coefficients $\nu,\kappa \in \mathcal{F} (K)$ that maintain the
diffusion like character of the residuals on which the formulas
(\ref{gm200}) and (\ref{gm202}) are motivated.

Under the action of the pullback of $\phi$, the residual equation
error $\,^*E^A = O_s(1)$ with respect to the scale parameter $\eta$.
It follows from the consistency theorem that we have a consistency
error $ |^*P^A + \,^*r^A | = O_s(\eta)$. Since the magnitudes of the
residuals $\,^*r^a$ and $\,^*r^{n+2}$ are $O_s(\eta)$ we can expect
that in application there will be some contamination of the
proposed residual model (\ref{gm200})-(\ref{gm202}) by unwanted
terms. It can also be expected that the free parameters that are
often contained in explicit formulas for the coefficients $\nu$
and $\kappa$ are not measurable and any fine tuning of these
parameters will not improve the order of magnitude estimate of the
consistency error.

Attempts have been made to generalize the empirical gradient models by way of
\begin{gather}
r^a  =  - \eta \sum_{b,c=1}^n \mathcal{V}_b ( \nu^{bc} \mathcal{V}_c v^a )
\quad 1 \leq a \leq n \label{gm320} \\
r^{n+2}  =  - \eta \sum_{b,c=1}^n \mathcal{V}_b ( \kappa^{bc} \mathcal{V}_c \omega )
 \label{gm330}
\end{gather}
where now  $\nu^{bc}, \kappa^{bc} \in \mathcal{F}(K)$ are second order
symmetric tensors. Choices of their functional dependence on the
filtered variables and their partial derivatives have been tried
containing additional free parameters that are adjusted to the
specific application being considered. However, repeating the
above procedure for these residual models reveals that the only
choices that can ensure a consistency error $O(\eta^2)$ are given
by
\begin{gather}
\nu^{bc}  \sim  - 2 v_c^b \label{gm340} \\
\kappa^{bc} \sim  - 2 v_c^b \label{gm350}
\end{gather}
The residual based on (\ref{gm320}) and (\ref{gm340}) is well
known in large eddy simulation and has been derived by series
expansion of the integral representation of the Gaussian filter
(\ref{me340}). Numerical experiments indicate that this model
correlates reasonably well with the fluid turbulence stresses and
strains inferred by direct numerical simualtions \cite{Pruett01}.
Unfortunately, computations of turbulent fluid flows using this
residual model are highly unstable and therefore the model is of
little use in numerical simulation.

While focus has been given here on spatial Gaussian filters  by
way of the spatial Laplacian differential operator (\ref{me140}),
it is important to keep in mind observations made in
\cite{Pruett01} that the residual models are very much dependent
on the type of filter being used. However, there appears no reason
why the spatial Laplacian differential operator could not be
replaced by any other elliptic operator on space/time. This will
widen the range of the type of filters that can be studied by the
methods presented here and could include spatial, spacio-temporal
or temporal only filters. While for these general elliptic
operators explicit representations for the consistency error as in
(\ref{me461}) may not be possible it can be expected that similar
order of magnitude estimates of the form (\ref{me430}) can be
obtained.

\begin{thebibliography}{00}

\bibitem{Ede90} D. G. B. Edelen, Alternatives to the Cartan-Kahler Theorem,
{\it Nonlinear Analysis, Theory, Methods and Applications} {\bf 15},
(1990), pp. 765-786.

\bibitem{Ede01} D. G. B. Edelen and J. Wang,
{\it Transformation Methods for Nonlinear Partial Differential Equations},
(World Scientific, 1992).

\bibitem{Lad85} O. A. Ladyzhenskaya, {\it The Boundary Value Problems
 of Mathematical Physics}, (Springer-Verlag, 1985).

\bibitem{Liu01} C. Liu, DNS/LES Perspective, In {\it DNS/LES Progress
and Challenges}, (Edited by C. Liu {\it et al.}), pp. i-vi,
(Greyden Press, 2001).

\bibitem{Pan04} G. Pantelis,
modelling nonlinear systems by scale considerations,
{\it Mathematical and Computer Modelling}, {\bf 7-8}, (2004), 797-807.

\bibitem{Pey83} R. Peyret and T. D. Taylor, {\it Computational Methods
for Fluid Flow}, (Springer-Verlag, (1983).

\bibitem{Pruett01} C. D. Pruett and N. A. Adams, \emph{A priori}
analysis of three subgrid-scale models for one-parameter families
of filters, {\it Phys. Fluids}, {\bf 12}, (2000), 1133.

\bibitem{Sob64} S. L. Sobolev, {\it Partial Differential Equations
of Mathematical Physics}, (Pergamon Press, 1964).

\bibitem{War01} F. W. Warner,
{\it Foundations of Differentiable Manifolds and Lie Groups},
(Springer Verlag, New Yok, 1971).

\end{thebibliography}

\section*{Addendum posted on July 2, 2007.}

Both equations \eqref{me105} and \eqref{me410} are meant to represent 
the expanded form of the statement $P^A|_{\eta=0^+}=0$.
To avoid confusion \eqref{me105} and \eqref{me410} should be written
\[
P^A( x,t,0, \widehat{v}^B , \frac{\partial \widehat{v}^B}{\partial x^i} , 
\frac{\partial^2 \widehat{v}^B}{\partial x^i \partial x^j} ) = 0.
\]
\smallskip

The statement following equation \eqref{me500} should read: \\
On $\Xi$ we can set $\mathcal{V}_\eta P^A = \mathcal{W} P^A$ and
 hence \eqref{me500} reduces to $\dots$
\smallskip

The statement following equations \eqref{gm220}--\eqref{gm222} should read:
 where we use the fact that 
 $\mathcal{V}_\eta (\nu \mathcal{L} v^a) = \mathcal{W} (\nu \mathcal{L} v^a)$ 
and $\mathcal{V}_\eta (\kappa \mathcal{L} \omega) 
= \mathcal{W} (\kappa \mathcal{L} \omega)$ on the point set $\Xi$.

\noindent End of Addendum.
\end{document}