SIAM Digital Library
 
 
 

SIAM J. on Scientific Computing

Year Range: 

All Journal content published prior to 1997 is part of LOCUS.

Search Issue | RSS Feeds RSS
Previous Issue Next Issue

2008

Volume 30, Issue 2, pp. 549-1103


A Level Set–Boundary Element Method for the Simulation of Underwater Bubble Dynamics

K. L. Tan, B. C. Khoo, and J. K. White

SIAM J. Sci. Comput. 30, pp. 549-571 (23 pages)

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
The front tracking–boundary element method, though commonly used to simulate underwater bubble dynamics, is not quite well suited for treating the complex topological changes that can occur. In this paper, we present a level set–boundary element method that overcomes this deficiency. This new technique combines the level set method to capture the bubble surface motions and the boundary element method to compute the surface velocities. We propose a level set-type potential function to capture the surface velocity potential, and we formulate a governing equation to evolve this function such that the functional value at the bubble surface is the surface velocity potential consistent with the unsteady Bernoulli equation. A convergence analysis performed using the Rayleigh bubble as reference establishes our implemented procedure to be second-order accurate. We also demonstrate the ease with which our method handles complex topological changes using axisymmetrical examples of bubble merging, bubble splitting, and bubble jet impact. The result for a highly complex problem comprising both jet impact and merging is also presented. Additionally, we develop a more computationally efficient technique that requires computations to be performed only within a narrow band containing the interface.

A Primal-Dual Active Set Algorithm for Three-Dimensional Contact Problems with Coulomb Friction

S. Hüeber, G. Stadler, and B. I. Wohlmuth

SIAM J. Sci. Comput. 30, pp. 572-596 (25 pages) | Cited 4 times

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
In this paper, efficient algorithms for contact problems with Tresca and Coulomb friction in three dimensions are presented and analyzed. The numerical approximation is based on mortar methods for nonconforming meshes with dual Lagrange multipliers. Using a nonsmooth complementarity function for the three-dimensional friction conditions, a primal-dual active set algorithm is derived. The method determines active contact and friction nodes and, at the same time, resolves the additional nonlinearity originating from sliding nodes. No regularization and no penalization are applied, and superlinear convergence can be observed locally. In combination with a multigrid method, it defines a robust and fast strategy for contact problems with Tresca or Coulomb friction. The efficiency and flexibility of the method is illustrated by several numerical examples.

Stochastic Lie Group Integrators

Simon J. A. Malham and Anke Wiese

SIAM J. Sci. Comput. 30, pp. 597-617 (21 pages) | Cited 3 times

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
We present Lie group integrators for nonlinear stochastic differential equations with noncommutative vector fields whose solution evolves on a smooth finite-dimensional manifold. Given a Lie group action that generates transport along the manifold, we pull back the stochastic flow on the manifold to the Lie group via the action, and subsequently pull back the flow to the corresponding Lie algebra via the exponential map. We construct an approximation to the stochastic flow in the Lie algebra via closed operations, and then push back to the Lie group and then to the manifold, thus ensuring that our approximation lies in the manifold. We call such schemes stochastic Munthe-Kaas methods after their deterministic counterparts. We also present stochastic Lie group integration schemes based on Castell–Gaines methods. These involve using an underlying ordinary differential integrator to approximate the flow generated by a truncated stochastic exponential Lie series. They become stochastic Lie group integrator schemes if we use Munthe-Kaas methods as the underlying ordinary differential integrator. Further, we show that some Castell–Gaines methods are uniformly more accurate than the corresponding stochastic Taylor schemes. Lastly we demonstrate our methods by simulating the dynamics of a free rigid body such as a satellite and an autonomous underwater vehicle both perturbed by two independent multiplicative stochastic noise processes.

Revisiting the Crowding Phenomenon in Schwarz–Christoffel Mapping

L. Banjai

