Synthesis of Polynomial Programs]Template-based Program Synthesis using Stellensätze

Amir Kafshdar Goharshady
Authors are ordered alphabetically.

Department of Computer Science and Engineering
Department of Mathematics
Hong Kong University of Science and Technology
Clear Water Bay, Kowloon
Hong Kong

S. Hitarth

Department of Computer Science and Engineering
Hong Kong University of Science and Technology
Clear Water Bay, Kowloon
Hong Kong

Fatemeh Mohammadi

Department of Mathematics and Department of Computer Science, KU Leuven, Belgium
KU Leuven

Harshit Jitendra Motwani

Department of Mathematics: Algebra and Geometry
Ghent University

§ INTRODUCTION

Program Synthesis Program synthesis is one of the most central research topics not only in the programming languages and verification communities, but in computer science as a whole <cit.>. It is often called the holy grail of computing <cit.>. Generally speaking, the goal in program synthesis is to automate the process of programming, i.e. given a specification, which is often a logical formula, automatically find a program (in a predefined language) that satisfies the specification.

Short History Various formulations of this problem have been studied since the 1950s (e.g. see <cit.> which was originally published in 1957). Similar problems were also considered in the context of constructive mathematics as far back as the 1930s <cit.>. Today, there is a rich body of literature that uses widely different techniques for synthesis. Some notable approaches to program synthesis include deductive synthesis using theorem provers <cit.>, counterexample-guided synthesis <cit.>, evolutionary algorithms <cit.>,  sketching/syntax-guided <cit.>, template-based/search-based <cit.>, proof-theoretic methods <cit.>, constraint-solving <cit.>, type-guided <cit.>, resource-guided <cit.>, and synthesis based on examples (input-output pairs) <cit.>. This is by no means a comprehensive list and there is no way we can do justice to all the work in this area. Moreover, most approaches, including ours, utilize more than one of the paradigms named above. Therefore, we refer to <cit.> for a more detailed discussion.

Search-based/Syntax-guided Synthesis While the initial classical formulations of program synthesis focused on automated generation of the whole program, given only a specification as a logical/semantic-based formula, it was later realized that scaling this approach is quite challenging and much efficiency can be gained by allowing the programmer to also provide syntactic templates that limit the search space of possible implementations <cit.>. Intuitively speaking, this would let the programmer focus on the high-level ideas behind the program, while relying on the synthesis engine to fill in the low-level details.

Various formalizations of similar ideas have appeared in the literature. For example, in <cit.>, the programmer has to specify a background theory, a correctness specification, which is semantic in nature, and a grammar which is used to syntactically limit the search space of possible implementations. In sketching <cit.>, the programmer provides the high-level structure or “sketch” of the implementation, while leaving low-level “holes” for the synthesis engine to fill. In template-based approaches, such as <cit.>, the programmer additionally provides a template for each hole, e.g. specifying that the hole should be filled by a linear/affine expression or that it should likely utilize a specific set of variables. These approaches have all led to highly successful synthesis engines. In particular, the SyGuS competition <cit.> is a great showcase of the remarkable year-on-year progress in program synthesis.

Completeness Given that most variants of program synthesis are undecidable, often due to Rice's theorem <cit.>, it is no surprise that the synthesis algorithms cannot be sound, complete and terminating at the same time. The undecidability also applies to almost all cases of syntax-guided synthesis, since they are more expressive than the verification of specific non-trivial semantic properties. As such, most synthesis engines, including the ones mentioned above, choose to sacrifice completeness and focus on efficient sound algorithms that can handle a wide variety of real-world use-cases. Completeness can only be achieved in limited circumstances, such as template-based synthesis of linear/affine programs with linear/affine invariants <cit.>. So, a natural question is to explore and characterize families of programs for which we can find efficient, sound and complete synthesis procedures. Of course, the difficulty lies in obtaining all three desired properties at the same time and not compromising on any of them.

Our Setting In this work, we consider template-based synthesis over polynomial programs, i.e. imperative programs with bounded real variables in which all expressions appearing in the assignments, loop guards, conditionals, and invariants (pre and post-conditions) are polynomials over program variables. See Section <ref> for a more formal statement. Following the sketching/template-based approach, we assume that the programmer provides a sketch of the program, together with one or more holes. The holes can appear both in the program itself or in the invariants and will be filled by polynomial expressions. The programmer can also provide a template for how each hole should be filled by specifying the program variables that can appear in the hole and the expected degree of the polynomial expression used for filling it. Note that this is without any loss of generality, since, in the absence of such information, our algorithm can simply iterate over all possible templates with a bounded degree.

@real: $i$, $s$, $n$;
@pre: $n \geq 1$;
$i = 0$;
$s = 0$;
@invariant: $\codetemplate{[\{i, s\},2]} \geq 0~\wedge~ \codetemplate{[\{i, s\},2]} \geq 0 ~\wedge~ \codetemplate{[\{i, n\}, 1]} \geq 0$
while($~\codetemplate{[\{i, n\}, 1]}~$)

$i$ = $i+1$;
$s$ = $\codetemplate{[\{i, s\},1]}$;

@post: $(n-1) \cdot n / 2 \leq s \leq n \cdot (n+1) / 2$;

An example polynomial program with holes.

As a simple example, consider the program in Figure <ref>, in which we have three bounded real variables $i, s,$ and $n$. The programmer has given us a sketch of the program, which includes the desired specification as a precondition and a postcondition. However, it also includes holes that should be filled by the synthesis engine. These holes are shown by dashed boxes. In this case, the programmer asks us to synthesize an invariant and a guard for the while loop, as well as an expression for the right-hand side of the assignment inside the while loop. Moreover, she has provided us with a template for each hole, specifying the variables that she expects would need to appear in that hole and the maximum degree of the polynomial that should be synthesized for filling the hole. For example, the while guard must be an affine (degree 1) polynomial over the variables $i$ and $n$. As we will see in the sequel, given this partial program as input, our approach is able to synthesize
the completed program in Figure <ref>, which basically sums up all the integers from $1$ to $\lfloor n \rfloor.$

@real: $i$, $s$, $n$;
@pre: $n \geq 1$;
$i = 0$;
$s = 0$;
@invariant: $s \geq i \cdot (i+1)/2 ~\wedge~ s \leq i \cdot (i+1)/2 ~\wedge~ i \leq n$
while($i \leq n - 1$)

$i$ = $i+1$;
$s$ = $s+i$;

@post: $(n-1) \cdot n / 2 \leq s \leq n \cdot (n+1) / 2$;

The synthesized (completed) version of the partial program in Figure <ref>.

Decidability and Hardness Note that the synthesis problem we consider here can be reduced, in polynomial time, to the decision problem of a formula in the first-order theory of the reals and is hence decidable. See Section <ref> for more details. However, the resulting formula is quite long and has a quantifier alternation. Thus, it is beyond the reach of modern SMT solvers and quantifier elimination methods, even for toy programs with three lines of code. In the same section, we also show that the problem is (strongly) NP-hard, so no strictly PTIME algorithms can be expected unless P=NP.  The NP-hardness holds even in the special case of linear/affine programs.

Our Contributions In this work, we consider the problem of template-based synthesis over polynomial imperative programs with bounded real variables. Our contributions are as follows:

* Complexity: We prove that the problem is decidable. This is achieved by a reduction to the first-order theory of the reals. We also provide a reduction showing that our problem is (strongly) NP-hard. Note that, as mentioned above, decision procedures for the first-order theory of the reals are notoriously slow and cannot handle even toy programs in practice.

* Practical Algorithm: We use classical theorems from polyhedral and real algebraic geometry, including Farkas' Lemma, Handelman's Theorem, Putinar's Positivstellensatz[Positivstellensatz is German for “positive locus theorem”. Its plural form is Positivstellensätze.], and the Real Nullstellensatz[Nullstellensatz is German for “zero locus theorem”.] to obtain a polynomial-time reduction to Quadratic Programming (QP). Since modern solvers, both SMT solvers and numerical ones, are quite efficient in handling real-world QP instances, this leads to a much more scalable approach to program synthesis in comparison to decision procedures for the first-order theory of the reals.

* Completeness: Crucially, we prove that our approach is not only sound but also semi-complete, i.e. it is complete as long as polynomials of sufficiently high degree are used in the templates. So, our approach achieves both soundness and completeness, while not suffering from the unscalability of quantifier elimination and decision procedures for the first-order theory of the reals.

* Experimental Results: We provide extensive experimental results, illustrating the efficiency of our approach and that it can successfully handle a variety of programs.

Novelty The main novelty of our approach is that it can handle a wide family of programs, i.e. polynomial imperative programs, without sacrificing completeness. It is quite rare in program synthesis to achieve soundness, semi-completeness and practical efficiency at the same time. To the best of our knowledge, our method of relying on theorems in real algebraic geometry is also novel in program synthesis. Moreover, we combine these theorems in non-trivial ways and obtain new mathematical results, which are of independent interest.

Limitations The main limitation of our approach is that it can only handle polynomial programs. This is to be expected, since we are using techniques from real algebraic geometry, which are highly dependent on having polynomial expressions/inequalities over real variables. Another limitation is that, although our approach is quite efficient in practice, its runtime is not polynomial in the size of the program or the number of variables. This is a natural consequence of our NP-hardness result and cannot be overcome unless P=NP.

Related Works Several previous works use ideas that are similar in spirit to our approach:

* Template-based Synthesis: The work <cit.> provides algorithms based on Farkas' Lemma for template-based synthesis of linear/affine programs. Since our work can handle polynomial programs of arbitrary degree, <cit.>'s results on linear/affine programs are subsumed by our approach. Moreover, while the use of Farkas' Lemma is shared among the two approaches, we also use theorems from real algebraic geometry, e.g. Handelman and Putinar, as well as non-trivial combinations of them, which are not considered in <cit.>. So, we have a different methodology and handle a wider family of programs. To the best of our knowledge, our algebro-geometric approach is novel in program synthesis.

* Invariant Generation: The work <cit.> provides template-based algorithms for automated generation of polynomial invariants over polynomial programs with real variables. Note that this is also a special case of our setting, in which the holes are limited to appear only in the invariants. Since we allow holes in both the program itself and the invariants, our setting is strictly more general. Moreover, while <cit.> also uses Positivstellensätze, it does not consider Nullstellensätze, and thus has a different mathematical basis. The work <cit.> considers the same problem as <cit.> and solves it using Gröbner basis computations. In comparison with both of these approaches, we handle a more general setting and have a different methodology.  There are also template-based approaches for generation of linear loop invariants, such as <cit.>. This is an even narrower special case of our setting. Finally, <cit.> considers invariant generation for probabilistic programs and uses a different Positivstellensatz, i.e. that of Stengle, to obtain an automated algorithm. This setting is incomparable to ours due to the presence of probabilistic behavior but absence of holes in the program itself.

* Termination analysis: Template-based algorithms that utilize Farkas' Lemma or Positivstellensätze have also been considered in the context of termination and reachability analysis <cit.>. However, this is a very different and orthogonal problem and, while there are clear similarities in the approaches, it is not possible to directly or experimentally compare them. Conceptually, the works in termination analysis assume that strong-enough invariants are given as part of the input[Such invariants can be generated by our approach or <cit.>.] and then reduce the problem of finding an affine/polynomial ranking function, or its probabilistic counterparts, to Linear Programming (LP). This is in contrast to our setting where the reduction is to Quadratic Programming (QP) and a reduction to LP is impossible due to the strong NP-hardness.

Organization In Section <ref>, we formalize our problem and define the syntax and semantics of our programs. Section <ref> is the core of the paper in which we provide our synthesis algorithm and prove that it is sound and semi-complete. We then study the complexity of the problem in Section <ref>, proving that it is decidable and strongly NP-hard. Finally, Section <ref> reports on an implementation of our algorithm and provides experimental results.

An imperative program is a sequence of instructions in a programming language that manipulate the memory to solve a problem. The task of writing the program is usually done by a programmer who first comes up with a high-level idea for the algorithm, and then they implement their algorithm by writing a program. Synthesis makes the later part easier; given logical specification for the program, we can synthesize a program that satisfies those specifications. Syntax-Guided Synthesis problem (SyGuS) module theory $T$ is a general synthesis problem of a function $f$ given semantic constraint as a formula built from symbols in theory $T$ and $f$, and syntactic constraint given as a (possibly infinite) set of expressions from $T$ specified using a context-free grammar.

In this paper, we focus on Syntax-Guided program synthesis problem module theory of Polynomial Arithmetic (over reals).

The treatment here generalises the methods for program synthesis used in <cit.> for linear programs to polynomial programs. Our main tools include the use of classical theorems from Real Algebraic Geoemtry <cit.> like Handelman's Theorem (Thoerem <ref>), Putinar’s Positivstellensatz (Theorem <ref>) and Real Nullstellensatz (Theorem <ref>). We provide with a sound and semi-complete synthesis algorithm for polynomial programs.

Outline of the paper: In Section <ref>, we defined the notions of programs, semantics and templates. We also give a mathematically precise formulation of the problem in terms of constraint pairs. In Section <ref>, we presented classical theorems from Real Algebraic Geometry required for design of the synthesis algorithm. We also presented Theorem <ref> which is a new mathematical contribution from this paper. In Section <ref>, we presented our synthesis algorithm along with proving it's soundness and semi-completeness. Finally in Section <ref> we presented a proof of our proposed concept by providing experimentation results on small programs.

§ TEMPLATE-BASED SYNTHESIS OF POLYNOMIAL PROGRAMS

In this work, we focus on imperative polynomial programs. Our programs can have real variables, conditionals and loops. Moreover, they include specifications as pre and post-conditions, as well as invariants for each loop, which are used to prove that the program satisfies the specification.

Variables and Valuations We assume that our programs have a finite set $\vars = \{v_1, \ldots, v_k\}$ of real-valued variables. A valuation is a function $\valuation: \vars \rightarrow \mathbb R$ that assigns a real value to each program variable.

Syntax We consider programs that can be generated from the following grammar:

\begin{align*}
P &:= (\Phi, C,\Phi) \\
C &:= \textbf{skip} \;\mid\; v \leftarrow \Pi \mid \texttt{ if } (\Phi)\; \{C\} \texttt{ else }\{ C\}  \mid\texttt{while } (\Phi, \Phi) \; \{C\} \mid  C; C\\
\Phi &:= \Pi \ge 0 \mid (\Phi \land \Phi) \mid \lnot \Phi\\
v & \in \vars \\
\Pi & \in \mathbb{R}[\vars]
\end{align*}

Intuitively, a program $P$ consists of a pre-condition, a sequence $C$ of commands, and a post-condition. We have a special command skip that does not do anything. We also allow assignments, conditionals and while loops. The header of each loop should include both a guard and a loop invariant. Moreover, the right-hand-side of every assignment is a polynomial expression $\Pi$ over program variables, and the guards, invariants and pre and post-conditions are boolean combinations $\Phi$ of polynomial inequalities over $\vars.$ For simplicity, we often write the programs as in Figure <ref> to make them more human-readable. We also use standard syntactic sugars, e.g. to allow other boolean operators.

Polynomial Assertions A polynomial assertion over a set $V$ of variables is a boolean combination of polynomial inequalities of the form $p \geq 0$ in which $p \in \mathbb{R}[V].$ In the sequel, without loss of generality, our algorithms assume that all polynomial assertions are in disjunctive normal form.

Polynomial Transition Systems We define the semantics of our programs using transition systems. A polynomial transition system (PTS) is a tuple $ (\vars,\locs,\loc_0,\theta_0,\loc_f, \theta_f, \transitions, \invariant)$ that consists of:

* A finite set $\vars$ of program variables,

* A finite set $\locs$ of locations,

* An initial location $\loc_0 \in \locs$,

* An initial polynomial assertion $\theta_0$,

* A final location $l_f \in \locs$,

* A final polynomial assertion $\theta_f,$

