SIAM Digital Library
 
 
 

SIAM J. on Scientific Computing

SECTION LISTING

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

Search Issue | RSS Feeds RSS
Next Issue

2010

Volume 32, Issue 1, pp. vii-475

* Special Issue: 2008 Copper Mountain Conference

back to top
RSS Feeds

Special Issue: 2008 Copper Mountain Conference

Ray Tuminaro

SIAM J. Sci. Comput. 32, pp. vii-vii ( pages)

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
The biennial Copper Mountain Conference on Iterative Methods was held April 6–11, 2008. This most recent meeting included both theoretical and applied presentations spanning a broad range of iterative/linear algebra topics such as saddle-point systems, eigensolvers, nonlinear solvers, inverse problems, and block solution methods. The diverse impact of iterative method research can be gauged by the broad range of application areas represented at the meeting. Some of these include computational fluid/gas flow, large-scale optimization, stochastic finite elements, porous media, image processing, magnetohydrodynamics, fluid-structure interaction, Markov processes, solar cells, and blood flow.
Howard Elman (University of Maryland) and Panayot Vassilevski (Center for Applied Scientific Computing at Lawrence Livermore National Laboratory) served as co-chairs. Management services were provided by Annette Anthony of Front Range Scientific Computations, Inc. The conference gratefully acknowledges support from Boeing, the Department of Energy, IBM, Lawrence Livermore National Laboratory, Los Alamos National Laboratory, 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 issue were open to the public. This was advertised in advance on web sites of both the SIAM Journal on Scientific Computing (SISC) and the Copper Mountain Conference.
The papers in this issue reflect leading-edge research in iterative methods. Many topics are covered including incomplete factorizations, convergence analysis, and multiple right hand side linear systems. Contributions also appear in areas such as optimization, inverse problems, signal restoration, and stochastic partital differential equations representing a broadening interest within these communities in iterative solvers. Additionally, this issue contains a healthy collection of papers in newer iterative solver topics associated with self-adapting methods, Markov chains, clustering, pagerank, magnetohydrodynamics, interface tracking, and quantum chromodynamics. Some of these submissions also appeared in the student paper competition, a unique event at the Copper Mountain conferences that has proven to be a very successful forum for presenting the latest results in the field. We have had increasingly strong participation from very talented recent Ph.D. holders and current Ph.D. students in this competition. This special issue clearly demonstrates the strong impact of iterative methods in various areas of applications and the need for continued research.
I would like to thank the guest editorial board, the SISC editorial board, and the anonymous referees for their hard work to help meet the strict deadlines. I would also like to acknowledge the SIAM staff for their efforts. Special appreciation is due to Mitch Chernoff (SIAM Publications Manager) and Colleen Hills (SIAM Editorial Assistant) for their work on this special issue.
We hope that this special issue will stimulate further progress in this exciting area of research.

General Constrained Energy Minimization Interpolation Mappings for AMG

Panayot S. Vassilevski

SIAM J. Sci. Comput. 32, pp. 1-13 (13 pages) | Cited 5 times

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
This report proposes a new class of interpolation procedures for use in algebraic multigrid (AMG). The procedure is general in that it applies to s.p.d. finite element discretization problems that include scalar and vector elliptic PDEs, and it can be adapted to the time-domain Maxwell (definite $H(\mathrm{curl})$) problems. It can also be viewed as an extension of previously proposed vector-fitting interpolation procedures and can be used in adaptive AMG cycles. We illustrate the performance of the new interpolation matrices on a number of test problems.

Towards Adaptive Smoothed Aggregation ($\alpha$SA) for Nonsymmetric Problems

M. Brezina, T. Manteuffel, S. MCormick, J. Ruge, and G. Sanders

