Showing papers in "Foundations of Computational Mathematics in 2008"

Random Sampling of Sparse Trigonometric Polynomials, II. Orthogonal Matching Pursuit versus Basis Pursuit

Abstract: We investigate the problem of reconstructing sparse multivariate trigonometric polynomials from few randomly taken samples by Basis Pursuit and greedy algorithms such as Orthogonal Matching Pursuit (OMP) and Thresholding. While recovery by Basis Pursuit has recently been studied by several authors, we provide theoretical results on the success probability of reconstruction via Thresholding and OMP for both a continuous and a discrete probability model for the sampling points. We present numerical experiments, which indicate that usually Basis Pursuit is significantly slower than greedy algorithms, while the recovery rates are very similar.

Online Gradient Descent Learning Algorithms

Abstract: This paper considers the least-square online gradient descent algorithm in a reproducing kernel Hilbert space (RKHS) without an explicit regularization term. We present a novel capacity independent approach to derive error bounds and convergence results for this algorithm. The essential element in our analysis is the interplay between the generalization error and a weighted cumulative error which we define in the paper. We show that, although the algorithm does not involve an explicit RKHS regularization term, choosing the step sizes appropriately can yield competitive error rates with those in the literature.

Semidefinite Characterization and Computation of Zero-Dimensional Real Radical Ideals

Abstract: For an ideal I⊆ℝ[x] given by a set of generators, a new semidefinite characterization of its real radical I(V ℝ(I)) is presented, provided it is zero-dimensional (even if I is not). Moreover, we propose an algorithm using numerical linear algebra and semidefinite optimization techniques, to compute all (finitely many) points of the real variety V ℝ(I) as well as a set of generators of the real radical ideal. The latter is obtained in the form of a border or Grobner basis. The algorithm is based on moment relaxations and, in contrast to other existing methods, it exploits the real algebraic nature of the problem right from the beginning and avoids the computation of complex components.

Symmetric Exponential Integrators with an Application to the Cubic Schrödinger Equation

Abstract: In this article, we derive and study symmetric exponential integrators. Numerical experiments are performed for the cubic Schrodinger equation and comparisons with classical exponential integrators and other geometric methods are also given. Some of the proposed methods preserve the L 2-norm and/or the energy of the system.

Convergence of the Magnus Series

Abstract: The Magnus series is an infinite series which arises in the study of linear ordinary differential equations. If the series converges, then the matrix exponential of the sum equals the fundamental solution of the differential equation. The question considered in this paper is: When does the series converge? The main result establishes a sufficient condition for convergence, which improves on several earlier results.

The Polynomial Method for Random Matrices

Abstract: We define a class of “algebraic” random matrices. These are random matrices for which the Stieltjes transform of the limiting eigenvalue distribution function is algebraic, i.e., it satisfies a (bivariate) polynomial equation. The Wigner and Wishart matrices whose limiting eigenvalue distributions are given by the semicircle law and the Marcenko–Pastur law are special cases. Algebraicity of a random matrix sequence is shown to act as a certificate of the computability of the limiting eigenvalue density function. The limiting moments of algebraic random matrix sequences, when they exist, are shown to satisfy a finite depth linear recursion so that they may often be efficiently enumerated in closed form. In this article, we develop the mathematics of the polynomial method which allows us to describe the class of algebraic matrices by its generators and map the constructive approach we employ when proving algebraicity into a software implementation that is available for download in the form of the RMTool random matrix “calculator” package. Our characterization of the closure of algebraic probability distributions under free additive and multiplicative convolution operations allows us to simultaneously establish a framework for computational (noncommutative) “free probability” theory. We hope that the tools developed allow researchers to finally harness the power of infinite random matrix theory.

The Beta-Jacobi Matrix Model, the CS Decomposition, and Generalized Singular Value Problems

