Search

In the present paper, we describe an adaptive modified rational global Lanczos algorithm for model-order reduction problems using multipoint moment matching-based methods. The major problem of these methods is the selection of some interpolation points. We first propose a modified rational global Lanczos process and then we derive Lanczos-like equations for the global case. Next, we propose adaptive techniques for choosing the interpolation points. Second-order dynamical systems are also considered in this paper, and the adaptive modified rational global Lanczos algorithm is applied to an equivalent state space model. Finally, some numerical examples will be given. [ABSTRACT FROM AUTHOR]

Based on the spectral divide‐and‐conquer algorithm by Nakatsukasa and Higham [SIAM J. Sci. Comput., 35(3):A1325–A1349, 2013], we propose a new algorithm for computing all the eigenvalues and eigenvectors of a symmetric banded matrix with small bandwidth, with the eigenvectors given implicitly as a product of orthonormal matrices stored in the so‐called hierarchically off‐diagonal low‐rank (HODLR) format. For this purpose, we combine our previous work on the fast computation of spectral projectors in the HODLR format, with a novel technique for extracting a basis for the range of such a HODLR matrix. Preliminary numerical experiments demonstrate that our algorithm exhibits quasi‐linear complexity for matrices that can be efficiently represented in the HODLR format throughout the divide‐and‐conquer algorithm, and allows for conveniently dealing with such large‐scale matrices. [ABSTRACT FROM AUTHOR]

In this paper, two accelerated divide-and-conquer (ADC) algorithms are proposed for the symmetric tridiagonal eigenvalue problem, which cost O(N2r) flops in the worst case, where N is the dimension of the matrix and r is a modest number depending on the distribution of eigenvalues. Both of these algorithms use hierarchically semiseparable (HSS) matrices to approximate some intermediate eigenvector matrices, which are Cauchy-like matrices and are off-diagonally low-rank. The difference of these two versions lies in using different HSS construction algorithms, one (denoted by ADC1) uses a structured low-rank approximation method and the other (ADC2) uses a randomized HSS construction algorithm. For the ADC2 algorithm, a method is proposed to estimate the off-diagonal rank. Numerous experiments have been carried out to show their stability and efficiency. These algorithms are implemented in parallel in a shared memory environment, and some parallel implementation details are included. Comparing the ADCs with highly optimized multithreaded libraries such as Intel MKL, we find that ADCs could be more than six times faster for some large matrices with few deflations. [ABSTRACT FROM AUTHOR]

Summary: Following the Perron theorem, the spectral radius of a primitive matrix is a simple eigenvalue. It is shown that for a primitive matrix A, there is a positive rank one matrix X such that B = A ∘ X, where ∘ denotes the Hadamard product of matrices, and such that the row (column) sums of matrix B are the same and equal to the Perron root. An iterative algorithm is presented to obtain matrix B without an explicit knowledge of X. The convergence rate of this algorithm is similar to that of the power method but it uses less computational load. A byproduct of the proposed algorithm is a new method for calculating the first eigenvector. [ABSTRACT FROM AUTHOR]

Summary: We present the augmented matrix for principal submatrix update (AMPS) algorithm, a finite element solution method that combines principal submatrix updates and Schur complement techniques, well‐suited for interactive simulations of deformation and cutting of finite element meshes. Our approach features real‐time solutions to the updated stiffness matrix systems to account for interactive changes in mesh connectivity and boundary conditions. Updates are accomplished by an augmented matrix formulation of the stiffness equations to maintain its consistency with changes to the underlying model without refactorization at each timestep. As changes accumulate over multiple simulation timesteps, the augmented solution algorithm enables tens or hundreds of updates per second. Acceleration schemes that exploit sparsity, memoization and parallelization lead to the updates being computed in real time. The complexity analysis and experimental results for this method demonstrate that it scales linearly with the number of nonzeros of the factors of the stiffness matrix. Results for cutting and deformation of three‐dimensional (3D) elastic models are reported for meshes with up to 50 000 nodes, and involve models of surgery for astigmatism and the brain. [ABSTRACT FROM AUTHOR]

This paper is concerned with the tight enclosure of matrix multiplication AB for two floating-point matrices A and B. The aim of this paper is to compute component-wise upper and lower bounds of the exact result C of the matrix multiplication AB by floating-point arithmetic. Namely, an interval matrix enclosing C is obtained. In this paper, new algorithms for enclosing C are proposed. The proposed algorithms are designed to mainly exploit the level 3 operations in BLAS. Although the proposed algorithms take around twice as much costs as a standard algorithm promoted by Oishi and Rump, the accuracy of the result by the proposed algorithms is better than that of the standard algorithm. At the end of this paper, we present numerical examples showing the efficiency of the proposed algorithms. Copyright © 2010 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

