Computability of the Solutions to Navier-Stokes Equations via Recursive Approximation

TL;DR

This paper investigates whether solutions to the 2D Navier-Stokes equations can be computed recursively, establishing a framework that makes the mild and strong solutions uniformly computable using a novel representation of divergence-free vector fields.

Contribution

It constructs a natural encoding for divergence-free vector fields on a 2D domain and proves the computability of both mild and strong solutions to the Navier-Stokes equations within this framework.

Findings

The representation makes the mild solution to the Stokes problem computable.

The strong local solution to the Navier-Stokes initial value problem is also computable.

The paper develops intricate estimates to establish these computability results.

Abstract

As one of the seven open problems in the addendum to their 1989 book "Computability in Analysis and Physics", Pour-El and Richards proposed ``... the recursion theoretic study of particular nonlinear problems of classical importance. Examples are the Navier-Stokes equation, the KdV equation, and the complex of problems associated with Feigenbaum's constant.'' In this paper, we approach the question of whether the Navier-Stokes Equation admits recursive solutions in the sense of Weihrauch's Type-2 Theory of Effectivity. A natural encoding (``representation'') is constructed for the space of divergence-free vector fields on 2-dimensional open square $\Omega = (-1, 1)^2$. This representation is shown to render first the mild solution to the Stokes Dirichlet problem and then a strong local solution to the nonlinear inhomogeneous incompressible Navier-Stokes initial value problem uniformly computable. Based on classical approaches, the proofs make use of many subtle and intricate estimates which are developed in the paper for establishing the computability results.

Equations

$$\partial_{t}\vec{u}-\bigtriangleup\vec{u}+(\vec{u}\cdot\nabla)\vec{u}+\nabla P=\vec{f},\quad\nabla\cdot\vec{u}=0,\quad\vec{u}(0)=\vec{a},\quad\vec{u}\big{|}_{\partial\Omega}\equiv\vec{0}$$

$${\partial }_t\vec{u}+\mathbb{A}\vec{u}+\mathbb{B}\vec{u}=\vec{g}(t>0),\vec{u}(0)=\vec{a}\in L_{2,0}^{\sigma }(\Omega )$$

$${\partial }_t\vec{v}+\mathbb{A}\vec{v}=\vec{0}\text{for }t\ge 0,\vec{v}(0)=\vec{a}\in L_{2,0}^{\sigma }(\Omega )$$

$$\vec{u}(t)=e^{-t\mathbb{A}}\vec{a}+{\int }_0^t{e}^{-(t-s)\mathbb{A}}\vec{g}(s)ds-{\int }_0^t{e}^{-(t-s)\mathbb{A}}\mathbb{B}\vec{u}(s)ds\text{for }t\ge 0$$

$${\vec{v}}_0(t)=e^{-t\mathbb{A}}\vec{a}+{\int }_0^t{e}^{-(t-s)\mathbb{A}}\vec{g}(s)ds,{\vec{v}}_{n+1}(t)={\vec{v}}_0(t)-{\int }_0^t{e}^{-(t-s)\mathbb{A}}\mathbb{B}{\vec{v}}_n(s)ds$$

$$\nabla P=\vec{f}-{\partial }_t\vec{u}+△\vec{u}-(\vec{u}\cdot \nabla )\vec{u}$$

$$\displaystyle\mathbb{S}\;:\;C\big{(}[0;\infty),L^{\sigma}_{2,0}(\Omega)\big{)}\times C\big{(}[0;\infty),$$

$$\displaystyle\;\mapsto\;\vec{v}_{n+1}\;\in\;C\big{(}[0;\infty),L^{\sigma}_{2,0}(\Omega)\big{)}$$

$$T:L^{\sigma}_{2,0}(\Omega)\times C\big{(}[0;\infty),L^{\sigma}_{2,0}(\Omega)\big{)}\to(0;\infty),\quad(\vec{a},\vec{f})\mapsto T(\vec{a},\vec{f})$$

$$\mathcal{S}:\subseteq L^{\sigma}_{2,0}(\Omega)\times C\big{(}[0;\infty),L^{\sigma}_{2,0}(\Omega)\big{)}\times[0;\infty)\to L^{\sigma}_{2,0}(\Omega)\!\times\!L_{2}(\Omega)$$

$$\vec{p}=(p_{1},p_{2})=\big{(}\sum_{i,j=0}^{N}a^{1}_{i,j}x^{i}y^{j},\sum_{i,j=0}^{N}a^{2}_{i,j}x^{i}y^{j}\big{)}$$

$$(i+1)a_{i+1,j}^1+(j+1)a_{i,j+1}^2=0,$$

$$(i+1)a_{i+1,N}^1=0,$$

$$(j+1)a_{N,j+1}^2=0,$$

$${\sum }_{i=0}^N{a}_{i,j}^1={\sum }_{i=0}^N{a}_{i,j}^2={\sum }_{i=0}^N(-1{)}^i{a}_{i,j}^1={\sum }_{i=0}^N(-1{)}^i{a}_{i,j}^2=0$$

$${\sum }_{j=0}^N{a}_{i,j}^1={\sum }_{j=0}^N{a}_{i,j}^2={\sum }_{j=0}^N(-1{)}^j{a}_{i,j}^1={\sum }_{j=0}^N(-1{)}^j{a}_{i,j}^2=0$$

$$\mathbb{T}_{k}p_{j}(x,y)\;=\;\left\{\begin{array}[]{ll}p_{j}(\frac{x}{1-2^{-k}},\frac{y}{1-2^{-k}}),&\;-1+2^{-k}\leq x,y\leq 1-2^{-k}\\ 0,&\mbox{otherwise}\end{array}\right.$$

$$\frac{\partial {\mathbb{T}}_k{p}_1}{\partial x}(x,y)+\frac{\partial {\mathbb{T}}_k{p}_2}{\partial y}(x,y)$$

$$=\frac{1}{1-2^{-k}}\left[\frac{\partial p_1}{\partial x^{\prime }}(x^{\prime },y^{\prime })+\frac{\partial p_2}{\partial y^{\prime }}(x^{\prime },y^{\prime })\right]=0$$

$$\gamma(\vec{x})\;:=\;\left\{\begin{array}[]{ll}\gamma_{0}\cdot\exp\Big{(}-\tfrac{1}{1-\|\vec{x}\|^{2}}\Big{)},&\text{if }1>\|\vec{x}\|:=\max\{|x_{1}|,|x_{2}|\}\\ 0,&\text{otherwise}\end{array}\right.$$

$$\big{(}f\ast g\big{)}(\vec{x})\;=\;\int_{\mathbb{R}^{2}}f(\vec{x}-\vec{y})\cdot g(\vec{y})\,d\vec{y}$$

$$\big{(}(n,k,\vec{p}),(n^{\prime},k^{\prime},\vec{p}^{\prime})\big{)}\;\mapsto\;\left(\int\big{|}(\gamma_{n}\ast\mathbb{T}_{k}\vec{p})(\vec{x})-(\gamma_{n^{\prime}}\ast\mathbb{T}_{k^{\prime}}\vec{p}^{\prime})(\vec{x})\big{|}^{2}\,d\vec{x}\right)^{1/2}$$

$$\phi\big{(}(r_{i,j})_{0\leq i,j\leq N},(s_{i,j})_{0\leq i,j\leq N}\big{)}=\left\{\begin{array}[]{ll}1,&\mbox{if (\ref{div}), (\ref{boundary-1}), and (\ref{boundary-2}) are satisfied}\\ &\text{ (with $r_{i,j}=a^{1}_{i,j}$ and $s_{i,j}=a^{2}_{i,j}$)}\\ 0,&\mbox{otherwise}\end{array}\right.$$

$$\mathbb{P}\vec{u}=(-{\partial }_y\phi ,{\partial }_x\phi )$$

$$△\phi ={\partial }_x{u}_2-{\partial }_y{u}_1\text{ in }\Omega ,\phi =0\text{ on }\partial \Omega$$

$$\{\sin (n\pi x)\cos (m\pi y){\}}_{n\ge 1,m\ge 0},\text{mbox}or\{\cos (n\pi x)\sin (m\pi y){\}}_{n\ge 0,m\ge 1}$$

$$u_i(x,y)=$$

$$=$$

$$u_{i,n,m}={\int }_0^1{\int }_0^1{u}_i(x,y)\cos (n\pi x)\sin (m\pi y)dxdy,\text{mbox}and$$

$${\tilde{u}}_{i,n,m}={\int }_0^1{\int }_0^1{u}_i(x,y)\sin (n\pi x)\cos (m\pi y)dxdy$$

$$\phi (x,y)={\sum }_{n,m=1}^{\infty }\frac{-n{u}_{2,n,m}+m{\tilde{u}}_{1,n,m}}{(n^2+m^2)\pi }\sin (n\pi x)\sin (m\pi y)$$

$$-{\partial }_y\phi ={\sum }_{n,m\ge 1}\frac{mn{u}_{2,n,m}-m^2{\tilde{u}}_{1,n,m}}{n^2+m^2}\sin (n\pi x)\cos (m\pi y)$$

$$\Vert (-{\partial }_y\phi ,{\partial }_x\phi )-(p,q){\Vert }_2,(p,q)\in \mathcal{P}.$$

Full text

The paper is dedicated to the memory of Professor Ker-I Ko.

1 Introduction

The (physical) Church-Turing Hypothesis [17] postulates that every physical phenomenon or effect can, at least in principle, be simulated by a sufficiently powerful digital computer up to any desired precision. Its validity had been challenged, though, in the sound setting of Recursive Analysis: with a computable $C^{1}$ initial condition to the Wave Equation leading to an incomputable solution [11, 13]. The controversy was later resolved by demonstrating that, in both physically [30, 1] and mathematically more appropriate Sobolev space settings, the solution is computable uniformly in the initial data [23]. Recall that functions $f$ in a Sobolev space are not defined pointwise but by local averages in the $L_{q}$ sense111We use $q\in[1,\infty]$ to denote the norm index, $P$ for the pressure field, $p$ for polynomials, $\mathcal{P}$ for a set of trimmed and mollified tuples of the latter, and $\mathbb{P}$ for the Helmholtz Projection. (in particular $q=2$ corresponding to energy) with derivatives understood in the distributional sense. This led to a series of investigations on the computability of linear and nonlinear partial differential equations [24, 25, 26].

The (incompressible) Navier-Stokes Equation

[TABLE]

describes the motion of a viscous incompressible fluid filling a rigid box $\overline{\Omega}$. The vector field \vec{u}=\vec{u}(\vec{x},t)=\big{(}u_{1},u_{2},\ldots,u_{d}\big{)} represents the velocity of the fluid and $P=P(\vec{x},t)$ is the scalar pressure with gradient $\nabla P$; $\bigtriangleup$ is the Laplace operator; $\nabla\cdot\vec{u}$ denotes componentwise divergence; $\vec{u}\cdot\nabla$ means, in Cartesian coordinates, $u_{1}\partial_{x_{1}}+u_{2}\partial_{x_{2}}+\ldots+u_{d}\partial_{x_{d}}$; and the function $\vec{a}=\vec{a}(\vec{x})$ with $\nabla\cdot\vec{a}=0$ provides the initial velocity and $\vec{f}$ is a given external force. Equation (1) thus constitutes a system of $d+1$ partial differential equations for $d+1$ functions.

The question of global existence and smoothness of its solutions, even in the homogeneous case $\vec{f}\equiv\vec{0}$, is one of the seven Millennium Prize Problems posted by the Clay Mathematics Institute at the beginning of the 21st century. Local existence has been established, though, in various $L_{q}$ settings [5]; and uniqueness of weak solutions in dimension 2, but not in dimension 3 [18, §V.1.5], [2, §V.1.3.1]. Nevertheless, numerical solution methods have been devised in abundance, often based on pointwise (or even uniform, rather than $L_{q}$) approximation and struggling with nonphysical artifacts [14]. In fact, the very last of seven open problems listed in the addendum to [12] asks for a “recursion theoretic study of $\dots$ the Navier-Stokes equation”. Moreover it has been suggested [16] that hydrodynamics could in principle be incomputable in the sense of allowing to simulate universal Turing computation and to thus ‘solve’ the Halting problem. And indeed recent progress towards (a negative answer to) the Millennium Problem [20] proceeds by simulating a computational process in the vorticity dynamics to construct a blowup in finite time for a PDE very similar to (1).

1.1 Overview

Using the sound framework of Recursive Analysis, we assert the computability of a local strong solution of (1) in the space $L^{\sigma}_{2,0}(\Omega)$ (see Section 2 for definition) from a given initial condition $\vec{a}\in L^{\sigma}_{2,0}(\Omega)$; moreover, the computation is uniform in the initial data. We follow a common strategy used in the classical existence proofs [3, 4, 5, 18, 21]:

1. i)

Eliminate the pressure $P$ by applying, to both sides of Equation (1), the Helmholtz projection \mathbb{P}:\big{(}L_{2}(\Omega)\big{)}^{2}\to L^{\sigma}_{2,0}(\Omega), thus arriving at the non-linear evolution equation

[TABLE]

where $L_{2}(\Omega)$ is the set of all square-integrable real-valued functions defined on $\Omega$, $\vec{g}=\mathbb{P}\vec{f}$, $\vec{u}=\mathbb{P}\vec{u}\in L^{\sigma}_{2,0}(\Omega)$, $\mathbb{A}=-\mathbb{P}\bigtriangleup$ denotes the Stokes operator, and $\mathbb{B}\vec{u}=\mathbb{P}\,(\vec{u}\cdot\nabla)\vec{u}$ is the nonlinearity. 2. ii)

Construct a mild solution $\vec{v}(t)\vec{a}=e^{-t\mathbb{A}}\vec{a}$ of the associated homogeneous linear equation

[TABLE] 3. iii)

Rewrite (2) using (ii) in an integral form [5, §2]

[TABLE]

and solve it by a limit/fixed-point argument using the following iteration scheme [5, Eq. (2.1)]:

[TABLE] 4. iv)