SIAM J. Sci. Comput. 30, pp. 618-636 (19 pages) | Cited 2 times

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
We address the problem of conformally mapping the unit disk to polygons with elongations. The elongations cause the derivative of the conformal map to be exponentially large in some regions. This crowding phenomenon creates difficulties in standard numerical methods for the computation of the conformal map. We make use of the Schwarz–Christoffel representation of the mapping and show that a simple change to the existing algorithms introduced by Trefethen [SIAM J. Sci. Statist. Comput., 1 (1980), pp. 82–102] makes it feasible to accurately compute conformal maps to polygons even in the presence of extreme crowding. For an efficient algorithm it is essential that a good initial guess for the parameters of the Schwarz–Christoffel map be available. A uniformly close initial guess can be obtained from the cross-ratios of certain quadrilaterals, as introduced in the CRDT algorithm of Driscoll and Vavasis [SIAM J. Sci. Comput., 19 (1998), pp. 1783–1803]. We present numerical experiments and compare our algorithms with the CRDT which has been particularly designed to combat crowding.

Continuation of Invariant Subspaces in Large Bifurcation Problems

David Bindel, James Demmel, and Mark Friedman

SIAM J. Sci. Comput. 30, pp. 637-656 (20 pages)

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
We summarize an algorithm for computing a smooth orthonormal basis for an invariant subspace of a parameter-dependent matrix, and describe how to extend it for numerical bifurcation analysis. We adapt the continued subspace to track behavior relevant to bifurcations, and use projection methods to deal with large problems. To test our ideas, we have integrated our code into MATCONT, a program for numerical continuation and bifurcation analysis.

Iterated Hard Shrinkage for Minimization Problems with Sparsity Constraints

Kristian Bredies and Dirk A. Lorenz

SIAM J. Sci. Comput. 30, pp. 657-683 (27 pages) | Cited 3 times

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
A new iterative algorithm for the solution of minimization problems in infinite-dimensional Hilbert spaces which involve sparsity constraints in form of $\ell^{p}$-penalties is proposed. In contrast to the well-known algorithm considered by Daubechies, Defrise, and De Mol, it uses hard instead of soft shrinkage. It is shown that the hard shrinkage algorithm is a special case of the generalized conditional gradient method. Convergence properties of the generalized conditional gradient method with quadratic discrepancy term are analyzed. This leads to strong convergence of the iterates with convergence rates $\mathcal{O}(n^{-1/2})$ and $\mathcal{O}(\lambda^n)$ for $p=1$ and $1 < p \leq 2$, respectively. Numerical experiments on image deblurring, backwards heat conduction, and inverse integration are given.

Multilevel Preconditioning of Two-dimensional Elliptic Problems Discretized by a Class of Discontinuous Galerkin Methods

Johannes K. Kraus and Satyendra K. Tomar

SIAM J. Sci. Comput. 30, pp. 684-706 (23 pages) | Cited 1 time

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
We present optimal-order preconditioners for certain discontinuous Galerkin (DG) finite element discretizations of elliptic boundary value problems. We consider two variants of hierarchical splittings where the underlying finite element meshes are geometrically nested. A specific assembling process is proposed which facilitates the local analysis of the angle between the related subspaces. By applying the corresponding two-level basis transformation recursively, a sequence of algebraic problems is generated that can be associated with a hierarchy of coarse versions of DG approximations of the original problem. New bounds for the constant $\gamma$ in the strengthened Cauchy–Bunyakowski–Schwarz inequality are derived. The presented numerical results support the theoretical analysis and demonstrate the potential of this approach.

Hierarchical Clustering of Massive, High Dimensional Data Sets by Exploiting Ultrametric Embedding

Fionn Murtagh, Geoff Downs, and Pedro Contreras

SIAM J. Sci. Comput. 30, pp. 707-730 (24 pages)

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
Coding of data, usually upstream of data analysis, has crucial implications for the data analysis results. By modifying the data coding—through use of less than full precision in data values—we can aid appreciably the effectiveness and efficiency of the hierarchical clustering. In our first application, this is used to lessen the quantity of data to be hierarchically clustered. The approach is a hybrid one, based on hashing and on the Ward minimum variance agglomerative criterion. In our second application, we derive a hierarchical clustering from relationships between sets of observations, rather than the traditional use of relationships between the observations themselves. This second application uses embedding in a Baire space, or longest common prefix ultrametric space. We compare this second approach, which is of $O(n \log n)$ complexity, to k-means.