SUMMARY In this paper, a rapid iterative algorithm is proposed to find robust approximations for the inverse of nonsingular matrices. The analysis of convergence reveals that this high-order method possesses eighth-order convergence. The interesting point is that, this rate is attained using less number of matrix-by-matrix multiplications in contrast to the existing methods of the same type in the literature. The extension of the method for finding Moore-Penrose inverse of singular or rectangular matrices is also presented. Numerical comparisons will be given to show the applicability, stability and consistency of the new scheme by paying special attention on the computational time. Copyright © 2013 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

In this paper, we study a class of weakly nonlinear complementarity problems arising from the discretization of free boundary problems. By reformulating the complementarity problems as implicit fixed-point equations based on splitting of the system matrices, we propose a class of modulus-based matrix splitting algorithms. We show their convergence by assuming that the system matrix is positive definite. Moreover, we give several kinds of typical practical choices of the modulus-based matrix splitting iteration methods based on the different splitting of the system matrix. Numerical experiments on two model problems are presented to illustrate the theoretical results and examine the numerical effectiveness of our modulus-based matrix splitting algorithms. Copyright © 2016 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

Among numerous iterative methods for solving the minimal nonnegative solution of an M-matrix algebraic Riccati equation, the structure-preserving doubling algorithm (SDA) stands out owing to its overall efficiency as well as accuracy. SDA is globally convergent and its convergence is quadratic, except for the critical case for which it converges linearly with the linear rate 1/2. In this paper, we first undertake a delineatory convergence analysis that reveals that the approximations by SDA can be decomposed into two components: the stable component that converges quadratically and the rank-one component that converges linearly with the linear rate 1/2. Our analysis also shows that as soon as the stable component is fully converged, the rank-one component can be accurately recovered. We then propose an efficient hybrid method, called the two-phase SDA, for which the SDA iteration is stopped as soon as it is determined that the stable component is fully converged. Therefore, this two-phase SDA saves those SDA iterative steps that previously have to have for the rank-one component to be computed accurately, and thus essentially, it can be regarded as a quadratically convergent method. Numerical results confirm our analysis and demonstrate the efficiency of the new two-phase SDA. Copyright © 2015 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

SUMMARY Novel memory-efficient Arnoldi algorithms for solving matrix polynomial eigenvalue problems are presented. More specifically, we consider the case of matrix polynomials expressed in the Chebyshev basis, which is often numerically more appropriate than the standard monomial basis for a larger degree d. The standard way of solving polynomial eigenvalue problems proceeds by linearization, which increases the problem size by a factor d. Consequently, the memory requirements of Krylov subspace methods applied to the linearization grow by this factor. In this paper, we develop two variants of the Arnoldi method that build the Krylov subspace basis implicitly, in a way that only vectors of length equal to the size of the original problem need to be stored. The proposed variants are generalizations of the so-called quadratic Arnoldi method and two-level orthogonal Arnoldi procedure methods, which have been developed for the monomial case. We also show how the typical ingredients of a full implementation of the Arnoldi method, including shift-and-invert and restarting, can be incorporated. Numerical experiments are presented for matrix polynomials up to degree 30 arising from the interpolation of nonlinear eigenvalue problems, which stem from boundary element discretizations of PDE eigenvalue problems. Copyright © 2013 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

SUMMARY The aim of the nuclear norm minimization problem is to find a matrix that minimizes the sum of its singular values and satisfies some constraints simultaneously. Such a problem has received more attention largely because it is closely related to the affine rank minimization problem, which appears in many control applications including controller design, realization theory, and model reduction. In this paper, we first propose an exact version alternating direction method for solving the nuclear norm minimization problem with linear equality constraints. At each iteration, the method involves a singular value thresholding and linear matrix equations which are solved exactly. Convergence of the proposed algorithm is followed directly. To broaden the capacity of solving larger problems, we solve approximately the subproblem by an iterative method with the Barzilai-Borwein steplength. Some extensions to the noisy problems and nuclear norm regularized least-square problems are also discussed. Numerical experiments and comparisons with the state-of-the-art method FPCA show that the proposed method is effective and promising. Copyright © 2011 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

