Keyword
DOI
Citation
Advanced Search

Search Content Type

article doi:

Citation:

Home
Journals
E-books
For Authors
- Journal Author Submissions
- Book Author Submissions
Subscription Info
Interactive Features
Journal Citations
Contact Us
- Feedback
- SIAM Website
Scitation

SIAM J. on Scientific Computing

SECTION LISTING

SPECIAL SECTION

BROWSE VOLUMES

Year Range:

2012

Volume 34

2011

Volume 33

Issue 6 | 2011 | pp. 3087-3561

Issue 5 | 2011 | pp. 2115-3086 †

Issue 4 | 2011 | pp. 1475-2114

Issue 3 | 2011 | pp. 1035-1474

Issue 2 | 2011 | pp. 455-1033

Issue 1 | 2011 | pp. 1-454

2010

Volume 32

2008

Volume 31

2007

Volume 30

2007

Volume 29

2006

Volume 28

2005

Volume 27

2004

Volume 26

2003

Volume 25

2002

Volume 24

2001

Volume 23

2000

Volume 22

1999

Volume 21

1998

Volume 20

1998

Volume 19

1997

Volume 18

1996

Volume 17

1995

Volume 16

1994

Volume 15

1993

Volume 14

1992

Volume 13

1991

Volume 12

1990

Volume 11

1989

Volume 10

1988

Volume 9

1987

Volume 8

1986

Volume 7

1985

Volume 6

1984

Volume 5

1983

Volume 4

1982

Volume 3

1981

Volume 2

1980

Volume 1

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

Select this Article		3D Composite Finite Elements for Elliptic Boundary Value Problems with Discontinuous Coefficients Tobias Preusser, Martin Rumpf, Stefan Sauter, and Lars Ole Schwen SIAM J. Sci. Comput. 33, pp. 2115-2143 (29 pages) Online Publication Date: September 01, 2011 Full Text: \| Download PDF
	+ -	Show Abstract For scalar and vector-valued elliptic boundary value problems with discontinuous coefficients across geometrically complicated interfaces, a composite finite element approach is developed. Composite basis functions are constructed, mimicking the expected jump condition for the solution at the interface in an approximate sense. The construction is based on a suitable local interpolation on the space of admissible functions. We study the order of approximation and the convergence properties of the method numerically. As applications, heat diffusion in an aluminum foam matrix filled with polymer and linear elasticity of microstructured materials, in particular, specimens of trabecular bone, are investigated. Furthermore, a numerical homogenization approach is developed for periodic structures and real material specimens which are not strictly periodic but are considered as statistical prototypes. Thereby, effective macroscopic material properties can be computed.

Select this Article		Methods for Pricing American Options under Regime Switching Y. Huang, P. A. Forsyth, and G. Labahn SIAM J. Sci. Comput. 33, pp. 2144-2168 (25 pages) \| Cited 1 time Online Publication Date: September 01, 2011 Full Text: \| Download PDF
	+ -	Show Abstract We analyze a number of techniques for pricing American options under a regime switching stochastic process. The techniques analyzed include both explicit and implicit discretizations with the focus being on methods which are unconditionally stable. In the case of implicit methods we also compare a number of iterative procedures for solving the associated nonlinear algebraic equations. Numerical tests indicate that a fixed point policy iteration, coupled with a direct control formulation, is a reliable general purpose method. Finally, we remark that we formulate the American problem as an abstract optimal control problem; hence our results are applicable to more general problems as well.

Select this Article		An Online Method for Interpolating Linear Parametric Reduced-Order Models David Amsallem and Charbel Farhat SIAM J. Sci. Comput. 33, pp. 2169-2198 (30 pages) Online Publication Date: September 01, 2011 Full Text: \| Download PDF
	+ -	Show Abstract A two-step online method is proposed for interpolating projection-based linear parametric reduced-order models (ROMs) in order to construct a new ROM for a new set of parameter values. The first step of this method transforms each precomputed ROM into a consistent set of generalized coordinates. The second step interpolates the associated linear operators on their appropriate matrix manifold. Real-time performance is achieved by precomputing inner products between the reduced-order bases underlying the precomputed ROMs. The proposed method is illustrated by applications in mechanical and aeronautical engineering. In particular, its robustness is demonstrated by its ability to handle the case where the sampled parameter set values exhibit a mode veering phenomenon.

Select this Article		Duality Based A Posteriori Error Estimation for Quasi-Periodic Solutions Using Time Averages M. Braack, E. Burman, and N. Taschenberger SIAM J. Sci. Comput. 33, pp. 2199-2216 (18 pages) Online Publication Date: September 06, 2011 Full Text: \| Download PDF
	+ -	Show Abstract We propose an a posteriori error estimation technique for the computation of average functionals of solutions for nonlinear time dependent problems based on duality techniques. The exact solution is assumed to have a periodic or quasi-periodic behavior favoring a fixed mesh strategy in time. We show how to circumvent the need of solving time dependent dual problems. The estimator consists of an averaged residual weighted by sensitivity factors coming from a stationary dual problem and an additional averaging error term coming from nonlinearities of the operator considered. In order to illustrate this technique the resulting adaptive algorithm is applied to several model problems: a linear scalar parabolic problem with known exact solution, the nonsteady Navier–Stokes equations with known exact solution, and finally to the well-known benchmark problem for Navier–Stokes (flow behind a cylinder) in order to verify the modeling assumptions.

Select this Article		Time Discretization Schemes for Poincaré Waves in Finite-Element Shallow-Water Models Daniel Y. Le Roux, Michel Dieme, and Abdou Sene SIAM J. Sci. Comput. 33, pp. 2217-2246 (30 pages) Online Publication Date: September 13, 2011 Full Text: \| Download PDF
	+ -	Show Abstract The finite-element spatial discretization of the linear shallow-water equations is examined in the context of several temporal discretization schemes. Three finite-element pairs are considered, namely, the $P^{}_{0}-P^{}_{1}$, $P^{NC}_{1}-P^{}_{1}$, and $RT^{}_{0}-P^{}_{0}$ schemes, and the backward and forward Euler, Crank–Nicolson, and second and third order Adams–Bashforth time stepping schemes are employed. A Fourier analysis is performed at the discrete level for the Poincaré waves, and it determines the stability limit of the schemes and the error in wave amplitude and phase that can be expected. Numerical solutions of test problems to simulate Poincaré waves illustrate the analytical results.