Recover the pressure $P$ from $\vec{u}$ by solving

[TABLE]

To make use of the above strategy for deriving an algorithm to compute the solution of (1), there are several difficulties which need to be dealt with. Firstly, a proper representation is needed for coding the solenoidals. The codes should be not only rich enough to capture the functional characters of these vector fields but also robust enough to retain the coded information under elementary function operations, in particular, integration. Secondly, since the Stokes operator $\mathbb{A}:\mbox{dom}(\mathbb{A})\to L^{\sigma}_{2,0}(\Omega)$ is neither continuous nor its graph dense in $(L^{\sigma}_{2,0}(\Omega))^{2}$, there is no convenient way to directly code $\mathbb{A}$ for computing the solution of the linear equation (3). The lack of computer-accessible information on $\mathbb{A}$ makes the computation of the solution $\vec{v}(t)\vec{a}=e^{-t\mathbb{A}}\vec{a}$ of (3) much more intricate. Thirdly, since the nonlinear operator $\mathbb{B}$ in the iteration (5) involves differentiation and multiplication, and a mere $L^{\sigma}_{2,0}$-code of $\vec{v}_{n}$ is not rich enough for carrying out these operations, it follows that there is a need for computationally derive a stronger code for $\vec{v}_{n}$ from any given $L^{\sigma}_{2,0}$-code of $\vec{a}$ so that $\mathbb{B}\vec{v}_{n}$ can be computed. This indicates that the iteration is to move back and forth among different spaces, and thus additional care must be taken in order to keep the computations flowing in and out without any glitches from one space to another. To overcome those difficulties arising in the recursion theoretic study of the Navier-Stokes equation, many estimates - subtle and intricate - are established in addition to the classical estimates.

The paper is organized as follows. Presuming familiarity with Weihrauch’s Type-2 Theory of Effectivity [22], Section 2 recalls the standard representation $\delta_{L_{2}}$ of $L_{2}(\Omega)$ and introduces a natural representation $\delta_{L_{2,0}^{\sigma}}$ of $L^{\sigma}_{2,0}(\Omega)$. Section 3 proves that the Helmholtz projection \mathbb{P}:\big{(}L_{2}(\Omega)\big{)}^{2}\to L^{\sigma}_{2,0}(\Omega) is \big{(}(\delta_{L_{2}})^{2},\delta_{L_{2,0}^{\sigma}}\big{)}-computable. Section 4 presents the proof that the solution to the linear homogeneous Dirichlet problem (3) is uniformly computable from the initial condition $\vec{a}$. Section 5 is devoted to show that the solution to the nonlinear Navier-Stokes problem (1) is uniformly computable locally. Subsection 5.1 recalls the Bessel (=fractional Sobolev) space $H_{2}^{s}(\Omega)\subseteq L_{2}(\Omega)$ of $s$-fold weakly differentiable square-integrable functions on $\Omega$ and its associated standard representation $\delta_{H_{2,0}^{s}}$, $s\geq 0$. For $s>1$ we assert differentiation $H_{2}^{s}(\Omega)\ni w\mapsto\partial_{x}w\in L_{2}(\Omega)$ to be \big{(}\delta_{H_{2,0}^{s}},\delta_{L_{2}}\big{)}-computable and multiplication $H_{2}^{s}(\Omega)\times L_{2}(\Omega)\ni(v,w)\mapsto vw\in L_{2}(\Omega)$ to be \big{(}\delta_{H_{2,0}^{s}}\times\delta_{L_{2}},\delta_{L_{2}}\big{)}-computable. Based on these preparations, Subsection 5.3 asserts that in the homogeneous case $\vec{g}\equiv\vec{0}$, the sequence, generated from the iteration map

[TABLE]

according to Equation (5), converges effectively uniformly on some positive (but not necessarily maximal) time interval $[0;T]$ whose length $T=T(\vec{a})>0$ is computable from the initial condition $\vec{a}$. Subsection 5.4 proves that the iteration map $\mathbb{S}$ is \big{(}[\rho\!\to\!\delta_{L_{2,0}^{\sigma}}]\times[\rho\!\to\!\delta_{L_{2,0}^{\sigma}}],[\rho\!\to\!\delta_{L_{2,0}^{\sigma}}]\big{)}-computable. We conclude in Subsection 5.5 with the final extensions regarding the inhomogeneity $\vec{f}$ and pressure $P$, thus establishing the main result of this work:

Theorem 1.1

There exists a \displaystyle\big{(}\delta_{L_{2,0}^{\sigma}}\!\!\times\!\big{[}\rho\!\to\!\delta_{L_{2,0}^{\sigma}}\big{]}\>,\rho\big{)}\text{-computable map }T,

[TABLE]

and a \big{(}\delta_{L_{2,0}^{\sigma}}\!\!\!\times\!\big{[}\rho\!\to\!\delta_{L_{2,0}^{\sigma}}\big{]}\times\rho\>,\delta_{L_{2,0}^{\sigma}}\big{)}\text{-computable partial map }\cal S,

[TABLE]

such that, for every $\vec{a}\in L^{\sigma}_{2,0}(\Omega)$ and \vec{f}\in C\big{(}[0;\infty),L^{\sigma}_{2,0}(\Omega)\big{)}, the function (\vec{u},P):\big{[}0;T(\vec{a},\vec{f})]\ni t\mapsto\mathcal{S}(\vec{a},\vec{f},t) constitutes a (strong local in time and weak global in space) solution to Equation (1).

Roughly speaking, a function is computable if it can be approximated by “computer-accessible” functions (such as rational numbers, polynomials with rational coefficients, and so forth) with arbitrary precision, where precision is given as an input; such a sequence of approximations is called an effective approximation. Thus in terms of effective approximations, the theorem states that the solution of Equation (1) can be effectively approximated locally in the time interval $[0,T(\vec{a},\vec{f})]$, where the time instance $T(\vec{a},\vec{f})$ is effectively approximable.

More precisely, in computable analysis, a map $F:X\to Y$ from a space $X$ with representation $\delta_{X}$ to a space $Y$ with representation $\delta_{Y}$ is said to be $(\delta_{X},\delta_{Y})$-computable if there exists a (Turing) algorithm (or any computer program) that computes a $\delta_{Y}$-name of $F(x)$ from any given $\delta_{X}$-name of $x$. A metric space $(X,d)$, equipped with a partial enumeration $\zeta:\subseteq\mathbb{N}\to X$ of some dense subset, gives rise to a canonical Cauchy representation $\delta_{\zeta}$ by encoding each $x\in X$ with a sequence $\bar{s}=(s_{0},s_{1},s_{2},\ldots)\in\mbox{dom}(\zeta)^{\omega}\subseteq\mathbb{N}^{\omega}$ such that d\big{(}x,\zeta(s_{k})\big{)}\leq 2^{-k} for all $k$ [22, §8.1]; in other words, $\{\zeta(s_{k})\}$ is an effective approximation of $x$. For example, approximating by (dyadic) rationals thus leads to the standard representation $\rho$ of $\mathbb{R}$; and for a fixed bounded $\Omega\subseteq\mathbb{R}^{2}$, the standard representation $\delta_{L_{2}}$ of $L_{2}(\Omega)$ encodes $f\in L_{2}(\Omega)$ by a sequence $\{p_{k}:k\in\mathbb{N}\}\subseteq\mathbb{Q}[\mathbb{R}^{2}]$ of $d$-variate polynomials with rational coefficients such that $\|f-p_{k}\|_{2}\leq 2^{-k}$, where $\|\cdot\|_{2}=\|\cdot\|_{L_{2}}$. Thus if both spaces $X$ and $Y$ admit Cauchy representations, then a function $f:X\to Y$ is computable if there is an algorithm that computes an effective approximation of $f(x)$ on any given effective approximation of $x$ as input. For represented spaces $(X,\delta_{X})$ and $(Y,\delta_{Y})$, $\delta_{X}\!\times\!\delta_{Y}$ denotes the canonical representation of the Cartesian product $X\!\times\!Y$. When $X$ and $Y$ are $\sigma$-compact metric spaces with respective canonical Cauchy representations $\delta_{X}$ and $\delta_{Y}$, $[\delta_{X}\!\to\!\delta_{Y}]$ denotes a canonical representation of the space $C(X,Y)$ of continuous total functions $f:X\to Y$, equipped with the compact-open topology [22, Theorem 3.2.11+Definition 3.3.13]. The representation $[\delta_{X}\!\to\!\delta_{Y}]$ supports type conversion in the following sense [22, Theorem 3.3.15]:

Fact 1.2

*On the one hand, the evaluation $(f,x)\mapsto f(x)$ is $([\delta_{X}\!\to\!\delta_{Y}]\!\times\!\delta_{X}\>,\>\delta_{Y})$-computable. On the other hand, a map $f:X\times Y\to Z$ is $(\delta_{X}\!\times\!\delta_{Y}\>,\>\delta_{Z})$-computable iff the map X\ni x\mapsto\big{(}y\mapsto f(x,y)\big{)}\in C(Y,Z) is $(\delta_{X}\>,\>[\delta_{Y}\!\to\!\delta_{Z}])$-computable.

We mention in passing that all spaces considered in this paper are equipped with a norm. Thus for any space $X$ considered below, a $\delta_{X}$-name of $f\in X$ is simply an effective approximation of $f$ despite the often cumbersome notations.

2 Representing Divergence-Free $L_{2}$ Functions on $\Omega$

Let us call a vector field $\vec{f}$ satisfying $\nabla\cdot\vec{f}=0$ in $\Omega$ divergence-free. A vector-valued function $\vec{p}$ is called a polynomial of degree $N$ if each of its components is a polynomial of degree less than or equal to $N$ with respect to each variable and at least one component is a polynomial of degree $N$. Let $L^{\sigma}_{2,0}(\Omega)$ — or $L^{\sigma}_{2,0}$ if the context is clear — be the closure in $L_{2}$-norm of the set $\{\vec{u}\in(C^{\infty}_{0}(\Omega))^{2}\,:\,\nabla\cdot\vec{u}=0\}$ of all smooth divergence-free functions with support of $\vec{u}$ and all of its partial derivatives contained in some compact subset of $\Omega$. Let $\mathbb{Q}[\mathbb{R}^{2}]$ be the set of all polynomials of two real variables with rational coefficients and $\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$ the subset of all $2$-tuples of such polynomials which are divergence-free in $\Omega$ and vanish on $\partial\Omega$. We note that the boundary value of a $L^{\sigma}_{2,0}(\Omega)$-function $\vec{u}$, $\vec{u}|_{\partial\Omega}$, is not defined unless $\vec{u}$ is (weakly) differentiable; if $\vec{u}$ is (weakly) differentiable, then $\vec{u}|_{\partial\Omega}=0$.

Notation 2.1

Hereafter we use $\|w\|_{2}$ for the $L_{2}$-norm $\|w\|_{L_{2}(\Omega)}$ if $w$ is real-valued, or for $\|w\|_{(L_{2}(\Omega))^{2}}$ if $w$ is vector-valued (in $\mathbb{R}^{2}$). We note that $\|\cdot\|_{L^{\sigma}_{2,0}(\Omega)}=\|\cdot\|_{(L_{2}(\Omega))^{2}}$. For any subset $A$ of $\mathbb{R}^{n}$, its closure is denoted as $\overline{A}$.

Proposition 1

1. a)

A polynomial tuple

[TABLE]

is divergence-free and boundary-free if and only if its coefficients satisfy the following system of linear equations with integer coefficients:

[TABLE]

and for all $0\leq i,j\leq N$,

[TABLE] 2. b)

$\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$* is dense in $L^{\sigma}_{2,0}(\Omega)$ w.r.t. $L_{2}$-norm.*

The proof of Proposition 1 is deferred to Appendix 0.A.

We may be tempted to use $\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$ as a set of names for coding/approximating the elements in the space $L^{\sigma}_{2,0}(\Omega)$. However, since the closure of $\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$ in $L_{2}$-norm contains $L^{\sigma}_{2,0}(\Omega)$ as a proper subspace, $\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$ is “too big” to be used as a set of codes for representing $L^{\sigma}_{2,0}(\Omega)$; one has to “trim” polynomials in $\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$ so that any convergent sequence of “trimmed” polynomials converges to a limit in $L^{\sigma}_{2,0}(\Omega)$. The trimming process is shown below. For each $k\in\mathbb{N}$ (where $\mathbb{N}$ is the set of all positive integers), let $\Omega_{k}=(-1+2^{-k};1-2^{-k})^{2}$. And for each $\vec{p}=(p_{1},p_{2})\in\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$, define $\mathbb{T}_{k}\vec{p}=(\mathbb{T}_{k}p_{1},\mathbb{T}_{k}p_{2})$, where

[TABLE]

$j=1,2$. Then $\mathbb{T}_{k}p_{j}$ and $\mathbb{T}_{k}\vec{p}$ have the following properties:

1. a)

$\mathbb{T}_{k}\vec{p}$ has compact support $\overline{\Omega}_{k}$ contained in $\Omega$. 2. b)