SUMMARY In this paper, we revisit the steady-state solution method for Markov Regenerative Processes (MRP) proposed in the work by German. This method solves the embedded Markov chain P of the MRP without storing the matrix P explicitly. We address three issues left open in German's Work: 1) the solution method is restricted to Power method; 2) it has been defined only for ergodic MRPs; and 3) no preconditioning is available to speed-up the computation. This paper discusses how to lift these limitations by extending the algorithm to preconditioned Krylov-subspace methods and by generalizing it to the non-ergodic case. An MRP-specific preconditioner is also proposed, which is built from a sparse approximation of the MRP matrix, computed via simulation. An experimental assessment of the proposed preconditioner is then provided. Copyright © 2011 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

This paper is motivated by some recent work of Fukuda, Ishiwata, Iwasaki, and Nakamura ( Inverse Problems 2009; :015007). We first design an algorithm for computing the eigenvalues of a specially structured matrix from the discrete Bogoyavlensky Lattice 2 (dBL2) system. A Lax representation for the dBL2 system is given in a matrix form. By considering the asymptotic behavior of dBL2 variables, some characteristic polynomials are then factorized. A new algorithm for computing the complex eigenvalues of a specially structured matrix is then introduced. Copyright © 2010 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

This paper concerns the LBM factorization of unsymmetric tridiagonal matrices, where L and M are unit lower triangular matrices and B is block diagonal with 1×1 and 2×2 blocks. In some applications, it is necessary to form this factorization without row or column interchanges while the tridiagonal matrix is formed. Bunch and Kaufman proposed a pivoting strategy without interchanges specifically for symmetric tridiagonal matrices, and more recently, Bunch and Marcia proposed pivoting strategies that are normwise backward stable for linear systems involving such matrices. In this paper, we extend these strategies to the unsymmetric tridiagonal case and demonstrate that the proposed methods both exhibit bounded growth factors and are normwise backward stable. Copyright © 2009 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

In this paper we consider an inverse problem for a damped vibration system from the noisy measured eigendata, where the mass, damping, and stiffness matrices are all symmetric positive-definite matrices with the mass matrix being diagonal and the damping and stiffness matrices being tridiagonal. To take into consideration the noise in the data, the problem is formulated as a convex optimization problem involving quadratic constraints on the unknown mass, damping, and stiffness parameters. Then we propose a smoothing Newton-type algorithm for the optimization problem, which improves a pre-existing estimate of a solution to the inverse problem. We show that the proposed method converges both globally and quadratically. Numerical examples are also given to demonstrate the efficiency of our method. Copyright © 2008 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

This paper extends the weighted low rank approximation (WLRA) approach to linearly structured matrices. In the case of Hankel matrices with a special block structure, an equivalent unconstrained optimization problem is derived and an algorithm for solving it is proposed. Copyright © 2005 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

The past 13 years have seen the development of many algorithms for approximating matrix functions in O(N) time, where N is the basis size. These O(N) algorithms rely on assumptions about the spatial locality of the matrix function; therefore their validity depends very much on the argument of the matrix function. In this article I carefully examine the validity of certain O(N) algorithms when applied to Hamiltonians of disordered systems. I focus on the prototypical disordered system, the Anderson model. I find that O(N) algorithms for the density matrix function can be used well below the Anderson transition (i.e. in the metallic phase;) they fail only when the coherence length becomes large. This paper also includes some experimental results about the Anderson model's behaviour across a range of disorders. Copyright © 2005 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

Many popular solution methods for large discrete ill-posed problems are based on Tikhonov regularization and compute a partial Lanczos bidiagonalization of the matrix. The computational effort required by these methods is not reduced significantly when the matrix of the discrete ill-posed problem, rather than being a general nonsymmetric matrix, is symmetric and possibly indefinite. This paper describes new methods, based on partial Lanczos tridiagonalization of the matrix, that exploit symmetry. Computed examples illustrate that one of these methods can require significantly less computational work than available structure-ignoring schemes. Copyright © 2004 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

In this paper, we develop a fast regularized structured total least squares (RSTLS) algorithm for solving the basic deconvolution problem. The algorithm is based on a particular implementation of the generalized Schur algorithm. We apply the new algorithm on a deconvolution problem arising in a medical application in renography. By means of this example, we show that regularization in the structured total least squares framework yields more accurate results compared to other estimators. Copyright © 2004 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

