Preconditioned iterative methods for a class of nonlinear eigenvalue problems

This paper proposes new iterative methods for the efficient computation of the smallest eigenvalue of symmetric nonlinear matrix eigenvalue problems of large order with a monotone dependence on the spectral parameter. Monotone nonlinear eigenvalue problems for differential equations have important applications in mechanics and physics. The discretization of these eigenvalue problems leads to nonlinear eigenvalue problems with very large sparse ill-conditioned matrices monotonically depending on the spectral parameter. To compute the smallest eigenvalue of large-scale matrix nonlinear eigenvalue problems, we suggest preconditioned iterative methods: preconditioned simple iteration method, preconditioned steepest descent method, and preconditioned conjugate gradient method. These methods use only matrix-vector multiplications, preconditioner-vector multiplications, linear operations with vectors, and inner products of vectors. We investigate the convergence and derive grid-independent error estimates for these methods. Numerical experiments demonstrate the practical effectiveness of the proposed methods for a model problem.     

The computation of bounds for the norm of the error in the conjugate gradient algorithm

In this paper we consider computing estimates of the norm of the error in the conjugate gradient (CG) algorithm. Formulas were given in a paper by Golub and Meurant (1997). Here, we first prove that these expressions are indeed upper and lower bounds for the A-norm of the error. Moreover, starting from these formulas, we investigate the computation of the l ₂-norm of the error. Finally, we define an adaptive algorithm where the approximations of the extreme eigenvalues that are needed to obtain upper bounds are computed when running CG leading to an improvement of the upper bounds for the norm of the error. Numerical experiments show the effectiveness of this algorithm. This revised version was published online in August 2006 with corrections to the Cover Date.     

On the quadrature error in operational quadrature methods for convolutions

Summary We analyze the quadrature error associated with operational quadrature methods for convolution equations. The assumptions are that the convolution kernel is inL ¹ and that its Laplace transform is analytic and bounded in an obtuse sector of the complex plane. Under these circumstances the Laplace transform has a slow variation property which admits a Fourier analysis of the quadrature error. We provide generalL ^p error estimates assuming suitable smoothness conditions on the function under convolution.     

Breakdowns and stagnation in iterative methods

One of the disadvantages of Krylov subspace iterative methods is the possibility of breakdown. This occurs when it is impossible to get the next approximation of the solution to the linear system of equationsAu=f. There are two different situations: lucky breakdown, when we have found the solution and hard breakdown, when the next Krylov subspace cannot be generated and/or the next approximate solution (iterate) cannot be computed. We show that some breakdowns depend on the chosen method of generating the basis vectors. Another undesirable feature of the iterative methods is stagnation. This occurs when the error does not change for several iterative steps. We investigate when iterative methods can stagnate and describe conditions which characterize stagnation. We show that in some cases stagnation can imply breakdown.     

A basic norm equivalence for the theory of multilevel methods

Summary Subspace decompositions of finite element spaces based onL ₂-like orthogonal projections play an important role for the construction and analysis of multigrid like iterative methods. Recently several authors have proved the equivalence of the associated discrete norms with theH ¹-norm. The present paper gives an elementary, self-contained derivation of this result which is based on the use ofK-functionals known from the theory of interpolation spaces.     

Schwarz alternating and iterative refinement methods for mixed formulations of elliptic problems,part II: Convergence theory

Summary In this paper we discuss bounds for the convergence rates of several domain decomposition algorithms to solve symmetric, indefinite linear systems arising from mixed finite element discretizations of elliptic problems. The algorithms include Schwarz methods and iterative refinement methods on locally refined grids. The implementation of Schwarz and iterative refinement algorithms have been discussed in part I. A discussion on the stability of mixed discretizations on locally refined grids is included and quantiative estimates for the convergence rates of some iterative refinement algorithms are also derived.Department of Mathematics, University of Wyoming, Laramie, WY 82071-3036. This work was supported in part by the National Science Foundation under Grant NSF-CCR-8903003, while the author was a graduate student at New York University, and in part by NSF Grant ASC 9003002, while the author was a Visiting, Assistant Researcher at UCLA.     

Estimating the largest singular values of large sparse matrices via modified moments