$\mathbb{T}_{k}\vec{p}$ is a polynomial with rational coefficients defined in $\Omega_{k}$. 3. c)

$\mathbb{T}_{k}\vec{p}$ is continuous on $\overline{\Omega}=[-1,1]^{2}$. 4. d)

$\mathbb{T}_{k}\vec{p}=0$ on $\partial\Omega_{k}$, for $\vec{p}$ vanishes on the boundary of $\Omega$. Thus $\mathbb{T}_{k}\vec{p}$ vanishes in the exterior region of $\Omega_{k}$ including its boundary $\partial\Omega_{k}$. 5. e)

$\mathbb{T}_{k}\vec{p}$ is divergence-free in $\Omega_{k}$ following the calculation below: for $(x,y)\in\Omega_{k}$, we have $(\frac{x}{1-2^{-k}},\frac{y}{1-2^{-k}})\in\Omega$ and

[TABLE]

for $\vec{p}$ is divergence-free in $\Omega$, where $x^{\prime}=\frac{x}{1-2^{-k}}$ and $y^{\prime}=\frac{y}{1-2^{-k}}$.

It follows from the discussion above that every $\mathbb{T}_{k}\vec{p}$ is a divergence-free polynomial of rational coefficients on $\Omega_{k}$ that vanishes in $[-1,1]^{2}\setminus\Omega_{k}$ and is continuous on $[-1,1]^{2}$. However, although the functions $\mathbb{T}_{k}\vec{p}$ are continuous on $[-1,1]^{2}$ and differentiable in $\Omega_{k}$, they can be non-differentiable along the boundary $\partial\Omega_{k}\subseteq\Omega$. To use these functions as names for coding elements in $L^{\sigma}_{2,0}(\Omega)$, it is desirable to smoothen them along the boundary $\partial\Omega_{k}$ so that they are differentiable in the entire $\Omega$. A standard technique for smoothing a function is to convolute it with a $C^{\infty}$ function. We use this technique to modify functions $\mathbb{T}_{k}\vec{p}$ so that they become divergence-free and differentiable on the entire region of $\Omega$. Let

[TABLE]

where $\gamma_{0}$ is a constant such that $\int_{\mathbb{R}^{2}}\gamma(\vec{x})\,d\vec{x}=1$ holds. The constant $\gamma_{0}$ is computable, since integration on continuous functions is computable [22, §6.4]. Let $\gamma_{k}(\vec{x})=2^{2k}\gamma(2^{k}\vec{x})$. Then, for all $k\in\mathbb{N}$, $\gamma_{k}$ is a $C^{\infty}$ function having support in the closed square $[-2^{-k},2^{-k}]^{2}$ and $\int_{\mathbb{R}^{2}}\gamma_{k}(\vec{x})\,d\vec{x}=1$. Recall that for differentiable functions $f,g:\mathbb{R}^{2}\to\mathbb{R}$ with compact support, the convolution $f\ast g$ is defined as follows:

[TABLE]

It is easy to see that for $n\geq k+1$ the support of $\gamma_{n}\ast\mathbb{T}_{k}\vec{p}:=(\gamma_{n}\ast\mathbb{T}_{k}p_{1},\gamma_{n}\ast\mathbb{T}_{k}p_{2})$ is contained in the closed square $[-1+2^{-(k+1)},1-2^{-(k+1)}]^{2}$. It is also known classically that $\gamma_{n}\ast\mathbb{T}_{k}\vec{p}$ is a $C^{\infty}$ function. Since $\gamma_{n}$ is a computable function and integration on compact domains is computable, the map $(n,k,\vec{p})\mapsto\gamma_{n}\ast\mathbb{T}_{k}\vec{p}$ is computable. Moreover the following metric is computable:

[TABLE]

Lemma 1

Every function $\gamma_{n}\ast\mathbb{T}_{k}\vec{p}$ is divergence-free in $\Omega$, where $n,k\in\mathbb{N}$, $n\geq k$, and $\vec{p}\in\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$.

Lemma 2

The set \mathcal{P}=\big{\{}\gamma_{n}\ast\mathbb{T}_{k}\vec{p}\>:\>n,k\in\mathbb{N},n\geq k+1,\vec{p}\in\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]\big{\}} is dense in $L^{\sigma}_{2,0}(\Omega)$.

See Appendices 0.B and 0.C for the proofs.

From Lemmas 1 and 2 it follows that $\mathcal{P}$ is a countable set that is dense in $L^{\sigma}_{2,0}(\Omega)$ (in $L_{2}$-norm) and every function in $\mathcal{P}$ is $C^{\infty}$, divergence-free on $\Omega$, and having a compact support contained in $\Omega$; in other words, $\mathcal{P}\subset\{\vec{u}\in C^{\infty}_{0}(\Omega)^{2}\,:\,\nabla\cdot\vec{u}=0\}$. Thus, L^{\sigma}_{2,0}(\Omega)=\mbox{the closure of \mathcal{P}$in$L_{2}-norm}. This fact indicates that the set $\mathcal{P}$ is qualified to serve as codes for representing $L^{\sigma}_{2,0}(\Omega)$.

Since the function $\phi:\bigcup_{N=0}^{\infty}\mathbb{Q}^{(N+1)^{2}}\times\mathbb{Q}^{(N+1)^{2}}\to\{0,1\}$, where

[TABLE]

is computable, there is a total computable function on $\mathbb{N}$ that enumerates $\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$. Then it follows from the definition of $\mathcal{P}$ that there is a total computable function $\alpha:\mathbb{N}\to\mathcal{P}$ that enumerates $\mathcal{P}$; thus, in view of the computable equation (13), \big{(}L^{\sigma}_{2,0}(\Omega),(\vec{u},\vec{v})\mapsto\|\vec{u}-\vec{v}\|_{2},\mathcal{P},\alpha\big{)} is a computable metric space. Let $\delta_{L_{2,0}^{\sigma}}:\mathbb{N}^{\omega}\to L^{\sigma}_{2,0}$ be the standard Cauchy representation of $L^{\sigma}_{2,0}$; that is, every function $\vec{u}\in L^{\sigma}_{2,0}(\Omega)$ is encoded by a sequence $\{\vec{p}_{k}:k\in\mathbb{N}\}\subseteq\mathcal{P}$, such that $\|\vec{u}-\vec{p}_{k}\|_{2}\leq 2^{-k}$. The sequence $\{\vec{p}_{k}\}_{k\in\mathbb{N}}$ is called a $\delta_{L_{2,0}^{\sigma}}$-name of $\vec{u}$, which is an effective approximation of $\vec{u}$ (in $L_{2}$-norm).

3 Computability of Helmholtz Projection

In this section, we show that the Helmholtz projection $\mathbb{P}$ is computable.

Proposition 2

The projection \mathbb{P}:\big{(}L_{2}(\Omega)\big{)}^{2}\to L^{\sigma}_{2,0}(\Omega) is \big{(}(\delta_{L_{2}})^{2},\delta_{L_{2,0}^{\sigma}}\big{)}-computable.

Proof

For simplicity let us set $\Omega=(0,1)^{2}$. The proof carries over to $\Omega=(-1,1)^{2}$ by a scaling on sine and cosine functions. We begin with two classical facts which are used in the proof:

• (i)

It follows from [6, pp.40]/[21] that for any \vec{u}=\big{(}u_{1},u_{2}\big{)}\in\big{(}L_{2}(\Omega)\big{)}^{2},

[TABLE]

where the scalar function $\varphi$ is the solution of the following boundary value problem:

[TABLE]

We note that $\mathbb{P}$ is a linear operator.

• (ii)

Each of $\{\sin(n\pi x)\sin(m\pi y)\}_{n,m\geq 1}$,

[TABLE]

is an orthogonal basis for $L^{2}(\Omega)$. Thus each \vec{u}=\big{(}u_{1},u_{2}\big{)} in \big{(}L_{2}(\Omega)\big{)}^{2}, $u_{i}$, $i=1$ or $2$, can be written in the following form:

[TABLE]

where

[TABLE]

[TABLE]

with the property that $\|u_{i}\|_{2}=\left(\sum_{n,m\geq 0}|u_{i,n,m}|^{2}\right)^{1/2}=\left(\sum_{n,m\geq 0}|\tilde{u}_{i,n,m}|^{2}\right)^{1/2}$. We note that the sequences $\{u_{i,n,m}\}$, $\{\tilde{u}_{i,n,m}\}$, and $\|u_{i}\|_{2}$ are computable from $\vec{u}$; cf. [22].

To prove that the projection is \big{(}(\delta_{L_{2}})^{2},\delta_{L_{2,0}^{\sigma}}\big{)}-computable, it suffices to show that there is an algorithm computing, given any accuracy $k\in\mathbb{N}$ and for any $\vec{u}\in(L^{2}(\Omega))^{2}$, a vector function $(p_{k},q_{k})\in\mathcal{P}$ such that $\|\mathbb{P}\vec{u}-(p_{k},q_{k})\|_{2}\leq 2^{-k}$. Let us fix $k$ and \vec{u}=\big{(}u_{1},u_{2}). Then a straightforward computation shows that the solution $\varphi$ of (15) can be explicitly written as

[TABLE]

It then follows that

[TABLE]

Similarly, we can obtain a formula for $\partial_{x}\varphi$. Since we have an explicit expression for $(-\partial_{y}\varphi,\partial_{x}\varphi)$, a search algorithm is usually a preferred choice for finding a $k$-approximation $(p_{k},q_{k})$ of $\mathbb{P}\vec{u}$ by successively computing the norms

[TABLE]

However, since $-\partial_{y}\varphi$ and $\partial_{x}\varphi$ are infinite series which involve limit processes, a truncating algorithm is needed so that one can compute approximations of the two limits before a search program can be executed. The truncating algorithm will find some $N(k,\vec{u})\in\mathbb{N}$ such that the $N(k,\vec{u})$-partial sum of $(-\partial_{y}\varphi,\partial_{x}\varphi)$ is a $2^{-(k+1)}$-approximation of the series; in other words, the algorithm chops off the infinite tails of the series within pre-assigned accuracy. The following estimate provides a basis for the desired truncating algorithm:

[TABLE]

A similar estimate applies to $\partial_{x}\varphi$. Since

[TABLE]

is computable, it follows that there is an algorithm computing $N(k,\vec{u})$ from $k$ and $\vec{u}$ such that the $N(k,\vec{u})$-partial sum of $(-\partial_{y}\varphi,\partial_{x}\varphi)$ is a $2^{-(k+1)}$-approximation of the series. Now we can search for $(p_{k},q_{k})\in\mathbb{P}$ that approximates the $N(k,\vec{u})$-partial sum in $L^{2}$-norm within the accuracy $2^{-(k+1)}$ as follows: enumerate $\mathcal{P}=\{\tilde{p}_{j}\}$, compute the $L^{2}$-norm of the difference between the $N(k,\vec{u})$-partial sum and $\tilde{p}_{j}$, halt the computation at $\tilde{p}_{j}$ when the $L^{2}$-norm is less that $2^{-(k+1)}$, and then set $(p_{k},q_{k})=\tilde{p}_{j}$. We note that each computation halts in finitely many steps. The search will succeed since $\mathbb{P}\vec{u}=(-\partial_{y}\varphi,\partial_{x}\varphi)\in L^{\sigma}_{2,0}$ and $\mathcal{P}$ is dense in $L^{\sigma}_{2,0}$. It is then clear that $\|\mathbb{P}\vec{u}-(p_{k},q_{k})\|_{2}\leq 2^{-k}$.

4 Computability of the linear problem

In this section, we show that the solution operator for the linear homogeneous equation (3) is uniformly computable from the initial data. We begin by recalling the Stokes operator and some of its classical properties. Let $\mathbb{A}=-\mathbb{P}\bigtriangleup$ be the Stokes operator as defined for instance in [3, §2] or [18, §III.2.1], where \mathbb{P}:\big{(}L_{2}(\Omega)\big{)}^{2}\to L^{\sigma}_{2,0} is the Helmholtz projection. It is known from the classical study that $\mathbb{A}$ is an unbounded but closed positively self-adjoint linear operator whose domain is dense in $L^{\sigma}_{2,0}$, and thus $-\mathbb{A}$ is the infinitesimal generator of an analytic semigroup; cf. [18, Theorem III.2.1.1] or [2, §IV.5.2]. In this case, the linear homogeneous equation (3) has the solution $\vec{u}(t)=e^{-\mathbb{A}t}\vec{a}$, where $\vec{u}(0)=\vec{a}$, $e^{-\mathbb{A}t}$ is the analytic semigroup generated by the infinitesimal generator $-\mathbb{A}$, and $\vec{u}(t)\in L^{\sigma}_{2,0}(\Omega)$ for $t\geq 0$. Furthermore, the following lemma shows that the solution $\vec{u}(t)$ decays in $L^{2}$-norm as time $t$ increases.

Lemma 3

For every $\vec{a}\in L^{\sigma}_{2,0}(\Omega)$ and $t\geq 0$,

[TABLE]

(Recall that $\|\cdot\|_{2}=\|\cdot\|_{L^{\sigma}_{2,0}(\Omega)}$; see Notation 2.1.)

Proof

