# Identifying the body force from partial observations of a 2D incompressible
velocity field

Aseel Farhat<EMAIL_ADDRESS>Florida State University
Department of Mathematics  Adam Larios<EMAIL_ADDRESS>University of Nebraska
Department of Mathematics  Vincent R. Martinez<EMAIL_ADDRESS>Hunter College, City University of New York
Department of Mathematics & Statistics
Graduate Center, City University of New York
Department of Mathematics  Jared P. Whitehead<EMAIL_ADDRESS>Brigham Young University
Department of Mathematics

###### Abstract

Using limited observations of the velocity field of the two-dimensional
Navier-Stokes equations, we successfully reconstruct the steady body force
that drives the flow. The number of observed data points is less than 10% of
the number of modes that describes the full flow field, indicating that the
method introduced here is capable of identifying complicated forcing
mechanisms from a relatively small collection of observations. In addition to
demonstrating the efficacy of this method on turbulent flow data generated by
simulations of the two-dimensional Navier-Stokes equations, we also rigorously
justify convergence of the derived algorithm. Beyond the practical
applicability of such an algorithm, the reliance of this method on the
dynamical evolution of the system yields physical insight into the turbulent
cascade.

††preprint: APS/123-QED

##  I Introduction

In an ancient dusty magic school, where wizards cast spells by moving their
hands in complex arcane patterns, a young student, Henri Pother, hides behind
a curtain watching a master wizard perform his most secret spell.
Unfortunately for Henri, he can only see through a few small holes in the
curtain. To make matters worse, the wizard is invisible! All Henri can see is
the motion of the dust in the air around the wizard’s invisible hands. Can
Henri reconstruct the movement of the wizard’s hands from only limited
observations of the air currents?

The above make-believe scenario is an analogue for the real-world setting of
weather modeling: We have equations that model weather, but we do not know how
these equations are forced (we don’t know the movement of the wizard’s hands).
On the other hand, we observe parts of the evolution on large-scales (the
motion of the dust particles from which one can systematically construct an
approximation of the surrounding air current velocity, for instance through
Particle Image Velocimetry [10]). The real question that we want to answer is
then: Can we systematically reconstruct the forcing on the system from only
sparse observations of the velocity? In the present paper, we show that the
answer is “yes,” at least in the context of the 2D incompressible Navier-
Stokes equations with periodic boundary conditions.111We consider periodic
boundary conditions here as a proof-of-concept, but there seems to be no
fundamental computational difficulty in adapting this algorithm to a wide
variety of practical settings, including the 3D case with physical boundary
conditions. This slightly more nuanced and computationally intensive setting
will be investigated in a future work. This is accomplished through an
algorithm that enforces exponential decay of state errors on observational
scales. We mathematically prove that that this algorithm reconstructs non-
potential driving forces that are external to the system, that is, state-
independent. We probe the practical efficacy of this algorithm in a scenario
of turbulent flow. The reference flow field is generated by a highly resolved
numerical solution of the 2D Navier-Stokes equation with a randomized low-mode
force to which the algorithm is agnostic. When state observations are made
across all directly-forced scales, we observe both model and state error to
achieve machine-precision convergence to the true values exponentially fast in
time.

Adequate control of fluid flow is a difficult problem that has been
extensively studied (see, e.g., [60, 61]) particularly in light of modern
advances in computing. Motivation for such studies is found in a number of
disciplines in the physical and engineering sciences, where the ability to
control either classical Newtonian fluids or complex non-Newtonian fluids is
of significant interest. In [6] a feedback control mechanism (motivated byt
the original work of [62]) was exploited to smoothly merge incoming partial
observations with a dynamic model given by a system of partial differential
equations that is simultaneously evolved forward in time, allowing for a
controlled estimate of the true full state of the unobserved system. This is
often referred to as the Azouani-Olson-Titi (AOT) algorithm or Continuous Data
Assimilation (CDA). The AOT approach has by now been studied in many fluid
systems, yielding both numerical evidence of its efficacy and rigorous
analysis to justify the convergence of the generated approximating state to
its true value. The current investigation exploits this feedback control
approach to produce a rigorously justified approach to the determination of a
steady, but apriori unknown, external driving force from partial observations
of the flow field in two-dimensional turbulence.

A common assumption found in many of the works that incorporate feedback
control in fluid systems is the existence and complete knowledge of a perfect
model for the dynamical evolution of the fluid itself, that is, the underlying
physical model is known exactly beforehand. Notable exceptions to this
framework are the studies performed in [44] and [23], where the true dynamic
model is known up to a set of unidentified scalar parameters. In [44], the
canonical Rayleigh-Bénard convection setting is investigated where the exact
value of the Prandtl number (the nondimensional quantity denoting the ratio of
kinematic viscosity to thermal diffusivity) is unknown. Observations extracted
from an identical setup, but with a different value of the Prandtl number are
then rigorously shown to drive the simulated solution toward the true state up
to a constant error that is Prandtl number-dependent. This analysis is then
supplemented with a battery of numerical simulations that verify the
mathematical conclusions and probe the efficacy of the method by exploring
situations outside of the regimes where convergence can mathematically be
guaranteed. In [23], a similar study is performed on the periodic 2D Navier-
Stokes equations with an unknown viscosity, but a further advancement is
introduced which exploits the numerical observation that the state error
relaxes to be relatively constant. This remaining state error is then used to
propose an approximate value to the true viscosity; this procedure can then be
iterated to produce a sequence of approximating values that reliably converge
to the true viscosity value up to numerical precision. In a series of recent
papers [24, 86, 85], the underlying mechanisms leading to the observed
convergence were identified in a mathematically salient way to supply an
analytical proof of this convergence. The current investigation takes this a
step further, providing both rigorous justification and numerical evidence for
an algorithm capable of recovering the entire inhomogeneous term in a fluid
system.

A variation to the approach for parameter recovery that was originally
introduced in [23] was proposed in [92]. There, a more principled perspective
based on enforcing a form of null-controllability was taken to develop a
continuously-updating parameter algorithm in contrast to the sequentially-
updating algorithm in [23]. This variation was numerically studied in the
context of the one-dimensional Kuramoto-Sivashinsky equation and was seen to
reliably infer multiple unknown parameters appearing in the system in a
concurrent fashion. The algorithm of interest in the present study is
primarily based on the one introduced in [92]. In contrast to [92], however,
the current article, firstly, negotiates a situation that possesses richer
dynamical behavior in two-dimensional turbulent flow, in addition to one that
contains significantly more parameters to infer in having to determine all
coefficients of the forcing function, the number of which may potentially be
very large, and secondly, supplies the first rigorous proof of convergence of
a parameter recovery algorithm that is based on the null-controllability
perspective from [92].

At this point, we remark on three particular recent works that directly relate
to the task of reconstructing driving forces in a system, namely [46, 5, 2]
and [85]. In the first set of references, the problem of inferring temperature
from velocity measurements [46, 5] or velocity from temperature measurements
[5, 2] in the context of Rayleigh-Bénard convection (RBC) is studied. In [46],
a downscaling algorithm for recovering temperature from velocity measurements
in the planar setting (2D) is analytically proven to converge [46] and its
efficacy was studied in several numerical tests in [5]. In [5, 2] the more
difficult case of inferring velocity from temperature measurements was also
studied numerically with mixed results; the former carried out numerical
experiments in the 2D setting and there identified a possible mechanism for
generating asynchrony between the approximating velocity field and the
reference velocity field through a Kelvin-Helmholtz instabilty, while the
latter carried out their experiments in the 3D setting and identified more
optimistic scenarios. In contrast to the present article, recovery of the
temperature in RBC may be viewed as a form of forcing recovery. However, one
notable and very important difference between [2] and the present article is
that in the setting of RBC, the driving force (temperature) present in the
velocity field is linearly coupled to the velocity, which the analysis in [46]
relies on in a crucial way. In the situation studied here, the driving force
is completely external to the system and no relation whatsoever between the
state and the driving force is posited apriori. Indeed, the success of the
algorithm studied here is independent of the existence of any relation between
the state and forcing.

In [85], the convergence properties are studied of an algorithm for inferring
$\mathbf{f}$ in (1) that is similar to the one treated in the current article.
However, in [85], the state variable is reconstructed by directly using the
observations as the low-mode approximation and subsequently generating the
unobserved high-modes with an incorrect forcing via a nudging-based method. A
correct forcing is then iteratively generated by exploiting the equation
itself as a nonlinear filtration of errors. In contrast, the algorithm studied
in the current article generates an approximation to the state variable by
enforcing exponential convergence of the state error at the observational
scale. As mentioned above, this idea was initially introduced in [92] in the
context of the one-dimensional Kuramoto-Sivashinsky equation. From this point
of view, we are able to supply a unifying perspective between the two
algorithms obtained by different derivations, that distinguishes each
algorithm within a family of variations of one common approach. We refer the
reader to Remark IV.1 for further details. Another notable remark is that the
work of [85] only provided theoretical justification of the algorithm studied
therein and did not provide any numerical evidence in that work. In this
article, we carry out both theoretical and numerical studies of the proposed
algorithm, and note that we were unable to successfully (at least up to the
level of machine precision) numerically implement the algorithm proposed in
[85]. Indeed, the current article complements the works referenced above
firstly, by proposing an algorithm that continuously processes observations
for reconstructing the force, secondly, by rigorously demonstrating its
convergence via mathematical proof, and thirdly, providing substantive
numerical evidence of its efficacy and robustness.

Although we focus on the idealized test-case of 2D isotropic turbulence, we
note that the methods presented here can be readily adapted to other dynamical
systems, such as the Magnetohydrodynamic (MHD) equations, the Boussinesq
equations of ocean flow, the setting of Rayleigh-Benard convection,
geophysical flows such as the surface quasi-geostrophic (SQG) equations or the
primitive equations of the ocean, pattern formation equations such as the
Cahn-Hilliard or Allen-Cahn equations, equations of flame fronts and crystal
growth, such as the Kuramoto-Sivashisnky equation, and many others. In
addition, it is straight-forward (at least, at the algorithmic level) to adapt
the methods presented here to the setting of physical boundaries. The 3D case
could also be considered, although with the usual caveats regarding the
analytical difficulties of the 3D Navier-Stokes equations (see, e.g., [17] for
some results on the AOT algorithm in 3D).