An Effective Fluid-Structure Interaction Formulation for Vascular Dynamics by Generalized Robin Conditions

F. Nobile and C. Vergara

SIAM J. Sci. Comput. 30, pp. 731-763 (33 pages) | Cited 5 times

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
In this work we focus on the modeling and numerical simulation of the fluid-structure interaction mechanism in vascular dynamics. We first propose a simple membrane model to describe the deformation of the arterial wall, which is derived from the Koiter shell equations and is applicable to an arbitrary geometry. Secondly, we consider a reformulation of the fluid-structure problem, in which the newly derived membrane model, thanks to its simplicity, is embedded into the fluid equations and will appear as a generalized Robin boundary condition. The original problem is then reduced to the solution of subsequent fluid equations defined on a moving domain and may be achieved with a fluid solver only. We also derive a stability estimate for the resulting numerical scheme. Finally, we propose new outflow absorbing boundary conditions, which are easy to implement and allow us to reduce significantly the spurious pressure wave reflections that typically appear in artificially truncated computational domains. We present several numerical results showing the effectiveness of the proposed approaches.

A Framework for Discrete Integral Transformations I—The Pseudopolar Fourier Transform

A. Averbuch, R. R. Coifman, D. L. Donoho, M. Israeli, and Y. Shkolnisky

SIAM J. Sci. Comput. 30, pp. 764-784 (21 pages) | Cited 6 times

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
The Fourier transform of a continuous function, evaluated at frequencies expressed in polar coordinates, is an important conceptual tool for understanding physical continuum phenomena. An analogous tool, suitable for computations on discrete grids, could be very useful; however, no exact analogue exists in the discrete case. In this paper we present the notion of pseudopolar grid (pp grid) and the pseudopolar Fourier transform (ppFT), which evaluates the discrete Fourier transform at points of the pp grid. The pp grid is a type of concentric-squares grid in which the radial density of squares is twice as high as usual. The pp grid consists of equally spaced samples along rays, where different rays are equally spaced in slope rather than angle. We develop a fast algorithm for the ppFT, with the same complexity order as the Cartesian fast Fourier transform; the algorithm is stable, invertible, requires only one-dimensional operations, and uses no approximate interpolations. We prove that the ppFT is invertible and develop two algorithms for its inversion: iterative and direct, both with complexity $O(n^{2}\log{n})$, where $n \times n$ is the size of the reconstructed image. The iterative algorithm applies conjugate gradients to the Gram operator of the ppFT. Since the transform is ill-conditioned, we introduce a preconditioner, which significantly accelerates the convergence. The direct inversion algorithm utilizes the special frequency domain structure of the transform in two steps. First, it resamples the pp grid to a Cartesian frequency grid and then recovers the image from the Cartesian frequency grid.

A Framework for Discrete Integral Transformations II—The 2D Discrete Radon Transform

A. Averbuch, R. R. Coifman, D. L. Donoho, M. Israeli, Y. Shkolnisky, and I. Sedelnikov

SIAM J. Sci. Comput. 30, pp. 785-803 (19 pages) | Cited 4 times

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
Although naturally at the heart of many fundamental physical computations, and potentially useful in many important image processing tasks, the Radon transform lacks a coherent discrete definition for two-dimensional (2D) discrete images which is algebraically exact, invertible, and rapidly computable. We define a notion of 2D discrete Radon transforms for 2D discrete images, which is based on summation along lines of absolute slope less than 1. Values at nongrid locations are defined using trigonometric interpolation on a zero-padded grid. Our definition is shown to be geometrically faithful: the summation avoids wrap-around effects. Our proposal uses a special collection of lines in $\mathbb{R}^{2}$ for which the transform is rapidly computable and invertible. We describe a fast algorithm using $O(N\log{N})$ operations, where $N =n^{2}$ is the number of pixels in the image. The fast algorithm exploits a discrete projection-slice theorem, which associates the discrete Radon transform with the pseudopolar Fourier transform [A. Averbuch et al., SIAM J. Sci. Comput., 30 (2008), pp. 764–784]. Our definition for discrete images converges to a natural continuous counterpart with increasing refinement.