Classically it is known that for any $\vec{a}\in L^{\sigma}_{2,0}(\Omega)$, $\vec{u}(t)=e^{-t\mathbb{A}}\vec{a}$ is in the domain of $\mathbb{A}$ for $t>0$. Thus if $\vec{a}=\vec{u}(0)$ itself is in the domain of $\mathbb{A}$, then so is $\vec{u}(t)$ for $t\geq 0$. Since $\mathbb{A}$ is positively self-adjoint, it follows that $\mathbb{A}^{\ast}=\mathbb{A}$ and $\langle\mathbb{A}\vec{u}(t),\vec{u}(t)\rangle:=\int_{\Omega}\mathbb{A}\vec{u}(t)(\vec{x})\cdot\vec{u}(t)(\vec{x})\,d\vec{x}>0$ for every $\vec{a}$ in the domain of $\mathbb{A}$ with $\vec{a}\not\equiv 0$ and $t\geq 0$. Now if we rewrite the equation (3) in the form of

[TABLE]

or equivalently $\frac{1}{2}\frac{d}{dt}\langle\vec{u},\vec{u}\rangle+\langle\mathbb{A}\vec{u},\vec{u}\rangle=0$, then $\frac{d}{dt}\langle\vec{u},\vec{u}\rangle\leq 0$ and consequently $\langle\vec{u},\vec{u}\rangle(t)\leq\langle\vec{u},\vec{u}\rangle(0)$; thus (17) holds true for $\vec{a}$ in the domain of $\mathbb{A}$. Since the domain of $\mathbb{A}$ is dense in $L^{\sigma}_{2,0}(\Omega)$, it follows that (17) holds true for all $\vec{a}\in L^{\sigma}_{2,0}(\Omega)$.     $\sqcup$

Proposition 3

For the linear homogenous equation (3), the solution operator \mathcal{S}:L^{\sigma}_{2,0}(\Omega)\to C\big{(}[0;\infty),L^{\sigma}_{2,0}(\Omega)\big{)}, $\vec{a}\mapsto(t\mapsto e^{-\mathbb{A}t}\vec{a})$, is $(\delta_{L_{2,0}^{\sigma}},[\rho\to\delta_{L_{2,0}^{\sigma}}])$-computable.

By the First Main Theorem of Pour-El and Richards [12, §II.3], the unbounded operator $\mathbb{A}$ does not preserve computability. In particular, the naive exponential series $\sum_{n}(-\mathbb{A}t)^{n}\vec{a}/n!$ does not establish Proposition 3.

Convention. For readability we will not notationally distinguish the spaces of vectors $\vec{u},\vec{a}$ and scalar functions $u,a$ in the proof below and the proof of Lemma  6.

Proof

We show how to compute a $\delta_{L_{2,0}^{\sigma}}$-name of $e^{-t\mathbb{A}}a$ on inputs $t\geq 0$ and $a\in L^{\sigma}_{2,0}(\Omega)$. Recall that a $\delta_{L_{2,0}^{\sigma}}$-name of $e^{-t\mathbb{A}}a$ is a sequence $\{q_{K}\}$, $q_{K}\in\mathcal{P}$, satisfying $\|e^{-t\mathbb{A}}a-q_{K}\|_{2}\leq 2^{-K}$ for all $K\in\mathbb{N}$. Again, for readability, we assume that $\Omega=(0,1)^{2}$.

We first consider the case where $a\in\mathcal{P}$ and $t>0$. The reason for us to start with functions in $\mathcal{P}$ is that these functions have stronger convergence property in the sense that, for any $a\in\mathcal{P}$, if $a=(a^{1},a^{2})$ is expressed in terms of the orthogonal basis $\{\sin(n\pi x)\sin(m\pi y)\}_{n,m\geq 1}$ for $L_{2}(\Omega)$: for $i=1,2$,

[TABLE]

where $a^{i}_{n,m}=\int_{0}^{1}\int_{0}^{1}a^{i}\sin(n\pi x)\sin(m\pi y)\,dx\,dy$, then the following series is convergent

[TABLE]

The inequality (19) holds true because functions in $\mathcal{P}$ are $C^{\infty}$. In fact, the series is not only convergent but its sum is also computable (from $a$) (see, for example, [28]).

Now let $K\in\mathbb{N}$ be any given precision. Since $-\mathbb{A}$ generates an analytic semigroup, it follows from [10, Section 2.5] that for $t>0$,

[TABLE]

where $\Gamma$ is the path composed from two rays $re^{i\beta}$ and $re^{-i\beta}$ with $0<r<\infty$ and $\beta=\frac{3\pi}{5}$. Thus we have an explicit expression for $e^{-t\mathbb{A}}a$, which involves a limit process – an infinite integral – indicating that a search algorithm is applicable for finding a desirable $K$-approximation provided that a finite approximation of $e^{-t\mathbb{A}}a$ can be computed by some truncating algorithm.

In the following, we construct such a truncating algorithm. We begin by writing the infinite integral in (20) as a sum of three integrals: two are finite and one infinite; the infinite one can be made arbitrarily small. Now for the details. Let $l$ be a positive integer to be determined; let $\Gamma_{1}$ be the path $re^{i\beta}$ with $0<r\leq l$; $\Gamma_{2}$ the path $re^{-i\beta}$ with $0<r\leq l$; and $\Gamma_{3}=\Gamma\setminus(\Gamma_{1}\cup\Gamma_{2})$. Since $a\in\mathcal{P}$, it follows that $-\mathbb{A}a=\mathbb{P}\bigtriangleup a=\bigtriangleup a$, which further implies that

[TABLE]

Note that for any $\lambda\in\Gamma$, $|\lambda+(n\pi)^{2}+(m\pi)^{2}|\neq 0$. From (20) and (4) we can write $e^{-t\mathbb{A}}a$ as a sum of three terms:

[TABLE]

where $\tilde{a}=(\lambda\mathbb{I}+\mathbb{A})^{-1}a$. The functions $\beta_{j}$, $j=1,2,3$, are in $L^{\sigma}_{2,0}(\Omega)$ as verified as follows: It follows from $a=(\lambda\mathbb{I}+\mathbb{A})\tilde{a}=(\lambda\mathbb{I}-\mathbb{P}\bigtriangleup)\tilde{a}$ and $\mathbb{P}\bigtriangleup\tilde{a}=\mathbb{P}(\bigtriangleup\tilde{a})\in L^{\sigma}_{2,0}(\Omega)$ that $\bigtriangledown(\mathbb{P}\bigtriangleup\tilde{a})=0$ and

[TABLE]

Since $\lambda\in\Gamma$, it follows that $\lambda\neq 0$; thus $\bigtriangledown\tilde{a}=0$. This shows that $\tilde{a}\in L^{\sigma}_{2,0}(\Omega)$. Then it follows from (22) that

[TABLE]

Hence $\beta_{j}\in L^{\sigma}_{2,0}(\Omega)$ for $1\leq j\leq 3$.

Next we show that $\beta_{1}$ and $\beta_{2}$ can be effectively approximated by finite sums while $\beta_{3}$ tend to zero effectively as $l\to\infty$. We start with $\beta_{3}$. Since $t>0$ and $\cos\beta=\cos\frac{3\pi}{5}<0$, it follows that

[TABLE]

effectively as $l\to\infty$. Thus there is some $l_{K}\in\mathbb{N}$, computable from $a$ and $t$, such that the following estimate is valid for $i=1,2$ when we take $l$ to be $l_{K}$:

[TABLE]

where $\beta_{3}=(\beta^{1}_{3},\beta^{2}_{3})$. Now let us set $l=l_{K}$ and estimate $\beta_{1}$. Since $\beta=\frac{3\pi}{5}<\frac{3\pi}{4}$, it follows that $\cos\beta<0$ and $|\cos\beta|<\sin\beta$. Consequently, for any $\lambda=re^{i\beta}$ on $\Gamma_{1}$, if $r\geq\frac{1}{\sin\beta}$, then $|re^{i\beta}+(n\pi)^{2}+(m\pi)^{2}|\geq r\sin\beta\geq 1$. On the other hand, if $0<r<\frac{1}{\sin\beta}$, then $r\cos\beta+(n\pi)^{2}+(m\pi)^{2}\geq\pi^{2}(n^{2}+m^{2})-r\sin\beta\geq 2\pi^{2}-1>1$, which implies that $|re^{i\beta}+(n\pi)^{2}+(m\pi)^{2}|\geq|r\cos\beta+(n\pi)^{2}+(m\pi)^{2}|>1$. Thus $|\lambda+(n\pi)^{2}+(m\pi)^{2}|\geq 1$ for every $\lambda\in\Gamma_{1}$. And so

[TABLE]

This estimate together with (19) implies that there exists a positive integer $k=k(t,a,K)$, computable from $t>0$, $a$ and $K$, such that

[TABLE]

Write $\beta_{1}(k)=(\beta^{1}_{1}(k),\beta^{2}_{1}(k))$ with

[TABLE]

$i=1,2$. Then

[TABLE]

Similarly, if we write $\beta_{2}(k)=(\beta^{1}_{2}(k),\beta^{2}_{2}(k))$ with

[TABLE]

then $\|\beta_{2}-\beta_{2}(k)\|_{2}\leq 2^{-(K+7)}$. The construction of the truncating algorithm is now complete; the algorithm outputs $\beta_{1}(k)+\beta_{2}(k)$ (uniformly) on the inputs $a\in\mathcal{P}$, $t>0$, and precision $K$; the output has the property that it is a finite sum involving a finite integral and $\|\beta_{1}(k)+\beta_{2}(k)-e^{-t\mathbb{A}}a\|_{2}\leq 2^{-(K+4)}$.

Now we are able to search for a desirable approximation in $\mathcal{P}$. Let us list $\mathcal{P}=\{\phi_{j}\,:\,j\in\mathbb{N}\}$ and compute $\|\phi_{j}-(\beta_{1}(k)+\beta_{2}(k))\|_{2}$. Halt the computation at $j=j(K)$ when

[TABLE]

The computation will halt since $\beta_{1},\beta_{2}\in L^{\sigma}_{2,0}(\Omega)$, $\|\beta_{1}-\beta_{1}(k)\|_{2}\leq 2^{-(K+7)}$, $\|\beta_{2}-\beta_{2}(k)\|_{2}\leq 2^{-(K+7)}$, and $\mathcal{P}$ is dense in $L^{\sigma}_{2,0}(\Omega)$ (in $L^{2}$-norm). Set $q_{K}=\phi_{j(K)}$. Then

[TABLE]

Next we consider the case where $a\in L^{\sigma}_{2,0}(\Omega)$ and $t>0$. In this case, the input $a$ is presented by (any) one of its $\delta_{L^{\sigma}_{2,0}}$-names, say $\{a_{k}\}$, where $a_{k}\in\mathcal{P}$. It is then clear from the estimate (17) and the discussion above that there is an algorithm that computes a $K$-approximation $p_{K}\in\mathcal{P}$ on inputs $t>0$, $a$ and precision $K$ such that $\|p_{K}-e^{-t\mathbb{A}}a\|_{2}\leq 2^{-K}$.

Finally we consider the case where $t\geq 0$ and $a\in L^{\sigma}_{2,0}(\Omega)$. Since $e^{-t\mathbb{A}}a=a$ for $t=0$ and we already derived an algorithm for computing $e^{-t\mathbb{A}}a$ for $t>0$, it suffices to show that $e^{-t\mathbb{A}}a\to a$ in $L^{2}$-norm effectively as $t\to 0$. Let $\{a_{k}\}$ be a $\delta_{L^{\sigma}_{2,0}}$-name of $a$. It follows from Theorem 6.13 of Section 2.6 [Paz83] that $\|e^{-t\mathbb{A}}a_{k}-a_{k}\|\leq Ct^{1/2}\|\mathbb{A}^{1/2}a_{k}\|$. Thus

[TABLE]

the right-hand side goes to 0 effectively as $t\to 0$.     $\sqcup$

We note that the computation of the approximations $q_{K}$ of $e^{-t\mathbb{A}}a$ does not require encoding $\mathbb{A}$. Let $W:L^{\sigma}_{2,0}(\Omega)\times[0,\infty)\to L^{\sigma}_{2,0}(\Omega)$, $(a,t)\mapsto e^{-t\mathbb{A}}a$. Then it follows from the previous Proposition and Fact 1.2 that $W$ is computable.

5 Extension to the nonlinear problem

We now proceed to the nonlinear problem (2) by solving its integral version (4) via the iteration scheme (5) but first restrict to the homogeneous case $\vec{g}\equiv\vec{0}$:

[TABLE]

Classically, it is known that for every initial condition $\vec{a}\in L^{\sigma}_{2,0}(\Omega)$ the sequence $\vec{u}_{m}=\vec{u}_{m}(t)$ converges near $t=0$ to a unique limit $\vec{u}$ solving (4) and thus (2). Since there is no explicit formula for the solution $\vec{u}$, the truncation/search type of algorithms such as those used in the proofs of Propositions 2 and 3 is no longer applicable for the nonlinear case. Instead, we use a method based on the fixed-point argument to establish the computability of $\vec{u}$. We shall show that the limit of the above sequence $\vec{u}_{m}=\vec{u}_{m}(t)$ has an effective approximation. The proof consists of two parts: first we study the rate of convergence and show that the sequence converges at a computable rate as $m\to\infty$ for $t\in[0;T]$ with some $T=T_{\vec{a}}>0$, where $T_{\vec{a}}$ is computable from $\vec{a}$; then we show that the sequence – as one entity – can be effectively approximated starting with the given $\vec{a}$. The precise statements of the two tasks are given in the following two propositions.

Proposition 4

There is a computable map $\mathbb{T}:L^{\sigma}_{2,0}(\Omega)\to(0,\infty)$, $\vec{a}\mapsto T_{\vec{a}}$, such that the sequence $\{{\vec{u}}_{m}\}$ converges effectively in $m$ and uniformly for $t\in[0;T_{\vec{a}}]$.