SIAM J. on Scientific Computing

SECTION LISTING

BROWSE VOLUMES

2011

Volume 33, Issue 5, pp. 2115-3086

† Special Section: 2010 Copper Mountain Conference

Your Recent History

Recently Viewed

Find articles by AUTHORNAME

Select this Article		Nonadiabatic Transitions Through Tilted Avoided Crossings Volker Betz and Benjamin D. Goddard SIAM J. Sci. Comput. 33, pp. 2247-2276 (30 pages) Online Publication Date: September 13, 2011 Full Text: \| Download PDF
	+ -	Show Abstract We investigate the transition of a quantum wave-packet through a one-dimensional avoided crossing of molecular energy levels when the energy levels at the crossing point are tilted. Using superadiabatic representations, and an approximation of the dynamics near the crossing region, we obtain an explicit formula for the transition wave function. Our results agree extremely well with high precision ab-initio calculations.

Select this Article		A Higher Order Scheme for a Tangentially Stabilized Plane Curve Shortening Flow with a Driving Force Martin Balažovjech and Karol Mikula SIAM J. Sci. Comput. 33, pp. 2277-2294 (18 pages) Online Publication Date: September 15, 2011 Full Text: \| Download PDF
	+ -	Show Abstract We introduce a new higher order scheme for computing a tangentially stabilized curve shortening flow with a driving force represented by an intrinsic partial differential equation for an evolving curve position vector. Our new scheme is a combination of the explicit forward Euler and the fully implicit backward Euler schemes. At any discrete time step, the solution is found efficiently using a few semi-implicit iterations. Basic properties of the new scheme are proved in the paper, and its precision is tested by comparing the results with known analytical solutions. For any choice of the time step, the new higher order scheme gives exact radius of evolving uniformly discretized circles in case of flow by curvature and in case of rotation by a constant tangential velocity. Such properties do not hold for other schemes solving flow by mean curvature like the classical explicit, semi-implicit, or fully implicit schemes. In general, the scheme is second order accurate, which is shown by comparing a numerically evolving encompassed area with known analytical expression. The behavior of the scheme is discussed on representative examples, and its advantages with respect to the balance between CPU time and precision are shown.

Select this Article		Tensor-Train Decomposition I. V. Oseledets SIAM J. Sci. Comput. 33, pp. 2295-2317 (23 pages) \| Cited 3 times Online Publication Date: September 22, 2011 Full Text: \| Download PDF
	+ -	Show Abstract A simple nonrecursive form of the tensor decomposition in $d$ dimensions is presented. It does not inherently suffer from the curse of dimensionality, it has asymptotically the same number of parameters as the canonical decomposition, but it is stable and its computation is based on low-rank approximation of auxiliary unfolding matrices. The new form gives a clear and convenient way to implement all basic operations efficiently. A fast rounding procedure is presented, as well as basic linear algebra operations. Examples showing the benefits of the decomposition are given, and the efficiency is demonstrated by the computation of the smallest eigenvalue of a 19-dimensional operator.

Select this Article		A General Framework for Deriving Integral Preserving Numerical Methods for PDEs Morten Dahlby and Brynjulf Owren SIAM J. Sci. Comput. 33, pp. 2318-2340 (23 pages) Online Publication Date: September 22, 2011 Full Text: \| Download PDF
	+ -	Show Abstract A general procedure for constructing conservative numerical integrators for time-dependent partial differential equations is presented. In particular, linearly implicit methods preserving a time discretized version of the invariant are developed for systems of partial differential equations with polynomial nonlinearities. The framework is rather general and allows for an arbitrary number of dependent and independent variables with derivatives of any order. It is proved formally that second order convergence is obtained. The procedure is applied to a test case, and numerical experiments are provided.

Select this Article		Fast Evaluation of Multiquadric RBF Sums by a Cartesian Treecode Robert Krasny and Lei Wang SIAM J. Sci. Comput. 33, pp. 2341-2355 (15 pages) \| Cited 1 time Online Publication Date: September 22, 2011 Full Text: \| Download PDF
	+ -	Show Abstract A treecode is presented for evaluating sums defined in terms of the multiquadric radial basis function (RBF), $\phi({\bf x}) = (\|{\bf x}\|^2+c^2)^{1/2}$, where ${\bf x} \in \mathbb{R}^3$ and $c \ge 0$. Given a set of $N$ nodes, evaluating an RBF sum directly requires CPU time that scales like $O(N^2)$. For a given level of accuracy, the treecode reduces the CPU time to $O(N\log N)$ using a far-field expansion of $\phi({\bf x})$. We consider two options for the far-field expansion: (1) a Laurent series previously used in applications of the Fast Multipole Method to multiquadric RBFs, and (2) a certain Taylor series previously used in treecode particle simulations, but not yet in the context of multiquadric RBFs. It is known that the Laurent series converges when the RBF parameter $c$ lies in an interval $0 \le c \le \bar{c}$, where $\bar{c}$ is proportional to the minimum node spacing, but here we show that the Taylor series converges uniformly for $c \ge 0$. We implement the treecode in Cartesian coordinates and use a recurrence relation to compute the Taylor coefficients. Numerical results exhibit the treecode error, CPU time, and memory usage in two test cases, random nodes in a cube and on the surface of a sphere. The treecode approach presented here is applicable to generalized multiquadrics in any dimension.

