Showing papers in "Siam Review in 2005"

SNOPT: An SQP Algorithm for Large-Scale Constrained Optimization

Abstract: Sequential quadratic programming (SQP) methods have proved highly effective for solving constrained optimization problems with smooth nonlinear functions in the objective and constraints. Here we consider problems with general inequality constraints (linear and nonlinear). We assume that first derivatives are available and that the constraint gradients are sparse. Second derivatives are assumed to be unavailable or too expensive to calculate. We discuss an SQP algorithm that uses a smooth augmented Lagrangian merit function and makes explicit provision for infeasibility in the original problem and the QP subproblems. The Hessian of the Lagrangian is approximated using a limited-memory quasi-Newton method. SNOPT is a particular implementation that uses a reduced-Hessian semidefinite QP solver (SQOPT) for the QP subproblems. It is designed for problems with many thousands of constraints and variables but is best suited for problems with a moderate number of degrees of freedom (say, up to 2000). Numerical results are given for most of the CUTEr and COPS test collections (about 1020 examples of all sizes up to 40000 constraints and variables, and up to 20000 degrees of freedom).

A Survey of Eigenvector Methods for Web Information Retrieval

Abstract: Web information retrieval is significantly more challenging than traditional well-controlled, small document collection information retrieval. One main difference between traditional information retrieval and Web information retrieval is the Web's hyperlink structure. This structure has been exploited by several of today's leading Web search engines, particularly Google and Teoma. In this survey paper, we focus on Web information retrieval methods that use eigenvector computations, presenting the three popular methods of HITS, PageRank, and SALSA.

Propagation, Observation, and Control of Waves Approximated by Finite Difference Methods

Abstract: This paper surveys several topics related to the observation and control of wave propagation phenomena modeled by finite difference methods. The main focus is on the property of observability, corresponding to the question of whether the total energy of solutions can be estimated from partial measurements on a subregion of the domain or boundary. The mathematically equivalent property of controllability corresponds to the question of whether wave propagation behavior can be controlled using forcing terms on that subregion, as is often desired in engineering applications. Observability/controllability of the continuous wave equation is well understood for the scalar linear constant coefficient case that is the focus of this paper. However, when the wave equation is discretized by finite difference methods, the control for the discretized model does not necessarily yield a good approximation to the control for the original continuous problem. In other words, the classical convergence (consistency + stability) property of a numerical scheme does not suffice to guarantee its suitability for providing good approximations to the controls that might be needed in applications. Observability/controllability may be lost under numerical discretization as the mesh size tends to zero due to the existence of high-frequency spurious solutions for which the group velocity vanishes. This phenomenon is analyzed and several remedies are suggested, including filtering, Tychonoff regularization, multigrid methods, and mixed finite element methods. We also briefly discuss these issues for the heat, beam, and Schrodinger equations to illustrate that diffusive and dispersive effects may help to retain the observability/controllability properties at the discrete level. We conclude with a list of open problems and future subjects for research.

What Color Is Your Jacobian? Graph Coloring for Computing Derivatives

Abstract: Graph coloring has been employed since the 1980s to efficiently compute sparse Jacobian and Hessian matrices using either finite differences or automatic differentiation. Several coloring problems occur in this context, depending on whether the matrix is a Jacobian or a Hessian, and on the specifics of the computational techniques employed. We consider eight variant vertex coloring problems here. This article begins with a gentle introduction to the problem of computing a sparse Jacobian, followed by an overview of the historical development of the research area. Then we present a unifying framework for the graph models of the variant matrix estimation problems. The framework is based upon the viewpoint that a partition of a matrix into structurally orthogonal groups of columns corresponds to distance-2 coloring an appropriate graph representation. The unified framework helps integrate earlier work and leads to fresh insights; enables the design of more efficient algorithms for many problems; leads to new algorithms for others; and eases the task of building graph models for new problems. We report computational results on two of the coloring problems to support our claims. Most of the methods for these problems treat a column or a row of a matrix as an atomic entity, and partition the columns or rows (or both). A brief review of methods that do not fit these criteria is provided. We also discuss results in discrete mathematics and theoretical computer science that intersect with the topics considered here.