Recall that a sequence $\{x_{m}\}$ in a metric space $(X,d)$ is effectively convergent if $d(x_{m},x_{m+1})\leq 2^{-m}$. In view of type conversion (Subsection 1.1), the following proposition asserts (ii):

Proposition 5

The map $\mathbb{S}:\mathbb{N}\times L^{\sigma}_{2,0}(\Omega)\times[0,\infty)\to L^{\sigma}_{2,0}(\Omega)$, $(m,\vec{a},t)\to\vec{u}_{m}(t)$ according to Equation (23), is \big{(}\nu\times\delta_{L_{2,0}^{\sigma}}\times\rho,\delta_{L_{2,0}^{\sigma}}\big{)}-computable.

The main difficulties in proving the two propositions are rooted in the nonlinearity of $\mathbb{B}$: the nonlinear operator $\mathbb{B}$ requires greater care in estimating the rate of convergence and demands richer codings for computation. Since information on $\mathbb{B}\vec{u}_{m}$ is required in order to compute $\vec{u}_{m+1}$, but $\mathbb{B}\vec{u}_{m}=\mathbb{P}\,(\vec{u}_{m}\cdot\nabla)\vec{u}_{m}$ involves both differentiation and multiplication, it follows that a $\delta_{L_{2,0}^{\sigma}}$-name of $\vec{u}_{m}$ may not contain enough information for computing $\mathbb{B}\vec{u}_{m}$. Moreover, since estimates of type $\|\mathbb{A}^{\alpha}\vec{u}_{m}(t)\|_{2}$, $0\leq\alpha\leq 1$, play a key role in proving Propositions 4 and 5, we need to computationally derive a richer code for $\vec{u}_{m}$ from a given $\delta_{L^{\sigma}_{2,0}}$-name of $\vec{u}_{m}$ in order to capture the fact that $\vec{u}_{m}$ is in the domain of $\mathbb{A}^{\alpha}$ for $t>0$.

5.1 Representing and Operating on Space $H^{s}_{2,0}(\Omega)$

We begin by recalling several definitions and facts. Let $\theta_{n,m}(x,y):=e^{i(nx+my)\pi}$, $n,m\geq 0$. Then, the sequence $\{\theta_{n,m}(x,y)\}_{n,m\geq 0}$ is a computable orthogonal basis of $L_{2}(\Omega)$. For any $s\geq 0$, $H_{2}^{s}(\Omega)$ is the set of all (generalized) functions $w(x,y)$ on $\Omega$ satisfying $\sum_{n,m\geq 0}(1+n^{2}+m^{2})^{s}|w_{n,m}|^{2}<\infty$, where $w_{n,m}=\int_{-1}^{1}\int_{-1}^{1}w(x,y)\theta_{n,m}(x,y)\,dx\,dy$. $H_{2}^{s}(\Omega)$ is a Banach space with a norm $\|w\|_{H_{2}^{s}}=(\sum_{n,m\geq 0}(1+n^{2}+m^{2})^{s}|w_{n,m}|^{2})^{1/2}$.

Let $D(\mathbb{A}^{\alpha})$ be the domain of $\mathbb{A}^{\alpha}$. Since

[TABLE]

and $D(\mathbb{A}^{\alpha})$, $0\leq\alpha\leq 1$, are the complex interpolation spaces of $L^{\sigma}_{2,0}(\Omega)$ and $D(\mathbb{A})$, we need to represent the subspace of $H^{s}_{2}(\Omega)$ in which the functions vanish on $\partial\Omega$. However, it is usually difficult to design a coding system for such subspaces. Fortunately, for $0\leq s<3/2$, it is known classically that

[TABLE]

where $H^{s}_{2,0}(\Omega)$ is the closure in $H^{s}_{2}$-norm of the set of all $C^{\infty}$-smooth functions defined on compact subsets of $\Omega$. For $H^{s}_{2,0}(\Omega)$, there is a canonical coding system

[TABLE]

(see (11) and (12) for the definitions of $\gamma_{n}$ and $\gamma_{n}\ast q$). Then every $w$ in $H^{s}_{2,0}(\Omega)$ can be encoded by a sequence $\{p_{k}\}\subset\mathcal{H}$ such that $\|p_{k}-w\|_{H^{s}_{2}}\leq 2^{-k}$; the sequence $\{p_{k}\}$, which are mollified polynomials with rational coefficients, is called a $\delta_{H^{s}_{2,0}}$-name of $w$. If $w=(w_{1},w_{2})\in H^{s}_{2,0}(\Omega)\times H^{s}_{2,0}(\Omega)$, a $\delta_{H^{s}_{2,0}}$-name of $w$ is a sequences $\{(p_{k},q_{k})\}$, $p_{k},q_{k}\in\mathcal{H}$, such that $(\|w_{1}-p_{k}\|^{2}_{H^{s}_{2,0}}+\|w_{2}-q_{k}\|^{2}_{H^{s}_{2,0}})^{1/2}\leq 2^{-k}$.

Notation 5.1

*We use $\|w\|_{H^{s}_{2}}$ to denote the $H^{s}_{2}$-norm of $w$ if $w$ is in $H^{s}_{2}(\Omega)$ or $H^{s}_{2}\times H^{s}_{2}$-norm of $w$ if $w$ is in $H^{s}_{2}(\Omega)\times H^{s}_{2}(\Omega)$. Also for readability we use $[\rho\to\delta_{H_{2,0}^{s}}]$ to denote the canonical representation of either C\big{(}[0;T];H_{2,0}^{s}(\Omega)\big{)} or C\big{(}[0;T];H_{2,0}^{s}(\Omega)\times H_{2,0}^{s}(\Omega)\big{)}.

Recall that C\big{(}[0;T];H_{2,0}^{s}(\Omega)\big{)} is the set of all continuous functions from the interval $[0;T]$ to $H_{2,0}^{s}(\Omega)$. A function u\in C\big{(}[0;T];H_{2,0}^{s}(\Omega)\big{)} is computable if there is a machine that computes a $\delta_{H_{2,0}^{s}}$-name of $u(t)$ when given a $\rho$-name of $t$ as input; and a map F:X\to C\big{(}[0;T];H_{2,0}^{s}(\Omega)\big{)} from a represented space $(X,\delta_{X})$ to C\big{(}[0;T];H_{2,0}^{s}(\Omega)\big{)} is computable if there is a machine that computes a $\delta_{H_{2,0}^{s}}$-name of $F(x)(t)$ when given a $\delta_{X}$-name of $x$ and a $\rho$-name of $t$. Let $X$ be either $L_{2}(\Omega)$, $L^{\sigma}_{2,0}(\Omega)$, $H_{2,0}^{s}(\Omega)$, or C\big{(}[0;T];H_{2,0}^{s}(\Omega)\big{)}. We remark again that a $\delta_{X}$-name of $f\in X$ is simply an effective approximation of $f$ because each space $X$ is equipped with a norm.

Lemma 4

For $s\geq 1$, differentiation $\partial_{x},\partial_{y}:H_{2,0}^{s}(\Omega)\to L_{2}(\Omega)$ is $(\delta_{H_{2,0}^{s}},\delta_{L_{2}})$-computable.

Proof

Let $\{p_{k}\}$ be a $\delta_{H_{2,0}^{s}}$-name of $w\in H^{s}_{2,0}(\Omega)$. Since $\partial_{x}(\gamma\ast q)=\gamma\ast\partial_{x}q$, the map $p_{k}\mapsto\partial_{x}p_{k}$ is computable; hence a polynomial $\tilde{p}$ in $\mathbb{Q}[\mathbb{R}^{2}]$ can be computed from $p_{k}$ such that $\max_{-1\leq x,y\leq 1}|\partial_{x}p_{k}-\tilde{p}_{k}|<2^{-k}$. Next let us express $w$ and $p_{k}$ in the orthogonal basis $\theta_{n,m}$: $w(x,y)=\sum_{n,m\geq 0}w_{n,m}e^{in\pi x}e^{im\pi y}$ and $p_{k}(x,y)=\sum_{n,m\geq 0}p_{k,n,m}e^{in\pi x}e^{im\pi y}$, where

[TABLE]

Since $s\geq 1$ and $\{p_{k}\}$ is a $\delta_{H_{2,0}^{s}}$-name of $w$, it follows that

[TABLE]

which further implies that

[TABLE]

Thus, by definition, $\{\tilde{p}_{k}\}$ is a $\delta_{L_{2}}$-name of $\partial_{x}w$.

It is known classically that every polygonal domain in $\mathbb{R}^{2}$ is Lipschitz (see, for example, [9]) and $H_{2}^{s}(U)$ is continuously embedded in $C(\overline{U})$ if $s>1$ and $U$ is a bounded Lipschitz domain, where $\overline{U}$ is the closure of $U$ in $\mathbb{R}^{2}$ and $C(\overline{U})$ is the set of all continuous functions on $\overline{U}$. Since $\Omega$ is a bounded polygonal domain, it follows that for any $s>1$, there is a constant $C_{s}>0$ such that $\|w\|_{C(\overline{\Omega})}\leq C_{s}\|w\|_{H_{2}^{s}(\Omega)}$, where $\|w\|_{C(\overline{\Omega})}=\|w\|_{\infty}=\max\{|w(x,y)|\,:\,(x,y)\in\overline{U}\}$.

Lemma 5

For $s>1$, multiplication $Mul:H_{2}^{s}(\Omega)\times L_{2}(\Omega)\to L_{2}(\Omega)$, $(v,w)\mapsto vw$, is $(\delta_{H_{2,0}^{s}}\times\delta_{L_{2}},\delta_{L_{2}})$-computable.

Proof

Assume that $\{p_{k}\}$ is a $\delta_{H_{2,0}^{s}}$-name of $v$ and $\{q_{k}\}$ is a $\delta_{L_{2}}$-name of $w$. For each $n\in\mathbb{N}$, pick $k(n)\in\mathbb{N}$ such that $C_{s}\|v\|_{H_{2}^{s}}\|w-q_{k(n)}\|_{2}\leq 2^{-(n+1)}$. Since $\|v\|_{H_{2}^{s}}$ is computable from $\{p_{k}\}$, the function $n\mapsto k(n)$ is computable from $\{p_{k}\}$ and $\{q_{k}\}$. Next pick $m(n)\in\mathbb{N}$ such that $\|q_{k(n)}\|_{C(\overline{\Omega})}\|v-p_{m(n)}\|_{H_{2}^{s}}\leq 2^{-(n+1)}$. It is clear that $m(n)$ is computable from $k(n)$, $\{q_{k}\}$, and $\{p_{k}\}$. The sequence $\{p_{m(n)}q_{k(n)}\}_{n}$ is then a $\delta_{L_{2}}$-name of $vw$, for it is a sequence of polynomials of rational coefficients and $\|vw-p_{m(n)}q_{k(n)}\|_{2}\leq\|v\|_{C(\overline{\Omega})}\|w-q_{k(n)}\|_{2}+\|q_{k(n)}\|_{C(\overline{\Omega})}\|v-p_{m(n)}\|_{H_{2}^{s}}\leq 2^{-n}$.

5.2 Some classical properties of fractional powers of $\mathbb{A}$

It is known that fractional powers of the Stokes operator $\mathbb{A}$ are well defined; cf. [10, Section 2.6]. In the following, we summarize some classical properties of the Stokes operator and its fractional powers; these properties will be used in later proofs.

Fact 5.2

Let $\mathbb{A}$ be the Stokes operator.

• (1)

For every $0\leq\alpha\leq 1$, let $D(\mathbb{A}^{\alpha})$ be the domain of $\mathbb{A}^{\alpha}$; this is a Banach space with the norm $\|\vec{u}\|_{D(\mathbb{A}^{\alpha})}:=\|\mathbb{A}^{\alpha}\vec{u}\|_{L^{\sigma}_{2,0}(\Omega)}=\|\mathbb{A}^{\alpha}\vec{u}\|_{2}$. In particular, $D(\mathbb{A}^{\alpha})$ is continuously embedded in $H_{2}^{2\alpha}$, that is, for every $\vec{u}\in D(\mathbb{A}^{\alpha})$,

[TABLE]

where $C$ is a constant independent of $\alpha$. Moreover, we have D(\mathbb{A}^{1/2})=L^{\sigma}_{2,0}(\Omega)\bigcap\{\vec{u}\in(H^{1}(\Omega))^{2};\mbox{\vec{u}=0$on$\partial\Omega}\}.

• (2)

For every nonnegative $\alpha$ the estimate

[TABLE]

is valid for all $\vec{u}\in L^{\sigma}_{2,0}(\Omega)$, where $C_{\alpha}$ is a constant depending only on $\alpha$. In particular, $C_{0}=1$. Moreover, the estimate implies implicitly that for every $\vec{u}\in L^{\sigma}_{2,0}(\Omega)$, $e^{-t\mathbb{A}}\vec{u}$ is in the domain of $\mathbb{A}$, and thus $e^{-t\mathbb{A}}\vec{u}$ vanishes on the boundary of $\Omega$ for $t>0$.

• (3)

If $\alpha\geq\beta>0$, then $D(\mathbb{A}^{\alpha})\subseteq D(\mathbb{A}^{\beta})$.

• (4)

For $0<\alpha<1$, if $\vec{u}\in D(\mathbb{A})$, then

[TABLE]

• (5)

