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

1998

Volume 19, Issue 6, pp. 1767-2110


MINRES and MINERR Are Better than SYMMLQ in Eigenpair Computations

Franciszek A. Dul

SIAM J. Sci. Comput. 19, pp. 1767-1782 (16 pages) | Cited 6 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
The SYMMLQ and MINRES methods for solving sparse, indefinite, symmetric systems of equations [SIAM J. Numer. Anal., 12 (1975), pp. 617--629] can be used inside the inverse or Rayleigh quotient iterations for computing the eigenpairs of large, sparse, generalized eigenproblems $A y = \lambda B y$. It is shown in this paper that the common opinion concerning the superiority of SYMMLQ over MINRES is not true in this context. Several results taken from quantum chemistry, acoustics, and structural mechanics are presented that show the superiority of MINRES and the new MINERR method over SYMMLQ, especially for ill-conditioned eigenproblems, where SYMMLQ may even diverge. Theoretical explanations of the observed phenomena are given. The use of MINRES or MINERR instead of SYMMLQ may thus lead to a significant increase in efficiency in large scale eigenvalue computations.

Numerical Conformal Mapping Using Cross-Ratios and Delaunay Triangulation

Tobin A. Driscoll and Stephen A. Vavasis

SIAM J. Sci. Comput. 19, pp. 1783-1803 (21 pages) | Cited 15 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
We propose a new algorithm for computing the Riemann mapping of the unit disk to a polygon, also known as the Schwarz--Christoffel transformation. The new algorithm, CRDT (for cross-ratios of the Delaunay triangulation), is based on cross-ratios of the prevertices, and also on cross-ratios of quadrilaterals in a Delaunay triangulation of the polygon. The CRDT algorithm produces an accurate representation of the Riemann mapping even in the presence of arbitrary long, thin regions in the polygon, unlike any previous conformal mapping algorithm. We believe that CRDT solves all difficulties with crowding and global convergence, although these facts depend on conjectures that we have so far not been able to prove. We demonstrate convergence with computational experiments. The Riemann mapping has applications in two-dimensional potential theory and mesh generation. We demonstrate CRDT on problems in long, thin regions in which no other known algorithm can perform comparably.

An Improved Fast Multipole Algorithm for Potential Fields

Tomasz Hrycak and Vladimir Rokhlin

SIAM J. Sci. Comput. 19, pp. 1804-1826 (23 pages) | Cited 27 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
A new version of the fast multipole method (FMM) for potential fields is presented. We introduce a new representation of potentials, in which most translation operators are diagonal. As a result, for double precision calculations in two dimensions we obtain an improvement of a factor of two to four in speed, compared to previously published algorithms; the improvement is expected to be much greater in three dimensions. The performance of the method is illustrated with several numerical examples.

A High-Order Iterative Implicit-Explicit Hybrid Scheme for Magnetohydrodynamics

Wenlong Dai and Paul R. Woodward

SIAM J. Sci. Comput. 19, pp. 1827-1846 (20 pages) | Cited 3 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
An iterative implicit-explicit hybrid scheme is proposed for one-dimensional ideal magnetohydrodynamical (MHD) equations. The scheme is accurate to second order in both space and time for all Courant numbers. Each kind of MHD wave, both compressible and incompressible, may be either implicitly or explicitly, or partially implicitly and partially explicitly treated depending on their associated Courant numbers in each numerical cell, and the scheme is able to smoothly switch between explicit and implicit calculations. The scheme is of Godunov type in both explicit and implicit regimes, in which the flux needed in a Godunov scheme is calculated from Riemann problems, is in a strictly conservation form, and is able to resolve MHD discontinuities. Only a single level of iterations is involved in the scheme, and the iterations solve both the implicit relations arising from upstream centered differences for all wave families and the nonlinearity of MHD equations. Multicolors proposed in the paper may significantly reduce the number of iterations required to reach a converged solution. Compared with the largest Courant number in a simulation, only a small number of iterations are needed in the scheme. The feature of the scheme has been shown through numerical examples. It is easy to vectorize the computer code of the scheme.