This paper proposes and studies the performance of a preconditioner suitable for solving a class of symmetric positive definite systems, Apx=b, which we call lower rank extracted systems (LRES), by the preconditioned conjugate gradient method. These systems correspond to integral equations with convolution kernels defined on a union of many line segments in contrast to only one line segment in the case of Toeplitz systems. The p × p matrix, Ap, is shown to be a principal submatrix of a larger N × N Toeplitz matrix, AN. The preconditioner is provided in terms of the inverse of a 2N × 2N circulant matrix constructed from the elements of AN. The preconditioner is shown to yield clustering in the spectrum of the preconditioned matrix similar to the clustering results for iterative algorithms used to solve Toeplitz systems. The analysis also demonstrates that the computational expense to solve LRE systems is reduced to O(N log N). Copyright © 2004 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

This paper analyses a novel method for constructing preconditioners for diagonally dominant symmetric positive-definite matrices. The method discussed here is based on a simple idea: we construct M by simply dropping offdiagonal non-zeros from A and modifying the diagonal elements to maintain a certain row-sum property. The preconditioners are extensions of Vaidya's augmented maximum-spanning-tree preconditioners. The preconditioners presented here were also mentioned by Vaidya in an unpublished manuscript, but without a complete analysis. The preconditioners that we present have only O(n+t2) nonzeros, where n is the dimension of the matrix and 1⩽t⩽n is a parameter that one can choose. Their construction is efficient and guarantees that the condition number of the preconditioned system is O(n2/t2) if the number of nonzeros per row in the matrix is bounded by a constant. We have developed an efficient algorithm to construct these preconditioners and we have implemented it. We used our implementation to solve a simple model problem; we show the combinatorial structure of the preconditioners and we present encouraging convergence results. [ABSTRACT FROM AUTHOR]

Additive Schwarz preconditioned GMRES is a powerful method for solving large sparse linear systems of equations on parallel computers. The algorithm is often implemented in the Euclidean norm, or the discrete l2 norm, however, the optimal convergence result is available only in the energy norm (or the equivalent Sobolev H1 norm). Very little progress has been made in the theoretical understanding of the l2 behaviour of this very successful algorithm. To add to the difficulty in developing a full l2 theory, in this note, we construct explicit examples and show that the optimal convergence of additive Schwarz preconditioned GMRES in l2 cannot be obtained using the existing GMRES theory. More precisely speaking, we show that the symmetric part of the preconditioned matrix, which plays a role in the Eisenstat–Elman–Schultz theory, has at least one negative eigenvalue, and we show that the condition number of the best possible eigenmatrix that diagonalizes the preconditioned matrix, key to the Saad–Schultz theory, is bounded from both above and below by constants multiplied by h-1/2. Here h is the finite element mesh size. The results presented in this paper are mostly negative, but we believe that the techniques used in our proofs may have wide applications in the further development of the l2 convergence theory and in other areas of domain decomposition methods. Copyright © 2002 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

Recently the authors proposed a simultaneous iteration algorithm for the computation of the partial derivatives of repeated eigenvalues and the corresponding eigenvectors of matrices depending on several real variables. This paper analyses the properties of that algorithm and extends it in several ways. The previous requirement that the repeated eigenvalue be dominant is relaxed, and the new generalized algorithm given here allows the simultaneous treatment of simple and repeated eigenvalues. Methods for accelerating convergence are examined. Numerical results support our theoretical analysis. Copyright © 2000 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

The purpose of this paper is to present an algorithm for matrix multiplication based on a formula discovered by Pan [7]. For matrices of order up to 10 000, the nearly optimum tuning of the algorithm results in a rather clear non-recursive one- or two-level structure with the operation count comparable to that of the Strassen algorithm [9]. The algorithm takes less workspace and has better numerical stability as compared to the Strassen algorithm, especially in Winograd's modification [2]. Moreover, its clearer and more flexible structure is potentially more suitable for efficient implementation on modern supercomputers. Copyright © 1999 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

SUMMARY The main purpose of this paper is the extension of the classical spectral direct and inverse analysis of Jacobi matrices for the non-self-adjoint setting. Matrices of this class appear in the context of non-Hermitian quantum mechanics. The reconstruction of a pseudo-Jacobi matrix from its spectrum and the spectra of two complementary principal matrices is investigated in the context of indefinite inner product spaces. An existence and uniqueness theorem is given, and a strikingly simple algorithm, based on the Euclidean division algorithm, to reconstruct the matrix from the spectral data is presented. A result of Friedland and Melkman stating a necessary and sufficient condition for a real sequence to be the spectrum of a non-negative Jacobi matrix is revisited and generalized. Namely, it is shown that a suitable set of prescribed eigenvalues defines a unique non-negative pseudo-Jacobi matrix, which is J-Hermitian for a fixed J. Copyright © 2012 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

