<figure><img src="image_1.png" /></figure>

HAL
open science

# Hilbert's Tenth Problem in Coq (Extended Version)

Dominique Larchey-Wendling, Yannick Forster

► To cite this version:

Dominique Larchey-Wendling, Yannick Forster. Hilbert's Tenth Problem in Coq (Extended Version). Logical Methods in Computer Science, 2022, 18 (1), 10.46298/lmcs-18(1:35)2022 . hal-03359505v2

HAL Id: hal-03359505

https://inria.hal.science/hal-03359505v2

Submitted on 10 Mar 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
HILBERT'S TENTH PROBLEM IN COQ
(EXTENDED VERSION *)

DOMINIQUE LARCHEY-WENDLING AND YANNICK FORSTER

a Université de Lorraine, CNRS, LORIA, Vandœuvre-lès-Nancy, France
*e-mail address:* dominique.larchey-wendling@loria.fr

b Saarland University, Saarland Informatics Campus, Saarbrücken, Germany
*e-mail address:* forster@cs.uni-saarland.de

ABSTRACT. We formalise the undecidability of solvability of Diophantine equations, i.e. polynomial equations over natural numbers, in Coq's constructive type theory. To do so, we give the first full mechanisation of the Davis-Putnam-Robinson-Matiyasevich theorem, stating that every recursively enumerable problem – in our case by a Minsky machine – is Diophantine. We obtain an elegant and comprehensible proof by using a synthetic approach to computability and by introducing Conway's FRACTRAN language as intermediate layer. Additionally, we prove the reverse direction and show that every Diophantine relation is recognisable by $\mu$-recursive functions and give a certified compiler from $\mu$-recursive functions to Minsky machines.

# 1. INTRODUCTION

Hilbert's tenth problem (H10) was posed by David Hilbert in 1900 as part of his famous 23 problems [Hil02] and asked for the "determination of the solvability of a Diophantine equation." A Diophantine equation<fnref n="1" /> is a polynomial equation over natural numbers (or, equivalently, integers) with constant exponents, e.g. $x^2 + 3z = yz + 2$. When Hilbert asked for "determination," he meant, in modern terms, a decision procedure, but computability theory was yet several decades short of being developed.

The first undecidable problems found by Church, Post and Turing were either native to mathematical logic or dependent on a fixed model of computation. H10, to the contrary, can be stated to every mathematician and its formulation is independent from a model of computation. Emil Post stated in 1944 that H10 "begs for an unsolvability proof" [Pos44]. From a computational perspective, it is clear that H10 is recursively enumerable (or *recognisable*), meaning there is an algorithm that halts on a Diophantine equation iff it is solvable.

Post's student Martin Davis conjectured that even the converse is true, i.e. that every recognisable set is also Diophantine. More precisely, he conjectured that if $A \subseteq \mathbb{N}^k$ is

<figure><img src="image_3.png" /></figure>

<footnote>Key words and phrases: Hilbert's tenth problem, Diophantine equations, undecidability, computability theory, reduction, Minsky machines, Fractran, Coq, type theory.</footnote>

<footnote>* extended version of [LWF19].</footnote>

<footnote n="1">¹Named after the Greek mathematician Diophantus of Alexandria, who started the study of polynomial equations in the third century.</footnote>
recognisable then $(a_1, \dots, a_k) \in A \leftrightarrow \exists x_1 \dots x_n, P(a_1, \dots, a_k, x_1, \dots, x_n) = 0$ holds for some polynomial $P$ in $k+n$ variables. He soon improved on a result by Gödel [Göd31] and gave a proof of his conjecture, however requiring up to one bounded universal quantification [Dav53]: $(a_1, \dots, a_k) \in A \leftrightarrow \exists z, \forall y < z, \exists x_1 \dots x_n, P(a_1, \dots, a_k, x_1, \dots, x_n, y, z) = 0$. Davis and Putnam [DP59] further improved on this, and showed that, provided a certain number-theoretic assumption holds, every recognisable set is *exponentially* Diophantine, meaning variables are also allowed to appear in exponents. Julia Robinson then in 1961 modified the original proof to circumvent the need for the assumption, resulting in the DPR theorem [DPR61], namely that every recognisable set is exponentially Diophantine. Due to another result from Robinson [Rob52], the gap now only consisted of proving that there is a Diophantine equation exhibiting exponential growth. In 1970, Yuri Matiyasevich showed that the Fibonacci sequence grows exponentially while being Diophantine, closing the gap and finishing the proof of the theorem nowadays called *DPRM theorem*, ultimately establishing that exponentiation is Diophantine itself [Mat70] (known as "Matiyasevich's theorem").

Even the most modern and simpler proofs of the DPRM theorem still require many preliminaries and complicated number-theoretic ideas, for an overview see [Mat16]. We formalise one such proof as part of our ongoing work on a library of undecidable problems [FLW18, FLWD$^+$20] in the proof assistant Coq [Tea21]. Since H10 is widely used as a seed [DR19, Gol81] for showing the undecidability of problems using *many-one reductions*, this will open further ways of extending the library. Given that our library already contains a formalisation of Minsky machines [FL19], we follow the approach of Jones and Matijasevič [JM84], who use register machines, being very well-suited since they already work on numbers. They encode full computations of register machines as Diophantine equations in one single, monolithic step. To make the proof more tractable for both mechanisation and explanation, we factor out an intermediate language, John Conway's FRACTRAN [Con87], which can simulate Minsky machines.

We first introduce three characterisations of Diophantine equations over natural numbers, namely *Diophantine logic* DIO_FORM (allowing to connect basic Diophantine equations with conjunction, disjunction and existential quantification), *elementary Diophantine constraints* DIO_ELEM (a finite set of constraints on variables, oftentimes used for reductions [DR19, Gol81]) and *single Diophantine equations* DIO_SINGLE, including parameters, as described above. H10 then asks about the solvability of single Diophantine equations with no parameters.

Technically, the reduction chain to establish the unsolvability of H10 starts at the halting problem for single-tape Turing machines *Halt*, reduced to the Post correspondence problem PCP in [FHS18]. In previous work [FL19] we have reduced PCP to a specialised halting problem for Minsky machines, which we use here in a slightly generalised form as MM. We then reduce Minsky machine halting to FRACTRAN termination. FRACTRAN is very natural to describe using polynomials, and the encoding does not rely on any complicated construction. The technical difficulty then only lies in the Diophantine encoding of the reflexive-transitive closure of a relation which follows from the direct elimination of bounded universal quantification, given that the proof in [Mat97] involves no detour via models of computation. In total, we obtain the following chain of reductions<fnref n="2" /> establishing the

<footnote n="2">²A *many-one reduction* $P \le Q$, later defined formally in this section, is a computable function from problem $P$ to problem $Q$ that maps instances of $P$ into instances of $Q$, pereserving both validity and invalidity.</footnote>
**PCP:** Post correspondence problem, see e.g. [FL19]. (matching)

**MM:** Given $n : \mathbb{N}$, a Minsky machine $P : \mathbb{L}I_n$ with $n$ registers, and $\vec{v} : \mathbb{N}^n$, does $(1, P)$ terminate from input state $(1, \vec{v})$? (termination)

**FRACTRAN:** Given a regular FRACTRAN program $Q : \mathbb{L}(\mathbb{N} \times \mathbb{N})$ and an input state $s$, does $Q$ terminate from input state $x$? (termination)

**DIO_FORM:** Given a Diophantine logic formula $A : \mathbb{D}_{\text{form}}$ and a valuation $\nu : \mathbb{V} \to \mathbb{N}$, does $[A]_{\nu}$ hold? (satisfaction)

**DIO_ELEM:** Given a list $l : \mathbb{L}\mathbb{D}_{\text{cstr}}$ of elementary Diophantine constraints and a valuation $\nu : \mathbb{V} \to \mathbb{N}$, does there exist $\varphi : U \to \mathbb{N}$ such that $\forall c \in l, [[c]]_{\nu}^{\varphi}$? (simultaneous satisfiability)

**DIO_SINGLE:** Given a single Diophantine equation $p \doteq q : \mathbb{D}_{\text{single}}(\mathbb{N}, \mathbb{N})$ and a valuation $\nu : \mathbb{N} \to \mathbb{N}$, does there exist $\varphi : \mathbb{N} \to \mathbb{N}$ s.t. $[[p]]_{\nu}^{\varphi} = [[q]]_{\nu}^{\varphi}$? (solvability)

**H10:** Given a single Diophantine equation $p \doteq q : \mathbb{D}_{\text{single}}(\mathbb{F}_n, \mathbb{F}_0)$ (over $\mathbb{N}$ with possibly $n$ variables but no parameters), does it have a solution in $\mathbb{N}$? (solvability)

**$\mu$-rec:** Given $n : \mathbb{N}$, an $n$-ary $\mu$-recursive function $f : \mathcal{A}_n$, and $\vec{v} : \mathbb{N}^n$, does $\vec{v}$ belong to the domain of $f$? (termination)

**H10<sub>Z</sub>:** Given a single Diophantine equation $p \doteq q : \mathbb{D}_{\text{single}}(\mathbb{F}_n, \mathbb{F}_0)$ (over $\mathbb{Z}$ with possibly $n$ variables but no parameters), does it have a solution in $\mathbb{Z}$? (solvability)

<caption>Figure 1: Summary description of some decision problems.</caption>

undecidability of H10 and it's many-one interreducibility with several decision problems:

Halt $\leq$ PCP $\leq$ MM $\leq$ FRACTRAN $\leq$ DIO_FORM $\leq$ DIO_ELEM $\leq$ DIO_SINGLE $\leq$ H10 $\leq$ $\mu$-rec $\leq$ MM

where Fig. 1 lists high-level descriptions of these problems. Furthermore, we prove that H10 $\leq$ H10<sub>Z</sub> via Lagrange's theorem. In the present paper, we focus on explaining this factorisation of the proof and give some details for the different stages. While we contribute Coq mechanisations of Matiyasevich's theorem and the elimination of bounded universal quantification, we treat them mainly as black-boxes and only elaborate on their challenging formalisation rather than the proofs themselves (see Section 2.3).

To the best of our knowledge, we are the first to give a *full verification* of the DPRM theorem and the undecidability of Hilbert's tenth problem in a proof assistant. We base the notion of recognisability in the DPRM theorem on Minsky machines.

When giving undecidability proofs via many-one reductions, it is critical to show that all reduction functions are actually computable. We could in theory verify the computability of all functions involved using an explicit model of computation. In pen-and-paper proofs, this approach is however almost never used, because implementing high-level mathematical transformations as provably correct low-level programs is a daunting task. Instead, we rely on a synthetic approach [FHS18, FKS19, FL19] based on the computability of all functions definable in Coq's constructive type theory, which is closer to the practice of pen-and-paper proofs. In this approach, a problem $P$ is considered undecidable if there is a reduction from an obviously undecidable problem, e.g. $Halt \leq P$.

The axiom-free Coq formalisation of all the results in this paper is available online [FLWD+21] and the main lemmas and theorems in the pdf version of the paper are hyper-linked with the html version of the source code at

https://github.com/uds-psl/coq-library-undecidability/tree/H10-LMCS-v1.1
Starting from our already existing library which included most of the Minsky machine code [FL19], the additional code for proving the undecidability of H10 and the DPRM theorem consists of about 8k loc including 3k loc for Matiyasevich's results alone, together with a 4k loc addition to our shared libraries; see Appendix A for more details. The paper itself can be read without in-depth knowledge of Coq or type theory.

**1.1. Contribution.** This paper is an extended journal version of a conference paper [LWF19] which, besides a full mechanisation of the DPRM theorem, contributed a novel refactoring of the proof via FRACTRAN improving the explainability of the DPRM theorem. Compared to the conference version, we contribute mechanised proofs showing that

• H10 reduces to solvability of Diophantine equations over integers, reduction obtained via a low-level implementation of Lagrange's theorem (Section 9);

• Diophantine relations are recognisable by $\mu$-recursive algorithms (Section 10);

• $\mu$-recursively recognisable relations are MM-recognisable (Section 10), thereby proving that all considered problems are in the same many-one reduction class;

• $\mu$-recursive algorithms can be simulated in the weak call-by-value $\lambda$-calculus (Section 11), thereby proving that all considered problems are in the same many-one reduction class as most problems in the Coq library of undecidable problems [FLWD$^+$20].

Apart from the new results, we have simplified the account of Diophantine logic considerably and expanded various explanations of proofs.

**1.2. Preliminaries.** For the text, while we cannot completely avoid it, we will try to minimize reliance on type theoretic language and notations. We write $\mathbb{P}$ for the (impredicative) type of propositions and **Type** for the (predicative hierarchy of) types of Coq. When $X$ and $Y$ are types, we write $X \to Y$ for functions from $X$ to $Y$.<fnref n="3" /> We write $x \times y$ or $x \cdot y$ for multiplication of natural numbers $x, y : \mathbb{N}$ and we will leave out the symbol where convenient. We write $\mathbb{L}X$ for the type of *lists* over $X$ and $l \dot{=} l'$ for the *concatenation* of two lists. We write $X^n$ for *vectors* $\vec{v}$ over type $X$ with length $n$, and $\mathbb{F}_n$ for the *finite type* with exactly $n$ elements. For $p : \mathbb{F}_n$, we write $\vec{v}_p$ for the $p$-th component of $\vec{v} : X^n$. Notations for lists are overloaded for vectors.

We write $\mathcal{O}X$ for the type of *options* over $X$ with elements $[x]$ where $x : X$ and $\emptyset$. We write $X + Y$ for the type-theoretic sum of types $X$ and $Y$, consisting of $\text{inl}~x$ for $x : X$ and $\text{inr}~y$ for $y : Y$. For a list $l : \mathbb{L}X$, $l_n : \mathcal{O}X$ denotes the $n$-th value in $l$ if it exists.

If $P : X \to \mathbb{P}$ is a predicate (on $X$) and $Q : Y \to \mathbb{P}$ is a predicate, we write $P \preceq Q$ if there is a function $f : X \to Y$ s.t. $\forall x : X, P x \leftrightarrow Q(f x)$, i.e. a *many-one reduction* from $P$ to $Q$. In the synthetic approach [FHS18, FKS19, FL19], the computability of the reduction $f$ is automatically ensured because $f$ is typeable in Coq without relying on any axiom.

<footnote n="3">³In the case where $X$ and $Y$ are moreover propositions, the type/proposition $X \to Y$ is understood as $X$ implies $Y$, inhabited with functions mapping proofs of $X$ into proofs of $Y$; the type theory of Coq fully implements the Curry-Howard isomorphism.</footnote>
## 2. DIOPHANTINE RELATIONS

Diophantine relations are composed of polynomials over natural numbers. There are several equivalent approaches to characterise these relations and oftentimes, the precise definition is omitted from papers. Basically, one can form equations between polynomial expressions and then combine these with conjunctions, disjunctions, and existential quantification.<fnref n="4" /> For instance, these operations are assumed as Diophantine producing operators in e.g. [JM84, Mat70, Mat97, Mat00]. Sometimes Diophantine relations are restricted to a single polynomial equation. Sometimes the exponentiation function $x, y \mapsto x^y$ is assumed as Diophantine [JM84]. To complicate the picture, Diophantine relations might equivalently range over $\mathbb{Z}$ (instead of $\mathbb{N}$) but expressions like $x^y$ implicitly assume that $y$ never gets a negative value.

Although seemingly diverging, these approaches are not contradictory because in the end, they characterise the same class of relations on natural numbers. However, mechanisation does not allow for such implicit assumptions. To give some mechanisable structure to some of these approaches, we propose three increasingly restricted characterisations of Diophantine relations: *Diophantine logic*, *elementary Diophantine constraints* and *single Diophantine equations*, between which we provide computable transformations in Sections 3 and 4. In Section 9, we also elaborate on the case of polynomials over $\mathbb{Z}$, i.e. we give an undecidability proof for Hilbert tenth problem over integers. But before turning to formal definitions, we motivate our approach for the automated analysis and recognition of Diophantine shapes.

### 2.1. Diophantine Shapes.
We introduce the meta-level notion of *Diophantine shape*. On purpose, this notion does not have a precise formal definition because it is a *dynamically evolving property* of meta-level expressions that upgrades itself as more and more closure results are proved about those shapes.

Of course, we cannot rely on some blurry notion to formally prove theorems about Diophantine relations. So at some point, we have to choose one (or several) formal representation(s) of Diophantine relations. Irrelevant to the actual syntax we finally pick up for the formal description of Diophantine relations, we call them *object-level representations*.

The problem we face is the following: how can we minimize the work we have to do to actually build an object-level witness corresponding to a given meta-level expression representing a Diophantine relation. Directly manipulating object-level syntax is far from the ideal way to proceed for a lazy Coq programmer,<fnref n="5" /> and indeed, this empirical lesson can be learned the painful way.

Let us illustrate this on the following complex example of meta-level expression

$$x = z^{x+k} \wedge \forall y, y < k \rightarrow \exists u v, u \equiv 2v [p] \wedge p = 2^y \wedge C_u^z \equiv x [v].$$

After a quick analysis of its structure, it appears that before being able to actually establish that this is a Diophantine relation, we would probably have to show that polynomials, identities, arithmetic congruences, the exponential function, binomial coefficients, conjunction, existential quantification and bounded universal quantification are all Diophantine admissible, hence to give *procedures* to derive object-level representations for all these meta-level constructions. This already amounts to significant work. But once this is done, we want to avoid both the hand-building of the object-level witness for the above expression, and the

<footnote n="4">⁴Universal quantification or negation are not accepted as is.</footnote>

<footnote n="5">⁵being able to defer repetitive work to computers is critical to the successful completion of mechanizations.</footnote>
proofs that it correctly reflects its semantics. We essentially require our framework to be able to automatically combine those procedures and build a provably correct witness for us.

To summarise, we aspire at the formal definition of an object-level representation and at the same time, at avoiding its direct manipulation. This is where the dynamic notion of Diophantine shape plays a central role. At first, there are very few basic Diophantine shapes, typically constants, variables, addition, multiplication, equality. But at some point, we would e.g. have a result stating that if the expressions *f* and *g* have a Diophantine shape, then so does the expression *f*⁹, that is the Diophantine admissibility of the exponential function, nowadays called Matiyasevich's theorem. Critically, shapes can be added dynamically as they are proved admissible as opposed to be carved in the stone of a given object-level syntax. We now describe how to do this in a successful way using some of the automation provided by Coq.

**2.2. Diophantine Logic.** We define the type $D_{\text{form}}$ of Diophantine formulæ for the abstract syntax of Diophantine logic. An atomic Diophantine logic formula is just expressing basic atomic identities between variables like $x_i \doteq x_j \dot{+} x_k$ or $x_i \doteq x_j \dot{\times} x_k$ and we combine those with binary disjunction, binary conjunction, and existential quantification.

$$A, B : D_{\text{form}} ::= x_i \doteq n \mid x_i \doteq x_j \mid x_i \doteq x_j \dot{+} x_k \mid x_i \doteq x_j \dot{\times} x_k \mid A \dot{\wedge} B \mid A \dot{\vee} B \mid \exists A$$

The letters *A*, *B* range over formulæ and *n*: $\mathbb{N}$ represents constant ranging over natural numbers. We use standard *De Bruijn syntax* with variables $x_0, x_1, \dots$ of type $V := \mathbb{N}$ for better readability. If we have $x_i: V$, we write $x_{1+i}$ for the next variable in $V$. As an example, the meta-level formula $\exists y, (y = 0 \wedge \exists z, y = z + k)$ would be represented as e.g. $\exists (x_0 \doteq 0 \wedge \exists (x_1 \doteq x_0 \dot{+} x_2))$, i.e. the variable $x_i$ refers to the *i*-th binder in the context. Notice that there is no type or syntactic distinction between variables and parameters in Diophantine logic. However some variables are bound in their context and others are free.

We provide a semantics for Diophantine logic. Given a valuation for variables $\nu: V \to \mathbb{N}$, we define the interpretation $[[x_i \doteq \dots]]_\nu: \mathbb{P}$ of atomic formulæ by

$$[[x_i \doteq n]]_\nu := \nu x_i = n \qquad [[x_i \doteq x_j \dot{+} x_k]]_\nu := \nu x_i = \nu x_j + \nu x_k$$

$$[[x_i \doteq x_j]]_\nu := \nu x_i = \nu x_j \qquad [[x_i \doteq x_j \dot{\times} x_k]]_\nu := \nu x_i = \nu x_j \times \nu x_k$$

and $[[A]]_\nu: \mathbb{P}$ for a compound formula $A: D_{\text{form}}$ by the following recursive equations:

$$[[A \dot{\wedge} B]]_\nu := [[A]]_\nu \wedge [[B]]_\nu \qquad [[A \dot{\vee} B]]_\nu := [[A]]_\nu \vee [[B]]_\nu \qquad [\dot{\exists} A]]_\nu := \exists n: \mathbb{N}, [[A]]_{n \cdot \nu}$$

where $n \cdot \nu: V \to \mathbb{N}$ is the standard De Bruijn extension<fnref n="6" /> of a valuation $\nu$ by $n$, defined by $n \cdot \nu(x_0) := n$ and $n \cdot \nu(x_{1+i}) := \nu x_i$.

We give a first object-level representation of Diophantine relations as members of type $(V \to \mathbb{N}) \to \mathbb{P}$ mapping valuations of variables to propositions. Moreover, they must be identical to $\lambda\nu.[[A]]_\nu$ for some Diophantine formula $A$, up to propositional extensionality. We give an informative content to this sub-type of $(V \to \mathbb{N}) \to \mathbb{P}$ to be able to do some computations with the witness $A: D_{\text{form}}$ of Diophantineineness, typically when moving to another formal representation like elementary Diophantine constraints in Section 3.