Abstract: We provide a solution to the β-Jacobi matrix model problem posed by Dumitriu and the first author. The random matrix distribution introduced here, called a matrix model, is related to the model of Killip and Nenciu, but the development is quite different. We start by introducing a new matrix decomposition and an algorithm for computing this decomposition. Then we run the algorithm on a Haar-distributed random matrix to produce the β-Jacobi matrix model. The Jacobi ensemble on ${\Bbb R}^{n}$, parametrized by β > 0, a > -1,and b > -1, is the probability distribution whose density is proportional to $\prod_{i}\lambda_{i}^{({\beta}/{2})(a+1)-1}(1-\lambda_{i})^{({\beta}/{2})(b+1)-1}\prod_{i 0. Observing a connection between Haar measure on the orthogonal (resp., unitary) group and pairs of real (resp., complex) Gaussian matrices, we find a direct connection between multivariate analysis of variance (MANOVA) and the new matrix model.

On the Hopf Algebraic Structure of Lie Group Integrators

Abstract: A commutative but not cocommutative graded Hopf algebra HN, based on ordered (planar) rooted trees, is studied. This Hopf algebra is a generalization of the Hopf algebraic structure of unordered rooted trees HC, developed by Butcher in his study of Runge-Kutta methods and later rediscovered by Connes and Moscovici in the context of noncommutative geometry and by Kreimer where it is used to describe renormalization in quantum field theory. It is shown that HN is naturally obtained from a universal object in a category of noncommutative derivations and, in particular, it forms a foundation for the study of numerical integrators based on noncommutative Lie group actions on a manifold. Recursive and nonrecursive definitions of the coproduct and the antipode are derived. The relationship between HN and four other Hopf algebras is discussed. The dual of HN is a Hopf algebra of Grossman and Larson based on ordered rooted trees. The Hopf algebra HC of Butcher, Connes, and Kreimer is identified as a proper Hopf subalgebra of HN using the image of a tree symmetrization operator. The Hopf algebraic structure of the shuffle algebra HSh is obtained from HN by a quotient construction. The Hopf algebra HP of ordered trees by Foissy differs from HN in the definition of the product (noncommutative concatenation for HP and shuffle for HN) and the definitions of the coproduct and the antipode, however, these are related through the tree symmetrization operator.

A Unifying Local Convergence Result for Newton's Method in Riemannian Manifolds

Abstract: We consider the problem of finding a singularity of a differentiable vector field X defined on a complete Riemannian manifold. We prove a unified result for theexistence and local uniqueness of the solution, and for the local convergence of a Riemannian version of Newton's method. Our approach relies on Kantorovich's majorant principle: under suitable conditions, we construct an auxiliary scalar equation φ(r) = 0 which dominates the original equation X(p) = 0 in the sense that the Riemannian-Newton method for the latter inherits several features of the real Newton method applied to the former. The majorant φ is derived from an adequate radial parametrization of a Lipschitz-type continuity property of the covariant derivative of X, a technique inspired by the previous work of Zabrejko and Nguen on Newton's method in Banach spaces. We show how different specializations of the main result recover Riemannian versions of Kantorovich's theorem and Smale's α-theorem, and, at least partially, the Euclidean self-concordant theory of Nesterov and Nemirovskii. In the specific case of analytic vector fields, we improve recent developments inthis area by Dedieu et al. (IMA J. Numer. Anal., Vol. 23, 2003, pp. 395-419). Some Riemannian techniques used here were previously introduced by Ferreira and Svaiter (J. Complexity, Vol. 18, 2002, pp. 304-329) in the context of Kantorovich's theorem for vector fields with Lipschitz continuous covariant derivatives.

On Smale's 17th Problem: A Probabilistic Positive Solution

Abstract: Smale's 17th Problem asks “Can a zero of n complex polynomial equations in n unknowns be found approximately, on the average [for a suitable probability measure on the space of inputs], in polynomial time with a uniform algorithm?” We present a uniform probabilistic algorithm for this problem and prove that its complexity is polynomial. We thus obtain a partial positive solution to Smale's 17th Problem.

Spectral Semi-discretisations of Weakly Non-linear Wave Equations over Long Times

Abstract: The long-time behaviour of spectral semi-discretisations of weakly non-linear wave equations is analysed. It is shown that the harmonic actions are approximately conserved for the semi-discretised system as well. This permits to prove that the energy of the wave equation along the interpolated semi-discrete solution remains well conserved over long times and close to the Hamiltonian of the semi-discrete equation. Although the momentum is no longer an exact invariant of the semi-discretisation, it is shown to be approximately conserved. All these results are obtained with the technique of modulated Fourier expansions.

Algorithms for Differential Invariants of Symmetry Groups of Differential Equations

Abstract: We develop new computational algorithms, based on the method of equivariant moving frames, for classifying the differential invariants of Lie symmetry pseudo-groups of differential equations and analyzing the structure of the induced differential invariant algebra. The Korteweg-deVries (KdV) and Kadomtsev-Petviashvili (KP) equations serve to illustrate examples. In particular, we deduce the first complete classification of the differential invariants and their syzygies of the KP symmetry pseudo-group.

On the Linear Stability of Splitting Methods

Abstract: A comprehensive linear stability analysis of splitting methods is carried out by means of a 2×2 matrix K(x) with polynomial entries (the stability matrix) and the stability polynomial p(x) (the trace of K(x) divided by two). An algorithm is provided for determining the coefficients of all possible time-reversible splitting schemes for a prescribed stability polynomial. It is shown that p(x) carries essentially all the information needed to construct processed splitting methods for numerically approximating the evolution of linear systems. By conveniently selecting the stability polynomial, new integrators with processing for linear equations are built which are orders of magnitude more efficient than other algorithms previously available.

Optimal Control and Geodesics on Quadratic Matrix Lie Groups

Abstract: The purpose of this paper is to extend the symmetric representation of the rigid body equations from the group SO (n) to other groups. These groups are matrix subgroups of the general linear group that are defined by a quadratic matrix identity. Their corresponding Lie algebras include several classical semisimple matrix Lie algebras. The approach is to start with an optimal control problem on these groups that generates geodesics for a left-invariant metric. Earlier work by Bloch, Crouch, Marsden, and Ratiu defines the symmetric representation of the rigid body equations, which is obtained by solving the same optimal control problem in the particular case of the Lie group SO (n). This paper generalizes this symmetric representation to a wider class of matrix groups satisfying a certain quadratic matrix identity. We consider the relationship between this symmetric representation of the generalized rigid body equations and the generalized rigid body equations themselves. A discretization of this symmetric representation is constructed making use of the symmetry, which in turn give rise to numerical algorithms to integrate the generalized rigid body equations for the given class of matrix Lie groups.

Computing the First Betti Number of a Semi-Algebraic Set

Abstract: In this paper we describe a singly exponential algorithm for computing the first Betti number of a given semi-algebraic set. Singly exponential algorithms for computing the zeroth Betti number, and the Euler–Poincare characteristic, were known before. No singly exponential algorithm was known for computing any of the individual Betti numbers other than the zeroth one. As a consequence we also obtain algorithms for computing semi-algebraic descriptions of the semi-algebraically connected components of any given real algebraic or semi-algebraic set in singly exponential time, which improves on the complexity of the previously published algorithms for this problem.

Point Counting in Families of Hyperelliptic Curves

Abstract: Let E Γ be a family of hyperelliptic curves defined by $Y^{2}=\bar {Q}(X,\ensuremath {\Gamma })$, where $\bar{Q}$ is defined over a small finite field of odd characteristic. Then with $\ensuremath {\bar {\gamma }}$ in an extension degree n field over this small field, we present a deterministic algorithm for computing the zeta function of the curve $E_{\ensuremath {\bar {\gamma }}}$ by using Dwork deformation in rigid cohomology. The time complexity of the algorithm is $\ensuremath {\mathcal {O}}(n^{2.667})$ and it needs $\ensuremath {\mathcal {O}}(n^{2.5})$ bits of memory. A slight adaptation requires only $\ensuremath {\mathcal {O}}(n^{2})$ space, but costs time $\ensuremath {\widetilde {\mathcal {O}}}(n^{3})$. An implementation of this last result turns out to be quite efficient for n big enough.

Meshfree Thinning of 3D Point Clouds

Abstract: An efficient data reduction scheme for the simplification of a surface given by a large set X of 3D point-samples is proposed. The data reduction relies on a recursive point removal algorithm, termed thinning, which outputs a data hierarchy of point-samples for multiresolution surface approximation. The thinning algorithm works with a point removal criterion, which measures the significances of the points in their local neighbourhoods, and which removes a least significant point at each step. For any point x in the current point set Y⊂X, its significance reflects the approximation quality of a local surface reconstructed from neighbouring points in Y. The local surface reconstruction is done over an estimated tangent plane at x by using radial basis functions. The approximation quality of the surface reconstruction around x is measured by using its maximal deviation from the given point-samples X in a local neighbourhood of x. The resulting thinning algorithm is meshfree, i.e., its performance is solely based upon the geometry of the input 3D surface point-samples, and so it does not require any further topological information, such as point connectivities. Computational details of the thinning algorithm and the required data structures for efficient implementation are explained and its complexity is discussed. Two examples are presented for illustration.

Modular Counting of Rational Points over Finite Fields

Abstract: Let Fq be the finite field of q elements, where q = ph. Let f(x) be a polynomial over Fq in n variables with m nonzero terms. Let N(f) denote the number of solutions of f(x) = 0 with coordinates in Fq. In this paper we give a deterministic algorithm which computes the reduction of N(f) modulo pb in O(n(8m)p(h+b)p) bit operations. This is singly exponential in each of the parameters {h, b, p}, answering affirmatively an open problem proposed by Gopalan, Guruswami, and Lipton.

Explicit Volume-Preserving Splitting Methods for Linear and Quadratic Divergence-Free Vector Fields

Abstract: We present new explicit volume-preserving methods based on splitting for polynomial divergence-free vector fields. The methods can be divided in two classes: methods that distinguish between the diagonal part and the off-diagonal part and methods that do not. For the methods in the first class it is possible to combine different treatments of the diagonal and off-diagonal parts, giving rise to a number of possible combinations.

Exceptional Sets and Fiber Products

Abstract: Exceptional sets where fibers have dimensions higher than the generic fiber dimension are of interest in mathematics and in application areas, such as the theory of overconstrained mechanisms.We show that fiber products promote such sets to become irreducible components, whereupon they can be found using techniques from numerical algebraic geometry for computing the irreducible decomposition. However, such a decomposition may contain components other than the exceptional loci we seek. We show that each irreducible component of the exceptional loci gives rise to a main component in a fiber product of sufficiently high order, and we give procedures for identifying these components. The methods are illustrated by finding the rulings of a general quadricin C3.

Primal Central Paths and Riemannian Distances for Convex Sets

Abstract: In this paper, we study the Riemannian length of the primal central path in a convex set computed with respect to the local metric defined by a self-concordant function. Despite some negative examples, in many important situations, the length of this path is quite close to the length of a geodesic curve. We show that in the case of a bounded convex set endowed with a ν-self-concordant barrier, the length of the central path is within a factor O(ν 1/4) of the length of the shortest geodesic curve.

Errata for Computing the Top Betti Numbers of Semialgebraic Sets Defined by Quadratic Inequalities in Polynomial Time

Abstract: For any l>0, we present an algorithm which takes as input a semi-algebraic set, S, defined by P 1≤0,…,P s ≤0, where each P i ∈R[X 1,…,X k ] has degree≤2, and computes the top l Betti numbers of S, b k−1(S),…,b k−l (S), in polynomial time. The complexity of the algorithm, stated more precisely, is $\sum_{i=0}^{\ell+2}{s\choose i}k^{2^{o(\min(\ell,s))}}$. For fixed l, the complexity of the algorithm can be expressed as $s^{\ell+2}+k^{2^{O(\ell)}}$, which is polynomial in the input parameters s and k. To our knowledge this is the first polynomial time algorithm for computing nontrivial topological invariants of semialgebraic sets in R k defined by polynomial inequalities, where the number of inequalities is not fixed and the polynomials are allowed to have degree greater than one. For fixed s, we obtain, by letting l=k, an algorithm for computing all the Betti numbers of S whose complexity is $k^{2^{O(S)}}$.

Coefficient Quantization in Banach Spaces

Abstract: Let (e i ) be a dictionary for a separable infinite-dimensional Banach space X. We consider the problem of approximation by linear combinations of dictionary elements with quantized coefficients drawn usually from a ‘finite alphabet’. We investigate several approximation properties of this type and connect them to the Banach space geometry of X. The existence of a total minimal system with one of these properties, namely the coefficient quantization property, is shown to be equivalent to X containing c 0. We also show that, for every e>0, the unit ball of every separable infinite-dimensional Banach space X contains a dictionary (x i ) such that the additive group generated by (x i ) is (3+e)−1-separated and 1/3-dense in X.

On Spherical Averages of Radial Basis Functions

Abstract: A radial basis function (RBF) has the general form $$s(x)=\sum_{k=1}^{n}a_{k}\phi(x-b_{k}),\quad x\in\mathbb{R}^{d},$$ where the coefficients a 1,…,a n are real numbers, the points, or centres, b 1,…,b n lie in ℝ d , and φ:ℝ d →ℝ is a radially symmetric function. Such approximants are highly useful and enjoy rich theoretical properties; see, for instance (Buhmann, Radial Basis Functions: Theory and Implementations, [2003]; Fasshauer, Meshfree Approximation Methods with Matlab, [2007]; Light and Cheney, A Course in Approximation Theory, [2000]; or Wendland, Scattered Data Approximation, [2004]). The important special case of polyharmonic splines results when φ is the fundamental solution of the iterated Laplacian operator, and this class includes the Euclidean norm φ(x)=‖x‖ when d is an odd positive integer, the thin plate spline φ(x)=‖x‖2log ‖x‖ when d is an even positive integer, and univariate splines. Now B-splines generate a compactly supported basis for univariate spline spaces, but an analyticity argument implies that a nontrivial polyharmonic spline generated by (1.1) cannot be compactly supported when d>1. However, a pioneering paper of Jackson (Constr. Approx. 4:243–264, [1988]) established that the spherical average of a radial basis function generated by the Euclidean norm can be compactly supported when the centres and coefficients satisfy certain moment conditions; Jackson then used this compactly supported spherical average to construct approximate identities, with which he was then able to derive some of the earliest uniform convergence results for a class of radial basis functions. Our work extends this earlier analysis, but our technique is entirely novel, and applies to all polyharmonic splines. Furthermore, we observe that the technique provides yet another way to generate compactly supported, radially symmetric, positive definite functions. Specifically, we find that the spherical averaging operator commutes with the Fourier transform operator, and we are then able to identify Fourier transforms of compactly supported functions using the Paley–Wiener theorem. Furthermore, the use of Haar measure on compact Lie groups would not have occurred without frequent exposure to Iserles’s study of geometric integration.