SIAM J. Sci. Comput. 32, pp. 14-39 (26 pages) | Cited 2 times

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
Applying smoothed aggregation (SA) multigrid to solve a nonsymmetric linear system, $A\mathbf{x} =\mathbf{b}$, is often impeded by the lack of a minimization principle that can be used as a basis for the coarse-grid correction process. This paper proposes a Petrov–Galerkin (PG) approach based on applying SA to either of two symmetric positive definite (SPD) matrices, $\sqrt{A^{t}A}$ or $\sqrt{AA^{t}}$. These matrices, however, are typically full and difficult to compute, so it is not computationally efficient to use them directly to form a coarse-grid correction. The proposed approach approximates these coarse-grid corrections by using SA to accurately approximate the right and left singular vectors of $A$ that correspond to the lowest singular value. These left and right singular vectors are used to construct the restriction and interpolation operators, respectively. A preliminary two-level convergence theory is presented, suggesting that more relaxation should be applied than for an SPD problem. Additionally, a nonsymmetric version of adaptive SA ($\alpha$SA) is given that automatically constructs SA multigrid hierarchies using a stationary relaxation process on all levels. Numerical results are reported for convection-diffusion problems in two dimensions with varying amounts of convection for constant, variable, and recirculating convection fields. The results suggest that the proposed approach is algorithmically scalable for problems coming from these nonsymmetric scalar PDEs (with the exception of recirculating flow). This paper serves as a first step for nonsymmetric $\alpha$SA. The long-term goal of this effort is to develop nonsymmetric $\alpha$SA for systems of PDEs, where the SA framework has proven to be well suited for adaptivity in SPD problems.

Smoothed Aggregation Multigrid for Markov Chains

H. De Sterck, T. A. Manteuffel, S. F. McCormick, K. Miller, J. Pearson, J. Ruge, and G. Sanders

SIAM J. Sci. Comput. 32, pp. 40-61 (22 pages) | Cited 4 times

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
A smoothed aggregation multigrid method is presented for the numerical calculation of the stationary probability vector of an irreducible sparse Markov chain. It is shown how smoothing the interpolation and restriction operators can dramatically increase the efficiency of aggregation multigrid methods for Markov chains that have been proposed in the literature. The proposed smoothing approach is inspired by smoothed aggregation multigrid for linear systems, supplemented with a new lumping technique that assures well-posedness of the coarse-level problems: the coarse-level operators are singular M-matrices on all levels, resulting in strictly positive coarse-level corrections on all levels. Numerical results show how these methods lead to nearly optimal multigrid efficiency for an extensive set of test problems, both when geometric and algebraic aggregation strategies are used.

Finite Element Approach to Clustering of Multidimensional Time Series

Illia Horenko

SIAM J. Sci. Comput. 32, pp. 62-83 (22 pages) | Cited 4 times

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
We present a new approach to clustering of time series based on a minimization of the averaged clustering functional. The proposed functional describes the mean distance between observation data and its representation in terms of $\mathbf{K}$ abstract models of a certain predefined class (not necessarily given by some probability distribution). For a fixed time series $x(t)$ this functional depends on $\mathbf{K}$ sets of model parameters $\Theta=(\theta_1,\dots,\theta_\mathbf{K})$ and $\mathbf{K}$ functions of cluster affiliations $\Gamma=(\gamma_1(t),\dots,\gamma_{\mathbf{K}}(t))$ (characterizing the affiliation of any element $x(t)$ of the analyzed time series to one of the $\mathbf{K}$ clusters defined by the considered model parameters). We demonstrate that for a fixed set of model parameters $\Theta$ the appropriate Tykhonov-type regularization of this functional with some regularization factor $\epsilon^2$ results in a minimization problem similar to a variational problem usually associated with one-dimensional nonhomogeneous partial differential equations. This analogy allows us to apply the finite element framework to the problem of time series analysis and to propose a numerical scheme for time series clustering. We investigate the conditions under which the proposed scheme allows a monotone improvement of the initial parameter guess with respect to the minimization of the discretized version of the regularized functional. We also discuss the interpretation of the regularization factor in the Markovian case and show its connection to metastability and exit times. The computational performance of the resulting method is investigated numerically on multidimensional test data and is applied to the analysis of multidimensional historical stock market data.

Adaptive Techniques for Improving the Performance of Incomplete Factorization Preconditioning

Anshul Gupta and Thomas George