* A finite set $\transitions$ of transitions. Each transition $\transition \in \transitions$ is a tuple $(\loc,\loc',\rho_{\transition})$, where $\loc,\loc'\in \locs$ are the pre and post locations, and $\rho_{\transition}$, called the transition relation, is a polynomial assertion over $\vars \cup \vars'$, where $\vars$ represents variables in the pre location and its primed version $\vars'$ represents the variables in the post location, i.e. after taking the transition. Specifically, we have $\vars' = \{v'~\mid~ v \in \vars \};$

* An invariant $\invariant$ that maps some of the locations $\loc \in \locs$ to a polynomial assertion $\invariant(\loc).$ Moreover, we have $\invariant(\loc_0) = \theta_0$ and $\invariant(\loc_f) = \theta_f.$ The set of locations for which an invariant is provided must be a cutset (see further below for a formal definition).

The translation from programs to transition systems is intuitive and straightforward. We can have one location $\ell \in \locs$ for each line of the program, and let $\ell_0$ correspond to the initial line and $\ell_f$ to the end of the program. Similarly, $\theta_0$ serves as the pre-condition and $\theta_f$ as the post-condition. Programs are more human-friendly and hence a nicer way to specify the input to the synthesis problem. In contrast, it is easier to define the algorithms using the more formal notion of PTS. As such, we will define our semantics based on a PTS instead of a program. Below, we fix a PTS $(\vars,\locs,\loc_0,\theta_0,\loc_f, \theta_f, \transitions, \invariant).$

The Polynomial Transition System (PTS) Corresponding to the Program in Figure <ref>.

Figure <ref> shows the PTS modeling the program of Figure <ref>. We have $\locs = \{3, 4, 6, 8, 9, 11\}$ and each location corresponds to a line of code. The invariants are bracketed and shown in green. Note that our algorithm works with polynomial inequalities, so $a = b$ is syntactic sugar for $a - b \geq 0 \wedge b - a \geq 0.$

Control Flow Graphs (CFGs) The control flow graph (CFG) of our PTS is a directed graph $G = (\locs, E)$ in which the locations serve as vertices and there is a directed edge from $\loc$ to $\loc'$ iff there is a transition $\transition \in \transitions$ of the form $\transition = (\loc, \loc', \rho_\transition).$

Cutsets, Cutpoints and Basic Paths A subset $\cutset \subseteq \locs$ of vertices of the CFG $G$ is called a cutset if

* $\loc_0, \loc_f \in \cutset$; and

* Every cycle in $G,$ be it a simple cycle or not, intersects $\cutset.$