Select this Article		A Well-Balanced Symplecticity-Preserving Gas-Kinetic Scheme for Hydrodynamic Equations under Gravitational Field Jun Luo, Kun Xu, and Na Liu SIAM J. Sci. Comput. 33, pp. 2356-2381 (26 pages) Online Publication Date: September 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract A well-balanced scheme for an isolated gravitational hydrodynamic system is defined as a scheme which exactly preserves an isothermal hydrostatic solution. In this paper, a well-balanced gas-kinetic symplecticity-preserving BGK (SP-BGK) scheme is developed. In the construction of the scheme, the gravitational potential is modeled as a piecewise constant function inside each cell with a potential jump at the cell interface. In the process of designing such a scheme, the energy conservation, Liouville's theorem, and the symplecticity-preserving property of a Hamiltonian flow play important roles in the description of particles penetration and reflection from a potential barrier. More importantly, the use of the symplecticity-preserving property is crucial in the evaluation of the moments of a postinteraction gas distribution function with a potential jump in terms of the moments of preinteraction distribution function. The SP-BGK method is the first well-balanced shock-capturing gas-kinetic scheme for the Navier–Stokes equation. A few theorems are proved for this scheme, which include the necessity to use an exact Maxwellian for keeping the isothermal hydrostatic state, the total mass and energy (the sum of kinetic, thermal, and gravitational ones) conservation, and the well-balanced property of the SP-BGK scheme to keep an isothermal hydrostatic state during the process of particle transport and collision. Many numerical examples are presented to validate the SP-BGK scheme.

Select this Article		Algorithms for Area Preserving Flows Catherine Kublik, Selim Esedoḡlu, and Jeffrey A. Fessler SIAM J. Sci. Comput. 33, pp. 2382-2401 (20 pages) Online Publication Date: September 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract We propose efficient and accurate algorithms for computing certain area preserving geometric motions of curves in the plane, such as area preserving motion by curvature. These schemes are based on a new class of diffusion generated motion algorithms using signed distance functions. In particular, they alternate two very simple and fast operations, namely convolution with the Gaussian kernel and construction of the distance function, to generate the desired geometric flow in an unconditionally stable manner. We present applications of these area preserving flows to large scale simulations of coarsening.

Select this Article		A Third Order Accurate Fast Marching Method for the Eikonal Equation in Two Dimensions Shahnawaz Ahmed, Stanley Bak, Joyce McLaughlin, and Daniel Renzi SIAM J. Sci. Comput. 33, pp. 2402-2420 (19 pages) Online Publication Date: September 29, 2011 Full Text: \| Download PDF
	+ -	Show Abstract In this paper, we develop a third order accurate fast marching method for the solution of the eikonal equation in two dimensions. There have been two obstacles to extending the fast marching method to higher orders of accuracy. The first obstacle is that using one-sided difference schemes is unstable for orders of accuracy higher than two. The second obstacle is that the points in the difference stencil are not available when the gradient is closely aligned with the grid. We overcome these obstacles by using a two-dimensional (2D) finite difference approximation to improve stability, and by locally rotating the grid 45 degrees (i.e., using derivatives along the diagonals) to ensure all the points needed in the difference stencil are available. We show that in smooth regions the full difference stencil is used for a suitably small enough grid size and that the difference scheme satisfies the von Neumann stability condition for the linearized eikonal equation. Our method reverts to first order accuracy near caustics without developing oscillations by using a simple switching scheme. The efficiency and high order of the method are demonstrated on a number of 2D test problems.

Select this Article		Stable and Spectrally Accurate Schemes for the Navier–Stokes Equations Jun Jia and Jie Liu SIAM J. Sci. Comput. 33, pp. 2421-2439 (19 pages) Online Publication Date: September 29, 2011 Full Text: \| Download PDF
	+ -	Show Abstract In this paper, we present an accurate, efficient and stable numerical method for the incompressible Navier–Stokes equations (NSEs). The method is based on (1) an equivalent pressure Poisson equation formulation of the NSE with proper pressure boundary conditions, which facilitates the design of high-order and stable numerical methods, and (2) the Krylov deferred correction (KDC) accelerated method of lines transpose ($\mbox{MoL}^{\it T}$), which is very stable, efficient, and of arbitrary order in time. Numerical tests with known exact solutions in three dimensions show that the new method is spectrally accurate in time, and a numerical order of convergence 9 was observed. Two-dimensional computational results of flow past a cylinder and flow in a bifurcated tube are also reported.

Select this Article		Schur Complement Preconditioners for Surface Integral-Equation Formulations of Dielectric Problems Solved with the Multilevel Fast Multipole Algorithm Tahi̇r Malas and Levent Gürel SIAM J. Sci. Comput. 33, pp. 2440-2467 (28 pages) Online Publication Date: October 04, 2011 Full Text: \| Download PDF
	+ -	Show Abstract Surface integral-equation methods accelerated with the multilevel fast multipole algorithm (MLFMA) provide a suitable mechanism for electromagnetic analysis of real-life dielectric problems. Unlike the perfect-electric-conductor case, discretizations of surface formulations of dielectric problems yield $2 \times 2$ partitioned linear systems. Among various surface formulations, the combined tangential formulation (CTF) is the closest to the category of first-kind integral equations, and hence it yields the most accurate results, particularly when the dielectric constant is high and/or the dielectric problem involves sharp edges and corners. However, matrix equations of CTF are highly ill-conditioned, and their iterative solutions require powerful preconditioners for convergence. Second-kind surface integral-equation formulations yield better conditioned systems, but their conditionings significantly degrade when real-life problems include high dielectric constants. In this paper, for the first time in the context of surface integral-equation methods of dielectric objects, we propose Schur complement preconditioners to increase their robustness and efficiency. First, we approximate the dense system matrix by a sparse near-field matrix, which is formed naturally by MLFMA. The Schur complement preconditioning requires approximate solutions of systems involving the (1,1) partition and the Schur complement. We approximate the inverse of the (1,1) partition with a sparse approximate inverse (SAI) based on the Frobenius norm minimization. For the Schur complement, we first approximate it via incomplete sparse matrix-matrix multiplications, and then we generate its approximate inverse with the same SAI technique. Numerical experiments on sphere, lens, and photonic crystal problems demonstrate the effectiveness of the proposed preconditioners. In particular, the results for the photonic crystal problem, which has both surface singularity and a high dielectric constant, shows that accurate CTF solutions for such problems can be obtained even faster than with second-kind integral equation formulations, with the acceleration provided by the proposed Schur complement preconditioners.