Lastly, we remark that the approaches described above for discovery of model
parameters while simultaneously reconstructing the true state of the system is
readily comparable to other data-driven model recovery techniques such as
SINDy (see [22], for example, or [40, 88] and the references therein). One
particular benefit of the feedback control mechanism exploited here for the
purpose of recovering the unknown forcing function in 2D turbulence, is that
the parameter estimation is performed on-the-fly in the fashion of continuous
data assimilation and without any need for post-processing. Indeed, the data
is inserted dynamically into the model in real-time, which is, in turn,
incrementally updated to obtain the correct parameters that characterize the
external forcing.

The rest of this article is organized as follows: Section II will establish
notation and the necessary preliminary results required by the analysis that
follows. Section III provides a formal derivation of the forcing update
algorithm. Section V provides the rigorous analysis that justifies the
algorithm in 2D. Section VI describes the numerical simulations that were
performed that demonstrate how well the algorithm performs in 2D, and Section
VII concludes with a takeaway message and potential for future work. A reader
that is less interested in the technical details of the mathematical analysis
may be well-served to focus on Sections III, VI and VII.

##  II The equations of motion and their mathematical setting

In this section, we provide a short description of the mathematical model of
interest and the functional setting in which we carry out our convergence
analysis. The development of the underlying algorithm in Section III can be
followed without some of the definitions provided here, but these definitions
and preliminary results are necessary for the rigorous analysis performed in
Section V that completely justifies the algorithm.

We consider an incompressible fluid in a periodic domain, $\Omega=[0,L]^{d}$,
of length $L>0$ (all of the analysis and simulations presented below assume
that $d=2$, but the heuristic derivation of the algorithm is independent of
the dimension $d$), subject to a time-independent external body force
$\mathbf{f}$, which is mean-free over $\Omega$. The evolution of the velocity
field is governed by the Navier-Stokes equations given by

$\partial_{t}\mathbf{u}+\mathbf{u}\cdot\nabla\mathbf{u}=-\nabla
p+\nu\nabla^{2}\mathbf{u}+\mathbf{f},\quad\nabla\cdot\mathbf{u}=0,$ (1)

where fluid density has been normalized to unity, the kinematic viscosity,
$\nu$, is known, and $\mathbf{u},p$ are assumed mean-free over $\Omega$. We
assume that the body force is unknown and that the velocity field is partially
observed. The main objective is to reconstruct the non-conservative component
of $\mathbf{f}$ at observational scales by leveraging the model jointly with
the observations. For convenience, let us therefore assume that
$\nabla\cdotp\mathbf{f}=0$.

An important non-dimensional quantity that serves as a proxy for the Reynolds
number when the Navier-Stokes equations is externally forced is the Grashof
number, which we will denote by $G$. In two-dimensions, $G$ is defined by

$\displaystyle G=\frac{\norm{\bf
f}_{2}}{(\kappa_{0}\nu)^{2}},\quad\text{where}\ \kappa_{0}:=\frac{2\pi}{L},$
(2)

where $\|\hskip 1.0pt\cdotp\|_{2}$ is the $L^{2}$ norm over $\Omega$. Note
that $\kappa_{0}$ is the smallest eigenvalue of the negative Laplacian
$-\Delta$, while in the context of three-dimensional turbulence, the Grashof
number is known to be related to the Reynolds number, Re, of the flow as
$G\sim\text{Re}^{2}$ (see [35]). For our purposes, we will assume that the
body force is twice weakly differentiable over $\Omega$, all of whose weak
partial derivatives are square-integrable over $\Omega$. It is known that the
corresponding initial value problem for the 2D NSE (1) is well-posed in this
setting and, moreover, the dynamics for the system possesses a finite-
dimensional global attractor when ${\bf f}$ is time-independent, see e.g.,
[32, 95, 102]; additional properties of the solution which are relevant to the
analysis performed later in the article is provided in Section VIII.

##  III Derivation of the algorithm and Heuristics for convergence

The method for parameter recovery studied in this paper makes crucial use of a
state-recovery algorithm originally introduced by Azouani, Olson, and Titi
(AOT) [6] in the context of continuous data assimilation for the 2D NSE
equations; the system associated to this method is introduced below in (4). In
the seminal work of AOT, a feedback control paradigm is introduced for
obtaining an approximation of the state using a sufficiently large, but finite
number of observations of the true solution, collected continuously in time.
This approximation asymptotically converges to the true solution corresponding
to the observed data at an exponential rate. This particular algorithm
differed from previously studied algorithms through the manner in which
observations were assimilated; in [6] observations were inserted into the
dynamical model as an external forcing term that enforced relaxation to the
true state of the system, whereas previous algorithms such as that studied in
[90] had inserted observations into the model by directly replacing dynamical
terms. The use of such schemes for numerical weather prediction has been and
continues to be extensively studied (see [1, 8, 100, 105] for instance). Since
its introduction to the mathematical fluid dynamics community, this algorithm
has been a topic of intense activity (see [3, 4, 5, 11, 14, 15, 16, 18, 19,
21, 23, 24, 26, 27, 73, 45, 46, 43, 47, 48, 49, 50, 54, 56, 57, 58, 64, 65,
66, 68, 69, 76, 75, 82, 78, 79, 80, 81, 38, 39, 84, 89]), and can be
considered a nonlinear complement to other data assimilation methods such as
variations on the Kalman filter (see e.g., [101]), variational methods,
particle filtering, and several other related methods (see [106, 94, 28, 42]
for just a few examples).

To describe our algorithm for parameter recovery, we shall denote the
observable state of the velocity by $I_{h}(\mathbf{u})$, where $I_{h}$ is an
(autonomous) bounded linear projection operator whose output is a suitable
interpolation of the observed data into the phase space of the system (1),
which consists only of divergence-free vector fields. In general, $I_{h}$ is
often referred to as an interpolant observable operator. In the context of
incompressible flows considered here, $I_{h}$ is understood to perform an
interpolation of the observed data, then orthogonally projects the result onto
solenoidal vector fields. Hence, $I_{h}^{2}=I_{h}$, $I_{h}\nabla r=0$, and
$\nabla\cdotp I_{h}\mathbf{v}=0$, for any sufficiently smooth scalar field $r$
and vector field $\mathbf{v}$. Here, $h>0$ quantifies the resolution of the
observational field in such a way that $\lim_{h\rightarrow
0}I_{h}(\mathbf{u})=\mathbf{u}$. Specific examples of the manner of
interpolation encoded in $I_{h}$ include Lagrangian interpolation of nodal
values or local averages which are distributed uniformly across the domain at
a mesh-size $h$, but also include large-scale filtering such as projection
onto finitely many Galerkin modes up to a cut-off frequency $1/h$.

We will assume that $I_{h}$ and its complementary projection $J_{h}=I-I_{h}$
satisfy the following boundedness properties:

$\displaystyle\|I_{h}(\partial^{\alpha}\phi)\|_{2}^{2}\leq{c_{m}^{2}}{h^{-2m}}\|\phi\|_{2}^{2},\quad\|J_{h}(\phi)\|_{2}^{2}\leq
C_{m}^{2}h^{2m}\sum_{|\alpha|=k}\|\partial^{\alpha}\phi\|_{2}^{2},$ (3)

for any multi-index $|\alpha|=m$ with $m\geq 0$, for some constants
$c_{m},C_{m}>0$, independent of $h$; such inequalities are satisfied in the
case where $I_{h}$ is given by spectral projection. In this special case,
$I_{h}$ also satisfies $I_{h}\partial^{\alpha}=\partial^{\alpha}I_{h}$, so
that $I_{h}(\nabla^{2}J_{h}(\mathbf{w}))=0$. For the remainder of the
manuscript, $I_{h}$ is therefore assumed to represent spectral projection,
that is, projection onto Fourier modes $|\mathbf{k}|\leq h^{-1}$, for
$\mathbf{k}\in(\kappa_{0}\mathbb{Z})^{2}$.

The unobserved scales of the flow are then dynamically approximated by
directly inserting the observations, $I_{h}(\mathbf{u})$, into (1) as a
feedback-control term, resulting in the system

$\partial_{t}\mathbf{v}+\mathbf{v}\cdot\nabla\mathbf{v}=-\nabla
q+\nu\nabla^{2}\mathbf{v}+\mathbf{g}-\mu I_{h}(\mathbf{v})+\mu
I_{h}(\mathbf{u}),\quad\nabla\cdotp\mathbf{v}=0,$ (4)

where $\mu$ denotes a tuning parameter, often referred to as the nudging
coefficient, and $\mathbf{g}$ is a divergence-free vector field chosen as a
putative approximation to the unknown forcing $\mathbf{f}$. A fundamental
property of the system is that in the absence of model error, i.e.,
$\mathbf{g}=\mathbf{f}$, then the approximating flow field, $\mathbf{v}$,
asymptotically synchronizes with the true flow field, $\mathbf{u}$, provided
that the true flow field is observed through sufficiently small scales, i.e.,
$h\ll 1$, and $\mu$ is appropriately tuned [7].

In the presence of model error, this fundamental property can be leveraged to
develop an ansatz for $\mathbf{f}$ by enforcing relaxation of the state error,
$\mathbf{w}=\mathbf{v}-\mathbf{u}$, on the observed scales. Indeed, the
evolution of $\mathbf{w}$ on the observed scales is governed by

$\displaystyle\partial_{t}I_{h}(\mathbf{w})=I_{h}\left(-\mathbf{v}\cdotp\nabla\mathbf{v}+\nu\nabla^{2}\mathbf{v}+\mathbf{g}-\partial_{t}\mathbf{u}\right)-\mu
I_{h}(\mathbf{w}).$

Assuming that $\mathbf{g}$ can be chosen instantaneously such that

$\displaystyle\mathbf{g}=\partial_{t}I_{h}(\mathbf{u})+I_{h}\left(\mathbf{v}\cdotp\nabla\mathbf{v}-\nu\nabla^{2}I_{h}\mathbf{u}-\nu\nabla^{2}J_{h}\mathbf{v}\right),$
(5)