We study Lanczos and polynomial algorithms with random start for estimating an eigenvector corresponding to the largest eigenvalue of an n × n large symmetric positive definite matrix. We analyze the two error criteria: the randomized error and the randomized residual error. For the randomized error, we prove that it is not possible to get distribution-free bounds, i.e., the bounds must depend on the distribution of eigenvalues of the matrix. We supply such bounds and show that they depend on the ratio of the two largest eigenvalues. For the randomized residual error, distribution-free bounds exist and are provided in the paper. We also provide asymptotic bounds, as well as bounds depending on the ratio of the two largest eigenvalues. The bounds for the Lanczos algorithm may be helpful in a practical implementation and termination of this algorithm. © 1998 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

A direct algorithm for the solution to the affine two-sided obstacle problem with an M-matrix is presented. The algorithm has the polynomial bounded computational complexity O( n) and is more efficient than those in ( Numer. Linear Algebra Appl. 2006; :543-551). Copyright © 2010 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

In this paper, two new matrix-form iterative methods are presented to solve the least-squares problem:\documentclass{article}\usepackage[mathscr]{euscript}\usepackage{amssymb,amsbsy,rotating}\footskip=0pc\pagestyle{empty}\begin{document}\[ \min\limits_{X\in {\mathscr{S}}_1,Y\in {\mathscr{S}}_2}\|AXB+CYD-E\| \]\end{document} and matrix nearness problem: \documentclass{article}\footskip=0pc\pagestyle{empty}\begin{document}\[ \min\limits_{[X: Y]\in S_{XY}}\|[X: Y]-[\widetilde{X}:\widetilde{Y}]\| \] \end{document} where matrices $A\in R^{p\times n_1},B\in R^{n_2\times q},C\in R^{p\times m_1},D\in R^{m_2\times q},E\in R^{p\times q},\widetilde{X}\in R^{n_1\times n_2}$ and $\widetilde{Y}\in R^{m_1\times m_2}$ are given; 𝒮1 and 𝒮2 are the set of constraint matrices, such as symmetric, skew symmetric, bisymmetric and centrosymmetric matrices sets and SXY is the solution pair set of the minimum residual problem. These new matrix-form iterative methods have also faster convergence rate and higher accuracy than the matrix-form iterative methods proposed by Peng and Peng (Numer. Linear Algebra Appl. 2006; 13: 473–485) for solving the linear matrix equation AXB+CYD=E. Paige's algorithms, which are based on the bidiagonalization procedure of Golub and Kahan, are used as the framework for deriving these new matrix-form iterative methods. Some numerical examples illustrate the efficiency of the new matrix-form iterative methods. Copyright © 2008 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

The paper presents fast algorithms for designing finite impulse response (FIR) notch filters. The aim is to design a digital FIR notch filter so that the magnitude of the filter has a deep notch at a specified frequency, and as the notch frequency changes, the filter coefficients should be able to track the notch fast in real time. The filter design problem is first converted into a convex optimization problem in the autocorrelation domain. The frequency response of the autocorrelation of the filter impulse response is compared with the desired filter response and the integral square error is minimized with respect to the unknown autocorrelation coefficients. Spectral factorization is used to calculate the coefficients of the filter. In the optimization process, the computational advantage is obtained by exploiting the structure of the Hessian matrix which consists of a Toeplitz plus a Hankel matrix. Two methods have been used for solving the Toeplitz-plus-Hankel system of equations. In the first method, the computational time is reduced by using Block–Levinson's recursion for solving the Toeplitz-plus-Hankel system of matrices. In the second method, the conjugate gradient method with different preconditioners is used to solve the system. Comparative studies demonstrate the computational advantages of the latter. Both these algorithms have been used to obtain the autocorrelation coefficients of notch filters with different orders. The original filter coefficients are found by spectral factorization and each of these filters have been tested for filtering synthetic as well as real-life signals. Copyright © 2007 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

In this paper, we consider the approximation of dense block Toeplitz-plus-Hankel matrices by sums of Kronecker products of Toeplitz-plus-Hankel matrices. We present an algorithm for efficiently computing the matrix approximation that requires the factorization of matrices of much smaller dimension than that of the original. The main results are described for block Toeplitz matrices with Toeplitz-plus-Hankel blocks, but the algorithms can be readily adjusted for other related structures that arise in image processing applications, such as block Toeplitz with Toeplitz blocks and block Toeplitz-plus-Hankel with Toeplitz-plus-Hankel blocks. Our work extends the techniques in Kamm and Nagy (SIAM J. Matrix Anal. Appl. 2000; 22:155–172) and Nagy et al. (SIAM J. Matrix Anal. Appl. 2004; 25:829–841), which consider similar matrices, but with the added restriction that the matrices have a banded/block-banded structure. We illustrate the effectiveness of our algorithm by using the output of the algorithm to construct preconditioners for systems from two different applications: diffuse optical tomography and atmospheric image deblurring. Copyright © 2007 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