Fast Iterative Schemes for Nonsymmetric Algebraic Riccati Equations Arising from Transport Theory

Zhong-Zhi Bai, Yong-Hua Gao, and Lin-Zhang Lu

SIAM J. Sci. Comput. 30, pp. 804-818 (15 pages)

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
For a class of nonsymmetric algebraic Riccati equations arising from transport theory, Lu [SIAM J. Matrix Anal. Appl., 26 (2005), pp. 679–685] reformulated them into a couple of fixed-point equations in a vector form and constructed an iterative method for finding their minimal positive solutions. In this paper, based on this fixed-point equation, we propose several fast and effective iterative methods for solving the above-mentioned special type of nonsymmetric algebraic Riccati equations. These iterative methods can be categorized into the class of nonlinear block relaxation iterations and can also be considered as accelerated variants of the known fixed-point iteration according to the principle of using the currently available information as promptly as possible. The monotone convergence theorems about these new iterative methods are established, and the sharper bounds about the solutions of the nonsymmetric algebraic Riccati equations are derived. Numerical implementations show that the new iterative methods are feasible and effective solvers for the nonsymmetric algebraic Riccati equations arising from transport theory.

FA-SART: A Frequency-Adaptive Algorithm in Pinhole SPECT Tomography

Vincent Israel-Jost

SIAM J. Sci. Comput. 30, pp. 819-836 (18 pages)

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
As a solution to the frequently reported problem of the high frequencies reconstructed sometimes many iterations after the low frequencies in tomography, we proposed in a recent paper [V. Israel-Jost, Ph. Choquet, and A. Constantinesco, Int. J. Biomed. Imaging, 2006, paper 34043.] a general adaptation that can apply to many algorithms: the incomplete backprojection. We stressed the fact that this type of frequency adapted (FA) algorithm works only when the linear problem of tomography to be solved involves deblurring, i.e., when the point spread functions (PSFs) associated with each voxel are large. Thus, depending on the desired level of detail, a parameter is set to achieve the reconstruction within a small number of iterations. The aim of this paper is to give a mathematical analysis of both the acceleration obtained and the convergence of an FA algorithm derived from the simultaneous algebraic reconstruction technique (SART) algorithm of Andersen and Kak, the so-called FA-SART algorithm.

Sharp Interface and Voltage Conservation in the Phase Field Method: Application to Cardiac Electrophysiology

Gregery T. Buzzard, Jeffrey J. Fox, and Fernando Siso-Nadal

SIAM J. Sci. Comput. 30, pp. 837-854 (18 pages) | Cited 1 time

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
We present a finite difference method for modeling the propagation of electrical waves in cardiac tissue using the cable equation with homogeneous Neumann boundary conditions. This method is novel in that it is derived by first discretizing the phase field method as described by Fenton et al. [Chaos, 15 (2005), p. 013502.] then taking a limit to recover a sharp interface. Our method provides for exact voltage conservation (up to floating point precision in application), can be used on arbitrary geometries, and can be implemented using only nearest neighbor coupling between grid points.

Deblurring Methods Using Antireflective Boundary Conditions

Martin Christiansen and Martin Hanke

SIAM J. Sci. Comput. 30, pp. 855-872 (18 pages) | Cited 2 times

Online Publication Date: February 14, 2008

Full Text: | Download PDF

Show Abstract
In this paper we consider the numerical solution of self-adjoint deblurring problems on bounded intervals. For these problems it has recently been shown that appropriate modeling of the solution near the boundary of the interval may significantly improve the numerical reconstructions. Among the alternatives the so-called antireflective boundary condition appears to be the best known choice. Here we develop an appropriate, i.e., stable and efficient implementation of this model in two steps, namely by (i) transforming the problem to homogeneous boundary values, and (ii) using the fast sine transform to solve the transformed problem. This approach allows us to incorporate regularization in a very straightforward way. Numerical reconstructions are superior qualitatively and quantitatively to those obtained with the reblurring method of Donatelli and Serra-Capizzano.