where $J_{h}=I-I_{h}$ denotes the complementary projection, it follows that
$\partial_{t}I_{h}(\mathbf{w})+\mu I_{h}(\mathbf{w})-\nu
I_{h}\nabla^{2}J_{h}(\mathbf{w})=0$. Note that this choice critically relies
on replacing the nudged velocity field $\mathbf{v}$ in the Laplacian term
($\nabla^{2}$) directly with the observed data $I_{h}\mathbf{u}$ combined with
the unobserved data from the nudged solution $J_{h}\mathbf{v}$. This is in
direct contrast to simply identifying the forcing update as a function of
$\mathbf{v}$. This ‘direct replacement’ strategy employed here is necessary
for the rigorous convergence of the forcing to the true value, but as noted
below does not appear to be necessary in practice.

By the orthogonality of $I_{h}(\mathbf{w})$ and $J_{h}(\mathbf{w})$, it
follows that the energy balance at observational scales one step forward in
time satisfies

$\displaystyle\|I_{h}(\mathbf{w})(t+\Delta t)\|_{2}=e^{-\mu\Delta
t}\|I_{h}(\mathbf{w}(t))\|_{2}.$ (6)

In particular, exponential relaxation of $I_{h}(\mathbf{w})$ is enforced. The
utility of this choice is readily seen; setting
$\mathbf{h}=\mathbf{g}-\mathbf{f}$ we find that

$\displaystyle I_{h}(\mathbf{h})$
$\displaystyle=I_{h}\left[J_{h}(\mathbf{w})\cdotp\nabla
J_{h}(\mathbf{w})+J_{h}(\mathbf{w})\cdotp\nabla\mathbf{u}+\mathbf{u}\cdotp\nabla
J_{h}(\mathbf{w})\right]$ $\displaystyle\quad-\nu
I_{h}(\nabla^{2}J_{h}\mathbf{w})+\mathcal{O}(I_{h}(\mathbf{w})).$ (7)

Under the assumption that $\mu$ is sufficiently large relative to the
observational density, $h$, and time-step, $\Delta t$, the resulting state
error one step forward, $\mathbf{w}(t+\Delta t)$ will satisfy the following
estimates (see the Appendix for details):

$\displaystyle\|\mathbf{w}(t+\Delta
t)\|_{2}\leq\mathcal{O}\left(\frac{\|\mathbf{h}(t)\|_{2}}{\mu}\right),\quad\|\nabla\mathbf{w}(t+\Delta
t)\|_{2}\leq\mathcal{O}\left(\frac{\|\mathbf{h}(t)\|_{2}}{\sqrt{\mu\nu}}\right).$
(8)

Further assuming that $I_{h}\mathbf{f}=\mathbf{f}$ (the forcing function lives
in the observation space even though it is unobservable itself) it follows
that $I_{h}\mathbf{h}=\mathbf{h}$. Then the model error in (8) is essentially
controlled by

$\displaystyle\|I_{h}\mathbf{h}(t+\Delta t)\|_{2}$
$\displaystyle\leq\mathcal{O}(\|I_{h}(\mathbf{w}(t+\Delta
t))\|_{2})+\mathcal{O}\left(\frac{\|\mathbf{h}(t)\|_{2}}{\sqrt{\mu\nu}}\right).$

Hence, upon invoking (6) for a sufficiently large time step $\Delta t$ and the
fact that $I_{h}(\mathbf{h})=\mathbf{h}$, one may then deduce that for $\mu$
chosen appropriately large

$\displaystyle\|\mathbf{h}(t+\Delta
t)\|_{2}\leq\frac{1}{2}\|\mathbf{h}(t)\|_{2}.$ (9)

Owing to (6), we see that (9) implies exponential convergence of $\mathbf{h}$
to zero upon further iteration.

##  IV Definition of the algorithm

The discussion above lends itself to an implementable algorithm: At stage $1$,
the algorithm is initialized with a pair
$(\mathbf{v}_{0}^{-1},\mathbf{g}^{0})$ (the superscripts are indices, not
exponentiation) that prescribes an initial velocity and a force for (4).
Integration of (4) forward-in-time produces
$\mathbf{v}^{0}(t)=\mathbf{v}^{0}(t;\mathbf{v}_{0}^{-1},\mathbf{g}^{0})$, for
all times $t\in I_{0}=[0,\infty)$. After a transient period that allows the
state error, $\mathbf{w}^{0}=\mathbf{v}^{0}-\mathbf{u}$ to establish a balance
with the model error, $\Delta\mathbf{g}^{0}=\mathbf{g}^{0}-\mathbf{f}$, a new
estimate of the forcing function modified from (5) is calculated at time
$t_{1}\gg 0$:

$\displaystyle\mathbf{g}^{1}$
$\displaystyle=\partial_{t}I_{h}(\mathbf{u})+I_{h}\left(\mathbf{v}^{0}\cdotp\nabla\mathbf{v}^{0}-\nu\nabla^{2}(I_{h}\mathbf{u}+J_{h}\mathbf{v}^{0})\right).$
(10)

This yields a new approximation, $\mathbf{g}^{1}$, to the force. Notably,
$\mathbf{g}^{1}=I_{h}(\mathbf{g}^{1})$. This concludes the initial cycle of
the algorithm. The process is then repeated: suppose that $\mathbf{g}^{\ell}$
have been produced in this way at stages $\ell=1,\dots,k$ at times
$t_{\ell}\gg t_{\ell-1}$, respectively, and that
$\mathbf{g}^{\ell}=I_{h}\mathbf{g}^{\ell}$. At stage $k$, (4) is re-
initialized at time $t=t_{k}$, with the pair
$(\mathbf{v}^{k-1}(t_{k}),\mathbf{g}^{k})$; this produces an approximating
velocity,
$\mathbf{v}^{k}(t)=\mathbf{v}^{k}(t;\mathbf{v}^{k-1}(t_{k}),\mathbf{g}^{k})$,
over the time interval $I_{k}=[t_{k},\infty)$. After a sufficiently long
interval of time for $\mathbf{w}^{k}=\mathbf{v}^{k}-\mathbf{u}$ to achieve
balance with the model error,
$\Delta\mathbf{g}^{k}=\mathbf{g}^{k}-\mathbf{f}$, the forcing function
estimate is updated according to

$\displaystyle\mathbf{g}^{k+1}=\partial_{t}I_{h}(\mathbf{u})+I_{h}\left(\mathbf{v}^{k}\cdotp\nabla\mathbf{v}^{k}-\nu\nabla^{2}(I_{h}\mathbf{u}+J_{h}\mathbf{v}^{k})\right),$
(11)

at time $t_{k+1}\gg t_{k}$ to produce $\mathbf{g}^{k+1}$; again
$\mathbf{g}^{k+1}=I_{h}\mathbf{g}^{k+1}$. In this way a sequence of
approximating forces is generated,
$\mathbf{g}^{1},\mathbf{g}^{2},\mathbf{g}^{3},\dots$, as observations are
assimilated.

######  Remark IV.1.

As mentioned in the introduction, the derivation of the algorithm presented in
the current article allows one to obtain many other algorithms. In particular,
the algorithm studied in [85] is obtained by replacing the $\mathbf{v}$ with
$I_{h}\mathbf{u}+(I-I_{h})\mathbf{v}$. All other variations of the algorithm
may be obtained by applying this substitution in some, rather than all of the
terms where $\mathbf{v}$ appears in (11). In this way, the point of view
introduced in [92] and subsequently adopted here allows for a more general
perspective.

One may notice in (11) that $I_{h}\mathbf{u}$ is being used directly in the
viscous term. At the moment, the analysis performed below is unable to replace
this specific term $I_{h}\mathbf{v}$. We believe this to be a technical point,
however, since the numerical tests performed below suggest that one should
nevertheless achieve convergence without this ‘direct replacement’ employed.
Nevertheless, establishing convergence of the version of the algorithm
corresponding to (5) remains an open problem.

######  Remark IV.2.

In the setting presented above, we have assumed that the velocity field is
mean-free. This assumption is typically justified by recalling the Galilean
invariance of the Navier-Stokes equations that allows one to shift the
reference frame of the fluid in order to equivalently view it with respect to
one that is mean-free. Coupled with the observation that the mean velocity is
preserved under the free-evolution of the Navier-Stokes equation, that is, in
the absence of external force, it then suffices to consider the scenario where
the initial fluid velocity is mean-free. In the presence of an external force,
the mean-free assumption remains valid provided that the external force is
also assumed to be mean-free. This assumption is not essential for the
analysis provided here, but does make the rigorous proof much simpler to
establish.

Our setting therefore implicitly assumes that 1) the force is restricted to
the observational subspace, i.e., $I_{h}\mathbf{f}=\mathbf{f}$, 2) that the
force is periodic in space, and 3) that the mean-force is zero, i.e.,
$\bar{\mathbf{f}}=|\Omega|^{-1}\int_{\Omega}\mathbf{f}(x)dx=0$. The third
assumption is not more restrictive than the assumption that one has access to
all large-scale motions of a fluid down to a particular length scale $h>0$.
Indeed, since $I_{h}$ is assumed to be the projection onto the subspace of
Fourier modes up to wave-number $|k|\leq h^{-1}$, where
$k\in(\kappa_{0}\mathbb{Z})^{2}$, and observations are made continuously in
time, if $\mathbf{f}$ happens to possess a non-zero mean, then one easily sees
that $\frac{d}{dt}I_{h>1}u=I_{h>1}\mathbf{f}=\hat{\mathbf{f}}(0)$. In
particular, in the setting we study, it is implicitly assumed that
$\hat{\mathbf{f}}(0)$ is accessible to the user from direct observations.

One particular scenario of interest in which the fluid possesses a nonzero
mean is that of shear flows with nonzero mean profile. Due to the remarks
above, this scenario is indeed within the purview of the setting studied here.
However, it is not interesting in the absence of boundaries, which is the main
setting of this article. The difficult problem of rigorously studying the
effects of boundaries, as they arise for instance from no-slip boundary
conditions, is relegated to future work. We simply point out that the choice
of considering the setting of 2D isotropic turbulence represents an important
first step towards understanding the problem of inferring unknown external
forces from low-mode observations in turbulent flows.

##  V Rigorous Convergence Analysis

The precise version of (8) asserts that the stage $k$ state error satisfies

$\displaystyle\|\mathbf{w}^{k}(t)\|_{2}\leq\mathcal{O}\left(\frac{\|\Delta\mathbf{g}^{k}\|_{2}}{\mu}\right),\quad\|\nabla\mathbf{w}^{k}(t)\|_{2}\leq\mathcal{O}\left(\frac{\|\Delta\mathbf{g}^{k}\|_{2}}{\sqrt{\mu\nu}}\right),$
(12)

