Exact Real Arithmetic with Perturbation Analysis and Proof of Correctness

In this article, we consider a simple representation for real numbers and propose top-down procedures to approximate various algebraic and transcendental operations with arbitrary precision. Detailed algorithms and proofs are provided to guarantee the correctness of the approximations. Moreover, we develop and apply a perturbation analysis method to show that our approximation procedures only recompute expressions when unavoidable. In the last decade, various theories have been developed and implemented to realize real computations with arbitrary precision. Proof of correctness for existing approaches typically consider basic algebraic operations, whereas detailed arguments about transcendental operations are not available. Another important observation is that in each approach some expressions might require iterative computations to guarantee the desired precision. However, no formal reasoning is provided to prove that such iterative calculations are essential in the approximation procedures. In our approximations of real functions, we explicitly relate the precision of the inputs to the guaranteed precision of the output, provide full proofs and a precise analysis of the necessity of iterations

Delta-Complete Decision Procedures for Satisfiability over the Reals

We introduce the notion of "\delta-complete decision procedures" for solving SMT problems over the real numbers, with the aim of handling a wide range of nonlinear functions including transcendental functions and solutions of Lipschitz-continuous ODEs. Given an SMT problem \varphi and a positive rational number \delta, a \delta-complete decision procedure determines either that \varphi is unsatisfiable, or that the "\delta-weakening" of \varphi is satisfiable. Here, the \delta-weakening of \varphi is a variant of \varphi that allows \delta-bounded numerical perturbations on \varphi. We prove the existence of \delta-complete decision procedures for bounded SMT over reals with functions mentioned above. For functions in Type 2 complexity class C, under mild assumptions, the bounded \delta-SMT problem is in NP^C. \delta-Complete decision procedures can exploit scalable numerical methods for handling nonlinearity, and we propose to use this notion as an ideal requirement for numerically-driven decision procedures. As a concrete example, we formally analyze the DPLL framework, which integrates Interval Constraint Propagation (ICP) in DPLL(T), and establish necessary and sufficient conditions for its \delta-completeness. We discuss practical applications of \delta-complete decision procedures for correctness-critical applications including formal verification and theorem proving.Comment: A shorter version appears in IJCAR 201

Efficient algorithms for computing the Euler-Poincar\'e characteristic of symmetric semi-algebraic sets

Let $\mathrm{R}$ be a real closed field and $\mathrm{D} \subset \mathrm{R}$ an ordered domain. We consider the algorithmic problem of computing the generalized Euler-Poincar\'e characteristic of real algebraic as well as semi-algebraic subsets of $\mathrm{R}^k$, which are defined by symmetric polynomials with coefficients in $\mathrm{D}$. We give algorithms for computing the generalized Euler-Poincar\'e characteristic of such sets, whose complexities measured by the number the number of arithmetic operations in $\mathrm{D}$, are polynomially bounded in terms of $k$ and the number of polynomials in the input, assuming that the degrees of the input polynomials are bounded by a constant. This is in contrast to the best complexity of the known algorithms for the same problems in the non-symmetric situation, which are singly exponential. This singly exponential complexity for the latter problem is unlikely to be improved because of hardness result ($\#\mathbf{P}$-hardness) coming from discrete complexity theory.Comment: 29 pages, 1 Figure. arXiv admin note: substantial text overlap with arXiv:1312.658

Complete Subdivision Algorithms, II: Isotopic Meshing of Singular Algebraic Curves

Given a real valued function f(X,Y), a box region B_0 in R^2 and a positive epsilon, we want to compute an epsilon-isotopic polygonal approximation to the restriction of the curve S=f^{-1}(0)={p in R^2: f(p)=0} to B_0. We focus on subdivision algorithms because of their adaptive complexity and ease of implementation. Plantinga and Vegter gave a numerical subdivision algorithm that is exact when the curve S is bounded and non-singular. They used a computational model that relied only on function evaluation and interval arithmetic. We generalize their algorithm to any bounded (but possibly non-simply connected) region that does not contain singularities of S. With this generalization as a subroutine, we provide a method to detect isolated algebraic singularities and their branching degree. This appears to be the first complete purely numerical method to compute isotopic approximations of algebraic curves with isolated singularities

Computing the Lambert W function in arbitrary-precision complex interval arithmetic

We describe an algorithm to evaluate all the complex branches of the Lambert W function with rigorous error bounds in interval arithmetic, which has been implemented in the Arb library. The classic 1996 paper on the Lambert W function by Corless et al. provides a thorough but partly heuristic numerical analysis which needs to be complemented with some explicit inequalities and practical observations about managing precision and branch cuts.Comment: 16 pages, 4 figure

A probabilistic approach to reducing the algebraic complexity of computing Delaunay triangulations

Computing Delaunay triangulations in $\mathbb{R}^d$ involves evaluating the so-called in\_sphere predicate that determines if a point $x$ lies inside, on or outside the sphere circumscribing $d+1$ points $p_0,\ldots ,p_d$. This predicate reduces to evaluating the sign of a multivariate polynomial of degree $d+2$ in the coordinates of the points $x, \, p_0,\, \ldots,\, p_d$. Despite much progress on exact geometric computing, the fact that the degree of the polynomial increases with $d$ makes the evaluation of the sign of such a polynomial problematic except in very low dimensions. In this paper, we propose a new approach that is based on the witness complex, a weak form of the Delaunay complex introduced by Carlsson and de Silva. The witness complex $\mathrm{Wit} (L,W)$ is defined from two sets $L$ and $W$ in some metric space $X$: a finite set of points $L$ on which the complex is built, and a set $W$ of witnesses that serves as an approximation of $X$. A fundamental result of de Silva states that $\mathrm{Wit}(L,W)=\mathrm{Del} (L)$ if $W=X=\mathbb{R}^d$. In this paper, we give conditions on $L$ that ensure that the witness complex and the Delaunay triangulation coincide when $W$ is a finite set, and we introduce a new perturbation scheme to compute a perturbed set $L'$ close to $L$ such that $\mathrm{Del} (L')= \mathrm{wit} (L', W)$. Our perturbation algorithm is a geometric application of the Moser-Tardos constructive proof of the Lov\'asz local lemma. The only numerical operations we use are (squared) distance comparisons (i.e., predicates of degree 2). The time-complexity of the algorithm is sublinear in $|W|$. Interestingly, although the algorithm does not compute any measure of simplex quality, a lower bound on the thickness of the output simplices can be guaranteed.Comment: 24 page