Optimal Control of Thermally Convected Fluid Flows

K. Ito and S. S. Ravindran

SIAM J. Sci. Comput. 19, pp. 1847-1869 (23 pages) | Cited 30 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
We examine the optimal control of stationary thermally convected fluid flows from the theoretical and numerical point of view. We use thermal convection as control mechanism; that is, control is effected through the temperature on part of the boundary. Control problems are formulated as constrained minimization problems. Existence of optimal control is given and a first-order necessary condition of optimality from which optimal solutions can be obtained is established. We develop numerical methods to solve the necessary condition of optimality and present computational results for control of cavity- and channel-type flows showing the feasibility of the proposed approach.

Robustness of an Elementwise Parallel Finite Element Method for Convection-Diffusion Problems

W. J. Layton, J. M. Maubach, and P. J. Rabier

SIAM J. Sci. Comput. 19, pp. 1870-1891 (22 pages) | Cited 4 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
We consider an elementwise data-parallel finite element procedure, recently proposed by Layton and Rabier [Appl. Math. Lett., 5 (1992), pp. 67--70], [J. Numer. Linear Algebra Appl., 2 (1995), pp. 363--394], applied to singularly perturbed convection-diffusion equations with possibly highly anisotropic diffusion. It is shown that the number of iterations required for the solution of the linear algebraic system is proportional to the inverse of the smallest grid element diameter, uniformly in the diffusion parameter and the degree of anisotropy. This is optimal, since the method can in some cases use only element matrices and load vectors and the algorithm requires only local communication on the physical mesh between adjacent elemental subdomains. Our analysis includes both the usual Galerkin formulation and the streamline upwind finite element formulations. The convergence result holds with conforming elements of any order or elements of arbitrary order from a new family of nonconforming elements.

Nonoscillatory Central Schemes for Multidimensional Hyperbolic Conservation Laws

Guang-Shan Jiang and Eitan Tadmor

SIAM J. Sci. Comput. 19, pp. 1892-1917 (26 pages) | Cited 107 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
We construct, analyze, and implement a new nonoscillatory high-resolution scheme for two-dimensional hyperbolic conservation laws. The scheme is a predictor-corrector method which consists of two steps: starting with given cell averages, we first predict pointvalues which are based on nonoscillatory piecewise-linear reconstructions from the given cell averages; at the second corrector step, we use staggered averaging, together with the predicted midvalues, to realize the evolution of these averages. This results in a second-order, nonoscillatory central scheme, a natural extension of the one-dimensional second-order central scheme of Nessyahu and Tadmor [J. Comput. Phys., 87 (1990), pp. 408--448].
As in the one-dimensional case, the main feature of our two-dimensional scheme is simplicity. In particular, this central scheme does not require the intricate and time-consuming (approximate) Riemann solvers which are essential for the high-resolution upwind schemes; in fact, even the computation of the exact Jacobians can be avoided. Moreover, the central scheme is "genuinely multidimensional" in the sense that it does not necessitate dimensional splitting.
We prove that the scheme satisfies the scalar maximum principle, and in the more general context of systems, our proof indicates that the scheme is positive (in the sense of Lax and Liu [CFD Journal, 5 (1996), pp. 1--24]). We demonstrate the application of our central scheme to several prototype two-dimensional Euler problems. Our numerical experiments include the resolution of shocks oblique to the computational grid; they show how our central scheme solves with high resolution the intricate wave interactions in the so-called double Mach reflection problem [J. Comput. Phys., 54 (1988), pp. 115--173] without following the characteristics; and finally we report on the accurate ray solutions of a weakly hyperbolic system [J. Comput. Appl. Math., 74 (1996), pp. 175--192], rays which otherwise are missed by the dimensional splitting approach. Thus, a considerable amount of simplicity and robustness is gained while achieving stability and high resolution.

Iterative Solution of Cyclically Reduced Systems Arising from Discretization of the Three-Dimensional Convection-Diffusion Equation

Chen Greif and James Varah