A Two-Dimensional Data Distribution Method for Parallel Sparse Matrix-Vector Multiplication

Abstract: A new method is presented for distributing data in sparse matrix-vector multiplication. The method is two-dimensional, tries to minimize the true communication volume, and also tries to spread the computation and communication work evenly over the processors. The method starts with a recursive bipartitioning of the sparse matrix, each time splitting a rectangular matrix into two parts with a nearly equal number of nonzeros. The communication volume caused by the split is minimized. After the matrix partitioning, the input and output vectors are partitioned with the objective of minimizing the maximum communication volume per processor. Experimental results of our implementation, Mondriaan, for a set of sparse test matrices show a reduction in communication volume compared to one-dimensional methods, and in general a good balance in the communication work. Experimental timings of an actual parallel sparse matrix-vector multiplication on an SGI Origin 3800 computer show that a sufficiently large reduction in communication volume leads to savings in execution time.

On the Finite Element Solution of the Pure Neumann Problem

Abstract: This paper considers the finite element approximation and algebraic solution of the pure Neumann problem. Our goal is to present a concise variational framework for the finite element solution of the Neumann problem that focuses on the interplay between the algebraic and variational problems. While many of the results that stem from our analysis are known by some experts, they are seldom derived in a rigorous fashion and remain part of numerical folklore. As a result, this knowledge is not accessible (or appreciated) by many practitioners---both novices and experts---in one source. Our paper contributes a simple, yet insightful link between the continuous and algebraic variational forms that will prove useful.

Reviving the Method of Particular Solutions

Abstract: Fox, Henrici, and Moler made famous a "method of particular solutions" for computing eigenvalues and eigenmodes of the Laplacian in planar regions such as polygons. We explain why their formulation of this method breaks down when applied to regions that are insufficiently simple and propose a modification that avoids these difficulties. The crucial changes are to introduce points in the interior of the region as well as on the boundary and to minimize a subspace angle rather than just a singular value or a determinant. Similar methods may be used to improve other "mesh-free" algorithms for a variety of computational problems.

The Effect of Dispersal Patterns on Stream Populations

Abstract: Individuals in streams are constantly subject to predominantly unidirectional flow. The question of how these populations can persist in upper stream reaches is known as the "drift paradox." We employ a general mechanistic movement-model framework and derive dispersal kernels for this situation. We derive thin- as well as fat-tailed kernels. We then introduce population dynamics and analyze the resulting integrodifferential equation. In particular, we study how the critical domain size and the invasion speed depend on the velocity of the stream flow. We give exact conditions under which a population can persist in a finite domain in the presence of stream flow, as well as conditions under which a population can spread against the direction of the flow. We find a critical stream velocity above which a population cannot persist in an arbitrarily large domain. At exactly the same stream velocity, the invasion speed against the flow becomes zero; for larger velocities, the population retreats with the flow.

Canonical Forms for Hermitian Matrix Pairs under Strict Equivalence and Congruence

Abstract: After a brief historical review and an account of the canonical forms attributed to Jordan and Kronecker, a systematic development is made of the simultaneous reduction of pairs of quadratic forms over the complex numbers and over the reals These reductions are by strict equivalence and by congruence, and essentially complete proofs are presented Some closely related results which can be derived from the canonical forms are also included They concern simultaneous diagonalization, a new criterion for the existence of positive definite linear combinations of a pair of Hermitian matrices, and the canonical structures of matrices which are self-adjoint in an indefinite inner product

Adaptive Smoothed Aggregation ($\alpha$SA) Multigrid