Select this Article		A Sinc Function Analogue of Chebfun Mark Richardson and Lloyd N. Trefethen SIAM J. Sci. Comput. 33, pp. 2519-2535 (17 pages) Online Publication Date: October 13, 2011 Full Text: \| Download PDF
	+ -	Show Abstract Chebfun is an established software system for computing with functions of a real variable, but its capabilities for handling functions with singularities are limited. Here an analogous system is described based on sinc function expansions instead of Chebyshev series. This experiment sheds light on the strengths and weaknesses of sinc function techniques. It also serves as a review of some of the main features of sinc methods, including construction, evaluation, zerofinding, optimization, integration, and differentiation.

Select this Article		Interpolating Irregularly Spaced Observations for Filtering Turbulent Complex Systems John Harlim SIAM J. Sci. Comput. 33, pp. 2620-2640 (21 pages) Online Publication Date: October 20, 2011 Full Text: \| Download PDF
	+ -	Show Abstract We present a numerically fast reduced filtering strategy, the Fourier domain Kalman filter with appropriate interpolations to account for irregularly spaced observations of complex turbulent signals. The design of such a reduced filter involves: (i) interpolating irregularly spaced observations to the model regularly spaced grid points, (ii) understanding under which situation the small scale oscillatory artifact from such interpolation won't degrade the filtered solutions, (iii) understanding when the interpolated covariance structure can be approximated by its diagonal terms when observations are corrupted by independent Gaussian noise, and (iv) applying a scalar Kalman filter formula on each Fourier component independently with an approximate diagonal interpolated covariance matrix. From the practical point of view, there is an emerging need to understand the effect of (i) toward the filtered solutions, for example, in utilizing the data produced from various interpolation techniques that merge multiple satellite measurements of atmospheric and ocean dynamical quantities. To understand point (iii) above and to see how many of nondiagonal terms are effectively ignored in (iv), we compute a ratio $\Lambda$ between the largest nondiagonal components and the smallest diagonal components of the interpolated covariance matrix. We find that for piecewise linear interpolation, this ratio, $\Lambda$, is always smaller than that of the alternative interpolation schemes such as trigonometric and cubic spline for any irregularly spaced observation networks. When observations are not so sparse with small noise, we find that the small scale oscillatory artifact in (ii) above is negligible when piecewise linear interpolation is used whereas for the other schemes such as the nearest neighbor, trigonometric, and cubic spline interpolation, the oscillatory artifact degrades the filtered solutions significantly. Finally, we also find that the reduced filtering strategy with piecewise linear interpolation produces more accurate filtered solutions than conventional approaches when observations are extremely irregularly spaced (such that the ratio $\Lambda$ is not so small) and very sparse.

Select this Article FREE		Special Section: 2010 Copper Mountain Conference Ray Tuminaro, Guest Editor , Michele Benzi, Xiao-Chuan Cai, Iain Duff, Howard Elman, Roland Freund, Kirk Jordan, Tim Kelley, David Keyes, Misha Kilmer, Sven Leyffer, Tom Manteuffel, Steve McCormick, David Silvester, Homer Walker, et al. SIAM J. Sci. Comput. 33, pp. 2685-2685 (1 page) Full Text: \| Download PDF
	+ -	Show Abstract The biennial Copper Mountain Conference on Iterative Methods was held April 4–9, 2010. This meeting included more than 140 presentations covering many scientific computing areas, such as uncertainty quantification, optimization, Markov chains, saddle-point systems, inverse problems, direct factorizations, Krylov methods, algebraic multigrid, software frameworks, and advanced computer architectures. A partial list of the represented application domains included acoustic scattering, material science, circuit simulation, neutron transport, image processing, magnetohydrodynamics, and porous media flow. Howard Elman (University of Maryland) and Ray Tuminaro (Sandia National Laboratories) served as co-chairs. Management services were provided by Annette Anthony of Front Range Scientific Computations, Inc. The Conference is gratefully for financial support from the Department of Energy, IBM, Lawrence Livermore National Laboratory, Los Alamos National Laboratory, Sandia National Laboratories, and the National Science Foundation. The Copper Mountain Conference is organized in cooperation with the Society for Industrial and Applied Mathematics (SIAM). Submissions to this special section were open to the public. This was advertised in advance on web sites of the SIAM Journal on Scientific Computing (SISC) and the Copper Mountain Conference. The collection of papers in this section all center on iterative methods. However, many diverse topics are covered spanning contributions with a theoretical orientation to those that are more applied or driven by a specific application area. Additionally, several papers are highly influenced by current computer architecture trends. Some of the submissions also appeared in the student paper competition, a unique event at the Copper Mountain conferences that has proven to be a successful forum for presenting the latest results in the field. We have had increasingly strong participation from talented recent Ph.D. holders and current Ph.D. students in this competition. I wish to thank the guest editorial board, the SIAM/SISC staff, and the referees for their hard work to help meet deadlines. Special acknowledgment is due to Mitch Chernoff (SIAM Publications Manager) and Heather Blythe (SIAM Senior Publications Coordinator) for their efforts on this special section. We hope that this special section will continue to foster new developments in this important area of research. The diversity of topics clearly demonstrates the significant impact of iterative methods across many scientific computing areas and the need for continued research.