SIAM J. Sci. Comput. 32, pp. 84-110 (27 pages)

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
Three techniques for improving the robustness and performance of incomplete factorization preconditioners for sparse systems with symmetric positive definite or mildly indefinite coefficient matrices are introduced. The primary contribution is two new block algorithms for incomplete factorization that result in an improvement in the performance of both the preconditioner generation and the iterative solution phases. One of the algorithms applies to matrices that have a natural block structure in their original form, and the other applies to matrices without natural dense blocks. Additionally, two relatively simple but highly effective heuristic strategies are introduced. These include selecting the solver based on the definiteness properties of the preconditioner and automatic selection and tuning of incomplete factorization parameters. All three techniques have adaptive components; i.e., the preconditioner-solver combination chooses parameters or algorithmic components based on the properties of the coefficient matrix and its incomplete factors. Two of the three techniques are applicable to incomplete LU factorization for unsymmetric systems as well.

Adaptive Time-Stepping for Incompressible Flow Part II: Navier–Stokes Equations

David A. Kay, Philip M. Gresho, David F. Griffiths, and David J. Silvester

SIAM J. Sci. Comput. 32, pp. 111-128 (18 pages)

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
We outline a new class of robust and efficient methods for solving the Navier–Stokes equations. We describe a general solution strategy that has two basic building blocks: an implicit time integrator using a stabilized trapezoid rule with an explicit Adams–Bashforth method for error control, and a robust Krylov subspace solver for the spatially discretized system. We present numerical experiments illustrating the potential of our approach.

Deflated and Restarted Symmetric Lanczos Methods for Eigenvalues and Linear Equations with Multiple Right-Hand Sides

Abdou M. Abdel-Rehim, Ronald B. Morgan, Dywayne A. Nicely, and Walter Wilcox

SIAM J. Sci. Comput. 32, pp. 129-149 (21 pages) | Cited 1 time

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
A deflated restarted Lanczos algorithm is given for both solving symmetric linear equations and computing eigenvalues and eigenvectors. The restarting limits the storage so that finding eigenvectors is practical. Meanwhile, the deflating from the presence of the eigenvectors allows the linear equations to generally have good convergence in spite of the restarting. Some reorthogonalization is necessary to control roundoff error, and several approaches are discussed. The eigenvectors generated while solving the linear equations can be used to help solve systems with multiple right-hand sides. Experiments are given with large matrices from quantum chromodynamics that have many right-hand sides.

Operator-Based Preconditioning of Stiff Hyperbolic Systems

Daniel R. Reynolds, Ravi Samtaney, and Carol S. Woodward

SIAM J. Sci. Comput. 32, pp. 150-170 (21 pages)

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
We introduce an operator-based scheme for preconditioning stiff components encountered in implicit methods for hyperbolic systems of PDEs posed on regular grids. The method is based on a directional splitting of the implicit operator, followed by a characteristic decomposition of the resulting directional parts. This approach allows for the solution of any number of characteristic components, from the entire system to only the fastest, stiffness-inducing waves. We apply the preconditioning method to stiff hyperbolic systems arising in magnetohydrodynamics and gas dynamics. We then present numerical results showing that this preconditioning scheme works well on problems where the underlying stiffness results from the interaction of fast transient waves with slowly-evolving dynamics, scales well to large problem sizes and numbers of processors, and allows for additional customization based on the specific problems under study.

An Iterative Method for Edge-Preserving MAP Estimation When Data-Noise Is Poisson

Johnathan M. Bardsley and John Goldes

SIAM J. Sci. Comput. 32, pp. 171-185 (15 pages)

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
In numerous applications of image processing, e.g., astronomical and medical imaging, data-noise is well modeled by a Poisson distribution. This motivates the use of the negative-log Poisson likelihood function for data fitting. However, difficulties arise when this likelihood is used due to the fact that it is nonquadratic and requires a nonnegativity constraint when minimized. Moreover, for ill-posed problems, which are our interest here, regularization is required. In this paper we present an edge-preserving, quadratic regularization function. The resulting regularized, negative-log Poisson likelihood minimization problem is solved with an iterative method that we present here. We prove that the method is convergent when applied to the problem of interest. Finally, in order to test the methodology, we apply the algorithm and edge-preserving regularization scheme to synthetically generated data.