Numerical Coupling of Electric Circuit Equations and Energy-Transport Models for Semiconductors

Markus Brunk and Ansgar Jüngel

SIAM J. Sci. Comput. 30, pp. 873-894 (22 pages) | Cited 2 times

Online Publication Date: February 22, 2008

Full Text: | Download PDF

Show Abstract
A coupled semiconductor-circuit model including thermal effects is proposed. The charged particle flow in the semiconductor devices is described by the energy-transport equations for the electrons and the drift-diffusion equations for the holes. The electric circuit is modeled by the network equations from modified nodal analysis. The coupling is realized by the node potentials providing the voltages applied to the semiconductor devices and the output device currents for the network model. The resulting partial differential-algebraic system is discretized in time by the 2-stage backward difference formula and in space by a mixed-hybrid finite-element method using Marini–Pietra elements. A static condensation procedure is applied to the coupled system reducing the number of unknowns. Numerical simulations of a one-dimensional $pn$-junction diode with time-dependent voltage and of a rectifier circuit show the heating of the electrons which allows one to identify hot spots in the devices. Moreover, the choice of the boundary conditions for the electron density and energy density is numerically discussed.

Dynamic Phase Boundaries for Compressible Fluids

T. Lu, Z. L. Xu, R. Samulyak, J. Glimm, and X. M. Ji

SIAM J. Sci. Comput. 30, pp. 895-915 (21 pages)

Online Publication Date: February 22, 2008

Full Text: | Download PDF

Show Abstract
We present an algorithm for the simulation of a generalized Riemann problem for phase transitions in compressible fluids. We model the transition as a tracked jump discontinuity. The emphasis here is on the coupling of the phase transition process to acoustic waves, which is required for the study of cavitation induced by strong rarefaction waves. The robustness of the proposed algorithm is verified by application to various physical regimes.

Consistency on Domain Boundaries for Linear PDAE Systems

Jens Neumann and Constantinos C. Pantelides

SIAM J. Sci. Comput. 30, pp. 916-936 (21 pages)

Online Publication Date: March 05, 2008

Full Text: | Download PDF

Show Abstract
This paper considers the determination of consistent variations of the dependent variables over the boundary of the domain of definition of a linear system of partial differential and algebraic equations (PDAEs). Particular emphasis is placed on the specifications (“boundary conditions”) imposed on different parts of the boundary and their consistency with the underlying PDAE system and with each other. Specifications imposed on overlapping parts of the boundary (e.g., faces with common edges, or edges with common vertices) often lead to inconsistencies (corner singularities) that are not trivial to detect, especially in PDAEs involving three or more dimensions. A symbolic/numerical algorithm is proposed for the analysis of PDAE systems defined over finite hyperrectangular domains of arbitrary dimensions.

Numerical Methods for Computing Nonlinear Eigenpairs: Part II. Non-Iso-Homogeneous Cases

Xudong Yao and Jianxin Zhou

SIAM J. Sci. Comput. 30, pp. 937-956 (20 pages)

Online Publication Date: March 05, 2008

Full Text: | Download PDF

Show Abstract
Standing (solitary)-wave/steady-state solutions in many nonlinear wave motions and Schrödinger flows lead to nonlinear eigenproblems. In [X. Yao and J. Zhou, SIAM J. Sci. Comput., 29 (2007), pp. 1355–1374], a Rayleigh-local minimax method is developed to solve iso-homogeneous eigenproblems. In this subsequent paper, a unified method in Banach spaces is developed for solving non-iso-homogeneous and even nonhomogeneous eigenproblems and applied to solve two models: the Gross–Pitaevskii problem in the Bose–Einstein condensate and the $p$-Laplacian problem in non-Newtonian flows/materials. First a new active Lagrange functional is formulated to establish a local minimax characterization. A local minimax method is then devised and implemented to solve the model problems. Numerical results are presented. Convergence results of the algorithm and an order of eigensolutions computed by the algorithm are also established.