SIAM J. Sci. Comput. 19, pp. 1918-1940 (23 pages) | Cited 6 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
We consider the system of equations arising from finite difference discretization of a three-dimensional convection-diffusion model problem. This system is typically nonsymmetric. We show that performing one step of cyclic reduction, followed by reordering of the unknowns, yields a system of equations for which the block Jacobi method generally converges faster than for the original system, using lexicographic ordering. The matrix representing the system of equations can be symmetrized for a large range of the coefficients of the underlying partial differential equation, and the associated iteration matrix has a smaller spectral radius than the one associated with the original system. In this sense, the three-dimensional problem is similar to the one-dimensional and two-dimensional problems, which have been studied by Elman and Golub. The process of reduction, the suggested orderings, and bounds on the spectral radii of the associated iteration matrices are presented, followed by a comparison of the reduced system with the full system and by details of the numerical experiments.

Preconditioned Mixed Spectral Element Methods for Elasticity and Stokes Problems

Luca F. Pavarino

SIAM J. Sci. Comput. 19, pp. 1941-1957 (17 pages) | Cited 3 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
This paper continues our study of efficient preconditioners for the indefinite systems obtained by discretizing the linear elasticity and Stokes problems with mixed spectral methods in three dimensions. Our study is extended here in two new directions, to include i) more general (multidomain) mixed spectral element methods and ii) a new class of block-triangular preconditioners. The stiffness matrices resulting from Stokes and mixed elasticity problems have the structure of saddle point problems with a penalty term, which is associated with the Poisson ratio for elasticity problems or with stabilization techniques for Stokes problems. The main result of this paper shows that the convergence rate of the proposed preconditioned iterative methods is independent of the penalty parameter $\nu$, the number of spectral elements N, and mildly dependent on the spectral degree n via the inf-sup constant. The number of iterations required by the block-triangular preconditioner is much smaller than the one required by the block-diagonal preconditioner. On the other hand, the block-triangular preconditioner is no longer symmetric and has a small extra cost. These results and several numerical experiments presented in the final section show that these algorithms are a practical and efficient strategy for the iterative solution of the indefinite systems arising from mixed spectral element discretizations of elliptic systems.

First-Order System Least Squares (FOSLS) for Convection-Diffusion Problems: Numerical Results

J. M. Fiard, T. A. Manteuffel, and S. F. McCormick

SIAM J. Sci. Comput. 19, pp. 1958-1979 (22 pages) | Cited 11 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
The focus of this paper is on planar linear convection-diffusion problems, to which we apply a special form of first-order system least squares (FOSLS [Cai et al., SIAM J. Numer. Anal., 31 (1994), pp. 1785--1799; Cai, Manteuffel, and McCormick, SIAM J. Numer. Anal., 34 (1997), pp. 425--454]). This we do by introducing the gradient of the primary variable, scaled by certain exponential functions. The convection-diffusion equation is then recast as a minimization principle for a functional corresponding to a sum of weighted L2 norms of the resulting first-order system. Discretization is accomplished by a Rayleigh--Ritz method based on standard finite element subspaces, and the resulting linear systems are solved by basic multigrid algorithms. The main goal here is to obtain optimal discretization accuracy and solver speed that is essentially uniform in the size of convection. Our results show that the FOSLS approach achieves this goal in general when the performance is measured either by the functional or by an equivalent weighted H1 norm. Included in our study is a multilevel adaptive refinement method based on locally uniform composite grids and local error estimates based on the functional itself.

A Wavelet-Optimized, Very High Order Adaptive Grid and Order Numerical Method

Leland Jameson

SIAM J. Sci. Comput. 19, pp. 1980-2013 (34 pages) | Cited 28 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
Wavelets detect information at different scales and at different locations throughout a computational domain. Furthermore, wavelets can detect the local polynomial content of computational data. Numerical methods are most efficient when the basis functions of the method are similar to the data present. By designing a numerical scheme in a completely adaptive manner around the data present in a computational domain, one can obtain optimal computational efficiency. This paper extends the numerical wavelet-optimized finite difference (WOFD) method to arbitrarily high order, so that one obtains, in effect, an adaptive grid and adaptive order numerical method which can achieve errors equivalent to errors obtained with a "spectrally accurate" numerical method.