The Cycle-Convergence of Restarted GMRES for Normal Matrices Is Sublinear

E. Vecharynski and J. Langou

SIAM J. Sci. Comput. 32, pp. 186-196 (11 pages)

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
We prove that the cycle-convergence of the restarted GMRES applied to a system of linear equations with a normal coefficient matrix is sublinear.

An Efficient Numerical Method for the Solution of the $L_2$ Optimal Mass Transfer Problem

Eldad Haber, Tauseef Rehman, and Allen Tannenbaum

SIAM J. Sci. Comput. 32, pp. 197-211 (15 pages)

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
In this paper we present a new computationally efficient numerical scheme for the minimizing flow approach for the computation of the optimal $L_2$ mass transport mapping. In contrast to the integration of a time dependent partial differential equation proposed in [S. Angenent, S. Haker, and A. Tannenbaum, SIAM J. Math. Anal., 35 (2003), pp. 61–97], we employ in the present work a direct variational method. The efficacy of the approach is demonstrated on both real and synthetic data.

Dynamically Adaptive Simulations with Minimal Memory Requirement—Solving the Shallow Water Equations Using Sierpinski Curves

Michael Bader, Christian Böck, Johannes Schwaiger, and Csaba Vigh

SIAM J. Sci. Comput. 32, pp. 212-228 (17 pages)

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
We present an approach to the numerical simulation of dynamically adaptive problems on recursively structured adaptive triangular grids. The intended application is the simulation of oceanic wave propagation (e.g., tsunami simulation) based on the shallow water equations. For the required 2D dynamically adaptive discretization, we adopt a grid generation process based on recursive bisection of triangles along marked edges. The recursive grid generation may be described via a respective refinement tree, which is sequentialized according to a Sierpinski space-filling curve. This allows a storage scheme for the adaptive grid that requires only a minimal amount of memory. Moreover, the sequentialization and, hence, the locality properties induced by the space-filling curve are retained throughout adaptive refinement and coarsening of the grid. Conforming adaptive refinement and coarsening, as well as time-stepping techniques for time-dependent systems of partial differential equations, are implemented using an inherently cache-efficient processing scheme, which is based on the use of stacks and stream-like data structures and a traversal of the adaptively refined grid along the Sierpinski curve. We demonstrate the computational efficiency of the approach on the solution of a simplified version of the shallow water equations, for which we use a discontinuous Galerkin discretization. Special attention is paid to the memory efficiency of the implementation.

First-Order System Least Squares for Incompressible Resistive Magnetohydrodynamics

J. H. Adler, T. A. Manteuffel, S. F. McCormick, and J. W. Ruge

SIAM J. Sci. Comput. 32, pp. 229-248 (20 pages) | Cited 2 times

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
Magnetohydrodynamics (MHD) is a fluid theory that describes plasma physics by treating the plasma as a fluid of charged particles. Hence, the equations that describe the plasma form a nonlinear system that couples Navier–Stokes equations with Maxwell's equations. This paper shows that the first-order system least squares (FOSLS) finite element method is a viable discretization for these large MHD systems. To solve this system, a nested-iteration–Newton–FOSLS–AMG approach is taken. Most of the work is done on the coarse grid, including most of the linearizations. We show that at most one Newton step and a few V-cycles are all that are needed on the finest grid. Here, we describe how the FOSLS method can be applied to incompressible resistive MHD and how it can be used to solve these MHD problems efficiently. A 3D steady state and a reduced 2D time-dependent test problem are studied. The latter equations can simulate a “large aspect-ratio” tokamak. The goal is to resolve as much physics from the test problems with the least amount of computational work. We show that this is achieved in a few dozen work units or fine grid residual evaluations.

Preconditioning Saddle-Point Systems with Applications in Optimization

H. Sue Dollar, Nicholas I. M. Gould, Martin Stoll, and Andrew J. Wathen