$\|\mathbb{A}^{-1/4}\mathbb{P}(\vec{u},\nabla)\vec{v}\|_{2}\leq M\|\mathbb{A}^{1/4}\vec{u}\|_{2}\|\mathbb{A}^{1/2}\vec{v}\|_{2}$* is valid for all $\vec{u},\vec{v}$ in the domain of $\mathbb{A}^{3/5}$, where $M$ is a constant independent of $\vec{u}$ and $\vec{v}$.*

Proof

See Lemmas 2.1, 2.2 and 2.3 in [3] for (1) and (2) except for $C_{0}=1$; $C_{0}=1$ is proved in Lemma 3. See Theorems 6.8 and 6.9 in Section 2.6 of [10] for (3) and (4); Lemma 3.2 in [3] for (5).

We record, without going into the details, that the constants $C$, $M$, and $C_{\alpha}$ ($0\leq\alpha\leq 1$) appeared in Fact 5.2 are in fact computable (some general discussions on the computability of Sobolev embedding constants and interpolation constants together with other constants in the PDE theory are forthcoming).

5.3 Proof of Proposition 4

In order to show that the iteration sequence is effectively convergent, we need to establish several estimates on various functions such as $\|\mathbb{A}^{\beta}u_{m}(t)\|_{2}$ and $\|\mathbb{A}^{\beta}(u_{m+1}(t)-u_{m}(t))\|_{2}$ for $\beta$ being some positive numbers. Subsequently, as a prerequisite, $u_{m}(t)$ must be in the domain of $\mathbb{A}^{\beta}$; thus the functions $u_{m}(t)$ are required to have higher smoothness than the given initial function $\vec{a}$ according to Fact 5.2-(1). This is indeed the case: For functions $u_{m}(t)$ obtained by the iteration (23), it is known classically that if $u_{m}(0)\in L_{2}(\Omega)$ then $u_{m}(t)\in H^{2\alpha}_{2}(\Omega)$ for $t>0$, where $0\leq\alpha\leq 1$. In other words, $u_{m}(t)$ undergoes a jump in smoothness from $t=0$ to $t>0$ (due to the integration). In the following lemma, we present an algorithmic version of this increase in smoothness.

Lemma 6

Let $\alpha=3/5$. Then for the iteration  (23)

[TABLE]

the mapping $\mathbb{S}_{H}:\mathbb{N}\times L^{\sigma}_{2,0}(\Omega)\times(0,\infty)\to H_{2,0}^{2\alpha}(\Omega)\times H_{2,0}^{2\alpha}(\Omega)$, $(m,\vec{a},t)\mapsto\vec{u}_{m}(t)$, is well-defined and $(\nu\times\delta_{L_{2,0}^{\sigma}}\times\rho,\delta_{H_{2,0}^{2\alpha}})$-computable.

We emphasize that the lemma holds true for $t>0$ only. Also the choice of $\alpha=3/5$ is somewhat arbitrary; in fact, $\alpha$ can be selected to be any rational number strictly between $\frac{1}{2}$ and $\frac{3}{4}$. The requirement $\alpha<\frac{3}{4}$ guarantees that $D(\mathbb{A}^{\alpha})\subset H^{2\alpha}_{2,0}(\Omega)\times H^{2\alpha}_{2,0}(\Omega)$ because $2\alpha<3/2$ (see (24)). The other condition $\alpha>\frac{1}{2}$ ensures that Lemma 5 can be applied for $2\alpha>1$.

Proof

We induct on $m$. Note that for any $t>0$ and any $a\in L^{\sigma}_{2,0}(\Omega)$, the estimates (25) and (26) imply that

[TABLE]

Combining this inequality with the following strengthened version of (19): for any $a\in\mathcal{P}$,

[TABLE]

(the inequality is valid since $a$ is $C^{\infty}$), a similar argument used to prove Proposition 3 works for $m=0$. Moreover, by type conversion (Fact 1.2), $a\in L^{\sigma}_{2,0}(\Omega)\mapsto u_{0}\in C((0,\infty),H^{6/5}_{2,0}(\Omega)\times H^{6/5}_{2,0}(\Omega))$ is $(\delta_{L^{\sigma}_{2,0}},[\rho\to\delta_{H^{6/5}_{2,0}}])$-computable.

Assume that $(j,a)\mapsto u_{j}$ is $(\nu,\delta_{L^{\sigma}_{2,0}},[\rho\to\delta_{H^{6/5}_{2,0}}])$-computable for $0\leq j\leq m$, where $a\in L^{\sigma}_{2,0}(\Omega)$, and $u_{j}\in C((0,\infty),(H^{6/5}_{2,0}(\Omega))^{2})$. We show how to compute a $\delta_{H_{2,0}^{6/5}}$-name for $u_{m+1}(t)=e^{-t\mathbb{A}}a-\int^{t}_{0}e^{-(t-s)\mathbb{A}}\mathbb{B}u_{m}(s)ds$ on inputs $m+1$, $a$ and $t>0$. Let us first look into the nonlinear term $\mathbb{B}u_{m}$. It is clear that $\mathbb{B}u_{m}(s)$ lies in $L^{\sigma}_{2,0}(\Omega)$ for $s>0$. Moreover, it follows from Lemmas 4 and 5, and Proposition 2 that the map $(u_{m},s)\mapsto\mathbb{B}u_{m}(s)$ is $([\rho\!\to\!\delta_{H_{2,0}^{2\alpha}}],\rho,\delta_{L_{2,0}^{\sigma}})$-computable for all $s\in(0,t]$. Now since $\mathbb{B}u_{m}(s)$ is in $L^{\sigma}_{2,0}(\Omega)$ for $s>0$, it follows from the case where $m=0$ that $(u_{m},s)\mapsto e^{-(t-s)\mathbb{A}}\mathbb{B}u_{m}(s)$ is $([\rho\!\to\!\delta_{H_{2,0}^{2\alpha}}],\rho,\delta_{H_{2,0}^{6/5}})$-computable for $0<s<t$.

Next let us consider the integral $\int^{t}_{0}e^{-(t-s)\mathbb{A}}\mathbb{B}\vec{u}_{m}(s)\,ds$; we wish to compute a $\delta_{H^{6/5}_{2,0}}$-name of the integral from $a$ and $t>0$. We make use of the following fact: For $\theta\geq 1$, the integration operator from $C([a,b];H^{\theta}_{2,0}(\Omega)\times H^{\theta}_{2,0}(\Omega))$ to $H^{\theta}_{2,0}(\Omega)\times H^{\theta}_{2,0}(\Omega)$, $F\mapsto\int^{b}_{a}F(t)(x)dt$, is computable from $a$, $b$, and $F$. This fact can be proved by a similar argument as the one used in the proof of Lemma 3.7 [24]. However, since the function $e^{-(t-s)\mathbb{A}}\mathbb{B}u_{m}(s)$ is not necessarily in $(H^{6/5}_{2}(\Omega))^{2}$ when $s=0$ or $s=t$, the stated fact cannot be directly applied to the given integral. To overcome the problem of possible singularities at the two endpoints, we use a sequence of closed subintervals $[t_{n},t-t_{n}]$ to approximate the open interval $(0,t)$, where $t_{n}=t/2^{n}$, $n\geq 1$. Then it follows from the stated fact and the induction hypotheses that a $\delta_{H^{6/5}_{2,0}}$-name, say $\{p_{n,K}\}$, of $u^{n}_{m+1}(t)=e^{-t\mathbb{A}}a-\int^{t-t_{n}}_{t_{n}}e^{-(t-s)\mathbb{A}}\mathbb{B}u_{m}(s)ds$ can be computed from inputs $n$, $u_{m}$, and $t>0$, which satisfies the condition that $\|u^{n}_{m+1}(t)-p_{n,K}\|_{H^{6/5}_{2}}\leq 2^{-K}$. Thus if we can show that the integrals $\int_{0}^{t_{n}}e^{-(t-s)\mathbb{A}}\mathbb{B}u_{m}(s)ds$ and $\int_{t-t_{n}}^{t}e^{-(t-s)\mathbb{A}}\mathbb{B}u_{m}(s)ds$ tend to zero effectively in $H_{2}^{6/5}\times H_{2}^{6/5}$-norm as $n\to\infty$, then we can effectively construct a $\delta_{H_{2,0}^{6/5}}$-name of $u_{m+1}(t)$ from $\{p_{n,K}\}_{n,K}$.

It remains to show that both sequences of integrals tend to zero effectively in $H_{2}^{6/5}\times H_{2}^{6/5}$-norm as $n\to\infty$. Since a similar argument works for both sequences, it suffices to show that the sequence $\mbox{Int}_{n}:=\int_{0}^{t_{n}}e^{-(t-s)\mathbb{A}}\mathbb{B}u_{m}(s)ds$ tends to zero effectively as $n\to\infty$. We are to make use of Fact 5.2-(1), (2), (5) for showing the effective convergence. The following two claims comprise the proof.

Claim I. Let $\beta=\frac{1}{2}$ or $\frac{1}{4}$. Then the map $(a,t,m,\beta)\mapsto M_{\beta,m}$ is computable, where $M_{\beta,m}$ is a positive number satisfying the condition

[TABLE]

(note that $M_{\beta,m}$ is independent of $s$).

Proof. Again we induct on $m$. For $m=0$, let $M_{\beta,0}=C_{\beta}\|a\|_{2}$, where $C_{\beta}$ is the constant in estimate (26) with $\alpha$ replaced by $\beta$ and $u$ by $a$. Then $M_{\beta,0}$ is computable from $a$ and $\beta$, and $\|\mathbb{A}^{\beta}u_{0}(s)\|_{2}\leq C_{\beta}s^{-\beta}\|a\|_{2}=M_{\beta,0}s^{-\beta}$ for any $s>0$. Assume that $M_{\beta,k}$, $0\leq k\leq m$, has been computed from $k,\beta,a$, and $t>0$. We show how to compute $M_{\beta,m+1}$. Since $u_{m+1}(s)$ has a singularity at $s=0$, it may not be in $H^{2\beta}_{2}(\Omega)\times H^{2\beta}_{2}(\Omega)$ at $s=0$ (recall that D(\mathbb{A}^{1/2})=L^{\sigma}_{2,0}(\Omega)\bigcap\{\vec{u}\in H^{1}_{2}(\Omega)\times H^{1}_{2}(\Omega):\mbox{\vec{u}=0$on$\partial\Omega}\}). Let us first compute a bound (in $L_{2}$-norm) for $\mathbb{A}^{\beta}\int^{s}_{\epsilon}e^{-(t-r)\mathbb{A}}\mathbb{B}u_{m}(r)dr$, where $0<\epsilon<s$. It follows from the induction hypothesis, Fact 5.2-(1), (2), (5), and Theorems 6.8 and 6.13 in [10] that

[TABLE]

Subsequently, we obtain that

[TABLE]

where $B(\frac{3}{4}-\beta,1/4)$ is the integral $\int_{0}^{1}(1-\theta)^{(\frac{3}{4}-\beta)-1}\theta^{\frac{1}{4}-1}d\theta$, which is the value of the Beta function $B(x,y)=\int^{1}_{0}(1-\theta)^{1-x}\theta^{1-y}d\theta$ at $x=\frac{3}{4}-\beta$ and $y=\frac{1}{4}$. It is clear that $B(\frac{3}{4}-\beta,1/4)$ is computable. Thus if we set

[TABLE]

then $M_{\beta,m+1}$ is computable and satisfies the condition that $\|\mathbb{A}^{\beta}u_{m+1}(s)\|_{2}\leq M_{\beta,m+1}s^{-\beta}$ for all $0<s<t$. The proof of Claim I is complete.

Claim II. $\left\|\int^{t_{n}}_{0}e^{-(t-s)\mathbb{A}}\mathbb{B}u_{m}(s)ds\right\|_{H_{2}^{6/5}}\to 0$ effectively as $n\to\infty$.

Proof. Once again, to avoid singularity of $u_{m}(s)$ at $s=0$, we begin with the following estimate: Let $0<\epsilon<t_{n}$. Then it follows from Fact 5.2-(1), (2), (5), (27), (30), and a similar calculation as performed in Claim I that

[TABLE]

which then implies that

[TABLE]

It is readily seen that $\int^{t_{n}}_{0}e^{-(t-s)\mathbb{A}}\mathbb{B}u_{m}(s)ds\|_{H_{2}^{6/5}}\to 0$ effectively as $n\to\infty$ (recall that $t_{n}=t/2^{n}$). The proof for the claim II, and thus for the lemma is now complete.

Remark 1

In our effort to compute an upper bound for $\|\mathbb{A}^{\beta}u_{m+1}(s)\|_{2}$, we start with the integral $\int^{s}_{\epsilon}e^{-(s-r)\mathbb{A}}\mathbb{B}u_{m}(r)dr$ because the integral might have a singularity at [math]; then we take the limit as $\epsilon\to 0$ to get the desired estimate (see computations of (5.3) and (5.3)). The limit exists because the bound, $C_{\beta+\frac{1}{4}}MM_{\frac{1}{4},m}M_{\frac{1}{2},m}B\left(\frac{3}{4}-\beta,\frac{1}{4}\right)$, is uniform in $r$ for $0<r\leq s$. In the rest of the paper, we will encounter several similar computations. In those later situations, we will derive the estimates starting with $\int^{t}_{0}$ instead of $\int^{t}_{\epsilon}$. There will be no loss in rigor because the integral is uniformly bounded with respect to the integrating variable, say $t$, for $t>0$.

Corollary 1