We describe a procedure for determining a few of the largest singular values of a large sparse matrix. The method by Golub and Kent which uses the method of modified moments for estimating the eigenvalues of operators used in iterative methods for the solution of linear systems of equations is appropriately modified in order to generate a sequence of bidiagonal matrices whose singular values approximate those of the original sparse matrix. A simple Lanczos recursion is proposed for determining the corresponding left and right singular vectors. The potential asynchronous computation of the bidiagonal matrices using modified moments with the iterations of an adapted Chebyshev semi-iterative (CSI) method is an attractive feature for parallel computers. Comparisons in efficiency and accuracy with an appropriate Lanczos algorithm (with selective re-orthogonalization) are presented on large sparse (rectangular) matrices arising from applications such as information retrieval and seismic reflection tomography. This procedure is essentially motivated by the theory of moments and Gauss quadrature.This author's work was supported by the National Science Foundation under grants NSF CCR-8717492 and CCR-910000N (NCSA), the U.S. Department of Energy under grant DOE DE-FG02-85ER25001, and the Air Force Office of Scientific Research under grant AFOSR-90-0044 while at the University of Illinois at Urbana-Champaign Center for Supercomputing Research and Development.This author's work was supported by the U.S. Army Research Office under grant DAAL03-90-G-0105, and the National Science Foundation under grant NSF DCR-8412314.     

AIR Tools — A MATLAB package of algebraic iterative reconstruction methods

We present a MATLAB package with implementations of several algebraic iterative reconstruction methods for discretizations of inverse problems. These so-called row action methods rely on semi-convergence for achieving the necessary regularization of the problem. Two classes of methods are implemented: Algebraic Reconstruction Techniques (ART) and Simultaneous Iterative Reconstruction Techniques (SIRT). In addition we provide a few simplified test problems from medical and seismic tomography. For each iterative method, a number of strategies are available for choosing the relaxation parameter and the stopping rule. The relaxation parameter can be fixed, or chosen adaptively in each iteration; in the former case we provide a new “training” algorithm that finds the optimal parameter for a given test problem. The stopping rules provided are the discrepancy principle, the monotone error rule, and the NCP criterion; for the first two methods “training” can be used to find the optimal discrepancy parameter.     

Localization of the Aronszajn-Slobodeckij norm and application to adaptive boundary element methods Part II. The three-dimensional case

Summary. In this paper we introduce new local a-posteriori error indicators for the Galerkin discretization of three-dimensional boundary integral equations. These error indicators are efficient and reliable for a wide class of integral operators, in particular for operators of negative order. They are based on local norms of the computable residual and can be used for controlling the adaptive refinement. The proofs of efficiency and reliability are based on the result that the Aronszajn-Slobodeckij norm (given by a double integral for a non-integer ) is localizable for certain functions. Neither inverse estimates nor saturation properties are needed. In this paper, we extend the two-dimensional results of a previous paper to the three-dimensional case. Received March 20, 2000 / Published online November 15, 2001     

Rounding error analysis of the classical Gram-Schmidt orthogonalization process

This paper provides two results on the numerical behavior of the classical Gram-Schmidt algorithm. The first result states that, provided the normal equations associated with the initial vectors are numerically nonsingular, the loss of orthogonality of the vectors computed by the classical Gram-Schmidt algorithm depends quadratically on the condition number of the initial vectors. The second result states that, provided the initial set of vectors has numerical full rank, two iterations of the classical Gram-Schmidt algorithm are enough for ensuring the orthogonality of the computed vectors to be close to the unit roundoff level.The work of the second author was supported in part by the US Department of Energy, Office of Basic Energy Science under LAB03-17 initiative, DOE contract No. DE-FG02-03ER25584, and in part by the TeraScale Optimal PDE Simulations (TOPS) SciDAC, DoE Contract No. DE-FC02-01ER25480The work of the third author was supported by the project 1ET400300415 within the National Program of Research ‘‘Information Society’’ and by the GA AS CR under grant No. IAA1030405.     

A sparse counterpart of Reichel and Gragg’s package QRUP

We describe how to maintain the triangular factor of a sparse QR factorization when columns are added and deleted and Q cannot be stored for sparsity reasons. The updating procedures could be thought of as a sparse counterpart of Reichel and Gragg’s package QRUP. They allow us to solve a sequence of sparse linear least squares subproblems in which each matrix B_k is an independent subset of the columns of a fixed matrix A, and B_k+1 is obtained by adding or deleting one column. Like Coleman and Hulbert [T. Coleman, L. Hulbert, A direct active set algorithm for large sparse quadratic programs with simple bounds, Math. Program. 45 (1989) 373-406], we adapt the sparse direct methodology of Björck and Oreborn of the late 1980s, but without forming A^TA, which may be only positive semidefinite. Our Matlab 5 implementation works with a suitable row and column numbering within a static triangular sparsity pattern that is computed in advance by symbolic factorization of A^TA and preserved with placeholders.     

Schwarz alternating and iterative refinement methods for mixed formulations of elliptic problems,part I: Algorithms and numerical results