Select this Article		Multigrid Smoothers for Ultraparallel Computing Allison H. Baker, Robert D. Falgout, Tzanio V. Kolev, and Ulrike Meier Yang SIAM J. Sci. Comput. 33, pp. 2864-2887 (24 pages) \| Cited 1 time Online Publication Date: October 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract This paper investigates the properties of smoothers in the context of algebraic multigrid (AMG) running on parallel computers with potentially millions of processors. The development of multigrid smoothers in this case is challenging, because some of the best relaxation schemes, such as the Gauss–Seidel (GS) algorithm, are inherently sequential. Based on the sharp two-grid multigrid theory from [R. D. Falgout and P. S. Vassilevski, SIAM J. Numer. Anal., 42 (2004), pp. 1669–1693] and [R. D. Falgout, P. S. Vassilevski, and L. T. Zikatanov, Numer. Linear Algebra Appl., 12 (2005), pp. 471–494] we characterize the smoothing properties of a number of practical candidates for parallel smoothers, including several $C$-$F$, polynomial, and hybrid schemes. We show, in particular, that the popular hybrid GS algorithm has multigrid smoothing properties which are independent of the number of processors in many practical applications, provided that the problem size per processor is large enough. This is encouraging news for the scalability of AMG on ultraparallel computers. We also introduce the more robust $\ell_1$ smoothers, which are always convergent and have already proven essential for the parallel solution of some electromagnetic problems [T. Kolev and P. Vassilevski, J. Comput. Math., 27 (2009), pp. 604–623].

Select this Article		A Fast Iterative Method for Solving the Eikonal Equation on Triangulated Surfaces Zhisong Fu, Won-Ki Jeong, Yongsheng Pan, Robert M. Kirby, and Ross T. Whitaker SIAM J. Sci. Comput. 33, pp. 2468-2488 (21 pages) Online Publication Date: October 06, 2011 Full Text: \| Download PDF
	+ -	Show Abstract This paper presents an efficient, fine-grained parallel algorithm for solving the Eikonal equation on triangular meshes. The Eikonal equation, and the broader class of Hamilton–Jacobi equations to which it belongs, have a wide range of applications from geometric optics and seismology to biological modeling and analysis of geometry and images. The ability to solve such equations accurately and efficiently provides new capabilities for exploring and visualizing parameter spaces and for solving inverse problems that rely on such equations in the forward model. Efficient solvers on state-of-the-art, parallel architectures require new algorithms that are not, in many cases, optimal, but are better suited to synchronous updates of the solution. In previous work [W. K. Jeong and R. T. Whitaker, SIAM J. Sci. Comput., 30 (2008), pp. 2512–2534], the authors proposed the fast iterative method (FIM) to efficiently solve the Eikonal equation on regular grids. In this paper we extend the fast iterative method to solve Eikonal equations efficiently on triangulated domains on the CPU and on parallel architectures, including graphics processors. We propose a new local update scheme that provides solutions of first-order accuracy for both architectures. We also propose a novel triangle-based update scheme and its corresponding data structure for efficient irregular data mapping to parallel single-instruction multiple-data (SIMD) processors. We provide detailed descriptions of the implementations on a single CPU, a multicore CPU with shared memory, and SIMD architectures with comparative results against state-of-the-art Eikonal solvers.

Select this Article		Interpolatory Projection Methods for Parameterized Model Reduction Ulrike Baur, Christopher Beattie, Peter Benner, and Serkan Gugercin SIAM J. Sci. Comput. 33, pp. 2489-2518 (30 pages) Online Publication Date: October 13, 2011 Full Text: \| Download PDF
	+ -	Show Abstract We provide a unifying projection-based framework for structure-preserving interpolatory model reduction of parameterized linear dynamical systems, i.e., systems having a structured dependence on parameters that we wish to retain in the reduced-order model. The parameter dependence may be linear or nonlinear and is retained in the reduced-order model. Moreover, we are able to give conditions under which the gradient and Hessian of the system response with respect to the system parameters is matched in the reduced-order model. We provide a systematic approach built on established interpolatory $\mathcal{H}_2$ optimal model reduction methods that will produce parameterized reduced-order models having high fidelity throughout a parameter range of interest. For single input/single output systems with parameters in the input/output maps, we provide reduced-order models that are optimal with respect to an $\mathcal{H}_2\otimes\mathcal{L}_2$ joint error measure. The capabilities of these approaches are illustrated by several numerical examples from technical applications.

Select this Article		Geometric Properties of the Icosahedral-Hexagonal Grid on the Two-Sphere Ning Wang and Jin-Luen Lee SIAM J. Sci. Comput. 33, pp. 2536-2559 (24 pages) Online Publication Date: October 18, 2011 Full Text: \| Download PDF
	+ -	Show Abstract An icosahedral-hexagonal grid on the two-sphere is created by dividing the faces of an icosahedron and projecting the vertices onto the sphere. This grid and its Voronoi tessellation have several desirable features for numerical simulations of physical processes on the sphere. While several methods to construct the icosahedral grid mesh have been proposed over the past decades, and empirical data have been collected to understand and help improve the grid, rarely have analytical analyses been done to investigate the basic geometric properties of the grid. In this paper, we present an analytical analysis of several geometric properties of the icosahedral grids based on two basic constructions: recursive and nonrecursive construction. We point out that these geometric properties can be improved with modified construction procedures. We demonstrate how these improvements impact the numerical integration of PDEs over the sphere.

Select this Article		A Comparison of the Dispersion and Dissipation Errors of Gauss and Gauss–Lobatto Discontinuous Galerkin Spectral Element Methods Gregor Gassner and David A. Kopriva SIAM J. Sci. Comput. 33, pp. 2560-2579 (20 pages) Online Publication Date: October 18, 2011 Full Text: \| Download PDF
	+ -	Show Abstract We examine the dispersion and dissipation properties of the Gauss and Gauss–Lobatto discontinuous Galerkin spectral element methods (DGSEMs), for linear wave propagation problems. We show that the inherent underintegration in the Gauss–Lobatto variant can be interpreted as a modal filtering of the highest polynomial mode. This in turn has a drastic impact on the dispersion and dissipation relations of the Gauss–Lobatto DGSEM compared to the Gauss variant. We show that the Gauss DGSEM is typically more accurate than the Gauss–Lobatto variant, needing fewer points per wavelength for a given accuracy while on the other hand being more restricted in the explicit time step choice. We show that the spectra of the DGSEM operators depend on the boundary conditions applied and that the ratio of the time step restrictions of the two schemes depends on the choice of boundary conditions.