For any $\vec{a}\in L^{\sigma}_{2,0}(\Omega)$ and $t>0$, let $\{u_{m}(t)\}$ be the sequence generated by the iteration scheme (23) based on $\vec{a}$. Then $u_{m}(t)\in Dom(\mathbb{A}^{3/5})\subset Dom(\mathbb{A}^{1/2})\subset Dom(\mathbb{A}^{1/4})$.

Proof

The corollary follows from Lemma 6 and Fact 5.2-(3).

Corollary 2

The map from $\mathcal{P}$ to $L_{2}(\Omega)$, $\vec{u}\mapsto\|\mathbb{A}^{\alpha}\vec{u}\|_{2}$, is $(\delta_{L_{2,0}^{\sigma}},\rho)$-computable, where $\alpha=1/8,1/4$, or $1/2$.

Proof

We prove the case when $\alpha=1/4$; the other two cases can be proved in exactly the same way. Since $\mathcal{P}$ is contained in the domain of $\mathbb{A}$, it follows from Theorem 6.9, Section 2.6 [10] that for every $\vec{u}\in\mathcal{P}$, $\mathbb{A}^{1/4}\vec{u}=\frac{\sin\pi/4}{\pi}\int^{\infty}_{0}t^{-3/4}\mathbb{A}(t\mathbb{I}+\mathbb{A})^{-1}\vec{u}dt$. By definition of $\mathcal{P}$, if $\vec{u}\in\mathcal{P}$, then $\vec{u}$ is $C^{\infty}$ with compact support in $\Omega$, and $\mathbb{A}\vec{u}=-\mathbb{P}\bigtriangleup\vec{u}=-\bigtriangleup\vec{u}$. Express each component of $\vec{u}=(u^{1},u^{2})$ in terms of the orthogonal basis $\{e^{in\pi x}e^{im\pi y}\}_{n,m}$ of $L_{2}(\Omega)$ in the form of $u^{i}=\sum_{n,m\geq 0}u^{i}_{n,m}e^{i\pi nx}e^{i\pi my}$, where $u^{i}_{n,m}=\int^{1}_{-1}\int_{-1}^{1}u^{1}(x,y)e^{i\pi nx}e^{i\pi my}dxdy$. Then a straightforward calculation shows that

[TABLE]

Since the integral is convergent and computable, it follows that $\mathbb{A}^{1/4}\vec{u}$ is computable from $\vec{u}$ and, consequently, $\|\mathbb{A}\vec{u}\|_{2}$ is computable.

Proof (Proof of Proposition 4)

For each $\vec{a}\in L^{\sigma}_{2,0}$, let $\{\vec{a}_{k}\}$, $\vec{a}_{k}=(a^{1}_{k},a^{2}_{k})\in\mathcal{P}$, be a $\delta_{L_{2,0}^{\sigma}}$-name of $\vec{a}$; i.e. $\|\vec{a}-\vec{a}_{k}\|_{2}\leq 2^{-k}$. Let $\widetilde{C}:=c_{1}MB_{1}$, where $M$ is the constant in Fact 5.2(4), $c_{1}=\max\{C_{1/4},C_{1/2},C_{3/4},1\}$, and

[TABLE]

with $B(a,b)=\int_{0}^{1}(1-t)^{a-1}t^{b-1}dt$, $a,b>0$, being the beta function. Then $M$ and $c_{1}$ are computable by assumption while $B_{1}$ is computable for the beta functions with rational parameters are computable. Note that $c_{1}B_{1}\geq 1$. Let

[TABLE]

Our goal is to compute a constant $\epsilon$, $0<\epsilon<1$, such that near $t=0$,

[TABLE]

where $L$ is a constant. Once this is accomplished, the proof is complete.

It follows from Corollary 1 that Fact 5.2-(5) holds true for all $u_{m}(t)$ and $v_{m}(t)$ with $t>0$. It is also known classically that

[TABLE]

(see, for example,[5]). The equality in the above estimate holds true because $\mathbb{A}^{-1/4}$ is a (bounded) linear operator. The estimate (33) indicates that, in order to achieve (32), there is a need in establishing some bounds on $\|\mathbb{A}^{\beta}u_{m}(t)\|_{2}$ and $\|\mathbb{A}^{\beta}v_{m}(t)\|_{2}$ which become ever smaller as $m$ gets larger uniformly for values of $t$ near zero. The desired estimates are developed in a series of claims beginning with the following one.

Claim 1. Let $\beta=\frac{1}{4}$ or $\frac{1}{2}$; let

[TABLE]

and

[TABLE]

Then there is a computable map from $L^{\sigma}_{2,0}(\Omega)$ to $(0,1)$, $\vec{a}\mapsto T_{\vec{a}}$, such that

[TABLE]

Proof. First we note that $t^{\beta}\|\mathbb{A}^{\beta}e^{-t\mathbb{A}}\vec{a}\|_{2}=0$ for any $\vec{a}\in L^{\sigma}_{2,0}(\Omega)$ if $t=0$; cf. Theorem 6.1 in [3]. Furthermore, it follows from (17) that the operator norm of $e^{-t\mathbb{A}}$, $\|e^{-t\mathbb{A}}\|_{op}$, is bound above by 1 for any $t>0$. Since $e^{-t\mathbb{A}}$ is the identity map on $L^{\sigma}_{2,0}(\Omega)$ when $t=0$, we conclude that $\max_{0\leq t\leq T}\|e^{-t\mathbb{A}}\|_{op}\leq 1$ for any $T>0$. Now for any $\vec{a}\in L^{\sigma}_{2,0}(\Omega)$, it follows from Fact  5.2-(2) and Theorems 6.8 and 6.13 in Section 2.6 of [10] ($\mathbb{A}^{\alpha}$ and $e^{-t\mathbb{A}}$ are interchangeable) that

[TABLE]

We note that although $\vec{a}$ is not necessarily in the domain of $\mathbb{A}$ but $\vec{a}_{k}\in\mathcal{P}$ and $\mathcal{P}$ is contained in the domain of $\mathbb{A}$; thus $\mathbb{A}^{\beta}\vec{a}_{k}$ is well defined. Furthermore, it follows from Corollary 2 that $\|\mathbb{A}^{\beta}\vec{a}_{k}\|_{2}$ is computable. Clearly one can compute a positive integer $\hat{k}$ such that

[TABLE]

then compute a positive number $T_{\vec{a}}$ such that

[TABLE]

The computations are performed on the inputs $\vec{a}$ and the constants $c_{1}$, $M$, and $B_{1}$. Consequently, $k^{\vec{a}}_{0}(T_{\vec{a}})<1/(8\widetilde{C})$. The proof of Claim 1 is complete.

We recall that, for a given $\vec{a}\in L^{\sigma}_{2,0}(\Omega)$, the iteration scheme (23) is based on the “seed” function $u_{0}(t)=e^{-t\mathbb{A}}\vec{a}$. Claim 1 asserts that the seed function has the property that $\max_{0\leq t\leq T_{\vec{a}}}t^{\beta}\|\mathbb{A}^{\beta}u_{0}(t)\|_{2}$ is bounded by $\tilde{K}^{\vec{a}}_{\beta,0}$, uniformly in $t$. We extend this property to the iteration sequence $\{u_{m}(t)\}$ in the next claim.

Claim 2. Let $\beta=\frac{1}{4}$ or $\frac{1}{2}$. Then there is a computable map $\mathbb{N}\times L^{\sigma}_{2,0}\to(0,\infty)$, $(m,\vec{a})\mapsto K^{\vec{a}}_{\beta,m}$, such that

[TABLE]

Proof. We induct on $m$. For $m=0$, let $K^{\vec{a}}_{\beta,0}=1/(8\widetilde{C})$. Then (34) follows from Claim 1. It is clear that $K^{\vec{a}}_{\beta,0}$ is computable.

For $m\geq 1$ and $t>0$, $K^{\vec{a}}_{\beta,m+1}$ is computed by the recursive formula:

[TABLE]

The recursive formula is derived similarly as that of (5.3). Since the upper bound is uniformly valid for all $0<t\leq T_{\vec{a}}$, it follows that it is also valid for $t=0$. The proof of Claim 2 is complete.

In the next claim, we show that the sequences $\{K^{\vec{a}}_{\beta,m}\}$, $\beta=1/4$ or $1/2$, are bounded above with an upper bound strictly less than $1/(2\widetilde{C})$.

Claim 3. Let $k^{\vec{a}}_{m}=\max\{K^{\vec{a}}_{\frac{1}{4},m},K^{\vec{a}}_{\frac{1}{2},m}\}$ and let $K=\frac{4k^{\vec{a}}_{0}(\sqrt{2}-1)}{\sqrt{2}}$. Then $k^{\vec{a}}_{m}\leq K<\frac{1}{2\widetilde{C}}$ for all $m\geq 1$.

Proof. It follows from Claim 2 that $k^{\vec{a}}_{0}=\frac{1}{8\widetilde{C}}$ and $k^{\vec{a}}_{m+1}\leq k^{\vec{a}}_{0}+\widetilde{C}(k^{\vec{a}}_{m})^{2}$ (recall that $\widetilde{C}=c_{1}MB_{1}$). To get a bound on $k^{\vec{a}}_{m}$, let’s write $k^{\vec{a}}_{m}=k^{\vec{a}}_{0}w_{m}$. Then $w_{m}$ satisfies the following inequality:

[TABLE]

which implies that

[TABLE]

Then a direct calculation shows that

[TABLE]

Thus

[TABLE]

And so if we pick $K=\frac{4k^{\vec{a}}_{0}(\sqrt{2}-1)}{\sqrt{2}}$, then $k^{\vec{a}}_{m}\leq K<\frac{1}{2\widetilde{C}}$ for all $m\geq 1$. The proof of Claim 3 is complete.

Next we present an upper bound for $t^{\alpha}\|\mathbb{A}^{\alpha}v_{m}(t)\|_{2}$, $m\geq 1$. Recall that $v_{m}(t)=u_{m+1}(t)-u_{m}(t)$.

Claim 4. For $t\in[0,T_{\vec{a}}]$, $0\leq\alpha<\frac{3}{4}$, and $m\geq 1$,

[TABLE]

Proof. First we observe that (36) is true for $t=0$. Next we assume that $0<t\leq T_{\vec{a}}$. Once again we induct on $m$. At $m=1$: We recall from the definition of $c_{1}$ and $B_{1}$ that $\frac{1}{2c_{1}B_{1}}\leq\frac{1}{2}$. Also it follows from (33), Claims 2 and 3 that $\|\mathbb{A}^{\frac{1}{2}}u_{1}(t)\|_{2}\leq K^{\vec{a}}_{\frac{1}{2},1}t^{-\frac{1}{2}}\leq Kt^{-\frac{1}{2}}$, $\|\mathbb{A}^{\frac{1}{4}}u_{0}(t)\|_{2}\leq Kt^{-\frac{1}{4}}$, $\|\mathbb{A}^{\frac{1}{4}}v_{0}(t)\|_{2}\leq 2Kt^{-\frac{1}{4}}$, and $\|\mathbb{A}^{\frac{1}{2}}v_{0}(t)\|_{2}\leq 2Kt^{-\frac{1}{2}}$. Making use of these inequalities we obtain the following estimate:

[TABLE]

Thus (36) is true for $m=1$.

Now assuming that (36) holds for all $1\leq j\leq m$, we show that (36) is also true for $m+1$. First it follows from Claims 2 and 3, and the induction hypothesis that for any $s\in(0,T_{\vec{a}})$,

[TABLE]

Similarly,

[TABLE]

Thus,

[TABLE]

These inequalities imply the desired estimate:

[TABLE]

The proof for Claim 4 is complete.

We now set $\alpha=0$, $\epsilon=2\widetilde{C}K$, and $L=2KC_{\frac{1}{4}}B\left(\frac{3}{4},\frac{1}{4}\right)$. Since $K<\frac{1}{2\widetilde{C}}$ by Claim 3, it follows that $0<\epsilon<1$ and

[TABLE]

Consequently, the iterated sequence $\{u_{m}(t)\}$ converges effectively to $u(t)$ and uniformly on $[0,T_{\vec{a}}]$.

We mention in passing the following fact that can be proved by similar computations of Claims 1 - 3: On input $(\vec{a},m,n)$, a positive number $T(\vec{a},m,n)$ can be computed such that $k^{\vec{a}}_{0}(T(\vec{a},m,n))<(8\tilde{C})^{-1}\cdot 2^{-n}$, $T(\vec{a},m,n+1)<T(\vec{a},m,n)$, and $\max_{0\leq t\leq T(\vec{a},m,n)}t^{\beta}\|\mathbb{A}^{\beta}u_{m}(t)\|_{2}\leq L^{\vec{a}}_{\beta,m}\cdot 2^{-n}$, where $L^{\vec{a}}_{\beta,m}$ is a constant independent of $t$ and $n$, and computable from $\vec{a}$ and $m$.

5.4 Proof of Proposition 5