for all $t\geq t_{k+1}$, where $t_{k+1}\gg t_{k}$, provided that $\mu
h^{2}\leq\nu$ and $\mu$ is chosen sufficiently large, depending only on the
maximum magnitude of the energy, enstrophy, and palenstrophy of $\mathbf{u}$;
a proof of the first inequality in (12) is provided in Section VIII, while the
proof of the second can be found in [85].

Now, from (1), (4), and (5), the model error at stage $k+1$ can be expanded
as:

$\displaystyle\Delta\mathbf{g}^{k+1}$
$\displaystyle=\left[I_{h}\left(\mathbf{w}^{k}\cdotp\nabla\mathbf{w}^{k}\right)+I_{h}\left(\mathbf{u}\cdotp\nabla\mathbf{w}^{k}\right)+I_{h}\left(\mathbf{w}^{k}\cdotp\nabla\mathbf{u}\right)\right]_{t=t_{k+1}}.$
(13)

Note that the fact that $\mathbf{f}$ is confined to the subspace spanned by
the observations, that is, $\mathbf{f}=I_{h}\mathbf{f}$, has been applied. It
is then clear from (13) that $I_{h}\Delta\mathbf{g}^{k+1}=\mathbf{g}^{k+1}$.
To assess the size of $\Delta\mathbf{g}^{k+1}$, observe that integration by
parts and the properties $I_{h}^{2}=I_{h}$,
$I_{h}\mathbf{g}^{k}=\mathbf{g}^{k}$, for all $k$, yield

$\displaystyle\langle
I_{h}\left(\mathbf{w}^{k}\cdotp\nabla\mathbf{w}^{k}\right),\Delta\mathbf{g}^{k+1}\rangle$
$\displaystyle=-\langle\mathbf{w}^{k}\cdotp\nabla\Delta\mathbf{g}^{k+1},\mathbf{w}^{k}\rangle,$
$\displaystyle\langle
I_{h}(\mathbf{u}\cdotp\nabla\mathbf{w}^{k}),\Delta\mathbf{g}^{k+1}\rangle$
$\displaystyle=-\langle\mathbf{u}\cdotp\nabla\Delta\mathbf{g}^{k+1},\mathbf{w}^{k}\rangle,$
$\displaystyle\langle
I_{h}(\mathbf{w}^{k}\cdotp\nabla\mathbf{u}),\Delta\mathbf{g}^{k+1}\rangle$
$\displaystyle=-\langle\mathbf{w}^{k}\cdotp\nabla\Delta\mathbf{g}^{k+1},\mathbf{u}\rangle,$

where brackets denote the $L^{2}$–inner product over $\Omega$. Now, upon
taking the $L^{2}$–inner product of (13) with $\Delta\mathbf{g}^{k+1}$,
integrating by parts, invoking the fact that $\Delta\mathbf{g}^{k+1}$ is
solenoidal, then applying the Cauchy-Schwarz inequality, Hölder’s inequality
(see the Appendix), and (3), yields

$\displaystyle\|\Delta\mathbf{g}^{k+1}\|_{2}^{2}$
$\displaystyle\leq\left(\|\mathbf{w}^{k}\|_{4}^{2}+2\|\mathbf{u}\|_{\infty}\|\mathbf{w}^{k}\|_{2}\right)\|\nabla\Delta\mathbf{g}^{k+1}\|_{2}.$

Let $R_{k}$ denote the supremum in time of
$\||\nabla|^{k}\mathbf{u}(t)\|_{2}$. By interpolation, the Cauchy-Schwarz
inequality, and (3) it follows that

$\displaystyle\|\Delta\mathbf{g}^{k+1}\|_{2}\leq\frac{c_{0}}{h}\left(\|\nabla\mathbf{w}^{k}\|_{2}\|\mathbf{w}^{k}\|_{2}+\sqrt{R_{2}R_{0}}\|\mathbf{w}^{k}\|_{2}\right),$

for some universal non-dimensional constant $c_{0}>0$ that depends on
$c_{1},c_{2}$ from (3) and the constants of interpolation; we refer to the
Appendix for additional details. Finally, (12) implies

$\displaystyle\|\Delta\mathbf{g}^{k+1}\|_{2}\leq
c_{0}\frac{\sqrt{R_{2}R_{0}}}{\mu
h}\left(\frac{\|\Delta\mathbf{g}^{k}\|_{2}}{\sqrt{\mu\nu
R_{2}R_{0}}}+1\right)\|\Delta\mathbf{g}^{k}\|_{2},$ (14)

where $c_{0}$ now denotes a larger, but still non-dimensional, constant.
Conditions ensuring the convergence of $\mathbf{g}^{k}$ to the true force,
$\mathbf{f}$, may then be given by

$\displaystyle 2c_{0}(R_{2}R_{0})^{1/2}h\leq\mu
h^{2}\leq\nu,\quad\mu\nu\geq\frac{1}{2}\frac{\|\Delta\mathbf{g}^{0}\|_{2}}{\sqrt{R_{2}R_{0}}}.$
(15)

Indeed, under (15), one may then deduce from (14) that

$\displaystyle\|\Delta\mathbf{g}^{k+1}\|_{2}\leq\frac{1}{2}\|\Delta\mathbf{g}^{k}\|_{2},\quad\text{for
all}\ k\geq 0.$ (16)

Observe that (15) constitutes a non-trivial set of pairs $(\mu,h)$ whenever
the observational density, $h^{-1}$, is sufficiently large in relation to the
viscosity, energy ($R_{0}$), and palenstrophy ($R_{2}$), of the flow, and the
nudging parameter, $\mu$, is tuned appropriately large. Note that the choice
of $\mu$ ultimately depends on the true forcing, $\mathbf{f}$, through the
values of $R_{0},~{}R_{2}$ and the initial error. In practice, this is often
approximately known, for instance, in terms of the Reynolds number of the
flow; we refer the reader to the Appendix for an additional comment on this
point.

##  VI Simulation Results

To demonstrate the utility of this algorithm in practice, we computationally
recover the forcing function for 2D forced incompressible Navier-Stokes.

  * •

The simulation was carried out in Matlab R2021a using a pseudo-spectral method
with explicit Euler time-stepping, respecting the CFL constraint and
2/3-dealiasing rule, and the mean-free condition enforced. The spatial
resolution was $2048^{2}$ (giving roughly $(2048\cdot 2/3)^{2}-1$ degrees of
freedom) on the domain $[-\pi,\pi)^{2}$ with time-step $\Delta t=0.0025$. The
reference $\mathbf{u}$ system was simulated at the stream-function level using
the Basdevant formula [9] to compute the nonlinearity efficiently. Namely, the
following equivalent form of the 2D Navier-Stokes equations were used:

$\psi_{t}=\nabla^{-2}[(\partial_{x}^{2}-\partial_{y}^{2})(u_{1}u_{2})+\partial_{xy}(u_{2}^{2}-u_{1}^{2})]+f^{\psi}$

where $\nabla^{\perp}:=\binom{-\partial_{y}}{\partial_{x}}$,
$\mathbf{u}=\binom{u_{1}}{u_{2}}=\nabla^{\perp}\psi$, and
$\mathbf{f}=\nabla^{\perp}f^{\psi}$. The equation for the data assimilation
$\mathbf{v}$ system was handled similarly.

  * •

The unobserved forcing function $f^{\psi}$ was determined by randomly (normal,
$N(0,1)$, distribution-seeded in Matlab with rng(0,'twister')) picking
amplitudes for wave-modes in Fourier space on an annulus with inner radius 16
and outer radius 64 (see FIG. 1 for an illustration of the forcing function).
Specifically, the forcing recovery algorithm was tasked to recover 11,672
(real-valued) unknowns of the forcing.

  * •

To ensure that the simulation is well-resolved, and also understand the
behavior at various wave lengths, we plot the (weighted) spectrum of various
functions. Namely, given a function $\rho$ where
$\rho(\mathbf{x})=\sum_{\mathbf{n}\in\mathbb{Z}^{2}}\widehat{\rho}_{\mathbf{n}}e^{i\mathbf{n}\cdot\mathbf{x}}$,
we denote its weighted spectrum by

$E_{s}[\rho](k):=(k+1)^{s}\Big{(}\sum_{\begin{subarray}{c}\mathbf{n}\in\mathbb{Z}^{2}\\\
k\leq|\mathbf{n}|<k+1\end{subarray}}|\widehat{\rho}_{\mathbf{n}}|^{2}\Big{)}^{1/2}$

where $s\in\\{0,1,2\\}$ is a weighting exponent for illustrative purposes, and
$k\in\mathbb{N}$ is referred to as the “wave number.” Roughly speaking, for
$\mathbf{u}=\nabla^{\perp}\psi$, $E_{1}(\psi)$ corresponds to the kinetic
energy spectrum (of $\mathbf{u}$) and $E_{2}(\psi)$ corresponds to the
enstrophy spectrum. Due to the mean-free condition, $E_{s}[h](0)=0$; hence it
is not plotted.

  * •

The viscosity was $\nu=10^{-4}$, and the Grashof number was $G=2.5\times
10^{6}$, leading to a solution with an energy spectrum resolved to machine
precision ($\approx 2.22\times 10^{-16}$) before the dealiasing cut-off at
$k=2*2048/3$.

  * •

Initial data for $\psi$ (i.e., for the $\mathbf{u}$ system) was given by a
simulation spun up from initially-zero data and run with the same forcing
until energy and enstrophy were roughly statistically stabilized; namely, time
$t=210$ (see FIG. 1). That is, no major qualitative differences were observed
in the spectrum for $t\geq 210$. It turns out that with these parameters, the
initial $\psi$ is fairly small; namely, $\|\psi\|_{L^{2}}\approx 0.0414346$.
However, this is comparable with initial data in other studies with similar
parameters, such as [58, 91].

(a) Visualization of reference force and initial data,
---
(bottom-left cut-away)
$\nabla^{2}f^{\psi}/\|\nabla^{2}f^{\psi}\|_{L^{\infty}}$,
(top-right cut-away) $\nabla^{2}\psi_{0}/\|\nabla^{2}\psi_{0}\|_{L^{\infty}}$.

(b) Weighted spectra of force $f^{\psi}$ and initial data $\psi_{0}$.