Summary We describe sequential and parallel algorithms based on the Schwarz alternating method for the solution of mixed finite element discretizations of elliptic problems using the Raviart-Thomas finite element spaces. These lead to symmetric indefinite linear systems and the algorithms have some similarities with the traditional block Gauss-Seidel or block Jacobi methods with overlapping blocks. The indefiniteness requires special treatment. The sub-blocks used in the algorithm correspond to problems on a coarse grid and some overlapping subdomains and is based on a similar partition used in an algorithm of Dryja and Widlund for standard elliptic problems. If there is sufficient overlap between the subdomains, the algorithm converges with a rate independent of the mesh size, the number of subdomains and discontinuities of the coefficients. Extensions of the above algorithms to the case of local grid refinement is also described. Convergence theory for these algorithms will be presented in a subsequent paper.This work was supported in part by the National Science Foundation under Grant NSF-CCR-8903003, while the author was a graduate student at New York University, and in part by the Army Research Office under Grant DAAL 03-91-G-0150, while the author was a Visiting Assistant Researcher at UCLA     

A priori error estimation for the dual mixed finite element method of the elastodynamic problem in a polygonal domain, II

In this paper we analyze a new dual mixed formulation of the elastodynamic system in polygonal domains by using an implicit scheme for the time discretization. After the analysis of stability of the fully discrete scheme, L^∞ in time, L² in space a priori error estimates for the approximation of the displacement, the strain, the pressure and the rotational are derived. Numerical tests are presented which confirm our theoretical results.     

Multiscale asymptotic expansion and finite element methods for the mixed boundary value problems of second order elliptic equation in perforated domains

In this paper, we will discuss the mixed boundary value problems for the second order elliptic equation with rapidly oscillating coefficients in perforated domains, and will present the higher-order multiscale asymptotic expansion of the solution for the problem, which will play an important role in the numerical computation . The convergence theorems and their rigorous proofs will be given. Finally a multiscale finite element method and some numerical results will be presented. This work is Supported by National Natural Science Foundation of China (grant # 10372108, # 90405016), and Special Funds for Major State Basic Research Projects( grant # TG2000067102)     

Existence and uniqueness of splittings for stationary iterative methods with applications to alternating methods

Summary. Given a nonsingular matrix , and a matrix of the same order, under certain very mild conditions, there is a unique splitting , such that . Moreover, all properties of the splitting are derived directly from the iteration matrix . These results do not hold when the matrix is singular. In this case, given a matrix and a splitting such that , there are infinitely many other splittings corresponding to the same matrices and , and different splittings can have different properties. For instance, when is nonnegative, some of these splittings can be regular splittings, while others can be only weak splittings. Analogous results hold in the symmetric positive semidefinite case. Given a singular matrix , not for all iteration matrices there is a splitting corresponding to them. Necessary and sufficient conditions for the existence of such splittings are examined. As an illustration of the theory developed, the convergence of certain alternating iterations is analyzed. Different cases where the matrix is monotone, singular, and positive (semi)definite are studied. Received September 5, 1995 / Revised version received April 3, 1996     

Convergence of block iterative methods for linear systems arising in the numerical solution of Euler equations

Summary We discuss block matrices of the formA=[A _ij ], whereA _ij is ak×k symmetric matrix,A _ij is positive definite andA _ij is negative semidefinite. These matrices are natural block-generalizations of Z-matrices and M-matrices. Matrices of this type arise in the numerical solution of Euler equations in fluid flow computations. We discuss properties of these matrices, in particular we prove convergence of block iterative methods for linear systems with such system matrices.     

Running error for the evaluation of rational Bézier surfaces

Error analysis of the usual method to evaluate rational Bézier surfaces is performed. The corresponding running error analysis is also carried out and the sharpness of our running error bounds is shown. We also modify the evaluation algorithm to include such error bounds without increasing significantly its computational cost.     

Convergence rates of some iterative methods for nonsymmetric algebraic Riccati equations arising in transport theory

We determine and compare the convergence rates of various fixed-point iterations for finding the minimal positive solution of a class of nonsymmetric algebraic Riccati equations arising in transport theory.     

Extensions of block-projections methods with relaxation parameters to inconsistent and rank-deficient least-squares problems

There exist many classes of block-projections algorithms for approximating solutions of linear least-squares problems. Generally, these methods generate sequences convergent to the minimal norm least-squares solution only for consistent problems. In the inconsistent case, which usually appears in practice because of some approximations or measurements, these sequences do no longer converge to a least-squares solution or they converge to the minimal norm solution of a “perturbed” problem. In the present paper, we overcome this difficulty by constructing extensions for almost all the above classes of block-projections methods. We prove that the sequences generated with these extensions always converge to a least-squares solution and, with a suitable initial approximation, to the minimal norm solution of the problem. Numerical experiments, described in the last section of the paper, confirm the theoretical results obtained.