Analytic one-dimensional maps and two-dimensional ordinary differential equations can robustly simulate Turing machines

Abstract

In this paper, we analyze the problem of finding the minimum dimension n such that an analytic map/ordinary differential equation over $R^{n}$ can simulate a Turing machine in a way that is robust to perturbations. We show that one-dimensional analytic maps are sufficient to robustly simulate Turing machines; but the minimum dimension for the analytic ordinary differential equations to robustly simulate Turing machines is two, under some reasonable assumptions. We also show that any Turing machine can be simulated by a two-dimensional $C^{\infty}$ ordinary differential equation on the compact sphere $S^{2}$ .

Keywords

Turing machines Turing universal analytic map Turing universal ordinary differential equation

1. Introduction

As it is well known (see e.g. [33]), a Turing machine is a mathematical model of computation that formalizes the notion of algorithm/computation over discrete structures such as the set of positive integers $N$ or integers $Z$ . In practice, Turing machines are computationally equivalent to a standard digital computer. Furthermore, following the work of Turing and others it is also known that some problems such as the Halting Problem or Hilbert’s 10th problem are noncomputable, i.e. there is no Turing machine (i.e. no algorithm) that solves those problems [29,34]. These remarkable results show that some problems are algorithmically unsolvable. Examples of other nontrivial behavior regarding Turing machines include, for instance, the existence of universal Turing machines, which can simulate the computation of any other Turing machine (see [34] or e.g. [33]), or the existence of self-reproducing Turing machines which output their own description [31].

Although the results above are considered classically over discrete structures (e.g. $N$ ), often they can be studied over continuous spaces such as $R^{n}$ . The idea is to simulate the computation of a Turing machine with a continuous map/flow. If a continuous system is able to simulate any Turing machine (or, equivalently, a universal Turing machine), then this system is usually referred to as Turing universal. A consequence of the noncomputability of the Halting Problem is that the long term behavior of Turing universal systems is highly complex (in a manner distinct from behaviors considered e.g. in chaos theory) and has some characteristics which are not computable (see e.g. [30]). However, in applications, it is often desirable that Turing universal systems are relatively simple and mathematically well-behaved so that they can be used in meaningful situations. For this reason, one might be interested in having properties such as low-dimensionality, reasonable smoothness, or robustness to perturbations for Turing universal systems.

In this paper, we investigate the problem of determining the lowest dimension n such that the analytic maps or analytic ODEs defined on $R^{n}$ can robustly simulate Turing machines.

It is well-known that piecewise affine and other types of maps and ODEs can simulate Turing machines on $R^{n}$ , see e.g. [2,3,8,9,20,21,23,24,26,30]). However, some authors have claimed (see e.g. [1,5,25,25,28]) that, when focusing on more physically realistic systems simulating Turing machines, one might desire other additional attributes such as robustness to noise (it is known that the addition of noise to some classes of systems reduces their computational power, see e.g. [1,28]) or smoothness of the dynamics since most classical physical systems are expressed with smooth (actually analytic) functions. In [25] the authors have shown that closed-form analytic maps are capable of simulating Turing machines with exponential slowdown in dimension one or in real time in dimension $⩾ 2$ . In [18] it was shown that, under a certain notion of robustness (see Theorems 1 and 3 below), the class of closed-form analytic maps on $R^{3}$ as well as the class of ODEs defined with analytic closed-form functions in $R^{6}$ can robustly simulate Turing machines. We recall that closed-form analytic maps are analytic functions which can be expressed in terms of elementary functions such as polynomials, trigonometric functions, exponential and logarithmic functions, their compositions, etc.

In the present paper, we show that (a) one-dimensional analytic maps can robustly simulate Turing machine on $R$ (Theorem 9) in real time (i.e. without the exponential slowdown of [25]; the simulation in [25] is not robust to perturbations either); (b) two-dimensional analytic ODEs can robustly simulate Turing machines (Theorem 11), which reduces the dimension of the ODEs from 6 as in [18] to 2; (c) none of one-dimensional autonomous analytic ODEs can simulate (robustly or not) a universal Turing machine under certain reasonable assumptions; and (d) Turing machines can be simulated on the two-dimensional unit ball – a compact subset of $R^{3}$ – albeit only with $C^{\infty}$ functions. As mentioned before, results (a) and (b) improve the results obtained in [25] and [18], respectively. Result (a) is proved by simulating Turing machines directly using the technique introduced in [18] rather than simulating three-counter Minsky machines (which are Turing universal) with Collatz functions and then embedding Collatz functions in closed-form functions as shown in [25]. The direct simulation removes the exponential slowdown and ensures the robustness to perturbations. To prove result (b), we developed a dimension-reduction scheme that preserves robustness with respect to perturbations. The dimension-reduction is not a trivial matter due to the fact that $R^{k}$ and $R^{n}$ are not homeomorphic if $k \neq n$ and, subsequently, there is no continuous bijection between the two spaces. Result (d) further extends the results in [18] – valid only on unbounded spaces – by showing that Turing machines can be simulated on the two-dimensional unit ball. The proof of result (d) makes use of the ideas from [14]. Similar to what is done in [14], we show that there is a $C^{\infty}$ ODE which can simulate Turing machines over the compact set $S^{n} = {x \in R^{n + 1} : ‖ x ‖ = 1}$ for $n ⩾ 2$ . In [14], a polynomial flow (defined with a polynomial ODE) having Turing universality on $S^{n}$ for $n ⩾ 17$ is first constructed and then this polynomial flow is used to construct a (Euler) partial differential equation which is Turing universal. And it is shown in [15] that there are Turing complete (stationary Euler) fluid-flows on a Riemannian 3-dimensional sphere. The difference between our result (d) and the results of [14,15] is that we use ODEs instead of partial differential equations (for the case of [15]), and result (d) is for $S^{2}$ instead of $S^{3}$ or $S^{17}$ as in [14,15], respectively.

The outline of the paper is as follows. In Section 2 we review the construction presented in [18] to create analytic maps on $R^{3}$ which can robustly simulate Turing machines. In Section 3 we present some auxiliary functions. Building on these results, we show in Section 4 that one-dimensional analytic maps can robustly simulate Turing machines. By iterating these maps with ODEs, we are able to show in Section 5 that two-dimensional ODEs can robustly simulate Turing machines. In Section 6, we construct a $C^{\infty}$ ODE that can simulate Turing machines over the compact set $S^{2}$ . Finally, in Section 7 we show that under reasonable hypothesis, no one-dimensional analytic ODE can simulate a universal Turing machine.

2. Simulating Turing machines in dimension three

In this section, we review several results from [18] which are useful for proving our main results.

We first recall some basic results from computability theory (see e.g. [33]). Given a finite set Σ (the alphabet), a word over Σ is a finite sequence $w = (w_{1}, \dots, w_{k}) \in Σ^{k}$ for some $k \in N_{0}$ (k is the length of the word), where $N_{0}$ is the set of all non-negative integers. Note that there is a special sequence, represented by ϵ, which denotes the word of length 0. As usual, for notational simplicity, we will denote the word $w = (w_{1}, \dots, w_{k})$ simply as $w = w_{1} ... w_{k}$ . The set of all words over Σ is denoted by $Σ^{*}$ . We also recall that a Turing machine is a discrete dynamical system defined by the iteration of a map, although it is usually viewed as a finite-state machine since this approach is often more convenient. More specifically, let Σ be an alphabet, and take some symbol $B \notin Σ$ , which is usually known as the blank symbol, and let Q be a finite set known as the set of states with some special elements $q_{0}, q_{h} \in Q$ , called the initial state and the final state, respectively. Then a Turing machine M is a map $F_{M} : Σ^{*} \times Σ^{*} \times Q \to Σ^{*} \times Σ^{*} \times Q$ that works as follows when viewed as a machine. It has a bi-infinite tape, divided into cells, and a head which is associated to some state of Q. Given some $(u, v, q) \in Σ^{*} \times Σ^{*} \times Q$ (the configuration of the Turing machine), where $u = u_{1} u_{2} \dots u_{n}$ and $v = v_{1} v_{2} \dots v_{p}$ , then the tape contents of the Turing machine at this configuration is $\begin{matrix} (1) & ... B B B v_{p} ... v_{2} v_{1} u_{1} u_{2} \dots u_{n} B B B ..., \end{matrix}$ while its associated state is q. In this case the Turing machine is also said to be reading symbol $v_{1}$ . Then, depending only on the value of the current state and of the symbol being read by the head, the machine simultaneously (i) updates its state, (ii) updates the symbol being read by the head and (iii) either moves the head one cell to the right, one cell to the left, or maintains the head on the same position.

A Turing machine M computes a function $f : Σ^{*} \to Σ^{*}$ as follows. Given a word w it starts its computation on the initial configuration $(w, ϵ, q_{0})$ , i.e. in the initial configuration the state is the initial state and the tape contains the input w only. Then M proceeds with the computation until it reaches the halting state $q_{h}$ . In this case we say that the Turing machine has halted with configuration $(u_{h}, v_{h}, q_{h}) \in Σ^{*} \times Σ^{*} \times Q$ . In this case its output will be $u_{h}$ , i.e. $u_{h} = f (w)$ . If M does not halt with input w, then $f (w)$ is undefined.

Given some Turing machine M as described above, let $k = 1 + # Σ$ and take an injective map $γ : Σ \to {0, 1, 2, \dots, k - 1}$ with $γ (B) = 0$ . Let $(u, v, q)$ be the current configuration of M and let us further assume that M has m states, represented by the numbers $1, \dots, m$ , and that if M reaches an halting configuration, then it moves to the same configuration (i.e. $F_{M} (u_{h}, v_{h}, q_{h}) = (u_{h}, v_{h}, q_{h})$ ). Take $\begin{matrix} (2) & \begin{matrix} y_{1} = γ (u_{1}) + γ (u_{2}) k + \dots + γ (u_{n}) k^{n - 1} \\ y_{2} = γ (v_{1}) + γ (v_{2}) k + \dots + γ (v_{p}) k^{p - 1} \end{matrix} \end{matrix}$ and suppose that q is the state associated to the current configuration. Then $(y_{1}, y_{2}, q) \in N^{3}$ encodes unambiguously the current configuration of M. Under these assumptions, the transition function of M can be encoded as a function ${\bar{f}}_{M} : N^{3} \to N^{3}$ . In [18] it was shown that ${\bar{f}}_{M}$ can be extended to a function $f_{M} : R^{3} \to R^{3}$ , which has the following properties: (i) it is capable of simulating M in the presence of perturbations; (ii) the function f is analytic, and each of its components can be expressed using only the following terms: variables, polynomial-time computable constants (see Remark 2 for a definition), +, −, ×, sin, cos, arctan. The precise statement of this result is given below, where $‖ f ‖ = {sup}_{x \in A} ‖ f (x) ‖$ for a function $f : A \subseteq R^{l} \to R^{j}$ , $‖ y ‖ = {max}_{1 ⩽ i ⩽ j} | y_{i} |$ for $y = (y_{1}, \dots, y_{j}) \in R^{j}$ , and $f^{[k]}$ denotes the kth iterate of the function $f : A \to A$ , which is defined as follows: $f^{[0]} (x) = x$ , $f^{[k + 1]} (x) = f^{[k]} (f (x))$ .

Theorem 1 ([18, p. 333]).

Let $ψ : N^{3} \to N^{3}$ be the transition function of a Turing machine M under the encoding described above, and let $0 < δ < ε < 1 / 2$ . Then ψ admits a globally analytic closed-form extension $f_{M} : R^{3} \to R^{3}$ such that the expression of each component of $f_{M}$ can be written using only the following terms: variables, polynomial-time computable constants, +, −, ×, sin, cos, arctan. Moreover, $f_{M}$ is robust to perturbations in the following sense: for all f such that $‖ f - f_{M} ‖ ⩽ δ$ , for all $j \in N$ , and for all ${\bar{x}}_{0} \in R^{3}$ satisfying $‖ {\bar{x}}_{0} - x_{0} ‖ ⩽ ε$ , where $x_{0} \in N^{3}$ represents a configuration according to the encoding described above, $\begin{matrix} ‖ f^{[j]} ({\bar{x}}_{0}) - ψ^{[j]} (x_{0}) ‖ ⩽ ε . \end{matrix}$

We note that the proof of this theorem is constructive and that $f_{M}$ can be obtained explicitly. A continuous-time version of Theorem 1 was also proved in [18].

Remark 2.
We note that we can define computable real constants and computable real functions using the approach of computable analysis. Hence the functions $f_{M}$ given by Theorem 1 are computable in the computable analysis sense. Computation in the setting of computable analysis is performed by computing with a discrete model over a (infinite) representation of a real number, a real function, etc., while computation as described in Theorem 1 uses an analytic function defined over a real space without making use of any discrete model. However, it can be shown that the apparent different approaches are in fact computationally equivalent [4,7,18]. We also recall that the notion of a computable real number or of a computable real function in the setting of computable analysis can be presented in several equivalent but different ways. For example, according to the approach presented in [22] (see also [10]), a number $c \in R$ is computable if there is a Turing machine M that, on input $n \in N$ , outputs (in finite time) a rational $q_{n}$ with the property that $| q_{n} - c | ⩽ 2^{- n}$ . If the Turing machine M runs in polynomial time (in n), then we say that c is computable in polynomial time. Similarly, a function $f : R \to R$ is computable if there is an oracle Turing machine M that computes $f (x)$ in the sense that, given as input $n \in N$ and any oracle $φ : N \to Q$ recording $x \in R$ (i.e. with the property that $| x - φ (n) | ⩽ 2^{- n}$ ), M outputs a rational number $q_{n}$ such that $| q_{n} - c | ⩽ 2^{- n}$ . A $C^{1}$ real function $f : R \to R$ is $C^{1}$ -computable if f and its derivative are both computable. These notions can be generalized to $R^{n}$ in a straightforward manner.
Theorem 3 ([18, p. 333]).

Let $ψ : N^{3} \to N^{3}$ be the transition function of a Turing machine M under the encoding described above; let $0 < ε ⩽ 1 / 4$ ; and let $0 ⩽ δ < 2 ε ⩽ 1 / 2$ . Then there exist

$η > 0$ satisfying $η < 1 / 2$ , which can be computed from ψ, ε, δ, and

an analytic closed-form function $g_{M} : R^{7} \to R^{6}$ which can be written using only the following terms: variables, polynomial-time computable constants, +, −, ×, sin, cos, arctan

such that the ODE