SIAM J. Sci. Comput. 32, pp. 249-270 (22 pages) | Cited 2 times

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
Saddle-point systems arise in many applications areas, in fact in any situation where an extremum principle arises with constraints. The Stokes problem describing slow viscous flow of an incompressible fluid is a classic example coming from PDEs and in the area of optimization such problems are ubiquitous. In this paper we present a framework into which many well-known methods for solving saddle-point systems fit. Based on this description we show how new approaches for the solution of saddle-point systems arising in optimization can be derived from the Bramble–Pasciak conjugate gradient approach widely used in PDEs and more recent generalizations thereof. In particular we derive a class of new solution methods based on the use of preconditioned conjugate gradients in nonstandard inner products and demonstrate how these can be understood through more standard machinery. We show connections to constraint preconditioning and give the results of numerical computations on a number of standard optimization test examples.

Optimal Solvers for PDE-Constrained Optimization

Tyrone Rees, H. Sue Dollar, and Andrew J. Wathen

SIAM J. Sci. Comput. 32, pp. 271-298 (28 pages) | Cited 4 times

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
Optimization problems with constraints which require the solution of a partial differential equation arise widely in many areas of the sciences and engineering, particularly in problems of design. The solution of such PDE-constrained optimization problems is usually a major computational task. Here we consider simple problems of this type: distributed control problems in which the 2- and 3-dimensional Poisson problem is the PDE. The large-dimensional linear systems which result from discretization and which need to be solved are of saddle-point type. We introduce two optimal preconditioners for these systems, which lead to convergence of symmetric Krylov subspace iterative methods in a number of iterations which does not increase with the dimension of the discrete problem. These preconditioners are block structured and involve standard multigrid cycles. The optimality of the preconditioned iterative solver is proved theoretically and verified computationally in several test cases. The theoretical proof indicates that these approaches may have much broader applicability for other PDEs.

Multilevel Approach For Signal Restoration Problems With Toeplitz Matrices

Malena I. Español and Misha E. Kilmer

SIAM J. Sci. Comput. 32, pp. 299-319 (21 pages)

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
We present a multilevel method for discrete ill-posed problems arising from the discretization of Fredholm integral equations of the first kind. In this method, we use the Haar wavelet transform to define restriction and prolongation operators within a multigrid-type iteration. The choice of the Haar wavelet operator has the advantage of preserving matrix structure, such as Toeplitz, between grids, which can be exploited to obtain faster solvers on each level where an edge-preserving Tikhonov regularization is applied. Finally, we present results that indicate the promise of this approach for restoration of signals and images with edges.

Marker Redistancing/Level Set Method for High-Fidelity Implicit Interface Tracking

Robert Nourgaliev, Samet Kadioglu, and Vincent Mousseau

SIAM J. Sci. Comput. 32, pp. 320-348 (29 pages)

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
A hybrid of the front tracking (FT) and the level set (LS) methods is introduced, combining advantages and removing drawbacks of both methods. The kinematics of the interface is treated in a Lagrangian (FT) manner, by tracking markers placed at the interface. The markers are not connected—instead, the interface topology is resolved in an Eulerian (LS) framework, by wrapping a signed distance function around Lagrangian markers each time the markers move. For accuracy and efficiency, we have developed a high-order “anchoring” algorithm and an implicit PDE-based redistancing. We have demonstrated that the method is 3rd-order accurate in space, near the markers, and therefore 1st-order convergent in curvature; this is in contrast to traditional PDE-based reinitialization algorithms, which tend to slightly relocate the zero level set and can be shown to be nonconvergent in curvature. The implicit pseudo-time discretization of the redistancing equation is implemented within the Jacobian-free Newton–Krylov (JFNK) framework combined with ILU(k) preconditioning. Due to the LS localization, the bandwidth of the Jacobian matrix is nearly constant, and the ILU preconditioning scales as $\sim N\log(\sqrt{N})$ in two dimensions, which implies efficiency and good scalability of the overall algorithm. We have demonstrated that the steady-state solutions in pseudo-time can be achieved very efficiently, with $\approx10$ iterations ($\mathrm{CFL}\approx10^4$), in contrast to the explicit redistancing which requires hundreds of iterations with $\mathrm{CFL}\leq1$.

An Inner-Outer Iteration for Computing PageRank

David F. Gleich, Andrew P. Gray, Chen Greif, and Tracy Lau