Assuming we have fixed a cutset $\cutset$, every location $c \in \cutset$ is called a cutpoint. A basic path is a finite sequence of transitions $\pi = \langle \transition_i \rangle_{i=1}^m = \langle (\loc_i, \loc'_i, \rho_{\transition_i})_{i=1}^m\rangle$ such that:

* $\loc_i = \loc'_{i-1}$ for all $i>1$.

* $\loc_1$ and $\loc'_m$ are cutpoints, but no internal location $\loc_i$ with $1 < i \leq m$ is a cutpoint.

We define the transition relation $\rho_\pi$ of the basic path $\pi$ as $\rho_{\tau_m} \circ \rho_{\tau_{m-1}} \circ \cdots \circ \rho_{\tau_1}.$ Informally, $\rho_\pi$ is the transition relation resulting from going through all the transitions of $\pi$ one-by-one and is hence equal to the sequential composition of $\rho_{\tau_i}$'s. If $\pi$ is the empty path, then we let $\rho_{\pi}$ be the identity relation $\bigwedge_{v \in \vars} v = v'.$

In the PTS of Figure <ref>, the cutpoints shown in red form a cutset $\cutset = \{3, 6, 11\}.$ In this example, our non-empty basic paths are $\langle \transition_1, \transition_2 \rangle, \langle \transition_4 \rangle,$ and $\langle \transition_3, \transition_5, \transition_6 \rangle.$

Remark on Invariants It is straightforward to see that it suffices to specify invariants in a cutset only and that the invariants can then be inductively expanded to every other line <cit.>. In the sequel, we assume that a cutset $\cutset$ is fixed and that the invariant $\invariant$ assigns a polynomial assertion $\invariant(\loc)$ to every cutpoint $\loc \in \cutset.$ This also fits our simple syntax above in which the invariant is only provided for the while loops and the start and end of the program (in the form of pre and post conditions). Nevertheless, our algorithmic approaches consider the PTS and can hence work with any arbitrary cutset $\cutset.$

States A state of the PTS is a pair $\state = (\loc, \valuation)$ in which $\loc \in \locs$ is a location and $\valuation: \vars \rightarrow \mathbb R$ is a valuation. Intuitively, a state specifies which line of the program we are in and what values are taken by our program variables.

Runs A run of our PTS is a sequence $R = \langle \state_i \rangle_{i=0}^m = \langle (\loc_i, \valuation_i) \rangle_{i=0}^m$ of states such that:

* $\val_0 \models \theta_0,$ i.e. the run starts with the initial location $\loc_0$ and an initial valuation that satisfies the precondition $\theta_0;$

* For every $i < m,$ there exists a transition $\transition_i \in \transitions$ that allows us to go from $\sigma_i$ to $\sigma_{i+1}.$ More formally:

* $\transition_i = (\loc_i, \loc_{i+1}, \rho_{\transition_i});$ and

* $(\val_i, \val_{i+1}') \models \rho_{\transition_i},$ i.e. if we consider $\val_i$ as a valuation on $\vars$ and $\val_{i+1}$ as a valuation on $\vars',$ then they jointly satisfy the transition condition $\rho_{\transition_i}.$

Valid Runs A run $R = \langle \state_i \rangle_{i=0}^m = \langle (\loc_i, \valuation_i) \rangle_{i=0}^m$ is called valid iff for every $i$ for which $\loc_i$ has an associated invariant $\invariant(\loc_i),$ we have $\valuation_i \models \invariant(\loc_i).$ Intuitively, a run is valid if it always satisfies the invariant. Specifically, if $\loc_m = \loc_f,$ then we must have $\val_m \models \theta_f$ in order for the run to be valid.

Valid Transition Systems A transition system (PTS) is called valid iff every run of the PTS is valid.

Inductively Valid Transition System A PTS is called inductively valid with respect to a cutset $\cutset$, iff it satisfies the following conditions for all cutpoints $\loc, \loc' \in \cutset$ and valuations $(\val, \val')$ over $\vars \cup \vars'$:

* Initiation: For every basic path $\pi$ from $\loc_0$ to $\loc,$ we have
$$\textstyle \val \models \theta_0 \wedge (\val, \val') \models \rho_\pi \Rightarrow \val' \models \invariant(\loc)'.$$
Here, $\invariant(\loc)'$ is the same as $\invariant(\loc),$ except that every program variable $v \in \vars$ is replaced with its primed version $v'.$

* Consecution: For every basic path $\pi$ from $\loc$ to $\loc',$ we have
$$\textstyle \val \models \invariant(\loc) \wedge (\val, \val') \models \rho_\pi \Rightarrow \val' \models \invariant(\loc')'.$$

* Finalization: For every basic path $\pi$ from $\loc$ to $\loc_f,$ we have
$$\textstyle \val \models \invariant(\loc) \wedge (\val, \val') \models \rho_\pi \Rightarrow \val' \models \theta'_f.$$

The intuition behind inductive validity is similar to that of inductive invariants <cit.>. Informally, a PTS is inductively valid if we can prove that it is valid by breaking each run into a sequence of basic paths and then show its validity by an inductive argument. It is straightforward to see that every inductively valid PTS is indeed valid <cit.>. Moreover, the finalization and initiation conditions above can be considered as special cases of consecution. We are including them separately to emphasize the point that validity requires the PTS to satisfy its specification (as specified by pre and post-conditions).

Non-determinism Note that, in order to simplify the presentation, we did not include non-determinism in the syntax of our programs. However, our transition systems have built-in support for non-determinism since one can have several outgoing transitions from the same node or, alternatively, a transition that does not uniquely specify the post-valuation $\valuation'$ in terms of the pre-valuation $\valuation.$ Therefore, all the results in the sequel are trivially extensible to the case of non-deterministic polynomial programs.

Our goal is to synthesize a valid PTS, given a sketch that includes holes in it. To do this, we first define variants of our programs and transition systems which can also include holes.

Sketches of Polynomial Programs
In this work, we assume that the programmer provides a sketch of the implementation. The sketch is a partial program, i.e. a program with holes, which is generated using the following grammar:

\begin{align*}
P &:= (\Phi, C,\Phi) \\
C &:= \textbf{skip} \;\mid\; v \leftarrow \Pi \mid v \leftarrow  \textcolor{red}{H} \mid \texttt{ if } (\Phi)\; \{C\} \texttt{ else }\{ C\}  \mid\texttt{while } (\Phi, \Phi) \; \{C\} \mid  C; C\\
\Phi &:= \Pi \ge 0 \mid \textcolor{red}{H } \ge 0 \mid (\Phi \land \Phi) \mid \lnot \Phi\\
\textcolor{red}{H} & := (\overline{\vars}, d)\\
v & \in \vars \\
\Pi & \in \mathbb{R}[\vars]\\
\overline{\vars} & \subseteq \vars\\
d & \in \mathbb N
\end{align*}

The main difference between this grammar and the previous one is the introduction of a new non-terminal $H$ which models a hole. Note that holes can appear in assignments or as part of the assertions $\Phi,$ be it in guards of if/while or the invariants or even pre and post-conditions. For each hole, the programmer specifies a set $\overline{\vars}$ of variables which she expects might need to appear in filling the hole, as well as a maximum degree $d \in \mathbb N,$ requiring that the hole should be filled by a polynomial of degree at most $d$. As an example, see the sketch in Figure <ref> in which the holes are shown in boxes. As usual, we do not follow the grammar exactly in order to make the sketch more human-readable. For example, we wrote the invariant outside the while loop, instead of having it as the second parameter in the while header.

Template and Symbolic Polynomials and Assertions Let $\vars$ be our finite set of program variables and $\tvars$ a separate finite set of template variables. We define a symbolic polynomial over $(\vars, \tvars)$ as a polynomial over $\vars$ whose every coefficient is itself a polynomial over $\tvars.$ More formally, every symbolic polynomial over $(\vars, \tvars)$ is of the form $p \in (\mathbb{R}[\tvars])[\vars].$ We say a symbolic polynomial $p$ is a template if every coefficient of $p$, when considered as a polynomial over $\tvars,$ has degree $0$ or $1$. Symbolic (resp. template) polynomial assertions are simply defined as boolean combinations of symbolic (resp. template) polynomial inequalities. As mentioned before, our algorithms assume without loss of generality that all assertions are in disjunctive normal form.

Sketches of Polynomial Transition Systems A sketch of a polynomial transition system (SPTS) is tuple $(\vars, \tvars, \locs,\loc_0,\theta_0,\loc_f, \theta_f, \transitions, \invariant)$ in which all parts serve the same purposes and are defined similarly as in a PTS. The only differences are:

* The assertions $\theta_0, \theta_f$ and invariants $\invariant(\loc)$ are now template polynomials over $(\vars, \tvars).$

* The transition relations $\rho_\transition$ are now template polynomials over $(\vars \cup \vars', \tvars).$

Moreover, we always assume $\tvars \cap \vars = \emptyset.$

The SPTS Corresponding to the Sketch Program of Figure <ref>. For brevity, we have not included equalities of the form $v' = v$ that ensure a variable $v$ does not change in a transition. See Figure <ref> for details of these equalities. The cutpoints are shown in red.

Figure <ref> shows a sketch polynomial transition system (SPTS) corresponding to the sketch program of Figure <ref>. To obtain this SPTS, we introduced 21 new template variables $\tvars = \{t_0, \ldots, t_{20}\}$ and replaced each hole with a template polynomial that corresponds to the programmer's desire for that hole. Note that the invariant of location $6$ and the transition relations for $\tau_3, \tau_4$ and $\tau_6$ are now symbolic, i.e. they are polynomials over the program variables $\{i, s, n\}$ in which every coefficient is itself a polynomial over $\tvars.$ In this case, all coefficients are simply variables in $\tvars.$ Our goal is to find suitable concrete real values for each $t_i$ such that the SPTS of Figure <ref> becomes a valid PTS, e.g. that of Figure <ref>.

We now have all the needed ingredients to formalize our Template-Based Synthesis problem over Polynomial programs (TBSP):

TBSP Given a sketch of a polynomial transition system $\Gamma = (\vars, \tvars, \locs,\loc_0,\theta_0,\loc_f, \theta_f, \transitions, \invariant)$, find a valuation $\tval: \tvars \rightarrow \mathbb R$ that assigns a real value to every template variable, such that $\Gamma^{\tval}$ becomes an inductively valid polynomial transition system. Here, $\Gamma^{\tval}$ is the PTS obtained by taking the sketch $\Gamma$ and replacing every appearance of each template variable $t \in \tvars$ by its corresponding value $\tval(t).$ Informally, our goal is to synthesize concrete values for each of the template variables so that the resulting program/PTS satisfies the programmer's desired specification.

Boundedness For technical reasons which will become evident in Section <ref>, we assume that our program variables $\vars$ have bounded values. More specifically, we consider a bound $M \in (0, \infty)$ and assume that no program variable ever takes a value larger than $M$ or smaller than $-M.$ We call this the boundedness assumption. The boundedness assumption is necessary for our completeness theorem (Theorem <ref>), but not for soundness. It is a merely theoretical assumption with no significant practical effect. In practice, we do not assume boundedness in our implementation in Section <ref> and can successfully handle all the cases without this extra assumption.

§ TEMPLATE-BASED SYNTHESIS OF POLYNOMIAL PROGRAMS

We consider program as a sequence of instructions, as in an imperative programming language. For instance, consider the following example program $P$:

@real: i,s,n;
@pre: n>=0;
i = 0;
while(i <= n, s >= i)
s = s+1;
i = i+1;

@post: s >= n;

The program variables in the above example program are $V = \{s,i,n\}$, and we say that the program is defined over the variables $V$. Note that we have the precondition and postcondition in our example program, which are polynomial assertion over the program variables. A program is said to be valid if every valuation of program variables that satisfies the pre condition, the valuation obtained after executing the program satisfies the post condition. The programmer sometimes may have a specification of a program via precondition and postcondition with partial implementation. We call such partial implementation as program template. The synthesis problem is to complete such partial implementation and generate a concrete valid program.

In this paper, we restrict ourselves to the so-called polynomial programs over a set of variables $V$ that can be generated from the following grammar:

\begin{align*}
P &:= (\Phi, E,\Phi) \\
E &:= v \leftarrow P(V) \mid \texttt{ if } (\Phi)\; \{E\} \texttt{ else }\{ E\}  \mid\texttt{while } (\Phi, \Phi) \; \{E\} \mid  E; E\\
\Phi &:= P(V) \ge 0 \mid (\Phi \land \Phi) \mid \lnot \Phi
\end{align*}

where $P(V)$ denotes any polynomial expression over the variables $V$ with real coefficients. The semantics of each instruction are as usual except for $while(\psi, \phi, E)$ which is a while loop with guard $\psi$ for the block of code $E$, and $\phi$ is a loop invariant. In general, the loop invariant $\phi$ should be strong enough for us to prove the correctness of the program with given precondition and postcondition.

Let the set of program variables be $V = \{x_1, x_2, \ldots, x_n\}$, $\nu \in \mathbb{R}^V$ be any valuation of $V$. Given a polynomial $P$ over $V$ and a valuation $\nu \in \mathbb{R}^V$, we denote $P(\nu)$ to be the evaluation of the polynomial $P$ on substituting each $x \in V$ with $\nu(x)$. On executing an instruction $E$, the valuation of the variables becomes $\nu' = E(\nu)$. We define the semantics of the program recursively in the Table <ref>.

Expression $E$     Final Valuation $\nu'$

$E_1;E_2$      $\nu' = E_2(E_1(\nu))$

$v \leftarrow  P(V)$     $\nu'(v) = P(\nu)$ and $\forall x \in V \setminus \{v\}. \nu'(x) = \nu(x)$

$\texttt{if } (\phi)\; \{E_1\} \texttt{ else }\{ E_2\}$     if $\nu \models \phi$ then $\nu' = E_1(\nu)$ else $\nu' = E_2(\nu)$

$\texttt{while } (guard, inv) \; \{E\}$     if $\nu \models \textit{guard}$ and $\nu \models \textit{inv}$, then $\nu' = E(\nu_1)$ where $\nu_1 = E(\nu)$ and $\nu_1 \models \textit{inv}$ else $\nu' = \nu$

Consider the following template $T$, for example, whose one realization is the program $P$ from Example <ref>. The synthesis algorithm will replace the holes with suitable polynomials, if they exist, such that the resulting program satisfies the constraints given by pre-conditions and post-conditions.
@real: i,s,n;
@function f1 = [(i),1]
@function f2 = [(i),1]
@pre: n>=0;
i = 0;
s = 0;
while(i <= n, s >= f1)
s = s+f2;
i = i+1;

@post: s >= n;

A template program basically specifies the intuition behind the final program and has holes at location where the program expects the synthesizer to fill in appropriate expressions. A hole could be thought of as a hint that programmer has provided to synthesis a correct program. In this example template program, we have two polynomial function symbols $F = \{f_1, f_2\}$. $[V',d]$, called a program hole, represents a symbolic polynomial $\textbf{P}$ over the variables  $V'$ with total degree $d \in \mathbb{N}$ and coefficients from the template variable set $U = \{a_1, a_2, a_3, a_4\}$. Given a template program, our goal is to synthesize values for all the template variables such that the resulting program is valid.

In the following we formalize what we mean by partial implementation, or template. A template of the polynomial program over the variables $V$ and set of polynomial function symbol $F$ is defined using the following grammar:

\begin{align*}\label{template-grammar1}
E &:= E' \mid \texttt{ if } (\Phi)\; \{E'\} \texttt{ else }\{ E\}  \mid\texttt{while } (\Phi, \Phi) \; \{E'\} \mid  E; E\\
X &:= f \leftarrow P[V] \mid f \leftarrow [V',d] \mid f \leftarrow f + f \mid  f \leftarrow f \times f \mid [\sigma](f) \\
\sigma &:= v \to P[V] \mid \sigma,\sigma\\
\Phi &:= P(V) \ge 0 \mid [V', d] \ge 0 \mid (\Phi \land \Phi) \mid \lnot \Phi\\
\end{align*}

The polynomial function symbols could be used by the programmer to specify either a concrete polynomial or a symbolic polynomial. The programmar can also evaluate a polynomial function symbol using a substitution map $\sigma: \mathbb{R}[V] \to \mathbb{R}[V]$. A substitution map that can be applied on polynomials and function symbols to obtain another polynomial expression.

Though this representation is programmer-friendly, to describe our algorithm succinctly we shall use an equivalent representation of programs and template programs. We provide with the notions of polynomial transition systems and symbolic polynomial transition systems which are respectively equivalent to polynomial programs and symbolic polynomial programs.

In the figure <ref>, we graphically describe the symbolic polynomial transition system $TS$ corresponding to the template program $T$ defined in Example <ref>. The location $l_0$ and $l_4$ are the inital and final location respectively. The set $C = \{l_2\}$ is the cutset.

The transition system for the template program defined in <ref>

We will first define few terms that are useful, and then define polynomial variants of linear transitions system and linear control-flow graph that were defined in <cit.>.

Given a set $V$ of program variables and $U$ of template variables, we define $\Psi(U, V)$ to be the set of finite boolean combination of polynomial inequalities over program variables with coefficients being real polynomials over template variables. More precisely, $\Psi(U,V)$ is the set of all possible finite boolean combinations of inequalities of the form $f \geq 0$ such that $f \in (\mathbb{R}[U])[V]$.

We can now define a polynomial transition system.

A polynomial transition system is a tuple $P = (V,L,l_0,\theta_0,\mathcal{T}, l_f, \theta_f)$ that consists of a set of variables $V$, a set of locations $L$, an initial location $l_0$, an initial assertion $\theta_0$, a final location $l_f$, a final assertion $\theta_f$ and set of transitions $\mathcal{T}$. Each transition $\tau \in \mathcal{T}$ is a tuple $(l,l',\rho_{\tau})$, where $l,l'\in L$ are the pre and post locations, and $\rho_{\tau}$ is an polynomial assertion over $V\cup V'$, where $V$ represents current-state variables and its primed version $V'$ represents the next-state variables.

Symbolic polynomial transition system over the program variables $V$ and has the same structure as the polynomial transition system except that its assertions are allowed to be symbolic polynomial assertion from $\Psi(U,V \cup V')$.

A symbolic polynomial transition system is a tuple $T = (V\cup U,L,l_0,\theta_0,\mathcal{T}, l_f, \theta_f)$ that consists of a set of variables $V$, a set of template variables $U$, a set of locations $L$, an initial location $l_0$, an initial symbolic polynomial assertion $\theta_0$, a final location $l_f$, a final symbolic polynomial assertion $\theta_f$ and set of transitions $\mathcal{T}$. Each transition $\tau \in \mathcal{T}$ is a tuple $(l,l',\rho_{\tau})$, where $l,l'\in L$ are the pre and post locations, and $\rho_{\tau}$ is a symbolic polynomial assertion over $V\cup V'$, where $V$ represents current-state variables and its primed version $V'$ represents the next-state variables. All the coefficients ofsymbolic polynomial assertions over $V\cup V'$ belong to the polynomial ring $\mathbb{R}[U]$.

A control-flow graph corresponds to the underlying structure of a given transition system.

Given a transition system $P = (V,L,l_0,\theta_0,\mathcal{T}, l_f, \theta_f)$, we define control-flow graph $\emph{CFG}(P)$ with locations $L$ as the vertex set and transitions $\mathcal{T}$ as the edge set. More precisely we have an edge from $l$ to $l'$ if and only if we have a transition $\tau = (l,l',\rho_{\tau}) \in \mathcal{T}$ for some $\rho_{\tau}$.

We shall now define the notion of cutset and basic path which would be help us in defining the inductive polynomial assertion map.

Given a transition system $P = (V,L,l_0,\theta_0,\mathcal{T}, l_f, \theta_f)$, a subset of vertices of $\emph{CFG}(P)$ is called a cutset $C$ if every cyclic path in $\emph{CFG}(P)$ passes through some vertex in $C$. An element of $C$ is called a cutpoint.

Given a transition system $P = (V,L,l_0,\theta_0,\mathcal{T}, l_f, \theta_f)$, a path $\pi$ between two cutpoints $l$ and $l'$ is called a basic path if it does not pass through any other cutpoint in $C$.

Given a transition system $P = (V,L,l_0,\theta_0,\mathcal{T}, l_f, \theta_f)$ with a cutset $C$ and an assertion $\eta_C(l)$ for each cutpoint $l$, we say that $\eta_C$ is an inductive polynomial assertion map for $C$ if it satisfies the following conditions for all cutpoints $l,l' \in C$:

* Initiation: For each basic path $\pi$ from $l_0$ to $l$, $\theta_0 \wedge \rho_{\pi} \models \eta_C(l)'$.

* Consecution: For each basic path $\pi$ from $l$ to $l'$, $\eta_C(l) \wedge \rho_{\pi} \models \eta_C(l')'$.

* Finalization: For each basic path $\pi$ from $l$ to $l_f$, $\eta_C(l) \wedge \rho_{\pi} \models \theta_f $

A inductive polynomial assertion map is said to be symbolic if it has one or more symbolic polynomial assertions. The definitions of control-flow graph, cutset and basic path also naturally extends to the symbolic polynomial transition system. Given a model $M: U \to \mathbb{R}$ of the template variables, let $M(T)$ and $M(\eta_C)$ denote the concrete polynomial transtion system and concrete inductive polynomial assertion map obtained by substituting the values of the template variables in the respective templates.

§.§ Problem Statement

Let $T$ be a symbolic polynomial transition system over the program variables $V = \{x_1,\dots,x_n\}$ and the set of template variables $U = \{c_1, \ldots, c_m\}$, $C$ be a cutset of $T$,  and $\eta_C$ be the symbolic polynomial inductive assertion map.

The synthesis problem for the tuple $(T, \eta_C)$ is to find a valuation $M: U \to \mathbb{R}$ of template variables such that $M(\eta_C)$ is valid inductive polynomial assertion map for the concrete transtion system $T$.

§ HARDNESS

?? TODO: FILL THIS IN

§ ALGORITHMS

§.§ Overview of our approach

The programmer first writes a template in the grammar defined at in the beginning of the section. In the pre-processing step, the algorithm first verifies the syntax of the program and substitute the program holes with symbolic polynomial of the specified degree. We call the program generated after the pre-processing step as symbolic polynomial program.

A constraint pair is a pair $(\alpha, \beta)$ where $\alpha = \{g_i \bowtie_i 0 \mid \bowtie_i \in \{\ge, >\}\}$ is a set of polynomial inequalities, and $\beta$ is an polynomial inequality of the form $f \bowtie 0 $ with $\bowtie \in \{\ge, > \}$. We can reduce the correctness of the symbolic polynomial program $T$ to validity of a set of constraint pairs $S$.

We will describe our approach to construct the constraint pairs with the following example.

We are given the symbolic polynomial transition system $TS$ as described in Example <ref>, the cutset $C = \{l_2\}$ and the symbolic polynomial assertion map $\eta_C$. Say, $\eta_C(l_2) = a_0s + a_1i \ge 0$. We will create three condition pairs corresponding to various restrictions on $\eta_C$ for it to be a valid inductive assertion map.

* Initiation: There is only one basic path from the initial location $l_0$ to $l_2$. We get the following constraint pair for the path $\pi_1 := l_0 \to l_1 \to l_2$:

$$S_1 := ([n_0 >= 0 \land i_0 == 0 \land s_0 == 0],[a_0s_0 + a_1i_0 \ge 0 ])$$

* Consecution: There is one basic path $\pi_2 := l_2 \to l_3 \to l_3$ from $l_2$ to $l_2$. The constraint pair for the basic path $\pi_2$ is:

$$S_2 := ([i_0 <= n_0 \land s_1 == s_0 + (a_0i_0) + a_1 \land i_1 == i_0 + 1], [a_0s_1 + a_1i_1 \ge 0])$$

* Finalization:
The path $\pi_3 := l_2 \to l_4$ is the only basic path from points of cutset $C$ to the final location $l_4$. For this path, we get the following constraint pair:

$$S_3 := ([a_0s_1 + a_1i_1 \ge 0, i_1 > n_0], [s_1 > n_0])$$

For each basic path $\pi$ from $l$ to $l_f$, $\eta_C(l) \wedge \rho_{\pi} \models \theta_f $

We get the three constraint pairs $\{S_1, S_2, S_3\}$.

The number of constraint pairs generated using our approach is at most quadratic in the number of locations if the given cutset contains all the locations, i.e, $C = L$. Moreover, if the cutset $C \subset L$, then the number of constraint pairs is linear in the to the number of basic paths between each pair of locations in the cutsets.

The proof follows directly from the construction.

§ OUR ALGORITHM FOR TBSP

In this section, we provide an automated, sound and semi-complete algorithm for the TBSP problem as defined in Section <ref>. We assume that the input to our algorithm is either a sketch polynomial transition system (SPTS) or a sketch of a polynomial program, containing holes. In the latter case, we would first process the input and translate it to an SPTS. We start by providing a high-level overview of our algorithm (Section <ref>) and an illustration over the example of Figure <ref>. Then, we present our mathematical toolkit, i.e. Farkas' Lemma, Handelman's Theorem, Putinar's Positivstellesatz and the Real Nullstellensatz in Section <ref> and manipulate them to be usable for our application. We use this toolkit to fill in the details of our algorithm in Section <ref>. Finally, we provide soundness and semi-completeness theorems in Section <ref>.

§.§ An Illustration and Overview of the Algorithm

Suppose that a sketch, i.e. a polynomial program with holes, is given as the input. As a running example, let us take Figure <ref>. Our algorithm fills in the holes using the following four steps:

Step 1: Creating Templates and Translation to SPTS For each hole in the sketch, the algorithm generates a template polynomial in which it introduces new template variables. Concretely, let the hole be of the form $H = (\overline{\vars}, d).$ This corresponds to a case in which the programmer desires the hole to be filled by a polynomial of degree at most $d$ over the set $\overline{\vars} \subseteq \vars.$ The algorithm first generates all monomials $m_i$ of degree at most $d$ over $\overline{\vars}.$ It then creates a template polynomial that includes all of these monomials and is of the form $p_H := \sum_{i} t_i \cdot m_i,$ in which each $t_i$ is a newly generated template variable. It replaces the hole $H$ by the template polynomial $p_H$ and then translates the program to an SPTS. The translation from programs to transition systems is a standard procedure, i.e. a simple parsing exercise, and hence we skip it for brevity. At the end of this step, we have an SPTS $\Gamma = (\vars, \tvars, \locs,\loc_0,\theta_0,\loc_f, \theta_f, \transitions, \invariant).$

Consider the sketch of Figure <ref>. The first hole in the invariant of line 5 should be filled by a polynomial of degree at most $2$ over the variables $\{i, s\}.$ Hence, the algorithm generates the following template for this hole:

\textstyle p_1 := t_0 + t_1 \cdot i + t_2 \cdot s + t_3 \cdot i^2 + t_4 \cdot i \cdot s + t_5 \cdot s^2.

Note that all the $t_i$'s are new template variables whose concrete values have to be synthesized later. It then establishes a similar template for the second hole in line 5:

\textstyle p_2 := t_6 + t_7 \cdot i + t_8 \cdot s + t_9 \cdot i^2 + t_{10} \cdot i \cdot s + t_{11} \cdot s^2.

The third hole in line 5 and the hole in line 6 should be affine expressions over $\{i, n\},$ so the algorithm generates:

\begin{matrix}
p_3 := t_{12} + t_{13} \cdot i + t_{14} \cdot n,\\
p_4 := t_{15} + t_{16} \cdot i + t_{17} \cdot n.
\end{matrix}

Finally, the last hole (line 9) should be filled by an affine expression over $\{i, s\}.$ So, we have:

\textstyle p_5 := t_{18} + t_{19} \cdot i + t_{20} \cdot s.

The algorithm first replaces each $i$-th hole with its corresponding template polynomial $p_i.$ It then translates the program to a transition system, leading to the SPTS of Figure <ref>.

Step 2: Generating Entailment Constraints Recall that our goal is to synthesize concrete real values for the template variables $t_i \in \tvars$ such that our SPTS becomes an inductively valid PTS. As such, it should satisfy the initiation, consecution and finalization constraints defined in Section <ref>. In this step, the algorithm considers every possible basic path $\pi$ and symbolically computes its transition relation and its corresponding initiation/consecution/finalization constraint. Then, the algorithm translates the constraint into an equivalent system of constraints in the following standard form:

\begin{equation} \label{eq:entailment-form}
f_1 \bowtie_1 0 ~\wedge~ f_2 \bowtie_2 0 ~\wedge~ \ldots \wedge f_r \bowtie_r 0 \Rightarrow f \bowtie 0
\end{equation}

in which $f$ and each $f_i$ are symbolic polynomials over $(\vars \cup \vars', \tvars)$ and each operator $\bowtie, \bowtie_i$ is either $\ge$ or $>$.

Consider the SPTS and cutset of Figure <ref>. As discussed in Example <ref>, the basic paths in this SPTS are
$\pi_1 = \langle \transition_1, \transition_2 \rangle, \pi_2 = \langle \transition_4 \rangle,$ and $\pi_3 = \langle \transition_3, \transition_5, \transition_6 \rangle.$
The algorithm symbolically computes the transition relation of each basic path:

\rho_{\pi_1} = \rho_{\tau_2} \circ \rho_{\tau_1} = \left[ n'=n ~\wedge~ i'=0 ~\wedge~ s'=0 \right]

\rho_{\pi_2} = \rho_{\tau_4} = \left[ t_{15} + t_{16} \cdot i + t_{17} \cdot n < 0 ~\wedge~ i'=i ~\wedge~ s'=s ~\wedge~ n'=n \right]

\rho_{\pi_3} = \rho_{\tau_6} \circ \rho_{\tau_5} \circ \rho_{\tau_3} = \left[ t_{15} + t_{16} \cdot i + t_{17} \cdot n \geq 0 ~\wedge~ i' = i+1 ~\wedge~ s'=t_{18} + t_{19} \cdot (i+1) + t_{20} \cdot s ~\wedge~ n' = n \right]

It then has to symbolically compute the inductive validity constraints. For $\pi_1$, we have the following initiation constraint:

$$\val \models \theta_0 \wedge (\val, \val') \models \rho_{\pi_1} \Rightarrow \val' \models \invariant(6)'.$$

The algorithm symbolically computes and expands this constraint and obtains:

n \geq 1 ~\wedge~ n' = n ~\wedge~ i' = 0 ~\wedge~ s' = 0 \Rightarrow

t_0 + t_1 \cdot i' + t_2 \cdot s' + t_3 \cdot i'^2 + t_4 \cdot i' \cdot s' + t_5 \cdot s'^2 \ge 0 ~\wedge~

t_6 + t_7 \cdot i' + t_8 \cdot s' + t_9 \cdot i'^2 + t_{10} \cdot i' \cdot s' + t_{11} \cdot s'^2 \geq 0 ~\wedge~

t_{12} + t_{13} \cdot i' + t_{14} \cdot n' \geq 0

Note that this constraint should hold for any valuation of $\vars \cup \vars'.$ At this point, since we know $n' = n$ and $i'=s'=0,$ we can further simplify the constraint to

n \geq 1 \Rightarrow t_0 \geq 0 ~\wedge~ t_6 \geq 0 ~\wedge~ t_{12} + t_{14} \cdot \textcolor{red}{n} \geq 0,

hence eliminating several variables and obtaining a much simpler but equivalent constraint. Our implementation in Section <ref> applies such sound heuristics for simplifying the constraints. However, the constraint is still not in the standard form of (<ref>) since there are several inequalities on the right-hand-side. The algorithm rewrites this constraint as the following equivalent set of constraints of the standard form (<ref>):

n \geq 1 \Rightarrow t_0 \geq 0;

n \geq 1 \Rightarrow t_6 \geq 0;

n \geq 1 \Rightarrow t_{12} + t_{14} \cdot n \geq 0.

Let us now consider $\pi_2.$ We need to generate a finalization constraint in this case. The algorithm symbolically computes

\invariant(6) ~\wedge~ \rho_{\pi_2} \Rightarrow \theta'_f

and then expands it and rewrites it as an equivalent set of standard constraints just as in the previous case. Finally, for $\pi_3,$ which starts and ends at the cutpoint $6$, the algorithm writes a consecution constraint:

\invariant(6) ~\wedge~ \rho_{\pi_3} \Rightarrow \invariant(6)'.

As above, the constraint will be simplified and translated to an equivalent set of standard constraints.

Step 3: Eliminating Program Variables and Reduction to QP We start by providing some intuition for this step. Let $\vars = \{v_1, \ldots, v_m\}$ be our program variables and $\tvars = \{t_1, \ldots, t_{m'}\}$ the template variables introduced in Step 1 above. By the end of the previous step, we have a set of standard constraints of form (<ref>). Let the constraints be $\kappa_1, \kappa_2, \ldots, \kappa_u.$ Our problem is now equivalent to finding concrete real values $t^*_i \in \mathbb R$ for each template variable $t_i \in \tvars$ such that the following formula is satisfied:

\begin{equation} \label{eq:firstorder}
\forall v_1, v_2, \ldots, v_m ~\;~\;~ (\kappa_1 ~\wedge~ \kappa_2 ~\wedge~ \ldots ~\wedge~ \kappa_u)
\end{equation}

or more formally

\begin{equation} \label{eq:firstorder22}
\textstyle \forall \val \in \mathbb{R}^{\vars} ~\;~\;~ \val \models \left( \bigwedge_{j=1}^u \kappa_j \left[ t_i \leftarrow t^*_i \right]_{i=1}^{m'} \right),
\end{equation}

where $\kappa_j \left[ t_i \leftarrow t^*_i \right]_{i=1}^{m'}$ is the formula $\kappa_j$ in which each $t_i$ is replaced by $t^*_i.$ This is not an easy problem in general (See Section <ref>), but the fact that all $\kappa_i$ constraints are in the standard form of (<ref>) can be exploited for a reduction to quadratic programming (QP). This step requires mathematical details that will be presented in Section <ref>. At this point, we consider a simple sound algorithm that captures the main intuition of this step and then fill in the details to make it semi-complete in Section <ref>.

Consider a standard constraint $\kappa_i$ of the following form:

f_1 \geq 0 ~\wedge~ f_2 \geq 0 ~\wedge~ \ldots \wedge f_r \geq 0 \Rightarrow f \geq  0

As mentioned above, this constraint should hold over all valuations. So, intuitively, we would like to guarantee that whenever all polynomials $f_1, f_2, \ldots, f_r$ are non-negative, then so is $f$. A simple sound solution is to try to write $f$ as a combination of the $f_i$'s in the following form:

\begin{equation} \label{eq:lhsrhs}
\textstyle f = \lambda_0 + \sum_{i=1}^r  \lambda_i \cdot f_i
\end{equation}

where each $\lambda_i$ is a new non-negative real variable. We write $\Lambda$ as a shorthand for $\{\lambda_0, \ldots, \lambda_r\}.$ This is clearly sound, since whenever the $f_i$'s are non-negative at a given valuation, then so is every combination of them with non-negative coefficients. Indeed, we can go one step further and let each $\lambda_i$ be an always-non-negative polynomial, instead of a real constant, and the exact same argument will go through. As we will see, this idea can lead to a semi-complete algorithm.

We can now present Step 3 of the algorithm. In this step, the algorithm considers each standard constraint $\kappa_i$ and symbolically computes the equation (<ref>) for it. Note that both sides of (<ref>) are symbolic polynomials over $(\vars, \tvars \cup \Lambda).$ In particular, they are polynomials over $\vars$. However, two polynomials are equal if and only if they have equal coefficients for every monomial. Hence, the algorithm equates the corresponding coefficients on the left and right-hand sides of (<ref>) and rewrites it as a system of quadratic equalities over $\tvars \cup \Lambda.$

Consider the standard constraint
n - 1 \geq 0 \Rightarrow t_{12} + t_{14} \cdot n \geq 0
$ that was generated in Step 2 above (Example <ref>). The algorithm symbolically computes (<ref>) and obtains

t_{12} + t_{14} \cdot n = \lambda_0 + \lambda_1 \cdot (n-1).

Since both sides of this equality are polynomials over $n$, in order for the equality to hold, they should have the same constant factor and the same coefficient for $n$, so the algorithm obtains the following constraints:

t_{12} = \lambda_0 - \lambda_1

t_{14} = \lambda_1

\lambda_0, \lambda_1 \geq 0

This is a quadratic programming (QP) instance in which none of the program variables, e.g. $n$, appear. The algorithm translates other standard constraints to QP in a similar manner and combines them conjunctively.

Step 4: Constraint Solving The algorithm passes the QP instance generated in the previous step to an SMT solver or a numerical solver. If the solver fails to find a solution, then the algorithm fails, too. Otherwise, the solver provides values for all $t_i$ and $\lambda_i$ variables. The algorithm plugs the values found for $t_i$'s back into the templates of Step 1, which successfully leads to an inductively valid PTS or equivalenty to a valid synthesized program that satisfies the desired specification.

Given the sketch of Figure <ref> as input, the synthesized program of Figure <ref> is obtained as one of the solutions of the QP instance in Step 4. However, the solution is not necessarily unique, and our algorithm might find a different solution based on the values found by the external solver.

§.§ Mathematical Toolkit

In Section <ref> above, we provided an overview of our algorithm. However, the details of Step 3, i.e. the reduction to QP, were not presented since they depend on certain mathematical prerequisites. In this section, we provide the mathematical tools and theorems that are crucial for this step of the algorithm. We first recall some notation and classical definitions. Then, we present several theorems from polyhedral and real algebraic geometry. Finally, we obtain tailor-made versions of these theorems in a format that can be used in Step 3 of our algorithm above. We refer to <cit.> for a more detailed treatment of these theorems.

Sums of Squares
A polynomial $h \in \mathbb{R}[x_1,\dots,x_n]$ is a sum of squares (SOS) iff there exist $k \geq 1$ and polynomials $g_1, \dots, g_k \in \mathbb{R}[x_1, \dots, x_n]$ such that $h = \sum_{i=1}^k g_i^2.$

Strong Positivity
Given a set $X \subseteq \mathbb{R}^n$ and a polynomial $g \in \mathbb{R}[x_1,\dots,x_n]$, we say that $g$ is strongly positive over X if $\inf_{x\in X} g(x) > 0$. We write this as $X \models g \gg 0.$

Notation Let $\Phi = \{g_1 \bowtie_1 0,\dots, g_k \bowtie_k 0\}$ be a set of polynomial inequalities where $\bowtie_i \in \{\geq,>\}$ and $g_i \in \mathbb{R}[x_1, \dots, x_n]$. We define $\sat(\Phi)$ as the set of all real valuations $\val$ over the variables $\{x_1,\dots,x_n\}$ that satisfy $\Phi.$ More formally, $\sat(\Phi) = \{ \val \in \mathbb{R}^n ~\vert~ \bigwedge_{i=1}^k  \val \models (g_i \bowtie_i 0)\}.$

Given a set $X \subseteq \mathbb{R}^n$, we define $\overline{X}$ to be the closure of $X$ with respect to the Euclidean topology of $\mathbb{R}^n$. For a set $\Phi$ of polynomial inequalities, we define $\overline{\Phi}$ as the system of polynomial inequalities obtained from $\Phi$ by replacing every strict ineqaulity with its non-strict counterpart.

We are now ready to present the main mathematical theorems that will be used in our work. Our presentation follows that of <cit.>, which also contains proofs of corollaries that are not proven here.

Consider a set $V = \{x_1, \ldots, x_r\}$ of real-valued variables and the following system $\Phi$ of equations over $V$:

\begin{equation*}
\Phi :=
\begin{cases}
a_{1,0} + a_{1,1} \cdot x_1 + \ldots + a_{1,r} \cdot x_r \ge 0\\
\hfil \vdots\\
a_{m,0} + a_{m,1} \cdot x_1 + \ldots + a_{m,r} \cdot x_r \ge 0\\
\end{cases}.
\end{equation*}

When $\Phi$ is satisfiable, it entails a linear inequality
$$\psi := c_0 + c_1 \cdot x_1 + \dots + c_r \cdot x_r \ge 0$$

if and only if $\psi$ can be written as non-negative linear combination of the inequalities in $\Phi$ and the trivial inequality $1 \ge 0,$ i.e. if there exist non-negative real numbers $y_0,\dots,y_m$ such that

\begin{matrix}
c_0 = y_0 + \sum_{i=1}^k y_i \cdot a_{i, 0};\\
c_1 = \sum_{i=1}^k y_i \cdot a_{i, 1};\\
\vdots\\
c_r = \sum_{i=1}^k y_i \cdot a_{i, r}.
\end{matrix}

Moreover, $\Phi$ is unsatisfiable if and only if $-1 \ge 0$ can be derived as above.

The importance of Farkas' Lemma for us is that if we have a standard constraint of form (<ref>) and if the constraint includes only linear/affine inequalities, then we can use this lemma in Step 3 of our algorithm to reduce the standard constraint to QP, just as we did in Section <ref>. Moreover, Farkas' Lemma guarantees that this approach is not only sound but also complete. A corner case that we have to consider is when $\Phi$ is itself unsatisfiable and thus $\Phi \Rightarrow \psi$ holds vacuously. Fortunately, Farkas' Lemma also provides a criterion for unsatisfiability. In practice, we work with the following corollary of Theorem <ref> which can also handle strict inequalities.

Consider a set $V = \{x_1, \ldots, x_r\}$ of real-valued variables and the following system of equations over $V$:

\begin{equation*}
\Phi :=
\begin{cases}
a_{1,0} + a_{1,1} \cdot x_1 + \ldots + a_{1,r} \cdot x_r \bowtie_1  0\\
\hfil \vdots\\
a_{m,0} + a_{m,1} \cdot x_1 + \ldots + a_{m,r} \cdot x_r \bowtie_m 0\\
\end{cases}
\end{equation*}

where $\bowtie_i \in \{>,\geq\}$ for all $1 \leq i \leq m$. When $\Phi$ is satisfiable, it entails a linear inequality
$$\psi := c_0 + c_1 \cdot x_1 + \dots + c_r \cdot x_r \bowtie 0$$

with $\bowtie \in \{>, \geq\}$, if and only if $\psi$ can be written as non-negative linear combination of inequalities in $\Phi$ and the trivial inequality $1 > 0$. Note that if $\psi$ is strict, then at least one of the strict inequalities should appear with a non-zero coefficient in the linear combination. Moreover, $\Phi$ is unsatisfiable if and only if either $-1 \ge 0$ or $0 > 0$ can be derived as above.

We now consider extensions of Farkas' Lemma which can help us handle non-linear standard constraints in Step 3. The first extension is Handelman's theorem, which can be applied when the inequalities on the left hand side of (<ref>) are linear/affine, but the right hand side is a polynomial of arbitrary degree. To present the theorem, we first need the concept of monoids.

Consider a set $V = \{x_1,\dots x_r \}$ of real-valued variables and the following system of linear inequalities over $V$:

\begin{equation*}
\Phi :=
\begin{cases}
a_{1,0} + a_{1,1} \cdot x_1 + \ldots + a_{1,r} \cdot x_r \bowtie_1  0\\
\hfil \vdots\\
a_{m,0} + a_{m,1} \cdot x_1 + \ldots + a_{m,r} \cdot x_r \bowtie_m 0\\
\end{cases}
\end{equation*}

where $\bowtie_i \in \{>, \geq \}$ for all $1 \leq i \leq m$. Let $g_i$ be the left hand side of the $i$-th inequality, i.e. $g_i(x_1,\dots,x_r) := a_{i,0} + a_{i, 1} \cdot x_1 + \dots a_{i,r} \cdot x_r$. The monoid of $\Phi$ is defined as:

\textstyle \monoid(\Phi) := \left\{ \prod_{i=1}^m g_i^{k_i} \mid m \in \mathbb{N} ~\wedge~\forall i ~\; k_i \in \mathbb{N} \cup \{0\} \right\}.

In other words, the monoid contains all polynomials that can be obtained as a multiplication of the $g_i$'s. Note that $1 \in \monoid(\Phi).$
We define $\monoid_d(\Phi)$ as the subset of polynomials in $\monoid(\Phi)$ of degree at most $d.$

Consider a set $V = \{x_1, \dots, x_r \}$ of real-valued variables and the following system of equations over $V$:

\begin{equation*}
\Phi :=
\begin{cases}
a_{1,0} + a_{1,1} \cdot x_1 + \ldots + a_{1,r} \cdot x_r \geq  0\\
\hfil \vdots\\
a_{m,0} + a_{m,1} \cdot x_1 + \ldots + a_{m,r} \cdot x_r \geq 0\\
\end{cases}.
\end{equation*}

If $\Phi$ is satisfiable, $\sat(\Phi)$ is compact, and $\Phi$ entails a polynomial inequality $g(x_1,\dots,x_r) > 0,$ then there exist non-negative real numbers $y_1,\dots y_u$ and polynomials $h_1,\dots,h_u \in \monoid(\Phi)$ such that:

\textstyle g = \sum_{i=1}^u y_i \cdot h_i.

The intuition here is that if every inequality in $\Phi$ holds, then all the LHS expressions in $\Phi$ are non-negative and hence any multiplication $h_i$ of them is also non-negative. As in the case of Farkas' Lemma, Handelman's theorem shows that this approach is not only sound but also complete. We also need a variant that can handle strict inequalities in $\Phi.$

Consider a set $V = \{x_1, \dots, x_r \}$ of real-valued variables and the following system of equations over $V$:

\begin{equation*}
\Phi :=
\begin{cases}
a_{1,0} + a_{1,1} \cdot x_1 + \ldots + a_{1,r} \cdot x_r \bowtie_1  0\\
\hfil \vdots\\
a_{m,0} + a_{m,1} \cdot x_1 + \ldots + a_{m,r} \cdot x_r \bowtie_m 0\\
\end{cases}
\end{equation*}

in which $\bowtie_i \in \{>,\geq\}$ for all $1 \leq i \leq m$. If $\Phi$ is satisfiable and $\sat(\Phi)$ is bounded, then $\Phi$ entails a strong polynomial inequality $g \gg 0 $ if and only if there exist constants $y_0 \in (0,\infty),$ and $y_1,\dots,y_u \in [0,\infty)$, and polynomials $h_1,\dots,h_u \in \monoid(\Phi)$ such that:

\textstyle g = y_0 + \sum_{i=1}^u y_i \cdot h_i.

Corollary <ref> above can handle a wider family of standard constraints than Corollary <ref>. However, it is also more expensive, since we now need to generate one new variable $y_i$ for every polynomial in $\monoid_d(\Phi)$ instead of $\Phi$ itself. Moreover, there is no bound $d$ in the theorem itself, so introducing $d$ would lead to semi-completeness instead of completeness, i.e. the approach would be complete only if a large enough value of $d$ is used. As such, in cases where both sides of (<ref>) are linear, Corollary <ref> is preferable. We now consider a more expressive theorem that can handle polynomials on boths sides of (<ref>).

Given a finite collection of polynomials $\{g,g_1,\dots,g_k \} \in \mathbb{R}[x_1,\dots,x_n]$, let $\Phi$ be the set of inequalities defined as

\Phi: \{g_1 \geq 0,  \dots , g_k \geq 0 \}.

If $\Phi$ entails the polynomial inequality $g > 0$ and there exist some $i$ such that $\sat(g_i \geq 0)$ is compact, then there exist polynomials $h_0, \dots, h_k \in \mathbb{R}[x_1,\dots,x_n]$ such that

\textstyle g = h_0 + \sum_{i=1}^{m} h_i \cdot g_i

and every $h_i$ is a sum of squares. Moreover, $\Phi$ is unsatisfiable if and only if $-1 > 0$ can be obtained as above, i.e. with $g = -1$.

As in the cases of Farkas and Handelman, we need a variant of of Theorem <ref> that can handle strict inequalities in $\Phi.$

Consider a finite collection of polynomials $\{g,g_1,\dots,g_k \} \in \mathbb{R}[x_1,\dots,x_n]$ and let

\Phi: \{g_1 \bowtie_1 0,  \dots , g_k \bowtie_k 0 \}

where $\bowtie_i \in \{>,\geq \}$ for all $1 \leq i \leq k$. Assume that there exist some $i$ such that $\sat(g_i \geq 0)$ is compact or equivalently $\sat(g_i \bowtie_i 0)$ is bounded. If $\Phi$ is satisfiable, then it entails the strong polynomial inequality $g \gg 0$, iff there exists a constant $y_0 \in (0, \infty)$ and  polynomials $h_0, \dots, h_k \in \mathbb{R}[x_1,\dots,x_n]$ such that

\textstyle g = y_0 + h_0 + \sum_{i=1}^{k} h_i \cdot g_i,

and every $h_i$ is a sum of squares.

Trying to use the corollary above for handling standard constraints of form (<ref>) in Step 3 of our algorithm leads to two problems: (i) we also need a criterion for unsatisfiability of $\Phi$ to handle the cases where $\Phi \Rightarrow \psi$ holds vacuously, and (ii) in our QP, we should somehow express the property that every $h_i$ is a sum of squares. We now show how each of these challenges can be handled. To handle (i), we need another classical theorem from real algebraic geometry.

Given polynomials $g,g_1,\dots,g_k \in \mathbb{R}[x_1,\dots,x_n]$, exactly one of the following two statements holds:

* There exists $x \in \mathbb{R}^n,$ such that $g_1(x) = \dots= g_k(x) = 0$, but $g(x) \neq 0$.

* There exists $\alpha \in \mathbb{N}\cup \{0\}$ and polynomials $h_1,\dots,h_k \in \mathbb{R}[x_1,\dots,x_n]$ such that $\sum_{i=1} ^k h_i \cdot g_i - h_0= g^{2 \cdot \alpha}$ and $h_0$ is a sum of squares.

We now combine the Real Nullstellensatz with Putinar's Postivstellensatz to obtain a criterion for unsatisfiability of $\Phi.$

Consider a finite collection of polynomials $\{g_1,\dots,g_k \} \in \mathbb{R}[x_1,\dots,x_n]$ and the following system of inequalities:

\[
\Phi: \{g_1 \bowtie_1 0,  \dots , g_k \bowtie_k 0 \}
\]

where $\bowtie_i \in \{>,\geq \}$ for all $1 \leq i \leq k$. $\Phi$ is unsatisfiable if and only if at least one of the following statements holds:

* There exist a constant $y_0 \in (0,\infty)$ and sum of square polynomials $h_0,\dots,h_k \in \mathbb{R}[x_1,\dots,x_n]$ such that

\begin{equation}\label{eqn:1}
\textstyle     -1 = y_0 + h_0 + \sum_{i=1}^k h_i \cdot g_i.
\end{equation}

* There exist $\alpha \in \mathbb{N} \cup \{0\}$ and $h_0, h_1,\dots,h_k \in \mathbb{R}[x_1,\dots,x_n,w_1,\dots,w_k]$, such that for some $1 \leq j \leq m $ with $\bowtie_j \in \{> \}$, we have

\begin{equation}\label{eqn:2}
\textstyle    w_j ^{4 \cdot \alpha} = \sum_{i=1} ^m h_i  \cdot (g_i - w_i ^2) - h_0.
\end{equation}

where $h_0$ is a sum of squares in $\mathbb{R}[x_1,\dots,x_n]$. Note that $w_1, \dots, w_k$ are new variables.

First we show that if any of the two equalities (<ref>) or (<ref>) holds then $\Phi$ is unsatisfiable. Suppose $\Phi$ is satisfiable and pick $\val \in \sat(\Phi)$. Then, the RHS of (<ref>) is positive at $\val,$ whereas the LHS is negative. So, (<ref>) cannot hold. Now define $\tilde{g_i}(x_1,\dots,x_n,w_1,\dots,w_k) :=  g_i(x_1,\dots,x_n) - w_i^2.$  Using this definition, we can rewrite (<ref>) as $w_j ^{4.\alpha} = \sum_{i=1} ^m h_i \cdot \tilde{g_i} - h_0$. Moreover, $g_j^{2 \cdot \alpha} = (\tilde{g_j} + w_j^2)^{2 \cdot \alpha}.$ Expanding the right hand side using the binomial theorem, we get $g_j^{2 \cdot \alpha} = w_j^{4 \cdot \alpha} + h_j'\cdot \tilde{g_j}$ for some $h_j' \in \mathbb{R}[x_1,\dots,x_n,w_1,\dots,w_k]$. Now, substituting $w_j^{4 \cdot \alpha}$ (<ref>), we get

\begin{equation}\label{eqn:3}
\textstyle g_j^{2 \cdot \alpha} = \sum_{i=1}^m h_i \cdot \tilde{g_i} - h_0 + h_j' \cdot \tilde{g_j}
\end{equation}

Let us extend $\val$, which is a valuation of $\{v_1,\dots ,v_n\}$ to a valuation $\val'$ over $\{v_1,\dots,v_n,w_1,\dots,w_k\}$ such that $\val' \models \tilde{g_i}(v_1,\dots,v_n,w_1,\dots,w_k) = 0$ for all $1 \leq i \leq k$. Note that such an extension is always possible. We get a contradiction by evaluating (<ref>) on $\val'$ as the LHS is positive, whereas the RHS is negative.

We now prove the other side. Suppose $\Phi$ is unsatisfiable. We have two possibilities: either $\overline{\Phi}$ is also unsatisfiable or $\overline{\Phi}$ is satisfiable. Suppose $\overline{\Phi}$ is unsatisfiable, then using Theorem <ref>, $\overline{\Phi}$ entails $-2 > 0$ and we can write $-2 = h_0 + \sum_{i=1}^k h_i \cdot g_i$ for sum of squares polynomials $h_0,\dots,h_k \in \mathbb{R}[x_1,\dots,x_n]$. Therefore,

\textstyle -1 = 1 + h_0 + \sum_{i=1}^k h_i \cdot g_i

which fits into (<ref>). Now we are left with the case when $\overline{\Phi}$ is satisfiable but $\Phi$ is unsatisfiable. We first reorder the inequalities in $\Phi$ such that the non-strict inequalities appear first in the order. Let $j$ be the smallest index for which $\Phi[1\dots j]$, i.e. the set of first $j$ inequalities in $\Phi$, is unsatisfiable. By definition, $\Phi[1\dots j-1]$ is satisfiable and hence $\overline{\sat(\Phi[1\dots j-1])} = \sat(\overline{\Phi}[1\dots j-1])$. We can rewrite $\Phi[1\dots j] = \Phi[1\dots j-1]\wedge (g_j > 0) $. As $\Phi[1\dots j] $  is unsatisfiable, we know that $\Phi[1\dots j-1] $ entails $g_j \leq 0$. More precisely, this means $\sat(\Phi[1\dots j-1]) \subseteq \sat(g_j \leq 0)$. On taking closures on both sides we get $\emph{SAT}( \overline{\Phi}[1\dots j-1]) \subseteq \emph{SAT}(g_j \leq 0)$. This implies that $\overline{\Phi}[1,\dots j-1]$ entails $g_j \leq  0$ and hence $\overline{\Phi}[1\dots j]$ entails $g_j = 0$. As defined above, $\tilde{g_i}(x_1,\dots,x_n,w_1,\dots,w_k) = g_i(x_1,\dots,x_n) - w_i^2$. Now, we will show that there is no valuation $\val^*$ over the variables $\{x_1,\dots, x_n,w_1,\dots,w_k\}$ such that for all $1\leq i \leq j$, $\tilde{g_i}(\val^*) = 0$ but $g_j(\val^*) \neq 0$. Suppose there exist such a valuation $\val^*$. Let us define $\val$ to be the restriction of $\val^*$ to $\{x_1,\dots,x_n\}$. For each $1 \leq i \leq j $, we get $g_i(\val) \geq  0$ as $\tilde{g_j}(\val^*) = 0$. We also get $g_j(\val) = g_j(\val^*) \neq 0 $. Hence, we get a contradiction with the previous result that $\overline{\Phi}[1\dots j] $ entails $g_j = 0$. Applying the Real Nullstellensatz (Theorem <ref>) to $\tilde{g_i}$'s and $g_j$, we have

\textstyle g_j^{2\cdot \alpha} = \sum_{i=1}^j \tilde{h_i} \cdot \tilde{g_i} - h_0

where $\alpha$ is a non-negative integer and $\tilde{h_i}$'s and $h_0$ are sums of sqaures in $\mathbb{R}[x_1,\dots,x_n,w_1,\dots,w_k]$. Using the definition of $g_j = \tilde{g_j} + w_j^2$ and the binomial theorem, we get $g_j^{2 \cdot \alpha} = w_j^{4 \cdot \alpha} + h_j' \cdot \tilde{g_j}$ for some $h_j' \in \mathbb{R}[x_1,\dots,x_n,w_1,\dots,w_k]$. Therefore, we finally get the following expression:

\textstyle w_j^{4 \cdot \alpha} = \sum_{i=1}^j \tilde{h_i} \cdot (g_i-w_i^2) - h_j' \cdot (g_j - w_j^2) - h_0

which fits into the format of (<ref>), hence completing the proof.

Finally, we provide the needed theorems to check that a certain polynomial $h$ is a sum of squares using QP.

Let $\vec{a}$ be the vector of all $\binom{n+d}{d}$ monomials of degree less than or equal to $d$ over the variables $\{x_1, \dots, x_n\}.$ A polynomial $p \in \mathbb{R}[x_1,\dots,x_n]$ of degree $2 \cdot d$ is a sum of squares if and only if there exist a positive semidefinite matrix $Q$ of order $\binom{n+d}{d}$ such that $p = a^T \cdot Q \cdot a.$

A symmetric square matrix $Q$ is positive semidefinite if and only if it has a Cholesky decomposition of the form $Q = L L^T$ where $L$ is a lower-triangular matrix with non-negative diagonal entries.

Based on Theorems <ref> and <ref>, a polynomial $p$ of degree $2 \cdot d$ is a sum of squares if and only if it can be written as
p = a^T \cdot L \cdot L^T \cdot a
such that diagonal entries of $L$ are non-negative. This representation provides us with a simple approach to generate a sum-of-squares polynomial of degree $2 \cdot d$ with symbolic coefficients and encoding them in QP. We first generate a lower triangular matrix $L$ of order $ \binom{n+d}{d}$ by creating one fresh variable for each entry in the lower triangle and adding the extra condition that the entries on the diagonal must be non-negative. Then, we symbolically compute $a^T \cdot L \cdot L^T \cdot a.$

§.§ Details of Step 3 of the Algorithm

We now have all the necessary ingredients to provide a variant of Step 3 of the algorithm that preserves completeness. We assume a positive integer $d$ is given as part of the input. This $d$ serves as an upper-bound on the degrees of polynomials/templates that we use in Handelman's theorem (the monoid) and the Stellensätze. Our approach is complete as long as a large enough $d$ is chosen.

Step 3: Eliminating Program Variables and Reduction to QP Recall that at the end of Step 2, we have a finite set of standard constraints of form (<ref>). In this step, the algorithm handles each standard constraint separately and reduces it to quadratic programming over template variables and newly-introduced variables, hence effectively eliminating the program variables and the quantification over them. Let
f_1 \bowtie_1 0 ~\wedge~ f_2 \bowtie_2 0 ~\wedge~ \ldots \wedge f_r \bowtie_r 0 \Rightarrow f \bowtie 0
$ be one of the standard constraints. The algorithm considers three cases: (i) if all the inequalities on both sides of the constraint are affine, then it applies Farkas' Lemma; (ii) if the LHS inequalities are affine but the RHS is a higher-degree polynomial, then the algorithm applies Handelman's theorem; and (iii) if the LHS contains higher-degree polynomials, the algorithm applies the Stellensätze and Theorem <ref>. Below, we define $\Phi : \{ f_1 \bowtie_1 0, f_2 \bowtie_2 0, \ldots , f_r \bowtie_r 0 \}$ and $\psi: f \bowtie 0.$

Step 3.(i). Applying Farkas' Lemma Assuming all the constraints in $\Phi$ and $\psi$ are affine, we can apply Corollary <ref>. Based on this corollary, we have to consider three cases disjunctively:

* $\Phi$ is satisfiable and entails $\psi$: The algorithm creates new template variables $y_0,\dots,y_r$ with the constraint $y_i \geq 0$ for every $i$. It then symbolically computes
f = y_0 + \sum_{i=1}^r y_i \cdot f_i.
$ The latter is a polynomial equality over $\vars.$ So, the algorithm equates the coefficients of corresponding monomials on both sides, hence reducing this case to QP. Additionally, if $\psi$ is a strict inequality, the algorithm adds the extra constraint $\sum_{\bowtie_i \in \{ > \}} y_i > 0.$

* $\Phi$ is unsatisfiable and $-1 \geq 0$ can be obtained: This is similar to the previous case, except that $-1$ should be written as $y_0 + \sum_{i=1}^r y_i \cdot f_i.$

* $\Phi$ is unsatisfiable and $0 > 0$ can be obtained: This is also similar to the last two cases. We have $0 = y_0 + \sum_{i=1}^r y_i \cdot f_i$ and $\sum_{\bowtie_i \in \{ > \}} y_i > 0.$

Note that the template variables $y$ are freshly generated in each case above. Also, we have to consider cases (2) and (3) because $\Phi$ is unsatisfiable in these cases and hence the constraint $\Phi \Rightarrow \psi$ always holds vacuously.

Step 3.(ii). Applying Handelman's Theorem Assuming all constraints in $\Phi$ are linear but $\psi$ is a higher-degree polynomial inequality, the algorithm applies Corollaries <ref> and <ref>. Again, we have to consider the same three cases as in Step 3.(i):

* $\Phi$ is satisfiable and entails $\psi$: We apply Corollary <ref>. The algorithm first symbolically computes $\monoid_d(\Phi) = \{h_1, h_2, \dots, h_u\}.$ It then generates new template variables $y_0, y_1, \dots, y_u$ and constrains them by setting $y_0, y_1, y_2, \dots, y_u \geq 0.$ If $\psi$ is a strict inequality, it further adds the constraint $y_0 > 0.$ It then symbolically computes the equality

\textstyle f = y_0 + \sum_{i=1}^u y_i \cdot h_i.

As usual, this is an equality whose both sides are polynomials over $\vars.$ So, the algorithm equates the coefficients of corresponding monomials on the LHS and RHS, which reduces this case to QP.

* Note that $\Phi$ consists of linear inequalities, so we can use Farkas' Lemma to check if $\Phi$ is unsatisfiable. As such, this step is the same as case (2) of Step 3.(i).

* This is the same as case (3) of Step 3.(i).

Step 3.(iii). Applying Stellensätze
If $\Phi$ includes polynomial inequalities of degree $2$ or larger, then we have to apply Corollary <ref> and Theorem <ref>. The algorithm considers three cases and combines them disjunctively:

* $\Phi$ is satisfiable and entails $\psi$: In this case, we apply Corollary <ref>. The algorithm generates template sum-of-squares polynomials $h_0, \ldots, h_r$ of degree $d$ and adds QP constraints that ensure each $h_i$ is a sum of squares (See the end of Section <ref>). It also generates a non-negative fresh variable $y_0$. If $\psi$ is strict, the algorithm adds the constraint $y_0 > 0.$ Finally, the algorithm symbolically computes

\textstyle	 f = y_0 + h_0 + \sum_{i=1}^r h_i \cdot f_i;

and equates the corresponding coefficients in the LHS and RHS to obtain QP constraints.

* $\Phi$ is unsatisfiable due to the first condition of Theorem <ref>: This case is handled similary to case (1) above, except that we have $-1 = y_0 + h_0 + \sum_{i=1}^r h_i \cdot f_i.$

* $\Phi$ is unsatisfiable due to the second condition of Theorem <ref>: The algorithm introduces $r$ new program variables $w_1, \dots, w_r.$ It then generates a sum-of-squares template polynomial $h$ over $\vars$ and arbitrary template polynomials $h_1, \dots, h_k$ over $\vars \cup \{w_1, \dots, w_r\}.$ All $h_i$'s are degree-$d$ templates. Finally, for every index $j$ that corresponds to a strict inequality, i.e. $\bowtie_j \in \{ > \},$ the algorithm symbolically computes

\textstyle	 w_j^{d'} = \sum_{i=1}^r h_i \cdot (f_i - w_i^2) - h_0,

in which $d'$ is the largest multiple of $4$ that does not exceed $d$. Note that this is an equality between two polynomials over $\vars \cup \{w_1, \dots, w_r\}.$ As before, the algorithm equates the coefficients of corresponding monomials on both sides of the equality and reduces it to QP.

After the algorithm runs Step 3 as above, all standard constraints generated in Step 2 will be reduced to QP and can hence be passed to an external solver in Step 4, as illustrated in Section <ref>.

Degree bounds We are using the same bound $d$ for the degree of the polynomials and templates in all cases above. This is not a requirement. One can fix different degree bounds for each case.

Handling Boundedness To achieve completeness, we need the boundedness assumption, i.e. that for every variable $v \in \vars,$ we always have $-M \leq v \leq M.$ To model this in the algorithm, we can add the boundedness inequalities to the left-hand side of every standard constraint. Additionally, we can get a concrete value for $M$ as part of the input, or treat $M$ symbolically, i.e. as a template variable, and let the QP solver synthesize a value for it.

§.§ Soundness and Completeness

We now prove that our algorithm is sound and semi-complete for TBSP. The soundness result needs no extra assumptions and can be obtained directly. The completeness on the other hand relies on several assumptions: (i) boundedness, (ii) having chosen a large enough degree bound $d$, and (iii) having invariants and post-conditions that, when written in DNF form, consist only of strongly positive polynomial inequalities of the form $g \gg 0.$

Given a sketch polynomial program or equivalently a sketch polynomial transition system (SPTS) of the form $(\vars, \tvars, \locs,\loc_0,\theta_0,\loc_f, \theta_f, \transitions, \invariant)$ together with a cutset $\cutset$ as input, every concrete polynomial transition system (PTS) synthesized by the algorithm above is inductively valid.

The standard constraints of form (<ref>) generated in Step 2 are equivalent to the initiation, consecution and finalization constraints in the definition of inductive validity. The reduction from standard constraints to QP in Step 3 is also sound since, in every case, it either writes the RHS of the standard constraint as a combination of the LHS polynomials, hence proving that it holds, or otherwise proves that the LHS is unsatisfiable and thus the standard constraint holds vacuously.

Consider a solvable sketch polynomial transition system (SPTS) of the form $(\vars, \tvars, \locs,\loc_0,\theta_0,\loc_f, \theta_f, \transitions, \invariant),$ together with a cutset $\cutset,$ that is given as input. Moreover, assume that:

* The boundedness assumption holds, i.e. there is a constant $M \in (0, \infty)$ such that for every $v \in \vars,$ we always have $-M \le v \le M.$

* Every invariant $\invariant(\loc)$ and post-condition $\theta_f,$ when written in disjunctive normal form, contains only strongly positive polynomial inequalities of the form $g \gg 0.$

Then, there exists a constant degree bound $d \in \mathbb N$, for which the algorithm above is guaranteed to successfully synthesize an inductively valid polynomial transition system (PTS).

Since our instance is solvable, there is a valuation $\tval$ for the template variables that yields an inductively valid PTS. We prove that for large enough $d$, the valuation $\tval$ is obtained by one of the solutions of the QP instance solved in Step 4 of the algorithm. The proof is pretty straightforward since we just have to check that all steps of our algorithm are complete. Step 2 is complete since it simply rewrites the inductive validity constraints as an equivalent set of standard constraints. For Step 3, we prove completeness of each case separately. Step 3.(i) is complete due to Farkas' Lemma (Corollary <ref>). In Step 3.(ii), if $\Phi$ is unsatisfiable, then the algorithm is complete based on Corollary <ref>. Otherwise, it is complete based on Corollary <ref>. However, since we are using a degree bound $d$ for the monoid, the completeness only holds if the chosen $d$ is large enough. Moreover, Corollary <ref> requires strong positivity of $g$, which corresponds to invariants and post-conditions in our use-case, and boundedness of $\sat(\Phi),$ which is a direct consequence of our boundedness assumption. Finally, Step 3.(iii) is complete due to Corollary <ref> and Theorem <ref>. These depend on $d$, strong positivity and boundedness in the exact same manner as in the case of Step 3.(ii).

Limitations of Completeness The main limitation in our completeness result is that it holds only if the degree bound $d$ chosen for the template polynomials is large enough. This is why we call it a semi-completeness theorem. In theory, it is possible to come up with adversarial instances in which the required degree is exponentially high <cit.>. In practice, we rarely, if ever, need to use a higher degree than that of the polynomials that are already part of the input SPTS. The second limitation is boundedness. This limitation cannot be lifted since both Handelman's Theorem and Putinar's Positivstellensatz assume compactness, which is equivalent to being closed and bounded in $\mathbb R^n$. Nevertheless, it does not have a significant practical effect and the algorithm remains sound even without this assumption. It is also noteworthy that the treatment of the linear/affine case using Farkas' Lemma requires neither boundedness nor any specific value of $d$ and is always complete. Finally, strong positivity in the invariants and post-conditions is needed because Putinar's Positivstellensatz and Theorem <ref> can only provide a sound and complete characterization of strongly positive polynomials over a bounded semi-algebraic set $X \subseteq R^n$ if we do not assume that $X$ itself is closed. In terms of the synthesis problem, this means that our algorithm is not guaranteed to be complete for inequalities of the form $f > 0$ in the invariants/post-conditions for which the value of $f$ in the runs of the program can get arbitrarily close to $0.$ However, this limitation is also not significant in practice because (i) our soundness does not depend on it, and (ii) $f + \epsilon \gg 0$ holds for any constant $\epsilon > 0.$ So, a small change (by any value $\epsilon > 0$) in the invariants/postconditions leads to an instance over which our completeness holds.

§ DECIDABILITY AND NP-HARDNESS

In this section, we study the decision variant of TBSP, in which we are given an SPTS of the form $\Gamma = (\vars, \tvars, \locs,\loc_0,\theta_0,\loc_f, \theta_f, \transitions, \invariant)$ and should decide whether there exists a valuation $\tval: \tvars \rightarrow \mathbb R$ for the template variables under which $\Gamma$ turns into an inductively valid PTS.

Reduction to the First-order Theory of the Reals Suppose $\vars = \{v_1, \dots, v_m\}$ and $\tvars = \{t_1, \dots, t_{m'}\}.$ Recall that Step 2 of our algorithm is sound and complete and reduces all the constraints to the form shown in Equation (<ref>). However, in the synthesis algorithm, we are looking for a valuation to the template variables $\tvars$ that ensures (<ref>). Hence, in the decision variant, we are interested in deciding the following first-order formula over the reals:

\begin{equation} \label{eq:tarski}
\exists t_1, t_2, \ldots, t_{m'} ~\;~\; \forall v_1, v_2, \ldots, v_m ~\;~\;~ (\kappa_1 ~\wedge~ \kappa_2 ~\wedge~ \ldots ~\wedge~ \kappa_u)
\end{equation}

It is well-known that the first-order theory of the reals is decidable and hence so is our problem. Specifically, one can apply Tarski's quantifier elimination method <cit.> to (<ref>). Note that while this proves decidability in theory, decision procedures and quantifier elimination on the first-order theory of the reals are notoriously unscalable and since (<ref>) contains a quantifier alternation, it is practically beyond the reach of all modern solvers, even for toy programs.

NP-hardness We prove our problem is strongly NP-hard, even for linear/affine transition systems, by providing a reduction from 3SAT. Consider a 3SAT instance
\bigwedge_{i=1}^m \bigvee_{j=1}^3 l_{i, j}
$ over the boolean variables $\{x_1, \ldots, x_n\}, $ i.e. each literal $l_{i, j}$ is either an $x_k$ or a $\neg x_k.$ We define the expression $e(l_{i,j})$ as $x_k$ if $l_{i, j} = x_k$ and $1-x_k$ if $l_{i, j} = \neg x_k.$ To reduce this 3SAT instance to our problem, we let $\vars = \{x_1, \ldots, x_n\}$ and $\tvars = \{t_1, \ldots, t_n\}.$ We then consider the sketch transition system corresponding to the following program:

@real: $x_1, x_2, \dots, x_n$
@pre: $x_1 = t_1 ~\wedge~ x_2 = t_2 ~\wedge~ \dots ~\wedge~ x_n = t_n$;
if ($x_1 > 0.5$)
$x_1$ = $1$;
$x_1$ = $0$;
if ($x_n > 0.5$)
$x_n$ = $1$;
$x_n$ = $0$;
@post: $\bigwedge_{i=1}^m e(l_{i, 1})+e(l_{i, 2})+e(l_{i, 3}) > 0.5$;

Note that the template variables are appearing in the pre-condition. It is easy to see that the synthesis instance has a solution if and only if the 3SAT instance is satisfiable.

§ IMPLEMENTATION AND EXPERIMENTAL RESULTS

In this section, we report on a prototype implementation of our algorithm and our experimental results.

Implementation We implemented our approach in Python and used SymPy <cit.> for symbolic computations. We also used the Z3 SMT solver <cit.> to handle the final QP instances. In all cases, we used quadratic templates for the holes and set the technical parameter (degree upper-bound) $d$ to $2.$

Experimental Setting All experimental results were obtained on an Intel i9-10980HK Processor (8 cores, 16 threads, 5.3 GHz, 16 MB Cache) with 32 GB of RAM running Microsoft Windows 10.

Heuristics We used the following simple heuristics to speed up the solution of the QP instance by the SMT solver:

* Simplifying Standard Constraints: We apply the simplification procedure illustrated in Example <ref>. Specifically, it is often the case that for some $v_1, v_2 \in \vars \cup \vars'$, the standard constraint contains $v_1 = v_2$ on its left-hand side. In these cases, we merge $v_1$ and $v_2$ into a single variable. Similarly, if the left-hand side of a standard constraint has $v_1 = c$ for a real constant $c$, then we simply replace every occurrence of $v_1$ with $c.$

* Strengthening the Constraints: Our algorithm introduces a large number of template variables, not only in Step 1, but also in Step 3, where each of the Stellensätze lead to the creation of new coefficient variables. Such new variables are created, for example, in the linear combinations appearing in Farkas/Handelman and the sum-of-squares template polynomials in Putinar. For most real-world cases, the vast majority of the template variables should be set to $0$ in solutions to the QP. However, their inclusion makes the QP much larger and more complicated and can even lead to the failure of other heuristics used by the SMT solver. Hence, we use a process which is intuitively the opposite of abstraction-refinement. We start with a strengthened version of our QP instance and repeatedly make it weaker, but never weaker than the original instance, until the SMT solver finds a solution. Specifically, we first add additional constraints of the form $t=0$ for every template variable $t$. This strengthens the constraints but makes the problem easier for the SMT solver, since it can now simply plug in the $0$ values and solve a lower-dimensional QP. Additionally, it would very likely lead to an unsatisfiable QP. However, Z3 can then provide us with an unsatisfiability core, i.e. a minimal subset of QP constraints that are unsatisfiable. If the core does not include any of our $t=0$ constraints, then we know that the original QP is also unsatisfiable. Otherwise, we simply remove all the $t=0$ constraints in the core and repeat the same process.

Experimental Results We now illustrate several programs that we synthesized using our implementation. In each case, we show the completed program, in which the synthesized parts are boxed. Due to space limitations, some examples are relegated to Appendix <ref>. Table <ref> provides a summary of the runtimes of our algorithm over these examples.

@prog: ClosestCubeRoot
@real: $a,x,s,r$;
@pre: $a \ge 0$;
$x = a$;
$r = 1$;
s = 3.25;
@invariant: $x \ge 0 \wedge -12 \cdot r^2 + 4 \cdot s = 1 \wedge 4 \cdot r^3 - 6 \cdot r^2 + 3 \cdot r + 4 \cdot x - 4 \cdot a = 1$
while ($x-s \ge 0$)

$x = \Hole{x-s}$;
$s = \Hole{3+s+6\cdot r}$;
$r = r + 1$;

@post: $4 \cdot r^3 + 6 \cdot r^2 + 3 \cdot r \geq 4 \cdot a \geq 4 \cdot r^3 - 6 \cdot r^2 + 3 \cdot r - 1$

@prog: ClosestSquareRoot
@real: $a,x,r$;
@pre: $a \ge 1$;
$x = 0.5 \cdot a$;
$r = 0$;
@invariant: $a = \Hole{2 \cdot x + r^2 - r} ~\wedge~ x \ge 0$
while ($x \ge r$)

$x = \Hole{x-r}$;
$r = r + 1$;

@post: $r^2-r \ge a-2 \cdot r ~\wedge~ r^2 - r \le a$

@prog: SquareRootFloor
@real: $a,su,t,n$;
@pre: $n \ge 0$;
$a = 0$;
$su = 1$;
$t = 1$;
@invariant: $a^2 \le n \wedge t = 2 \cdot a + 1 \wedge su = (a+1)^2$
while ($su \le n$)

$a = a + 1$;
$t = \Hole{2+t}$;
$su = \Hole{su + t}$;

@post((a * a) <= n);

@prog: SquareRootApproximation
@real: $a,err,r,q,p$;
@pre: $a \ge 1 \wedge err \ge 0$;
$r = a -1$;
$q = 1$;
$p = 0.5$;
@invariant: $p \ge 0 \wedge r \ge 0 \wedge a = q^2 + 2 \cdot r \cdot p$
while ($2 \cdot p \cdot r \ge err$)

if ($2 \cdot r - 2 \cdot q - p \ge 0$)

$r = \Hole{2\cdot r -2\cdot q -p}$;
$q = p + \Hole{q}$;
$p = p / 2$;

$r = \Hole{2 \cdot r}$;
$p = p \cdot 0.5$;

@post:$q^2 \geq a - err \wedge q^2 \leq a$

@prog: ConsecutiveCubes
@real: $N,n,x,y,z,s$;
@pre: $1 \ge 0$;
$n = 0$;
$x = 0$;
$y = 1$;
$z = 6$;
$s = 0$;
@invariant: $z = 6 \cdot n \wedge y = 3\cdot n^2+3\cdot n+1 \wedge x=n^3$
while ($n \le N$)

if ($x = n^3$)

$s = s + x$;

$n = n + 1$;
$x = \Hole{y+x}$;
$y = \Hole{y+z}$;
$z = \Hole{6+z}$;

@post: $1 \ge 0$;

* PositivityEnforcement: We are synthesizing an expression such that adding it to the given polynomial makes it positive for all possible inputs.

* SquareCompletion: We are synthesizing a number such that adding it to the given polynomial makes it always positive.

* SlidingBody:
Consider the following Physics problem: An object is placed on an inclined plane of given length $l$. The inclination is given in terms of the sine and cosine of its angle. The coefficients of static and kinetic friction are also known. We want to find whether the body will slide and, if so, find the time it will take for it to reach the bottom of the plane.

A physics student can easily write a program template as shown in the example for solving this problem using the principles and equations from classical mechanics. For example, the student can use an if block to decide if the body will move or not. Then, inside the block, they can use an equation of motion that holds for the various parameters involved in the problem. Now, to find the time it takes to slide down, they can leave a program hole in the template.

* ArchimedesPrinciple:
A cuboid with dimensions ($l$, $b$, $h$) of density $\rho$ is placed in a fluid of density $\rho'$. We want to determine whether the body will float, and if it floats, find the percentage of the volume of the body that will be outside the fluid.

* ClosestCubeRoot<cit.>: Given a positive real number $n$, we need to synthesize a program that computes the closest integer cube root of $n$.

* SquareRootFloor<cit.>: Given a positive real number $n$, we need to synthesize a program that computes the floor of the square root of $n$.

* SquareRootApproximation<cit.>: Given a positive real number $n$ and an error threshold $\epsilon$, we need to synthesize a program that computes $\sqrt n$.

* ConsecutiveCubes<cit.>: We want to synthesize a program that computes the sum of first $n$ cubes.

Example     Step 2     Step 3     Step 4 (Z3)

ArchimedesPrinciple     0.10     0.27     0.22

ClosestCubeRoot     0.13     12.48     3.20

ClosestSquareRoot     0.09     2.21     1.31

ConsecutiveCubes     0.14     35.00     20.92

PolynomialApproximation     0.09     0.10     0.05

PositivityEnforcement     0.08     0.08     0.03

SlidingBody     0.10     0.30     0.22

SquareCompletion     0.07     0.08     0.04

SquareRootApproximation     0.12     14.87     5.54

SquareRootFloor     0.11     8.09     2.31

The runtimes of the different steps of our algorithm in seconds.

In summary, our algorithm is able to handle various polynomial synthesis problems, including examples from physics and  <cit.> in a few seconds. Hence, not only is the algorithm semi-complete in theory, but it is also applicable in practice.

§ ACKNOWLEDGMENTS

This section was removed to preserve anonymity.

§.§ Heuristics using Counter Example Guided Abstraction Refinement (CEGAR)

After we reduce the constraint pairs to system of polynomial inequations and equations over template variables in step <ref> of the synthesis algorithm, we then pass the system to solvers in step <ref> of the algorithm. Generally for examples arising from the real world use cases we get a large number of polynomial inequations and equations. In practice, these solvers only work efficiently for very small system of equations. In this section we describe the heuristics we have developed to handle relatively larger examples.

The idea of the heuristics is as follows: Given an instance $I$ of the system of polynomial inequations, we put constraints on some of the template variables $l$ to be equal to zero to get a simplified instance $I'$ of the problem. If our solvers solve the instance $I'$, then we also get a solution for the original instance $I$. Otherwise we relax some of the constraints on template variables and run the entire process again.

Now we formally describe the above idea in the following two steps.

Abstraction: We find an abstraction of the instance $I$ of the system of polynomial inequations by imposing additional constraints on a set of template variables $l \in L' \subseteq L$ to $0$. Note that even though it seems that we have added more constraints to the system, we are actually dropping the dimension of the search space for solutions. As we have assigned $0$ to all the template variables $l \in L'$, the solver only has to search for solutions in the template variables $U \cup L \setminus L'$. Therefore the instance $I'$ of the problem which we get after assignment is relatively easier to solve than the original instance $I$ of the problem. Now, if the solver returns a solution $M: (L \cup U \setminus L') \to \mathbb{R}$ to the instance $I'$ then we also get a model $M: L \cup U \cup L' \to \mathbb{R}$ as valuation of variables in $L'$ is fixed to zero in $I$ to obtain $I'$. Hence, this method is sound.

Suppose the solver returns that the instance $I'$ of the problem is unsatisfiable then there are two possibilites, i.e., either the original instance of the problem $I$ is also unsatisfiable or we need to relax some constraints on template variables $L'$. In other words, we need to refine the instance $I'$.

Refinement: To refine the instance $I'$ of the problem, we need to select the template variables $l \in L'$ for which we need to relax the constraints. We now find the minimal set of constraints $C$ which make the instance $I'$ unsatisfiable. Suppose there is a subset $K$ of $L'$ such that for each $l \in K$ we have that $l = 0$ is contained in the set of constraints $C$. We then remove all variables contained in $K$ from $L'$ and obtain a refinement of instance $I'$ of the problem. If such a susbet $K$ of $L'$ doesn't exists then the original instance $I$ is unsatisfiable.

Note that on every refinement step, we remove at least one variable from $L'$. In the worst case scenario, we have to set $L' = \emptyset$, and this will refine the instance $I'$ to the original instance $I$, and this shows that the method is complete.

In our experiments, we observed that the heuristics were able to solve many instances which were timing out in solvers otherwise.

#1     #1#1#1     #1     #1     #1     #1#1          #1#1

[Abate, David, Kesseli, Kroening, and
PolgreenAbate et al2018]
authorpersonAlessandro Abate, personCristina
David, personPascal Kesseli, personDaniel Kroening,
and personElizabeth Polgreen. year2018.
Counterexample guided inductive synthesis modulo
theories. In booktitleCAV. pages270–288.

[Agrawal, Chatterjee, and NovotnỳAgrawal
et al2017]
authorpersonSheshansh Agrawal,
personKrishnendu Chatterjee, and personPetr
Novotnỳ. year2017.
Lexicographic ranking supermartingales: an
efficient approach to termination of probabilistic programs. In

[Alur, Bodik, Juniwal, Martin, Raghothaman,
Seshia, Singh, Solar-Lezama, Torlak, and UdupaAlur et al2013]
authorpersonRajeev Alur, personRastislav
Bodik, personGarvit Juniwal, personMilo MK Martin,
personMukund Raghothaman, personSanjit A Seshia,
personRishabh Singh, personArmando Solar-Lezama,
personEmina Torlak, and personAbhishek Udupa.
Syntax-Guided Synthesis.
journalFMCAD (year2013).

[Alur, Fisman, Padhi, Reynolds, Singh, and
UdupaAlur et al2019]
authorpersonRajeev Alur, personDana Fisman,
personSaswat Padhi, personAndrew Reynolds,
personRishabh Singh, and personAbhishek Udupa.
titleSyGuS-Comp: Syntax-guided Synthesis Competition.

[Alur, Singh, Fisman, and Solar-LezamaAlur
et al2018]
authorpersonRajeev Alur, personRishabh
Singh, personDana Fisman, and personArmando
Solar-Lezama. year2018.
Search-based program synthesis.
journalCommun. ACM volume61,
number12 (year2018), pages84–93.

[Asadi, Chatterjee, Fu, Goharshady, and
MahdaviAsadi et al2021]
authorpersonAli Asadi, personKrishnendu
Chatterjee, personHongfei Fu, personAmir Kafshdar
Goharshady, and personMohammad Mahdavi.
Polynomial reachability witnesses via
Stellensätze. In booktitlePLDI.

[Blekherman, Parrilo, and ThomasBlekherman
et al2012]
authorpersonGrigoriy Blekherman,
personPablo A Parrilo, and personRekha R Thomas.
booktitleSemidefinite optimization and convex
algebraic geometry.

[Bochnak, Coste, and RoyBochnak
et al2013]
authorpersonJacek Bochnak, personMichel
Coste, and personMarie-Françoise Roy.
booktitleReal algebraic geometry.
Vol. volume36.
publisherSpringer Science & Business Media.

[Chatterjee, Fu, and GoharshadyChatterjee
et al2016]
authorpersonKrishnendu Chatterjee,
personHongfei Fu, and personAmir Kafshdar
Goharshady. year2016.
Termination analysis of probabilistic programs
through Positivstellensatz’s. In booktitleCAV.

[Chatterjee, Fu, Goharshady, and
GoharshadyChatterjee et al2020b]
authorpersonKrishnendu Chatterjee,
personHongfei Fu, personAmir Kafshdar Goharshady,
and personEhsan Kafshdar Goharshady.
Polynomial invariant generation for
non-deterministic recursive programs. In booktitlePLDI.

[Chatterjee, Fu, and NovotnỳChatterjee
et al2020a]
authorpersonKrishnendu Chatterjee,
personHongfei Fu, and personPetr Novotnỳ.
Termination analysis of probabilistic programs with
journalFoundations of Probabilistic Programming
(year2020), pages221–258.

authorpersonAlonzo Church.
Application of recursive arithmetic to the problem
of circuit synthesis.
journalJournal of Symbolic Logic
volume28, number4 (year1963).

[Colón, Sankaranarayanan, and
SipmaColón et al2003]
authorpersonMichael A Colón,
personSriram Sankaranarayanan, and personHenny B
Sipma. year2003.
Linear invariant generation using non-linear
constraint solving. In booktitleCAV.

authorpersonJulius Farkas.
Theory of simple inequalities.
journalJournal for pure and applied mathematics
(Crelles Journal) volume1902, number124
(year1902), pages1–27.

[Fedyukovich, Prabhu, Madhukar, and
GuptaFedyukovich et al2019]
authorpersonGrigory Fedyukovich,
personSumanth Prabhu, personKumar Madhukar, and
personAarti Gupta. year2019.
Quantified invariants via syntax-guided synthesis.
In booktitleCAV. pages259–277.

[Feng, Zhang, Jansen, Zhan, and XiaFeng
et al2017]
authorpersonYijun Feng, personLijun Zhang,
personDavid N Jansen, personNaijun Zhan, and
personBican Xia. year2017.
Finding polynomial loop invariants for
probabilistic programs. In booktitleATVA.

authorpersonAmir Goharshady.
titleParameterized and algebro-geometric advances in
static program analysis.
thesistypePh.D. Dissertation. schoolInstitute
of Science and Technology Austria.

authorpersonCordell Green.
Application of theorem proving to problem solving.
In booktitleReadings in Artificial Intelligence.

authorpersonSumit Gulwani.
Automating string processing in spreadsheets using
input-output examples. In booktitlePOPL.

[Gulwani, Harris, and SinghGulwani
et al2012]
authorpersonSumit Gulwani, personWilliam R
Harris, and personRishabh Singh.
Spreadsheet data manipulation using examples.
journalCommun. ACM volume55,
number8 (year2012), pages97–105.

[Gulwani, Jha, Tiwari, and VenkatesanGulwani
et al2011]
authorpersonSumit Gulwani, personSusmit
Jha, personAshish Tiwari, and personRamarathnam
Venkatesan. year2011.
Synthesis of loop-free programs. In
booktitlePLDI. pages62–73.

[Gulwani, Polozov, and SinghGulwani
et al2017]
authorpersonSumit Gulwani, personOleksandr
Polozov, and personRishabh Singh.
Program synthesis.
journalFoundations and Trends in Programming
Languages volume4, number1-2
(year2017), pages1–119.

[Guo, James, Justo, Zhou, Wang, Jhala, and
PolikarpovaGuo et al2019]
authorpersonZheng Guo, personMichael James,
personDavid Justo, personJiaxiao Zhou,
personZiteng Wang, personRanjit Jhala, and
personNadia Polikarpova. year2019.
Program synthesis by type-guided abstraction
refinement. In booktitlePOPL. pages1–28.

authorpersonDavid Handelman.
Representing polynomials by positive linear
functions on compact convex polyhedra.
journalPacific J. Math. volume132,
number1 (year1988), pages35–62.

authorpersonRobin Hartshorne.
booktitleAlgebraic geometry.
Vol. volume52.
publisherSpringer Science & Business Media.

[Knoth, Wang, Polikarpova, and HoffmannKnoth
et al2019]
authorpersonTristan Knoth, personDi Wang,
personNadia Polikarpova, and personJan Hoffmann.
Resource-guided program synthesis. In
booktitlePLDI. pages253–268.

authorpersonAndrej Kolmogoroff.
Zur deutung der intuitionistischen logik.
journalMathematische Zeitschrift
volume35, number1 (year1932),

authorpersonJohn R Koza.
Genetic programming as a means for programming
computers by natural selection.
journalStatistics and computing
volume4, number2 (year1994),

authorpersonKrzysztof Krawiec.
booktitleBehavioral program synthesis with genetic
programming. Vol. volume618.

[Liu, Fu, Yu, Song, and LiLiu
et al2022]
authorpersonHongming Liu, personHongfei Fu,
personZhiyong Yu, personJiaxin Song, and
personGuoqiang Li. year2022.
Location-by-Location Linear Invariant Generation
with Farkas' Lemma. In booktitleOOPSLA.

[Mai and MagronMai and Magron2022]
authorpersonNgoc Hoang Anh Mai and
personVictor Magron. year2022.
On the complexity of Putinar–Vasilescu's
journalJournal of Complexity
(year2022), pages101663.

[Manna and WaldingerManna and
authorpersonZohar Manna and
personRichard J Waldinger. year1971.
Toward automatic program synthesis.
journalCommun. ACM volume14,
number3 (year1971), pages151–165.

[Meurer, Smith, Paprocki, Čertík,
Kirpichev, Rocklin, Kumar, Ivanov, Moore, Singh, et alMeurer
et al2017]
authorpersonAaron Meurer,
personChristopher P Smith, personMateusz Paprocki,
personOndřej Čertík, personSergey B
Kirpichev, personMatthew Rocklin, personAMiT Kumar,
personSergiu Ivanov, personJason K Moore,
personSartaj Singh, et al
SymPy: symbolic computing in Python.
journalPeerJ Computer Science
volume3 (year2017), pagese103.

[Moura and BjørnerMoura and
authorpersonLeonardo de Moura and
personNikolaj Bjørner. year2008.
Z3: An efficient SMT solver. In
booktitleTACAS. pages337–340.

[Neumann, Ouaknine, and WorrellNeumann
et al2020]
authorpersonEike Neumann, personJoël
Ouaknine, and personJames Worrell.
On ranking function synthesis and termination for
polynomial programs. In booktitleCONCUR.

[Pnueli and RosnerPnueli and Rosner1989]
authorpersonAmir Pnueli and personRoni
Rosner. year1989.
On the synthesis of a reactive module. In
booktitlePOPL. pages179–190.

[Polikarpova, Kuraj, and
Solar-LezamaPolikarpova et al2016]
authorpersonNadia Polikarpova, personIvan
Kuraj, and personArmando Solar-Lezama.
Program synthesis from polymorphic refinement
types. In booktitlePLDI. pages522–538.

[Preiner, Niemetz, and BierePreiner
et al2017]
authorpersonMathias Preiner, personAina
Niemetz, and personArmin Biere.
Counterexample-guided model synthesis. In
booktitleTACAS. pages264–280.

authorpersonMihai Putinar.
Positive polynomials on compact semi-algebraic
journalIndiana University Mathematics Journal
volume42, number3 (year1993),

[Reynolds, Deters, Kuncak, Tinelli, and
BarrettReynolds et al2015]
authorpersonAndrew Reynolds, personMorgan
Deters, personViktor Kuncak, personCesare Tinelli,
and personClark Barrett. year2015.
Counterexample-guided quantifier instantiation for
synthesis in SMT. In booktitleCAV.

authorpersonHenry Gordon Rice.
Classes of recursively enumerable sets and their
decision problems.
journalTransactions of the American Mathematical
society volume74, number2 (year1953),

authorpersonEnric Rodríguez-Carbonell.
titleSome programs that need polynomial invariants in
order to be verified.
howpublished<https://www.cs.upc.edu/ erodri/webpage/polynomial_invariants/list.html>.

[Sankaranarayanan, Sipma, and
MannaSankaranarayanan et al2004a]
authorpersonSriram Sankaranarayanan,
personHenny B Sipma, and personZohar Manna.
Constraint-based linear-relations analysis. In
booktitleSAS. pages53–68.

[Sankaranarayanan, Sipma, and
MannaSankaranarayanan et al2004b]
authorpersonSriram Sankaranarayanan,
personHenny B Sipma, and personZohar Manna.
Non-linear loop invariant generation using
Gröbner bases. In booktitlePOPL.

[Shaw, Swartout, and GreenShaw
et al1975]
authorpersonDavid E Shaw, personWilliam R
Swartout, and personC Cordell Green.
Inferring LISP Programs From Examples. In
booktitleIJCAI. pages260–267.

[Si, Lee, Zhang, Albarghouthi, Koutris, and
NaikSi et al2018]
authorpersonXujie Si, personWoosuk Lee,
personRichard Zhang, personAws Albarghouthi,
personParaschos Koutris, and personMayur Naik.
Syntax-guided synthesis of datalog programs. In
booktitleESEC/FSE. pages515–527.

authorpersonDavid Canfield Smith.
booktitlePygmalion: a creative programming
publisherStanford University.

[Sobania, Schweim, and RothlaufSobania
et al2022]
authorpersonDominik Sobania, personDirk
Schweim, and personFranz Rothlauf.
A Comprehensive Survey on Program Synthesis with
Evolutionary Algorithms.
journalIEEE Transactions on Evolutionary
Computation (year2022).

authorpersonArmando Solar-Lezama.
booktitleProgram synthesis by sketching.
publisherUniversity of California, Berkeley.

authorpersonArmando Solar-Lezama.
The sketching approach to program synthesis. In
booktitleAPLAS. pages4–13.

[Srivastava, Gulwani, and FosterSrivastava
et al2010]
authorpersonSaurabh Srivastava, personSumit
Gulwani, and personJeffrey S Foster.
From program verification to program synthesis. In
booktitlePOPL. pages313–326.

[Srivastava, Gulwani, and FosterSrivastava
et al2013]
authorpersonSaurabh Srivastava, personSumit
Gulwani, and personJeffrey S Foster.
Template-based program verification and program
journalInternational Journal on Software Tools for
Technology Transfer volume15, number5
(year2013), pages497–518.

[Takisaka, Oyabu, Urabe, and HasuoTakisaka
et al2018]
authorpersonToru Takisaka, personYuichiro
Oyabu, personNatsuki Urabe, and personIchiro
Hasuo. year2018.
Ranking and repulsing supermartingales for
reachability in probabilistic programs. In booktitleATVA.

authorpersonAlfred Tarski.
A Decision Method for Elementary Algebra and
journalJournal of Symbolic Logic
volume14, number3 (year1949).

authorpersonAnne S Troelstra.
Aspects of constructive mathematics.
In booktitleStudies in Logic and the Foundations of
Mathematics. Vol. volume90. pages973–1052.

authorpersonDavid S Watkins.
booktitleFundamentals of matrix computations.
Vol. volume64.
publisherJohn Wiley & Sons.

§ COMPARISON WITH ROSETTE AND SKETCH

Benchmark     2c|Our Approach     2c|Sketch     2c|Rosette     2c|Direct
1c|Result     1c|Time     1c|Result     1c|Time     1c|Result     1c|Time     1c|Result     1c|Time
ArchimedesPrinciple         0.59         0.33         0.37         0.25
ClosestCubeRoot         15.81         0.46         >43200         >43200
ClosestSquareRoot         3.61         969.17         0.15         >43200
ConsecutiveCubes         56.06         778.04         >43200         >43200
PolynomialApproximation         0.24         0.45         0.12         0.14
PositivityEnforcement         0.19         0.30         >43200         0.09
SlidingBody         0.62         0.47         < 0.01         0.13
SquareCompletion         0.19         0.29         0.10         0.09
SquareRootApproximation          20.53         0.39         >43200         >43200
Cohendiv         26.11         6.297         0.11         0.29
MannadivCarre         4.07         >43200         0.25         316.5
MannadivCube         4.03         0.33         >43200         92.9
MannadivInd         5.25         0.408         0.11         0.18
Wensley         52.75         0.44         >43200         >43200
Euclidex2         29.68         >43200         0.34         >43200
LCM1         10.99         >43200         >43200          87.29
Fermat1         17.95         >43200         0.17         117.14
Fermat2         9.30         >43200         >43200         0.31
Petter3         2.51         0.493         >43200         124.51
Petter5         1.84         2.72         >43200         117.33

Comparison of our approach with previous synthesis methods. All times are in seconds.

Table <ref> provides a comparison between our approach and the state-of-the-art synthesis tools Sketch and Rosette, as well as a direct encoding to SMT that was passed to Z3. The rest of this document contains the details of programs output by these tools. All experiments were run on the same machine as in our original submission with a time limit of 12 hours per instance. Note that Sketch is not sound for real variables and wrong answers are expected. However, in most cases Sketch's output is just setting everything to 0, which is wrong even for integer variables.

§ MORE EXAMPLES FROM SECTION <REF>

@prog: PositivityEnforcement
@real: $x,n$;
@pre: $x \le 100$;
$x$ = $n^2 -5 \cdot n$;
$n$ = $\Hole{5 \cdot n + 1}$;
$x$ = $x + n$;
@post: $x \geq 0$;

@prog: SquareCompletion
@real: $x,n$;
@pre: $x \le 100$;
$x$ = $n^2 - 5\cdot n$;
$n$ = $\Hole{\frac{25}{4}}$;
$x$ = $x + n$;
@post: $x \geq 0$;

@prog: PolynomialApproximation
@real $t,x,e,g_1,g_2,n$;
@macro $f_1 = x^3$;
@macro $f_2 = x^3+x$;
@macro $f_3 = \Hole{5+x+x^3}$;
@pre: $1 \geq 0$;
$t = 0$;
$e = 5$;
if ($x \le t$)

$g_1 = f_3 - f_1$;
$n = g_1$;

$g_2 = f_3 - f_2$;
$n = g_2$;

@post: $n \le e$;

@prog: SlidingBody
@real: $move,e,t,g,nus,nuk,sin,cos,l$;
@pre: $t \geq 0 \wedge nus \geq nuk$;
$g = 9.8$;
$nuk = 0.2$;
$nus = 0.3$;
$sin = 0.5$;
$cos = 0.86$;
$l = 10$;
if ($sin \geq nus \cdot cos$)

$move = 1$;
$t = \Hole{\frac{21}{2}-8.0056060865 \cdot move}$;
$e = t \cdot t \cdot g \cdot sin - g \cdot nuk \cdot cos - 2 \cdot l$;

$move = 0$;
$e = 0$;
$t = \Hole{\frac{21}{2}-8.0056060865 \cdot move}$;

@post: $e = 0 \wedge t \geq 0;$

@prog: ArchimedesPrinciple
@real: $fl,m,l,b,h,\rho,e,x$;
@pre: $x \geq 0 \wedge m \geq 0$;
$l = 10$;
$b = 5$;
$h = 2$;
$\rho = 2$;
if ( $m \leq \rho \cdot l \cdot  b \cdot h$)

$fl = 1$;
$x = \Hole{\frac{fl}{200}}$;
$e = m - \rho \cdot l \cdot b \cdot h \cdot x$;

$e = 0$;
$fl = 0$;
$x = \Hole{\frac{fl}{200}}$;

@post: $e = 0 \wedge 0 \le x \le 1$;

@prog: Cohendiv
@real: $x,y,q,r,dd$;
@pre: $x \geq 0 \wedge y \geq 1$;
$ q = 0$;
$ r = x$;
@invariant: $x = q \cdot y + r \wedge r \geq 0$
while ($r \geq y$)

$d = 1$;
$dd = y$;
@invariant: $dd = y \cdot d \wedge x = q\cdot y + r $
$\wedge r \geq 0 \wedge r \geq y \cdot d$
while($r \geq 2 \cdot dd$)

$d = 2 \cdot d$;
$dd = \Hole{2 \cdot dd}$;

$r = r - dd$;
$q = \Hole{q + d}$;

@post: $1 \geq 0$;

@prog: MannadivCarre
@real $n,y,x,t$;
@pre: $n \geq 0$;
$y = n$;
$x = 0$;
$t = 0$;
@invariant: $x^2 + 2\cdot t + y - n = 0$
while($y\cdot (y-1) \geq 0$)

if($t = x$)

$y = \Hole{y - 1}$;
$t = 0$;
$x = x + 1$;

$y = y - 2$;
$t = t + 1$;

@post: $1 \geq 0$;

@prog: MannadivCube
@real $n,y,x,t$;
@pre: $n \geq 0$;
$y = n$;
$x = 0$;
$t = 0$;
@invariant: $x^3 + 3\cdot t + y - n = 0$
while($y \geq 0$)

if($t = x*x$)

$y = \Hole{y - 3\cdot x - 1}$;
$t = 0$;
$x = x + 1$;

$y = y - 3$;
$t = t + 1$;

@post: $1 \geq 0$;

@prog: MannadivInd
@real $x1,x2,y1,y2,y3$;
@pre: $x1 \geq 0 \wedge x2 \geq 0$;
$y1 = 0$;
$y2 = 0$;
$y3 = x1$;
@invariant: $y1\cdot x2  + y2 + y3 = x1$
while($y3 \geq 0$)

if($y2 + 1 = x2$)

$y1 = y1 + 1$;
$y2 = \Hole{0}$;
$y3 = y3 - 1$;

$y2 = \Hole{y2 + 1}$;
$y3 = y3 - 1$;

@post: $1 \geq 0$;

@prog: Wensley
@real: $P,Q,E,a,b,d,y$;
@pre: $Q \geq P \wedge P \geq 0 \wedge E \geq 0$;
$ a = 0$;
$ b = Q/2$;
$ d = 1$;
$ y = 0$;
@invariant: $ a = Q\cdot y \wedge b = Q\cdot d/2$
$\wedge P - d\cdot Q \leq y \cdot Q \wedge y\cdot Q \leq P \wedge Q \geq 0$
while ($d \geq E$)

if($P \leq a + b$)

$b = \Hole{b/2}$;
$d = d/2$;

$a = a + b$;
$y = \Hole{y + d/2}$;
$ b = b/2$;
$ d = d/2$;

@post: $P \geq y \cdot Q \wedge y\cdot Q \geq P - E \cdot Q$;

@prog: Euclidex2
@real $a,b,p,q,r,s,x,y$;
@pre: $x \geq 0 \wedge y \geq 0$;
$a = x$;
$b = y$;
$p = 1$;
$q = 0$;
$r = 0$;
$s = 1$;
@invariant: $a = y \cdot r + x \cdot p \wedge b = x \cdot q + y \cdot s$
while($1 \geq 0$)

if($a \geq b$)

$a = \Hole{a-b}$;
$p = \Hole{p-q}$;
$r = r - s$;

$b = \Hole{b-a}$;
$q = q - p$;
$s = \Hole{s-r}$;

@post: $1 \geq 0$;

@prog: LCM1
@real $a,b,x,y,u,v$;
@pre: $a \geq 0 \wedge b \geq 0$;
$x = a$;
$y = b$;
$u = b$;
$v = 0$;
@invariant: $x\cdot u + y \cdot v = a\cdot b$
while($1 \geq 0$)

@invariant: $x\cdot u + y \cdot v = a\cdot b$
while($1 \geq 0$)

$x = \Hole{x-y}$;
$v = v + u$;

@invariant: $x\cdot u + y \cdot v = a\cdot b$
while($1 \geq 0$)

$y = y - x$;
$u = \Hole{u+v}$;

@post: $1 \geq 0$;

@prog: Fermat1
@real $N,R,u,v,r$;
@pre: $N \geq 1 \wedge (R-1)^2  \leq N-1 \wedge N \leq R^2$;
$u = 2\cdot R + 1$;
$v = 1$;
$r = R^2 - N$;
@invariant: $N\geq 1 \wedge 4 \cdot (N + r) = u^2- v^2 - 2 \cdot u + 2 \cdot v$
while($1\geq 0$)

@invariant: $N\geq 1 \wedge 4 \cdot (N + r) = u^2- v^2 - 2 \cdot u + 2 \cdot v$
while($1 \geq 0$)

$r = r - v$;
$v = \Hole{v + 2}$;

@invariant: $N\geq 1 \wedge 4 \cdot (N + r) = u^2- v^2 - 2 \cdot u + 2 \cdot v$
while($1 \geq 0$)

$r = \Hole{r+u}$;
$u = u + 2$;

@post: $1 \geq 0$;

@prog: Fermat2
@real $N,R,u,v,r$;
@pre: $N \geq 1 \wedge (R-1)^2  \leq N-1 \wedge N \leq R^2$;
$u = 2\cdot R + 1$;
$v = 1$;
$r = R^2 - N$;
@invariant: $N\geq 1 \wedge 4 \cdot (N + r) = u^2- v^2 - 2 \cdot u + 2 \cdot v$
while($r\geq 1$)

if($r \geq 0$)

$r = r - v$;
$v = v + 2$;

$r = \Hole{r+u}$;
$u = u + 2$;

@post: $1 \geq 0$;

@prog: Petter3
@real $n,k,N,s$;
@pre: $N \geq 0$;
$s = 0$;
$n = 0$;
@invariant: $s = 0.25 \cdot n^2 \cdot (n+1)^2 \wedge n \leq N + 1$
while($n \leq N$)

$s = s + \Hole{nˆ3}$;
$n = n + 1$;

@post: $1 \geq 0$;

@prog: Petter5
@real $n,k,N,s$;
@pre: $N \geq 0$;
$s = 0$;
$n = 0$;
@invariant: $s = n^6 \wedge n \leq N + 1$
while($n \leq N$)

$s = s + \Hole{6 \cdot n^5 + 15 \cdot n^4 + 20 \cdot n^3 + 15 \cdot n^2 + 6 \cdot n + 1}$;
$n = n + 1$;

@post: $1 \geq 0$;

§.§ Examples used in Section <ref>

These examples are motivated from the benchmarks provided in <cit.>.