Select this Article		An Algorithm for the Principal Component Analysis of Large Data Sets Nathan Halko, Per-Gunnar Martinsson, Yoel Shkolnisky, and Mark Tygert SIAM J. Sci. Comput. 33, pp. 2580-2594 (15 pages) \| Cited 1 time Online Publication Date: October 18, 2011 Full Text: \| Download PDF
	+ -	Show Abstract Recently popularized randomized methods for principal component analysis (PCA) efficiently and reliably produce nearly optimal accuracy—even on parallel processors—unlike the classical (deterministic) alternatives. We adapt one of these randomized methods for use with data sets that are too large to be stored in random-access memory (RAM). (The traditional terminology is that our procedure works efficiently out-of-core.) We illustrate the performance of the algorithm via several numerical examples. For example, we report on the PCA of a data set stored on disk that is so large that less than a hundredth of it can fit in our computer's RAM.

Select this Article		Nonsymmetric Preconditioner Updates in Newton–Krylov Methods for Nonlinear Systems Stefania Bellavia, Daniele Bertaccini, and Benedetta Morini SIAM J. Sci. Comput. 33, pp. 2595-2619 (25 pages) Online Publication Date: October 20, 2011 Full Text: \| Download PDF
	+ -	Show Abstract Newton–Krylov methods, a combination of Newton-like methods and Krylov subspace methods for solving the Newton equations, often need adequate preconditioning in order to be successful. Approximations of the Jacobian matrices are required to form preconditioners, and this step is very often the dominant cost of Newton–Krylov methods. Therefore, working with preconditioners may destroy the “Jacobian-free” (or matrix-free) setting where the single Jacobian-vector product can be provided without forming and storing the element of the true Jacobian. In this paper, we propose and analyze a preconditioning technique for sequences of nonsymmetric Jacobian matrices based on the update of an earlier preconditioner. The proposed strategy can be implemented in a matrix-free manner. Numerical experiments on popular test problems confirm the effectiveness of the approach in comparison with the standard ILU-preconditioned Newton–Krylov approaches.

Select this Article		Optimal Control of Stochastic Flow over a Backward-Facing Step Using Reduced-Order Modeling Max Gunzburger and Ju Ming SIAM J. Sci. Comput. 33, pp. 2641-2663 (23 pages) Online Publication Date: October 20, 2011 Full Text: \| Download PDF
	+ -	Show Abstract We explore a stochastic optimal control problem for incompressible flow in a channel by using reduced-order modeling (ROM) based on proper orthogonal decomposition (POD). The control acts on a portion of the boundary and has a stochastic perturbation. The objective of control is to minimize the enstrophy, i.e., the $L^2$ norm of the vorticity. Owing to the high computational requirements, the classical Monte Carlo method is not an appropriate option. The POD basis is constructed in the usual manner from a set of “snapshots” obtained by sampling the solution of the Navier–Stokes system at several time instants and over several independent realizations of the random data. Then, the stochastic ROM is constructed via a Galerkin method with respect to the POD basis. The optimal control problem is solved using the stochastic ROM and compared to results obtained using a full finite element discretization that employs many more degrees of freedom. These comparisons show the effectiveness of the use of the stochastic ROM for obtaining accurate optimal solutions of the stochastic control problem at low cost.

Select this Article		Superconvergence of Discontinuous Galerkin Solutions for Delay Differential Equations of Pantograph Type Qiumei Huang, Hehu Xie, and Hermann Brunner SIAM J. Sci. Comput. 33, pp. 2664-2684 (21 pages) Online Publication Date: October 25, 2011 Full Text: \| Download PDF
	+ -	Show Abstract This paper is concerned with the superconvergence of the discontinuous Galerkin solutions for delay differential equations with proportional delays vanishing at $t = 0$. Two types of superconvergence are analyzed here. The first is based on interpolation postprocessing to improve the global convergence order by finding the superconvergence points of discontinuous Galerkin solutions. The second type follows from the integral iteration which just requires a local integration procedure applied to the discontinuous Galerkin solution, thus increasing the order of convergence. The theoretical results are illustrated by a broad range of numerical examples.

Select this Article		Estimates of the Norm of the Error in Solving Linear Systems with FOM and GMRES Gérard Meurant SIAM J. Sci. Comput. 33, pp. 2686-2705 (20 pages) Online Publication Date: October 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract We provide formulas for the norm of the error when solving nonsymmetric linear systems with the full orthogonalization method (FOM) and the generalized minimum residual method (GMRES) as well as relations between the error norm and the residual norm. From these formulas we are able to compute estimates of the norm of the error during the iterations. Since stopping criteria based on the norm of the residual may sometimes be misleading, such estimates could lead to a more robust way to stop the iterations. Numerical experiments show that the proposed norm estimates work nicely on difficult linear systems.