**Definition 2.1.** The class of Diophantine relations $D_{\text{rel}}: ((V \to \mathbb{N}) \to \mathbb{P}) \to \text{Type is the informative sub-type defined for } R: (V \to \mathbb{N}) \to \mathbb{P}$ by

$$D_{\text{rel}} R := \sum A : D_{\text{form}}, (\forall \nu, [[A]]_{\nu} \leftrightarrow R \nu).$$

<footnote n="6">⁶The notation $n \cdot \nu$ emphasizes that the value $n$ is pushed ahead of the infinite sequence $\nu x_0; \nu x_1; \nu x_2; \dots$.</footnote>
Note that $\Sigma$ denotes type-theoretic dependent pairs. Hence an inhabitant $w$ of $D_{\text{rel}} R$ is a (dependent) pair $(A, H_A)$ where $A = \pi_1(w)$ is a Diophantine formula and $H_A = \pi_2(w)$ a proof that $[A]_{\cdot}$ and $R$ are extensionally equivalent.<fnref n="7" /> The informal notion of Diophantine shape will correspond to the dynamically growing collection of established closure properties of the class $D_{\text{rel}}$ of Diophantine relations. Definition 2.1 of the sub-type $D_{\text{rel}}$ already entails that $D_{\text{rel}}$ is closed under conjunction, disjunction, existential quantification and renaming.

**Proposition 2.2.** Let $R, S : (V \to \mathbb{N}) \to \mathbb{P}$ be relations, $T : \mathbb{N} \to (V \to \mathbb{N}) \to \mathbb{P}$ be a relation with a singled out variable, and $\rho : V \to V$ be a renaming function. We have the maps:

1. $D_{\text{rel}} R \to D_{\text{rel}} S \to D_{\text{rel}}(\lambda\nu.R\nu \wedge S\nu);$

4. $(\forall\nu, S\nu \leftrightarrow R\nu) \to D_{\text{rel}} R \to D_{\text{rel}} S;$

2. $D_{\text{rel}} R \to D_{\text{rel}} S \to D_{\text{rel}}(\lambda\nu.R\nu \lor S\nu);$

5. $D_{\text{rel}} R \to D_{\text{rel}}(\lambda\nu.R(\nu \circ \rho)).$

3. $D_{\text{rel}}(\lambda\nu.T(\nu x_0)(\lambda x_i.\nu x_{1+i})) \to D_{\text{rel}}(\lambda\nu.\exists u, T u \nu);$

Understood as Diophantine shapes, maps number 1–3 recognise the logical connectives of conjunction, disjunction and existential quantification as newly allowed shapes. Map number 5 allows renaming (free) variables hence Diophantine shapes are closed under renaming.<fnref n="8" /> Map number 4 provides a way to replace the goal $D_{\text{rel}} S$ with $D_{\text{rel}} R$ once a proof that they are logically equivalent is established. Hence, if $S$ cannot be analysed because it does not currently have a Diophantine shape, it can still be replaced by an equivalent relation $R$, hopefully better behaved; see e.g. the proof of Proposition 2.4.

Working with Diophantine relations already gives a satisfying implementation of Diophantine shapes but it is sometimes more convenient to manipulate Diophantine functions instead of relations so we define the following specialization.

**Definition 2.3.** The class of Diophantine functions $D_{\text{fun}} : ((V \to \mathbb{N}) \to \mathbb{N}) \to \text{Type is the informative sub-type defined for } f : (V \to \mathbb{N}) \to \mathbb{N}$ by $D_{\text{fun}} f := D_{\text{rel}}(\lambda\nu.\nu x_0 = f(\lambda x_i.\nu x_{1+i}))$.

We extend Diophantine shapes with polynomials expressions and equations between them. To illustrate the mechanics behind Diophantine shape recognition, for once we give a detailed account of the steps that are usually performed automatically in this framework.

**Proposition 2.4.** Let $x_i : V, n : \mathbb{N}, \rho : V \to V$ and $f,g : (V \to \mathbb{N}) \to \mathbb{N}$. We have the maps:

1. $D_{\text{fun}}(\lambda\nu.\nu x_i);$

4. $D_{\text{fun}} f \to D_{\text{fun}} g \to D_{\text{fun}}(\lambda\nu.f\nu + g\nu);$

2. $D_{\text{fun}}(\lambda\nu.n);$

5. $D_{\text{fun}} f \to D_{\text{fun}} g \to D_{\text{fun}}(\lambda\nu.f\nu \times g\nu);$

3. $D_{\text{fun}} f \to D_{\text{fun}}(\lambda\nu.f(\nu \circ \rho));$

6. $D_{\text{fun}} f \to D_{\text{fun}} g \to D_{\text{rel}}(\lambda\nu.f\nu = g\nu).$

*Proof.* Items 1, 2 and map 3 are for projections, constants and renaming functions respectively. The non-trivial cases are for $+$, $\times$ and $=$. We cover the cases of $+$ and then $=$ in details to illustrate how the representations of Diophantine relations/functions behave in proof scripts. In particular, we prove the results reasoning backwards (as is usually done in Coq), applying established theorems to convert a given proof goal into (hopefully) simpler proof goals.

For the goal $D_{\text{fun}}(\lambda\nu.f\nu + g\nu)$, unfolding the assumptions $D_{\text{fun}} f$ and $D_{\text{fun}} g$ we have

$$D_{\text{rel}}(\lambda\nu.\nu x_0 = f(\lambda x_i.\nu x_{1+i})) \quad \text{and} \quad D_{\text{rel}}(\lambda\nu.\nu x_0 = g(\lambda x_i.\nu x_{1+i})) \quad (2.1)$$

and we want to establish $D_{\text{fun}}(\lambda\nu.f\nu + g\nu)$, i.e.

$$D_{\text{rel}}(\lambda\nu.\nu x_0 = f(\lambda x_i.\nu x_{1+i}) + g(\lambda x_i.\nu x_{1+i})). \quad (2.2)$$

<footnote n="7">⁷For the efficiency of computations, we usually hide the purely logical part $H_A$ into an *opaque proof term*.</footnote>

<footnote n="8">⁸$\nu \circ \rho : V \to \mathbb{N}$ denotes the composition $\lambda x_i.\nu x_{\rho_i}$ of the valuation $\nu : V \to \mathbb{N}$ with the renaming $\rho : V \to V$.</footnote>
By map 4 of Proposition 2.2, we replace Eq. (2.2) with the (obviously) equivalent goal

$$ \mathbb{D}_{\text{rel}}(\lambda\nu.\exists a\exists b, \nu x_0 = a + b \land a = f(\lambda x_i.\nu x_{1+i}) \land b = g(\lambda x_i.\nu x_{1+i})) \quad (2.3) $$

and we then apply map 3 of Proposition 2.2 twice to get the goal

$$ \mathbb{D}_{\text{rel}}(\lambda\nu.\nu x_2 = \nu x_1 + \nu x_0 \land \nu x_1 = f(\lambda x_i.\nu x_{3+i}) \land \nu x_0 = g(\lambda x_i.\nu x_{3+i})). $$

We now apply twice map 1 of Proposition 2.2 and we get the three following sub-goals:

$$ \mathbb{D}_{\text{rel}}(\lambda\nu.\nu x_2 = \nu x_1 + \nu x_0) \quad \mathbb{D}_{\text{rel}}(\lambda\nu.\nu x_1 = f(\lambda x_i.\nu x_{3+i})) \quad \mathbb{D}_{\text{rel}}(\lambda\nu.\nu x_0 = g(\lambda x_i.\nu x_{3+i})) $$

(1) For the first sub-goal, we use the formula $x_2 \doteq x_1 \dotplus x_0 : \mathbb{D}_{\text{form}}$ as object-level witness;

(2) for the second sub-goal, we consider the renaming function $\rho_1 : V \to V$ defined by $\rho_1(x_0) := x_1$ and $\rho_1(x_{1+i}) := x_{3+i}$ and derive $\mathbb{D}_{\text{rel}}(\lambda\nu.\nu x_1 = f(\lambda x_i.\nu x_{3+i}))$ by applying map 5 of Proposition 2.2 to the hypothesis $\mathbb{D}_{\text{rel}}(\lambda\nu.\nu x_0 = f(\lambda x_i.\nu x_{1+i}))$ in Eqs. (2.1);

(3) the third and last sub-goal $\mathbb{D}_{\text{rel}}(\lambda\nu.\nu x_0 = g(\lambda x_i.\nu x_{3+i}))$ is solved similarly with the renaming function $\rho_0 : V \to V$ defined by $\rho_0(x_0) := x_0$ and $\rho_0(x_{1+i}) := x_{3+i}$.

We now deal with map 6, hence with the goal $\mathbb{D}_{\text{rel}}(\lambda\nu.f\nu = g\nu)$ under the same previous assumptions $\mathbb{D}_{\text{fun}} f$ and $\mathbb{D}_{\text{fun}} g$, i.e. Eqs. (2.1). We proceed in a somewhat less detailed explanation. We replace the goal by $\mathbb{D}_{\text{rel}}(\lambda\nu.\exists a\exists b, a = b \land a = f\nu \land b = g\nu)$ which is equivalent and then, after applying the maps of Proposition 2.2, we get three sub-goals $\mathbb{D}_{\text{rel}}(\lambda\nu.\nu x_1 = \nu x_0)$, $\mathbb{D}_{\text{rel}}(\lambda\nu.\nu x_1 = f(\lambda x_i.\nu x_{2+i}))$ and $\mathbb{D}_{\text{rel}}(\lambda\nu.\nu x_0 = g(\lambda x_i.\nu x_{2+i}))$. In turn, the first sub-goal corresponds to the witness $x_1 \doteq x_0 : \mathbb{D}_{\text{form}}$, while the second and third sub-goals follow from Eqs. (2.1) respectively using straightforward renaming functions. $\square$

On paper these proofs look somehow complicated by the need to infer the renaming functions but from a *mechanisation point of view*, Coq's unification algorithm automatically solves such goals. Provided we populate the *hint database* with enough admissible shapes, we can automate the analysis of the meta-level syntax to establish Diophatineness and reflect a meta-level expression of Diophantine shape into the corresponding object-level witness of type $\mathbb{D}_{\text{form}}$ together with the proof that it is an appropriate witness, hence packed into the types $\mathbb{D}_{\text{rel}} R$ for relational expressions or $\mathbb{D}_{\text{fun}} f$ for functional expressions.

With Propositions 2.2 and 2.4, we populate the hint database for relations with the shapes conjunction, disjunction, existential quantification, renaming and identity between two functional expressions, and for functions, we add the shapes of projections, constants, addition, multiplication and renaming. In our implementation, the analysis of Diophantine shapes is performed by the automatic *dio auto* tactic. With such an automated approach, the remaining (and sometimes difficult) work occurs when we apply map 4 of Proposition 2.2, that is, we have to find an equivalent expression of Diophantine shape, like in Eq. (2.3) and to prove it is indeed equivalent to Eq. (2.2), which, unlike that specific example, might be non-trivial; see e.g. the discussion in Section 5.

**Proposition 2.5.** Let $f, g : (V \to \mathbb{N}) \to \mathbb{N}$. We have the maps:

1. $\mathbb{D}_{\text{rel}}(\lambda\nu True)$;

3. $\mathbb{D}_{\text{fun}} f \to \mathbb{D}_{\text{fun}} g \to \mathbb{D}_{\text{rel}}(\lambda\nu.f\nu \le g\nu)$;

2. $\mathbb{D}_{\text{rel}}(\lambda\nu False)$;

4. $\mathbb{D}_{\text{fun}} f \to \mathbb{D}_{\text{fun}} g \to \mathbb{D}_{\text{rel}}(\lambda\nu.f\nu < g\nu)$;

5. $\mathbb{D}_{\text{fun}} f \to \mathbb{D}_{\text{fun}} g \to \mathbb{D}_{\text{rel}}(\lambda\nu.f\nu \ne g\nu)$.
*Proof.* For e.g. $\mathbb{D}_{\text{rel}}(\lambda\nu. f \nu < g \nu)$, we shift to the equivalent $\mathbb{D}_{\text{rel}}(\lambda\nu.\exists a, 1 + a + f \nu = g \nu)$ using map 4 of Proposition 2.2 and finish the proof calling *dio auto*<fnref n="9" /> □

Again, we populate the hint database with the new shapes of Proposition 2.5. We follow up with the slightly more complex example of the “does not divide” relation defined by $u \nmid v := \neg(\exists k, v = k \times u)$. At this point, this expression cannot be recognized as a Diophantine shape because it contains a negation.

**Proposition 2.6.** $\forall fg : (V \to \mathbb{N}) \to \mathbb{N}$, $\mathbb{D}_{\text{fun}} f \to \mathbb{D}_{\text{fun}} g \to \mathbb{D}_{\text{rel}}(\lambda\nu.f\nu \nmid g\nu)$.

*Proof.* However, using Euclidean division, we (easily) prove the equivalence

$$u \nmid v \leftrightarrow (u = 0 \land v \neq 0 \lor \exists a b, v = a \times u + b \land 0 < b < u)$$

and this new expression can now be recognised as a Diophantine shape. Using this equivalence in combination with map 4 of Proposition 2.2, we replace the goal $\mathbb{D}_{\text{rel}}(\lambda\nu.f\nu \nmid g\nu)$ with

$$\mathbb{D}_{\text{rel}}(\lambda\nu.f\nu = 0 \land g\nu \neq 0 \lor \exists a b, g\nu = a \times f\nu + b \land 0 < b < f\nu)$$

and then let the magic of *dio auto* unfold. □

Again, once established, we can add the map $\mathbb{D}_{\text{fun}} f \to \mathbb{D}_{\text{fun}} g \to \mathbb{D}_{\text{rel}}(\lambda\nu.f\nu \nmid g\nu)$ in the Diophantine hint database so that later encountered proof goals $\mathbb{D}_{\text{rel}}(\lambda\nu.f\nu \nmid g\nu)$ can be immediately solved by *dio auto*.

In this above described approach, the recovery of the object-level witness A of Definition 2.1 from meta-level syntax is automatic and hidden by the use of the *dio auto* tactic associated with the ever growing hint database. This allows us to proceed as in e.g. Matiyasevich papers where he usually transforms a relation into an equivalent Diophantine shape, accumulating more and more Diophantine shapes on the way. Instead of having to manipulate object-level witnesses by hand, obfuscating sometimes simple to understand proofs, we use Diophantine shapes as the cornerstone of the faithful implementation of existing pen and paper scripts.

**2.3. Exponentiation and Bounded Universal Quantification.** For now, we introduce the elimination of the exponential relation and then of bounded universal quantification as black boxes expressed in the framework of Diophantine shapes, i.e. new closure properties of the classes $\mathbb{D}_{\text{fun}}/\mathbb{D}_{\text{rel}}$.

While we do contribute implementations for both of these hard results, on purpose, we choose to avoid the detailed presentation of these mechanised proofs for several reasons:

* first of all, there are already fully detailed pen and paper accounts of these results and we implemented two of these somewhat faithfully;

* then, in our modular approach, the proof of these admissibility results can be ignored without hindering the understanding of the overall structure of the main results, e.g. H10;

* finally, already the above cited pen and paper proofs assume some not so standard results in arithmetic like e.g. Lucas's theorem, and we favoured giving accounts of those assumed theorems instead of simply reproducing the rest of the existing arguments.

Hence, for the moment, we postpone remarks and discussions about the Diophantine-ness of the exponential function and the Diophantine admissibility of bounded universal quantification to Section 5.

<footnote n="9">⁹Notice that the actual implemented proofs might differ slightly because we sometimes optimize the shape of expressions for smaller witnesses, especially for these basic shapes which pop up over and over again.</footnote>
**Theorem 2.7 (Exponential).** $\forall fg : (V \to N) \to N, D_{\text{fun}} f \to D_{\text{fun}} g \to D_{\text{fun}} (\lambda\nu.(f\nu)^{g\nu})$.

To prove it, one needs a meta-level Diophantine shape for the exponential relation, the *proof of which is nothing short of extraordinary*. This landmark result is due to Matiyasevich [Mat70], but we have implemented the shorter and more up-to-date proof of [Mat00].

**Theorem 2.8 (Bounded Univ. Quantification).** For $f : (V \to N) \to N$ and $T : N \to (V \to N) \to P$, we have $D_{\text{fun}} f \to D_{\text{rel}}(\lambda\nu.T(\nu x_0) (\lambda x_i.\nu x_{1+i})) \to D_{\text{rel}}(\lambda\nu.\forall u, u < f\nu \to Tu\nu)$.

This map can be compared with map 3 of Proposition 2.2 and allows to recognise bounded universal quantification as a legitimate Diophantine shape. We have implemented the direct proof of Matiyasevich [Mat97] which does not involve a detour through a model of computation. Notice that the bound $f\nu$ in $\forall u, u < f\nu \to \dots$ is not assumed constant otherwise the elimination of the quantifier would proceed as a simple reduction to a finitary conjunction.

**2.4. Reflexive-Transitive Closure is Diophantine.** With these tools – elimination of the exponential relation and of bounded universal quantification – we can show that the reflexive and transitive closure of a Diophantine binary relation is itself Diophantine. We assume a binary relation $R : N \to N \to P$ over natural numbers. The Diophantineness of $R$ can be formalised by assuming that e.g. $\lambda\nu.R(\nu x_1) (\nu x_0)$ is a Diophantine relation. We show that the $i$-th iterate of $R$ is Diophantine (where $i$ is non-constant).

**Lemma 2.9.** For any binary relation $R : N \to N \to P$ and $f,g,h : (V \to N) \to N$, we have

$$D_{\text{fun}} f \to D_{\text{fun}} g \to D_{\text{fun}} h \to D_{\text{rel}}(\lambda\nu.R(\nu x_1) (\nu x_0)) \to D_{\text{rel}}(\lambda\nu.R^{h\nu}(f\nu) (g\nu)).$$

*Proof.* Using Euclidean division, we define the `is.`*digit* `c q n d` predicate stating that $d$ is the $n$-th digit of the base $q$ development of number $c$, as a Diophantine sentence:

$$\texttt{is\_digit c q n d} := d < q \land \exists a b, c = (aq + d)q^n + b \land b < q^n.$$

The Diophantineness of this follows from previous Diophantine shapes, including the exponential (Theorem 2.7). Then we define the `is_seq` $R c q i$ predicate stating that the first $i+1$ digits of $c$ in base $q$ form an $R$-chain, again with a Diophantine shape, established using hypothesis $D_{\text{rel}}(\lambda\nu.R(\nu x_1) (\nu x_0))$ and the Diophantine admissibility of bounded universal quantification (Theorem 2.8):

$$\texttt{is\_seq R c q i} := \forall n, n < i \to \exists u v, \texttt{is\_digit c q n u} \land \texttt{is\_digit c q (1+n) v} \land R u v$$

Then we encode $R^i u v$ by stating that there exists a (large enough) $q$ and a number $c$ such that the first $i+1$ digits of $c$ in base $q$ form an $R$-chain starting at $u$ and ending at $v$:

$$R^i u v \leftrightarrow \exists q c, \texttt{is\_seq R c q i} \land \texttt{is\_digit c q 0} \land \texttt{is\_digit c q i v}$$

and this expression is accepted as a Diophantine shape by `dio auto`. Then assuming Diophantineness of $f, g$ and $h$, we easily derive that $\lambda\nu.R^{h\nu}(f\nu) (g\nu)$ is Diophantine. $\square$

We fill in Lemma 2.9 in the Diophantine hint database and we derive the Diophantineness of the reflexive-transitive closure as a direct consequence of the equivalence $R^* u v \leftrightarrow \exists i, R^i u v$.

**Theorem 2.10 (RT-closure).** For any binary relation $R : N \to N \to P$, we have the map
$$\forall fg : (V \to N) \to N, D_{\text{fun}} f \to D_{\text{fun}} g \to D_{\text{rel}}(\lambda\nu.R(\nu x_1) (\nu x_0)) \to D_{\text{rel}}(\lambda\nu.R^*(f\nu) (g\nu)).$$
### 3. ELEMENTARY DIOPHANTINE CONSTRAINTS

We now shift to another, seemingly less expressive, object-level representation of Diophantine relations. Elementary Diophantine constraints are very simple equations where only one instance of either $\dagger$ or $\dot{x}$ is allowed. Schematically, starting from Diophantine logic, we remove disjunction and existential quantification and encode conjunctions into the structure of a list. We give a direct proof that any Diophantine logic formula is semantically equivalent to the simultaneous satisfiability of a list of elementary Diophantine constraints.

Starting from two copies of $\mathbb{N}$, one called $\mathbb{U}$ with $u, v, w$ ranging over $\mathbb{U}$ for existentially quantified variables, and another one $\mathbb{V} = \{x_0, x_1, \dots\}$ for parameters, we define the type of elementary Diophantine constraints by:<fnref n="10" />

$$ c : \mathbb{D}_{\text{cstr}} ::= u \doteq n \mid u \doteq v \mid u \doteq x_i \mid u \doteq v \dagger w \mid u \doteq v \dot{x} w \quad \text{where } n : \mathbb{N} $$

Notice that these constraints do not have a “real” inductive structure, they are flat and of size (number of symbols) either 3 or 5. Given two interpretations, $\varphi : \mathbb{U} \to \mathbb{N}$ for variables and $\nu : \mathbb{V} \to \mathbb{N}$ for parameters, it is trivial to define the semantics $[[c]]_\nu^\varphi : \mathbb{P}$ of a single constraint $c$ of type $\mathbb{D}_{\text{cstr}}$:

$$ \begin{align*} [[u \doteq n]]_\nu^\varphi &:= \varphi u = n & [[u \doteq v]]_\nu^\varphi &:= \varphi u = \varphi v & [[u \doteq v \dagger w]]_\nu^\varphi &:= \varphi u = \varphi v + \varphi w \\ [[u \doteq x_i]]_\nu^\varphi &:= \varphi u = \nu x_i & [[u \doteq v \dot{x} w]]_\nu^\varphi &:= \varphi u = \varphi v \times \varphi w \end{align*} $$

Given a list $l : \mathbb{L}\mathbb{D}_{\text{cstr}}$ of constraints, we write $[[l]]_\nu^\varphi$ when all the constraints in $l$ are simultaneously satisfied, i.e. $[[l]]_\nu^\varphi := \forall c, c \in l \to [[c]]_\nu^\varphi$. We show the following result:

**Theorem 3.1.** For any Diophantine formula $A : \mathbb{D}_{\text{form}}$ one can compute a list of elementary Diophantine constraints $l : \mathbb{L}\mathbb{D}_{\text{cstr}}$ such that $\forall \nu : \mathbb{V} \to \mathbb{N}, [[A]]_\nu \leftrightarrow \exists \varphi : \mathbb{U} \to \mathbb{N}, [[l]]_\nu^\varphi$.

Put in other terms, for any given interpretation $\nu$ of parameters, $[[A]]_\nu$ holds if and only if the constraints in $l$ are simultaneously satisfiable. Any Diophantine logic formula is equivalent to the satisfiability of the conjunction of finitely many elementary Diophantine constraints. The proof of Theorem 3.1 spans the rest of this section. We will strengthen the result a bit to be able to get an easy argument by induction on $A$.

**Definition 3.2.** Given a relation $R : (\mathbb{V} \to \mathbb{N}) \to \mathbb{P}$ and an interval $[u_a, u_{a+n}] \subseteq \mathbb{U}$, an elementary representation of $R$ in $[u_a, u_{a+n}]$ is given by:

(1) a list $\mathcal{E} : \mathbb{L}\mathbb{D}_{\text{cstr}}$ of constraints and a reference variable $\mathfrak{r} : \mathbb{U}$;

(2) proofs that $\mathfrak{r}$ and the (existentially quantified) variables of $\mathcal{E}$ belong to $[u_a, u_{a+n}]$;

(3) a proof that the constraints in $\mathcal{E}$ are always (simultaneously) satisfiable, i.e. $\forall \nu \exists \varphi [[\mathcal{E}]_\nu^\varphi]$;

(4) a proof that the list $(\mathfrak{r} \doteq 0) :: \mathcal{E}$ is equivalent to $R$, i.e. $\forall \nu, R \nu \leftrightarrow (\exists \varphi, \varphi \mathfrak{r} = 0 \wedge [[\mathcal{E}]_\nu^\varphi])$.

It is obvious that an elementary representation of $\lambda\nu.[[A]]_\nu$ in any interval $[u_a, u_{a+n}]$ is enough to prove Theorem 3.1 because of item 4 of Definition 3.2. But actually, computing such a representation is simpler than proving Theorem 3.1 directly.<fnref n="11" /> Below, the size of a Diophantine formula $A$, denoted $|A|$, is defined as the number of nodes of its syntactic tree.

**Lemma 3.3.** For any $a : \mathbb{N}$ and any $A : \mathbb{D}_{\text{form}}$, one can compute $n \le 8|A|$ and an elementary representation of the relation $\lambda\nu.[[A]]_\nu$ in $[u_a, u_{a+n}]$.

<footnote n="10">¹⁰The equation $u \doteq v$ is redundant because it could be replaced with $z \doteq 0 \dot{\wedge} u \doteq z \dagger v$, for some fresh $z$. However we keep $u \doteq v$ in the syntax because this simplifies arguments when parameters $x_i$ are mapped to existential variables $v$ in the proof of Lemma 3.3, the type $\mathbb{D}_{\text{cstr}}$ being thus closed under this transformation.</footnote>

<footnote n="11">¹¹Proving Theorem 3.1 directly involves renamings of existential variables and might produce exponential blow-up in the number of constraints when handled naively.</footnote>
*Proof.* We show the result by structural induction on A.

* If A is e.g. $x_i \doteq x_j \dot{\times} x_k$, we get the representation with $n := 8$ and the pair $(\mathcal{E}, \mathbf{r})$ with

$$ \mathcal{E} := \left[ \begin{array}{l} u_{a+7} \doteq u_a \doteq u_{a+1} ; \quad u_{a+6} \doteq u_a \doteq u_{a+2} ; \quad u_{a+6} \doteq u_{a+1} \doteq u_{a+5} ; \\ u_{a+5} \doteq u_{a+3} \dot{\times} u_{a+4} ; \quad u_{a+4} \doteq x_k ; \quad u_{a+3} \doteq x_j ; \quad u_{a+2} \doteq x_i \end{array} \right] \qquad \mathbf{r} := u_{a+7} $$

Property (2) is obviously satisfied. Property (3), i.e. the satisfiability of $\mathcal{E}$ whatever the values of $x_i$, $x_j$ and $x_k$, is simple to establish: indeed, the values of $u_{a+2} = x_i$, $u_{a+3} = x_j$, $u_{a+4} = x_k$ and $u_{a+5} = x_j x_k$ are uniquely determined. Pick $u_a := x_j x_k$ and $u_{a+1} := x_i$ and then again $u_{a+6} = u_{a+7} = x_i + x_j x_k$ are both uniquely determined. This assignment of variables satisfies all the constraints in $\mathcal{E}$.

For property (4), let us now add the extra constraint $u_{a+7} \doteq 0$ and consider a valuation satisfying the constraints in $(u_{a+7} \doteq 0) :: \mathcal{E}$. We must have $u_{a+7} = u_a + u_{a+1} = 0$, hence $u_a = u_{a+1} = 0$ which entails $u_{a+6} = u_{a+2} = u_{a+5}$. But $u_{a+2} = x_i$ and $u_{a+5} = u_{a+3} u_{a+4} = x_j x_k$. Hence $x_i \doteq x_j \dot{\times} x_k$ is satisfied.

Conversely, parameter values satisfying $x_i \doteq x_j \dot{\times} x_k$ can be extended to a valuation of variables satisfying $(u_{a+7} \doteq 0) :: \mathcal{E}$ in a unique way with $u_{a+7}, u_a, u_{a+1} := 0, u_{a+6}, u_{a+5}, u_{a+2} := x_i, u_{a+3} := x_j$ and $u_{a+4} := x_k$;

* We proceed similarly for the other atomic cases $x_i \doteq n$, $x_i \doteq x_j$ and $x_i \doteq x_j + x_k$;

* When A is $B \dot{\wedge} C$, we get a representation in $[u_a, u_{a+n_A}]$ by induction. Hence, let $(\mathcal{E}_B, \mathbf{r}_B)$ be the representation of B in $[u_a, u_{a+n_B}]$. Then, inductively again, let $(\mathcal{E}_C, \mathbf{r}_C)$ be a representation of C in $[u_{a+n_B}, u_{a+n_B+n_C}]$. We define $\mathbf{r}_A := u_{a+n_A+n_B}$ and $\mathcal{E}_A := (\mathbf{r}_A \doteq \mathbf{r}_B + \mathbf{r}_C) :: \mathcal{E}_B \supset \mathcal{E}_C$ and then $(\mathcal{E}_A, \mathbf{r}_A)$ represents $A = B \dot{\wedge} C$ in $[u_a, u_{a+1+n_B+n_C};^{12}$

* The case of $B \dot{\wedge} C$ is similar: simply replace $\mathbf{r}_A \doteq \mathbf{r}_B + \mathbf{r}_C$ with $\mathbf{r}_A \doteq \mathbf{r}_B \dot{\times} \mathbf{r}_C$;

* We finish with the case when A is $\exists B$. Let $(\mathcal{E}_B, \mathbf{r}_B)$ be a representation of B in $[u_a, u_{a+n_B}]$. Let $\sigma$ be the substitution mapping parameters in V and defined by $\sigma(x_0) := u_{a+n_B}$ and $\sigma(x_{1+i}) := x_i$; existential variables in U are left unmodified by this substitution. Then $(\sigma(\mathcal{E}_B), \mathbf{r}_B)$ is a representation of $A = \exists B$ in $[u_a, u_{a+1+n_B}]$.

This concludes the recursive construction of a representation of $\lambda\nu.[[A]]_\nu$. $\square$

## 4. SINGLE DIOPHANTINE EQUATIONS

We now give our last and most naive object-level representation of Diophantine relations as a single polynomial equation. In this section, we show how a list of elementary Diophantine constraints can be simulated by a single identity between two Diophantine polynomials. We use the following well known convexity identity to achieve the reduction.

**Proposition 4.1.** Let $(p_1, q_1), \dots, (p_n, q_n)$ be a sequence of pairs in $\mathbb{N} \times \mathbb{N}$. Then

$$ \sum_{i=1}^{n} 2p_i q_i = \sum_{i=1}^{n} p_i^2 + q_i^2 \iff p_1 = q_1 \wedge \dots \wedge p_n = q_n. $$

There are many possible justifications for this equivalence which, if ranging over $\mathbb{Z}$, would more conventionally be written as $\sum_i (p_i - q_i)^2 = 0 \iff \forall i, p_i = q_i$. This form does not hold

<footnote n="12">¹²Since the intervals $[u_a, u_{a+n_B}]$ and $[u_{a+n_B}, u_{a+n_B+n_C}]$ are built disjoint, there is no difficulty in merging valuations whereas this usually involves renamings when existential variables are not carefully chosen.</footnote>
directly for $\mathbb{N}$ however, hence the reformulation without using the minus/− operator.<fnref n="13" /> See Appendix D for an elementary justification of the result.

Similarly to elementary Diophantine constraints, we define Diophantine polynomials distinguishing the types of $U$ of bound variables and $V$ of parameters (or free variables) but the types $U$ and $V$ are not fixed copies of $\mathbb{N}$ anymore, but type parameters of arbitrary value.

**Definition 4.2.** The type of Diophantine polynomials $\mathbb{D}_{\text{poly}}(U, V)$ and the type of single Diophantine equations $\mathbb{D}_{\text{single}}(U, V)$ are defined by:

$$p, q : \mathbb{D}_{\text{poly}}(U, V) ::= u : U \mid x : V \mid n : \mathbb{N} \mid p \dot{=} q \mid p \dot{=} q \qquad E : \mathbb{D}_{\text{single}}(U, V) ::= p \dot{=} q.$$

For $\varphi : U \to \mathbb{N}$ and $\nu : V \to \mathbb{N}$ we define the semantic interpretations of polynomials $[[p]]_{\nu}^{\varphi} : \mathbb{N}$ and single Diophantine equations $[[p] \doteqdot q]]_{\nu}^{\varphi} : \mathbb{P}$ in the obvious way.

**Theorem 4.3.** For any list $l : \mathbb{L}\mathbb{D}_{\text{cstr}}$ of elementary Diophantine constraints, one can compute a single Diophantine equation $E : \mathbb{D}_{\text{single}}(\mathbb{N}, \mathbb{N})$ such that $\forall \nu \forall \varphi, [[E]]_{\nu}^{\varphi} \leftrightarrow [[l]]_{\nu}^{\varphi}$.

*Proof.* We write $l = [p_1 \doteqdot q_1; \dots; p_n \doteqdot q_n]$ and then use Proposition 4.1. In the code, we moreover show that the size of $E$ is linear in the length of $l$. If needed, one could also show that the degree of the polynomial is less than 4. $\square$

**Corollary 4.4.** Let $R : (\mathbb{N} \to \mathbb{N}) \to \mathbb{P}$. Assuming $\mathbb{D}_{\text{rel}} R$, one can compute a single Diophantine equation $p \doteqdot q : \mathbb{D}_{\text{single}}(\mathbb{N}, \mathbb{N})$ such that $\forall \nu, R \nu \leftrightarrow \exists \varphi, [[p]]_{\nu}^{\varphi} = [[q]]_{\nu}^{\varphi}$.

*Proof.* Direct combination of Definition 2.1 and Theorems 3.1 and 4.3. In the formalisation, we also show that the size of the obtained single Diophantine equation is linearly bounded by the size of the witness formula contained in the proof of $\mathbb{D}_{\text{rel}} R$. $\square$

We have shown that the automation we designed to recognise relations of Diophantine shape entail that these relations are also definable by satisfiability of a single equation between Diophantine polynomials, so these tools are sound w.r.t. a formally restrictive characterisation of Diophantineness. One could argue that the above existential quantifier $\exists \varphi$ encodes infinitely many existential quantifiers but it can easily be replaced by finitely many existential quantifiers over the bound variables that actually occur in $p$ or $q$.

**Proposition 4.5.** For any single Diophantine equation $p \doteqdot q : \mathbb{D}_{\text{single}}(\mathbb{N}, V)$, one can compute $n : \mathbb{N}$ and a new single Diophantine equation $p' \doteqdot q' : \mathbb{D}_{\text{single}}(\mathbb{F}_n, V)$ such that for any $\nu : V \to \mathbb{N}$, $(\exists \varphi : \mathbb{N} \to \mathbb{N}, [[p]]_{\nu}^{\varphi} = [[q]]_{\nu}^{\varphi}) \leftrightarrow (\exists \varphi : \mathbb{F}_n \to \mathbb{N}, [[p']_{\nu}^{\varphi} = [[q']_{\nu}^{\varphi}})$.

*Proof.* We pick $n$ greater that the number of bound variables which occur in either $p$ or $q$. This subset of $\mathbb{N}$ can be faithfully embedded into the finite type $\mathbb{F}_n$ and we use such a renaming to compute $(p', q')$. Remark that the size of $(p', q')$ is the same as that of $(p, q)$. $\square$

By Corollary 4.4 and Proposition 4.5, we see that a Diophantine logic formula $A : \mathbb{D}_{\text{form}}$ potentially containing inner existential quantifiers and representing the Diophantine relation $\lambda\nu. [[A]]_{\nu}$ can effectively be reduced to a single Diophantine equation $p' \doteqdot q' : \mathbb{D}_{\text{single}}(\mathbb{F}_n, V)$ such that $[[A]]_{\nu} \leftrightarrow \exists \varphi : \mathbb{F}_n \to \mathbb{N}, [[p']_{\nu}^{\varphi} = [[q']_{\nu}^{\varphi}}$. Because $\mathbb{F}_n$ is the finite type of $n$ elements, the (higher order) existential quantifier $\exists \varphi$ simply encodes $n$ successive (first order) existential quantifiers. The existential quantifiers that occur deep inside $A$ are not erased by the reduction, they are moved at the outer level and to be ultimately understood as solvability for some polynomial equation of which the parameters match the free variables of $A$.

<footnote n="13">¹³Remember that $0 - 1 = 0$ holds in $\mathbb{N}$.</footnote>
## 5. REMARKS ON THE IMPLEMENTATION OF MATIYASEVICH'S THEOREMS

*Matiyasevich's theorem* stating that there is a Diophantine description of the exponential relation $x = y^z$ is a masterpiece which concluded the line of work by Davis, Putnam and Robinson, starting at Davis's conjecture in 1953. Already in 1952, Julia Robinson discovered that in order to show the exponential relation Diophantine, it suffices to find a single binary Diophantine relation exhibiting exponential growth [Rob52], a so-called Robinson predicate, i.e. a predicate $J(u, v)$ in two variables s.t. $J(u, v)$ implies $v < u^u$ and for every $k$ there are $u, v$ with $J(u, v)$ and $v > u^k$. Robinson's insight meant the only thing missing to prove what is nowadays called the DPRM theorem, was a *single* polynomial equation capturing *any* freely chosen Robinson predicate. Similar to other famous hard problems of mathematics, the question is easy to state, but from the start of the study of Diophantine equations to the late 60s, no such relation was known, rendering the problem one of the most baffling questions for mathematicians and computer scientists alike.

In 1970, Yuri Matiyasevich [Mat70] discovered that $v = \text{fib}_{2u}$ is both a Robinson predicate and Diophantine. Here $(\text{fib}_n)_{n \in \mathbb{N}}$ is the well known Fibonacci sequence defined by the second order recurrence relation $\text{fib}_0 = 0$, $\text{fib}_1 = 1$ and $\text{fib}_{n+2} = \text{fib}_{n+1} + \text{fib}_n$. Combined with previous results, this concluded the multi-decades effort to establish the Diophantineness of all recursively enumerable predicates, implying a negative solution to Hilbert's tenth problem. That proof which included the original proof of Matiyasevich [Mat70] was later simplified. For instance, exploiting similar ideas but in the easier context of the solutions of another second order equation – namely Pell's equations $x^2 - (a^2 - 1)y^2 = 1$ with parameter $a > 1$, – Martin Davis [Dav73] gave a standalone proof of the DPRM-theorem where recursively enumerable predicates are characterised by a variant of $\mu$-recursive functions. In that paper, Davis also provided a proof of the admissibility of bounded universal quantification using the Chinese remainder theorem to encode finite sequences of numbers. There exists more recent and simpler proofs of this admissibility result as well, see e.g. [Mat97].

Before we discuss the mechanisation of the Diophantineness of both the exponential relation and of bounded universal quantification, we want to remark on the difficulty of mechanising the former proof. Both on its own and as a stepping stone towards the negative solution to Hilbert's tenth problem, it is clear that Matiyasevich's theorem was an extremely difficult question which required superior intellectual resources to be solved. The mechanisation of a modernised form of the proof, although not trivial, cannot be compared to the difficulty of finding a solution. In particular, the modern proof relies on very mature background theories, lowering the number of possible design choices for the mechanisation. Moreover, very detailed pen and paper accounts of the proof are available, which can be followed closely.

An aspect that is more challenging in mechanisation than on paper are proofs regarding the computability of certain functions. Since paper proofs oftentimes rely on a vague notion of algorithm, most of the reasoning about these algorithms is hand-waved away by computer scientists, relying on the implicit understanding of what is an algorithm. By using a synthetic approach to computability [FHS18, FKS19, FL19], we make the notion of an algorithm precise and thus enable mechanisation, at the same time circumventing the verification of low-level programs.

### 5.1. Exponential is Diophantine (Theorem 2.7)
For our mechanised proof, we rely on a more recent account of Matiyasevich's theorem from [Mat00], which, among the many
options we considered, seemed the shortest. The proof employs the equation $x^2 - bxy + y^2 = 1$ for $b \ge 2$, also called Pell's equation in [Dav73]. We use the second order recurrence relation $\alpha_b(-1) = -1$, $\alpha_b(0) = 0$ and $\alpha_b(n+2) = b\alpha_b(n+1) - \alpha_b(n)$ to describe the set of solutions of Pell's equation by $\{(\alpha_b(n), \alpha_b(n+1)) | n \in \mathbb{N}\}$. The recurrence can be characterised by the following square $2 \times 2$ matrix equation:

$$A_b(n) = (B_b)^n \quad \text{with} \quad A_b(n) := \begin{pmatrix} \alpha_b(n+1) & -\alpha_b(n) \\ \alpha_b(n) & -\alpha_b(n-1) \end{pmatrix} \quad \text{and} \quad B_b := \begin{pmatrix} b & -1 \\ 1 & 0 \end{pmatrix}$$

Then, studying the properties of the sequence $n \mapsto \alpha_b(n)$ in $\mathbb{N}$ or $\mathbb{Z}$, one can show that $\alpha_2(n) = n$ for all $n$ and $n \mapsto \alpha_b(n)$ grows exponentially for $b \ge 3$. Studying the properties of the same sequence in $\mathbb{Z}/p\mathbb{Z}$ (for varying values of the modulus $p$), one can for instance show that $n = \alpha_2(n) \equiv \alpha_b(n) [b-2]$, which relates $n$ and $\alpha_b(n)$ modulo $(b-2)$. With various intricate but elementary results;<fnref n="14" /> such as e.g. $\alpha_b(k) | \alpha_b(m) \leftrightarrow k | m$ and $\alpha_b^2(k) | \alpha_b(m) \leftrightarrow k\alpha_b(k) | m$ (both for $b \ge 2$ and any $k, m \in \mathbb{N}$), one can show that $a, b, c \mapsto 3 < b \land a = \alpha_b(c)$ has a Diophantine representation. In our formalisation, we get a Diophantine logic formula of size 1445 as a witness (see dio_rel_alpha_size).

Once $\alpha_b(n)$ is proven Diophantine, one can recover the exponential relation $x, y, z \mapsto x = y^z$ using the eigenvalue $\lambda$ of the matrix $B_b$ which satisfies $\lambda^2 - b\lambda - 1 = 0$. By wisely choosing $m := bq - q^2 - 1$, one gets $\lambda \equiv q[m]$ and thus, using the corresponding eigenvector, one derives $q\alpha_b(n) - \alpha_b(n-1) \equiv q^n[m]$. For a large enough value of $m$, hence a large enough value<fnref n="15" /> of $b$, this gives a Diophantine representation of $q^n$. In our code, we get a Diophantine logic formula of size 4903 as a witness (see dio fun expo example_size).

The main libraries which are needed to solve Pell's equation and characterise its solutions are linear algebra (or at least square $2 \times 2$ matrices) over commutative rings such as $\mathbb{Z}$ and $\mathbb{Z}/p\mathbb{Z}$, a good library for modular algebra ($\mathbb{Z}/p\mathbb{Z}$), and the binomial theorem over rings. Without the help of the Coq ring tactic, such a development would be extremely painful. These libraries are then used again to derive the Diophantine encoding of the exponential.

## 5.2. Admissibility of Bounded Universal Quantification (Theorem 2.8).
As hinted earlier, we provide an implementation of the algorithm for the elimination of bounded universal quantification described in [Mat97]. It does not involve the use of a model of computation, hence does not create a chicken-and-egg problem when used for the proof of the DPRM theorem. The technique of [Mat97] uses the exponential function and thus Theorem 2.7 (a lot), and a combination of arithmetic and bitwise operations over $\mathbb{N}$ through base 2 and base $2^q$ representations of natural numbers.