z^{'} = g_{M} (t, z)

robustly simulates M in the following sense: for all g satisfying

‖ g - g_{M} ‖ ⩽ δ < 1 / 2

and for every

x_{0} \in N^{3}

that encodes a configuration according to the encoding described above, if

{\bar{x}}_{0}, {\bar{y}}_{0} \in R^{3}

satisfy the conditions

‖ {\bar{x}}_{0} - x_{0} ‖ ⩽ ε

and

‖ {\bar{y}}_{0} - x_{0} ‖ ⩽ ε

, then the solution

z (t)

\begin{matrix} z^{'} = g (t, z), z (0) = ({\bar{x}}_{0}, {\bar{y}}_{0}) \end{matrix}

satisfies, for all

j \in N

and for all

t \in [j, j + 1 / 2]

\begin{matrix} (3) & ‖ z_{2} (t) - ψ^{[j]} (x_{0}) ‖ ⩽ η, \end{matrix}

where

z \equiv (z_{1}, z_{2})

, with

z_{1} \in R^{3}

and

z_{2} \in R^{3}

3. Some useful auxiliary functions

In this section, we present several functions and results which are needed in subsequent sections. In a first reading, the reader can skim through the results without going into the proofs.

The function Υ presented in the next lemma can be seen as a generalization of the error-correcting function $l_{2}$ from [18, Lemma 9]. More specifically, the function Ψ can be viewed as a function that improves the accuracy of approximations within distance $⩽ 1 / 5$ of an integer (the function $l_{2}$ of [18, Lemma 9] has a similar property, but only works for the integers 0 and 1), where the correction factor is bounded by $e^{- y}$ , and $y > 0$ is the second argument of Ψ.

Lemma 4.
Let $Ψ : R^{2} \to R$ be given by $Ψ (x, y) = x - \frac{1}{2 π} arcsin (sin (2 π x) (1 - e^{- y - 2}))$ . Then $| Ψ (x, y) - k | < e^{- y} | x - k |$ whenever $| x - k | ⩽ 1 / 5$ for some $k \in Z$ and $y ⩾ 0$ .
Proof.
Let $\overline{Ψ} : R \to R$ be defined by $\overline{Ψ} (x) = x - \frac{1}{2 π} arcsin (sin (2 π x))$ . We note that, since $sin (2 π x)$ has period 1, $\overline{Ψ} (x) = k$ if $x \in [k - 1 / 4, k + 1 / 4]$ for some $k \in Z$ . However, although it is continuous, the function $\overline{Ψ}$ is not analytic, since it is well-known that if an analytic function is constant in a non-empty interval, e.g. $[3 / 4, 5 / 4]$ , then it should be constant everywhere on the real line $R$ , which is not the case for $\overline{Ψ}$ . The problem is that, although the composition of analytic functions yields again an analytic function, the derivative of $arcsin y$ is not defined when $y = - 1$ or $y = 1$ and thus $\overline{Ψ}$ is not analytic when $x = k - 1 / 4$ or $x = k + 1 / 4$ for some $k \in Z$ . Note, however, that arcsin is analytic in $(- 1, 1)$ . Hence, we can multiply $sin (2 π x)$ by a value $1 - e^{- y}$ (or $1 - e^{- y - 2}$ , which will be more convenient later on), which is slightly less than 1, to ensure that the resulting function Ψ is analytic, since in this way we guarantee that $- 1 < sin (2 π x) (1 - e^{- y}) < 1$ for any $x \in R$ and $y ⩾ 0$ . Next we notice that, by the mean value theorem $\begin{array}{l} | arcsin a - arcsin b | & ⩽ | a - b | max_{x \in [a, b]} \frac{1}{\sqrt{1 - x^{2}}} \\ = | a - b | max (\frac{1}{\sqrt{1 - a^{2}}}, \frac{1}{\sqrt{1 - b^{2}}}) . \end{array}$ Let us now take $g (x) = 7 x - sin (2 π x)$ . We note that $g (0) = 0$ and that $g^{'} (x) = 7 - 2 π cos (2 π x) > 0$ . Hence we conclude that g strictly increases in $[0, 1 / 5]$ , which implies $7 | x | ⩾ | sin (2 π x) |$ when $x \in [- 1 / 5, 1 / 5]$ . This implies that for $x \in [k - 1 / 5, k + 1 / 5]$ , where $k \in Z$ is arbitrary, we have $\begin{array}{l} | Ψ (x, y) - k | & = | x - \frac{1}{2 π} arcsin (sin (2 π x) (1 - e^{- y - 2})) - (x - \frac{1}{2 π} arcsin (sin (2 π x))) | \\ = \frac{1}{2 π} | arcsin (sin (2 π x)) - arcsin (sin (2 π x) (1 - e^{- y - 2})) | \\ ⩽ \frac{1}{2 π} \frac{1}{\sqrt{1 - {sin}^{2} (2 π / 5)}} | sin (2 π x) (1 - (1 - e^{- y - 2})) | \\ = \frac{1}{2 π} \frac{1}{\sqrt{1 - {sin}^{2} (2 π / 5)}} | sin (2 π (x - k)) | e^{- y - 2} \\ < | sin (2 π (x - k)) | e^{- y} e^{- 2} \\ ⩽ 7 | x - k | e^{- y} e^{- 2} \\ ⩽ | x - k | e^{- y} . \end{array}$ □

We now present another error-correcting function $σ : R \to R$ which was first presented in [18, Proposition 5]. This function is a uniform contraction around integers. Unlike Ψ, one cannot prescribe the amount of error reduction around each integer with a single application of the map σ. On the other hand its use is not restricted to a $1 / 5$ -neighborhood of integers and can be used on larger neighborhoods. This last property will be handy later on.
Lemma 5 ([18]).

Let $σ : R \to R$ be the function defined by $σ (x) = x - 0.2 sin (2 π x)$ . Let $ε \in [0, 1 / 2)$ . Then there is some contracting factor $λ_{ε} \in (0, 1)$ such that, $\forall δ \in [- ε, ε]$ , $\forall n \in Z$ , $| σ (n + δ) - n | < λ_{ε} δ$ .

The constants $λ_{ε}$ can usually be explicitly obtained. For example, as shown in [18], we can take $λ_{1 / 4} = 0.4 π - 1 \approx 0.2566371$ .

It is well known that there are bijective functions from $N^{2}$ to $N$ . An example (see e.g. [32, pp. 26–27]) is the dovetailing pairing map $I : N^{2} \to N$ defined by the formula $\begin{matrix} (4) & I (x, y) = \frac{{(x + y)}^{2} + 3 x + y}{2} . \end{matrix}$ Using this map we can obtain a bijective map $I_{k} : N^{k} \to N$ , for $k ⩾ 2$ , by defining $I_{k}$ recursively: $I_{2} (x_{1}, x_{2}) = I (x_{1}, x_{2})$ ; $I_{k + 1} (x_{1}, \dots, x_{k}, x_{k + 1}) = I_{2} (I_{k} (x_{1}, \dots, x_{k}), x_{k + 1})$ . We now show that the maps $I_{k}$ can be extended to $R^{k}$ robustly around the integers. Since each $I_{k}$ is a (multivariate) polynomial, to achieve this objective we have to analyze how the error is propagated via the application of a polynomial map. The following lemma is from [6], and can be viewed as an extension of a similar result proved in [18, Lemma 11] for the case of monomials. For multivariate polynomials, the multi-index notation is used for compactness as follows: a monomial $x_{1}^{α_{1}} \dots x_{k}^{α_{k}}$ is represented by $x^{α}$ , where $x = (x_{1}, \dots, x_{k})$ , $α = (α_{1}, \dots, α_{k})$ , $| α | = α_{1} + \dots + α_{k}$ is the degree of the monomial, and the degree of a (multivariate) polynomial is the maximum degree of all the monomials which appear in its expression.

Lemma 6 ([6, Lemma 4]).

Let $P : R^{k} \to R$ be a multivariate polynomial of degree k and let $x, y \in R^{k}$ be such that $‖ x ‖, ‖ y ‖ ⩽ M$ for some $M ⩾ 0$ . Then $\begin{matrix} | P (x) - P (y) | ⩽ k M^{k - 1} Σ P ‖ x - y ‖ \end{matrix}$ where $Σ P$ denotes the sum of the absolute values of the coefficients of P.

Now we are ready to state the result that shows the existence of robust analytic extensions of $I_{k}$ for each k.

Proposition 7.
For each $k \in N$ , $k ⩾ 2$ , there exists an analytic function $Υ_{k} : R^{k} \to R$ with the following properties:
If $x \in N^{k}$ , then $Υ_{k} (x) = I_{k} (x)$ ;

For any $x \in R^{k}$ , if there is some $y \in N^{k}$ such that $‖ x - y ‖ ⩽ 1 / 5$ , then $| Υ_{k} (x) - I_{k} (y) | ⩽ ‖ x - y ‖ ⩽ 1 / 5$ .

Proof.
We start with the case $k = 2$ . Since $\begin{matrix} I_{2} (x_{1}, x_{2}) = \frac{{(x_{1} + x_{2})}^{2} + 3 x_{1} + x_{2}}{2} = \frac{x_{1}^{2} + 2 x_{1} x_{2} + x_{2}^{2} + 3 x_{1} + x_{2}}{2}, \end{matrix}$ it is clear that the function $I_{2}$ is well-defined for all $x_{1}, x_{2} \in R$ . In the remaining of this proof, we assume that $I_{2}$ is defined over $R^{2}$ . Since $I_{2}$ is a polynomial of degree 2, by Lemma 6 we conclude that $\begin{matrix} | I_{2} (x) - I_{2} (y) | ⩽ 8 max (‖ x ‖, ‖ y ‖) ‖ x - y ‖ . \end{matrix}$ Since $‖ x ‖ ⩽ 1 + ‖ x ‖^{2}$ for all $x \in R^{2}$ , it follows that $\begin{matrix} (5) & | I_{2} (x) - I_{2} (y) | ⩽ 8 (2 + ‖ x ‖^{2} + ‖ y ‖^{2}) ‖ x - y ‖ . \end{matrix}$ Set $\begin{matrix} Υ_{2} (x_{1}, x_{2}) = I_{2} (Ψ (x_{1}, 32 (1 + ‖ x ‖_{2}^{2})), Ψ (x_{2}, 32 (1 + ‖ x ‖_{2}^{2}))), \end{matrix}$ where $‖ x ‖_{2}^{2} = x_{1}^{2} + x_{2}^{2}$ . As a composition of analytic functions, $Υ_{2} : R^{2} \to R$ is clearly analytic. If $x_{1}, x_{2} \in N$ , it is trivial to verify that $Υ_{2} (x_{1}, x_{2}) = I_{2} (x_{1}, x_{2})$ , which implies property 1. For property 2, let us assume that $‖ x - y ‖ ⩽ ε ⩽ 1 / 5$ for some $y \in N^{k}$ . Using Lemma 4 and the inequality $e^{- l} < (1 / l)$ for all $l ⩾ 1$ , we obtain the following estimate, where $ε = ‖ x - y ‖ ⩽ 1 / 5$ : $\begin{array}{l} ‖ (Ψ (x_{1}, 32 (1 + ‖ x ‖_{2}^{2})), Ψ (x_{2}, 32 (1 + ‖ x ‖_{2}^{2}))) - (y_{1}, y_{2}) ‖ \\ ⩽ ε e^{- 32 (1 + ‖ x ‖_{2}^{2})} \\ (6) & ⩽ \frac{ε}{32 (1 + ‖ x ‖_{2}^{2})} . \end{array}$ Recall that, on $R^{2}$ , $‖ \cdot ‖_{2}$ denotes the Euclidean norm while $‖ \cdot ‖$ denotes the maximum norm. Since $‖ x - y ‖ ⩽ \frac{1}{5}$ and $‖ y ‖_{2} ⩽ \sqrt{2} ‖ y ‖$ , it follows that $‖ x - y ‖_{2} ⩽ \sqrt{2} ‖ x - y ‖ ⩽ 1 / 2$ , which further implies that $\begin{array}{l} ‖ y ‖_{2}^{2} & = ‖ y - x + x ‖_{2}^{2} \\ ⩽ {(‖ x - y ‖_{2} + ‖ x ‖_{2})}^{2} \\ ⩽ {(\frac{1}{2} + ‖ x ‖_{2})}^{2} \\ = \frac{1}{4} + ‖ x ‖_{2} + ‖ x ‖_{2}^{2} \\ ⩽ \frac{1}{4} + 1 + ‖ x ‖_{2}^{2} + ‖ x ‖_{2}^{2} < 2 + 2 ‖ x ‖_{2}^{2} . \end{array}$ Then it follows from this inequality, (5), and (6) that $\begin{array}{l} | Υ_{2} (x) - I_{2} (y) | \\ = | Υ_{2} (x) - Υ_{2} (y) | \\ ⩽ 8 max (‖ x ‖, ‖ y ‖) ‖ (Ψ (x_{1}, 32 (1 + ‖ x ‖_{2}^{2})), Ψ (x_{2}, 32 (1 + ‖ x ‖_{2}^{2}))) - (y_{1}, y_{2}) ‖ \\ ⩽ 8 (2 + ‖ x ‖^{2} + ‖ y ‖^{2}) \frac{ε}{32 (1 + ‖ x ‖_{2}^{2})} \\ ⩽ (2 + ‖ x ‖_{2}^{2} + ‖ y ‖_{2}^{2}) \frac{ε}{(2 + 2 ‖ x ‖_{2}^{2} + (2 + 2 ‖ x ‖_{2}^{2}))} (‖ \cdot ‖ ⩽ ‖ \cdot ‖_{2}) \\ ⩽ (2 + ‖ x ‖_{2}^{2} + ‖ y ‖_{2}^{2}) \frac{ε}{(2 + ‖ x ‖_{2}^{2} + ‖ y ‖_{2}^{2})} = ε \end{array}$ which proves property 2 for $k = 2$ .

For the case where $k > 2$ , the result is obtained inductively by setting $\begin{matrix} Υ_{k + 1} (x_{1}, \dots, x_{k}, x_{k + 1}) = Υ_{2} (Υ_{k} (x_{1}, \dots, x_{k}), x_{k + 1}) . \end{matrix}$ Property 1 is immediate; property 2 follows from the estimate below: $\begin{array}{l} ‖ Υ_{k + 1} (x_{1}, \dots, x_{k}, x_{k + 1}) - I_{k + 1} (y_{1}, \dots, y_{k}, y_{k + 1}) ‖ \\ = ‖ Υ_{k + 1} (x_{1}, \dots, x_{k}, x_{k + 1}) - Υ_{k + 1} (y_{1}, \dots, y_{k}, y_{k + 1}) ‖ \\ = ‖ Υ_{2} (Υ_{k} (x_{1}, \dots, x_{k}), x_{k + 1}) - Υ_{2} (Υ_{k} (y_{1}, \dots, y_{k}), y_{k + 1}) ‖ \\ ⩽ ‖ (Υ_{k} (x_{1}, \dots, x_{k}), x_{k + 1}) - (Υ_{k} (y_{1}, \dots, y_{k}), y_{k + 1}) ‖ \\ ⩽ max (‖ Υ_{k} (x_{1}, \dots, x_{k}) - Υ_{k} (y_{1}, \dots, y_{k}) ‖, ‖ x_{k + 1} - y_{k + 1} ‖) \\ ⩽ max (‖ x_{1} - y_{1} ‖, ‖ x_{2} - y_{2} ‖, \dots, ‖ x_{k + 1} - y_{k + 1} ‖) \\ ⩽ ‖ x - y ‖ . \end{array}$ □

As shown above, $I_{k} : N^{k} \to N$ provides a bijection between $N^{k}$ and $N$ that can be robustly extended to an analytic function $Υ_{k} : R^{k} \to R$ . We now show that the inverse function of $I_{k}$ can also be robustly extended to an analytic function from $R$ to $R^{k}$ . We write $I_{k}^{- 1} (z) = (J_{k, 1} (z), \dots, J_{k, k} (z))$ , where $J_{k, 1}, \dots, J_{k, k} : N \to N$ . The following result shows that $I_{k}^{- 1} : {N \to N}^{k}$ can be robustly extended to an analytic function from $R$ to $R^{k}$ . Proposition 8.
For each $k \in N$ , $k ⩾ 2$ , and for each $1 ⩽ i ⩽ k$ , there exists an analytic function $Ω_{k, i} : R \to R$ with the following property: for any $x \in R$ , if there is some $n \in N$ such that $| x - n | ⩽ 1 / 5$ , then $| Ω_{k, i} (x) - J_{k, i} (n) | ⩽ 1 / 5$ .
Proof.
First we prove the result when $k = 2$ . Let us assume that $i = 1$ (the case where $i = 2$ is similar). Then, by definition, $J_{2, 1} (I_{2} (x_{1}, x_{2})) = x_{1}$ . Since $I_{2} (x_{1}, x_{2}) ⩾ x_{i}$ for $i = 1, 2$ or, more generally, $I_{k} (x_{1}, \dots, x_{k}) ⩾ x_{i}$ for all $i = 1, \dots, k$ (see e.g. [32, p. 27]), we have the following algorithm to compute $J_{2, 1}$ , given some input $x \in N$ :
For all $i = 1, \dots, x$

For all $j = 1, \dots, x$

If $I_{2} (i, j) = x$ , then output i

Next j

Next i

This algorithm always stops with the correct result. Hence $J_{2, 1}$ can be computed by a one tape Turing machine M. Furthermore, using well-known techniques, we can assume that M has the following properties: (i) the tape alphabet of the Turing machine is ${B, 1}$ where B denotes the blank symbol; (ii) the input alphabet is ${1}$ ; (iii) each input $z \in N$ and the respective output of the computation is represented in unary, i.e. by a sequence of z 1’s; and (iv) $J_{2, 1}$ is computed by M in time $P (n) = P (x)$ , where P is a polynomial which can be explicitly obtained and which is assumed to be an increasing function. Regarding condition (i), we notice that there are universal Turing machines which only use the alphabet ${B, 1}$ and hence we do not lose computational computational power with respect to Turing machines using more symbols. For example, if we have a Turing machine $M_{1}$ which tape alphabet has $k > 3$ symbols (including the blank symbol B), then we can create a Turing machine $M_{2}$ with tape alphabet ${B, 0, 1}$ which simulates $M_{1}$ by taking some fixed $l \in N$ satisfying $l ⩾ {log}_{2} (k)$ such that each symbol of $M_{1}$ is represented by distinct strings (blocks) of ${0, 1}^{*}$ of length l, with the exception of the blank symbol of $M_{1}$ which is represented by a block of l blank symbols in $M_{2}$ . By its turn, $M_{2}$ can be simulated by a Turing machine $M_{3}$ with tape alphabet ${B, 1}$ by coding each symbol of $M_{2}$ as a string of length 2 in $M_{3}$ , e.g. by coding 0, 1, and B as $1 B$ , 11, and $B B$ , respectively. We note that regarding condition (iv), the expression of P depends on the exact implementation details of M, but we prefer to omit the exact description of M and of P for brevity (these can be obtained as usual, although the procedure is a bit tedious and hence not of much interest for this proof).

Let $g_{M} : R^{7} \to R^{6}$ be the function given by Theorem 3 such that $y^{'} = g_{M} (y)$ simulates M (taking $ε = 1 / 5$ in that theorem), and let $η < 1 / 2$ be the associated constant such that (3) holds. (The value of δ in Theorem 3 is not really relevant; one may simply take $δ = 1 / 5$ ).) Let $l \in N$ be chosen such that $σ^{[l]} (η) ⩽ 1 / 5$ (see Lemma 5). Given an input $x \in N$ for the Turing machine M, let us assume that this input is encoded in unary (i.e. x is represented by a sequence of x 1’s) when processed by M. We can then transform this unary coding of x into another integer value $φ_{1} (x)$ via the coding (2), where $γ (B) = 0$ and $γ (1) = 1$ , which can then be used to create an initial condition for $y^{'} = g_{M} (y)$ such that this IVP simulates M with input x. Note that although $x \in N$ and $φ_{1} (x) \in N$ , we do not necessarily have $x = φ_{1} (x)$ . Since initially the tape will be empty, with the exception of the input, and M will be on its initial state, which we assume to be the state 1 (we can assume, without loss of generality, that the states of M correspond to the elements of ${1, 2, \dots, m}$ ), then the initial configuration of M will be coded as $(φ_{1} (x), 0, 1)$ . Let $Φ (t, φ_{1} (x))$ denote the solution of $y^{'} = g_{M} (y)$ with initial condition associated to the configuration $(φ_{1} (x), 0, 1)$ and let $π_{i}^{k} : R^{k} \to R$ be the projection $π_{i}^{k} (x_{1}, \dots, x_{k}) = x_{i}$ for $1 ⩽ i ⩽ k$ . Note that Φ is analytic and that M computes $J_{2, 1}$ . We will use these facts to create the function $Ω_{2, 1}$ to be defined as $Ω_{2, 1} (x) = φ_{2} \circ π_{4}^{6} \circ σ^{[l]} \circ Φ (P (x + 1), φ_{1} (x))$ for some analytic functions $φ_{1}$ , $φ_{2}$ yet to be defined, which essentially translate the value of $x \in N$ into the coding (2) of its unary representation (case of $φ_{1}$ ) and, reciprocally, converts the coding of the unary representation back to the number encoded by this representation (note again that $x \in N$ may not be equal to the number $\bar{x} \in N$ encoding the symbolic representation – unary, binary, etc. – of x given by (2)). Since $| x - n | ⩽ 1 / 5$ and P is assumed to be increasing, it follows that $x + 1 > n ⩾ 0$ and $P (x + 1) ⩾ P (n)$ . Hence, if $| x - n | ⩽ 1 / 5$ implies that $| φ_{1} (x) - φ_{1} (n) | ⩽ 1 / 5$ , we get that $Φ (P (x + 1), φ_{1} (x))$ will return the coding of the output of M with the input encoding the number $x \in N$ (note that although the relation (3) is in general valid only in intervals of the format $[j, j + 1 / 2]$ with $j \in N$ , but since we have assumed that the image of an halting configuration is itself, it follows from the results of [18] that (3) is valid for all times $[j + 1 / 2, j + 1]$ after the Turing machine has halted. See also Remark 15). We now only have to define the functions $φ_{1}$ and $φ_{2}$ . Let us now first turn our attention to $φ_{1}$ . Note that given some $n \in N$ , the number $2^{n} - 1$ will represent n in unary when using the coding (2) (taking $k = 2$ , since $γ (B) = 0$ by definition, and by taking $γ (1) = 1$ ). Hence it makes sense to take $φ_{1} (n) = 2^{n} - 1$ . However, we cannot take $φ_{1} (x)$ to be $2^{x} - 1$ , because in that case we cannot ensure that $| x - n | ⩽ 1 / 5$ implies $| φ_{1} (x) - φ_{1} (n) | ⩽ 1 / 5$ . To avoid this problem, we improve the accuracy of x using the function Ψ from Lemma 4, obtaining an improved estimate $\bar{x}$ satisfying $| 2^{\bar{x}} - 2^{n} | ⩽ 1 / 5$ . We now determine the accuracy improvement needed to achieve this objective. Note that the exponential function $2^{x}$ is strictly increasing and thus, by the mean value theorem, we have $\begin{matrix} | 2^{\bar{x}} - 2^{n} | ⩽ 2^{max (\bar{x}, n)} ln 2 | x - n | < 2^{x + 1} | \bar{x} - n | . \end{matrix}$ Hence, if we have $| \bar{x} - n | ⩽ 2^{- (x + 4)}$ , we get $| 2^{\bar{x}} - 2^{n} | ⩽ 1 / 5$ . This is achieved if $\bar{x} = Ψ (x, x + 2)$ , due to Lemma 4 and from the property that $| x - n | ⩽ 1 / 5 < 2^{- 2}$ . Hence we can take $φ_{1} (x) = 2^{Ψ (x, x + 2)} - 1$ .

We now proceed with a similar reasoning for $φ_{2}$ . We first note that if $n \in N$ codes the exact output of M according to (2), and thus represents in unary some number $i \in N$ , we will have $n = 2^{i} - 1$ as we have already seen. This implies that $i = {log}_{2} (n + 1)$ . Now we have to analyze again the effect of replacing n by some real value x satisfying $| x - n | ⩽ 1 / 5$ . By the mean value theorem, we have (note also that $n ⩾ 0$ since $n \in N$ ) $\begin{array}{l} | {log}_{2} (x + 1) - {log}_{2} (n + 1) | & ⩽ \frac{1}{ln 2 (min (x, n) + 1)} | (x + 1) - (n + 1) | \\ ⩽ \frac{1}{ln 2 (n + 4 / 5)} | x - n | \\ ⩽ \frac{5}{4 ln 2} | x - n | \\ < 2 | x - n | . \end{array}$ Therefore, to ensure that $| x - n | ⩽ 1 / 5$ implies that $| {log}_{2} (x + 1) - {log}_{2} (n + 1) | ⩽ 1 / 5$ , it is enough to take $φ_{2} (x) = {log}_{2} (Ψ (x, 2) + 1)$ (using Lemma 4) or $φ_{2} (x) = {log}_{2} (σ (x) + 1)$ (using Lemma 5 and noting, as mentioned in [18, Remark 6], that we can take $λ_{1 / 4} = 0.4 π - 1 \approx 0.2566371$ ).

Proceeding similarly for the case of $J_{2, 2}$ , we conclude that $Ω_{2, 2} (x) = φ_{2} \circ π_{4}^{6} \circ σ^{[l]} \circ Φ_{M^{'}} (P_{M^{'}} (x + 1), φ_{1} (x))$ , where $M^{'}$ is a TM machine computing $J_{2, 2}$ which is similar to M, with the difference that in Step 3 of the pseudo-algorithm above we take “If $I_{2} (i, j) = x$ , then output j”. The results for $k > 2$ follow inductively. □

4. Analytic one-dimensional maps robustly simulate Turing machines

We now present one of the main results of this paper. Let M be a one-tape Turing machine and let $(y_{1}, y_{2}, q) \in N^{3}$ be the encoding of a configuration as given in Section 2 and (2). In what follows each configuration $(y_{1}, y_{2}, q)$ is encoded in the single value $\begin{matrix} c = C (y_{1}, y_{2}, q) = I_{3} (y_{1}, y_{2}, q) \in N . \end{matrix}$ Thus we can consider that, under this new encoding, the transition function of a Turing machine is a map $ψ : N \to N$ .

Theorem 9.
Let $ψ : N \to N$ be the transition function of a Turing machine M, under the encoding described above, and let $0 ⩽ δ < 1 / 5$ . Then there is an analytic function $g_{M} : R \to R$ that robustly simulates M in the following sense: for all g such that $‖ g - g_{M} ‖ ⩽ δ$ , and for all ${\bar{x}}_{0} \in R$ satisfying $| {\bar{x}}_{0} - x_{0} | ⩽ 1 / 5$ , where $x_{0} \in N$ represents some configuration, one has for all $j \in N$ $\begin{matrix} (7) & | g^{[j]} ({\bar{x}}_{0}) - ψ^{[j]} (x_{0}) | ⩽ 1 / 5 . \end{matrix}$
Proof.
Most of the work to prove this theorem was already done in Section 3. The idea is to compose the robust simulation of Turing machines of Theorem 1 from $R^{3}$ to $R^{3}$ with the map $Υ_{3}$ (which allows to go from dimension 3 to 1) and its inverse (which allows to go back from dimension 1 to 3). Since all involved functions are robust to perturbations, the same can be said to their composition. Hence, the composition will provide a robust simulation of Turing machines from $R$ to $R$ .

We now present the details. Let us first define a function ${\bar{g}}_{M} : R \to R$ that robustly simulates M in a weaker sense that just the input can be perturbed and not ${\bar{g}}_{M}$ itself. More specifically, let us define a function ${\bar{g}}_{M} : R \to R$ with the following property: for all ${\bar{x}}_{0} \in R$ satisfying $| {\bar{x}}_{0} - x_{0} | ⩽ 1 / 5$ , where $x_{0} \in N$ represents some configuration, one has for all $j \in N$ $\begin{matrix} (8) & | {\bar{g}}_{M}^{[j]} ({\bar{x}}_{0}) - ψ^{[j]} (x_{0}) | ⩽ 1 / 5 . \end{matrix}$ To achieve this purpose, let $f_{M}$ be the corresponding 3-dimensional map simulating M obtained via Theorem 1 with $ε = 1 / 5$ . Then we take $\begin{matrix} {\bar{g}}_{M} (x) = Υ_{3} \circ f_{M} (Ω_{3, 1} (x), Ω_{3, 2} (x), Ω_{3, 3} (x)) . \end{matrix}$ It then follows from Theorem 1, Proposition 7, and Proposition 8 that property (8) is satisfied. We now only have to take care of the perturbations to ${\bar{g}}_{M}$ . Let σ be the function defined in Lemma 5. Let $j \in N$ be some integer such that $0 < λ_{1 / 4}^{j} / 5 < 1 / 5 - δ$ and take $\begin{matrix} g_{M} (x) = σ^{[j]} \circ {\bar{g}}_{M} (x) . \end{matrix}$ Then, by property (8) and Lemma 5, if ${\bar{x}}_{0} \in R$ satisfies $| {\bar{x}}_{0} - x_{0} | ⩽ 1 / 5$ , where $x_{0} \in N$ represents some configuration, one has $\begin{matrix} | g_{M} ({\bar{x}}_{0}) - ψ (x_{0}) | ⩽ 1 / 5 - δ . \end{matrix}$ If $‖ g - g_{M} ‖ ⩽ δ$ , then we conclude that $\begin{array}{l} | g ({\bar{x}}_{0}) - ψ (x_{0}) | & ⩽ | g ({\bar{x}}_{0}) - g_{M} ({\bar{x}}_{0}) | + | g_{M} ({\bar{x}}_{0}) - ψ (x_{0}) | \\ ⩽ δ + (1 / 5 - δ) \\ ⩽ 1 / 5 . \end{array}$ By using this last inequality and by iterating g and ψ, we conclude that the property (7) holds. □
Remark 10.
In the statement of Theorem 9, we could have picked some fixed $ε > 0$ satisfying $δ < ε ⩽ 1 / 5$ and, instead of assuming that $| {\bar{x}}_{0} - x_{0} | ⩽ 1 / 5$ , we could have assumed that $| {\bar{x}}_{0} - x_{0} | ⩽ ε$ and required that $| g^{[j]} ({\bar{x}}_{0}) - ψ^{[j]} (x_{0}) | ⩽ ε$ for condition (7). To see this it would be enough to compose $g_{M}$ with $σ^{[l]}$ , where σ is given by Lemma 5 and $l \in N$ is such that $λ_{1 / 4}^{l} / 4 ⩽ ε - δ$ , with $λ_{1 / 4} = 0.4 π - 1 \approx 0.2566371$ .

5. Analytic two-dimensional ODEs can robustly simulate Turing machines

In this section we construct an analytic two-dimensional ODE that robustly simulates Turing machines in the sense of Theorem 3. To prove this result, we simulate the iteration of the one-dimensional analytic function provided by Theorem 9 using a two-dimensional analytic ODE. Then it follows from Theorem 9 that this ODE will simulate a TM. The approach is similar to that used in [18]. More precisely, the following theorem is proved in this section.

Theorem 11.
Let $ψ : N \to N$ be the transition function of a Turing machine M, under the encoding described in Section 4 , and let $0 ⩽ δ < 2 / 5$ . Then there exist:
$η > 0$ satisfying $η < 2 / 5 < 1 / 2$ , which can be computed from δ; and

an analytic function $g_{M} : R^{3} \to R^{2}$
such that the ODE $z^{'} = g_{M} (t, z)$ robustly simulates M in the following sense: for all g satisfying $‖ g - g_{M} ‖ ⩽ δ < 2 / 5$ and for all $x_{0} \in N$ which encodes a configuration according to the encoding described above, if ${\bar{x}}_{0}, {\bar{y}}_{0} \in R$ satisfy the conditions $‖ {\bar{x}}_{0} - x_{0} ‖ ⩽ 1 / 5$ and $‖ {\bar{y}}_{0} - x_{0} ‖ ⩽ 1 / 5$ , then the solution $z (t)$ of $\begin{matrix} z^{'} = g (t, z), z (0) = ({\bar{x}}_{0}, {\bar{y}}_{0}) \end{matrix}$ satisfies, for all $j \in N_{0}$ and for all $t \in [j, j + 1 / 2]$ , $\begin{matrix} ‖ z_{2} (t) - ψ^{[j]} (x_{0}) ‖ ⩽ η, \end{matrix}$ where $z (t) \equiv (z_{1} (t), z_{2} (t)) \in R^{2}$ .

The remaining of this section is devoted to the proof of Theorem 11. We first present the main ideas in [18] for simulating the iteration of a map defined over the integers, which admits a robust analytic real extension, using an analytic ODE. The basic idea is to simulate the iteration of a map f (in our case the transition function $g_{M}$ of a Turing machine M given by Theorem 9) using a two dimensional ODE with solution $z_{1}$ , $z_{2}$ . In each time interval $[k, k + 1]$ , $k \in N$ , one of the components will be fixed on the first half unit interval $[k, k + 1 / 2]$ – serving as a “memory” of the result $f^{[k]} (x_{0})$ obtained most recently – to allow the other component to be updated from $f^{[k]} (x_{0})$ to its correct value $f^{[k + 1]} (x_{0})$ . In the second half unit interval $[k + 1 / 2, k + 1]$ the roles of the components are reversed. Hence, at time $t = k + 1$ both components will have the value $f^{[k + 1]} (x_{0})$ (see Fig. 1). Since this procedure repeats itself for subsequent time intervals, one will be able to simulate a Turing machine with an ODE. The fact that we will be using analytic functions poses additional challenges. The challenges are to be overcome by using the robustness to perturbation of the function $g_{M}$ from Theorem 9.

Figure 1.
Iterating the exponential function $2^{x}$ with an ODE.

We begin with a construction that uses an analytic ODE to approximate a value $b \in R$ in a finite (given) amount of time with some (given) accuracy. This construction will be needed to derive the ODE simulating the iteration of the map. Consider the following basic ODE $\begin{matrix} (9) & y^{'} = c {(b - y)}^{3} ϕ (t), \end{matrix}$ which was already studied in [11,13], [18, Section 7]. The ODE can be easily solved by separating variables, which gives rise to the following result.
Lemma 12 ([18]).

Consider a point $b \in R$ (the target), some $γ > 0$ (the targeting error), time instants $t_{0}$ (departure time) and $t_{1}$ (arrival time), with $t_{1} > t_{0}$ , and a function $ϕ : R \to R$ with the property that $ϕ (t) ⩾ 0$ for all $t ⩾ t_{0}$ and $\int_{t_{0}}^{t_{1}} ϕ (t) d t > 0$ . Then the IVP defined by ( 9 ) (the targeting equation) with the initial condition $y (t_{0}) = y_{0}$ and $\begin{matrix} (10) & c ⩾ \frac{1}{2 γ^{2} \int_{t_{0}}^{t_{1}} ϕ (t) d t} \end{matrix}$ has the property that $| y (t) - b | < γ$ for $t ⩾ t_{1}$ , independently of the initial condition $y_{0} \in R$ .

However, since we wish the ODE simulating Turing machines to be robust to perturbations, we have to analyze a perturbed version of (9).

Lemma 13 ([18]).

Consider a point $b \in R$ (the target), some $γ > 0$ (the targeting error), time instants $t_{0}$ (departure time) and $t_{1}$ (arrival time), with $t_{1} > t_{0}$ , and a function $ϕ : R \to R$ with the property that $ϕ (t) ⩾ 0$ for all $t ⩾ t_{0}$ and $\int_{t_{0}}^{t_{1}} ϕ (t) d t > 0$ . Let $ρ, δ ⩾ 0$ and let $\bar{b}$ , $E : R \to R$ be functions with the property that $| \overline{b} (t) - b | ⩽ ρ$ and $| E (t) | ⩽ δ$ for all $t ⩾ t_{0}$ . Then the IVP defined by $\begin{matrix} (11) & z^{'} = c {(\overline{b} (t) - z)}^{3} ϕ (t) + E (t), \end{matrix}$ with the initial condition $z (t_{0}) = {\bar{z}}_{0}$ , where c satisfies ( 10 ), has the property that $| z (t_{1}) - b | < ρ + γ + δ (t_{1} - t_{0})$ , independently of the initial condition ${\bar{z}}_{0} \in R$ .

Proceeding along the lines of the argument presented in [18], we now show how the map given by Theorem 9 can be iterated by a 2-dimensional ODE. Although our objective is to obtain an analytic function $g_{M}$ defining an ODE $z^{'} = g_{M} (t, z)$ which simulates a given Turing machine M, in a first step we iterate the map given by Theorem 9 by using a non-analytic ODE.

Following the approach provided in [12, p. 37], let $θ : R \to R$ be the $C^{\infty}$ function defined by $\begin{matrix} (12) & θ (x) = 0 if x ⩽ 0, θ (x) = e^{- \frac{1}{x}} if x ⩾ 0 . \end{matrix}$ Next we define the $C^{\infty}$ function $v : R \to R$ given by $\begin{matrix} (13) & v (0) = 0, v^{'} (x) = \bar{c} θ (- sin 2 π x), \end{matrix}$ where $\begin{matrix} \bar{c} = {(\int_{0}^{1} θ (- sin 2 π x) d x)}^{- 1} = {(\int_{1 / 2}^{1} e^{\frac{1}{sin 2 π x}} d x)}^{- 1} \end{matrix}$ The function v has the property that $v (x) = n$ , whenever $x \in [0, n + 1 / 2]$ . We now get the following lemma.

Lemma 14.
The $C^{\infty}$ function $r : R \to R$ defined by $r (x) = v (x + 1 / 4)$ has the property that $r (x) = n$ , whenever $x \in [n - 1 / 4, n + 1 / 4]$ , for all integers n.

Note that the function r can be seen as a function that returns the integer part of a real number around a $1 / 4$ -vicinity of an integer.

We now consider the ODE $\begin{matrix} (14) & \{\begin{matrix} z_{1}^{'} = \tilde{c} {(\tilde{f} (r (z_{2})) - z_{1})}^{3} θ (sin 2 π t), \\ z_{2}^{'} = \tilde{c} {(r (z_{1}) - z_{2})}^{3} θ (- sin 2 π t), \end{matrix} \end{matrix}$ where $\tilde{f} : R \to R$ is an extension of the function $f : N \to N$ , $z_{1} (0) = z_{2} (0) = x_{0} \in N$ and $\tilde{c}$ is a constant yet to be defined. We will next show that (14) iterates the map f near integers. Its behavior is depicted in Fig. 1 when iterating the exponential function $2^{x}$ . Suppose that $t \in [0, 1 / 2]$ . Then $z_{2}^{'} (t) = 0$ and thus $z_{2} (t) = x_{0}$ and $r (z_{2}) = x_{0}$ . In this manner, the first equation of (14) behaves like the targeting equation (9), where $b = \tilde{f} (r (z_{2})) = f (x_{0})$ , $t_{0} = 0$ , $t_{1} = 1 / 2$ , and $ϕ (t) = θ_{j} (sin 2 π t)$ and $\tilde{c}$ has to satisfy the condition (10) for c. Now note that $sin 2 π t ⩾ 1 / \sqrt{2}$ and thus $- 1 / sin (2 π t) ⩾ - \sqrt{2}$ when $t \in [1 / 8, 3 / 8]$ , which implies that $\begin{array}{l} \int_{0}^{1 / 2} e^{- \frac{1}{sin 2 π t}} d t & ⩾ \int_{1 / 4}^{2 / 8} e^{- \frac{1}{sin 2 π t}} d t \\ ⩾ \frac{1}{4} e^{- \sqrt{2}} . \end{array}$ This implies that $\begin{array}{l} \frac{1}{2 γ^{2} \int_{0}^{1 / 2} θ (sin 2 π t) d t} & = \frac{1}{2 γ^{2} \int_{0}^{1 / 2} e^{- \frac{1}{sin 2 π t}} d t} \\ ⩽ \frac{2 e^{\sqrt{2}}}{γ^{2}} . \end{array}$ Therefore, due to Lemma 12 and (10), if we take $\tilde{c} ⩾ 2 e^{\sqrt{2}} / γ^{2}$ the first equation of (14) becomes a targetting equation on the time interval $[0, 1 / 2]$ associated to a targetting error γ. In particular, if we pick $γ = 1 / 5$ , then we can pick $\tilde{c} = 206$ such that (10) holds and thus $\begin{matrix} | b - z_{1} (1 / 2) | = | f (x_{0}) - z_{1} (1 / 2) | ⩽ 1 / 5 . \end{matrix}$ On the time interval $[1 / 2, 1]$ , the roles of $z_{1}$ and $z_{2}$ are switched: we will have $z_{1}^{'} (t) = 0$ which implies that $z_{1} (t) = z_{1} (1 / 2)$ for all $t \in [1 / 2, 1]$ . We thus conclude that $r (z_{1} (t)) = f (x_{0})$ for $t \in [1 / 2, 1]$ and therefore the second equation of (14) behaves like the targeting equation (9) with $b = r (z_{1} (t)) = f (x_{0})$ , $t_{0} = 1 / 2$ , $t_{1} = 1$ , and $ϕ (t) = θ_{j} (- sin 2 π t)$ and $\tilde{c} = 206$ . Again, using Lemma 12 and similar arguments as in the previous case, we conclude that $\begin{matrix} | b - z_{2} (1) | = | f (x_{0}) - z_{2} (1) | ⩽ 1 / 5 . \end{matrix}$ In the following time interval $[1, 3 / 2]$ , the cycle repeats itself and we have $z_{2}^{'} (t) = 0$ and thus $\tilde{f} (r (z_{2})) = f (f (x_{0})) = f^{[2]} (x_{0})$ . Using a similar reasoning, we conclude that for all $k \in N_{0}$ we have $z_{2} (t) = f^{[k]} (x_{0})$ for all $t \in [k, k + 1 / 2]$ , (assuming $f^{[0]} (x) = x$ ) and $z_{1} (t) = f^{[k + 1]} (x_{0})$ for all $t \in [k + 1 / 2, k + 1]$ .

We thus have shown how to iterate (the extension of) a discrete function $f : N \to N$ with an ODE. Nonetheless, the ODE is still not analytic as required. Remark that if the ODE is analytic, then $z_{1}^{'}$ and $z_{2}^{'}$ cannot be 0 in half-unit intervals, since it is well-known that if an analytic function ( $z_{1}^{'}$ and $z_{2}^{'}$ in our case) takes the value zero in a non-empty interval, then this function has to be identically equal to 0 on all its domain. Therefore, instead of requiring that $z_{1}^{'}$ and $z_{2}^{'}$ take the value 0 in alternating half-unit intervals, we require that these functions take values very close to zero. Since the map $g_{M}$ given by Theorem 9 is robust to perturbations on its input, this will ensure that the whole simulation of $g_{M}$ with a two-dimensional ODE can still be performed, even if $z_{1}$ and $z_{2}$ are not strictly constant in the half-intervals $[k + 1 / 2, k + 1]$ and $[k, k + 1 / 2]$ , respectively. However, we still have to ensure that $| z_{1}^{'} (t) |$ and $| z_{2}^{'} (t) |$ are sufficiently small in the half-unit intervals of interests to guarantee that the iteration can be carried faithfully. To better understand how this can be achieved, we have to analyze the effects of introducing perturbations in (14), with the help of Lemma 13 since now the “targets”will be slightly perturbed.

Proceeding as in [18], the non-analytic function $θ_{j} (sin 2 π t)$ in the first equation of (14) is replaced by an analytic periodic function with period 1 which is close to zero when $t \in [1 / 2, 1]$ . As shown in [18], this can be done by considering the function s defined by $\begin{matrix} (15) & s (t) = \frac{1}{2} ({sin}^{2} (2 π t) + sin (2 π t)) . \end{matrix}$ which ranges between 0 and 1 in $[0, 1 / 2]$ (and, in particular, between $4 / 5$ and 1 when $x \in [0.17, 0.33]$ ), and between $- \frac{1}{8}$ and 0 on the time interval $[1 / 2, 1]$ . Then we take the analytic function $ϕ : R^{2} \to [0, 1]$ defined by $\begin{matrix} ϕ (t, y) = Ψ (s (t), y), \end{matrix}$ we conclude that $\int_{0}^{1 / 2} ϕ (t, y) d t > 4 / 5 \times (0.33 - 0.17) = 0.128 > 0$ (assuming that $y ⩾ 5$ ) and $| ϕ (t, y) | < e^{- y} / 8$ for all $t \in [1 / 2, 1]$ (i.e. y allows us to provide an error bound for $z_{1}^{'} (t)$ in the time interval $[1 / 2, 1]$ ). Since ϕ has period 1 on t, we conclude that $\int_{k}^{k + 1 / 2} ϕ (t, y) d t > 0.128 > 0$ and $| ϕ (t, y) | < e^{- y} / 8$ for all $t \in [k + 1 / 2, k + 1]$ , where $k \in N$ is arbitrary and $y ⩾ 5$ . Therefore ϕ satisfies the assumptions of the function ϕ in Lemma 13 on the time interval $[0, 1 / 2]$ .

We can now proceed with the main construction that simulates the iteration of the map given by Theorem 9 with an analytic ODE. Take $γ > 0$ to be a value such that $2 γ + δ / 2 ⩽ 1 / 5$ , and let $g_{M}$ be the map given by Theorem 9 (use as value for δ in the statement of the Theorem 9 the value $η / 2 < 1 / 5$ , where $η = (γ + δ) / 2 + 1 / 5 < 2 / 5$ ). Consider the ODE $z^{'} = h_{M} (t, z)$ given by $\begin{array}{l} (16) & z_{1}^{'} = c_{1} {(z_{1} - σ^{[l]} \circ g_{M} \circ σ^{[l]} (z_{2}))}^{3} ϕ_{1} (t, z_{1}, z_{2}), \\ z_{2}^{'} = c_{2} {(z_{2} - σ^{[l]} (z_{1}))}^{3} ϕ_{2} (t, z_{1}, z_{2}), \end{array}$ where $z (t) = (z_{1} (t), z_{2} (t)) \in R^{2}$ , with initial conditions $z_{1} (0) = {\bar{x}}_{0}$ , $z_{2} (0) = {\bar{y}}_{0}$ , where ${\bar{x}}_{0}, {\bar{y}}_{0} \in R$ are approximations of some initial configuration $x_{0}$ satisfying $| {\bar{x}}_{0} - x_{0} | ⩽ 1 / 5$ and $| {\bar{y}}_{0} - x_{0} | ⩽ 1 / 5$ , $l \in N$ is such that $σ^{[l]} (η) ⩽ γ < 1 / 5$ (see Lemma 5), and $\begin{array}{l} ϕ_{1} (t, z_{1}, z_{2}) = ϕ (t, \frac{c_{1}}{γ} {(z_{1} - σ^{[l]} \circ g_{M} \circ σ^{[l]} (z_{2}))}^{4} + \frac{c_{1}}{γ} + 10), \\ (17) & ϕ_{2} (t, z_{1}, z_{2}) = ϕ (- t, \frac{c_{2}}{γ} {(z_{2} - σ^{[l]} (z_{1}))}^{4} + \frac{c_{2}}{γ} + 10), \end{array}$ and $c_{1}$ , $c_{2}$ are constants associated to a targeting error of value γ (i.e. they are chosen such as the constant c in (10) and by noting that $\int_{0}^{1 / 2} ϕ (t, y) > 0.128$ , we can take $c_{1}, c_{2} = 4 γ^{- 2}$ ). Note that, since ϕ satisfies the assumptions of function ϕ in Lemma 13, the same will happen to $ϕ_{1}$ and $ϕ_{2}$ . The ODE (16) simulates the Turing machine M similarly to the ODE (14), by iterating $g_{M}$ . However (17) has a more complicated expression, since it has to deal with the fact that $z_{1}^{'} (t)$ and $z_{2}^{'} (t)$ are not exactly zero in half-unit intervals.

Since we want that the ODE $z^{'} = h_{M} (t, z)$ to robustly simulate M, let us assume that the right hand-side of the equations in (16) can be subject to an error of absolute value not exceeding δ. To start the analysis of the behavior of (16), let us first consider the time interval $[0, 1 / 2]$ . Since $| x |^{3} ⩽ x^{4} + 1$ for all $x \in R$ , we conclude that $ϕ_{2}$ is less than $min (γ {(c_{2} | z_{2} - σ^{[n_{2}]} (z_{1}) |^{3})}^{- 1}, 1 / 10)$ . This implies, together with the assumption that $z_{2}^{'}$ in (17) is perturbed by an amount not exceeding δ, that $| z_{2}^{'} (t) | ⩽ γ + δ$ for all $t \in [0, 1 / 2]$ which implies that $| z_{2} (t) - z_{2} (0) | ⩽ (γ + δ) / 2$ for all $t \in [0, 1 / 2]$ . Since the initial condition ${\bar{x}}_{0}$ satisfies $| {\bar{x}}_{0} - x_{0} | ⩽ 1 / 5$ where $x_{0} \in N$ , we conclude that $\begin{matrix} | z_{2} (t) - x_{0} | ⩽ \frac{γ + δ}{2} + 1 / 5 = η < \frac{2}{5} for all t \in [0, 1 / 2] \end{matrix}$ which implies that $\begin{matrix} | σ^{[l]} (z_{2} (t)) - x_{0} | ⩽ γ < \frac{1}{5} for all t \in [0, 1 / 2] . \end{matrix}$ Due to Theorem 9, we conclude that $\begin{matrix} | g_{M} \circ σ^{[l]} (z_{2} (t)) - ψ (x_{0}) | ⩽ \frac{1}{5} for all t \in [0, 1 / 2] . \end{matrix}$ Since $σ^{[l]} (η) ⩽ γ$ and $η > 1 / 5$ , we get that $| σ^{[l]} \circ g_{M} \circ σ^{[l]} (z_{2} (t)) - g_{M} (x_{0}) | < γ$ for all $t \in [0, 1 / 2]$ . Then the behavior of $z_{1}$ is given by Lemma 13 and $\begin{matrix} (18) & | z_{1} (\frac{1}{2}) - ψ (x_{0}) | < 2 γ + δ / 2 ⩽ \frac{1}{5} . \end{matrix}$

In the next half-unit interval $[1 / 2, 1]$ the roles of $z_{1}$ and $z_{2}$ are switched as before and one concludes, using similar arguments to those used in the time interval $[0, 1 / 2]$ that $| z_{1}^{'} (t) | ⩽ γ + δ$ for $t \in [1 / 2, 1]$ . Therefore $| z_{1} (t) - z_{1} (1 / 2) | ⩽ (γ + δ) / 2$ for all $t \in [1 / 2, 1]$ . This inequality together with (18) yields that $\begin{matrix} | z_{1} (t) - ψ (x_{0}) | < 1 / 5 + (γ + δ) / 2 = η for all t \in [1 / 2, 1] . \end{matrix}$ Since $σ^{[l]} (η) ⩽ γ$ , we conclude that $| σ^{[l]} (z_{1} (t)) - ψ (x_{0}) | ⩽ γ$ and Lemma 13 gives us $\begin{matrix} (19) & | z_{2} (1) - ψ (x_{0}) | < 2 γ + δ / 2 ⩽ \frac{1}{5} . \end{matrix}$ On the time interval $[1, 3 / 2]$ , the roles of $z_{1}$ and $z_{2}$ are again switched. There we will have $| z_{2}^{'} (t) | ⩽ γ + δ$ for all $t \in [1, 3 / 2]$ which implies that $| z_{2} (t) - z_{2} (1) | ⩽ (γ + δ) / 2$ for all $t \in [1 / 2, 1]$ . This inequality together with (19) gives us $\begin{matrix} | z_{2} (t) - ψ (x_{0}) | < 1 / 5 + (γ + δ) / 2 = η for all t \in [1, 3 / 2] . \end{matrix}$ We thus conclude that this analysis can be repeated in subsequent time intervals and thus that for all $j \in N$ , if $t \in [j, j + \frac{1}{2}]$ then $| z_{2} (t) - ψ^{[j]} (x_{0}) | ⩽ 1 / 5$ . This concludes the proof.
Remark 15.
From the proof of Theorem 11, we conclude that if the Turing machine halts after $n_{0}$ steps with configuration $c_{h}$ and if we assume that $ψ (c_{h}) = c_{h}$ , then for any $n \in N$ satisfying $n ⩾ n_{0} + 1$ , we have $\begin{matrix} | z_{2} (n) - c_{h} | ⩽ \frac{1}{5} . \end{matrix}$ Furthermore, on the half-unit time interval $[n, n + 1 / 2]$ , we will get that $\begin{matrix} | σ^{[l]} (z_{1} (t)) - c_{h} | < γ < 1 / 5 \end{matrix}$ for $t \in [n, n + 1 / 2]$ . Assuming that $δ = 0$ in Theorem 11, we then conclude from (16) that $z_{2} (t)$ starts its trajectory from a value $1 / 5$ -near to $c_{h}$ and will monotically converge to a value $σ^{[l]} (z_{1} (t))$ which, although may change with time, will never leave an $1 / 5$ -vicinity of $c_{h}$ for $t \in [n, n + 1 / 2]$ . This implies that $\begin{matrix} (20) & | z_{2} (t) - c_{h} | ⩽ \frac{1}{5} \end{matrix}$ for all $t \in [n, n + 1 / 2]$ . Since (20) holds on $[n, n + 1 / 2]$ as we have seen from the proof of Theorem 11, we conclude that (20) holds for $t ⩾ n_{0} + 1 / 2$ .
6. Simulating Turing machines over compact sets

So far all our simulations of Turing machines are carried out over the non-compact set $R^{n}$ . In this section we show that it is possible to simulate Turing machines with ODEs on the 2-dimensional compact sphere $S^{2} = {x \in R^{3} : ‖ x ‖ = 1}$ . The technique used in the construction is based on a similar result proved in [14]. It is shown in [14] that there is a computable polynomial vector field on the sphere $S^{n} = {x \in R^{n + 1} : ‖ x ‖ = 1}$ simulating a universal Turing machine; in other words, the polynomial vector field is Turing complete, where $n ⩾ 17$ is arbitrary. The underlying idea of the construction presented in [14] is to map a vector field defined in $R^{n}$ to $S^{n}$ using the stereographic projection, and to remove the singularity at the north pole by using a suitable reparametrization on the polynomial vector fields. We show in this section that the dimension can be lowered from $n ⩾ 17$ to $n = 2$ at the cost of using a non-polynomial $C^{\infty}$ vector field. This is done by using the stereographic projection map applied to the function $g_{M}$ of Theorem 11. However, since $g_{M}$ is not polynomial, a different technique from that of [14] will have to be used to remove the singularity at the north pole, and that technique only works for $C^{\infty}$ vector fields. It would be interesting to know if the simulation of a Turing machine on $S^{2}$ can be achieved with an analytic vector field.

The notion introduced in the following definition will be used throughout this section.

Definition 16.
Let $f : R^{n + k} \to R^{j}$ , where $n, k, j \in N_{0}$ , with $j ⩾ 1$ . Given an expression $f (x_{1}, \dots, x_{n}, y_{1}, \dots, y_{k})$ , we say that f is bounded by a constant on the variables $x_{1}, \dots, x_{n}$ and bounded by a polynomial on the variables $y_{1}, \dots, y_{k}$ if $\begin{matrix} ‖ f (x_{1}, \dots, x_{n}, y_{1}, \dots, y_{k}) ‖ ⩽ p (y_{1}, \dots, y_{k}) \end{matrix}$ for all $(x_{1}, \dots, x_{n}, y_{1}, \dots, y_{k}) \in R^{n + k}$ in the domain of f.

For example $g (t, x) = x sin (2 π t)$ is bounded by a constant on t and polynomially bounded on x and $h (t) = sin (2 π t)$ is bounded by a constant on t (to simplify notation, in this latter case where h is not polynomially bounded on any other variables, we will just say that h is bounded by a constant). Some functions which are bounded by a constant include sin, cos, arctan.

The next result shows when a Turing universal vector field can be defined on $S^{n}$ .
Theorem 17.
Let $\begin{matrix} (21) & y^{'} = f (t, y) \end{matrix}$ be a $C^{1}$ -computable ODE simulating a Turing machine M, where $f : R^{n + 1} \to R^{n}$ is a $C^{\infty}$ function with the property that both f and all its partial derivatives $\frac{\partial^{| k |} f}{\partial t^{k_{0}} \partial y_{1}^{k_{1}} \dots \partial y_{n}^{k_{n}}}$ are bounded by a constant on t and are polynomially bounded on the variables $x_{1}, \dots, x_{n}$ , assuming that the argument of f is $(t, x) \in R^{n + 1}$ , where $k_{0}, k_{1}, \dots, k_{n} \in N_{0}$ , $k = (k_{0}, k_{1}, \dots, k_{n})$ , and $| k | = k_{0} + k_{1} + \dots + k_{n}$ . Then from f one can compute a $C^{\infty}$ vector field F defined on $S^{n}$ that also simulates M.
Proof.
We recall that the (inverse) stereographic projection $φ : R^{n} \to S^{n} \subseteq R^{n + 1}$ is given by $\begin{matrix} φ (x_{1}, x_{2}, \dots, x_{n}) = (\frac{r^{2} - 1}{1 + r^{2}}, \frac{2 x_{1}}{1 + r^{2}}, \frac{2 x_{2}}{1 + r^{2}}, \dots, \frac{2 x_{n}}{1 + r^{2}}) \end{matrix}$ where $r^{2} = x_{1}^{2} + x_{2}^{2} + \dots + x_{n}^{2}$ . Suppose that $f (t, x) = (f_{1} (t, x), \dots, f_{n} (t, x))$ . Then we can write the vector field defined by f on $R^{n}$ as $\begin{matrix} f = \sum_{i = 1}^{n} f_{i} \frac{\partial}{\partial x_{i}} . \end{matrix}$ We recall that if $g : M \to N$ is a $C^{1}$ map between two manifolds $M = R^{n}$ and $N = R^{k}$ , then for each $p \in M$ the map g induces a linear map $g_{} : T_{p} M \to T_{g (p)} N$ from the tangent space $T_{p} M$ of M at p to the tangent space of N at $g (p)$ . We also recall that $(\partial / \partial x_{1} |_{p}, \dots, \partial / \partial x_{n} |_{p})$ forms a basis for $T_{p} M$ and, similarly, if $({\bar{x}}_{1}, \dots, {\bar{x}}_{k})$ are coordinates for $N = R^{k}$ , then $(\partial / \partial {\bar{x}}_{1} |_{g (p)}, \dots, \partial / \partial {\bar{x}}_{k} |_{g (p)})$ forms a basis for $T_{g (p)} N$ . Moreover, the matrix that (locally) defines the linear map $g_{}$ , relative to the bases $(\partial / \partial x_{1} |_{p}, \dots, \partial / \partial x_{n} |_{p})$ and $(\partial / \partial {\bar{x}}_{1} |_{g (p)}, \dots, \partial / \partial {\bar{x}}_{k} |_{g (p)})$ , is the Jacobian of g. In the case of the map φ, and if we take $(y_{0}, y_{1}, \dots, y_{n})$ as coordinates for $R^{n + 1}$ , we obtain the following (note that the variable t in the expression of f can be seen as a fixed parameter): $\begin{array}{l} φ_{} (f) & = \sum_{i = 1}^{n} f_{i} φ_{} (\frac{\partial}{\partial x_{i}}) \\ = \sum_{i = 1}^{n} f_{i} \sum_{j = 0}^{n} \frac{\partial φ_{j}}{\partial x_{i}} \frac{\partial}{\partial y_{j}} \\ (22) & = \sum_{i = 1}^{n} f_{i} \cdot ((1 - y_{0}) y_{i} \frac{\partial}{\partial y_{0}} + (1 - y_{0} - y_{i}^{2}) \frac{\partial}{\partial y_{i}} - \sum_{\begin{array}{c} j = 0 \\ j \notin {0, i} \end{array}}^{n} y_{i} y_{j} \frac{\partial}{\partial y_{j}}), \end{array}$ where $f_{i}$ is evaluated at $(t, φ^{- 1} (y_{0}, y_{1}, \dots, y_{n})) = (t, \frac{y_{1}}{1 - y_{0}}, \dots, \frac{y_{n}}{1 - y_{0}})$ . This implies that $φ_{} (f)$ is a vector field of class $C^{\infty}$ , except at the north pole $y_{N P} = (1, 0, \dots, 0)$ of $S^{n}$ where it is not defined. Let us now consider the following ODE $\begin{array}{l} (23) & τ^{'} = K (\bar{x}) \\ {\bar{x}}^{'} = f (τ, \bar{x}) K (\bar{x}) \end{array}$ where $K : R^{n} \to R$ is such that $K (x) > 0$ for any $x \in R^{n}$ and $τ (0) = 0$ , $\bar{x} (0) = x_{0}$ . Since $τ^{'} (t) > 0$ , τ is strictly increasing and thus τ admits an inverse $τ^{- 1}$ . Next we note that if $\tilde{x} = x \circ τ$ , where x is a solution of (21) and $τ (0) = 0$ , we have that $\begin{array}{l} {\tilde{x}}^{'} (t) & = {(x (τ (t)))}^{'} \\ = f (τ (t), x (τ (t))) τ^{'} (t) \\ (24) & = f (τ (t), \tilde{x} (t)) K (\bar{x} (t)) . \end{array}$ It is not difficult to see that if $(τ, \bar{x})$ is a solution for (23), then $\bar{x}$ is also a solution to (24). Since the solution of the ODE (23) is unique, by the Picard-Lindelöf theorem, we conclude that $\bar{x} (t) = \tilde{x} (t) = x (τ (t))$ . Furthermore, since $τ^{'} (t) > 0$ for any $t \in R$ , we conclude that any solution curve of (21) with initial condition $y (0) = x_{0}$ also provides a solution curve for the last n components of the solution of (23) with initial condition $τ (0) = 0$ , $\bar{x} (0) = x_{0}$ , up to some time reparametrization, and vice versa.

Thus, by taking $K (x) = e^{- \frac{2}{1 + r^{2}}} = e^{- \frac{2}{1 + x_{1}^{2} + \dots + x_{n}^{2}}}$ , we conclude that the solution curves of $\begin{array}{l} τ^{'} = e^{- \frac{2}{1 + r^{2}}} \\ (25) & x^{'} = e^{- \frac{2}{1 + r^{2}}} f (τ, x) = h (t, x) \end{array}$ and of (21) are the same, up to a time reparametrization τ given by the ODE $τ^{'} = e^{- \frac{2}{1 + r^{2}}} > 0$ . Note that the right-hand side of (25) is formed by $C^{1}$ -computable functions, which means that the solution $(τ, x)$ to (25) is also computable, since the solution to a $C^{1}$ -computable system is also computable [16,19]. Hence, when simulating Turing machines, if the result of the nth step of the computation of the Turing machine being simulated by (21) can be read in the time interval $[a_{n}, b_{n}]$ , then the result of the nth step when the simulation is performed by (25) can be read on the time interval $[τ^{- 1} (a_{n}), τ^{- 1} (b_{n})]$ . Note that $τ^{- 1}$ can be computed from τ and hence from f. Indeed, we know that the derivative of $τ^{- 1}$ is given by $\begin{array}{l} {(τ^{- 1} (a))}^{'} & = \frac{1}{τ^{'} (τ^{- 1} (a))} \\ = \frac{1}{K (\bar{x} (τ^{- 1} (a)))} \\ = \frac{1}{K (x (τ \circ τ^{- 1} (a)))} \\ \frac{1}{K (x (a))} . \end{array}$ Hence $τ^{- 1}$ can be obtained as the solution of the initial-value problem (IVP) defined by ${(τ^{- 1} (t))}^{'} = 1 / K (x (t))$ , $τ^{- 1} (0) = 0$ . Since the right-hand side of the ODE defining this IVP is computable from x and hence from f, we conclude that $τ^{- 1}$ is computable from f [16,19]. In particular, if f is computable, then so is $τ^{- 1}$ .

We now have (we can again assume that t is a fixed parameter for h, where h is given by (25)) $\begin{matrix} φ_{} (h) = e^{- (1 - y_{0})} φ_{} (f) . \end{matrix}$ From (22) we conclude that (note that $| y_{i} | ⩽ 1$ ) $\begin{array}{l} ‖ φ_{} (h) (y) ‖ = & ‖ e^{- (1 - y_{0})} φ_{} (f) (y) ‖ \\ = & ‖ e^{- (1 - y_{0})} ‖ ‖ φ_{} (f) (y) ‖ \\ = & ‖ e^{- (1 - y_{0})} ‖ ‖ \sum_{i = 1}^{n} f_{i} (t, \frac{y_{1}}{1 - y_{0}}, \dots, \frac{y_{n}}{1 - y_{0}}) \cdot \\ \cdot ((1 - y_{0}) y_{i} \frac{\partial}{\partial y_{0}} + (1 - y_{0} - y_{i}^{2}) \frac{\partial}{\partial y_{i}} - \sum_{\begin{array}{c} j = 0 \\ j \notin {0, i} \end{array}}^{n} y_{i} y_{j} \frac{\partial}{\partial y_{j}}) ‖ \\ ⩽ & e^{- (1 - y_{0})} \cdot 6 n \cdot p (\frac{y_{1}}{1 - y_{0}}, \dots, \frac{y_{n}}{1 - y_{0}}) . \end{array}$ This latter result implies that $\begin{matrix} (26) & lim_{\begin{array}{c} y \to y_{N P} \\ y \in S^{n} - {y_{N P}} \end{array}} φ_{} (h) (y) = 0 . \end{matrix}$ Similarly, if we assume that $\begin{matrix} φ_{} (h) (y) = e^{- (1 - y_{0})} \sum_{i = 0}^{n} s_{i} (t, y) \frac{\partial}{\partial y_{i}} \end{matrix}$ where $s_{i}$ are functions, from (22) and from the assumption that all the partial derivatives of f of the form $\frac{\partial^{| k |} f}{\partial t^{k_{0}} \partial x_{1}^{k_{1}} \dots \partial x_{n}^{k_{n}}}$ are polynomially bounded on $x_{1}, \dots, x_{n}$ , we can conclude that each partial derivative of $s_{i}$ is polynomially bounded by some polynomial $\tilde{p} (\frac{y_{1}}{1 - y_{0}}, \dots, \frac{y_{n}}{1 - y_{0}})$ , which implies that $\begin{matrix} lim_{\begin{array}{c} y \to y_{N P} \\ y \in S^{n} - {y_{N P}} \end{array}} \frac{\partial^{| k |} φ_{} (h) (y)}{\partial t^{k_{0}} \partial y_{1}^{k_{1}} \dots \partial y_{n}^{k_{n}}} = 0 . \end{matrix}$ Therefore we can extend $φ_{} (h) (y)$ to a $C^{\infty}$ vector field $\tilde{g}$ defined on the entire sphere $S^{n}$ if we assume that the value of $\tilde{g}$ and of its partial derivatives is 0 at the north pole $y_{N P} = (1, 0, \dots, 0)$ of $S^{n}$ . We thus have defined a $C^{\infty}$ vector field $\tilde{g}$ on the entire sphere $S^{n}$ , and $\tilde{g}$ is Turing universal. □

With the help of the above theorem, we may hope we could just make use of the vector field $g_{M}$ of Theorem 11 as the vector field f in Theorem 17 to prove that there is a Turing universal vector field in $S^{2}$ . However, there are several problems with this approach as listed below: (1) the approach requires that $g_{M} (t, x_{1}, x_{2})$ and all of its partial derivatives are bounded by a constant on t and are polynomially bounded on $x_{1}$ and $x_{2}$ . A major problem in this respect is that $g_{M}$ uses in its expression the function Ψ defined in Lemma 4 that does not necessarily have polynomially bounded derivatives because the function arcsin is used in the definition of Ψ and the derivative of arcsin is not polynomially bounded as its argument approaches $- 1$ or 1. (2) The expression of $g_{M}$ relies on the expression of the function ${\bar{g}}_{M}$ given by Theorem 3. Thus one must show that ${\bar{g}}_{M} (t, x_{1}, x_{2}, x_{3}, y_{1}, y_{2}, y_{3})$ is bounded polynomially on $x_{1}$ , $x_{2}$ , $x_{3}$ , $y_{1}$ , $y_{2}$ , $y_{3}$ . (3) The argument of the functions $Ω_{k, i} : R \to R$ from Proposition 8 is provided as the initial condition of an ODE. It is not straightforward to analyze the dependence of $Ω_{k, i}$ and of its derivatives on its argument.

In the following we present the solutions to the listed problems. First we note that in Theorem 17 the vector field is no longer required to be analytic (it only has to be $C^{\infty}$ ). Therefore we can substitute the function Ψ defined in Lemma 4 by the function $r : R \to R$ defined in Lemma 14. We recall that the function r has the property that $r (x) = n$ , whenever $x \in [n - 1 / 4, n + 1 / 4]$ , for all integers n. Thus if we take $Ψ (x, y) = r (x)$ , the properties stated for Ψ in Lemma 4 remain true. Therefore, we can replace Ψ with r when defining the vector field $g_{M}$ of Theorem 11. In this case, the properties stated in Theorem 11 remain true with $g_{M}$ being a $C^{\infty}$ function rather than an analytic function. However, we gain the advantage that r and its derivatives are polynomially bounded as we shall show now. Indeed, it follows from its definition that $| r (x) | ⩽ | x | + 1 ⩽ x^{2} + 2$ (recall that $| x | ⩽ x^{2} + 1$ ). For the derivatives of θ, we note that if $a (x) = e^{- 1 / x}$ , then it is readily seen by induction that for each $n \in N_{0}$ there is a polynomial $P_{n}$ such that $a^{(n)} (x) = P_{n} (1 / x) e^{- 1 / x}$ , with $a^{(0)} (x) = a (x)$ (and thus $P_{0} = 1$ ). This implies that ${lim}_{x \to 0^{+}} a^{(n)} (x) = 0$ and ${lim}_{x \to + \infty} a^{(n)} (x) = K_{n}$ , where $K_{n}$ is the constant term of $P_{n}$ . Therefore, by definition of the limit, there exists $b_{n}, ε_{n} > 0$ such that $| a^{(n)} (x) | ⩽ 1$ whenever $x \in (0, ε_{n}]$ and $| a^{(n)} (x) | ⩽ K_{n} + 1$ whenever $x ⩾ b_{n}$ . Now set $M_{n} = {max}_{x \in [ε, M]} | a^{(n)} (x) | \in R$ . Then $| a^{(n)} (x) | ⩽ max (1, K_{n} + 1, M_{n})$ for all $x \in (0, + \infty)$ , which implies that θ in (12) as well as its derivatives are polynomially bounded on $[0, + \infty)$ . Consequently, the function v from (13) as well as the derivatives of v are polynomially bounded following Lemmas 18 and 19 to be presented in a moment. Lemma 14 together with Lemmas 18 and 19 then imply that r as well as its derivatives are polynomially bounded.

Some notations are in order for the statements of the next two lemmas. Given multi-indexes $α = (α_{1}, \dots, α_{n})$ , $β = (β_{1}, \dots, β_{n}) \in N_{0}^{n}$ and $x = (x_{1}, \dots, x_{n}) \in R^{n}$ , let $\begin{array}{l} | α | = α_{1} + \dots + α_{n} \\ α + β = (α_{1} + β_{1}, \dots, α_{n} + β_{n}) \\ α! = (α_{1}!) \cdot (α_{2}!) \cdot \dots \cdot (α_{n}!) \\ β ⩽ α iff β_{1} ⩽ α_{1}, \dots, β_{n} ⩽ α_{n} \\ (\binom{α}{β}) = (\binom{α_{1}}{β_{1}}) \dots (\binom{α_{n}}{β_{n}}) = \frac{α!}{β! (α - β)!} for β ⩽ α \\ D_{x}^{α} = \frac{\partial^{| α |}}{\partial x_{1}^{α_{1}} \dots \partial x_{n}^{α_{n}}} (D_{x}^{0} is the identity operator) \\ x^{α} = x_{1}^{α_{1}} \dots x_{n}^{α_{n}} \end{array}$
Lemma 18.
Suppose that $f, g : R^{n} \to R$ are $C^{\infty}$ functions which are bounded by a constant on the variables $x_{1}, \dots, x_{i}$ and are bounded by a polynomial on the variables $x_{i + 1}, \dots, x_{n}$ , as well as all their partial derivatives. Then:
The $C^{\infty}$ function $f \pm g : R^{n} \to R$ defined by $(f + g) (x) = f (x) + g (x)$ is bounded by a constant on the variables $x_{1}, \dots, x_{i}$ and bounded by a polynomial on the variables $x_{i + 1}, \dots, x_{n}$ , as well as all its partial derivatives.

The $C^{\infty}$ function $f \times g : R^{n} \to R$ defined by $(f \times g) (x) = f (x) \cdot g (x)$ is bounded by a constant on the variables $x_{1}, \dots, x_{i}$ and bounded by a polynomial on the variables $x_{i + 1}, \dots, x_{n}$ , as well as all its partial derivatives.

Proof.
Let $α = (α_{1}, \dots, α_{n})$ be a multi-index. For point 1, we note that $\begin{matrix} \frac{\partial^{| α |}}{\partial x^{α}} (f \pm g) (x) = \frac{\partial^{| α |} f (x)}{\partial x^{α}} \pm \frac{\partial^{| α |} g (x)}{\partial x^{α}} \end{matrix}$ and thus the result follows immediately from the assumption.

For the product $f \times g$ , the claim follows directly from the general Leibniz rule for multivariate functions (see e.g. [17, Proof of Lemma 2.6]): $\begin{matrix} \frac{\partial^{| α |}}{\partial x^{α}} (f \times g) (x) = \sum_{0 ⩽ β ⩽ α} (\binom{α}{β}) \frac{\partial^{| β |} f (x)}{\partial x^{β}} \frac{\partial^{| α - β |} g (x)}{\partial x^{α - β}} . \end{matrix}$ □
Lemma 19.
Suppose that $f : R^{j} \to R^{n}$ and $g : R^{k} \to R^{j}$ are $C^{\infty}$ functions with the following properties:
f and its partial derivatives are polynomially bounded on its arguments $z_{1}, \dots, z_{j}$ ;

g and its partial derivatives are bounded by a constant on the variables $x_{1}, \dots, x_{i}$ and bounded by a polynomial on the variables $x_{i + 1}, \dots, x_{k}$ .

Then $f \circ g$ is a $C^{\infty}$ function with the property that $f \circ g$ as well as all its partial derivatives are bounded by a constant on the variables $x_{1}, \dots, x_{i}$ and by a polynomial on the variables $x_{i + 1}, \dots, x_{k}$ .
Proof.
To prove this theorem, we will use a multivariate version of the Faà di Bruno formula which allows us to compute the higher order partial derivatives of the composition of multivariate functions f and g. Some notational matter is in order first. Let $α = (α_{1}, \dots, α_{n}) \in N_{0}^{n}$ be a multi-index. Then following the approach of [27], let us assume that $χ^{0} = 1$ regardless of whether χ is a number or a differential operator. We say that a multi-index α can be decomposed into s parts $p_{1}, \dots, p_{s} \in N_{0}^{n}$ with multiplicities $m_{1}, \dots, m_{s} \in N_{0}^{d}$ if the decomposition $α = | m_{1} | p_{1} + \dots + | m_{s} | p_{s}$ holds and all parts are different. In this case the total multiplicity is defined as $m = m_{1} + \dots + m_{s}$ . The list $(s, p, m)$ is called a d-decomposition, or simply just a decomposition, of α. Then, assuming that $z = f \circ g (x)$ and $y = g (x)$ , we have $\begin{matrix} \frac{\partial^{| α |}}{\partial x^{α}} f \circ g (x) = α! \sum_{(s, p, m) \in D} \frac{\partial^{| m |}}{\partial y^{m}} f (g (x)) \prod_{k = 1}^{s} \frac{1}{m_{k}!} (\frac{1}{p_{k}!} \frac{\partial^{| p_{k} |}}{\partial x^{p_{k}}} g (x)) f \circ g (x) \end{matrix}$ where $D$ is the set of all decompositions of α. From this later formula and the hypothesis on f, g we conclude that $f \circ g$ is a $C^{\infty}$ function with the property that $f \circ g$ as well as all its partial derivatives are bounded by a constant on the variables $x_{1}, \dots, x_{i}$ and are bounded by a polynomial on the variables $x_{i + 1}, \dots, x_{k}$ . □

Since the function $f_{M}$ of Theorem 1 can be written using only the following terms: variables, polynomial-time computable constants, +, −, ×, sin, cos, arctan, we conclude from Lemmas 18 and 19 that $f_{M}$ as well as all its partial derivatives are bounded by a polynomial. Furthermore, as the function $g_{M}$ from Theorem 9 is obtained using the functions $Υ_{3}$ , $f_{M}$ , $Ω_{3, 1}$ , $Ω_{3, 2}$ , $Ω_{3, 3}$ , σ, again by Lemmas 18 and 19 it suffices to show that $Υ_{3}$ , $Ω_{3, 1}$ , $Ω_{3, 2}$ , $Ω_{3, 3}$ , σ as well as their partial derivatives are (or can be made) bounded by a polynomial. The case of the function σ in Lemma 5 is immediate from Lemmas 18 and 19. The case of $Υ_{3}$ can be treated by replacing $Ψ (x, y)$ in the expression of $Υ_{3}$ given in the proof of Proposition 7 by the $C^{\infty}$ function r, as explained above. Then it follows immediately that $Υ_{3}$ maintains its properties given by Proposition 7, except analiticity ( $Υ_{3}$ will only be $C^{\infty}$ ) and, meanwhile, Lemmas 18 and 19 imply that $Υ_{3}$ and its partial derivatives are polynomially bounded.

The situation for $Ω_{3, 1}$ , $Ω_{3, 2}$ , $Ω_{3, 3}$ is more subtle, as the arguments to these functions are passed as initial conditions of an ODE. What we are going to do is to create new $C^{\infty}$ functions ${\overline{Ω}}_{3, 1}$ , ${\overline{Ω}}_{3, 2}$ , ${\overline{Ω}}_{3, 3}$ such that these new functions maintain the useful properties of $Ω_{3, 1}$ , $Ω_{3, 2}$ , $Ω_{3, 3}$ on the one hand and, on the other hand, the new functions together with their partial derivatives are polynomially bounded. Once this is done, the old functions $Ω_{3, 1}$ , $Ω_{3, 2}$ , $Ω_{3, 3}$ can then be replaced by the new functions ${\overline{Ω}}_{3, 1}$ , ${\overline{Ω}}_{3, 2}$ , ${\overline{Ω}}_{3, 3}$ in the expression of $g_{M}$ in the proof of Theorem 9. The function $g_{M}$ as well as its partial derivatives are now polynomially bounded on all variables except the time t. But since the time variable only appears inside the functions $ϕ_{1}$ and $ϕ_{2}$ defined by (17) in the format of $sin (2 π t)$ as defined in (15), and $sin (2 π t)$ and its derivatives are obviously bounded by a constant, it follows that the theorem below holds true.
Theorem 20.
Let $ψ : N \to N$ be the transition function of a Turing machine M, under the encoding described in Section 4 . Then there exist $0 < η < 2 / 5$ and a $C^{\infty}$ function $g_{M} : R^{3} \to R^{2}$ such that the ODE $z^{'} = g_{M} (t, z)$ simulates M in the following sense: for all $x_{0} \in N$ which encodes a configuration according to the encoding described above, if ${\bar{x}}_{0}, {\bar{y}}_{0} \in R$ satisfy the conditions $‖ {\bar{x}}_{0} - x_{0} ‖ ⩽ 1 / 5$ and $‖ {\bar{y}}_{0} - x_{0} ‖ ⩽ 1 / 5$ , then the solution $z (t)$ of $\begin{matrix} z^{'} = g_{M} (t, z), z (0) = ({\bar{x}}_{0}, {\bar{y}}_{0}) \end{matrix}$ satisfies, for all $j \in N_{0}$ and for all $t \in [j, j + 1 / 2]$ , $\begin{matrix} ‖ z_{2} (t) - ψ^{[j]} (x_{0}) ‖ ⩽ η, \end{matrix}$ where $z (t) \equiv (z_{1} (t), z_{2} (t)) \in R^{2}$ . Furthermore $g_{M} (t, z_{1}, z_{2})$ and its partial derivatives are polynomially bounded on $z_{1}$ and $z_{2}$ and bounded by a constant on t.
Theorem 21.
Let M be a Turing machine. Then one can compute from f a $C^{\infty}$ vector field F defined on $S^{2}$ which also simulates M.
Proof.
Immediate from Theorems 17 and 20. □

It remains to show that the new $C^{\infty}$ functions ${\overline{Ω}}_{3, 1}$ , ${\overline{Ω}}_{3, 2}$ , ${\overline{Ω}}_{3, 3}$ can be constructed such that they retain the useful properties of $Ω_{3, 1}$ , $Ω_{3, 2}$ , $Ω_{3, 3}$ , but these new functions and their partial derivatives are polynomially bounded. We begin with a preliminary lemma.
Lemma 22.
Let $f : R^{n + 1} \to R^{n}$ be a $C^{\infty}$ function such that f and its partial derivatives are polynomially bounded. Suppose that $x_{f} : R \to R$ is the first coordinate of a solution x of the IVP $\begin{matrix} x^{'} = f (t, x) . \end{matrix}$ If x is polynomially bounded, then all the derivatives of $x_{f}$ are polynomially bounded.
Proof.
We show the result by induction on the order of the derivative by showing that $x^{(k)} = f^{k} (t, x (t))$ , where $f^{k}$ is a $C^{\infty}$ function which is polynomially bounded, as well as all its partial derivatives. Then we will conclude that $x^{(k)}$ is polynomially bounded as well as all its partial derivatives by Lemma 19. The base case $\begin{matrix} x^{'} (t) = f (t, x (t)) \end{matrix}$ is trivial. Let us now assume that $x^{(k)} = f^{k} (t, x (t))$ , where $f^{k}$ is polynomially bounded as well as all its partial derivatives. Then $\begin{array}{l} x_{i}^{(k + 1)} (t) & = {(f_{i}^{k} (t, x (t)))}^{'} \\ = \frac{\partial f_{i}^{k}}{\partial t} (t, x (t)) + \sum_{j = 1}^{n} \frac{\partial f_{i}^{k}}{\partial y_{j}} (t, x (t)) \frac{d x_{j}}{d t} (t) \\ = \frac{\partial f_{i}^{k}}{\partial t} (t, x (t)) + \sum_{j = 1}^{n} \frac{\partial f_{i}^{k}}{\partial y_{j}} (t, x (t)) f_{j}^{k} (t, x (t)) \\ = f_{i}^{k + 1} (t, x (t)) . \end{array}$ By taking $f^{k + 1} = (f_{1}^{k + 1}, \dots, f_{n}^{k + 1})$ , we conclude that $x^{(k + 1)} = f^{k + 1} (t, x (t))$ and by Lemmas 18 and 19 we conclude that $f^{k + 1}$ is polynomially bounded as well as all its partial derivatives, thus showing the result. □

To define new $C^{\infty}$ functions ${\overline{Ω}}_{3, 1}$ , ${\overline{Ω}}_{3, 2}$ , ${\overline{Ω}}_{3, 3}$ as in Proposition 8, we recall that the key point of the proof of Proposition 8 was to consider the bijection $I : N^{2} \to N$ given by (4) and to obtain real extensions of the components $J_{2, 1}$ and $J_{2, 2}$ which form the inverse function of I, i.e. $I^{- 1} (z) = (J_{2, 1} (z), J_{2, 2} (z))$ . Then the result would follow inductively when obtaining extensions $Ω_{k, i}$ from $J_{k, i}$ for $k > 2$ . Subsequently, the only required modification is to obtain suitable real extensions ${\overline{Ω}}_{2, 1}$ , ${\overline{Ω}}_{2, 2}$ of $J_{2, 1}$ , $J_{2, 2}$ , respectively. We shall demand that if $n \in N_{0}$ , then $| {\overline{Ω}}_{2, i} (z) - J_{2, i} (n) | ⩽ 1 / 5$ whenever $| z - n | ⩽ 1 / 4$ for $i = 1, 2$ , so that ${\overline{Ω}}_{2, i}$ has the same properties as of $Ω_{2, i}$ regarding Proposition 8, except that ${\overline{Ω}}_{2, i}$ is $C^{\infty}$ instead of analytic. We also require that ${\overline{Ω}}_{2, i}$ and its derivatives are polynomially bounded. To obtain ${\overline{Ω}}_{2, i}$ , we first construct a $C^{\infty}$ function ${\tilde{Ω}}_{2, i}$ with the property that $| {\tilde{Ω}}_{2, i} (z) - J_{2, i} (n) | ⩽ 1 / 4$ whenever $z \in [n, n + 1 / 2]$ for $i = 1, 2$ , and ${\tilde{Ω}}_{2, i}$ as well as its partial derivatives are polynomially bounded. By setting ${\overline{Ω}}_{2, i} (z) = σ \circ {\tilde{Ω}}_{2, i} (z + 1 / 4)$ , where σ is given by Lemma 5, we conclude from the above and from Lemmas 18 and 19 that ${\overline{Ω}}_{2, i}$ has the desired properties. Before defining ${\tilde{Ω}}_{2, 1}$ and ${\tilde{Ω}}_{2, 2}$ , we note that I is obtained by dovetailing by enumerating the pairs in the diagonals below from the (one element) diagonal starting on $(0, 0)$ and then moving to the next diagonal

In other words, we have $I (0, 0) = 0$ , $I (0, 1) = 1$ , $I (1, 0) = 2$ , $I (0, 2) = 3$ , and so on. We note that the sum of the coordinates in each diagonal is constant. From here we see that the graphs for $J_{2, 1}$ and $J_{2, 2}$ , provided in Figs 2 and 3, respectively, have certain regularities which will be explored to obtain ${\tilde{Ω}}_{2, 1}$ and ${\tilde{Ω}}_{2, 2}$ .

Figure 2.
Graph of the function $J_{2, 1}$ . Since the function is discrete, the image are only the blue points (the red line is given as a visualization helper).

Let us start with the case of ${\tilde{Ω}}_{2, 1}$ . We first analyze the behavior of $J_{2, 1}$ . Let us suppose that the argument z of $J_{2, 1} (z)$ codes a pair $(0, n)$ at the start of the diagonal with sum n. Then $J_{2, 1} (z) = 0$ , $J_{2, 1} (z + 1) = 1$ , …, $J_{2, 1} (z + n) = n$ , $J_{2, 1} (z + n + 1) = 0$ , $J_{2, 1} (z + n + 2) = 1$ , …Thus to simulate $J_{2, 1}$ we need to track the sum s of the diagonal, and increase it by one when $J_{2, 1} (z)$ reaches the value of s. On the next value $z + 1$ , we will have that $J_{2, 1} (z + 1)$ will take the value 0, and each time its argument z increases by one, $J_{2, 1} (z)$ also increases by one, until it reaches the (new) value of the sum of the diagonal and the cycle repeats itself. We will simulate this behavior with ODEs. Before showing how this can be done, consider the auxiliary function $ξ : R \to R$ defined as $\begin{matrix} ξ^{'} (x) = \{\begin{matrix} 0 & if x ⩽ 1 / 4 \\ c_{ξ} θ (- (x - 1 / 4) (x - 3 / 4)) & if 1 / 4 < x < 3 / 4 \\ 0 & if x ⩾ 3 / 4 \end{matrix} \end{matrix}$ where $ξ (x) = 0$ and $c_{ξ} = {(\int_{1 / 4}^{3 / 4} θ (- (x - 1 / 4) (x - 3 / 4)) d x)}^{- 1}$ . It is not difficult to see that $θ (- (x - 1 / 4) (x - 3 / 4)) > 0$ when $1 / 4 < x < 3 / 4$ . Hence we have that $ξ (x) = 0$ whenever $x ⩽ 1 / 4$ , $0 < ξ (x) < 1$ when $1 / 4 < x < 3 / 4$ and $ξ (x) = 1$ when $x ⩾ 3 / 4$ . Furthermore ξ is $C^{\infty}$ and ξ and all its derivatives are polynomially bounded due to Lemmas 18 and 19.

Figure 3.
Graph of the function $J_{2, 2}$ . Since the function is discrete, the image are only the blue points (the red line is given as a visualization helper).

Let us now present the ODE which will define ${\tilde{Ω}}_{2, 1}$ $\begin{matrix} (27) & \{\begin{matrix} x_{1}^{'} = \tilde{c} {(ξ (r (s_{2}) - r (x_{2})) (1 + r (x_{2})) - x_{1})}^{3} θ (sin 2 π t) \\ x_{2}^{'} = \tilde{c} {(r (x_{1}) - x_{1})}^{3} θ (- sin 2 π t) \\ s_{1}^{'} = \tilde{c} {(r (s_{2}) + ξ (r (x_{2}) + 1 - r (s_{2})) - s_{1})}^{3} θ (sin 2 π t) \\ s_{2}^{'} = \tilde{c} {(r (s_{1}) - s_{1})}^{3} θ (- sin 2 π t) \end{matrix} \end{matrix}$ with $x_{1} (0) = x_{2} (0) = s_{1} (0) = s_{2} (0)$ and $\tilde{c} = 206$ . The behavior of the ODE (27) is similar to the one of (14). The variable updates are done on alternating time intervals. The variable $x_{2}$ stores the value of the function $J_{2, 1}$ on time intervals with the format $[k, k + 1 / 2]$ , i.e. we will have $| z_{2} (t) - J_{2, 1} | ⩽ 1 / 4$ whenever $t \in [k, k + 1 / 2]$ , with $k \in N_{0}$ . The variable $s_{2}$ will give the current sum of the diagonal on time intervals with the format $[k, k + 1 / 2]$ . We first update the variables $x_{1}$ and $s_{1}$ on time intervals $[k, k + 1 / 2]$ to be able to use the “memorized”values of $x_{2}$ and $s_{2}$ when updating $x_{1}$ and $s_{1}$ . We note that $x_{1}$ must be increased by one unit from its previous value (stored on $x_{2}$ ) until it reaches the value of the sum of the diagonal, which is stored in $s_{2}$ . On that moment we will have $ξ (r (s_{2}) - r (x_{2})) = 0$ on the equation for $x_{1}^{'}$ (if the value of $x_{2}$ is less than the value of the sum of the diagonal stored in $s_{2}$ , then $ξ (r (s_{2}) - r (x_{2})) = 1$ and $x_{1}$ is incremented by one) and $x_{1}$ will be reset to the value 0 starting the cycle again. The analysis for $s_{1}$ is similar: its value will be essentially constant as long as $ξ (r (x_{2}) + 1 - r (s_{2})) = 0$ , which happens when $r (s_{2}) - 1 ⩾ r (x_{2})$ . When $r (x_{2}) = r (s_{2})$ , we will have $ξ (r (x_{2}) + 1 - r (s_{2})) = 1$ and $s_{1}$ will be incremented by one from its previous value. We can see that the ODE (27) behaves as desired and that we can take ${\tilde{Ω}}_{2, 1} (t) = x_{2} (t)$ . Since the right-hand sides of (27) are polynomially bounded as well as their derivatives and $(x_{1}, x_{2}, s_{1}, s_{2})$ is also polynomially bounded, implying by Lemma 22 that ${\tilde{Ω}}_{2, 1}$ and all its derivatives are polynomially bounded.

The case for $J_{2, 2}$ is similar. Let us suppose that the argument z of $J_{2, 2} (z)$ codes a pair $(0, n)$ at the start of the diagonal with sum n. Then $J_{2, 2} (z) = n, J_{2, 2} (z + 1) = n - 1, \dots, J_{2, 2} (z + n) = 0, J_{2, 2} (z + n + 1) = n + 1, J_{2, 2} (z + n + 2) = n, \dots$ Thus to simulate $J_{2, 2}$ we need again to track the sum s of the diagonal, but now we need to increase it by one when $J_{2, 1} (z)$ reaches the value 0. On the next value $z + 1$ , we will have that $J_{2, 1} (z + 1)$ will take the value $s + 1$ , and each time its argument z increases by one, $J_{2, 1} (z)$ decreases by one, until it reaches 0 and the cycle repeats itself. This behavior can be simulated in a similar way to (27) by the following ODE $\begin{matrix} (28) & \{\begin{matrix} x_{1}^{'} = \tilde{c} {(ξ (x_{2}) (r (x_{2}) - 1) + ξ (1 - x_{2}) (1 + r (s_{2})) - x_{1})}^{3} θ (sin 2 π t) \\ x_{2}^{'} = \tilde{c} {(r (x_{1}) - x_{1})}^{3} θ (- sin 2 π t) \\ s_{1}^{'} = \tilde{c} (ξ (x_{2}) {(r (s_{2}) + ξ (1 - x_{2}) (r (s_{2}) + 1) - s_{1})}^{3} θ (sin 2 π t) \\ s_{2}^{'} = \tilde{c} {(r (s_{1}) - s_{1})}^{3} θ (- sin 2 π t) \end{matrix} \end{matrix}$

We can do an analysis to (28) similar to the one of (27) to conclude that we can take ${\tilde{Ω}}_{2, 2} (t) = x_{2} (t)$ on (28). This concludes the proof of Theorem 20.
7. Can one-dimensional ODEs simulate Turing machines?

As we have seen in the previous section, analytic two-dimensional ODEs can robustly simulate Turing machines. But what about one-dimensional ODEs? In this section we show that no one-dimensional autonomous ODE can simulate a universal Turing machine under some reasonable conditions.

First let us give a more precise meaning to the notion of an ODE simulating a Turing machine. Let M be a Turing machine. Since ODEs are defined on $R^{k}$ , to simulate the Turing machine M with an ODE we first need to encode a configuration of M as a point of $R^{k}$ . However, since the coding of a configuration might not be unique, as it happens in the previous sections, we map each configuration to a set of possible encodings of that configuration. Hence we have to consider a map χ which maps configurations of M into non-empty subsets of $R^{k}$ . Then given a configuration $c_{M}$ , any point of $χ (c_{M})$ is assumed to represent the configuration $c_{M}$ . In this manner we can consider the case when $c_{M}$ is represented by a single point in $R^{k}$ (when $χ (c_{M})$ is a singleton) or when $c_{M}$ is represented by several points of $R^{k}$ . For example, in Theorem 9 we have assumed that any point in a $1 / 5$ -vicinity of $I_{3} (y_{1}, y_{2}, q)$ , where $y_{1}$ and $y_{2}$ are given by (2) represents the configuration $c_{M}$ which is encoded by $I_{3} (y_{1}, y_{2}, q)$ , i.e. $\begin{matrix} χ (c_{M}) = {x \in R : | x - I_{3} (y_{1}, y_{2}, q) | ⩽ 1 / 5} . \end{matrix}$ Note that it makes sense to assume that if $c_{M}$ and $c_{M}^{'}$ are distinct configurations, then $χ (c_{M}) \cap χ (c_{M}^{'}) = \emptyset$ . However, this assumption may be too weak, since even if $χ (c_{M}) \cap χ (c_{M}^{'}) = \emptyset$ nothing prevents e.g. that $χ (c_{M})$ and $χ (c_{M}^{'})$ are fractal (e.g. Cantor-like) sets which are intermingled and thus very hard to separate in practice. To avoid such undesirable instances we impose a natural separation-condition on χ so that $χ (c_{M})$ and $χ (c_{M}^{'})$ are separated by disjoint open subsets of $R^{k}$ . More precisely, let ${c_{i}}_{i \in N}$ denote all configurations of a given Turing machine M (recall that a Turing machine has at most countably many configurations). Then we assume that there are two computable maps $a : {N \to Q}^{k}$ , $r : N \to Q$ such that for all $i \in N$ , $χ (c_{i}) \subseteq \overline{B (a (i), r (i))} = {x \in R^{k} : ‖ x - a (i) ‖ ⩽ r (i)}$ and, moreover, if $i \neq j$ then $\overline{B (a (i), r (i))} \cap \overline{B (a (j), r (j))} = \emptyset$ .

Now we say that the ODE $\begin{matrix} (29) & y^{'} = f (y), \end{matrix}$ where $f : R^{k} \to R^{k}$ , simulates a Turing machine with the coding χ if given an arbitrary configuration $c_{0}$ of M and some point $y_{0} \in χ (c_{0})$ one has that the solution y to (29) with initial condition $y (0) = y_{0}$ satisfies $y (n) \in χ (ψ^{[n]} (c_{0}))$ for all $n \in N$ , where ψ is the transition function of M. In the following, we show that no one-dimensional ODE can simulate a universal Turing machine under the separation-condition.

Theorem 23.

Let M be a universal Turing machine. Then no ODE $y^{'} = f (y)$ can simulate M in the sense explained above, where $f : R \to R$ is a computable function with only isolated zeros.

Proof.

Let M be a universal Turing machine. We may assume that it has only one halting state and it cleans its tape immediately before it halts. This implies that M has only one halting configuration $c_{k}$ , $k \in N$ ; moreover, the problem of deciding whether M halts on input w, $w \in Σ^{*} = {0, 1}^{*}$ , is undecidable.

Let us assume, by contradiction, that there is an ODE (29) which simulates M, where $f : R \to R$ is a computable function. Then f must admit a zero in $B_{k} = \overline{B (a (k), r (k))} = [a_{k} - r_{k}, a_{k} + r_{k}]$ where $a_{k} = a (k)$ and $r_{k} = r (k)$ . Assume otherwise that this is not the case. Then since computable functions are continuous, it must be either $f (x) < 0$ for all $x \in B_{k}$ or $f (x) > 0$ for all $x \in B_{k}$ . Moreover, since $B_{k}$ is compact, it follows that ${min}_{x \in B_{k}} | f (x) | = δ > 0$ . As a result, any solution starting on $B_{k}$ must leave it in time $⩽ 2 r_{k} / δ$ and never return to $B_{k}$ afterwards (note that a solution of (29) is a continuous function which must move continuously along the real line). But this is impossible because $ψ^{[n]} (c_{k}) = c_{k}$ for all $n \in N$ and (29) simulates M. Hence $B_{k}$ must contain at least one zero $x_{k}$ , which is computable because f is computable and the zeros of f are isolated (it is well-known that isolated zeros of computable functions are computable. See e.g. [10, Theorem 7.8]).

Let w be some input with the property that M halts on w, and suppose that the initial configuration associated to w is $c_{i_{w}}$ . Then $c_{i_{w}} \neq c_{k}$ (note that in an universal Turing machine the initial state cannot be an halting state). Hence, if $y_{0} \in χ (c_{i_{w}})$ , then $y_{0} \notin B_{k}$ . Let us assume without loss of generality that $y_{0} < a_{k} - r_{k}$ . We note that the solution of the IVP (29) with $y (0) = y_{0}$ must reach $B_{k}$ because M halts on input w. This implies that $f (x) > 0$ for all $x \in [y_{0}, a_{k} - r_{k})$ . There are two cases to be considered:

$f (x) > 0$ for all $x \in (- \infty, a_{k} - r_{k})$ . In this case, let $I = (- \infty, a_{k} + r_{k}]$ .

$f (\bar{x}) = 0$ for some $\bar{x} \in (- \infty, a_{k} - r_{k})$ . In this case, we must have $\bar{x} < y_{0}$ , which implies that $d = min {| x - y_{0} | : x \in [\bar{x}, y_{0}] and f (x) = 0} > 0$ . By continuity of f there is some $\tilde{x} \in [\bar{x}, y_{0}]$ such that $f (\tilde{x}) = 0$ and $| \tilde{x} - y_{0} | = d$ . In this case, let $I = (\tilde{x}, a_{k} + r_{k}]$ (note that $\tilde{x}$ is computable because it is an isolated zero of a computable function).

If there is a word w such that M halts on input w with the property that there is some $z \in χ (c_{i_{w}})$ satisfying $a_{k} < z$ , we repeat the above procedure on the half line $[a_{k}, + \infty)$ obtaining $I : = I \cup [a_{k} + r_{k}, + \infty)$ provided that $f (x) < 0$ for all $x \in (a_{k} + r_{k}, + \infty)$ or $I : = I \cup [a_{k} + r_{k}, \tilde{x})$ provided that $f (x) = 0$ for some $x \in (a_{k} + r_{k}, + \infty)$ , where $\tilde{x}$ is obtained similarly as in the previous case.

From the arguments above, we conclude that M halts on word w iff $χ (c_{i_{w}}) \subseteq I$ . We will use this fact to show that the halting problem is decidable, a contradiction. Let us assume, without loss of generality that $I = ({\tilde{x}}_{1}, {\tilde{x}}_{2})$ (the cases where one or more extremities of I are unbounded is dealt with similarly). Suppose first that for all $i \in N$ , ${{\tilde{x}}_{1}, {\tilde{x}}_{2}} \cap B_{i} = \emptyset$ . Then to decide whether M halts on input w proceed as follows. Given the initial configuration $c_{i_{w}}$ associated to input w, test whether $a_{i_{w}} \in I$ by testing whether ${\tilde{x}}_{1} < a_{i_{w}} < {\tilde{x}}_{2}$ . Notice that, since ${{\tilde{x}}_{1}, {\tilde{x}}_{2}} \cap B_{i_{w}} = \emptyset$ , this test can be done in finite time. If the test suceeds, then accept w otherwise reject it.

Let us now suppose that ${\tilde{x}}_{1} \in B_{j_{1}}$ and ${\tilde{x}}_{2} \in B_{j_{2}}$ (the cases where (i) ${\tilde{x}}_{1} \in B_{j_{1}}$ and for all $i \in N$ , ${\tilde{x}}_{2} \notin B_{i}$ or (ii) ${\tilde{x}}_{2} \in B_{j_{2},}$ and for all $i \in N$ , ${\tilde{x}}_{1} \notin B_{i}$ can be treated similarly). Then we can take $j_{1}$ and $j_{2}$ and their respective output for the halting problem as constants used by the algorithm and, if $i_{w} = j_{1}$ or $i_{w} = j_{2}$ , we can output the correct result for these cases. More specifically, given an input w, compute the initial configuration $c_{i_{w}}$ associated to this input and compute its index $i_{w}$ . Then test if $i_{w} = j_{1}$ . If the test $i_{w} = j_{1}$ succeeds, then accept (reject) if M halts (does not halt, respectively) starting on configuration $c_{j_{1}}$ . Otherwise, test if $i_{w} = j_{2}$ . If the test $i_{w} = j_{2}$ succeeds, then accept (reject) if M halts (does not halt, respectively) starting on configuration $c_{j_{2}}$ . If both tests fail, test whether $a_{i_{w}} \in I$ by testing whether ${\tilde{x}}_{1} < a_{i_{w}} < {\tilde{x}}_{2}$ . We note that in this case it must be ${{\tilde{x}}_{1}, {\tilde{x}}_{2}} \cap B_{i_{w}} = \emptyset$ , and thus this test can be done in finite time by checking if $a_{j_{1}} < a_{i_{w}} < a_{j_{2}}$ . If the test succeeds, then accept w; otherwise reject it.

In other words, if the ODE (29) simulates M, then we can decide the halting problem, a contradiction. □

Footnotes

Acknowledgements

D. Graça wishes to thank Marco Mackaaij and Nenad Manojlović for helpful discussions. D. Graça was partially funded by FCT/MCTES through national funds and when applicable co-funded EU funds under the project UIDB/50008/2020. This project has received funding from the European Union’s Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No 731143.

References

Asarin and

Bouajjani , Perturbed Turing machines and hybrid systems, in: Proc. 16th Annual IEEE Symposium on Logic in Computer Science, 2001, pp. 269–278.

Bournez , Achilles and the tortoise climbing up the hyper-arithmetical hierarchy, Theoretical Computer Science 210(1) (1999), 21–71. doi:10.1016/S0304-3975(98)00096-6.

Bournez and

M.L.

Campagnolo , A survey on continuous time computations, in: New Computational Paradigms. Changing Conceptions of What Is Computable,

S.B.

Cooper ,

Löwe and

Sorbi , eds, Springer, New York, 2008, pp. 383–423.

Bournez ,

M.L.

Campagnolo ,

D.S.

Graça and

Hainry , Polynomial differential equations compute all real computable functions on computable compact intervals, Journal of Complexity 23(3) (2007), 317–335. doi:10.1016/j.jco.2006.12.005.

Bournez ,

D.S.

Graça and

Hainry , Computation with perturbed dynamical systems, Journal of Computer and System Sciences 79(5) (2013), 714–724. doi:10.1016/j.jcss.2013.01.025.

Bournez ,

D.S.

Graça and

Pouly , On the complexity of solving initial value problems, in: Proc. 37h International Symposium on Symbolic and Algebraic Computation (ISSAC 2012), 2012, arXiv:1202.4407.

Bournez ,

D.S.

Graça and

Pouly , Polynomial time corresponds to solutions of polynomial ordinary differential equations of polynomial length, Journal of the ACM 64(6) (2017), 38:1–38:76. doi:10.1145/3127496.

Bournez and

Pouly , A survey on analog models of computation, in: Handbook of Computability and Complexity in Analysis, Springer International Publishing, Cham, 2021, pp. 173–226. doi:10.1007/978-3-030-59234-9_6.

M.S.

Branicky , Universal computation and other capabilities of hybrid and continuous dynamical systems, Theoretical Computer Science 138(1) (1995), 67–100. doi:10.1016/0304-3975(94)00147-B.

10.

Brattka ,

Hertling and

Weihrauch , A tutorial on computable analysis, in: New Computational Paradigms: Changing Conceptions of What Is Computable,

S.B.

Cooper ,

Löwe and

Sorbi , eds, Springer, 2008, pp. 425–491. doi:10.1007/978-0-387-68546-5_18.

11.

Braverman , Computational complexity of Euclidean sets: Hyperbolic Julia sets are poly-time computable, in: Proc. 6th Workshop on Computability and Complexity in Analysis (CCA 2004),

Brattka ,

Staiger and

Weihrauch , eds, Electronic Notes in Theoretical Computer Science, Vol. 120, Elsevier, 2005, pp. 17–30.

12.

M.L.

Campagnolo , Computational Complexity of Real Valued Recursive Functions and Analog Circuits, PhD thesis, Instituto Superior Técnico/Universidade Técnica de Lisboa, 2002.

13.

M.L.

Campagnolo ,

Moore and

J.F.

Costa , Iteration, inequalities, and differentiability in analog computers, Journal of Complexity 16(4) (2000), 642–660. doi:10.1006/jcom.2000.0559.

14.

Cardona ,

Miranda and

Peralta-Salas , Turing Universality of the Incompressible Euler Equations and a Conjecture of Moore, International Mathematics Research Notices (2021), rnab233. doi:10.1093/imrn/rnab233.

15.

Cardona ,

Miranda ,

Peralta-Salas and

Presas , Constructing Turing complete Euler flows in dimension 3, Proceedings of the National Academy of Sciences 118(19) (2021). doi:10.1073/pnas.2026818118.

16.

Collins and

D.S.

Graça , Effective computability of solutions of ordinary differential equations the thousand monkeys approach, in: Proc. 5th International Conference on Computability and Complexity in Analysis (CCA 2008),

Brattka ,

Dillhage ,

Grubba and

Klutsch , eds, Electronic Notes in Theoretical Computer Science, Vol. 221, Elsevier, 2008, pp. 103–114.

17.

G.M.

Constantine and

T.H.

Savits , A multivariate Faa di Bruno formula with applications, Trans. Amer. Math. Soc. 348(2) (1996). doi:10.1090/S0002-9947-96-01501-2.

18.

D.S.

Graça ,

M.L.

Campagnolo and

Buescu , Computability with polynomial differential equations, Advances in Applied Mathematics 40(3) (2008), 330–349. doi:10.1016/j.aam.2007.02.003.

19.

D.S.

Graça ,

Zhong and

Buescu , Computability, noncomputability and undecidability of maximal intervals of IVPs, Transactions of the American Mathematical Society 361(6) (2009), 2913–2927. doi:10.1090/S0002-9947-09-04929-0.

20.

Kawamura , Lipschitz continuous ordinary differential equations are polynomial-space complete, Computational Complexity 19(2) (2010), 305–332. doi:10.1007/s00037-010-0286-0.

21.

Kawamura ,

Ota ,

Rösnick and

Ziegler , Computational complexity of smooth differential equations, Logical Methods in Computer Science 10(1:6) (2014), 1–15.

22.

K.-I.

Ko , Complexity Theory of Real Functions, Birkhäuser, 1991.

23.

Koiran , A family of universal recurrent networks, Theoretical Computer Science 168(2) (1996), 473–480. doi:10.1016/S0304-3975(96)00088-6.

24.

Koiran ,

Cosnard and

Garzon , Computability with low-dimensional dynamical systems, Theoretical Computer Science 132 (1994), 113–128. doi:10.1016/0304-3975(94)90229-1.

25.

Koiran and

Moore , Closed-form analytic maps in one and two dimensions can simulate universal Turing machines, Theoretical Computer Science 210(1) (1999), 217–223. doi:10.1016/S0304-3975(98)00117-0.

26.

Kurganskyy and

Potapov , Computation in one-dimensional piecewise maps and planar pseudo-billiard systems, in: Proc. 4th International Conference on Unconventional Computation (UC 2005), Lecture Notes in Computer Science, Vol. 3699, Springer, 2005, pp. 169–175. doi:10.1007/11560319_16.

27.

T.-W.

Ma , Higher chain formula proved by combinatorics, The Electronic Journal of Combinatorics 16(1) (2009).

28.

Maass and

Orponen , On the effect of analog noise in discrete-time analog computations, Neural Computation 10(5) (1998), 1071–1095. doi:10.1162/089976698300017359.

29.

Matiyasevich , Hilbert’s 10th Problem, The MIT Press, 1993.

30.

Moore , Generalized shifts: Unpredictability and undecidability in dynamical systems, Nonlinearity 4(2) (1991), 199–230. doi:10.1088/0951-7715/4/2/002.

31.

Moschovakis , Kleene’s amazing second recursion theorem, The Bulletin of Symbolic Logic 16 (2010), 189–239. doi:10.2178/bsl/1286889124.

32.

Odifreddi , Classical Recursion Theory, Vol. 1, Elsevier, 1989.

33.

Sipser , Introduction to the Theory of Computation, 3rd edn, Cengage Learning, 2012.

34.

A.M.

Turing , On computable numbers, with an application to the entscheidungsproblem, Proceedings of the London Mathematical Society s2–42(1) (1937), 230–265. doi:10.1112/plms/s2-42.1.230.