Space-Time Continuous Analysis of Waveform Relaxation for the Heat Equation

Martin J. Gander and Andrew M. Stuart

SIAM J. Sci. Comput. 19, pp. 2014-2031 (18 pages) | Cited 14 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
Waveform relaxation algorithms for partial differential equations (PDEs) are traditionally obtained by discretizing the PDE in space and then splitting the discrete operator using matrix splittings. For the semidiscrete heat equation one can show linear convergence on unbounded time intervals and superlinear convergence on bounded time intervals by this approach. However, the bounds depend in general on the mesh parameter and convergence rates deteriorate as one refines the mesh.
Motivated by the original development of waveform relaxation in circuit simulation, where the circuits are split in the physical domain into subcircuits, we split the PDE by using overlapping domain decomposition. We prove linear convergence of the algorithm in the continuous case on an infinite time interval, at a rate depending on the size of the overlap. This result remains valid after discretization in space and the convergence rates are robust with respect to mesh refinement. The algorithm is in the class of waveform relaxation algorithms based on overlapping multisplittings. Our analysis quantifies the empirical observation by Jeltsch and Pohl [SIAM J. Sci. Comput., 16 (1995), pp. 40--49] that the convergence rate of a multisplitting algorithm depends on the overlap.
Numerical results are presented which support the convergence theory.

Asymptotic-Induced Domain Decomposition Methods for Kinetic and Drift Diffusion Semiconductor Equations

Axel Klar

SIAM J. Sci. Comput. 19, pp. 2032-2050 (19 pages) | Cited 8 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
This paper deals with domain decomposition methods for kinetic and drift diffusion semiconductor equations. In particular accurate coupling conditions at the interface between the kinetic and drift diffusion domain are given. The cases of slight and strong nonequilibrium situations at the interface are considered and numerical examples are shown.

A Solution-Based Triangular and Tetrahedral Mesh Quality Indicator

M. Berzins

SIAM J. Sci. Comput. 19, pp. 2051-2060 (10 pages) | Cited 12 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
A new mesh quality measure for triangular and tetrahedral meshes is presented. This mesh quality measure is based on both geometrical and solution information and is derived by considering the error when linear triangular and tetrahedral elements are used to approximate a quadratic function. The new measure is shown to be related to existing measures of mesh quality but with the advantage that local solution information in the form of scaled derivatives along edges is taken into account.

The Perfectly Matched Layer in Curvilinear Coordinates

Francis Collino and Peter Monk

SIAM J. Sci. Comput. 19, pp. 2061-2090 (30 pages) | Cited 42 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
In 1994 Bérenger showed how to construct a perfectly matched absorbing layer for the Maxwell system in rectilinear coordinates. This layer absorbs waves of any wavelength and any frequency without reflection and thus can be used to artificially terminate the domain of scattering calculations. In this paper we show how to derive and implement the Bérenger layer in curvilinear coordinates (in two space dimensions). We prove that an infinite layer of this type can be used to solve time harmonic scattering problems. We also show that the truncated Bérenger problem has a solution except at a discrete set of exceptional frequencies (which might be empty). Finally numerical results show that the curvilinear layer can produce accurate solutions in the time and frequency domain.

Geometric Mesh Partitioning: Implementation and Experiments

John R. Gilbert, Gary L. Miller, and Shang-Hua Teng

SIAM J. Sci. Comput. 19, pp. 2091-2110 (20 pages) | Cited 13 times

Online Publication Date: July 25, 2006

Full Text: | Download PDF

Show Abstract
We investigate a method of dividing an irregular mesh into equal-sized pieces with few interconnecting edges. The method's novel feature is that it exploits the geometric coordinates of the mesh vertices. It is based on theoretical work of Miller, Teng, Thurston, and Vavasis, who showed that certain classes of "well-shaped" finite-element meshes have good separators. The geometric method is quite simple to implement: we describe a \sc Matlab code for it in some detail. The method is also quite efficient and effective: we compare it with some other methods, including spectral bisection.
Close

close