The Diophantine admissibility of bitwise operations over $\mathbb{N}$ is based on the relation stating that every bit of $a$ is lower or equal than the corresponding bit in $b$ and denoted $a \le b$. The equation $a \le b \leftrightarrow C_b^a$ is odd<fnref n="16" /> gives a Diophantine representation for $a \le b$ and then bitwise operators are derived from $\le$ in combination with regular addition $+$, in particular, the digit by digit AND operation called “projection”. To obtain that $a \le b$ holds if and only if $C_b^a = 1$ [2], we prove Lucas’s theorem [Luc78] which allows for the computation of the

<footnote n="14">¹⁴by elementary we certainly do not mean either simple or obvious, but we mean that they only involve standard tools from modular and linear algebra.</footnote>

<footnote n="15">¹⁵the largeness of which is secured using $\alpha$ itself again, but with other input values. But this works only in the case where $n > 0$ and $q > 0$. The cases where $n = 0$ (and hence $q^n = 1$) or $q = 0$ and $n > 0$ (and hence $q^n = 0$) are trivial and treated separately.</footnote>

<footnote n="16">¹⁶where $C_b^a$ denotes the binomial coefficient with the usual convention that $C_b^a = 0$ when $a > b$.</footnote>
binomial coefficient in base $p$. It states that $C_b^a \equiv C_{b_n}^{a_n} \times \cdots \times C_{b_0}^{a_0} [p]$ holds when $p$ is prime and $a = a_n p^n + \cdots + a_0$ and $b = b_n p^n + \cdots + b_0$ are the respective base $p$ representations of $a$ and $b$; see Appendix B for an elementary combinatorial proof of Lucas's theorem.

A Diophantine representation of the binomial coefficient can be obtained via e.g. the binomial theorem: $C_n^k$ is the $k$-th digit of the development of $(1+q)^n = \sum_{i=0}^n C_n^i q^i$ in base $q=2^{n+1}$. This gives a Diophantine representation using the `is.`digit relation of Lemma 2.9.

The rest of the admissibility proof for bounded universal quantification $\forall i, i < n \rightarrow A$ is a very nice encoding of vectors of natural numbers of type $\mathbb{N}^n$ into natural numbers $\mathbb{N}$ such that regular addition $+$ (resp. multiplication $\times$) somehow performs parallel/simultaneous additions (resp. multiplications) on the encoded vectors. More precisely, a vector $(a_1, \dots, a_n) \in [0, 2^q - 1]^n$ of natural numbers is encoded as the “cipher” $a_1 r^2 + a_2 r^4 + a_3 r^8 + \dots + a_n r^{2^n}$ with $r = 2^{4q}$. In these *sparse ciphers*, only the digits occurring at $r^{2^i}$ are non-zero. We remark that none of the parameters, including $n$ or $q$, are constant in the encoding.

Besides the low-level inductive proof of Lucas's theorem presented in Appendix B, the essential library for the removal of bounded universal quantification consists of tools to manipulate the type $\mathbb{N}$ simultaneously and smoothly both as (a) usual natural numbers and (b) sparse base $r = 2^{4q}$ encodings of vectors of natural numbers in $[0, 2^q - 1]$. Notice that $r$ is defined as $r = 2^{2q}$ in [Mat97] but we favour the alternative choice $r = 2^{4q}$ which allows for an easier soundness proof for vector multiplication because there is no need to manage for digit overflows (see Appendix C).

A significant step in the Diophantine encoding of $+$ and $\times$ on $\mathbb{N}^n$ is the Diophantine encoding of $u = \sum_{i=1}^n r^{2^i}$ and $u' = \sum_{i=2}^{n+1} r^{2^i}$ as the ciphers of the constant vectors $[1; \dots; 1] \in \mathbb{N}^n$ and $[0; 1; \dots; 1] \in \mathbb{N}^{n+1}$ respectively, obtained by masking $u^2$ with $w = \sum_{i=0}^{2^{n+1}} r^i$ and $2w$.

Finally, it should be noted that prior to the elimination of the quantifier in $\forall i, i < n \rightarrow A$, the Diophantine formula $A$ is first normalised into a conjunction of elementary constraints using Theorem 3.1, and then the elimination is performed on that list of elementary constraints, encoding e.g. $v_0 \doteq v_1 \dotplus v_2$ and $v_0 \doteq v_1 \dottimes v_2$ with their respective sparse cipher counterparts.

## 6. MINSKY MACHINES REDUCE TO FRACTRAN

In previous work, we have reduced the halting problem for Turing machines to PCP [FHS18] and on to a specialised halting problem for Minsky machines [FL19] in Coq. The specialised halting problem asked whether a machine on a given input halts in a configuration with all registers containing zeros. In order to define Minsky machine recognisability, we consider a general halting problem which allows *any* final configuration, final meaning that computation cannot further proceed. The adaptation of the formal proofs reducing PCP via binary stack machines to Minsky machines is quite straightforward and reuses the certified compiler for low-level languages defined in [FL19].

We first show that one can remove self loops from Minsky machines, i.e. instructions which jump to their own location, using the compositional reasoning techniques developed in [FL19]. We then formalise the FRACTRAN language [Con87] and show how the halting problem for Minsky machines can be encoded into the halting problem for FRACTRAN programs. While the verification of Minsky machines can be complex and needs preliminary thoughts on compositional reasoning, the translation from Minsky machines to FRACTRAN is elementary and needs no heavy machinery.
**6.1. Minsky Machines.** We employ Minsky machines [Min67] with instructions $\iota : I_n ::= INC (\alpha : \mathbb{F}_n) \mid DEC (\alpha : \mathbb{F}_n) (p : \mathbb{N})$. A Minsky machine with $n$ registers is a sequence of consecutively indexed instructions $s : \iota_0; \dots s + k : \iota_k$; represented as a pair $(s : \mathbb{N}, [\iota_0; \dots ; \iota_k] : \mathbb{L} I_n)$. Its state $(i, \vec{v})$ is a program counter (PC) value $i : \mathbb{N}$ and a vector of values for registers $\vec{v} : \mathbb{N}^n$. INC $\alpha$ increases the value of register $\alpha$ and the PC by one. DEC $\alpha$ $p$ decreases the value of register $\alpha$ by one if that is possible and increases the PC, or, if the register is already 0, jumps to PC value $p$. Given a Minsky machine $(s, P)$, we write $(s, P) \llbracket_M (i_1, \vec{v}_1) \succcurlyeq^k (i_2, \vec{v}_2)$ when $(s, P)$ transforms state $(i_1, \vec{v}_1)$ into $(i_2, \vec{v}_2)$ in $k$ steps of computation. For $(s, P)$ to do a step in state $(i, \vec{v})$ the instruction at label $i$ in $(s, P)$ is considered. We define out $i$ $(s, P) := i < s \lor |P| + s \le i$ to characterize when label $i$ is outside of the code of $(s, P)$. In that case (and only that case), no computation step can occur. We define the halting problem for Minsky Machines as

$$ \mathbf{MM}(n : \mathbb{N}, P : \mathbb{L} I_n, \vec{v} : \mathbb{N}^n) := (1, P) \llbracket_M (1, \vec{v}) \downarrow \\ \text{where } \begin{cases} (s, P) \llbracket_M (i, \vec{v}) \rightsquigarrow (j, \vec{w}) := \exists k, (s, P) \llbracket_M (i, \vec{v}) \succcurlyeq^k (j, \vec{w}) \land \text{out } j (s, P) \\ (s, P) \llbracket_M (i, \vec{v}) \downarrow := \exists j \vec{w}, (s, P) \llbracket_M (i, \vec{v}) \rightsquigarrow (j, \vec{w}) \end{cases} $$

meaning that the machine $(1, P)$ has a terminating computation starting at state $(1, \vec{v})$, the value of the final state being irrelevant. Notice that since $\mathbf{MM}(n, P, \vec{v})$ considers only Minsky machines starting at PC value 1, the PC value 0 is always outside of their code (i.e. out 0 $(1, P)$ holds), hence computations can be halted by jumping there. We refer to [FL19] for a more in-depth formal description of those counter machines. Note that the halting problem defined there as $\mathbf{MM}_0(n : \mathbb{N}, P : \mathbb{L} I_n, \vec{v} : \mathbb{N}^n) := (1, P) \llbracket_M (1, \vec{v}) \rightsquigarrow (0, \vec{0})$ is more specific than the problem $\mathbf{MM}$ above defined but both are proved undecidable in our library.

We say that a machine *has a self loop* if it contains an instruction of the form $i : \mathbf{DEC} \ \alpha \ i$, i.e. jumps to itself in case the register $\alpha$ has value 0, leading necessarily to non-termination (in case the PC reaches value $i$). For every machine $P$ with self loops, we can construct an equivalent machine $Q$ using one additional register $\alpha_0$ with constant value 0, which has the same behaviour but no self loops. Since the effect of a self loop $i : \mathbf{DEC} \ \alpha \ i$ is either decrement and move to the next instruction at $i + 1$ if $\alpha > 0$ or else enter in a forever loop at $i$, it is easily simulated by a jump to a length-2 cycle, i.e. replacing $i : \mathbf{DEC} \ \alpha \ i$ with $i : \mathbf{DEC} \ \alpha \ j$ and adding $j : \mathbf{DEC} \ \alpha_0 (j + 1)$; $j + 1 : \mathbf{DEC} \ \alpha_0 \ j$ somewhere near the end of the program.

**Theorem 6.1.** Given a Minsky machine $P$ with $n$ registers one can compute a machine $Q$ with $1 + n$ registers and no self loops s.t. for any $\vec{v}$,

$$ (1, P) \llbracket_M (1, \vec{v}) \downarrow \leftrightarrow (1, Q) \llbracket_M (1, 0 \coloneqq \vec{v}) \downarrow. $$

*Proof.* We explain how any Minsky machine $(1, P)$ with $n$ registers can be transformed into an equivalent one that uses an extra 0 valued spare register $\alpha_0 = 0 \in \mathbb{F}_{1+n}$ and avoids self loops. Let $k$ be the length of $P$ and let $P'$ be the Minsky machine with $1 + n$ registers defined by performing a 1-1 replacement of instructions of $(1, P)$:

* instructions of the form $i : INC \ \alpha$ are replaced by $i : INC (1 + \alpha)$;

* self loops $i : DEC \ \alpha \ i$ are replaced by $i : DEC (1 + \alpha) (2 + k)$;

* proper inside jumps $i : DEC \ \alpha \ j$ for $i \neq j$ and $1 \le j \le k$ are replaced by $i : DEC (1 + \alpha) j$;

* and outside jumps $i : DEC \ \alpha \ j$ for $j = 0 \lor k < j$ are replaced by $i : DEC (1 + \alpha) 0$.

Then we define $Q := P' \oplus [DEC \ \alpha_0 \ 0; DEC \ \alpha_0 \ (3 + k); DEC \ \alpha_0 \ (2 + k)]$. Notice that $P'$ is immediately followed $DEC \ \alpha_0 \ 0$, i.e. by an unconditional jump to 0 (because $\alpha_0$ has value 0), and that $(1, Q)$ ends with the length-2 cycle composed of $2 + k : DEC \ \alpha_0 \ (3 + k); 3 + k : $
DEC $\alpha_0 (2 + k)$. We show that $(1, Q)$ is a program without self loops (obvious) that satisfies the required simulation equivalence. Indeed, self loops are replaced by jumps to the length-2 cycle that uses the unmodified register $\alpha_0$ to loop forever. One should just be careful that the outside jumps of $(1, P)$ do not accidentally fall into that cycle and this is why we redirect them all to PC value 0, which halts the computation because out 0 $(1, Q)$ holds. □

A predicate $R : \mathbb{N}^n \to \mathbb{P}$ is MM-recognisable if there exist $m : \mathbb{N}$ and a Minsky machine $P : \mathbb{L} |_{n+m}$ of $(n+m)$ registers such that for any $\vec{v} : \mathbb{N}^n$ we have $R \vec{v} \leftrightarrow (1, P) \//_M (1, \vec{v} \oplus \vec{0}) \downarrow$. The last $m$ registers serve as spare registers during the computation. Notice that not allowing for spare registers would make e.g. the empty predicate un-recognisable.<fnref n="17" /> It is possible to limit the number of (spare) registers but that question is not essential in our development.

## 6.2. The FRACTRAN language.
We formalise the language FRACTRAN, introduced as a universal programming language for arithmetic by Conway [Con87]. A FRACTRAN program $Q$ consists of a list of positive fractions $[p_1/q_1; \dots; p_n/q_n]$. The current state of a FRACTRAN program is just a natural number $s$. The first fraction $p_i/q_i$ in $Q$ such that $s \cdot (p_i/q_i)$ is still integral determines the successor state, which then is $s \cdot (p_i/q_i)$. If there is no such fraction in $Q$, the program terminates.

We make this precise inductively for $Q$ being a list of fractions $p/q : \mathbb{N} \times \mathbb{N}$:

$$ \frac{q \cdot y = p \cdot x}{(p/q :: Q) \//_F x \succ y} \qquad \frac{q \nmid p \cdot x}{(p/q :: Q) \//_F x \succ y} \qquad Q \//_F x \succ y $$

i.e. at state $x$ the first fraction $p/q$ in $Q$ where $q$ divides $p \cdot x$ is used, and $x$ is multiplied by $p$ and divided by $q$. For instance, the FRACTRAN program $[5/7; 2/1]$ runs forever when starting from state 7, producing the sequence $5 = 7 \cdot (5/7)$, $10 = 5 \cdot (2/1)$, $20 = 10 \cdot (2/1)$ ...<fnref n="18" />

We say that a FRACTRAN program $Q = [p_1/q_1; \dots; p_n/q_n]$ is regular if none of its denominators is 0, i.e. if $q_1 \neq 0, \dots, q_n \neq 0$. For a FRACTRAN program $Q : \mathbb{L} (\mathbb{N} \times \mathbb{N})$ and $s : \mathbb{N}$, we define the decision problem as the question “does $Q$ halt when starting from $s$”:

$$ \text{FRACTRAN}(Q, s) := Q \//_F s \downarrow \quad \text{with} \quad Q \//_F s \downarrow := \exists x, Q \//_F s \succ^* x \land \forall y, \neg Q \//_F x \succ y $$

Following [Con87], we now show how (regular) FRACTRAN halting can be used to simulate Minsky machines halting. The idea is to use a simple Gödel encoding of the states of a Minsky machine. We first fix two infinite sequences of prime numbers $p_0, p_1, \dots$ and $q_0, q_1, \dots$ all distinct from each other. We define the encoding of $n$-register Minsky machine states as $(i, \vec{v}) := p_i q_0^{x_0} \cdots q_{n-1}^{x_{n-1}}$ where $\vec{v} = [x_0, \dots, x_{n-1}]$:

* To simulate the step semantics of Minsky machines for $i : \text{INC} \alpha$, we divide the encoded state by $p_i$ and multiply by $p_{i+1}$ for the change in PC value, and increment the register $\alpha$ by multiplying with $q_\alpha$, hence we add the fraction $p_{i+1}q_\alpha/p_i$;

* To simulate $i : \text{DEC} \alpha j$ when $\vec{v}_\alpha = 1 + n$ we divide by $p_i$, multiply by $p_{i+1}$ and decrease register $\alpha$ by dividing by $q_\alpha$, hence we add the fraction $p_{i+1}/p_i q_\alpha$;

* To simulate $i : \text{DEC} \alpha j$ when $\vec{v}_\alpha = 0$ we divide by $p_i$ and multiply by $p_j$. To make sure that this is only executed when the previous rule does not apply, we add the fraction $p_j/p_i$ after the fraction $p_{i+1}/p_i q_\alpha$.

<footnote n="17">¹⁷For any Minsky machine $(1, P)$, if it starts on large enough register values, for instance if they are all greater than the length of $P$, then no jump can occur and the machine terminates after its last instruction executes. Such unfortunate behavior can be circumvented with a 0-valued spare register.</footnote>

<footnote n="18">¹⁸No FRACTRAN program can ever stop when it contains a fraction having an integer value like 2/1.</footnote>
In short, we define the encoding of labelled instructions and then programs as

$$ \overline{(i, \text{INC } \alpha)} := [\mathfrak{p}_{i+1} \mathfrak{q}_{\alpha} / \mathfrak{p}_{i}] \qquad \overline{(i, \text{DEC } \alpha j)} := [\mathfrak{p}_{i+1} / \mathfrak{p}_{i} \mathfrak{q}_{\alpha}; \mathfrak{p}_{j} / \mathfrak{p}_{i}] \qquad \overline{(i, [\iota_{0}; \dots; \iota_{k}])} := \overline{(i, \iota_{0})} + \dots + \overline{(i + k, \iota_{k})}. $$

Notice that we only produce regular programs and that a self loop like $i : \text{DEC } \alpha i$, jumping on itself when $\vec{v}_{\alpha} = 0$, will generate the fraction $\mathfrak{p}_{i}/\mathfrak{p}_{i}$ potentially capturing any state $(j, \vec{v})$ even when $j \neq i$. So this encoding does not work on Minsky machines containing self loops because the corresponding FRACTRAN program would never terminate, even when the PC never reaches self loops.

**Lemma 6.2.** If $(1, P)$ has no self loops then $(1, P) \bowtie_{M} (1, \vec{v}) \downarrow \leftrightarrow \overline{(1, P)} \bowtie_{F} \overline{(1, \vec{v})} \downarrow$.

*Proof.* Let $(i, P)$ be a Minsky machine with no self loops. We show that the simulation of $(i, P)$ by $(\overline{i, P})$ is 1-1, i.e. each step is simulated by one step. We first show the forward simulation, i.e. that $(i, P) \bowtie_{M} (i_1, \vec{v}_1) \succ (i_2, \vec{v}_2)$ entails $(\overline{i, P}) \bowtie_{F} (\overline{i_1, \vec{v}_1}) \succ (\overline{i_2, \vec{v}_2})$, by case analysis. Conversely we show that if $(i, P) \bowtie_{F} (i_1, \vec{v}_1) \succ st$ holds then $st = (i_2, \vec{v}_2)$ for some $(i_2, \vec{v}_2)$ such that $(i, P) \bowtie_{M} (i_1, \vec{v}_1) \succ (i_2, \vec{v}_2)$. Backward simulation involves the totality of MM one step semantics and the determinism of regular FRACTRAN one step semantics combined with the forward simulation.

Using these two simulation results, the desired equivalence follows by induction on the length of terminating computations. $\square$

**Theorem 6.3.** For any n-register Minsky machine P one can compute a regular FRACTRAN program Q s.t. $(1, P) \bowtie_{M} (1, [x_1; \dots; x_n]) \downarrow \leftrightarrow Q \bowtie_{F} \mathfrak{p}_1 \mathfrak{q}_1^{x_1} \cdots \mathfrak{q}_n^{x_n} \downarrow$ holds for any $x_1, \dots, x_n$.

*Proof.* Using Theorem 6.1, we first compute a Minsky machine $(1, P_1)$ equivalent to $(1, P)$ but with one extra 0-valued spare register and no self loops. Then we apply Lemma 6.2 to $(1, P_1)$ and let $Q := \overline{(1, P_1)}$. The program Q is obviously regular and given $\vec{v} = [x_1; \dots; x_n]$, the encoding of the starting state $(1, 0 :: \vec{v})$ for $(1, P_1)$ is $\mathfrak{p}_1 \mathfrak{q}_0^0 \mathfrak{q}_1^{x_1} \cdots \mathfrak{q}_n^{x_n}$ hence the result. $\square$

This gives us a formal constructive proof that (regular) FRACTRAN is Turing complete as a model of computation and is consequently undecidable.

**Corollary 6.4.** *Halt reduces to FRACTRAN.*

*Proof.* Theorem 6.3 gives us a reduction from MM to FRACTRAN which can be combined with the reduction of Halt to PCP from [FHS18] and a slight modification of PCP to MM from [FL19]. $\square$

# 7. DIOPHANTINE ENCODING OF FRACTRAN

We show that a single step of FRACTRAN computation is a Diophantine relation.

**Lemma 7.1.** For any FRACTRAN program $Q : \mathbb{L}(\mathbb{N} \times \mathbb{N})$, one can compute a map

$$ \forall fg : (V \to \mathbb{N}) \to \mathbb{N}, \mathbb{D}_{\text{fun}} f \to \mathbb{D}_{\text{fun}} g \to \mathbb{D}_{\text{rel}} (\lambda \nu. Q \bowtie_{F} f \nu \succ g \nu). $$

*Proof.* The map is built by induction on $Q$. If $Q = []$, then we show $[] \bowtie_{F} f \nu \succ g \nu \leftrightarrow \text{False}$, and thus $\mathbb{D}_{\text{rel}}(\lambda \nu.Q \bowtie_{F} f \nu \succ g \nu)$ by map 4 of Proposition 2.2 followed by dio auto. If $Q$ is a composed list $Q = p/q :: Q'$, then we show the equivalence

$$ ((p/q :: Q') \bowtie_{F} f \nu \succ g \nu) \leftrightarrow (q \cdot (g \nu) = p \cdot (f \nu)) \lor (q \nmid (p \cdot (f \nu))) \land (Q' \bowtie_{F} f \nu \succ g \nu) $$
and we derive $\mathbb{D}_{\text{rel}}(\lambda\nu.Q \//_F f\nu \succ g\nu)$ by map 4 of Proposition 2.2 followed by **dio auto**, the induction hypothesis being used locally as a hint for the tactic. □

In addition, the “Q has terminated at $x$” predicate is Diophantine for any FRACTRAN program $Q$. The proof is similar to the previous one:

**Lemma 7.2.** For any FRACTRAN program $Q : \mathbb{L}(\mathbb{N} \times \mathbb{N})$, one can compute a map

$$ \forall f : (V \to \mathbb{N}) \to \mathbb{N}, \mathbb{D}_{\text{fun}} f \to \mathbb{D}_{\text{rel}}(\lambda\nu.\forall y, \neg Q \//_F f\nu \succ y). $$

*Proof.* The map $\forall f, \mathbb{D}_{\text{fun}} f \to \mathbb{D}_{\text{rel}}(\lambda\nu.\forall y, \neg Q \//_F f\nu \succ y)$ is built by induction on $Q$. If $Q = []$, then we show $(\forall y, \neg [] \//_F f\nu \succ y) \leftrightarrow \text{True}$, and thus $\mathbb{D}_{\text{rel}}(\lambda\nu.\forall y, \neg Q \//_F f\nu \succ y)$ by map 4 of Proposition 2.2 followed by **dio auto**. If $Q = p/q :: Q'$, then we show the equivalence

$$ (\forall y, \neg Q \//_F f\nu \succ y) \leftrightarrow (q \nmid (p \cdot (f\nu))) \wedge (\forall y, \neg Q' \//_F f\nu \succ y) $$

and we get $\mathbb{D}_{\text{rel}}(\lambda\nu.\forall y, \neg Q \//_F f\nu \succ y)$ by map 4 of Proposition 2.2 followed by **dio auto**, the induction hypothesis being used as a hint again. □

We can now deduce a core result of the paper which states that FRACTRAN programs have Diophantine termination predicates.

**Theorem 7.3.** If $Q : \mathbb{L}(\mathbb{N} \times \mathbb{N})$ is a FRACTRAN program then one can compute a map

$$ \forall f : (V \to \mathbb{N}) \to \mathbb{N}, \mathbb{D}_{\text{fun}} f \to \mathbb{D}_{\text{rel}}(\lambda\nu.Q \//_F f\nu \downarrow). $$

*Proof.* By definition we have $Q \//_F f\nu \downarrow \leftrightarrow \exists x (Q \//_F f\nu \succ * x \wedge \forall y, \neg Q \//_F x \succ y)$ and hence we obtain the claim using Theorem 2.10 together with Lemma 7.1 and Lemma 7.2. □

We conclude with the undecidability of Hilbert's tenth problem by a reduction chain starting from the Halting problem for single tape Turing machines:

**Theorem 7.4 (Hilbert's tenth problem).** We have the following reduction chain

$$ \text{Halt} \preceq \text{PCP} \preceq \text{MM} \preceq \text{FRACTRAN} \preceq \text{DIO\_FORM} \preceq \text{DIO\_ELEM} \preceq \text{DIO\_SINGLE} \preceq \text{H10} $$

and as a consequence, H10 is undecidable.

*Proof.* The proof combines previous results as Theorems 6.3 and 7.3 and Corollary 4.4. □

## 8. THE DAVIS-PUTNAM-ROBINSON-MATIYASEVICH THEOREM

We give a proof of an instance of the DPRM theorem stating that recursively enumerable predicates are Diophantine<fnref n="19" /> Here we assume that the informal notion of “recursive enumerability” (justified by Church’s thesis) can be characterised by Minsky machines recognisability as defined in Section 6.1.

**Proposition 8.1.** The Gödel encoding is Diophantine, i.e. we have $\mathbb{D}_{\text{fun}}(\lambda\nu.q_1^{\nu x_0} \cdots q_n^{\nu x_{n-1}})$.

*Proof.* By induction on $n : \mathbb{N}$ using Proposition 2.4 and Theorem 2.7. Notice that the $q_i$'s are hard-coded in the Diophantine representation, which means we do not need to encode the algorithm that actually computes them, which would otherwise be very painful. □

<footnote n="19">¹⁹By instance, we mean that the DPRM is an open theorem bound to be extended for any newly proposed Turing complete model of computation.</footnote>
**Lemma 8.2.** For any FRACTRAN program Q we have $\mathbb{D}_{\text{rel}}(\lambda\nu.Q \//_{F} \mathfrak{p}_1\mathfrak{q}_1^{\nu x_0} \cdots \mathfrak{q}_n^{\nu x_{n-1}} \downarrow)$.

*Proof.* By Theorem 7.3, we only have to show $\mathbb{D}_{\text{fun}} f$ for $f \nu := \mathfrak{p}_1\mathfrak{q}_1^{\nu x_0} \cdots \mathfrak{q}_n^{\nu x_{n-1}}$. This follows from Propositions 2.4 and 8.1. $\square$

To simplify the notation $[[p]]_{\vec{v}}^{\vec{w}}$ below, we abusively identify the vector $\vec{v} : \mathbb{N}^n$ (resp. $\vec{w} : \mathbb{N}^m$) with the valuation $\lambda(i : \mathbb{F}_n).\vec{v}_i$ (resp. $\lambda(j : \mathbb{F}_m).\vec{w}_j$) that accesses the components of the vector $\vec{v}$ (resp. $\vec{w}$).

**Theorem 8.3 (DPRM).** Any MM-recognisable relation $R : \mathbb{N}^n \to \mathbb{P}$ is Diophantine: one can compute a single Diophantine equation $p \doteq q : \mathbb{D}_{\text{single}}(\mathbb{F}_m, \mathbb{F}_n)$ with $n$ parameters and $m$ variables s.t. $\forall \vec{v} : \mathbb{N}^n, R \vec{v} \leftrightarrow \exists \vec{w} : \mathbb{N}^m, [[p]]_{\vec{v}}^{\vec{w}} = [[q]]_{\vec{v}}^{\vec{w}}$.

*Proof.* By definition, $R : \mathbb{N}^n \to \mathbb{P}$ is recognised by some Minsky machine $P$ with $(n+m)$ registers, i.e. $R \vec{v} \leftrightarrow (1, P) \//_M (1, \vec{v} \equiv \vec{0}) \downarrow$. By Theorem 6.3, we compute a FRACTRAN program $Q$ s.t. $(1, P) \//_M (1, [v_1; \dots; v_n; w_1; \dots; w_m]) \downarrow \leftrightarrow Q \//_F \mathfrak{p}_1\mathfrak{q}_1^{v_1} \cdots \mathfrak{q}_n^{v_n} \mathfrak{q}_{n+1}^{w_1} \cdots \mathfrak{q}_{n+m}^{w_m} \downarrow$.

Hence we deduce $R [v_1; \dots; v_n] \leftrightarrow Q \//_F \mathfrak{p}_1\mathfrak{q}_1^{v_1} \cdots \mathfrak{q}_n^{v_n} \downarrow$. As a consequence, the relation $\lambda\nu.R [\nu x_0; \dots; \nu x_{n-1}]$ is Diophantine by Lemma 8.2. By Corollary 4.4, there is a Diophantine equation $p \doteq q : \mathbb{D}_{\text{single}}(\mathbb{N}, \mathbb{V})$ such that $R [\nu x_0; \dots; \nu x_{n-1}] \leftrightarrow \exists \varphi, [[p]]_{\nu}^{\varphi} = [[q]]_{\nu}^{\varphi}$. Notice that the value $\nu x_i$ of any parameter of $p \doteq q$ greater than $x_n$ does not influence its solvability.

Now let $m$ be an upper bound of the number of (existentially quantified) variables in $p \doteq q$. We injectively map those variables in $\mathbb{F}_m$ and we project the parameters of $p \doteq q$ onto $\mathbb{F}_n$ by replacing every parameter greater than $x_n$ with the 0 constant. We get a Diophantine equation $p' \doteq q' : \mathbb{D}_{\text{single}}(\mathbb{F}_m, \mathbb{F}_n)$ of which the solvability at $\vec{v}$ is equivalent to $R \vec{v}$. $\square$

## 9. HILBERT'S TENTH PROBLEM OVER INTEGERS

In our formalisation, polynomials are defined over natural numbers, that is both constants and solutions come from $\mathbb{N}$. The standard way to extend the undecidability of H10 to a formalisation based on integers is via Lagrange's theorem, stating that an integer is positive if and only if it is the sum of four squares.

Similar to Definition 4.2, we define polynomials over integers:

**Definition 9.1.** The type of Diophantine polynomials $\mathbb{D}_{\text{poly}}^{\mathbb{Z}}(\mathbb{U}, \mathbb{V})$ over $\mathbb{Z}$ is defined by:

$$p, q : \mathbb{D}_{\text{poly}}^{\mathbb{Z}}(\mathbb{U}, \mathbb{V}) ::= u : \mathbb{U} \mid x : \mathbb{V} \mid z : \mathbb{Z} \mid p \dot{+} q \mid p \dot{\times} q.$$

The interpretation of a polynomial $p : \mathbb{D}_{\text{poly}}^{\mathbb{Z}}(\mathbb{U}, \mathbb{V})$ in $\mathbb{Z}$ given $\varphi : \mathbb{U} \to \mathbb{Z}, \nu : \mathbb{V} \to \mathbb{Z}$ and denoted $[[p]]_{\nu}^{\varphi}$, is defined in the obvious way. Again, if $p : \mathbb{D}_{\text{poly}}^{\mathbb{Z}}(\mathbb{F}_m, \mathbb{F}_n)$ we abusively write $[[p]]_{\vec{v}}^{\vec{w}}$ when $\vec{v} : \mathbb{N}^n$ and $\vec{w} : \mathbb{N}^m$. We can then define

$$\text{H10}_{\mathbb{Z}}(n, p : \mathbb{D}_{\text{poly}}^{\mathbb{Z}}(\mathbb{F}_n, \mathbb{F}_0)) := \exists \vec{w} : \mathbb{N}^n, [[p]]_{\downarrow}^{\vec{w}} = 0$$

that is “does the polynomial equation $p(x_1, \dots, x_n) = 0$ in (at most) $n$ variables have a solution in $\mathbb{Z}$.” We first outline a proof of Lagrange’s theorem and then reduce H10 to H10$_{\mathbb{Z}}$.
**9.1. Lagrange's theorem.** The proof we have implemented roughly follows the "classical proof" in Wikipedia's account of the theorem. Their use of the "classical" qualifier should be understood as typical/standard, and certainly not as opposed to constructive/intuitionistic. The below proof perfectly fits in our constructive setting.

**Proposition 9.2** (Euler, 1748). Let us assume the two equations $n = a_1^2 + b_1^2 + c_1^2 + d_1^2$ and $m = a_2^2 + b_2^2 + c_2^2 + d_2^2$ hold in $\mathbb{Z}$. Let us define the four relative integers:

$$a := a_1a_2 + b_1b_2 + c_1c_2 + d_1d_2$$

$$b := a_1b_2 - b_1a_2 + d_1c_2 - c_1d_2$$

$$c := a_1c_2 - c_1a_2 + b_1d_2 - d_1b_2$$

$$d := a_1d_2 - d_1a_2 + c_1b_2 - b_1c_2.$$

Then the identity $nm = a^2 + b^2 + c^2 + d^2$ holds.

*Proof.* This holds in any commutative ring and the proof just calls the **ring tactic**. □

Because of Proposition 9.2, "being the sum of four squares" is a multiplicative property. Hence to show that it holds for every natural number, it is enough to establish it for primes.

**Theorem 9.3** (Prime induction). Let $P: \mathbb{N} \to \mathbb{P}$ be a predicate. To establish $\forall n: \mathbb{N}, Pn$, it is enough to prove these four induction steps:

$$P0 \quad P1 \quad \forall ab: \mathbb{N}, Pa \to Pb \to P(ab) \quad \forall p: \mathbb{N}, \text{prime } p \to Pp.$$

*Proof.* This is one possible form of the fundamental theorem of arithmetic. For the proof, first show by strong induction on $n: \mathbb{N}$ that one can discriminate whether $n < 2$ or compute a prime factor of $n$, including the possibility that $n$ itself is prime. Then, to prove the prime induction principle, proceed by strong induction again. □

**Lemma 9.4.** If $p: \mathbb{N}$ is prime then $\exists n ab: \mathbb{N}, np = 1 + a^2 + b^2 \wedge 0 < n < p$.

*Proof.* Let us first rule out the case $p = 2$ which has an obvious solution. So let us write $p = 2m + 1$ because all the other primes are odd. In the field $\mathbb{Z}/p\mathbb{Z}$, let us study the modular equation $a^2 \equiv -(1+b^2)[p]$. Because $\mathbb{Z}/p\mathbb{Z}$ has no zero divisor, the map $a \mapsto a^2$ from $[0, m] \to \mathbb{Z}/p\mathbb{Z}$ is injective. As a consequence, so is the map $b \mapsto -(1+b^2)$ from $[0, m] \to \mathbb{Z}/p\mathbb{Z}$. Hence none of the two lists $[a^2: \mathbb{Z}/p\mathbb{Z} | a \in [0, m]]$ and $[-(1+b^2): \mathbb{Z}/p\mathbb{Z} | b \in [0, m]]$ contain a duplicate. Since they have a combined length of $2(m+1) > 2m+1$ and their concatenation is contained in a list of length $p = 2m+1$ enumerating $\mathbb{Z}/p\mathbb{Z}$, by the pigeon hole principle, they must intersect and this gives $a, b \in [0, m]$ such that $a^2 \equiv -(1+b^2)[p]$. This in turn gives $n: \mathbb{N}$ such that $np = 1 + a^2 + b^2$. Given that $2a < p$ and $2b < p$, we deduce $0 < n < p$. □

Lagrange's theorem gives a Diophantine characterisation of those relative integers which are positive as the sum of four squares.

**Theorem 9.5** (Lagrange, 1770). For any relative integer $z: \mathbb{Z}$, $z$ is positive if and only if there exists $a, b, c, d: \mathbb{Z}$ such that $z = a^2 + b^2 + c^2 + d^2$.

*Proof.* It is enough to show that any natural number $n: \mathbb{N}$ is the sum of four squares, i.e. there exists $a, b, c, d: \mathbb{Z}$ such that $n = a^2 + b^2 + c^2 + d^2$. By Euler's four-square identity (see Proposition 9.2) and the principle of prime induction (see Theorem 9.3), we only need to show the property for prime numbers.

We fix a prime number $p$ and define the predicate "$mp$ is the sum of four squares" as

$$P(m: \mathbb{N}) := \exists a b c d: \mathbb{Z}, mp = a^2 + b^2 + c^2 + d^2.$$
We want to show that $P1$ holds. By Lemma 9.4, we know that $Pn$ holds for some $0 < n < p$. We are going to decrease this value of $n$ until it reaches 1, i.e. to establish $P1$, it is enough to give a proof of

$$ \forall m : \mathbb{N}, 1 < m < p \rightarrow Pm \rightarrow \exists r : \mathbb{N}, 1 \le r < m \wedge Pr $$

and then finish the argument by strong induction (implementing “infinite descent” here).

So let us assume $1 < m < p$ and $mp = x_1^2+x_2^2+x_3^2+x_4^2$. For each $i = 1, \dots, 4$ we compute a “small” representative of $x_i$ in $\mathbb{Z}/m\mathbb{Z}$, i.e. $y_i : \mathbb{Z}$ such that $x_i \equiv y_i [m]$ and $4y_i^2 \le m^2$. Then we have $mp = x_1^2+x_2^2+x_3^2+x_4^2 \equiv y_1^2+y_2^2+y_3^2+y_4^2 [m]$ hence $y_1^2+y_2^2+y_3^2+y_4^2 \equiv 0 [m]$ and we get $r : \mathbb{N}$ such that $y_1^2+y_2^2+y_3^2+y_4^2 = rm$. Because $y_1^2+y_2^2+y_3^2+y_4^2 \le m^2$, we deduce $r \le m$. Now we rule out the cases $r = 0$ and $r = m$:

* if $r = 0$ then $y_1 = \dots = y_4 = 0$ and thus $x_i^2 \equiv 0 [m^2]$ for $i = 1, \dots, 4$. As a consequence, $mp = x_1^2+x_2^2+x_3^2+x_4^2 \equiv 0 [m^2]$ and thus $m$ divides $p$, contradicting the primality of $p$;

* if $r = m$ then $y_1^2+y_2^2+y_3^2+y_4^2 = m^2$ and since $4y_i^2 \le m^2$ holds for $i = 1, \dots, 4$, we deduce $4y_1^2 = \dots = 4y_4^2 = m^2$ and thus $m = 2q$ and $y_i = \pm q$ for $i = 1, \dots, 4$. As a consequence $x_i^2 \equiv q^2 [m^2]$ and thus $x_1^2+x_2^2+x_3^2+x_4^2 \equiv 4q^2 [m^2]$ hence $mp \equiv 0 [m^2]$, in contradiction with the primality of $p$ again.

So we have $1 \le r < m$, $rm = y_1^2+y_2^2+y_3^2+y_4^2$ and $x_i \equiv y_i [m]$ for $i = 1, \dots, 4$. We also have $mp = x_1^2+x_2^2+x_3^2+x_4^2$ and using Euler's four squares identity we get $a, b, c, d : \mathbb{Z}$ such that $(mp)(rm) = a^2+b^2+c^2+d^2$. Using $x_i \equiv y_i [m]$ and the values of $a, b, c$ and $d$ as defined in Proposition 9.2, we show that $a \equiv 0 [m]$, $b \equiv 0 [m]$, $c \equiv 0 [m]$ and $d \equiv 0 [m]$, i.e. $m$ divides $a, b, c$ and $d$. Hence $rp = (a/m)^2+(b/m)^2+(c/m)^2+(d/m)^2$ which establishes $Pr$ for $1 \le r < m$, as required. $\square$

**9.2. H10 reduces to H10<sub>Z</sub>.** Let $f : \mathbb{D}_{\text{poly}}(\mathbb{F}_n, \mathbb{F}_0) \rightarrow \mathbb{D}_{\text{poly}}^{\mathbb{Z}}(\mathbb{F}_{4n}, \mathbb{F}_0)$ be the function replacing every variable $i : \mathbb{F}_n$ in a polynomial by the polynomial expression

$$ (4i) \dot{\times} (4i) \dot{\times} (4i+1) \dot{\times} (4i+1) \dot{\times} (4i+2) \dot{\times} (4i+2) \dot{\times} (4i+3) \dot{\times} (4i+3). $$

The function $f$ can now be used to define the reduction from $H10$ to $H10_Z$:

**Corollary 9.6.** $H10 \le H10_Z$.

*Proof.* Given polynomials $p, q : \mathbb{D}_{\text{poly}}(\mathbb{F}_n, \mathbb{F}_0)$ the reduction function returns the polynomial $f(p) + (-1) \dot{\times} f(q) : \mathbb{D}_{\text{poly}}^{\mathbb{Z}}(\mathbb{F}_{4n}, \mathbb{F}_0)$.

If $H10(n, p, q)$ holds, i.e. if $[[p]]_v^w = [[q]]_v^w$ for $\vec{w} = [w_1; \dots; w_n]$ we know via Lagrange's Theorem 9.5 that there exist $a_1, b_1, c_1, d_1, \dots, a_n, b_n, c_n, d_n$ s.t. $w_i = a_i^2+b_i^2+c_i^2+d_i^2$ and thus that $[[f(p)]]_v^w = [[p]]_v^w$ for $\vec{v} : \mathbb{Z}^{4n} := [a_1; b_1; c_1; d_1; \dots; a_n; b_n; c_n; d_n]$ as well as $[[f(q)]]_v^w = [[q]]_v^w$. Thus, $[[f(p) \dot{\times} (-1) \dot{\times} f(q)]]_v^w = 0$ and $H10_Z(4n, f(p) \dot{\times} (-1) \dot{\times} f(q))$.

For the other direction, let $\vec{v} : \mathbb{Z}^{4n}$ be given s.t. $[[f(p) \dot{\times} (-1) \dot{\times} f(q)]]_v^w = 0$. We thus know that $[[f(p)]]_v^w = [[f(q)]]_v^w$. Because $\vec{v} : \mathbb{Z}^{4n}$, it has the form $\vec{v} = [a_1; b_1; c_1; d_1; \dots; a_n; b_n; c_n; d_n]$. Then all elements of $\vec{w} := [a_1^2+b_1^2+c_1^2+d_1^2; \dots; a_n^2+b_n^2+c_n^2+d_n^2]$ are natural numbers and the equations $[[f(p)]]_v^w = [[p]]_v^w$ and $[[f(q)]]_v^w = [[q]]_v^w$ hold. Thus $[[p]]_v^w = [[q]]_v^w$, i.e. $H10(n, p, q)$. $\square$
<figure><img src="image_4.png" /><figcaption>Figure 2: Relational semantics for μ-recursive algorithms.</figcaption></figure>

## 10. μ-RECURSIVE ALGORITHMS

In order to show that MM, FRACTRAN, and H10 are in the same many-one reduction class (i.e. interreducible via many-one reductions), we introduce μ-recursive algorithms as intermediate layer. Programming in this well-known model of computation resembles functional programming (in a first-order language) and we will use it first for a reduction H10 $\preceq$ μ-rec. Afterwards, we explain a compiler to Minsky machines, yielding a reduction from μ-rec to MM. The next section will then connect μ-recursive algorithms to the weak call-by-value λ-calculus.

**10.1. μ-Recursive Recognisability.** We define a type $\mathcal{A}_k$ of μ-recursive algorithms representing μ-recursive partial functions in $N^k \to N$.

<figure><img src="image_5.png" /></figure>

We represent these μ-recursive partial functions using the standard relational semantics $[f] : N^k \to N \to P$ of the μ-recursive algorithm $f : \mathcal{A}_k$, formalised in Figure 2 as a fixpoint definition. Intuitively, $cst_n$ represents the constant $n$ of arity 0, zero is the constant 0-function of arity 1, succ the successor function of arity 1, projection $pr_{j_p}$ returns the $p$-th argument, $comp f g$ where $g$ is a $k$-vector of functions of arity $i$ first applies each element of $g$ to the $i$ inputs and then $f$ to the resulting $k$ numbers. $rec f g$ is performing primitive recursion on the first argument. If the argument is 0, $f$ is used. If the argument is $1+n$, $g$ is applied to $n$, the recursive call and the rest of the arguments. Finally, minimisation $\mu f$ performs unbounded search returning the smallest number $x$ s.t. $f$ on $x$ returns 0 and $f$ terminates on a non-zero value for every $y < x$.

Following [Lar17], in Fig. 3 we also formalise a cost aware big-step evaluation predicate $[f; v] -[c] x$, where $f : \mathcal{A}_k, v : N^k, x : N$, and where $c : N$ denotes the cost of a computation, the cost being tailored toward the naive step-indexed evaluator to be defined in the next section.<fnref n="20" /> We can directly relate the relational semantics with the big-step semantics:

**Lemma 10.1.** $\forall k (f : \mathcal{A}_k) (\vec{v} : N^k) (x : N), [[f]] \vec{v} x \leftrightarrow \exists c, [f; \vec{v}] -[c] x$.

<footnote n="20">²⁰In the Coq code, the predicate has one additional auxiliary argument to ease the correctness proof for the step-indexed interpreter which is only internal to the proof and thus omitted here.</footnote>
<figure><img src="image_6.png" /><figcaption>Figure 3: Cost aware big-step semantics for μ-recursive algorithms.</figcaption></figure>

With these characterizations of the semantics of μ-recursive algorithms, we can define μ-recursive halting as follows:

$$ \mu\text{-rec}(k : \mathbb{N}, f : \mathcal{A}_k, \vec{v} : \mathbb{N}^k) := \exists x, \llbracket f \rrbracket \vec{v} x $$

A function $f : \mathcal{A}_k$ is called *total* if it terminates on any input, i.e. $\forall \vec{v} : \mathbb{N}^k \exists x, \llbracket f \rrbracket \vec{v} x$.

**Definition 10.2.** We call a relation $k$-ary relation $R$ over natural numbers $\mu$-recursively recognisable if there is a $\mu$-recursive algorithm $f : \mathcal{A}_k$ s.t. for any $\vec{v} : \mathbb{N}^k$ the equivalence $R \vec{v} \leftrightarrow \exists x, \llbracket f \rrbracket \vec{v} x$ holds.

**10.2. Compiling μ-recursive algorithms into Minsky machines.** We describe how to compile a μ-recursive algorithm $f : \mathcal{A}_k$ into a Minsky machine. To avoid solving complicated inequality constraints over the bounded numbers in $\mathbb{F}_m$,<fnref n="21" /> we here work with Minsky machines where registers are indexed with $\mathbb{N}$ instead of $\mathbb{F}_m$. Hence the number of registers is not bounded a priori but of course, for a given Minsky machine, the number of register that actually occur in the code is bounded. And indeed, in the end of the process in Theorem 10.5, we project to Minsky machines with registers in $\mathbb{F}_m$ (for some $m$) after the μ-recursive algorithm has been fully compiled. We will not enter into the details on how particular μ-recursive operators are implemented but instead focus on the global invariant of the compiler.

For a μ-recursive algorithm $f : \mathcal{A}_k$ with $k$ inputs and one output, we produce a Minsky machine $P : I_{\mathbb{N}}$ with the following structure for registers, described by the four extra parameters $i, p, o$ and $m$, all of type $\mathbb{N}$:

* $i$ : is the PC-index of the first instruction of $P$;

* $p$ : the $k$ inputs for $f$ are to be read in the $k$ registers $\{p, \dots, p+k-1\}$ of $P$;

* $o$ : the output of $f$ is to be written in the register $o$ of $P$;

<footnote n="21">²¹Unlike $\mathbb{N}$ or $\mathbb{Z}$, the inductive type $\mathbb{F}_m$ representing the numbers bounded by $m$ is not equipped with powerful linear constraints solving tactics such as *lia* or *omega*. Working with bounded numbers here would require from us to manipulate embeddings of e.g. the disjoint sum $\mathbb{F}_i + \mathbb{F}_j$ into $\mathbb{F}_k$ (for some large enough $k$), which would also impact types that depend on $i, j$ and $k$. For those who experienced it, this should instantly trigger the bad memory of *setoid hell*. Using the above mentioned tools, it is much simpler to work with intervals of $\mathbb{N}$ instead.</footnote>
m : P can use spare registers above m and assume their initial value is 0.

Moreover, we require that all registers except the output register o are returned to their initial value when the computation is terminated. In particular, spare registers return to their initial value 0.

**Definition 10.3.** We say that P properly compiles f under the constraints i, p, o, m and we write ra_compiled k (f : A_k) (i p o m : N) (P : l_N) if for any v : N^k and p : N -> N such that (1) ∀u ≥ m, ρ_u = 0 and (2) ∀q : F_k, ρ_q+p = v_q, we have

**soundness:** ∀x : N, [[f]] v x → (i, P) //_M (i, p) >_* (|P| + i, p{o ← x});

**completeness:** (i, P) //_M (i, p) ↓ → ∃x, [[f]] v x.

In these N-indexed register machines, the state of a machine cannot be described by a finite vector v : N^k but is instead represented by an environment p : N → N mapping register indices (in N) to the values they contain (in N also, but with a different meaning). The two conditions (1) and (2) state that p is null above register m and p contains v starting at register p. And under these conditions, (i, P) should simulate a terminating computation of f on v outputting x (soundness), while conversely, whenever (i, P) terminates starting from (i, p), this entails that f terminates on v (completeness).

We can now construct a certified compiler from μ-recursive algorithms to N-indexed Minsky machines. Literally, Theorem 10.4 below states that whenever the output register o does not belong to the input registers {p, ..., p + k - 1} or the spare registers {m, m + 1, ...} and that spare registers have indices above input registers, then one can properly compile f at PC value i.

**Theorem 10.4.** Given a μ-recursive algorithm f : A_k, we can build a term:

$$ \forall ipom, o < m \to \neg(p \le o < k + p) \to k + p \le m \to \sum P, \text{ra\_compiled } k \text{ f ipom } P. $$

The proof proceeds by structural induction on f : A_k. The details are quite involved and not exposed in here: we just wanted to state the above invariant, which is closed under each μ-recursive constructor. The sub-Minsky machines are built and composed using the compositional techniques already presented in [FL19]. Notice that the statement of ra_compiled takes the output value x in [[f]] v x into account while in the below Theorem 10.5, we only care about termination. Because the termination of e.g. the μ-recursive composition comp f g depends not just on the termination of g (for all p : F_k) but on the actual output values w_p of [[g_p]] v w_p, it is necessary to be more precise in the stated invariant used in the inductive construction of the compiler.

**Theorem 10.5.** Given a μ-recursive algorithm f : A_k, one can compute n : N and a list of Minsky machine instructions P : l_{F_{k+1+n}}, such that for any v : N^k,

$$ (\exists x, \llbracket f \rrbracket \vec{v} x) \leftrightarrow (1, P) \mathrel{/\mathrel{\parallel}}_M (1, \vec{v} + \vec{0}) \downarrow. $$

*Proof.* By Theorem 10.4, we compile f into P' : l_N under the constraints i := 1, p := 0, o := k and m := 1 + k. Then we compute n such that k + 1 + n is a strict upper bound of all the register indices that occur in P'. Now we can map P' : l_N to P : l_{F_{k+1+n}} while preserving its semantics. Notice that contrary to Theorem 10.4, the resulting output value of f or (1, P) can now be disregarded since we only care about termination. □

**Corollary 10.6.** μ-rec ⊆ MM.
**10.3. Diophantine relations are μ-recursively recognisable.** We encode n-ary Diophantine relations using μ-recursive algorithms. Following the DPRM Theorem 8.3, we consider a n-ary relation $R : \mathbb{N}^n \to \mathbb{P}$ to be Diophantine if there is a single Diophantine equation $p \doteq q : \mathbb{D}_{\text{single}}(\mathbb{F}_m, \mathbb{F}_n)$ with n parameters and m (existential) variables such that

$$ \forall \vec{v} : \mathbb{N}^n, R \vec{v} \leftrightarrow \exists \vec{w} : \mathbb{N}^m, [[p]]_{\vec{v}}^{\vec{w}} = [[q]]_{\vec{v}}^{\vec{w}}. $$

Again we abusively confuse vectors with maps using the equivalence $X^n \simeq \mathbb{F}_n \to X$. We will prove that every Diophantine relation in this sense is μ-recursively recognisable (see Definition 10.2). The proof idea is relatively straightforward.

First, we implement $\text{eval} : \mathbb{D}_{\text{poly}}(\mathbb{F}_m, \mathbb{F}_n) \to \mathcal{A}_{m+n}$ evaluating any polynomial:

**Lemma 10.7.** Given any Diophantine polynomial $p : \mathbb{D}_{\text{poly}}(\mathbb{F}_m, \mathbb{F}_n)$, one can compute a μ-recursive algorithm $\text{eval}_p : \mathcal{A}_{m+n}$ s.t. $[[\text{eval}_p]](\vec{w} + \vec{v}) [[p]]_{\vec{v}}^{\vec{w}}$.

*Proof.* The implementation relies on implementations of addition and multiplication on natural numbers, which are relatively straightforward for μ-recursive algorithms. □

Secondly, we implement a bijection between $\mathbb{N}$ and $\mathbb{N}^m$ for any $m$:

**Lemma 10.8.** There are functions $\text{pr}_m : \mathbb{N} \to \mathbb{N}^m$ and $\text{inj}_m : \mathbb{N}^m \to \mathbb{N}$ s.t. $\text{pr}_m(\text{inj}_m \vec{v}) = \vec{v}$. Furthermore, given $i : \mathbb{F}_m$, one can compute a μ-recursive algorithm $\text{project}_i : \mathcal{A}_1$ s.t. $[[\text{project}_i]](x ::[]) (\text{pr}_m x)_i$, i.e. all the $m$ components of $\text{pr}_m$ are μ-recursive.

Thirdly, given an equation $p \doteq q$, we compute a test algorithm $\text{test}_{p,q} : \mathcal{A}_{1+n}$, always terminating and returning 0 iff $p \doteq q$ is satisfied when decoding the first argument as a vector $\vec{w}$ using $\text{pr}_m$:

**Lemma 10.9.** Given any $p \doteq q : \mathbb{D}_{\text{single}}(\mathbb{F}_m, \mathbb{F}_n)$, one can compute a μ-recursive algorithm $\text{test}_{p,q} : \mathcal{A}_{1+n}$ s.t. $[[\text{test}_{p,q}]](x :: \vec{v}) 0 \leftrightarrow [[p]]_{\vec{v}}^{\vec{w}} = [[q]]_{\vec{v}}^{\vec{w}}$ where $\vec{w} := \text{pr}_m(x)$. Furthermore, $\text{test}_{p,q}$ is primitive recursive hence total.

To finish, given a vector $\vec{v} : \mathbb{N}^n$ describing the values of parameters, we use minimisation $\mu$ to search for a number $x$ such that $\vec{w} := \text{pr}_m(x)$ is a solution for $p \doteq q$:

**Theorem 10.10.** Given $p \doteq q : \mathbb{D}_{\text{single}}(\mathbb{F}_m, \mathbb{F}_n)$, one can compute a μ-recursive algorithm $\text{find}_{p,q} : \mathcal{A}_n$ s.t. for any $\vec{v} : \mathbb{N}^n$, we have the equivalence

$$ (\exists x : \mathbb{N}, [[\text{find}_{p,q}]] \vec{v} x) \leftrightarrow \exists \vec{w} : \mathbb{N}^m, [[p]]_{\vec{v}}^{\vec{w}} = [[q]]_{\vec{v}}^{\vec{w}}. $$

Notice that while $\text{find}_{p,q}$ cannot be primitive recursive (it does not terminate when the equation $p \doteq q$ is unsolvable), it is however implemented as a single unbounded minimization applied to an otherwise primitive recursive algorithm.

**Corollary 10.11.** $\text{H10} \preceq \mu-\text{rec}$.

*Proof.* Direct application of Theorem 10.10 with no parameters, i.e. $n := 0$. □
## 11. THE WEAK CALL-BY-VALUE λ-CALCULUS L

Until now, we have shown that the problems MM, FRACTRAN, H10, and μ-rec are all interreducible w.r.t. many-one reductions. We contribute our proofs to the Coq library of undecidability proofs [FLWD+20], which contains several other well-known reduction proofs. Amongst them is a chain of reductions published in related work and discussed in Section 12, establishing that the halting problem of Turing machines reduces to MM and thus to all other problems we consider.

It is also possible to prove that all considered problems are in fact interreducible to the halting problem of Turing machines. We demonstrate one technique to do so, based on the weak call-by-value λ-calculus L, which is already shown interreducible with the halting problem of Turing machines in [FKSW21].

We briefly introduce L in this section and then reduce μ-rec to halting in L. Syntactically, L is just the untyped λ-calculus, De Bruijn style:

$$s, t, u : L ::= n \mid st \mid \lambda s \quad \text{where } n : \mathbb{N}$$

We define a weak call-by-value evaluation predicate $s \triangleright t$, coinciding with other definitions of reduction for closed terms, where $s[t/0]$ means “replace De Bruijn variable 0 in $s$ by $t$”.

$$\frac{\lambda s \triangleright \lambda s}{st \triangleright u}$$

The halting problem for L can then be defined as $\text{WCBV}(s : L) := \exists t, s \triangleright t$.

It is possible to encode natural numbers, vectors over natural numbers, and other data types into L using Scott's encoding. We will denote with $\bar{\cdot}$ such encoding functions. For further details on encodings and L we refer to [FS17], since they do not actually matter to understand the reduction.

In general, a predicate reduces to the L-halting problem WCBV if it is L-recognisable:

**Definition 11.1.** We say that a function $f : X_1 \to \dots \to X_n \to Y$ is L-computable if there exists a closed term $t$ s.t. $\forall x_1 : X_1, \dots, \forall x_n : X_n. t \bar{x}_1 \dots \bar{x}_n \triangleright f x_1 \dots x_n$.

Given a type X encodable in L and a predicate P over X we say that P is L-recognisable if there is an L-computable function $f : X \to \mathbb{N} \to \mathbb{B}$ s.t. $P x \leftrightarrow \exists n, f x n = \text{true}$.

**Theorem 11.2.** Let X be encodable in L. For a predicate P over X we have that $P \preceq \text{WCBV}$ if and only if P is L-recognisable.

To instantiate the theorem to μ-rec one would like to give a step-indexed evaluation function for μ-recursive algorithms, i.e. for $f : A_k$ a function $[[f]] : \mathbb{N} \to \mathbb{N}^k \to \mathbb{O} \mathbb{N}$ s.t. $[f; \vec{v}] - [c] x \leftrightarrow [[f]]c \vec{v} = [x]$ and then show that it is computable in L. The second step, i.e. proving that this function is L-computable, is in principle fully automatic using the certifying extraction framework from Forster and Kunze [FK19], which is implemented using tools from the MetaCoq projet [SAB+20]. The framework supports extracting non-dependent, non-mutual, non-nested, simply-typed Coq functions to L automatically, and also generates a proof of correctness of the extract. However, recall that the type $A_k$ of μ-recursive functions is both heavily dependent and used nested constructions: It has a type index in $k$, makes use of finite types $\mathbb{F}_k$, and has a nested use of the dependent vector type, i.e. contains subterms of type $(A_k)^i$. Thus, a direct step-indexed interpreter will not be extractable. Instead, we use a general technique and implement a step-indexed interpreter working on the syntactic skeleton of μ-recursive functions.
A syntactic skeleton for a type *I* mirrors the constructors of *I*, but without any dependent types. If *I* furthermore has nested applications of types $N_1, \dots, N_n$, the syntactic skeleton also contains the constructors of $N_1, \dots, N_n$. That is, if *I* has *i* constructors and $N_1, \dots, N_n$ have $m_1, \dots, m_n$ constructors respectively, the syntactic skeleton of *I* has $i + m_1 + \dots + m_n$ constructors. The syntactic skeleton $\mathcal{A}'$ of the type $\mathcal{A}$ is defined as follows:

$$
\begin{array}{ccc}
\begin{array}{c} n : \mathbb{N} \\ \overline{\text{cst}_n : \mathcal{A}'} \end{array} & \begin{array}{c} \overline{\text{zero} : \mathcal{A}'} \\ \overline{\text{succ} : \mathcal{A}'} \end{array} & \begin{array}{c} j : \mathbb{N} \\ \overline{\text{prj}_j : \mathcal{A}'} \end{array} \\
\begin{array}{c} \overline{\text{comp f g} : \mathcal{A}'} \end{array} & \begin{array}{c} \overline{\text{rec f g} : \mathcal{A}'} \end{array} & \begin{array}{c} \overline{\text{cons f g} : \mathcal{A}'} \end{array}
\end{array}
\qquad
\begin{array}{c} \overline{\text{nil} : \mathcal{A}'} \\ \overline{\mu f : \mathcal{A}'} \end{array}
$$

Note that we re-use constructor names. In the syntactic skeleton we have $j : \mathbb{N}$ instead of $j : \mathbb{F}_k$ for $\text{prj}_j$, $g : \mathcal{A}'$ instead of $G : (\mathcal{A}_k)^i$ for $\text{comp f g}$, and the constructors $\text{cons}$ and $\text{nil}$ as the constructors of the vector type are added.

It is straightforward to implement mutually recursive functions $\text{erase} : \mathcal{A}_k \to \mathcal{A}'$ and $\text{erase}': (\mathcal{A}_k)^i \to \mathcal{A}'$ which are essentially the identity. We call a skeleton $s : \mathcal{A}' \text{ valid}$ if it corresponds to a $\mu$-recursive function or a vector, i.e. if $\exists k : \mathbb{N}$, $(\exists f : \mathcal{A}_k, s = \text{erase } f) \lor (\exists i : \mathbb{N}, \exists l : (\mathcal{A}_k)^i, s = \text{erase}' l)$.

For $\mathcal{A}'$ we can now define a step-indexed evaluation function

$$ \llbracket \cdot \rrbracket : \mathcal{A}' \to \mathbb{N} \to \mathbb{N} \to \mathbb{L}\mathbb{N} \to \mathcal{O}(\mathbb{N} + \mathbb{L}\mathbb{N}). $$

A call $\llbracket f \rrbracket_c^m l$ uses $c$ as step-index, $m$ as auxiliary counter to implement unbounded search, and $l$ as input. If $\llbracket f \rrbracket_c^m l = \lfloor \text{inl } v \rfloor$, then $v$ is the value of the evaluation. If $\llbracket f \rrbracket_c^m l = \lfloor \text{inr } l' \rfloor$, then $f$ encoded a list of functions via $\text{cons}$ and $\text{nil}$ which pointwise evaluated to $l' : \mathbb{L}\mathbb{N}$. If $\llbracket f \rrbracket_c^m l = \emptyset$ either if usual the step-index $c$ was not big enough, the function $f$ does not terminate on input $l$, or $f$ is not a valid skeleton.

$$
\begin{align*}
\llbracket \text{cst}_n \rrbracket_{1+c}^m l &:= \lfloor \text{inl } n \rfloor \\
\llbracket \text{zero} \rrbracket_{1+c}^m l &:= \lfloor \text{inl } 0 \rfloor \\
\llbracket \text{succ} \rrbracket_{1+c}^m (x :: l) &:= \lfloor \text{inl } (1+x) \rfloor \\
\llbracket \text{prj}_p \rrbracket_{1+c}^m l &:= \lfloor \text{inl } x \rfloor \quad \text{for } l_p = \lfloor x \rfloor \\
\llbracket \text{comp f g} \rrbracket_{1+c}^m l &:= \lfloor \text{inl } x \rfloor \quad \text{for } \llbracket g \rrbracket_c^m l = \lfloor \text{inr } l' \rfloor \text{ and } \llbracket f \rrbracket_c^m l' = \lfloor \text{inl } x \rfloor \\
\llbracket \text{rec f g} \rrbracket_{1+c}^m (0 :: l) &:= \lfloor \text{inl } x \rfloor \quad \text{for } \llbracket f \rrbracket_c^m l = \lfloor \text{inl } x \rfloor \\
\llbracket \text{rec f g} \rrbracket_{1+c}^m ((1+n) :: l) &:= \lfloor \text{inl } x \rfloor \quad \text{for } \llbracket \text{rec f g} \rrbracket_c^m (n :: l) = \lfloor \text{inl } y \rfloor \text{ and } \llbracket g \rrbracket_c^m (n :: y :: l) = \lfloor \text{inl } x \rfloor \\
\llbracket \text{nil} \rrbracket_{1+c}^m l &:= \lfloor \text{inr } \llbracket \cdot \rrbracket \rfloor \\
\llbracket \text{cons f g} \rrbracket_{1+c}^m l &:= \lfloor \text{inr } (x :: l') \rfloor \quad \text{for } \llbracket f \rrbracket_c^m l = \lfloor \text{inl } x \rfloor \text{ and } \llbracket g \rrbracket_c^m l = \lfloor \text{inr } l' \rfloor \\
\llbracket \mu f \rrbracket_{1+c}^m l &:= \lfloor \text{inl } m \rfloor \quad \text{for } \llbracket f \rrbracket_c^0 (m :: l) = \lfloor \text{inl } 0 \rfloor \\
\llbracket \mu f \rrbracket_{1+c}^m l &:= \lfloor \text{inl } x \rfloor \quad \text{for } \llbracket f \rrbracket_c^0 (m :: l) = \lfloor \text{inl } (1+y) \rfloor \text{ and } \llbracket \mu f \rrbracket_c^{1+m} l = \lfloor \text{inl } x \rfloor \\
\llbracket f \rrbracket_c^m l &:= \emptyset \quad \text{in all other cases}
\end{align*}
$$

**Lemma 11.3.** For all $f : \mathcal{A}_k, \vec{v} : \mathbb{N}^k, c : \mathbb{N}$, and $n : \mathbb{N}$ we have

$$ [f; \vec{v}] - [c] n \leftrightarrow \llbracket f \rrbracket_c^0(v21 \vec{v}) = \lfloor \text{inl } n \rfloor $$

Note that on the right hand side of this equivalence, the vector $\vec{v}$ is just mapped to its underlying list $v21 \vec{v}$.

**Theorem 11.4.** $\mu$-rec $\leq$ WCBV.
*Proof.* We use Theorem 11.2 and have to provide a function $F$ s.t. $\mu$-rec($k$ : $\mathbb{N}$, $f$ : $A_k, \vec{v}$ : $\mathbb{N}^k$) $\leftrightarrow \exists c, F_c(k, f, \vec{v}) = \text{true}$.

We define $F(k : \mathbb{N}, f : A_k, \vec{v} : \mathbb{N}^k) c := \text{true}$ if $[[f]]_c^0(\mathbf{v} \mathbb{1} \vec{v}) = \lfloor \text{inl } n \rfloor$ for some $n$ and $F(k, f, v) c := \text{false}$ otherwise.

Now $F$ is part of the Coq fragment supported by the extraction framework in [FK19], since the step-indexed evaluator $[[\cdot]]$ does not use any dependent, nested, mutual, or non-simple types. Thus $F$ is L-computable.

Furthermore we have as wanted that

$$ \begin{align*} \mu\text{-rec}(k, f, \vec{v}) &\leftrightarrow \exists n, [[f]] \vec{v} n \\ &\leftrightarrow \exists n c, [f; \vec{v}]-[c] n \\ &\leftrightarrow \exists n c, [[f]]_c^0(\mathbf{v} \mathbb{1} \vec{v}) = \lfloor \text{inl } n \rfloor \\ &\leftrightarrow \exists c, F(k, f, \vec{v}) c = \text{true} \end{align*} \quad \square $$

## 12. RELATED WORK

Regarding formalisations of Hilbert's tenth problem, there are various unfinished and preliminary results in different proof assistants: Carneiro [Car18] formalises Matiyasevich's theorem (Diophantineness of exponentiation) in Lean, but does not consider computational models or the DPRM theorem. Pak formalises results regarding Pell's equation [Pak17] and proves that Diophantine sets are closed under union and intersection [Pak18], both as parts of the Mizar Mathematical Library. Stock et al. [S+18, BDP+19] report on an unfinished formalisation of the DPRM theorem in Isabelle based on [Mat00]. They cover some parts of the proof, but acknowledge for important missing results like Lucas's or "Kummer's theorem" and a "formalisation of a register machine". Moreover, none of the cited reports considers the computability of the reductions involved or the verification of a universal machine in the chosen model of computation yet, one of them being a necessary proof goal for an actual undecidability result in the classical meta-theories of Isabelle/HOL and Mizar.

Regarding undecidability proofs in type theory, Forster, Heiter, and Smolka [FHS18] reduce the halting problem of Turing machines to PCP. Forster and Larchey-Wendling [FL19] reduce PCP to provability in linear logic via the halting problem of Minsky machines, which we build on. Forster, Kirst and Smolka develop the notion of synthetic undecidability in Coq and prove the undecidability of various notions in first-order logic [FKS19]. Spies and Forster mechanise the undecidability proof of second-order unification by reduction from H10 [SF20] originally shown by Goldfarb [Gol81]. Forster, Kunze, and Wuttke reduce the halting problem of multi-tape Turing machines to single-tape Turing machines [FKW20]. Dudenhefner and Rehof [DR19] mechanise a recently simplified undecidability proof for System F inhabitation.

Regarding formalisations of $\mu$-recursive functions, Larchey-Wendling [Lar17] shows that every total $\mu$-recursive function can directly be computed in Coq and Carneiro [Car19] mechanises standard computability theory based on $\mu$-recursive functions. Xu, Zhang, and Urban [XZU13] mechanise $\mu$-recursive functions and Turing machines in Isabelle and prove their computational equivalence. Their proof uses Abacus machines as intermediate layer, which are similar to our Minsky machines.
ACKNOWLEDGEMENTS

We would like to thank Gert Smolka, Dominik Kirst, and Simon Spies for helpful discussion regarding the presentation. The first author was partially supported by the TICAMORE project (ANR grant 16-CE91-0002).

REFERENCES

[BDP<sup>+</sup>19] Jonas Bayer, Marco David, Abhik Pal, Benedikt Stock, and Dierk Schleicher. The DPRM Theorem in Isabelle (Short Paper). In *10th International Conference on Interactive Theorem Proving (ITP 2019)*, volume 141 of *Leibniz International Proceedings in Informatics (LIPIcs)*, pages 33:1–33:7, Dagstuhl, Germany, 2019. Schloss Dagstuhl–Leibniz-Zentrum fuer Informatik. doi:10.4230/LIPIcs.ITP.2019.33.

[Car18] Mario Carneiro. A Lean formalization of Matiyasevič's theorem, 2018. arXiv: 1802.01795.

[Car19] Mario Carneiro. Formalizing Computability Theory via Partial Recursive Functions. In *10th International Conference on Interactive Theorem Proving (ITP 2019)*, volume 141 of *Leibniz International Proceedings in Informatics (LIPIcs)*, pages 12:1–12:17, Dagstuhl, Germany, 2019. Schloss Dagstuhl–Leibniz-Zentrum fuer Informatik. doi:10.4230/LIPIcs.ITP.2019.12.

[Con87] John H. Conway. *FRACTRAN: A Simple Universal Programming Language for Arithmetic*, pages 4–26. Springer New York, New York, NY, 1987. doi: 10.1007/978-1-4612-4808-8_2.

[Dav53] Martin Davis. Arithmetical problems and recursively enumerable predicates. *The Journal of Symbolic Logic*, 18(1):33–41, 1953. doi:10.2307/2266325.

[Dav73] Martin Davis. Hilbert's Tenth Problem is Unsolvable. *The American Mathematical Monthly*, 80(3):233–269, 1973. URL: 10.1080/00029890.1973.11993265.

[DP59] Martin Davis and Hilary Putnam. *A computational proof procedure; Axioms for number theory; Research on Hilbert's Tenth Problem*. Air Force Office of Scientific Research, Air Research and Development, 1959.

[DPR61] Martin Davis, Hilary Putnam, and Julia Robinson. The decision problem for exponential Diophantine equations. *Annals of Mathematics*, pages 425–436, 1961. doi:10.2307/1970289.

[DR19] Andrej Dudenhefner and Jakob Rehof. A Simpler Undecidability Proof for System F Inhabitation. In Peter Dybjer, José Espírito Santo, and Luís Pinto, editors, *24th International Conference on Types for Proofs and Programs (TYPES 2018)*, volume 130 of *Leibniz International Proceedings in Informatics (LIPIcs)*, pages 2:1–2:11, Dagstuhl, Germany, 2019. Schloss Dagstuhl–Leibniz-Zentrum fuer Informatik. doi:10.4230/LIPIcs.TYPES.2018.2.

[FHS18] Yannick Forster, Edith Heiter, and Gert Smolka. Verification of PCP-Related Computational Reductions in Coq. In *ITP 2018*, pages 253–269. Springer, 2018. doi:10.1007/978-3-319-94821-8_15.

[FK19] Yannick Forster and Fabian Kunze. A Certifying Extraction with Time Bounds from Coq to Call-By-Value Lambda Calculus. In *10th International Conference on Interactive Theorem Proving (ITP 2019)*, volume 141 of *Leibniz International Proceedings in Informatics (LIPIcs)*, pages 17:1–17:19, Dagstuhl, Germany, 2019.
Schloss Dagstuhl-Leibniz-Zentrum fuer Informatik. doi:10.4230/LIPIcs.ITP. 2019.17.

[FKS19] Yannick Forster, Dominik Kirst, and Gert Smolka. On synthetic undecidability in coq, with an application to the entscheidungsproblem. In *Proceedings of the 8th ACM SIGPLAN International Conference on Certified Programs and Proofs*, CPP 2019, page 38–51, New York, NY, USA, 2019. Association for Computing Machinery. doi:10.1145/3293880.3294091.

[FKSW21] Yannick Forster, Fabian Kunze, Gert Smolka, and Maximilian Wuttke. A Mechanised Proof of the Time Invariance Thesis for the Weak Call-By-Value λ-Calculus. 193:19:1–19:20, 2021. doi:10.4230/LIPIcs.ITP.2021.19.

[FKW20] Yannick Forster, Fabian Kunze, and Maximilian Wuttke. Verified Programming of Turing Machines in Coq. In *Proceedings of the 9th ACM SIGPLAN International Conference on Certified Programs and Proofs. ACM*, 2020. doi:10.1145/3372885.3373816.

[FL19] Yannick Forster and Dominique Larchey-Wendling. Certified Undecidability of Intuitionistic Linear Logic via Binary Stack Machines and Minsky Machines. In *CPP 2019*, pages 104–117. ACM, 2019. doi:10.1145/3293880.3294096.

[FLW18] Yannick Forster and Dominique Larchey-Wendling. Towards a library of formalised undecidable problems in Coq: The undecidability of intuitionistic linear logic. *Workshop on Syntax and Semantics of Low-level Languages, Oxford, 2018*.

[FLWD$^{+20}$] Yannick Forster, Dominique Larchey-Wendling, Andrej Dudenhefner, Edith Heiter, Dominik Kirst, Fabian Kunze, Gert Smolka, Simon Spies, Dominik Wehr, and Maximilian Wuttke. A Coq library of undecidable problems. In *The Sixth International Workshop on Coq for Programming Languages (CoqPL 2020)*, 2020. URL: https://github.com/uds-psl/coq-library-undecidability.

[FLWD$^{+21}$] Yannick Forster, Dominique Larchey-Wendling, Andrej Dudenhefner, Lennard Gäher, Edith Heiter, Marc Hermes, Dominik Kirst, Fabian Kunze, Maximilian Wuttke, Gert Smolka, Simon Spies, and Dominik Wehr. The Coq Library of Undecidability Proofs: Hilbert's Tenth problem in Coq (LMCS) v1.1, May 2021. doi:10.5281/zenodo.4827453.

[FS17] Yannick Forster and Gert Smolka. Weak call-by-value lambda calculus as a model of computation in coq. In *Interactive Theorem Proving*, pages 189–206, Cham, 2017. Springer International Publishing. doi:10.1007/978-3-319-66107-0_13.

[Göd31] Kurt Gödel. Über formal unentscheidbare Sätze der Principia Mathematica und verwandter Systeme I. *Monatshefte für mathematik und physik*, 38(1):173–198, 1931. doi:10.1063/1.3051400.

[Gol81] Warren D. Goldfarb. The undecidability of the secondorder unification problem. *Theoretical Computer Science*, 13:225–230, 1981. doi:10.1016/0304-3975(81) 90040-2.

[Hil02] David Hilbert. Mathematical problems. *Bulletin of the American Mathematical Society*, 8(10):437–479, 1902. doi:10.1090/S0002-9904-1902-00923-3.

[JM84] J. P. Jones and Y. V. Matijasevič. Register Machine Proof of the Theorem on Exponential Diophantine Representation of Enumerable Sets. *J. Symb. Log.*, 49(3):818–829, 1984. doi:10.2307/2274135.
[Lar17] Dominique Larchey-Wendling. Typing Total Recursive Functions in Coq. In *ITP 2017*, pages 371–388. Springer, 2017. doi:10.1007/978-3-319-66107-0_24.

[Luc78] Edouard Lucas. Théorie des Fonctions Numériques Simplement Périodiques. [Continued]. *American Journal of Mathematics*, 1(3):197–240, 1878. doi:10.2307/2369308.

[LWF19] Dominique Larchey-Wendling and Yannick Forster. Hilbert's Tenth Problem in Coq. In *4th International Conference on Formal Structures for Computation and Deduction (FSCD 2019)*, volume 131 of *Leibniz International Proceedings in Informatics (LIPIcs)*, pages 27:1–27:20, Dagstuhl, Germany, 2019. Schloss Dagstuhl–Leibniz-Zentrum fuer Informatik. doi:10.4230/LIPIcs.FSCD.2019.27.

[Mat70] Yuri V. Matijasevič. Enumerable sets are Diophantine. In *Soviet Mathematics: Doklady*, volume 11, pages 354–357, 1970.

[Mat97] Yuri V. Matiyasevich. A new technique for obtaining Diophantine representations via elimination of bounded universal quantifiers. *J. Math. Sci.*, 87(1):3228–3233, 1997.

[Mat00] Yuri V. Matiyasevich. On Hilbert's Tenth Problem. Expository Lectures 1, Pacific Institute for the Mathematical Sciences, University of Calgary, February 2000. URL: http://www.mathtube.org/sites/default/files/lecture-notes/Matiyasevich.pdf.

[Mat16] Yuri V. Matiyasevich. Martin davis and hilbert's tenth problem. In *Martin Davis on Computability, Computational Logic, and Mathematical Foundations*, volume 10 of *Outstanding Contributions to Logic*, pages 35–54. Springer, 2016. doi:10.1007/978-3-319-41842-1_2.

[Min67] Marvin L. Minsky. *Computation: finite and infinite machines*. Prentice-Hall, Inc., 1967.

[Pak17] Karol Pak. The Matiyasevich Theorem. Preliminaries. *Formalized Mathematics*, 25(4):315–322, 2017. doi:10.1515/forma-2017-0029.

[Pak18] Karol Pak. Diophantine sets. Preliminaries. *Formalized Mathematics*, 26(1):81–90, 2018. doi:10.2478/forma-2018-0007.

[Pos44] Emil L. Post. Recursively enumerable sets of positive integers and their decision problems. *bulletin of the American Mathematical Society*, 50(5):284–316, 1944. doi:10.1090/S0002-9904-1944-08111-1.

[Rob52] Julia Robinson. Existential definability in arithmetic. *Transactions of the American Mathematical Society*, 72(3):437–449, 1952. doi:10.1090/S0002-9947-1952-0048374-2.

[S+18] Benedikt Stock et al. Hilbert meets Isabelle: Formalisation of the DPRM theorem in Isabelle. *Isabelle Workshop 2018*, 2018. doi:10.29007/3q4s.

[SAB+20] Matthieu Sozeau, Abhishek Anand, Simon Boulier, Cyril Cohen, Yannick Forster, Fabian Kunze, Gregory Malecha, Nicolas Tabareau, and Théo Winterhalter. The MetaCoq project. *J. Autom. Reason.*, 64(5):947–999, 2020. doi:10.1007/s10817-019-09540-0.

[SF20] Simon Spies and Yannick Forster. Undecidability of higher-order unification formalised in Coq. In *Proceedings of the 9th ACM SIGPLAN International Conference on Certified Programs and Proofs*, pages 143–157, 2020. doi:10.1145/3372885.3373832.
[SM19] Matthieu Sozeau and Cyprien Mangin. Equations reloaded: High-level dependently-typed functional programming and proving in coq. *Proc. ACM Program. Lang.*, 3(ICFP), July 2019. doi:10.1145/3341690.

[tea] The Mathematical Components team. Mathematical Components. URL: https://math-comp.github.io.

[Tea21] The Coq Development Team. The Coq Proof Assistant, January 2021. doi: 10.5281/zenodo.1003420.

[XZU13] Jian Xu, Xingyuan Zhang, and Christian Urban. Mechanising turing machines and computability theory in isabelle/hol. In *Interactive Theorem Proving*, pages 147–162, Berlin, Heidelberg, 2013. Springer Berlin Heidelberg. doi: 10.1007/978-3-642-39634-2_13.

## APPENDIX A. SOME REMARKS ABOUT THE COQ CODE CONTENTS

The file names below are hyper-linked to the corresponding files in the following specific release of the Coq library of undecidability proofs:

https://github.com/uds-psl/coq-library-undecidability/tree/H10-LMCS-v1.1

To help at understanding the H10 code from a high-level point of view, we provide two additional Coq source code files:

* the file *summary.v* which gives specific pointers to the main problems contained in the reduction chain described in this paper;

* the file *standalone.v* which gives a (combined) reduction from the binary BPCP to H10 but where the statements of those two problems have been rewritten to minimize dependencies.

Additionally, we here give a detailed overview of the structure of the code corresponding to the results presented in this paper, and which was contributed to our Coq library of undecidability proofs. The following *lines of code (loc)* measurements combine both definitions and proof scripts but do not account for comments. Notice that there are more files in the whole library than those needed to actually cover H10, but here, we only present the latter. In total, we contribute 16k loc to our undecidability project, 4k being additions to its shared libraries as extensions of the Coq standard library.

Concerning the multi-purpose shared libraries in Shared/Libs/DLW/Utils:

* we implemented finitary sums/products (over monoids) up to the binomial theorem (Newton) over non-commutative rings in *sums.v* and *binomial.v* for a total of 550 loc;

* we implemented bitwise operations over $\mathbb{N}$, both a lists of bits in *bool_list.v* and Peano *nat* in *bool_nat.v* for a total of 1700 loc;

* we implemented many results about Euclidean division and Bézout's identity in *gcd.v*, prime numbers and their unboundedness in *prime.v*, and base *p* representations in *power_decomp.v* for a total of 1200 loc;

* we implemented miscellaneous libraries for the reification of *boundedaderivation.v* (120 loc), the Pigeon Hole Principle in *php.v* (350 loc) and iterations of binary relations in *rel_iter.v* (230 loc).

Concerning the libraries for Minsky machines and FRACTRAN programs:
* by a slight update to the existing code [FL19], we proved in mm_comp.v that MM-termination (on any state) is undecidable (10 loc). Both the pre-existing result (undecidability of MM-termination on the zero state) and the new result derive from the correctness of the compiler of binary stack machines into Minsky machines;

* we implemented the removal of self loops in Minsky machines in mm_noemer.v (340 loc);

* we construct two infinite sequences of primes $p_i$ and $q_i$ in prime_seq.v (240 loc);

* FRACTRAN definitions and basic results occur in fractran_utils.v (310 loc) and the verified compiler from Minsky machines to FRACTRAN occurs in MM_FRACTRAN.v (300 loc);

Concerning the libraries for proving Matiyasevich's theorems:

* we implemented a library for modular arithmetic ($\mathbb{Z}/p\mathbb{Z}$) in Zp.v (920 loc);

* we implemented a library for 2 × 2-matrix computation including exponentiation and determinants in matrix.v (210 loc);

* we implemented an elementary proof of Lucas's theorem in luca.v (330 loc);

* we implemented the “classical proof” of Lagrange's theorem in lagrange.v (520 loc);

* the solution $\alpha_b(n)$ of Pell's equation and its (modular) arithmetic properties up to a proof of its Diophantineness are in alpha.v (1150 loc);

* from $\alpha_b(n)$, we implement the meta-level Diophantine encoding of the exponential in expo_diophantine.v (150 loc);

* we implement the sparse ciphers used in the Diophantine elimination of bounded universal quantification in cipher.v (1450 loc).

Concerning the object-level Diophantine libraries:

* the definition of Diophantine logic and basic results is in dio_logic.v (540 loc);

* the definition of elementary Diophantine constraints and the reduction from Diophantine logic is in dio_elem.v (440 loc);

* the definition of single Diophantine equations and the reduction from elementary Diophantine constraints is in dio_single.v (350 loc);

* we implement the object-level Diophantine encoding of the exponential relation in dio_expo.v (60 loc); but all the work is done in the previously mentioned libraries;

* the object-level Diophantine encoding of bounded universal quantification spans over dio_binary.v, dio_cipher.v and dio_bounded.v (460 loc);

* we derive the object-level Diophantine encoding of the reflexive-transitive closure in dio.rtClosure.v (40 loc);

* we implement the object-level Diophantine encoding of the FRACTRAN termination predicate in fractran_dio.v (80 loc).

Concerning $\mu$-recursive algorithms, reducing H10 and reduced to Minsky machines:

* building on the pre-existing developments corresponding to [Lar17] (1000 loc);

* the $\mu$-recursive algorithm that searches for a solution to a single Diophantine equation in ra_utils.v, ra_dio_pol.v, recomp.v and ra_recomp.v (1370 loc);

* extensions of the Minsky machines library for N-indexed registers in env.v, mme_def.v and mme_utils.v (600 loc);

* the certified compiler from $\mu$-recursive algorithms to Minsky machines in ra_mm.v and ra_mm_env.v (1250 loc).

Concerning the reduction from $\mu$-recursive algorithms to L:

* the step-indexed evaluator and the reduction are in MuRec.v (450 loc)
* the framework used to extract the step-indexed evaluator is in the directory L and consists of about 2500 lines of code. More details can be found in [FK19].

To finish, the main undecidability results and the DPRM:

* the undecidability of Minsky machines is in MM\_undec.v (20 loc);

* the reduction from MM to FRACTRAN is in MM\_FRACTRAN.v (90 loc);

* the Diophantine encoding of FRACTRAN termination is in FRACTRAN\_DIO.v (70 loc);

* the whole reduction chain leading to the undecidability of H10 is in H10\_undec.v (60 loc);

* the reduction from H10 to H10<sub>Z</sub> is in H10Z.v (200 loc)

* and the DPRM theorem is in DPRM.v (170 loc).

## APPENDIX B. LUCAS'S THEOREM

Lucas's theorem allows for the computation of the binomial coefficient $C_m^n$ modulo a prime number $p$ using the base $p$ expansions of $m$ and $n$. There are various proofs of this theorem but most of them involve high-level concepts like *generating functions* or *group action* and we choose instead to implement a low-level combinatorial proof of the theorem. Such a proof provides specific combinatorial insights into Lucas's results.

While a high-level proof could rightfully be considered as more beautiful, it would also likely require a lot of additional library code. We could of course rely on an external library such as MathComp [tea] but this would generate a new external dependency, additional to the already involved Equations [SM19] and MetaCoq [SAB<sup>+</sup>20] libraries, increasing the likelihood of synchronization issues between all these developments. An additional benefit of a low-level proof is that it could reasonably be imported in alternate developments (e.g. [BDP<sup>+</sup>19]) in other proof assistants (e.g. Isabelle/HOL) without the need to convert a substantial helper library.

Before we enter this low-level proof of Lucas's theorem, we must give a light-weight, working and formal definition of binomial coefficients. For this, we use Pascal's identity as a ground for a fixpoint definition:

$$C_m^0 := 1 \quad \text{and} \quad C_0^{1+n} := 0 \quad \text{and} \quad C_{1+m}^{1+n} := C_m^n + C_m^{1+n}$$

where we use the compact notation $C_m^n$ for a more compact typesetting of the upcoming equations. Starting from Pascal's definition, one can derive the following identities

$$n! (m-n)! C_m^n = m! \quad \text{for } n \le m \quad \text{and} \quad C_m^n = 0 \quad \text{for } m < n$$

and we call the leftmost one as the binomial identity<sup>22</sup>

**Lemma B.1.** Let $p$ be a prime number and let us consider $M = mp + m_0$ and $N = np + n_0$ with $m_0, n_0 < p$. Then the identity $C_M^N = C_m^n C_{m_0}^{n_0} [p]$ holds.

*Proof.* Let us first notice that the identity is trivial (in $\mathbb{N}$ already) when $n > m$ or $n = m \land n_0 > m_0$ because both sides evaluate to 0. So below, we only consider the cases $n \le m \land n_0 \le m_0$ or $n < m \land m_0 < n_0$ which cover the remaining part of the domain.

Before we split those two cases, let us define $\phi: \mathbb{N} \to \mathbb{N}$ and $\psi: \mathbb{N} \to \mathbb{N}$ by

$$\phi_r^i := (ip + 1) \cdots (ip + r) \quad \text{and} \quad \psi_i := \phi_{p-1}^0 \cdots \phi_{p-1}^{i-1}.$$

<footnote n="22"><sup>22</sup>We do not enter the details of these inductive proofs which are standard exercises in basic arithmetic.</footnote>
Notice that $\phi_r^i$ generalize the factorial function as $r! = \phi_r^0$ but we will always use $\phi_r^i$ with $r < p$. Projected on $\mathbb{Z}/p\mathbb{Z}$, both $\phi$ and $\psi$ simplify to factorial as it is easy to show$^{23}$

$$ \phi_r^i \equiv r! [p] \quad \text{and} \quad \psi_i \equiv (p-1)! [p]. \qquad (B.1) $$

Moreover, $\phi_r^i$ is invertible in $\mathbb{Z}/p\mathbb{Z}$ for any $r < p$, and in particular $(p-1)! = \phi_{p-1}^0$ and $\psi_i$ are also invertible in $\mathbb{Z}/p\mathbb{Z}$. Using $\phi$, $\psi$ and the semiring structure of $\mathbb{N}$, we establish the identity

$$ (ip + r)! = \phi_{p-1}^0 \cdot (1p + 0) \cdots \phi_{p-1}^{i-1} \cdot (ip + 0) \cdot \phi_r^i \psi_i \qquad (B.2) $$

In the case $n \le m$ and $n_0 \le m_0$, we have $M - N = (m - n)p + (m_0 - n_0)$ and we rewrite the binomial equation $N!(M - N)!C_M^N = M!$ in $\mathbb{N}$ as

$$ (np + n_0)! ((m - n)p + (m_0 - n_0))! C_M^N = (mp + m_0)! $$

and then rewrite again with Eq. (B.2) and the binomial equation $m! = n!(m-n)!C_m^n$ as

$$ n! p^n \phi_{n_0}^n \psi_n \cdot (m - n)! p^{m-n} \phi_{m_0 - n_0}^{m-n} \psi_{m-n} \cdot C_M^N = n! (m - n)! C_m^n \cdot p^m \phi_{m_0}^m \psi_m $$

In $\mathbb{N}$, we simplify by $p^m p^{m-n} = p^m$ and $n!(m - n)!$ and get

$$ \phi_{n_0}^n \psi_n \phi_{m_0 - n_0}^{m-n} \psi_{m-n} C_M^N = C_m^n \phi_{m_0}^m \psi_m $$

Now switching to $\mathbb{Z}/p\mathbb{Z}$ using Eqs. (B.1), we derive:

$$ n_0! (p-1)!^n (m_0 - n_0)! (p-1)!^{m-n} C_M^N = C_m^n m_0! (p-1)!^m [p] $$

Because $(p-1)!$ is invertible in $\mathbb{Z}/p\mathbb{Z}$, we can simplify by $(p-1)!^n(p-1)!^{m-n} = (p-1)!^m$. Then we rewrite using the binomial equation $m_0! = n_0!(m_0 - n_0)!C_{m_0}^n$ and deduce

$$ n_0! (m_0 - n_0)! C_M^N = C_m^n m_0! = C_m^n n_0! (m_0 - n_0)! C_{m_0}^{n_0} [p] $$

Finally, as both $n_0 < p$ and $m_0 - n_0 < p$ hold, then $n_0!(m_0 - n_0)!$ is invertible in $\mathbb{Z}/p\mathbb{Z}$ and, simplifying, we get $C_M^N = C_m^n C_{m_0}^{n_0} [p]$ as required.

Then we consider the alternative case where $n < m$ and $m_0 < n_0 < r$. In this case we have $M - N = (m - (n+1))p + (p - (n_0 - m_0))$ and again we develop the binomial equation $N!(M - N)!C_M^N = M!$ in $\mathbb{N}$ as

$$ (np + n_0)! ((m - (n + 1))p + (p - (n_0 - m_0)))! C_M^N = (mp + m_0)! $$

and thus, using Eq. (B.2) and the binomial equation $m! = n!(m-n)!C_m^n$, rewriting it into

$$ n! p^n \phi_{n_0}^n \psi_n (m - (n + 1))! p^{m-(n+1)} \phi_{p-(n_0-m_0)}^{m-(n+1)} \psi_{m-(n+1)} C_M^N = n! (m - n)! C_m^n p^m \phi_{m_0}^m \psi_m $$

that we then simplify by $n!(m - (n + 1))!$ and $p^n p^{m-(n+1)} = p^{m-1}$ to get

$$ \phi_{n_0}^n \psi_n \phi_{p-(n_0-m_0)}^{m-(n+1)} \psi_{m-(n+1)} C_M^N = (m-n) C_m^n p \phi_{m_0}^m \psi_m $$

Notice that $p$ appears (at least once) on the right-hand side, thus switching to $\mathbb{Z}/p\mathbb{Z}$, we get:

$$ n_0! (p-1)!^n (p - (n_0 - m_0))! (p-1)!^{m-(n+1)} C_M^N = 0 [p] $$

Now, as $n_0 < p$ and $p - (n_0 - m_0) < p$, the factorials $n_0!$ and $(p - (n_0 - m_0))!$ are both invertible in $\mathbb{Z}/p\mathbb{Z}$. This is also the case of any power of $(p-1)!$. We deduce that $C_M^N = 0$ [p] which leads us to conclude $C_M^N = C_m^n C_{m_0}^{n_0} [p]$ as $C_{m_0}^{n_0} = 0$ trivially follows from $m_0 < n_0$. $\square$

<footnote>$^{23}$Notice that Wilson's theorem establishes $(p-1)! \equiv -1$ [p] but this equation is not needed in this proof.</footnote>
**Theorem B.2 (Lucas [Luc78]).** Whenever $u_0, \dots, u_n, v_1, \dots, v_n < p$, the following identity holds:

$$ \begin{pmatrix} u_n p^n + \dots + u_0 \\ v_n p^n + \dots + v_0 \end{pmatrix} \equiv \begin{pmatrix} u_n \\ v_n \end{pmatrix} \dots \begin{pmatrix} u_0 \\ v_0 \end{pmatrix} [p]. $$

*Proof.* By induction on the number $n$ of coefficients using Lemma B.1. □

**Definition B.3.** Let $a = a_n 2^n + \dots + a_0$ and $b = b_n 2^n + \dots + b_0$ representations in base 2. We denote $a \prec b$ if $a_i \le b_i$ holds for any $i \le n$.

Notice that the definition of $a \prec b$ is irrelevant to which base 2 representations of $a$ and $b$ are picked up, i.e. it is not influenced by trailing zeros that might appear in front of (the representations of) $a$ or $b$.

**Corollary B.4.** For any $a, b: \mathbb{N}$ we have $a \prec b \leftrightarrow C_b^a = 1$ [2].

## APPENDIX C. AVOIDING OVERFLOWS IN THE PROOF OF THEOREM 2.8

The section explains why we slightly modified the original proof of the elimination of bounded universal quantification [Mat97] to avoid overflows when dealing with parallel multiplications, in search for a Diophantine encoding of the simultaneous equations $a_1 = b_1 c_1, \dots, a_n = b_n c_n$, which can also be viewed as the identity between the vectors $\vec{a} = [a_1; \dots; a_n]$ and the scalar/dot product $\vec{b} \cdot \vec{c} = [b_1 c_1; \dots; b_n c_n]$.

Following Matiyasevich [Mat97], the component of those vectors are encoded sparsely within the digits in base $r$ of the corresponding so called ciphers (see below for examples). Using the natural product of ciphers, he can somehow recover the identity between the components of $\vec{a}$ and those of $\vec{b} \cdot \vec{c}$. The basis $r$ for the encoding is chosen large enough w.r.t. the component of $\vec{a}, \vec{b}$ and $\vec{c}$ to avoid digit overflow during the natural product of ciphers. However, this product generates artefacts on in-between digits that need to be masked out and, with Matiyasevich's choice for $r$ as $r := 2^{2q}$, some overflow could occur within these artefacts, a subtlety which we speculate he might have missed. This overflow pointlessly complicates the correctness proof of masking. In our Coq code, we propose to avoid any digit overflow by increasing $r$ to $r := 2^{4q}$, but, as we realized later on, any power of 2 greater than $2^{2q+1}$ would work as well.

Let us switch to the technical details explaining how the overflow occurs and how it can be avoided, with little impact on the rest of the argument. For the moment, we describe what Matiyasevich does in [Mat97] so we stick with his definition of $r$ as $r := 2^{2q}$ and consider the encoding of vectors of $D_q^n$ as ciphers, written in base $r$, where $q \ge 1$ and $D_q := \{0, \dots, 2^q - 1\}$ is the set of allowed digits for vectors. He considers five vectors and their respective ciphers:

<table><thead><tr><th>vector/length</th><th>cipher</th></tr></thead><tbody><tr><td>$\vec{u} = [1; \dots; 1] \in D_q^n$</td><td>$u = 1 \cdot r^{21} + \dots + 1 \cdot r^{2n}$</td></tr><tr><td>$\vec{u}' = [0; 1; \dots; 1] \in D_q^{n+1}$</td><td>$u' = 0 \cdot r^{21} + 1 \cdot r^{22} + \dots + 1 \cdot r^{2n+1}$</td></tr><tr><td>$\vec{a} = [a_1; \dots; a_n] \in D_q^n$</td><td>$a = a_1 r^{21} + \dots + a_n r^{2n}$</td></tr><tr><td>$\vec{b} = [b_1; \dots; b_n] \in D_q^n$</td><td>$b = b_1 r^{21} + \dots + b_n r^{2n}$</td></tr><tr><td>$\vec{c} = [c_1; \dots; c_n] \in D_q^n$</td><td>$c = c_1 r^{21} + \dots + c_n r^{2n}$.</td></tr></tbody></table>

Each summand is an actual digit of the cipher in base $r$. Notice that $q$ is to be chosen large enough so that the digits in $D_q$ cover the components of the vectors $\vec{a}, \vec{b}$ and $\vec{c}$. Also,
the ciphers have a sparse representation in base $r$, i.e. only the digits which correspond to powers $r^{2^i}$ are possibly non-zero. Finally, most of the digits of base $r$ are unused, i.e. only $2^q$ of them out of $r = 2^{2q}$ many, to avoid overflows when adding or multiplying digits. But remark that which overflows are avoided (or not) is the purpose of the whole discussion here.

The problem solved in [Mat97] is to find, depending on $q/r$ and $n$, Diophantine representations for (the ciphers of) $\vec{u}$ and $\vec{u}'$, and most importantly, for the identities $\vec{a} = \vec{b} + \vec{c}$ (the sum) and $\vec{a} = \vec{b} \cdot \vec{c}$ (the scalar product), implementing some kind of parallel addition and multiplication. We skip over the Diophantine encodings of $\vec{u}$ and $\vec{u}'$, of which the description spans about one fourth of [Mat97], to directly consider the parallel sum and multiplication. For $\vec{a} = \vec{b} + \vec{c}$, the encoding is as straightforward as the addition of ciphers $a = b + c$, because the basis $r$ is large enough to avoid overflows of digit additions, i.e. $2(2^q - 1) < r = 2^{2q}$, and no artefacts appear at powers different from $r^{2^i}$.

The case of parallel multiplication is more complicated. To recover the identity between $\vec{a}$ and the scalar product $\vec{b} \cdot \vec{c} = [b_1 c_1; \dots; b_n c_n]$, Matiyasevich also uses the natural product of the ciphers but it generates artefacts on in-between digits that need to be masked out. Indeed, considering the identities

$$au = \sum_{i=1}^{n} a_i r^{2^i} \times \sum_{i=1}^{n} r^{2^i} = \sum_{i=1}^{n} a_i r^{2^i+1} + \sum_{1 \le i < j \le n} (a_i + a_j) r^{2^i+2^j} \quad (C.1)$$

and

$$bc = \sum_{i=1}^{n} b_i r^{2^i} \times \sum_{i=1}^{n} c_i r^{2^i} = \sum_{i=1}^{n} b_i c_i r^{2^i+1} + \sum_{1 \le i < j \le n} (b_i c_j + b_j c_i) r^{2^i+2^j} \quad (C.2)$$

by filtering out the two parts $\sum_{1 \le i < j \le n} \dots$ on the right hand side of the $+$ sign, and then identifying the two values, we may recover the identity $\vec{a} = \vec{b} \cdot \vec{c}$.

For this, he uses $(r-1)u' = (r-1)r^{2^2} + \dots + (r-1)r^{2^{n+1}}$ as a mask. This is the meaning of Equation (40) of page 3232 in [Mat97], (slightly) reformulated here as

$$au \& (r-1)u' = bc \& (r-1)u' \quad (C.3)$$

where $\&$ is the digit by digit AND operator.<fnref n="24" /> However, to show the equivalence between Equation (C.3) and the identity $\vec{a} = \vec{b} \cdot \vec{c}$, we have to make sure that the masking operator $x \mapsto x \& (r-1)u'$ effectively erases the two rightmost sums $\sum_{1 \le i < j \le n} \dots$ in Equations (C.1) and (C.2), while keeping the left sums intact. Proving this property is not even considered worthy of an explanation in [Mat97] but it turned out to be not that obvious.

We first remark that the digit by digit AND operator & does not commute with $+$ (in general), hence computing the mask can be tricky. However, & commutes with the digit by digit OR operator. Moreover OR and $+$ coincide whenever the binary digits of the added numbers do not overlap (there are never two 1 on the same binary digit power in each number). As $r$ is a power of 2, this generalises from base 2 to base $r$, and provided the below summations are base $r$ representations, one can commute the sum/+ with $\&$ and apply masking component by component, e.g. one can show the identity

$$\left( \sum_{i \in I} x_i r^{\varphi(i)} \right) \& \left( \sum_{i \in I} m_i r^{\varphi(i)} \right) = \sum_{i \in I} (x_i \& m_i) r^{\varphi(i)}$$

<footnote n="24">²⁴The Diophantine encoding of the AND operator is itself derived from Lucas's theorem.</footnote>
as soon as $x_i < r$ and $m_i < r$ holds for any $i \in I$, and $\phi: I \to \mathbb{N}$ injectively maps the finite type $I$. But this implies that we have to check that Equations (C.1) and (C.2) contain such genuine base $r$ representations.

In that spirit, if $i < j$ holds then the powers $r^{2^k}$ and $r^{2^i+2^j}$ are never the same power of $r$. So these represent two distinct digit positions in base $r$ and hence, in Equations (C.1) and (C.2), those digits from the left of $+$, and those from the right of $+$ do not overlap/interfere. This also motivates the use of the mask $(r-1)u' = (r-1)r^{2^2} + \cdots + (r-1)r^{2^{n+1}}$ which indeed filters out digits at any $r^{2^i+2^j}$ position while leaving those at $r^{2^i}$ intact.

The way Equations (C.1) and (C.2) are written unfortunately suggests that these are base $r$ representations of both au and bc. However, formally, we have to verify that the two rightmost sums $\sum_{1 \le i < j \le n} \dots$ are proper base $r$ representations. This is not the case for (C.2) and we here give a counter example. On the one hand, it is true that $2^i + 2^j = 2^{i'} + 2^{j'}$ implies $(i, j) = (i', j')$ when $i < j$ and $i' < j'$, ensuring that digits at those powers of $r$ do not accidentally overlap. On the other hand, it possible that the value $b_i c_j + b_j c_i$ overflows the digit range of base $r$, even though considerable room was initially reserved to avoid that situation. For instance, when $b_i = c_j = b_j = c_i = 2^q - 1$ (the maximum allowed digit in $D_q$), then we get $b_i c_j + b_j c_i = 2^{2q+1} - 2^{q+2} + 2 \ge r + 2$ as soon as $q > 1$. So this component at digit position $r^{2^i+2^j}$ overflows and we would have to consider the spill out at digit position $r^{2^i+2^j+1}$ and show that it does not overlap with the other digits or other spill outs. This is a property which we think holds true but the proof of it would add significant formal complexity.

On the other hand, with our proposed alternative choice of $r := 2^{4q}$, keeping the same set $D_q = \{0, \dots, 2^q - 1\}$ for the allowed digits in vectors, the overflow does not occur any more, i.e. no spill out, and Equations (C.1) and (C.2) become genuine base $r$ representations, allowing smooth component by component masking.

Because there is an obvious way out of a tricky overflow management problem, we think that it is possible that Matiyasevich simply did not notice the eventuality of an overflow in the to be masked parts of Equation (C.2). We speculate this because this overflow does not occur in the parts which are left intact by the mask, at the left of the $+$ sign. With this remark, we do not imply that Equation (40) of [Mat97] is improper in any way. However, the formal proof of its equivalence with the identity $\vec{a} = \vec{b} \cdot \vec{c}$ is really more complicated when overflows occur and we did not try to prove it in the case where $r := 2^{2q}$, this situation being straightforward to avoid with a change (of otherwise low impact) in the value of $r$.

## APPENDIX D. PROOF OF PROPOSITION 4.1

Let $(p_1, q_1), \dots, (p_n, q_n)$ be a sequence of pairs in $\mathbb{N} \times \mathbb{N}$. We establish the equivalence:

$$ \sum_{i=1}^{n} 2p_i q_i = \sum_{i=1}^{n} p_i^2 + q_i^2 \iff p_1 = q_1 \wedge \cdots \wedge p_n = q_n. $$

*Proof.* We give an elementary arithmetic justification of the result, proof which involves none of the high-level tools of mathematical analysis. We first show the two following statements

$$ 2ab \le a^2 + b^2 \quad \text{and} \quad 2ab = a^2 + b^2 \iff a = b \quad \text{for any } a, b : \mathbb{N} \quad (D.1) $$

Assuming w.l.o.g. that $a \le b$, we can write $b = a + \delta$ with $\delta \in \mathbb{N}$ and then, for $\bowtie \in \{\le, =\}$ we have $2ab \bowtie a^2 + b^2 \leftrightarrow 2a^2 + 2a\delta \bowtie a^2 + a^2 + 2a\delta + \delta^2 \leftrightarrow 0 \bowtie \delta^2$ hence the desired result.
From the left inequality (D.1), we easily generalize by induction on $n$ and obtain the following inequality:

$$ \sum_{i=1}^{n} 2p_i q_i \le \sum_{i=1}^{n} p_i^2 + q_i^2 \qquad (D.2) $$

Now we can proceed with the proof of the main stated equivalence. The *if case* is obvious so we only describe the *only if case*. Hence we show that $\sum_{i=1}^{n} 2p_i q_i = \sum_{i=1}^{n} p_i^2 + q_i^2$ entails $p_i = q_i$ for all $i \in [1, n]$. We proceed by induction on $n$ again. The base case $n = 0$ is trivial. For the inductive step $1 + n$, let us assume

$$ \sum_{i=1}^{n} 2p_i q_i + 2p_{n+1} q_{n+1} = \sum_{i=1}^{n} (p_i^2 + q_i^2) + p_{n+1}^2 + q_{n+1}^2. \qquad (D.3) $$

By the left inequality of (D.1) and inequality (D.2), we have both $2p_{n+1}q_{n+1} \le p_{n+1}^2 + q_{n+1}^2$ and $\sum_{i=1}^{n} 2p_i q_i \le \sum_{i=1}^{n} p_i^2 + q_i^2$. The only possibility for the identity in hypothesis (D.3) to hold is that both inequalities are in fact identities, hence both $2p_{n+1}q_{n+1} = p_{n+1}^2 + q_{n+1}^2$ and $\sum_{i=1}^{n} 2p_i q_i = \sum_{i=1}^{n} p_i^2 + q_i^2$ hold. From this we derive $p_{n+1} = q_{n+1}$ by the equivalence on the right of (D.1) and $p_1 = q_1, \dots, p_n = q_n$ by the induction hypothesis. $\square$

<figure><img src="image_7.png" /></figure>