compare MPDATA simulations of decaying turbulence
in a triply-periodic cube with the pseudospectral results
of Herring and Kerr [3]. For two values of physical
viscosity corresponding to direct numerical simulations
(DNS) where all dynamical scales are resolved, the MP-
DATA and the pseudospectral solutions compare closely
in all integral measures.
We also compare simulations using zero viscosity. The
pseudospectral simulation, in this Euler-equation limit
of Navier-Stokes, shows an enstrophy blowup at finite
time. The MPDATA simulation tracks the pseudospec-
tral results for a while, but does not show any blowup
of enstrophy. The pseudospectral result is unphysical
- all physical flows exhibit viscous dissipation at some
finite length scale. In contrast, enstrophy in the MP-
DATA simulation remains uniformly bounded, and the
solution appears physically reasonable. However, the
question remains whether the MPDATA simulation is
the result of a well-posed physical problem, and if so,
what this problem is.
To address this question - the central issue of the pa-
per - we analyze the zero-viscosity MPDATA results
by combining theoretical arguments with computational
experiments. In section III. we derive a relation between
the numerical energy spectra at different resolutions. In
section IV. we validate this relation computationally. A
theoretical consequence is the existence of an asymp-
totic spectrum in the continuum limit. For zero explicit
viscosity, we interpret the asymptotic spectrum as the
high Reynolds number (Re) limit of viscous flows [4].
We elaborate, and then summarize our conclusions in
section V.
II. MPDATA
The simulations presented in this paper all employed
the nonoscillatory fluid solver MPDATA. We would
emphasize that that implicit turbulence modeling ap-
pears to be a property of all nonoscillatory finite vol-
ume (NFV) schemes, however MPDATA does have some
unique features. MPDATA was developed originally for
applications in meteorology. Since it is likely unfamiliar
to the aeronautics community, we take this opportunity
to summarize some of its features. The interested reader
can find a comprehensive description of MPDATA in
[15], including both the underlying concepts and the
details of implementation.
MPDATA in its basic form is sign preserving, but
not monotonicity preserving. For meteorological appli-
cations, we have found that sign preservation is often
sufficient and leads to a less diffusive solution. MP-
DATA is fully second-order accurate and conservative.
A variety of options have been documented [15] that
extend MPDATA to full monotonicity preservation, to
third-order accuracy, and to fields that do not preserve
sign (of which the most important is momentum).

Unlike most NFV methods, MPDATA is not based
on the idea of flux limiting. Instead it is based more
directly on upwinding. In practical terms, the algo-
rithm consists of a series of donor cell steps; the first
step provides a first-order accurate solution while sub-
sequent steps compensate higher-order errors as identi-
fied from a modified equation analysis. One important
consequence of this approach is that MPDATA is fully
multidimensional - i.e., has no spatial splitting errors
- which implies significantly reduced mesh dependence.
MPDATA is a full fluid solver. In analyzing the trun-
cation error of approximations to the momentum equa-
tion, one finds error terms that depend on the interac-
tion of the advection with the forcing terms, including
the pressure gradient. In implementations of NFV al-
gorithms that treat advection separately from the forc-
ings, this error is uncompensated, reducing the order of
accuracy of the solution and potentially leading to os-
cillations and even instability (see [12]). In MPDATA
we compensate this error effectively by integrating the
forcing terms along a flow trajectory rather than at a
point.
MPDATA is implemented in the 3D program EULAG
for simulating rotating, stratified flows in complex ge-
ometries, [14] [16]. The name EULAG alludes to the ca-
pability to solve the fluid equations in either an Eulerian
(flux form) or a Lagrangian (advective form) framework.
The latter uses a semi-Lagrangian algorithm in which
an MPDATA type scheme performs as an interpolation
routine. However all the simulations in this paper use
the Eulerian framework. EULAG can be run for incom-
pressible or anelastic fluids; in either case, we solve an
elliptic equation for pressure using a preconditioned gen-
eralized conjugate residual solver [13]. EULAG is fully
parallelized using message-passing and runs efficiently
on a variety of platforms.
III. Theory
In this section, we summarize the theoretical results
of [7] and describe an extension to support our analysis.
111.1 Background
In [7], the authors describe a rationale for implicit tur-
bulence modeling. Their analysis begins by deriving the
modified equation for MPDATA applied to 1D Burg-
ers' equation. Among the third-order truncation terms,
there appears a nonlinear dispersive term of the form
Azexuxx. The authors then construct the governing
equation for a finite volume of Burgers' fluid. These
equations are derived from the point equations, but are
different due to the nonlinearity of the latter, a fact
that has long been appreciated by theorists and mod-
elers studying turbulence. What is unexpected is that
a straightforward and justifiable derivation of the finite
volume equations leads directly to nonlinearly disper-

2