We now come to the proof of Proposition 5. We need to show that the map $\mathbb{S}:\mathbb{N}\times L^{\sigma}_{2,0}\times[0,\infty)\to L^{\sigma}_{2,0}$, $(m,\vec{a},t)\mapsto\vec{u}_{m}(t)$, is $(\nu\times\delta_{L_{2,0}^{\sigma}}\times\rho,\delta_{L_{2,0}^{\sigma}})$-computable. By a similar argument as we used for proving Lemma 6, we are able to compute $\vec{u}_{m}(t)$ on the input $(m,\vec{a},t)$, where $m\in\mathbb{N}$, $\vec{a}\in L^{\sigma}_{2,0}(\Omega)$, and $t>0$. We note that $\vec{u}_{m}(0)=\vec{u}_{0}(0)=\vec{a}$ for all $m\in\mathbb{N}$. Thus, to complete the proof, it suffices to show that there is a modulus function $\eta:\mathbb{N}\times\mathbb{N}\to\mathbb{N}$, computable from $\vec{a}$, such that $\|\vec{u}_{m+1}(t)-\vec{a}\|_{2}\leq 2^{-k}$ whenever $0<t<2^{-\eta(m+1,k)}$. Now for the details. Given $\vec{a}$ and $k$. Refereeing to the last paragraph of the previous subsection and Fact  5.2-(2), (5), we obtain the following estimate: for $0<t<T(\vec{a},m,n)$

[TABLE]

Thus if $\|e^{-t\mathbb{A}}\vec{a}-\vec{a}\|_{2}\leq 2^{-(k+1)}$ and

[TABLE]

then

[TABLE]

Since $e^{-t\mathbb{A}}\vec{a}$ is computable in $t$ by Proposition 3 and $\vec{a}=e^{-0\mathbb{A}}\vec{a}$, there is a computable function $\theta_{1}:\mathbb{N}\to\mathbb{N}$ such that $\|e^{-t\mathbb{A}}\vec{a}-\vec{a}\|_{2}\leq 2^{-(k+1)}$ whenever $0<t<2^{-\theta_{1}(k)}$. Let $\theta_{2}:\mathbb{N}\times\mathbb{N}\to\mathbb{N}$ be a computable function satisfying $C_{1/4}ML^{\vec{a}}_{1/4,m}L^{\vec{a}}_{1/2,m}B(3/4,1/4)\cdot 2^{-2\theta_{2}(m,k)}\leq 2^{-(k+1)}$. Let $\eta(m+1,k)$ be a positive integer such that $2^{-\eta(m+1,k)}\leq\min\{2^{-\theta_{1}(k)},T(\vec{a},m,\theta_{2}(m,k))\}$. Then $\eta$ is the desired modulus function. The proof of Proposition 5 is complete.

Propositions 4 and 5 show that the solution $\vec{u}$ of the integral equation (4) is an effective limit of the computable iterated sequence $\{{\vec{u}}_{m}\}$ starting with ${\vec{u}}_{0}=\vec{a}$ on $[0,T_{\vec{a}}]$; consequently, $\vec{u}$ itself is also computable. Thus we obtain the desired preliminary result:

Theorem 5.3

There is a computable map $T:L^{\sigma}_{2,0}(\Omega)\to(0,\infty)$, $\vec{a}\mapsto T(\vec{a})$, such that $\vec{u}(t)$, the solution of the integral equation (4), is computable in $L^{\sigma}_{2,0}$ from $\vec{a}$ and $t$ for $\vec{a}\in L^{\sigma}_{2,0}$ and $t\in[0;T(\vec{a})]$.

5.5 The Inhomogeneous Case and Pressure

It is known [5, Theorem 2.3] that, also in the presence of an inhomogeneity \vec{g}\in C\big{(}[0;T],L^{\sigma}_{2,0}(\Omega)\big{)}, the iterate sequence (5) converges to a unique solution $\vec{u}$ of Equation (2) near $t=0$. Similarly to (the proofs of) Propositions 5, 4, and [24, Lemma 3.7], this solution is seen to be computable. Moreover, $\vec{g}=\mathbb{P}\vec{f}$ is computable from \vec{f}\in\big{(}L_{2}(\Omega)\big{)}^{2} according to Proposition 2. Finally the right-hand side of Equation (6) equals

[TABLE]

which, by the definition of $\mathbb{P}$ projecting onto the solenoidal subspace, is conservative (=rotation-free/a pure divergence). Hence the path integral $\int_{\vec{0}}^{\vec{x}}\vec{h}(\vec{y})\cdot d\vec{\gamma}(\vec{y})$ does not depend on the chosen path from $\vec{0}$ to $\vec{x}$ and well-defines $P(\vec{x})$. This concludes our proof of Theorem 1.1.

Appendix 0.A Proof of Proposition 1

(a)  For a divergence-free and boundary-free polynomial, its coefficients must satisfy a system of linear equations. In the following, we derive explicitly this system of linear equations in the 2-dimensional case, i.e. $\Omega=(-1,1)^{2}$. Let \vec{p}=(p_{1},p_{2})=\big{(}\sum_{i,j=0}^{N}a^{1}_{i,j}x^{i}y^{j},\sum_{i,j=0}^{N}a^{2}_{i,j}x^{i}y^{j}\big{)} be a divergence-free and boundary-free polynomial of real coefficients. (If the degree of $p_{1}$ or $p_{2}$ is less than $N$, then zeros are placed for the coefficients of missing terms.) Then, by definition,

[TABLE]

which implies that all coefficients in $\nabla\cdot\vec{p}$ must be zero; or equivalently, Equation (a)) holds true. Turning to the boundary conditions, along the line $x=1$, since

[TABLE]

is identically zero, it follows that $\sum_{i=0}^{N}a^{1}_{i,j}=\sum_{i=0}^{N}a^{2}_{i,j}=0$ for $0\leq j\leq N$. There are similar types of restrictions on the coefficients of $\vec{p}$ along the lines $x=-1$, $y=1$, and $y=-1$. In summary, $\vec{p}$ vanishes on $\partial\Omega$ if and only if for all $0\leq j,i\leq N$, both (8) and (9) hold true.

In the 3-dimensional case, a similar calculation shows that a polynomial triple \vec{p}(x,y,z)=\big{(}p_{1}(x,y,z),p_{2}(x,y,z),p_{3}(x,y,z)\big{)} is divergence-free and boundary-free if and only if its coefficients satisfies a system of linear equations with integer coefficients.

(b)  In [8] it is shown that for any real number $s\geq 3$ and for any function $\vec{w}\in\mathcal{N}_{div}^{s}\cap H_{2,0}^{1,\sigma}(\Omega)^{d}$, the following holds:

[TABLE]

where $\Omega=(-1,1)^{d}$,

[TABLE]

$\mathcal{P}_{N}$ is the set of all $d$-tuples of real polynomials with $d$ variables and degree less than or equal to $N$ with respect to each variable, $H_{2,0}^{1,\sigma}(\Omega)$ is the closure in $H_{2}^{1}(\Omega)$ of $C^{\infty}_{0}(\Omega)$, and $C$ is a constant independent of $N$. This estimate implies that every function $\vec{w}\in L^{\sigma}_{2,0}$ can be approximated with arbitrary precision by divergence-free and boundary-free real polynomials as follows: for any $n\in\mathbb{N}$, since $\{\vec{u}\in C^{\infty}_{0}(\Omega)^{d}\,:\,\nabla\cdot\vec{u}=0\}$ is dense in $L^{\sigma}_{2,0}$, there is a divergence-free $C^{\infty}$ function $\vec{u}$ with compact support in $\Omega$ such that $\|\vec{w}-\vec{u}\|_{L_{2}}\leq 2^{-(n+1)}$. Then it follows from the above inequality that there exists a positive integer $N$ and a divergence-free and boundary-free polynomial $\vec{p}$ of degree $N$ with real coefficients such that $\|\vec{u}-\vec{p}\|_{L_{2}}\leq\|\vec{u}-\vec{p}\|_{H^{3}(\Omega)^{d}}\leq 2^{-(n+1)}$. Consequently, $\|\vec{w}-\vec{p}\|_{L_{2}}\leq\|\vec{w}-\vec{u}\|_{L_{2}}+\|\vec{u}-\vec{p}\|_{L_{2}}\leq 2^{-n}$.

It remains to show that $\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$, the divergence-free and boundary-free polynomial tuples with rational coefficients, is dense (in $L_{2}$-norm) in the set of all polynomial tuples with real coefficients which are divergence-free on $\Omega$ and boundary-free on $\partial\Omega$. To this end we note that, according to part (a), the divergence-free and boundary-free polynomials can be characterized, independent of their coefficient field, in terms of a homogeneous system of linear equations with integer coefficients. Then it follows from the lemma below that the set of the rational solutions of this system is dense in the set of its real solutions. And since $\Omega$ is bounded (=relatively compact), the approximations to its coefficients of a polynomial yields (actually uniform) the approximations to the polynomial itself:

[TABLE]

Lemma 7

Let $A\in\mathbb{Q}^{m\times n}$ be a rational matrix. Then the set $\operatorname{kernel}_{IQ}(A):=\{\vec{x}\in\mathbb{Q}^{n}:A\cdot\vec{x}=\vec{0}\}$ of rational solutions to the homogeneous system of linear equations given by $A$ is dense in the set $\operatorname{kernel}_{\mathbb{R}}(A)$ of real solutions.

Proof

For d:=\operatorname{dim}\big{(}\operatorname{kernel}_{\mathbb{R}}(A)\big{)}, Gaussian Elimination yields a basis $B=(\vec{b}^{1},\ldots,\vec{b}^{d})$ of $\operatorname{kernel}(A)$; in fact it holds $B\in\mathbb{Q}^{n\times d}$ and

[TABLE]

for every field $\mathbb{F}\supseteq\mathbb{Q}$: Observe that the elementary row operations Gaussian Elimination employs to transform $A$ into echelon form containing said basis $B$ consist only of arithmetic (=field) operations! (We deliberately do not require $B$ to be orthonormal; cf. [29, §3].) Now $\operatorname{image}_{\mathbb{Q}}(B)$ is obviously dense in $\operatorname{image}_{\mathbb{R}}(B)$.

Appendix 0.B Proof of Lemma 1

Note that $\gamma_{n}\ast\mathbb{T}_{k}\vec{p}=(\gamma_{n}\ast\mathbb{T}_{k}p_{1},\gamma_{n}\ast\mathbb{T}_{k}p_{2})$. For each $\vec{p}\in\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$ and $n\geq k$, since

[TABLE]

for $\mathbb{T}_{k}p_{1}=0$ in the exterior region of $\Omega_{k}$ including its boundary $\partial\Omega_{k}$. Note that $\frac{\partial\gamma_{n}}{\partial s}$ is continuous on $\mathbb{R}^{2}$; $\frac{\partial\gamma_{n}}{\partial s}(x-s,y-t)\cdot\mathbb{T}_{k}p_{1}(s,t)$ is continuous on $[-1,1]^{2}$ for any given $x,y\in\mathbb{R}$; $\frac{\partial\mathbb{T}_{k}p_{1}}{\partial s}(s,t)$ is continuous in $(-1+2^{-n},1-2^{-n})$ and $\mathbb{T}_{k}p_{1}$ is continuous on $[-1+2^{-n},1-2^{-n}]$ for any given $t\in[-1;1]$. Thus, we can apply the integration by parts formula to the integral

[TABLE]

as follows:

[TABLE]

Then it follows from (0.B) and (0.B) that for any $(x,y)\in\Omega$,

[TABLE]

A similar calculation yields that for any $(x,y)\in\Omega$,

[TABLE]

Thus, for any $(x,y)\in\Omega$ and $n\geq k$,

[TABLE]

for $\mathbb{T}_{k}\vec{p}=(\mathbb{T}_{k}p_{1},\mathbb{T}_{k}p_{2})$ is divergence-free on $\Omega_{k}$. This proves that for any $\vec{p}\in\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$ and $n\geq k$, $\gamma_{n}\ast\mathbb{T}_{k}\vec{p}$ is divergence-free on $\Omega$.

Appendix 0.C Proof of Lemma 2

Since for each $\vec{p}\in\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$ and $k\in\mathbb{N}$, $\gamma_{n}\ast\mathbb{T}_{k}\vec{p}\to\mathbb{T}_{k}\vec{p}$ effectively and uniformly on $\Omega_{k}$ as $n\to\infty$, it suffices to show that $\{\mathbb{T}_{k}\vec{p}:k\in\mathbb{N},\vec{p}\in\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]\}$ is dense in $L^{\sigma}_{2,0}(\Omega)$. On the other hand, since $\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$ is dense in $L^{\sigma}_{2,0}(\Omega)$, we only need to show that for each $\vec{p}\in\mathbb{Q}^{\sigma}_{0}[\mathbb{R}^{2}]$ and $m\in\mathbb{N}$, there is a $k\in\mathbb{N}$ such that $2^{-m}\geq\|\vec{p}-\mathbb{T}_{k}\vec{p}\|_{\infty}=\max\{|p_{1}(\vec{x})-\mathbb{T}_{k}p_{1}(\vec{x})|,|p_{2}(\vec{x})-\mathbb{T}_{k}p_{2}(\vec{x})|:\vec{x}\in\bar{\Omega}\}$.

Since $p_{i}$ is uniformly continuous on $\bar{\Omega}$, there exists a $k\in\mathbb{N}$ such that $|p_{i}(x,y)-p_{i}(x^{\prime},y^{\prime})|\leq 2^{-m}$ whenever $|x-x^{\prime}|,|y-y^{\prime}|\leq 2^{-k+1}$, and, in particular, for $x^{\prime}=\frac{x}{1-2^{-k}}$ and $y^{\prime}=\frac{y}{1-2^{-k}}$. Also, since $p_{i}(x,y)=0$ for $(x,y)\in\partial\Omega$, $|p_{i}(x,y)|\leq 2^{-m-1}$ for all $(x,y)\in\Omega\setminus\Omega_{k}$. This establishes $|p_{i}(x,y)-\mathbb{T}_{k}p_{i}(x,y)|\leq 2^{-m}$ on $\bar{\Omega}$.