Figure 1:  Initial data and reference force. In (a), the normalized initial
vorticity $-\nabla^{2}\psi_{0}/\|\nabla^{2}\psi_{0}\|_{L^{\infty}}$ is shown
rather than the initial stream function $\psi_{0}$, to highlight details.
Similarly, $-\nabla^{2}f^{\psi}/\|\nabla^{2}f^{\psi}\|_{L^{\infty}}$ is shown
in (a) rather than $f^{\psi}$. For reference, we computed that
$\|\nabla^{2}\psi_{0}\|_{L^{\infty}}\approx 2.8256$ and
$\|\nabla^{2}f^{\psi}\|_{L^{\infty}}\approx 0.7443$. In (b), the weighted
spectra of the initial data and stream function are shown.

  * •

Initial data and forcing for the $\mathbf{v}$ system were identically zero.

  * •

The interpolation operator, $I_{h}$, was taken to be a projection onto the
“observed” Fourier modes of the stream function, that is, those on wave-modes
$k$, $0<|\mathbf{k}|_{\infty}\leq 64$. In particular, only $\approx 0.893\%$
of the unknowns in the solution were observed. The nudging parameter was taken
to be $\mu=1.9/\Delta t$.

  * •

The time-derivative that appears in the ansatz forcing function (11) was
computed via a forward-Euler finite difference approximation of
$\partial_{t}I_{h}(\mathbf{v})$.

  * •

The nonlinearity of the nudged equation was computed the same way as in the
solution of the ‘truth’ system (using $\mathbf{v}$ instead of $\mathbf{u}$),
and the same can be said for the forcing estimate in (11).

  * •

When we refer to the “$L^{2}$ Error” in figures or the text (say, the $L^{2}$
error between $\psi$ and $\psi_{DA}$), we specifically mean the following:

$\displaystyle\|\psi-\psi_{DA}\|_{L^{2}}:=\left(\int_{[-\pi,\pi)^{2}}|\psi(\mathbf{x})-\psi_{DA}(\mathbf{x})|^{2}\,d\mathbf{x}\right)^{1/2}$
(17)

In particular, no normalization factor is used. Of course, in simulations,
this is computed discretely (via Parseval’s identity), which we do via a
Riemann sum.

The forcing estimation was implemented in two different ways:

  1. 1.

As outlined here, a direct replacement strategy was used for the Laplacian
term, i.e. (11) was used. This is the method justified by the rigorous
analysis.

  2. 2.

Rather than replace the Laplacian term with a combination of the low modes
from the observed truth $I_{h}\mathbf{u}$ and the high modes of the nudged
system $J_{h}\mathbf{v}$, we also simply used the low modes of the nudged
system, i.e. replacing (11) with

$\mathbf{g}^{k+1}=\partial_{t}I_{h}(\mathbf{u})+I_{h}(\mathbf{v}^{k}\cdot\nabla\mathbf{v}^{k}-\nu\nabla^{2}\mathbf{v}^{k}).$
(18)

Although only the first case is rigorously justified, we found no qualitative
computational difference between these two schemes, i.e. convergence of both
the state and approximated forcing function were nearly identical between the
two approaches. This is demonstrated in FIG 3. In this Figure, the direct
replacement (DR) simulations are those where the forcing update is given by
(11) whereas the exact Laplacian refers to (18). Note that there is very
minimal difference between the convergence rates in either case.

(a) Observation of all forcing modes, i.e., $|\mathbf{k}|_{\infty}\leq 64$.
Solution and forcing converge exponentially fast.

(b) Observation of only forcing modes $\mathbf{k}$ with
$|\mathbf{k}|_{\infty}\leq 60$. Convergence is halted around error $\approx
10^{-4}$.

Figure 2:  Convergence plots for the exact Laplacian algorithm with force
updated continuously (i.e., at each time step) in various norms. Comparing (a)
and (b), it appears that observing all modes contained in the force gives
convergence to machine precision, while not observing only a small number of
modes leads to lack of convergence.

In addition to variations in the forcing update, we also updated the forcing
continuously (or at every discrete time step in the simulation) or at
uniformly-spaced time intervals. In all cases, the force and the solution
converged exponentially fast with very little distinguishing differences
amongst the different simulations, as shown in FIG. 3. We anticipate this is
because the relaxation time scale required between forcing updates is akin to
$1/\mu$ which for the selected values of $\mu$ is on the same order as the
actual time step of the numerical algorithm so that ‘continuous’ updating of
the forcing is on the requisite time scale anyway. We also observed the time-
evolution of the convergence of the energy spectrum of the errors of the
stream function and the forcing function (FIG. 4). These figures demonstrate
that the state (streamfunction in this case) converges on the observable
scales (below the cutoff of $|\mathbf{k}|_{\infty}=64$) almost
instantaneously, but the unobserved scales of the state, and the full forcing
function (which is everything in the observed range) converge much slower.

Figure 3:  (log-linear plot) $L^{2}$ errors vs. time. “DR” refers to the
direct replacement algorithm, while “EX” refers to the exact Laplacian
algorithm. “-c” indicates that the force was updated continuously (i.e., at
each time step), while “-d” indicates discrete force updates, every $0.25$
time units.

(a) Spectral error in stream function: $E_{0}[\psi-\psi_{DA}](k)$ vs. $k$.

(b) Spectral error in forcing: $E_{0}[f^{\psi}-f^{\psi}_{DA}](k)$ vs. $k$.

(c) Relative spectral error in solution:
$E_{0}[\psi-\psi_{DA}](k)/E_{0}[\psi](k)$ vs. $k$.

(d) Relative spectral error in forcing:
$E_{0}[f^{\psi}-f^{\psi}_{DA}](k)/E_{0}[f^{\psi}](k)$ vs. $k$.

Figure 4:  (log-log plot) Spectra333The reason for the difference in the
$x$-axes is that the force is generated from only a projection onto the
observational modes (namely, the modes $\mathbf{k}=(k_{1},k_{1})$ for which
$|\mathbf{k}|_{\infty}:=\max(|k_{1}|,|k_{2}|)\leq 64$). In two dimensions,
this means that for
$|\mathbf{k}|_{2}=\sqrt{|k_{1}|^{2}+|k_{2}|^{2}}>64\sqrt{2}\approx 90.5$, the
$\mathbf{k}^{\text{th}}$ modes will be zero. Since $\log(0)$ is undefined,
these do not show up on a log-log plot. of the error in the stream function
(a) and the error of the forcing function (b) at times $t=0.0,0.5,1.0$, and
also $t=2,4,6\ldots,40$ using direct replacement scheme. Colors move from blue
($t=0.0$) to red ($t=40.0)$. Vertical dashed line is the observational wave-
number cut-off at $|\mathbf{k}|_{\infty}=64$. Horizontal dotted line is
machine precision ($\epsilon\approx 2.22\times 10^{-16}$). (c) Relative error
in spectrum spectrum at times $t=0.0,0.5,1.0$, and also $t=2,4,6\ldots,80$
using direct replacement scheme.(d) Relative error in forcing spectrum at
times $t=0.0,0.5,1.0$, and also $t=2,4,6\ldots,80$ using direct replacement
scheme.

##  VII Conclusions

In conclusion, we have developed a novel algorithm for recovering an unknown
large scale forcing function in the 2D Navier-Stokes equations that can be
implemented on the fly. The algorithm is both rigorously and numerically
justified, simple to implement, and does not require an ensemble of
simulations to generate the desired forcing function. In particular, the
algorithm here will be applicable to a host of other situations where
dissipative mechanisms and driving forces are simultaneously present; this is
typically the case in many hydrodynamic settings, especially ones arising in
the context of climate and weather modeling. We emphasize that the rigorous
result provided here is restricted to two dimensions only due to the current
lack of a global-in-time bound on solutions for the gradient of the velocity
field for the three dimensional Navier-Stokes equations. Under suitable
assumptions on the regularity of solutions to the 3D equations, we are
confident that the same algorithm presented here will recover the forcing in
that setting. Hence, we do not anticipate that the recovery of the forcing
critically relies on the inverse cascade in 2D turbulence which tends to
aggregate large scale structures [20]. Further investigations (both numerical
and analytical) are required to confirm these conjectures however.

###### Acknowledgements.

We wish to acknowledge the ADAPT group and J. Murri and B. Pachev in
particular who invigorated discussions surrounding these results developed
here. A.F. was partially supported by NSF grant DMS-2206493. A.L. was
supported in part by NSF grants DMS-2206762 and CMMI-1953346, and USGS grant
G23AS00157 number GRANT13798170. V.M. was partially supported by NSF grants
DMS-2206491 and DMS-2213363. J.P.W. was partially supported by NSF grant
DMS-2206762, as well as the Simons Foundation travel grant under 586788.

##  VIII Appendix

We provide some additional supporting details for the convergence analysis
performed after (13). In particular we utilize several inequalities in the
rigorous analysis which are described below. First we note that we will make
use of the $L^{p}$ norm on functions defined for functions
$\phi:\Omega\rightarrow\mathbb{R}$ on the domain $\Omega=[0,L]^{2}$ as

$\|\phi\|_{p}=\left(\int_{\Omega}|\phi(\mathbf{x})|^{p}d\mathbf{x}\right)^{1/p},$

with $p=\infty$ corresponding to the supremum norm over the entire domain.

This allows us to describe the following inequalities defined for all
functions $\phi,\psi,\chi:\Omega\rightarrow\mathbb{R}$ (extension to vector-
valued functions is immediate):

  * •

Cauchy-Schwarz inequality:

$\int_{\Omega}|\phi(\mathbf{x})\psi(\mathbf{x})|d\mathbf{x}\leq\|\phi\|_{2}\|\psi\|_{2}.$
(19)

  * •

Agmon interpolation inequality:

$\|\phi\|_{\infty}\leq c\|\nabla\phi\|_{2}^{1/2}\|\nabla^{2}\phi\|_{2}^{1/2},$
(20)

where $c$ is a universal constant.

  * •

Ladyzhenskaya interpolation inequality:

$\|\phi\|_{4}\leq c\|\nabla\phi\|_{2}^{1/2}\|\phi\|_{2}^{1/2},$ (21)

where once again $c$ is a universal constant.

  * •

Generalized Hölder inequality:

$\|\phi\psi\chi\|_{r}\leq\|\phi\|_{p}\|\psi\|_{q}\|\chi\|_{s},$ (22)

where $\frac{1}{s}+\frac{1}{q}+\frac{1}{p}=\frac{1}{r}$.

  * •