SIAM J. Sci. Comput. 32, pp. 349-371 (23 pages)

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
We present a new iterative scheme for PageRank computation. The algorithm is applied to the linear system formulation of the problem, using inner-outer stationary iterations. It is simple, can be easily implemented and parallelized, and requires minimal storage overhead. Our convergence analysis shows that the algorithm is effective for a crude inner tolerance and is not sensitive to the choice of the parameters involved. The same idea can be used as a preconditioning technique for nonstationary schemes. Numerical examples featuring matrices of dimensions exceeding 100,000,000 in sequential and parallel environments demonstrate the merits of our technique. Our code is available online for viewing and testing, along with several large scale examples.

Iterative Solvers for the Stochastic Finite Element Method

Eveline Rosseel and Stefan Vandewalle

SIAM J. Sci. Comput. 32, pp. 372-397 (26 pages) | Cited 2 times

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
This paper presents an overview and comparison of iterative solvers for linear stochastic partial differential equations (PDEs). A stochastic Galerkin finite element discretization is applied to transform the PDE into a coupled set of deterministic PDEs. Specialized solvers are required to solve the very high-dimensional systems that result after a finite element discretization of the resulting set. This paper discusses one-level iterative methods, based on matrix splitting techniques; multigrid methods, which apply a coarsening in the spatial dimension; and multilevel methods, which make use of the hierarchical structure of the stochastic discretization. Also Krylov solvers with suitable preconditioning are addressed. A local Fourier analysis provides quantitative convergence properties. The efficiency and robustness of the methods are illustrated on two nontrivial numerical problems. The multigrid solver with block smoother yields the most robust convergence properties, though a cheaper point smoother performs as well in most cases. Multilevel methods based on coarsening the stochastic dimension perform in general poorly due to a large computational cost per iteration. Moderate size problems can be solved very quickly by a Krylov method with a mean-based preconditioner. For larger spatial and stochastic discretizations, however, this approach suffers from its nonoptimal convergence properties.

Least-Squares Finite Element Methods for Quantum Electrodynamics

J. Brannick, C. Ketelsen, T. Manteuffel, and S. McCormick

SIAM J. Sci. Comput. 32, pp. 398-417 (20 pages)

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
A significant amount of the computational time in large Monte Carlo simulations of lattice field theory is spent inverting the discrete Dirac operator. Unfortunately, traditional covariant finite difference discretizations of the Dirac operator present serious challenges for standard iterative methods. For interesting physical parameters, the discretized operator is large and ill-conditioned and has random coefficients. More recently, adaptive algebraic multigrid (AMG) methods have been shown to be effective preconditioners for Wilson's discretization [J. Brannick et al., Lect. Notes Comput. Sci. Eng., 55 (2006), pp. 499–506], [J. Brannick et al., Phys. Rev. Lett., 100 (2008), pp. 041601–041604] of the Dirac equation. This paper presents an alternate discretization of the two-dimensional Dirac operator of quantum electrodynamics (QED) based on least-squares finite elements. The discretization is systematically developed and, physical properties of the resulting matrix system are discussed. Finally, numerical experiments are presented that demonstrate the effectiveness of adaptive smoothed aggregation ($\alpha$SA) multigrid as a preconditioner for the discrete field equations.

A Fully Implicit Domain Decomposition Algorithm for Shallow Water Equations on the Cubed-Sphere

Chao Yang, Jianwen Cao, and Xiao-Chuan Cai

SIAM J. Sci. Comput. 32, pp. 418-438 (21 pages)

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
Popular approaches for solving shallow water equations (SWEs) for climate modeling are explicit and semiimplicit methods, and both have certain constraints on the time step size. In this paper, we propose and study a fully implicit method which imposes no limit on the time step size but requires the solution of a large sparse nonlinear system of equations at every time step. The focus of the paper is a parallel, fully coupled, Newton–Krylov–RAS algorithm with a Jacobian matrix explicitly calculated on a weakly nonmatching cubed-sphere mesh. Here, RAS is a restricted additive Schwarz method. We show numerically that with such a preconditioned, nonlinearly implicit method the time step size is no longer constrained by the CFL condition, and we report superlinear speedup of the algorithm on machines with thousands of processors for problems with smooth and nonsmooth solutions.