This paper presents a new algebraic extension of the Rayleigh quotient multigrid (RQMG) minimization algorithm to compute the smallest eigenpairs of a symmetric positive definite pencil (A, M). Earlier versions of RQMG minimize the Rayleigh quotient over a hierarchy of geometric grids. We replace the geometric mesh information with the algebraic information defined by an algebraic multigrid preconditioner. At each level, we minimize the Rayleigh quotient with a block preconditioned algorithm. Numerical experiments illustrate the efficiency of this new algorithm to compute several eigenpairs. Copyright © 2007 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

In this paper, we describe and analyse several block matrix iterative algorithms for solving a saddle point linear system arising from the discretization of a linear-quadratic elliptic control problem with Neumann boundary conditions. To ensure that the problem is well posed, a regularization term with a parameter α is included. The first algorithm reduces the saddle point system to a symmetric positive definite Schur complement system for the control variable and employs conjugate gradient (CG) acceleration, however, double iteration is required (except in special cases). A preconditioner yielding a rate of convergence independent of the mesh size h is described for Ω ⊂ R2 or R3, and a preconditioner independent of h and α when Ω ⊂ R2. Next, two algorithms avoiding double iteration are described using an augmented Lagrangian formulation. One of these algorithms solves the augmented saddle point system employing MINRES acceleration, while the other solves a symmetric positive definite reformulation of the augmented saddle point system employing CG acceleration. For both algorithms, a symmetric positive definite preconditioner is described yielding a rate of convergence independent of h. In addition to the above algorithms, two heuristic algorithms are described, one a projected CG algorithm, and the other an indefinite block matrix preconditioner employing GMRES acceleration. Rigorous convergence results, however, are not known for the heuristic algorithms. Copyright © 2007 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

In this paper, two direct algorithms for solving the two-sided obstacle problem with an M-matrix are presented. The algorithms are well defined and have polynomial computational complexity. Copyright © 2006 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

The anti-reflective boundary condition for image restoration was recently introduced as a mathematically desirable alternative to other boundary conditions presently represented in the literature. It has been shown that, given a centrally symmetric point spread function (PSF), this boundary condition gives rise to a structured blurring matrix, a submatrix of which can be diagonalized by the discrete sine transform (DST), leading to an O(n2 log n) solution algorithm for an image of size n × n. In this paper, we obtain a Kronecker product approximation of the general structured blurring matrix that arises under this boundary condition, regardless of symmetry properties of the PSF. We then demonstrate the usefulness and efficiency of our approximation in an SVD-based restoration algorithm, the computational cost of which would otherwise be prohibitive. Copyright © 2005 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

The solution of large sparse linear systems of equations is one of the cornerstones of scientific computation. In many applications it is important to be able to solve these systems as rapidly as possible. One approach for very large problems is to reorder the system matrix to bordered block diagonal form and then to solve the block system using a coarse-grained parallel approach. In this paper, we consider the problem of efficiently ordering unsymmetric systems to singly bordered block diagonal form. Algorithms such as the MONET algorithm of Hu et al. (Comput. Chem. Eng. 23 (2000) 1631) that depend upon computing a representation of AAT can be prohibitively expensive when a single (or small number of) matrix factorizations are required. We therefore work with the graph of AT + A (or BT + B, where B is a row permutation of A) and propose new reordering algorithms that use only vertex separators and wide separators of this graph. Numerical experiments demonstrate that our methods are efficient and can produce bordered forms that are competitive with those obtained using MONET. Copyright © 2005 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

The LBLT factorization of Bunch for solving linear systems involving a symmetric indefinite tridiagonal matrix T is a stable, efficient method. It computes a unit lower triangular matrix L and a block 1 × 1 and 2 × 2 matrix B such that T=LBLT. Choosing the pivot size requires knowing a priori the largest element σ of T in magnitude. In some applications, it is required to factor T as it is formed without necessarily knowing σ. In this paper, we present a modification of the Bunch algorithm that can satisfy this requirement. We demonstrate that this modification exhibits the same bound on the growth factor as the Bunch algorithm and is likewise normwise backward stable. Copyright © 2005 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