Abstract: Substantial effort has been focused over the last two decades on developing multilevel iterative methods capable of solving the large linear systems encountered in engineering practice. These systems often arise from discretizing partial differential equations over unstructured meshes, and the particular parameters or geometry of the physical problem being discretized may be unavailable to the solver. Algebraic multigrid (AMG) and multilevel domain decomposition methods of algebraic type have been of particular interest in this context because of their promises of optimal performance without the need for explicit knowledge of the problem geometry. These methods construct a hierarchy of coarse problems based on the linear system itself and on certain assumptions about the smooth components of the error. For smoothed aggregation (SA) multigrid methods applied to discretizations of elliptic problems, these assumptions typically consist of knowledge of the near-kernel or near-nullspace of the weak form. This paper introduces an extension of the SA method in which good convergence properties are achieved in situations where explicit knowledge of the near-kernel components is unavailable. This extension is accomplished in an adaptive process that uses the method itself to determine near-kernel components and adjusts the coarsening processes accordingly.

A Simply Stabilized Running Model

Abstract: The spring-loaded inverted pendulum (SLIP), or monopedal hopper, is an archetypal model for running in numerous animal species. Although locomotion is generally considered a complex task requiring sophisticated control strategies to account for coordination and stability, we show that stable gaits can be found in the SLIP with both linear and "air" springs, controlled by a simple fixed-leg reset policy. We first derive touchdown-to-touchdown Poincare maps under the common assumption of negligible gravitational effects during the stance phase. We subsequently include and assess these effects and briefly consider coupling to pitching motions. We investigate the domains of attraction of symmetric periodic gaits and bifurcations from the branches of stable gaits in terms of nondimensional parameters.

Bouncing Ball Modes and Quantum Chaos

Abstract: Quantum ergodicity for classically chaotic systems has been studied extensively both theoretically and experimentally in mathematics and physics. Despite this long tradition we are able to prove a new rigorous result using only elementary calculus. In the case of the famous Bunimovich stadium shown in Figure 1, we prove that the wave functions have to spread to any neighborhood of the wings.

A Fast and Stable Solution Method for the Radiative Transfer Problem

Abstract: Radiative transfer theory considers radiation in turbid media and is used in a wide range of applications. This paper outlines a problem formulation and a solution method for the radiative transfer problem in multilayer scattering and absorbing media using discrete ordinate model geometry. A selection of different steps is brought together. The main contribution here is the synthesis of these steps, all of which have been used in different areas, but never all together in one method. First, all necessary steps to get a numerically stable solution procedure are treated, and then methods are introduced to increase the speed by a factor of several thousand. This includes methods for handling strongly forward-scattering media. The method is shown to be unconditionally stable, though the problem was previously considered numerically intractable.

On the Occurrence of Superlinear Convergence of Exact and Inexact Krylov Subspace Methods

Abstract: Krylov subspace methods often exhibit superlinear convergence. We present a general analytic model which describes this superlinear convergence, when it occurs. We take an invariant subspace approach, so that our results apply also to inexact methods, and to nondiagonalizable matrices. Thus, we provide a unified treatment of the superlinear convergence of GMRES, conjugate gradients, block versions of these, and inexact subspace methods. Numerical experiments illustrate the bounds obtained.

Convergence of Polynomial Restart Krylov Methods for Eigenvalue Computations

Abstract: Krylov subspace methods have led to reliable and effective tools for resolving large-scale, non-Hermitian eigenvalue problems. Since practical considerations often limit the dimension of the approximating Krylov subspace, modern algorithms attempt to identify and condense significant components from the current subspace, encode them into a polynomial filter, and then restart the Krylov process with a suitably refined starting vector. In effect, polynomial filters dynamically steer low-dimensional Krylov spaces toward a desired invariant subspace through their action on the starting vector. The spectral complexity of nonnormal matrices makes convergence of these methods difficult to analyze, and these effects are further complicated by the polynomial filter process. The principal object of study in this paper is the angle an approximating Krylov subspace forms with a desired invariant subspace. Convergence analysis is posed in a geometric framework that is robust to eigenvalue ill-conditioning, yet remains relatively uncluttered. The bounds described here suggest that the sensitivity of desired eigenvalues exerts little influence on convergence, provided the associated invariant subspace is well-conditioned; ill-conditioning of unwanted eigenvalues plays an essential role. This framework also gives insight into the design of effective polynomial filters. Numerical examples illustrate the subtleties that arise when restarting non-Hermitian iterations.