Advantages of Nonlinear-Programming-Based Methodologies for Inequality Path-Constrained Optimal Control Problems—A Numerical Study

Shivakumar Kameswaran and Lorenz T. Biegler

SIAM J. Sci. Comput. 30, pp. 957-981 (25 pages) | Cited 1 time

Online Publication Date: March 05, 2008

Full Text: | Download PDF

Show Abstract
We address some advantages of nonlinear programming (NLP)-based methods for inequality path-constrained optimal control problems. The analysis is carried out on the differential algebraic equation (DAE)-constrained optimization problem presented in [J. T. Betts and S. L. Campbell, Discretize Then Optimize, Technical report M&CT-TECH-03-01, The Boeing Company, Chicago, IL, 2003], which relates to a one-dimensional transient heat conduction problem. This optimal control problem possesses the essential characteristics of a typical inequality path-constrained optimal control problem. In a DAE-constrained optimization setting, the index of the path constraint can be arbitrarily increased, and this facilitates the study of the effect of the index of a path constraint on the solution obtained by using an NLP-based methodology. We show that the direct transcription approach leads to a problem whose limiting behavior does not satisfy the linear independence constraint qualification but satisfies the weaker Mangasarian–Fromowitz constraint qualification. The direct transcription approach leads to a convex quadratic programming problem with the result that the control obtained is the unique global minimizer. We also contrast this with the classical indirect approach where we show that it is difficult to numerically integrate the ODEs over the constrained arc because of stability and error control reasons. These observations explain some of the results of Betts and Campbell. This supports the fact that NLP-based methodologies have additional flexibility with respect to constraint qualifications, and this can be put to use in the case of inequality path-constrained optimal control problems to obtain well-defined solutions.

Symmetric Permutations for I-matrices to Delay and Avoid Small Pivots During Factorization

Jan Mayer

SIAM J. Sci. Comput. 30, pp. 982-996 (15 pages)

Online Publication Date: March 05, 2008

Full Text: | Download PDF

Show Abstract
In this article, we present several new permutations for I-matrices making these more suitable for incomplete LU-factorization preconditioners used in solving linear systems by iterative methods. A general matrix can be transformed by row permutation as well as row and columns scaling into an I-matrix, i.e., a matrix having elements of modulus 1 on the diagonal and elements of modulus of no more than 1 elsewhere. Reordering rows and columns by the same permutation clearly preserves I-matrices. In this article, we consider such reordering techniques which make the permuted matrix more suitable for an incomplete LU-factorization preconditioner than the original I-matrix. We use a multilevel ILUC, an incomplete LU-factorization preconditioner using Crout's implementation of Gaussian elimination without pivoting to test these reorderings. The combination of I-matrix preprocessing with the various algorithms presented here and the multilevel incomplete LU-factorizations forms a powerful preconditioning method for unsymmetric, highly indefinite problems. The C++ code has been made available in the software package ILU++.

S-ROCK: Chebyshev Methods for Stiff Stochastic Differential Equations

Assyr Abdulle and Stephane Cirilli

SIAM J. Sci. Comput. 30, pp. 997-1014 (18 pages) | Cited 2 times

Online Publication Date: March 05, 2008

Full Text: | Download PDF

Show Abstract
We present and analyze a new class of numerical methods for the solution of stiff stochastic differential equations (SDEs). These methods, called S-ROCK (for stochastic orthogonal Runge–Kutta Chebyshev), are explicit and of strong order 1 and possess large stability domains in the mean-square sense. For mean-square stable stiff SDEs, they are much more efficient than the standard explicit methods proposed so far for stochastic problems and give significant speed improvement. The explicitness of the S-ROCK methods allows one to handle large systems without linear algebra problems usually encountered with implicit methods. Numerical results and comparisons with existing methods are reported.

Adaptive, Fast, and Oblivious Convolution in Evolution Equations with Memory