Structured matrix algebras ℒ and a generalized BFGS-type iterative scheme have been recently investigated to introduce low-complexity quasi-Newton methods, named ℒQN, for solving general (non-structured) minimization problems. In this paper we introduce the ℒkQN methods, which exploit ad hoc algebras at each step. Since the structure of the updated matrices can be modified at each iteration, the new methods can better fit the Hessian matrix, thereby improving the rate of convergence of the algorithm. Copyright © 2005 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

The QR algorithm is one of the classical methods to compute the eigendecomposition of a matrix. If it is applied on a dense n × n matrix, this algorithm requires O(n3) operations per iteration step. To reduce this complexity for a symmetric matrix to O(n), the original matrix is first reduced to tridiagonal form using orthogonal similarity transformations. In the report (Report TW360, May 2003) a reduction from a symmetric matrix into a similar semiseparable one is described. In this paper a QR algorithm to compute the eigenvalues of semiseparable matrices is designed where each iteration step requires O(n) operations. Hence, combined with the reduction to semiseparable form, the eigenvalues of symmetric matrices can be computed via intermediate semiseparable matrices, instead of tridiagonal ones. The eigenvectors of the intermediate semiseparable matrix will be computed by applying inverse iteration to this matrix. This will be achieved by using an O(n) system solver, for semiseparable matrices. A combination of the previous steps leads to an algorithm for computing the eigenvalue decompositions of semiseparable matrices. Combined with the reduction of a symmetric matrix towards semiseparable form, this algorithm can also be used to calculate the eigenvalue decomposition of symmetric matrices. The presented algorithm has the same order of complexity as the tridiagonal approach, but has larger lower order terms. Numerical experiments illustrate the complexity and the numerical accuracy of the proposed method. Copyright © 2005 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

In this paper, we propose and analyse preconditioning strategies for non-Hermitian Toeplitz linear systems. These techniques used in connection with the GMRES algorithm allow to solve in an optimal way the corresponding linear systems. Copyright © 2004 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

For decades considerable efforts have been exerted to resolve the inverse eigenvalue problem for non-negative matrices. Yet fundamental issues such as the theory of existence and the practice of computation remain open. Recently, it has been proved that, given an arbitrary (n–1)-tuple ℒ = (λ2,...,λn) ∈ ℂn–1 whose components are closed under complex conjugation, there exists a unique positive real number ℛ(ℒ), called the minimal realizable spectral radius of ℒ, such that the set {λ1,...,λn} is precisely the spectrum of a certain n × n non-negative matrix with λ1 as its spectral radius if and only if λ1 ⩾ ℛ(ℒ). Employing any existing necessary conditions as a mode of checking criteria, this paper proposes a simple bisection procedure to approximate the location of ℛ(ℒ). As an immediate application, it offers a quick numerical way to check whether a given n-tuple could be the spectrum of a certain non-negative matrix. Copyright © 2004 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

This paper extends the weighted low rank approximation (WLRA) approach towards linearly structured matrices. In the case of Hankel matrices an equivalent unconstrained optimization problem is derived and an algorithm for solving it is proposed. The correctness of the latter algorithm is verified on a benchmark problem. Finally the statistical accuracy and numerical efficiency of the proposed algorithm is compared with that of STLNB, a previously proposed algorithm for solving Hankel WLRA problems. Copyright © 2004 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

In a recent work, the author introduced a robust multilevel incomplete factorization algorithm using spanning trees of matrix graphs (Proceedings of the 1999 International Conference on Preconditioning Techniques for Large Sparse Matrix Problems in Industrial Applications, Hubert H. Humphrey Center, University of Minnesota, 1999, 251–257). Based on this idea linear and non-linear algebraic multilevel iteration (AMLI) methods are investigated in the present paper. In both cases, the preconditioner is constructed recursively from the coarsest to finer and finer levels. The considered W-cycles only need diagonal solvers on all levels and additionally evaluate a second-degree matrix polynomial (linear case), or, perform ν inner GCG-type iterations (non-linear case) on every other level. This involves the same type of preconditioner for the corresponding Schur complement. The non-linear variant has the additional benefit of being free from any method parameters to be estimated. Based on the same type of approximation property similar convergence rates are obtained for linear and non-linear AMLI, even for a very small number ν of inner iterations, e.g. ν =2,3. The presented methods are robust with respect to anisotropy and discontinuities in the coefficients of the PDEs and can also be applied to unstructured-grid problems. Copyright © 2002 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