Product Eigenvalue Problems

Abstract: Many eigenvalue problems are most naturally viewed as product eigenvalue problems. The eigenvalues of a matrix $A$ are wanted, but $A$ is not given explicitly. Instead it is presented as a product of several factors: $A = A_{k}A_{k-1}\cdots A_{1}$. Usually more accurate results are obtained by working with the factors rather than forming $A$ explicitly. For example, if we want eigenvalues/vectors of $B^{T}B$, it is better to work directly with $B$ and not compute the product. The intent of this paper is to demonstrate that the product eigenvalue problem is a powerful unifying concept. Diverse examples of eigenvalue problems are discussed and formulated as product eigenvalue problems. For all but a couple of these examples it is shown that the standard algorithms for solving them are instances of a generic $GR$ algorithm applied to a related cyclic matrix.

Dominance of Cyclic Solutions and Challenges in the Scheduling of Robotic Cells

Abstract: We consider the problem of scheduling operations in bufferless robotic cells that produce identical parts. Maximizing the long-term average throughput of parts is an important problem in both theory and practice. We define an appropriate state space required to analyze this problem and show that cyclic schedules which repeat a fixed sequence of robot moves indefinitely are the only ones that need to be considered. For the different classes of robotic cells studied in the literature, we discuss the current state of knowledge with respect to cyclic schedules. Finally, we discuss the importance of two fundamental open problems concerning optimal cyclic schedules, special cases for which these problems have been solved, and attempts to solve the general case.

Hyperasymptotics and the Linear Boundary Layer Problem: Why Asymptotic Series Diverge

Abstract: The simplest problem with boundary layers, $\epsilon^{2} u_{xx} - u = - f(x)$, is used to illustrate (i) why the perturbation series in powers of $\epsilon$ is asymptotic but divergent, (ii) why the optimally truncated expansion is ``superasymptotic'' in the sense that that error is proportional to $\exp(- \mbox{[constant]} / \epsilon)$, and (iii) how to obtain an improved "hyperasymptotic" approximation.

Modeling Basketball Free Throws

Abstract: This paper presents a mathematical model for basketball free throws. It is intended to be a supplement to an existing calculus course and could easily be used as a basis for a calculus project. Students will learn how to apply calculus to model an interesting real-world problem, from problem identification all the way through to interpretation and verification. Along the way we will introduce topics such as optimization (univariate and multiobjective), numerical methods, and differential equations.

k Workers in a Circular Warehouse: A Random Walk on a Circle, without Passing

Abstract: We consider the problem of stochastic flow of multiple particles traveling on a closed loop, with a constraint that particles move without passing. We use a Markov chain description that reduces the problem to a generalized random walk on a hyperplane (with boundaries). By expressing positions via a moving reference frame, the geometry of the no-passing criteria is greatly simplified, with the resultant condition expressible as the coordinate system planes which bound the first orthant. To determine state transition probabilities, we decompose transitions into independent events and construct a digraph representation in which calculating transition probability is reduced to a shortest path determination on the digraph. The resultant decomposition digraph is self-converse, and we exploit that property to establish the necessary symmetries to find the stationary density for the process.

Iterated Impact Dynamics of N -Beads on a Ring