Computing and Deflating Eigenvalues While Solving Multiple Right-Hand Side Linear Systems with an Application to Quantum Chromodynamics

Andreas Stathopoulos and Konstantinos Orginos

SIAM J. Sci. Comput. 32, pp. 439-462 (24 pages) | Cited 1 time

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
We present a new algorithm that computes eigenvalues and eigenvectors of a Hermitian positive definite matrix while solving a linear system of equations with conjugate gradient (CG). Traditionally, all the CG iteration vectors could be saved and recombined through the eigenvectors of the tridiagonal projection matrix, which is equivalent theoretically to unrestarted Lanczos. Our algorithm capitalizes on the iteration vectors produced by CG to update only a small window of vectors that approximate the eigenvectors. While this window is restarted in a locally optimal way, the CG algorithm for the linear system is unaffected. Yet, in all our experiments, if the window has more than a properly chosen but small number of vectors, it converges to the required eigenvectors at a rate identical to unrestarted Lanczos. After the solution of the linear system, eigenvectors that have not accurately converged can be improved in an incremental fashion by solving additional linear systems. In this case, eigenvectors identified in earlier systems can be used to deflate, and thus accelerate, the convergence of subsequent systems. We have used this algorithm with excellent results in lattice quantum chromodynamics applications, where hundreds of right-hand sides may be needed. Specifically, about 70 eigenvectors are obtained to full accuracy after solving 24 right-hand sides. Deflating these from the large number of subsequent right-hand sides removes the dreaded critical slowdown, where the conditioning of the matrix increases as the quark mass reaches a critical value. Our experiments show almost a constant number of iterations for our method, regardless of quark mass, and speedups of 8 over original CG for light quark masses.

The Iterative Solver RISOLV with Application to the Exterior Helmholtz Problem

Hillel Tal-Ezer and Eli Turkel

SIAM J. Sci. Comput. 32, pp. 463-475 (13 pages)

Online Publication Date: February 05, 2010

Full Text: | Download PDF

Show Abstract
The innermost computational kernel of many large-scale scientific applications is often a large set of linear equations of the form $Ax=b$ which typically consumes a significant portion of the overall computational time required by the simulation. The traditional approach for solving this problem is to use direct methods. This approach is often preferred in industry because direct solvers are robust and effective for moderate size problems. However, direct methods can consume a huge amount of memory, and CPU time, in large-scale cases. In these cases, iterative techniques are the only viable alternative. Unfortunately, iterative methods lack the robustness of direct methods. The situation is especially difficult when the matrix is nonsymmetric. A lot of research has been devoted to trying to develop a robust iterative algorithm for nonsymmetric systems. The present paper describes a new robust and efficient algorithm aimed at solving iteratively nonsymmetric linear systems. It is based on looking for an approximation to the “optimal” polynomial $P_m(z)$ which satisfies $||P_m(z)||_{\infty}=\min_{Q\in\Pi_m}||Q(z)||_{\infty}$, $z\in D$, where $\Pi_m$ is the set of all polynomials of degree $m$ which satisfies $Q_m(0)=1$ and $D$ is a domain in the complex plane which includes all the eigenvalues of $A$. The resulting algorithm is an efficient one, especially in the case where we have a set of linear systems which share the same matrix $A$. We present several applications, including the exterior Helmholtz problem, which leads to a large indefinite, nonsymmetric, and complex system.
Close

Your Recent History Recent History Help

Recently Viewed

  1. Computational Solution of Blow-Up Problems for Semilinear Parabolic PDEs on Unbounded Domains
  2. Extrapolated Implicit-Explicit Time Stepping
  3. GPU-Based Volume Reconstruction from Very Few Arbitrarily Aligned X-Ray Images
  4. Multivariate Regression and Machine Learning with Sums of Separable Functions
  5. Construction of Data-Sparse $\mathcal{H}^2$-Matrices by Hierarchical Compression
© SIAM By using SIAM Publications Online you agree to abide by the Terms and Conditions of Use. Banner art adapted from a figure by Hinke M. Osinga and Bernd Krauskopf (University of Bristol, UK).
close