Select this Article		Algebraic Multigrid for Linear Systems Obtained by Explicit Element Reduction Thomas A. Brunner and Tzanio V. Kolev SIAM J. Sci. Comput. 33, pp. 2706-2731 (26 pages) Online Publication Date: October 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract We consider sparse linear systems, where a set of “interior” degrees of freedom have been eliminated in order to reduce the problem size. This elimination is assumed to be local, so the “interior” principal submatrix is block-diagonal, and the resulting Schur complement is still sparse. For it to be beneficial, the elimination process should lead to reduced memory requirements, but equally importantly, it should also result in an algebraic problem that can be solved efficiently. In this paper we propose a general element reduction approach and show how the elimination process can exploit a particular “subzonal” discretization to maintain the sparsity of the Schur complement. We also investigate algebraic multigrid (AMG) solution algorithms applied to the reduced problem, and we discuss the influence of the local elimination on solver-related properties of the matrix, such as near-nullspace preservation and the availability of stable subspace decompositions. We focus on BoomerAMG, a parallel variant of classical Ruge–Stüben AMG, applied to scalar diffusion problems [V. Henson and U. Yang, Appl. Numer. Math., 41 (2002), pp. 155–177], and the auxiliary-space Maxwell solver (AMS) for electromagnetic diffusion applications [T. Kolev and P. Vassilevski, J. Comput. Math., 27 (2009), pp. 604–623]. In the electromagnetic case, we establish algebraically a reduced version of the Hiptmair–Xu decomposition from [R. Hiptmair and J. Xu, SIAM J. Numer. Anal., 45 (2007), pp. 2483–2509] and consider a modification of the reduction process that targets the singular problems arising in simulations with pure void (zero conductivity) regions. For scalar diffusion problems, our particular stencil analysis shows that the reduction has a positive effect on meshes with stretched elements. We present a number of two-dimensional, three-dimensional, and axisymmetric numerical experiments, which demonstrate that the combination of an appropriately chosen local elimination with the use of the BoomerAMG and AMS solvers can lead to significant improvements in the overall solution time.

Select this Article		Peano—A Traversal and Storage Scheme for Octree-Like Adaptive Cartesian Multiscale Grids Tobias Weinzierl and Miriam Mehl SIAM J. Sci. Comput. 33, pp. 2732-2760 (29 pages) Online Publication Date: October 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract Almost all approaches to solving partial differential equations (PDEs) are based upon a spatial discretization of the computational domain—a grid. This paper presents an algorithm to generate, store, and traverse a hierarchy of $d$-dimensional Cartesian grids represented by a $(k=3)$-spacetree, a generalization of the well-known octree concept, and it also shows the correctness of the approach. These grids may change their adaptive structure throughout the traversal. The algorithm uses $2d+4$ stacks as data structures for both cells and vertices, and the storage requirements for the pure grid reduce to one bit per vertex for both the complete grid connectivity structure and the multilevel grid relations. Since the traversal algorithm uses only stacks, the algorithm's cache hit rate is continually higher than 99.9 percent, and the runtime per vertex remains almost constant; i.e., it does not depend on the overall number of vertices or the adaptivity pattern. We use the algorithmic approach as the fundamental concept for a mesh management for $d$-dimensional PDEs and for a matrix-free PDE solver represented by a compact discrete $3^d$-point operator. In the latter case, one can implement a Jacobi smoother, a Krylov solver, or a geometric multigrid scheme within the presented traversal scheme which inherits the low memory requirements and the good memory access characteristics directly.

Select this Article		Analysis of Augmented Lagrangian-Based Preconditioners for the Steady Incompressible Navier–Stokes Equations Michele Benzi and Zhen Wang SIAM J. Sci. Comput. 33, pp. 2761-2784 (24 pages) Online Publication Date: October 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract We analyze a class of modified augmented Lagrangian-based preconditioners for both stable and stabilized finite element discretizations of the steady incompressible Navier–Stokes equations. We study the eigenvalues of the preconditioned matrices obtained from Picard linearization, and we devise a simple and effective method for the choice of the augmentation parameter $\gamma$ based on Fourier analysis. Numerical experiments on a wide range of model problems demonstrate the robustness of these preconditioners, yielding fast convergence independent of mesh size and only mildly dependent on viscosity on both uniform and stretched grids. Good results are also obtained on linear systems arising from Newton linearization. We also show that performing inexact preconditioner solves with an algebraic multigrid algorithm results in excellent scalability. Comparisons of the modified augmented Lagrangian preconditioners with other state-of-the-art techniques show the competitiveness of our approach.

Select this Article		Iterative Methods for Neutron Transport Eigenvalue Problems Fynn Scheben and Ivan G. Graham SIAM J. Sci. Comput. 33, pp. 2785-2804 (20 pages) Online Publication Date: October 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract We discuss iterative methods for computing criticality in nuclear reactors. In general this requires the solution of a generalized eigenvalue problem for an unsymmetric integro-differential operator in six independent variables, modeling transport, scattering, and fission, where the dependent variable is the neutron angular flux. In engineering practice this problem is often solved iteratively, using some variant of the inverse power method. Because of the high dimension, matrix representations for the operators are often not available and the inner solves needed for the eigenvalue iteration are implemented by matrix-free inner iterations. This leads to technically complicated inexact iterative methods, for which there appears to be no published rigorous convergence theory. For the monoenergetic homogeneous model case with isotropic scattering and vacuum boundary conditions, we show that, before discretization, the general nonsymmetric eigenproblem for the angular flux is equivalent to a certain related eigenproblem for the scalar flux, involving a symmetric positive definite weakly singular integral operator (in space only). This correspondence to a symmetric problem (in a space of reduced dimension) permits us to give a convergence theory for inexact inverse iteration and related methods. In particular this theory provides rather precise criteria on how accurate the inner solves need to be in order for the whole iterative method to converge. We also give examples of discretizations which have a corresponding symmetric finite-dimensional reduced form. The theory is illustrated with numerical examples for several test problems of physical relevance, using GMRES as the inner solver.

Select this Article		A Second Order Discretization of Maxwell's Equations in the Quasi-Static Regime on OcTree Grids Lior Horesh and Eldad Haber SIAM J. Sci. Comput. 33, pp. 2805-2822 (18 pages) Online Publication Date: October 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract In this study we consider adaptive mesh refinement for the solution of Maxwell's equations in the quasi-static or diffusion regime. We propose a new finite volume OcTree discretization for the problem and show how to construct second order stencils on Yee grids, extending the known first order discretization stencils. We then develop an effective preconditioner to the problem. We show that our preconditioner performs well for discontinuous conductivities as well as for a wide range of frequencies.