Abstract: When N-beads slide along a frictionless hoop, their collision sequence gives rise to a dynamical system that can be studied via matrix products. It is of general interest to understand the distribution of velocities and the corresponding eigenvalue spectrum that a given collision sequence can produce. We formulate the problem for general N and state some basic theorems regarding the eigenvalues of the collision matrices and their products. The case of three beads of masses m1, m2, m3 is studied in detail. We exploit the fact that each collision sequence can be viewed as a billiard trajectory in a right triangle with non-standard reflection rules. Existence of families of periodic orbits are proven, and orbits that densely fill the triangle are computed. Eigenvalue distributions and position and velocity histograms are computed as a function of the restitution coefficient, both periodic and dense collision sequences are discussed, and a series of conjectures based on computational evidence are formulated. Comparisons are made between the eigenvalue distributions and autocorrelation matrices associated with dense trajectories generated from a chaotic collision sequence and spectra from matrix sequences generated from random orderings, and we describe how the three-bead system could be used as the basis for a random number generating algorithm that is computationally efficient.

Survey and Review

Abstract: Most readers of SIAM Review routinely use public-key cryptography, perhaps without realizing it. Whenever you visit a webpage starting with https and send information, it is encoded so that the information you send, e.g., your credit card number, is transmitted securely. This can be done rapidly using standard cryptography techniques, provided the computers on each end have the same secret key to use for encoding and decoding the information. But how do the computers exchange this key in the first place? For this, a more computationally intensive public-key system is used that allows the encoded message---in this case, the secret key---to be sent in such a way that only the intended recipient can decode it. This can be done using information the web server makes publicly available to anyone who wants to send a secret message. The basis for public-key cryptosystems is a "one-way function" f(x) that can be made publicly available in a form that is easy to evaluate by anyone, or by any web browser, but that is essentially impossible to invert except by someone who possesses additional information. To send a message x the sender computes f(x) and sends this value. Presumably only the intended recipient has the means to invert the function and determine x. Many public-key cryptosystems have been proposed in the past few decades, mostly based on number theoretic or combinatorial problems that are computationally intractable. One of the first and best known is the RSA system, based loosely on the fact that it is easy to multiply two huge primes together but computationally infeasible to factor the product. Other cryptosystems are based on the "discrete log problem": given elements y and g of a finite field, determine an integer x such that y = gx. Neal Koblitz and Alfred J. Menezes, the authors of the following Survey and Review paper, are number theorists who have worked extensively on an approach known as "elliptic curve cryptography." This is based on a discrete log problem in which the finite field is replaced by a group action on integer solutions to certain polynomial equations in x and y. The "product" of two points that lie on the solution curve in the x-y plane is a third point on the curve, defined by a simple geometric construction. The paper provides an introduction to these and other public-key cryptosystems, as well as some insight into the practical aspects of designing secure cryptosystems. As in most branches of applied mathematics, there is often a large gap between what can be rigorously proved and what is needed in practice. In the case of cryptography it can be very difficult to prove that a proposed system is in fact secure. Even if the underlying problem can be shown to be "hard" in some technical sense, e.g., NP-hard, there may be other ways discovered to decipher a given message without solving the general case of the problem, or "side-channel" attacks may be devised that exploit vulnerabilities that have little to do with the underlying mathematical theory. Turning elegant mathematical ideas into viable security systems for the real world is a fascinating challenge in applying mathematics.

A Software Package for Lie Algebraic Computations

Abstract: The paper presents a computer algebra package that facilitates Lie algebraic symbolic computations required in the solution of a variety of problems, such as the solution of right-invariant differential equations evolving on Lie groups. Lie theory is a powerful tool, helpful in the analysis and design of modern nonlinear control laws, nonlinear filters, and the study of particle dynamics. The practical application of Lie theory often results in highly complex symbolic expressions that are difficult to handle efficiently without the aid of a computer software tool. The aim of the package is to facilitate and encourage further research relying on Lie algebraic computations.

Critical Percolation on a Bethe Lattice Revisited

Abstract: We introduce the model of independent percolation on general graphs with emphasis on the Bethe lattice, for which we prove the existence of a phase transition with the precise calculation of the critical value, and derive the value of the critical exponents $\gamma$ and $\b$.