Grönwall inequality:

Suppose that $f(t),g(t):[0,T]\rightarrow\mathbb{R}$ are continuously
differentiable, and satisfy

$\frac{df}{dt}\leq\alpha f(t)+g(t),$ (23)

for some constant $\alpha$. Then

$f(t)\leq f(0)e^{\alpha t}+\int_{0}^{t}e^{\alpha(t-s)}g(s)ds.$ (24)

The rigorous analysis also makes use of a priori bounds on the velocity field
corresponding to (1). To state these bounds, suppose
$\mathbf{f},\nabla\mathbf{f}\in L^{2}(\Omega)$ such that $\mathbf{f}$ is
divergence-free. Let $\kappa_{0}$ and $G$ be defined by (2) and denote a shape
factor of the force by

$\sigma_{1}:=\frac{\kappa_{0}\|\nabla\mathbf{f}\|_{2}}{\|\mathbf{f}\|_{2}}.$

Strong solutions of (1) satisfy

$\displaystyle\begin{split}\|\mathbf{u}(t)\|_{2}^{2}&\leq\|\mathbf{u}_{0}\|_{2}^{2}e^{-\kappa_{0}^{2}\nu
t}+\kappa_{0}^{2}\nu^{2}G^{2}<R_{0}\\\
\|\nabla\mathbf{u}(t)\|_{2}^{2}&\leq\|\nabla\mathbf{u}_{0}\|_{2}^{2}e^{-\kappa_{0}^{2}\nu
t}+\nu^{2}G^{2}<R_{1}\\\
\|\Delta\mathbf{u}(t)\|_{2}^{2}&\leq\|\Delta\mathbf{u}_{0}\|_{2}^{2}e^{-\kappa_{0}^{2}\nu
t}+c\kappa_{0}\nu^{2}(\sigma_{1}+G)^{2}G^{2}<R_{2},\end{split}$ (25)

for all $t\geq 0$. The first two bounds are classical and may be found in [33,
102, 51], whereas the third can be found in, for instance, [85]; we have taken
the $R_{j}$ to be strict upper bounds on these estimates. Note that when $t$
is sufficiently large, one may assume that $R_{0},R_{1},R_{2}$ are independent
of the initial velocity; this is a safe assumption in the analysis performed
above since we must always wait long enough before the first update to ensure
that $\mathbf{u}(t),\Delta\mathbf{u}(t)$ henceforth remain bounded by $R_{0}$
and $R_{2}$, respectively.

Lastly, the crucial ingredient used to close all of the estimates and arrive
at (14) is provided by the state error estimates (12). The bound stated for
$\nabla\mathbf{w}^{k}$ was originally proven in [85]. We state its precise
form here: Let $\mathbf{u}_{0},\mathbf{v}_{0}$ be divergence-free velocity
fields, which are mean-free over $\Omega$, and whose derivatives up to second-
order are square-integrable over $\Omega$. Then there exist universal
constants $\overline{c}_{0},\underline{c}_{0}\geq 1$ such that if $\mu,N$
satisfy

$\displaystyle\underline{c}_{0}\nu\left(\sigma_{1}+G\right)G^{2}\leq{\mu
h^{2}}\leq\overline{c}_{0},$

then for each $k\geq 1$, there exists $t_{k}>t_{k-1}>0$ (where $t_{0}:=0$),
such that

$\displaystyle\|\nabla\mathbf{w}^{k}(t)\|_{2}$
$\displaystyle\leq\left(\frac{2C_{1}\kappa_{0}^{2}\nu}{\mu}\right)^{1/2}\frac{\|\Delta\mathbf{g}^{k}(t_{k})\|_{2}}{\kappa_{0}\nu},$

for all $t\in[t_{k},\infty)$, for some universal constant $C_{1}\geq 1$
independent of $k$.

On the other hand, the bound for $\mathbf{w}^{k}$ was not required in [85].
For the sake of completeness, we provide the details of this estimate here.
Observe that the system governing the evolution of
$\mathbf{w}^{k}=\mathbf{v}^{k}-\mathbf{u}$ is given by

$\displaystyle\partial_{t}\mathbf{w}^{k}-$
$\displaystyle\nu\nabla^{2}\mathbf{w}^{k}+\mathbf{w}^{k}\cdotp\nabla\mathbf{w}^{k}+\mathbf{w}^{k}\cdotp\nabla\mathbf{u}+\mathbf{u}\cdotp\nabla\mathbf{w}^{k}=-\nabla
r^{k}+\Delta\mathbf{g}^{k}-\mu
I_{h}\mathbf{w}^{k},\quad\nabla\cdotp\mathbf{w}^{k}=0$

where $r^{k}=q^{k}-p$, where $q^{k},p$ denote the pressures corresponding to
$\mathbf{v}^{k},\mathbf{u}$, respectively. Then the corresponding energy
balance is given by

$\displaystyle\frac{1}{2}$
$\displaystyle\frac{d}{dt}\|\mathbf{w}^{k}\|_{2}^{2}+\nu\|\nabla\mathbf{w}^{k}\|_{2}^{2}+\mu\|\mathbf{w}^{k}\|_{2}^{2}=-\langle\mathbf{w}^{k}\cdotp\nabla\mathbf{u},\mathbf{w}^{k}\rangle+\langle\Delta\mathbf{g}^{k},\mathbf{w}^{k}\rangle+\mu\|J_{h}\mathbf{w}^{k}\|_{2}^{2}.$

The first of these terms on the right-hand side can be estimated with Hölder’s
inequality, Ladyzhenskaya’s inequality, the Cauchy-Schwarz inequality, and
(25) to obtain

$\displaystyle|\langle\mathbf{w}^{k}\cdotp\nabla\mathbf{u},\mathbf{w}^{k}\rangle|\leq\|\nabla\mathbf{u}\|_{2}\|\mathbf{w}^{k}\|_{4}^{2}\leq
cR_{1}\|\nabla\mathbf{w}^{k}\|_{2}\|\mathbf{w}^{k}\|_{2}\leq
c\frac{R_{1}^{2}}{\mu}\|\nabla\mathbf{w}^{k}\|_{2}^{2}+\frac{\mu}{4}\|\mathbf{w}^{k}\|_{2}^{2}.$

The second term can be estimated with the Cauchy-Schwarz inequality to obtain

$\displaystyle|\langle\Delta\mathbf{g}^{k},\mathbf{w}^{k}\rangle|\leq\|\Delta\mathbf{g}^{k}\|_{2}\|\mathbf{w}^{k}\|_{2}\leq\frac{1}{\mu}\|\Delta\mathbf{g}^{k}\|_{2}^{2}+\frac{\mu}{4}\|\mathbf{w}^{k}\|_{2}^{2}.$
(26)

Lastly, the third term can be estimated using the Cauchy-Schwarz inequality
and (3) to obtain

$\displaystyle\mu\|J_{h}\mathbf{w}^{k}\|_{2}^{2}\leq C_{1}^{2}\mu
h^{2}\|\nabla\mathbf{w}^{k}\|_{2}^{2}.$

Now let us assume that $\mu,h$ satisfy

$2C_{1}^{2}\mu h^{2}\leq\nu,\quad c\mu\nu\geq R_{1}^{2}.$

We may then combine the above estimates to arrive at

$\displaystyle\frac{d}{dt}\|\mathbf{w}^{k}\|_{2}^{2}+\mu\|\mathbf{w}^{k}\|_{2}^{2}\leq\frac{1}{\mu}\|\Delta\mathbf{g}^{k}\|_{2}^{2},$

from which (12) follows upon taking $t$ sufficiently large to evaluate
$\mathbf{g}^{k}$, and finally an application of Grönwall’s inequality.