This paper defines topology relations of elements treated as overlapping lists of nodes. In particular, the element topology makes use of element faces, element vertices and boundary faces which coincide with the actual (geometrical) faces, vertices and boundary faces in the case of true finite elements. The element topology is used in an agglomeration algorithm to produce agglomerated elements (a non-overlapping partition of the original elements) and their topology is then constructed, thus allowing for recursion. The main part of the algorithms is based on operations on Boolean sparse matrices and the implementation of the algorithms can utilize any available (parallel) sparse matrix format. Applications of the sparse matrix element topology to AMGe (algebraic multigrid for finite element problems), including elementwise constructions of coarse non-linear finite element operators are outlined. An algorithm to generate a block nested dissection ordering of the nodes for generally unstructured finite element meshes is given as well. The coarsening of the element topology is illustrated on a number of fine-grid unstructured triangular meshes. Published in 2002 by John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

This paper gives fast O(n2) algorithms for the factorization of positive-definite and indefinite Hankel matrices. The algorithms are based on the concept of displacement structure and are valid for the more general class of Hankel-like matrices. The positive-definite algorithm is proven to be backward stable. The indefinite algorithm uses a look-ahead step that is naturally suggested by displacement approach. Our error analysis suggests a new criterion for the size of the look-ahead step and our numerical experiments suggest that the use of the new criterion allows us to ensure numerical stability in practice. Copyright © 2001 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

This paper suggests a new algorithm, called IBBLU, for constructing incomplete LU-factorizations of general non-singular sparse matrices. This algorithm is based on the block bordering idea and uses factorized sparse approximate inverses to approximate the inverses of pivot blocks. The triangular factors L and U are represented in explicit-implicit block form, which enhances the flop performance of the preconditioning. The algorithm suggested is theoretically justified for M- and H-matrices, and its competitiveness with the best available algebraic preconditioning methods in both symmetric and unsymmetric cases is demonstrated numerically. Copyright © 2000 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

We consider algebraic multilevel preconditioning methods based on the recursive use of a 2 × 2 block incomplete factorization procedure in which the Schur complement is approximated by a coarse grid matrix. As is well known, for discrete second-order elliptic PDEs, optimal convergence properties are proved for such basic two-level schemes under mild assumptions on the PDE coefficients, but their recursive use in a simple V-cycle algorithm does not generally lead to optimal order convergence. In the present paper, we analyse the combination of these techniques with a smoothing procedure much the same as the one used in standard multigrid algorithms, except that smoothing is not required on the finest grid. Theoretical results prove optimal convergence properties for the V-cycle under an assumption similar to the ‘approximation property’ of the classical multigrid convergence theory. On the other hand, numerical experiments made on both 2D and 3D problems show that the condition number is close to that of the two-level method. Further, the method appears robust in the presence of discontinuity and anisotropy, even when the material interfaces are not aligned with the coarse grid. Copyright © 1999 John Wiley & Sons, Ltd. [ABSTRACT FROM AUTHOR]

The implicit Cholesky algorithm has been developed by several authors during the last 10 years but under different names. We identify the algorithm with a special version of Rutishauser's LR algorithm. Intermediate quantities in the transformation furnish several attractive approximations to the smallest singular value. The paper extols the advantages of using shifts with the algorithm. The nonorthogonal transformations improve accuracy. [ABSTRACT FROM AUTHOR]

In this paper, we demonstrate the scalability of the Reduced System Conjugate Gradient (RSCG) algorithm on distributed memory architectures. We present speed-up results obtained on the Intel iPSC/860 that compare one-, two-, and three-dimensional decompositions of the domain lot both positive definite and positive semidefinite test problems. We develop a model for the RSCG algorithm to analyze computational and communication costs. The model is validated using experimental data and then used to examine and predict behavior of the RSCG algorithm as a function of architecture parameters including communication latency and transmission times and memory access costs. [ABSTRACT FROM AUTHOR]

With the growing demands from disciplinary and interdisciplinary fields of science and engineering for the numerical solution of the nonsymmetric eigenvalue problem, competitive new techniques have been developed for solving the problem. In this paper we examine the state of the art of the algorithmic techniques and the software scene for the problem. Some current developments are also outlined. [ABSTRACT FROM AUTHOR]

The OR algorithm is a basic algorithm for computing the eigenvalues of dense matrices. For efficiency reasons it is prerequisite that the algorithm is applied only after the original matrix has been reduced to a matrix of a particular shape, most notably Hessenberg and tridiagonal, which is preserved during the iterative process. In certain circumstances a reduction to another matrix shape may be advantageous. In this paper, we identify which zero patterns of symmetric matrices are preserved under the OR algorithm. [ABSTRACT FROM AUTHOR]