María López-Fernández, Christian Lubich, and Achim Schädle

SIAM J. Sci. Comput. 30, pp. 1015-1037 (23 pages) | Cited 1 time

Online Publication Date: March 05, 2008

Full Text: | Download PDF

Show Abstract
To approximate convolutions which occur in evolution equations with memory terms, a variable-step-size algorithm is presented for which advancing $N$ steps requires only $O(N\log N)$ operations and $O(\log N)$ active memory, in place of $O(N^2)$ operations and $O(N)$ memory for a direct implementation. A basic feature of the fast algorithm is the reduction, via contour integral representations, to differential equations which are solved numerically with adaptive step sizes. Rather than the kernel itself, its Laplace transform is used in the algorithm. The algorithm is illustrated on three examples: a blowup example originating from a Schrödinger equation with concentrated nonlinearity, chemical reactions with inhibited diffusion, and viscoelasticity with a fractional order constitutive law.

Balanced Truncation Model Reduction for a Class of Descriptor Systems with Application to the Oseen Equations

Matthias Heinkenschloss, Danny C. Sorensen, and Kai Sun

SIAM J. Sci. Comput. 30, pp. 1038-1063 (26 pages) | Cited 1 time

Online Publication Date: March 05, 2008

Full Text: | Download PDF

Show Abstract
We discuss the computation of balanced truncation model reduction for a class of descriptor systems which include the semidiscrete Oseen equations with time-independent advection and the linearized Navier–Stokes equations, linearized around a steady state. The purpose of this paper is twofold. First, we show how to apply standard balanced truncation model reduction techniques, which apply to dynamical systems given by ordinary differential equations, to this class of descriptor systems. This is accomplished by eliminating the algebraic equation using a projection. The second objective of this paper is to demonstrate how the important class of ADI/Smith-type methods for the approximate computation of reduced order models using balanced truncation can be applied without explicitly computing the aforementioned projection. Instead, we utilize the solution of saddle point problems. We demonstrate the effectiveness of the technique in the computation of reduced order models for semidiscrete Oseen equations.

A Condition Number Analysis of a Line-Surface Intersection Algorithm

Gun Srijuntongsiri and Stephen A. Vavasis

SIAM J. Sci. Comput. 30, pp. 1064-1081 (18 pages) | Cited 1 time

Online Publication Date: March 07, 2008

Full Text: | Download PDF

Show Abstract
We propose an algorithm based on Newton's method and subdivision for finding all zeros of a polynomial system in a bounded region of the plane. This algorithm can be used to find the intersections between a line and a surface, which has applications in graphics and computer-aided geometric design. The algorithm can operate on polynomials represented in any basis that satisfies a few conditions. The power basis, the Bernstein basis, and the first-kind Chebyshev basis are among those compatible with the algorithm. The main novelty of our algorithm is an analysis showing that its running is bounded only in terms of the condition number of the polynomial's zeros and a constant depending on the polynomial basis.

Analysis of Aggregation-Based Multigrid

Adrian C. Muresan and Yvan Notay

SIAM J. Sci. Comput. 30, pp. 1082-1103 (22 pages) | Cited 4 times

Online Publication Date: March 07, 2008

Full Text: | Download PDF

Show Abstract
Aggregation-based multigrid with standard piecewise constant like prolongation is investigated. Unknowns are aggregated either by pairs or by quadruplets; in the latter case the grouping may be either linewise or boxwise. A Fourier analysis is developed for a model two-dimensional anisotropic problem. Most of the results are stated for an arbitrary smoother (which fits with the Fourier analysis framework). It turns out that the convergence factor of two-grid schemes can be bounded independently of the grid size. With a sensible choice of the (linewise or boxwise) coarsening, the bound is also uniform with respect to the anisotropy ratio, without requiring a specialized smoother. The bound is too large to guarantee optimal convergence properties with the V-cycle or the standard W-cycle, but a W-cycle scheme accelerated by the recursive use of the conjugate gradient method exhibits near grid independent convergence.
Close

close