Select this Article		Domain-Decomposition-Type Methods for Computing the Diagonal of a Matrix Inverse Jok M. Tang and Yousef Saad SIAM J. Sci. Comput. 33, pp. 2823-2847 (25 pages) Online Publication Date: October 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract This paper presents two methods based on domain decomposition concepts for determining the diagonal of the inverse of specific matrices. The first uses a divide-and-conquer principle and the Sherman–Morrison–Woodbury formula and assumes that the matrix can be decomposed into a $2 \times 2$ block-diagonal matrix and a low-rank matrix. The second method is a standard domain decomposition approach in which local solves are combined with a global correction. Both methods can be successfully combined with iterative solvers and sparse approximation techniques. The efficiency of the methods usually depends on the specific implementation, which should be fine-tuned for different test problems. Preliminary results for some two-dimensional (2D) problems are reported to illustrate the proposed methods.

Select this Article		Efficient Parallel Nonnegative Least Squares on Multicore Architectures Yuancheng Luo and Ramani Duraiswami SIAM J. Sci. Comput. 33, pp. 2848-2863 (16 pages) Online Publication Date: October 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract We parallelize a version of the active-set iterative algorithm derived from the original works of Lawson and Hanson [Solving Least Squares Problems, Prentice-Hall, 1974] on multicore architectures. This algorithm requires the solution of an unconstrained least squares problem in every step of the iteration for a matrix composed of the passive columns of the original system matrix. To achieve improved performance, we use parallelizable procedures to efficiently update and downdate the $QR$ factorization of the matrix at each iteration, to account for inserted and removed columns. We use a reordering strategy of the columns in the decomposition to reduce computation and memory access costs. We consider graphics processing units (GPUs) as a new mode for efficient parallel computations and compare our implementations to that of multicore CPUs. Both synthetic and nonsynthetic data are used in the experiments.

Select this Article		Algebraic Multigrid for High-Order Hierarchical $H(curl)$ Finite Elements James H. Lai and Luke N. Olson SIAM J. Sci. Comput. 33, pp. 2888-2902 (15 pages) Online Publication Date: October 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract Classic multigrid methods are often not directly applicable to nonelliptic problems such as curl-type partial differential equations (PDEs). Curl-curl PDEs require specialized smoothers that are compatible with the gradient-like (near) null space. Moreover, recent developments have focused on replicating the grad-curl-div de Rham complex in a multilevel hierarchy through smoothed aggregation based algebraic multigrid. These approaches have been successful for Nédélec finite elements (i.e., $H(curl)$ edge elements), but do not extend naturally to high-order representations. In this paper we consider hierarchical high-order Whitney elements for the curl-curl eddy current problem and devise a scalable multilevel approach. Our method generates a hierarchy similar to $p$-type multigrid, which results in a coarse level that is amenable to further coarsening through the established process of a multilevel complex. The natural hierarchy of the elements induces an effective interpolation operator and motivates the construction of a compatible gradient smoothing process. We detail the multilevel solver for a hierarchical $H(curl)$ basis and present numerical results in support of the method.

Select this Article		Preconditioning Iterative Methods for the Optimal Control of the Stokes Equations Tyrone Rees and Andrew J. Wathen SIAM J. Sci. Comput. 33, pp. 2903-2926 (24 pages) Online Publication Date: October 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract Solving problems regarding the optimal control of partial differential equations (PDEs)—also known as PDE-constrained optimization—is a frontier area of numerical analysis. Of particular interest is the problem of flow control, where one would like to effect some desired flow by exerting, for example, an external force. The bottleneck in many current algorithms is the solution of the optimality system—a system of equations in saddle point form that is usually very large and ill conditioned. In this paper we describe two preconditioners—a block diagonal preconditioner for the minimal residual method and a block lower-triangular preconditioner for a nonstandard conjugate gradient method—which can be effective when applied to such problems where the PDEs are the Stokes equations. We consider only distributed control here, although we believe other problems could be treated in the same way. We give numerical results, and we compare these with those obtained by solving the equivalent forward problem using similar techniques.

Select this Article		On-the-Fly Adaptive Smoothed Aggregation Multigrid for Markov Chains Eran Treister and Irad Yavneh SIAM J. Sci. Comput. 33, pp. 2927-2949 (23 pages) Online Publication Date: October 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract A new adaptive algebraic multigrid scheme is developed for the solution of Markov chains, where the hierarchy of operators is adapted on-the-fly in a setup process that is interlaced with the solution process. The setup process feeds the solution process with improved operators, while the solution process provides the adaptive setup process with better approximations on which to base further-improved operators. The approach is demonstrated using Petrov–Galerkin smoothed aggregation where only the prolongation operator is smoothed, while the restriction remains of low order. Results show that the on-the-fly adaptive scheme can improve the performance of multigrid solvers that require extensive setup computations, in both serial and parallel environments.

Select this Article		LSMR: An Iterative Algorithm for Sparse Least-Squares Problems David Chin-Lung Fong and Michael Saunders SIAM J. Sci. Comput. 33, pp. 2950-2971 (22 pages) Online Publication Date: October 27, 2011 Full Text: \| Download PDF
	+ -	Show Abstract An iterative method LSMR is presented for solving linear systems $Ax=b$ and least-squares problems $\min \\|Ax-b\\|_2$, with $A$ being sparse or a fast linear operator. LSMR is based on the Golub–Kahan bidiagonalization process. It is analytically equivalent to the MINRES method applied to the normal equation $A^T\! Ax = A^T\! b$, so that the quantities $\\|A^T\! r_k\\|$ are monotonically decreasing (where $r_k = b - Ax_k$ is the residual for the current iterate $x_k$). We observe in practice that $\\|r_k\\|$ also decreases monotonically, so that compared to LSQR (for which only $\\|r_k\\|$ is monotonic) it is safer to terminate LSMR early. We also report some experiments with reorthogonalization.