## References

  * [1] H. D. Abarbanel, S. Shirman, D. Breen, N. Kadakia, D. Rey, E. Armstrong, and D. Margoliash.  A unifying view of synchronization for data assimilation in complex nonlinear networks.  Chaos, 27(12):126802, 2017.
  * [2] L. Agasthya, P. Clark Di Leoni, and L. Biferale.  Reconstructing rayleigh–bénard flows out of temperature-only measurements using nudging.  Physics of Fluids, 34(1), 01 2022.  015128\.
  * [3] D. A. Albanez and M. J. Benvenutti.  Continuous data assimilation algorithm for simplified bardina model.  Evolution Equations & Control Theory, 7(1):33, 2018.
  * [4] D. A. Albanez, H. J. Nussenzveig Lopes, and E. S. Titi.  Continuous data assimilation for the three-dimensional Navier-Stokes-$\alpha$ model.  Asymptotic Analysis, 97(1-2):139–164, 2016.
  * [5] M. Altaf, E. Titi, T. Gebrael, O. Knio, L. Zhao, M. McCabe, and I. Hoteit.  Downscaling the 2D bénard convection equations using continuous data assimilation.  Computational Geosciences, 21(3):393–410, 2017.
  * [6] A. Azouani, E. Olson, and E. Titi.  Continuous data assimilation using general interpolant observables.  J. Nonlinear Sci., 24(2):277–304, 2014.
  * [7] A. Azouani, E. Olson, and E. S. Titi.  Continuous data assimilation using general interpolant observables.  Journal of Nonlinear Science, 24(2):277–304, 2014.
  * [8] S.-J. Baek.  Synchronization in chaotic systems: Coupling of chaotic maps, data assimilation and weather forecasting.  University of Maryland, College Park, 2007.
  * [9] C. Basdevant.  Technical improvements for direct numerical simulation of homogeneous three-dimensional turbulence.  Journal of Computational Physics, 50(2):209–214, 1983.
  * [10] C. Basdevant.  Digital partial image velocimetry.  Experiments in fluids, 10(2):181–193, 1991.
  * [11] H. Bessaih, E. Olson, and E. Titi.  Continuous data assimilation with stochastically noisy data.  J. Nonlinear Sci., 28(0):729–753, 2015.
  * [12] A. Biswas, K. Brown, and V. Martinez.  Higher-order synchronization and a refined paradigm for global interpolant observables.  Ann. Appl. Math., 38(3):pp. 1–60, 2022.
  * [13] A. Biswas, C. Foias, C. Mondaini, and E. Titi.  Downscaling data assimilation algorithm with applications to statistical solutions of the Navier-Stokes equations.  Annales de l’Institut Henri Poincaré C, Analyse non linéaire, 36(2):295–326, 2019.
  * [14] A. Biswas, C. Foias, C. F. Mondaini, and E. S. Titi.  Downscaling data assimilation algorithm with applications to statistical solutions of the Navier–Stokes equations.  In Annales de l’Institut Henri Poincaré C, Analyse non linéaire, volume 36, pages 295–326. Elsevier, 2019.
  * [15] A. Biswas, J. Hudson, A. Larios, and Y. Pei.  Continuous data assimilation for the 2D magnetohydrodynamic equations using one component of the velocity and magnetic fields.  Asymptot. Anal., 108(1-2):1–43, 2018.
  * [16] A. Biswas and V. R. Martinez.  Higher-order synchronization for a data assimilation algorithm for the 2D Navier–Stokes equations.  Nonlinear Analysis: Real World Applications, 35:132–157, 2017.
  * [17] A. Biswas and R. Price.  Continuous data assimilation for the three-dimensional Navier–Stokes equations.  SIAM Journal on Mathematical Analysis, 53(6):6697–6723, 2021.
  * [18] J. Blocher, V. R. Martinez, and E. Olson.  Data assimilation using noisy time-averaged measurements.  Physica D: Nonlinear Phenomena, 376:49–59, 2018.
  * [19] D. Blömker, K. Law, A. M. Stuart, and K. C. Zygalakis.  Accuracy and stability of the continuous-time 3DVAR filter for the Navier–Stokes equation.  Nonlinearity, 26(8):2193, 2013.
  * [20] G. Boffetta, R. E. Ecke, et al.  Two-dimensional turbulence.  Annual review of fluid mechanics, 44(1):427–451, 2012.
  * [21] C. E. Brett, K. F. Lam, K. Law, D. McCormick, M. R. Scott, and A. Stuart.  Accuracy and stability of filters for dissipative PDEs.  Physica D: Nonlinear Phenomena, 245(1):34–45, 2013.
  * [22] S. L. Brunton, J. L. Proctor, and J. N. Kutz.  Discovering governing equations from data by sparse identification of nonlinear dynamical systems.  Proceedings of the national academy of sciences, 113(15):3932–3937, 2016.
  * [23] E. Carlson, J. Hudson, and A. Larios.  Parameter recovery for the 2 dimensional Navier–Stokes equations via continuous data assimilation.  SIAM J. Sci. Comput., 42(1):A250–A270, 2020.
  * [24] E. Carlson, J. Hudson, A. Larios, V. R. Martinez, E. Ng, and J. Whitehead.  Dynamically learning the parameters of a chaotic system using partial observations.  Discrete Contin Dyn Syst Ser A, 42(8):3809–3839, 2022.
  * [25] E. Carlson, J. Hudson, A. Larios, V. R. Martinez, E. Ng, and J. P. Whitehead.  Dynamically learning the parameters of a chaotic system using partial observations.  Discrete and Continuous Dynamical Systems, 42(8):3809–3839, 2022\.
  * [26] E. Carlson, A. Larios, and E. S. Titi.  Super-exponential convergence rate of a nonlinear continuous data assimilation algorithm: The 2D Navier–Stokes equations paradigm.  2023\.  (submitted) arXiv:2304.01128.
  * [27] E. Carlson, L. Van Roekel, M. Petersen, H. C. Godinez, and A. Larios.  CDA algorithm implemented in MPAS-O to improve eddy effects in a mesoscale simulation.  2023\.  (submitted).
  * [28] A. Carrassi, M. Bocquet, L. Bertino, and G. Evensen.  Data assimilation in the geosciences: An overview of methods, issues, and perspectives.  Wiley Interdisciplinary Reviews: Climate Change, 9(5):e535, 2018\.
  * [29] E. Celik and E. Olson.  Data assimilation using time-delay nudging in the presence of gaussian noise.  arXiv:2108.03631, 2022.  (preprint).
  * [30] E. Celik, E. Olson, and E. S. Titi.  Spectral filtering of interpolant observables for a discrete-in-time downscaling data assimilation algorithm.  SIAM J. Appl. Dyn. Syst., 18(2):1118–1142, 2019.
  * [31] B. Cockburn, D. Jones, and E. Titi.  Determining degrees of freedom for nonlinear dissipative equations.  C.R. Acad. Sci. Paris Sér. I Math, 321:563–568, 1995.
  * [32] P. Constantin and C. Foias.  Navier-Stokes equations.  Chicago Lectures in Mathematics. University of Chicago Press, Chicago, IL, 1988.
  * [33] P. Constantin and C. Foias.  Navier-Stokes equations.  Chicago Lectures in Mathematics. University of Chicago Press, Chicago, IL, 1988.
  * [34] R. Dascaliuc, C. Foias, and M. S. Jolly.  Relations between energy and enstrophy on the global attractor of the 2-D Navier–Stokes equations.  J. Dynam. Differential Equations, 17(4):643–736, 2005.
  * [35] R. Dascaliuc, C. Foias, and M. S. Jolly.  On the asymptotic behavior of average energy and enstrophy in 3D turbulent flows.  Phys. D, 238(7):725–736, 2009.
  * [36] P. A. Davidson.  Turbulence: An Introduction for Scientists and Engineers.  Oxford University Press, 2015.
  * [37] S. Desamsetti, H. P. Dasari, S. Langodan, E. S. Titi, O. Knio, and I. Hoteit.  Efficient dynamical downscaling of general circulation models using continuous data assimilation.  Quarterly Journal of the Royal Meteorological Society, 145(724):3175–3194, 2019.
  * [38] P. C. Di Leoni, A. Mazzino, and L. Biferale.  Inferring flow parameters and turbulent configuration with physics-informed data assimilation and spectral nudging.  Physical Review Fluids, 3(10):104604, 2018.
  * [39] P. C. Di Leoni, A. Mazzino, and L. Biferale.  Synchronization to big data: Nudging the Navier-Stokes equations for data assimilation of turbulent flows.  Physical Review X, 10(1):011023, 2020.
  * [40] F. Djeumou and U. Topcu.  Learning to reach, swim, walk and fly in one trial: Data-driven control with scarce data and side information.  In Learning for Dynamics and Control Conference, pages 453–466. PMLR, 2022.
  * [41] V. Eswaran and S. B. Pope.  An examination of forcing in direct numerical simulations of turbulence.  Computers & Fluids, 16(3):257–278, 1988.
  * [42] G. Evensen, F. C. Vossepoel, and P. J. van Leeuwen.  Data assimilation fundamentals: A unified formulation of the state and parameter estimation problem.  Springer Nature, 2022.
  * [43] A. Farhat, L. E., and E. Titi.  A data assimilation algorithm: the paradigm of the 3D Leray-$\alpha$ model of turbulence.  Nonlinear Partial Differential Equations Arising from Geometry and Physics. London Mathematical Society. Lecture Notes Series. Cambridge University Press, 2017.
  * [44] A. Farhat, N. E. Glatt-Holtz, V. R. Martinez, S. A. McQuarrie, and J. P. Whitehead.  Data assimilation in large Prandtl Rayleigh-Bénard convection from thermal measurements.  SIAM J. Appl. Dyn. Syst., 19(1):510–540, 2020.
  * [45] A. Farhat, H. Johnston, M. Jolly, and E. Titi.  Assimilation of nearly turbulent Rayleigh-Bénard flow through vorticity or local circulation measurements: A computational study.  J. Sci. Comput., pages 1–15, 2018.
  * [46] A. Farhat, M. S. Jolly, and E. S. Titi.  Continuous data assimilation for the 2D bénard convection through velocity measurements alone.  Physica D: Nonlinear Phenomena, 303:59–66, 2015.
  * [47] A. Farhat, E. Lunasin, and E. Titi.  Abridged continuous data assimilation for the 2D Navier-Stokes equations utilizing measurements of only one component of the velocity field.  J. Math. Fluid Mech., 18(1):1–23, 2016.
  * [48] A. Farhat, E. Lunasin, and E. Titi.  Data assimilation algorithm for 3D Bénard convection in porous media employing only temperature measurements.  J. Math. Anal. Appl., 438(1):492–506, 2016.
  * [49] A. Farhat, E. Lunasin, and E. Titi.  On the Charney conjecture of data assimilation employing temperature measurements alone: The paradigm of 3D planetary geostrophic model.  Math. of Climate and Wea. Forecasting, 2:61–74, 2016.
  * [50] A. Farhat, E. Lunasin, and E. Titi.  Continuous data assimilation for a 2D Bénard convection system through horizontal velocity measurements alone.  J. Nonlinear Sci., 27:1065–1087, 2017.
  * [51] C. Foias, O. Manley, R. Rosa, and R. Temam.  Navier-Stokes equations and turbulence, volume 83 of Encyclopedia of Mathematics and its Applications.  Cambridge University Press, Cambridge, 2001.
  * [52] C. Foias, O. Manley, R. Temam, and Y. Treve.  Number of modes governing two-dimensional viscous, incompressible flows.  Physical Review Letters, 50(14):1031, 1983.
  * [53] C. Foias, C. Mondaini, and E. Titi.  A discrete data assimilation scheme for the solutions of the 2D Navier-Stokes equations and their statistics.  SIAM J. Appl. Dyn. Syst., 15(4):2019–2142, 2016.
  * [54] C. Foias, C. F. Mondaini, and E. S. Titi.  A discrete data assimilation scheme for the solutions of the two-dimensional Navier–Stokes equations and their statistics.  SIAM Journal on Applied Dynamical Systems, 15(4):2109–2142, 2016\.
  * [55] C. Foias and G. Prodi.  Sur le comportement global des solutions non-stationnaires des équations de Navier-Stokes en dimension $2$.  Rend. Sem. Mat. Univ. Padova, 39:1–34, 1967.
  * [56] T. Franz, A. Larios, and C. Victor.  The bleeps, the sweeps, and the creeps: Convergence rates for dynamic observer patterns via data assimilation for the 2D Navier-Stokes equations.  Comput. Methods Appl. Mech. Engrg., 392:Paper No. 114673, 19, 2022\.
  * [57] M. Gardner, A. Larios, L. G. Rebholz, D. Vargun, and C. Zerfas.  Continuous data assimilation applied to a velocity-vorticity formulation of the 2D Navier–Stokes equations.  Electron. Res. Arch., 29(3):2223–2247, 2021.
  * [58] M. Gesho, E. Olson, and E. S. Titi.  A computational study of a data assimilation algorithm for the two-dimensional Navier-Stokes equations.  Communications in Computational Physics, 19(4):1094–1110, 2016\.
  * [59] C. H. Gibson.  Turbulence in the ocean, atmosphere, galaxy, and universe.  Journal of Cosmology, 17(16):7485–7523, 2011.
  * [60] M. D. Gunzburger.  Perspectives in flow control and optimization.  SIAM, 2002.
  * [61] M. D. Gunzburger.  Flow control, volume 68.  Springer Science & Business Media, 2012.
  * [62] J. E. Hoke and R. A. Anthes.  The initialization of numerical models by a dynamic-initialization technique.  Monthly Weather Review, 104(12):1551–1556, 1976.
  * [63] H. Ibdah, C. Mondaini, and E. Titi.  Fully discrete numerical schemes of a data assimilation algorithm: uniform-in-time error estimates.  IMA Journal of Numerical Analysis, 40(4):2584–2625, 11 2019.
  * [64] H. A. Ibdah, C. F. Mondaini, and E. S. Titi.  Fully discrete numerical schemes of a data assimilation algorithm: uniform-in-time error estimates.  IMA Journal of Numerical Analysis, 40(4):2584–2625, 2020.
  * [65] M. Jolly, V. Martinez, E. Olson, and E. Titi.  Continuous data assimilation with blurred-in-time measurements of the surface quasi-geostrophic equation.  Chin. Ann. Math. Ser. B, 40:721––764, 2019.
  * [66] M. Jolly, V. Martinez, and E. Titi.  A data assimilation algorithm for the 2D subcritical surface quasi-geostrophic equation.  Adv. Nonlinear Stud., 35:167–192, 2017.
  * [67] M. S. Jolly and A. Pakzad.  Data assimilation with higher order finite element interpolants.  arXiv preprint arXiv:2108.03631, 2021.
  * [68] M. S. Jolly, T. Sadigov, and E. S. Titi.  A determining form for the damped driven nonlinear schrödinger equation—fourier modes case.  Journal of Differential Equations, 258(8):2711–2744, 2015.
  * [69] M. S. Jolly, T. Sadigov, and E. S. Titi.  Determining form and data assimilation algorithm for weakly damped and driven Korteweg–deVries equation—Fourier modes case.  Nonlinear Analysis: Real World Applications, 36:287–317, 2017.
  * [70] D. Jones and E. Titi.  Determining finite volume elements for the 2D Navier-Stokes equations.  Phys. D, 60:165–174, 1992.
  * [71] D. Jones and E. Titi.  On the number of determining nodes for the 2D Navier-Stokes equations.  J. Math. Anal., 168:72–88, 1992.
  * [72] D. Jones and E. Titi.  Upper bounds on the number of determining modes, nodes and volume elements for the Navier-Stokes equations.  Indiana Math. J., 42:875–887, 1993.
  * [73] O. Knio, S. Desamsetti, H. Dasari, S. Langodan, I. Hoteit, and E. Titi.  Efficient dynamical downscaling of general circulation models using continuous data assimilation.  In APS Division of Fluid Dynamics Meeting Abstracts, pages C20–009, 2019.
  * [74] A. N. Kolmogorov.  The local structure of turbulence in incompressible viscous fluid for very large Reynolds numbers.  Cr Acad. Sci. URSS, 30:301–305, 1941.
  * [75] A. Larios and Y. Pei.  Approximate continuous data assimilation of the 2D Navier–Stokes equations via the Voigt-regularization with observable data.  Evol. Equ. Control Theory, 9(3):733–751, 2020.
  * [76] A. Larios and Y. Pei.  Nonlinear continuous data assimilation.  Evol. Equ. Control Theory, 2023.  (accepted, in press).
  * [77] A. Larios, Y. Pei, and L. Rebholz.  Global well-posedness of the velocity-vorticity-Voigt model of the 3D Navier–Stokes equations.  J. Differential Equations, 266(5):2435–2465, 2019.
  * [78] A. Larios, M. R. Petersen, and C. Victor.  Application of continuous data assimilation in high-resolution ocean modeling.  (submitted) arXiv:2308.02705, 2023.
  * [79] A. Larios, L. G. Rebholz, and C. Zerfas.  Global in time stability and accuracy of IMEX-FEM data assimilation schemes for Navier–Stokes equations.  Comput. Methods Appl. Mech. Engrg., 345:1077–1093, 2019.
  * [80] A. Larios and C. Victor.  Continuous data assimilation with a moving cluster of data points for a reaction diffusion equation: A computational study.  Commun. Comput. Phys., 29(4):1273–1298, 2021.
  * [81] A. Larios and C. Victor.  Continuous data assimilation for the 3D and higher-dimensional Navier–Stokes equations with higher-order fractional diffusion.  (submitted) arXiv:2307.00096, 2023.
  * [82] A. Larios and C. Victor.  The second-best way to do sparse-in-time continuous data assimilation: Improving convergence rates for the 2D and 3D Navier–Stokes equations.  (submitted) arXiv:2303.03495, 2023.
  * [83] C. Marchioro.  An example of absence of turbulence for any Reynolds number.  Comm. Math. Phys., 105(1):99–106, 1986.
  * [84] P. A. Markowich, E. S. Titi, and S. Trabelsi.  Continuous data assimilation for the three-dimensional Brinkman–Forchheimer-extended Darcy model.  Nonlinearity, 29(4):1292, 2016.
  * [85] V. Martinez.  On the reconstruction of unknown driving forces from low-mode observations in the 2D Navier-Stokes equations.  arXiv:2208.00541v1, page 15, 31 Jul 2022.
  * [86] V. R. Martinez.  Convergence analysis of a viscosity parameter recovery algorithm for the 2D Navier–Stokes equations.  Nonlinearity, 35(5):2241, 2022.
  * [87] F. R. Menter.  Two-equation eddy-viscosity turbulence models for engineering applications.  AIAA journal, 32(8):1598–1605, 1994.
  * [88] R. Mojgani, A. Chattopadhyay, and P. Hassanzadeh.  Discovery of interpretable structural model errors by combining bayesian sparse regression and data assimilation: A chaotic kuramoto–sivashinsky test case.  Chaos: An Interdisciplinary Journal of Nonlinear Science, 32(6):061105, 2022.
  * [89] C. Mondaini and E. Titi.  Uniform-in-time error estimates for the postprocessing Galerkin method applied to a data assimilation algorithm.  SIAM J. Numer. Anal., 56(1):78–110, 2018.
  * [90] E. Olson and E. Titi.  Determining modes for continuous data assimilation in 2D turbulence.  J. Stat. Phys., 113(5-6):799–840, 2003.
  * [91] E. Olson and E. S. Titi.  Determining modes and Grashof number in 2D turbulence: a numerical case study.  Theor. Comp. Fluid Dyn., 22(5):327–339, 2008.
  * [92] B. Pachev, J. P. Whitehead, and S. A. McQuarrie.  Concurrent multi-parameter learning demonstrated on the kuramoto-sivashinsky equation.  SIAM Scientific Computing, 2022.
  * [93] S. B. Pope.  Turbulent Flows.  Cambridge University Press, 2000.
  * [94] J. Poterjoy.  A localized particle filter for high-dimensional nonlinear systems.  Monthly Weather Review, 144(1):59–76, 2016.
  * [95] J. Robinson.  Infinite-Dimensional Dynamical Systems: An Introduction to Dissipative Parabolic PDEs and the Theory of Global Attractors. Cambridge texts in applied mathematics.  Appl. Mech. Rev., 56(4):B54–B55, 2003.
  * [96] C. Schneide, A. Pandey, K. Padberg-Gehle, and J. Schumacher.  Probing turbulent superstructures in Rayleigh-Bénard convection by Lagrangian trajectory clusters.  Physical Review Fluids, 3(11):113501, 2018.
  * [97] J. Schumacher, J. D. Scheel, D. Krasnov, D. A. Donzis, V. Yakhot, and K. R. Sreenivasan.  Small-scale universality in fluid turbulence.  Proceedings of the National Academy of Sciences, 111(30):10961–10965, 2014.
  * [98] K. M. Smith, C. Caulfield, and J. Taylor.  Turbulence in forced stratified shear flows.  Journal of Fluid Mechanics, 910, 2021.
  * [99] R. J. Stevens, A. Blass, X. Zhu, R. Verzicco, and D. Lohse.  Turbulent thermal superstructures in Rayleigh-Bénard convection.  Physical Review Fluids, 3(4):041501, 2018.
  * [100] I. G. Szendro, M. A. Rodríguez, and J. M. López.  On the problem of data assimilation by means of synchronization.  Journal of Geophysical Research: Atmospheres, 114(D20), 2009.
  * [101] P. Tandeo, P. Ailliot, M. Bocquet, A. Carrassi, T. Miyoshi, M. Pulido, and Y. Zhen.  A review of innovation-based methods to jointly estimate model and observation error covariance matrices in ensemble data assimilation.  Monthly Weather Review, 148(10):3973–3994, 2020.
  * [102] R. Temam.  Infinite-Dimensional Dynamical Systems in Mechanics and Physics, volume 68 of Applied Mathematical Sciences.  Springer-Verlag, New York, 1997.
  * [103] R. Temam.  Infinite-Dimensional Dynamical Systems in Mechanics and Physics, volume 68 of Applied Mathematical Sciences.  Springer-Verlag, New York, 1997.
  * [104] D. Ting.  Basics of Engineering Turbulence.  Academic Press, 2016.
  * [105] S.-C. Yang, D. Baker, H. Li, K. Cordes, M. Huff, G. Nagpal, E. Okereke, J. Villafane, E. Kalnay, and G. S. Duane.  Data assimilation as synchronization of truth and model: Experiments with the three-variable lorenz system.  Journal of the atmospheric sciences, 63(9):2340–2354, 2006.
  * [106] S.-C. Yang, M. Corazza, A. Carrassi, E. Kalnay, and T. Miyoshi.  Comparison of local ensemble transform kalman filter, 3DVAR, and 4DVAR in a quasigeostrophic model.  Monthly weather review, 137(2):693–709, 2009.
  * [107] P. Yeung and K. Ravikumar.  Advancing understanding of turbulence through extreme-scale computation: Intermittency and simulations at large problem sizes.  Physical Review Fluids, 5(11):110517, 2020.