Toeplitz Linear Systems Research Articles

Fast Iterative Method with a Second-Order Implicit Difference Scheme for Time-Space Fractional Convection–Diffusion Equation

In this paper we want to propose practical numerical methods to solve a class of initial-boundary problem of time-space fractional convection-diffusion equations (TSFCDEs). To start with, an implicit difference method based on two-sided weighted shifted Gr\"{u}nwald formulae is proposed with a discussion of the stability and convergence. We construct an implicit difference scheme (IDS) and show that it converges with second order accuracy in both time and space. Then, we develop fast solution methods for handling the resulting system of linear equation with the Toeplitz matrix. The fast Krylov subspace solvers with suitable circulant preconditioners are designed to deal with the resulting Toeplitz linear systems. Each time level of these methods reduces the memory requirement of the proposed implicit difference scheme from $\mathcal{O}(N^2)$ to $\mathcal{O}(N)$ and the computational complexity from $O(N^3)$ to $O(N\log N)$ in each iterative step, where $N$ is the number of grid nodes. Extensive numerical example runs show the utility of these methods over the traditional direct solvers of the implicit difference methods, in terms of computational cost and memory requirements.

Convergence rate of GMRES on tridiagonal block Toeplitz linear systems

Iterative methods such as generalized minimal residual (GMRES) method are used to solve large sparse linear systems. This paper is considered the GMRES method for solving tridiagonal block Toeplitz linear systems with diagonal blocks, and establishes upper bounds for GMRES residuals. The coefficient matrix becomes an m-tridiagonal Toeplitz matrix, and tridiagonal toeplitz systems are subcases of these systems. Also, we show that the GMRES method on linear system computes the exact solution in at most N steps.

Fast Toeplitz linear system inversion for solving two-dimensional acoustic inverse problem

Abstract The coefficient inverse problem for the acoustic equation is considered. We propose the method for reconstructing the density based on the N-approximation by the finite system of one-dimensional problems and the two-dimensional M. G. Krein approach. The two-dimensional analogue of the M. G. Krein approach is applied to reduce the non-linear inverse problem to a family of linear integral equations. We consider the fast algorithm for solving the relevant linear system, based on using the block-Toeplitz structure of the matrix. The algorithm applied to the M. G. Krein equation allows to obtain the solution of the whole family of the integral equations by solving only one linear system. Results of numerical calculations are presented.

Essential spectral equivalence via multiple step preconditioning and applications to ill conditioned Toeplitz matrices

In this note, we study the fast solution of Toeplitz linear systems with coefficient matrix Tn(f), where the generating function f is nonnegative and has a unique zero at zero of any real positive order θ. As preconditioner we choose a matrix τn(f) belonging to the so-called τ algebra, which is diagonalized by the sine transform associated to the discrete Laplacian. In previous works, the spectral equivalence of the matrix sequences {τn(f)}n and {Tn(f)}n was proven under the assumption that the order of the zero is equal to 2: in other words the preconditioned matrix sequence {τn−1(f)Tn(f)}n has eigenvalues, which are uniformly away from zero and from infinity. Here we prove a partial generalization of the above result when θ<2. Furthermore, by making use of multiple step preconditioning, we show that the matrix sequences {τn(f)}n and {Tn(f)}n are essentially spectrally equivalent for every θ>2, i.e., for every θ>2, there exist mθ and a positive interval [αθ,βθ] such that all the eigenvalues of {τn−1(f)Tn(f)}n belong to this interval, except at most mθ outliers larger than βθ: while the essential bound from above is proven, the bound from below is only observed numerically. Such a nice property, already known only when θ is an even positive integer greater than 2, is coupled with the fact that the preconditioned sequence has an eigenvalue cluster at one, so that the convergence rate of the associated preconditioned conjugate gradient method is optimal. As a conclusion we discuss possible generalizations and we present selected numerical experiments.

Preconditioned HSS method for large multilevel block Toeplitz linear systems via the notion of matrix‐valued symbol

SummaryWe perform a spectral analysis of the preconditioned Hermitian/skew‐Hermitian splitting (PHSS) method applied to multilevel block Toeplitz linear systems in which the coefficient matrix Tn(f) is associated with a Lebesgue integrable matrix‐valued function f. When the preconditioner is chosen as a Hermitian positive definite multilevel block Toeplitz matrix Tn(g), the resulting sequence of PHSS iteration matrices Mn belongs to the generalized locally Toeplitz class. In this case, we are able to compute the symbol ϕ(f,g) describing the asymptotic eigenvalue distribution of Mnwhen n→∞ and the matrix size diverges. By minimizing the infinity norm of the spectral radius of the symbol ϕ(f,g), we are also able to identify effective PHSS preconditioners Tn(g) for the matrix Tn(f). A number of numerical experiments are presented and commented, showing that the theoretical results are confirmed and that the spectral analysis leads to efficient PHSS methods. Copyright © 2015 John Wiley &amp; Sons, Ltd.

An inexact shift‐and‐invert Arnoldi algorithm for Toeplitz matrix exponential

SummaryWe revisit the shift‐and‐invert Arnoldi method proposed in [S. Lee, H. Pang, and H. Sun.Shift‐invert Arnoldi approximation to the Toeplitz matrix exponential, SIAM J. Sci. Comput., 32: 774–792, 2010] for numerical approximation to the product of Toeplitz matrix exponential with a vector. In this approach, one has to solve two large‐scale Toeplitz linear systems in advance. However, if the desired accuracy is high, the cost will be prohibitive. Therefore, it is interesting to investigate how to solve the Toeplitz systems inexactly in this method. The contribution of this paper is in three regards. First, we give a new stability analysis on the Gohberg–Semencul formula (GSF) and define the GSF condition number of a Toeplitz matrix. It is shown that when the size of the Toeplitz matrix is large, our result is sharper than the one given in [M. Gutknecht and M. Hochbruck.The stability of inversion formulas for Toeplitz matrices, Linear Algebra Appl., 223/224: 307–324, 1995]. Second, we establish a relation between the error of Toeplitz systems and the residual of Toeplitz matrix exponential. We show that if the GSF condition number of the Toeplitz matrix is medium‐sized, then the Toeplitz systems can be solved in a low accuracy. Third, based on this relationship, we present a practical stopping criterion for relaxing the accuracy of the Toeplitz systems and propose an inexact shift‐and‐invert Arnoldi algorithm for the Toeplitz matrix exponential problem. Numerical experiments illustrate the numerical behavior of the new algorithm and show the effectiveness of our theoretical results. Copyright © 2015 John Wiley &amp; Sons, Ltd.

The upper and lower bounds for generalized minimal residual method on a tridiagonal Toeplitz linear system

The generalized minimal residual (GMRES) method is widely used to solve a linear system . This paper establishes upper and lower bounds for GMRES residuals for solving an tridiagonal Toeplitz linear system. For normal matrix A, this problem has been studied previously by Li [Convergence of CG and GMRES on a tridiagonal Toeplitz linear system, BIT 47(3) (2007), 577–599.]. Also, Li and Zhang [The rate of convergence of GMRES on a tridiagonal Toeplitz linear system, Numer. Math. 112 (2009), pp. 267–293.] for non-symmetric matrix A, presented upper bound for GMRES residuals. In fact, our main goal in this paper is to find the upper and lower bounds for GMRES residuals on normal tridiagonal Toeplitz linear systems, and lower bounds for residuals of GMRES on solving non-normal tridiagonal Toeplitz linear systems.

A Preconditioned MINRES Method for Nonsymmetric Toeplitz Matrices

Circulant preconditioning for symmetric Toeplitz linear systems is well established; theoretical guarantees of fast convergence for the conjugate gradient method are descriptive of the convergence seen in computations. This has led to robust and highly efficient solvers based on use of the fast Fourier transform exactly as originally envisaged in [G. Strang, Stud. Appl. Math., 74 (1986), pp. 171--176]. For nonsymmetric systems, the lack of generally descriptive convergence theory for most iterative methods of Krylov type has provided a barrier to such a comprehensive guarantee, though several methods have been proposed and some analysis of performance with the normal equations is available. In this paper, by the simple device of reordering, we rigorously establish a circulant preconditioned short recurrence Krylov subspace iterative method of minimum residual type for nonsymmetric (and possibly highly nonnormal) Toeplitz systems. Convergence estimates similar to those in the symmetric case are established.

Transformations of matrix structures work again

Matrices with the structures of Toeplitz, Hankel, Vandermonde and Cauchy types are omnipresent in modern computations in Sciences, Engineering, and Signal and Image Processing. These four matrix classes have distinct features, but in our 1990 paper in Mathematics of Computation we showed that Vandermonde and Hankel multipliers can be applied to transform each class into the others, and then we demonstrated how by using these transforms we can readily extend any successful matrix inversion algorithm from one of these classes to all the others. The power of this approach was widely recognized later, when novel numerically stable algorithms solved nonsingular Toeplitz linear systems of equations in quadratic (versus classical cubic) arithmetic time based on transforming Toeplitz into Cauchy matrix structures. More recent papers combined the same transformation with a link of the Cauchy matrices to the Hierarchical Semiseparable matrix structure, which is a specialization of matrix representations employed by the Fast Multipole Method. This produced numerically stable algorithms that approximated the solution of a nonsingular Toeplitz linear system of equations in nearly linear arithmetic time. We first revisit the successful method of structure transformation, covering it comprehensively. Then we analyze the latter approximation algorithms for Toeplitz linear systems and extend them to approximate multiplication of Vandermonde and Cauchy matrices by a vector and to approximate solution of Vandermonde and Cauchy linear systems of equations provided that they are nonsingular and well-conditioned. We decrease the arithmetic cost of the known numerical approximation algorithms for these tasks from quadratic to nearly linear, and similarly for the computations with the matrices having structures that generalize the structures of Vandermonde and Cauchy matrices and for polynomial and rational evaluation and interpolation. We also accelerate a little further the known numerical approximation algorithms for a nonsingular Toeplitz or Toeplitz-like linear system by employing distinct transformations of matrix structures, and we briefly comment on some natural research challenges, particularly some promising applications of our techniques to high precision computations.

Spectral behavior of preconditioned non-Hermitian multilevel block Toeplitz matrices with matrix-valued symbol

This note is devoted to preconditioning strategies for non-Hermitian multilevel block Toeplitz linear systems associated with a multivariate Lebesgue integrable matrix-valued symbol. In particular, we consider special preconditioned matrices, where the preconditioner has a band multilevel block Toeplitz structure, and we complement known results on the localization of the spectrum with global distribution results for the eigenvalues of the preconditioned matrices. In this respect, our main result is as follows. Let Ik≔(-π,π)k, let Ms be the linear space of complex s×s matrices, and let f,g:Ik→Ms be functions whose components fij,gij:Ik→C,i,j=1,…,s, belong to L∞. Consider the matrices Tn-1(g)Tn(f), where n≔(n1,…,nk) varies in Nk and Tn(f),Tn(g) are the multilevel block Toeplitz matrices of size n1⋯nks generated by f,g. Then {Tn-1(g)Tn(f)}n∈Nk∼λg-1f, i.e. the family of matrices {Tn-1(g)Tn(f)}n∈Nk has a global (asymptotic) spectral distribution described by the function g-1f, provided g possesses certain properties (which ensure in particular the invertibility of Tn(g) for all n) and the following topological conditions are met: the essential range of g-1f, defined as the union of the essential ranges of the eigenvalue functions λj(g-1f),j=1,…,s, does not disconnect the complex plane and has empty interior. This result generalizes the one obtained by Donatelli, Neytcheva, Serra-Capizzano in a previous work, concerning the non-preconditioned case g=1. The last part of this note is devoted to numerical experiments, which confirm the theoretical analysis and suggest the choice of optimal GMRES preconditioning techniques to be used for the considered linear systems.

An algorithm for solving nonsymmetric penta-diagonal Toeplitz linear systems

Fast algorithms for solving symmetric penta-diagonal Toeplitz linear systems have been proposed in McNally (2010) [7], Nemani (2010) [8] and Xu (2010) [11]. In this paper, we discuss the general nonsymmetric problem and propose an algorithm for solving nonsymmetric penta-diagonal Toeplitz linear systems. Numerical examples are given to illustrate the efficiency of our algorithm.

A parallel linear solver for multilevel Toeplitz systems with possibly several right-hand sides

Abstract A Toeplitz matrix has constant diagonals; a multilevel Toeplitz matrix is defined recursively with respect to the levels by replacing the matrix elements with Toeplitz blocks. Multilevel Toeplitz linear systems appear in a wide range of applications in science and engineering. This paper discusses an MPI implementation for solving such a linear system by using the conjugate gradient algorithm. The implementation techniques can be generalized to other iterative Krylov methods besides conjugate gradient. These techniques include the use of an arbitrary dimensional process grid for handling the multilevel Toeplitz structure, a communication-hiding approach for performing matrix–vector multiplications, the incorporation of multilevel circulant preconditioning for accelerating convergence, an efficient orthogonalization manager for detecting linear dependence in block iterations, and an algorithmic rearrangement to eliminate all-reduce synchronizations. The combined use of these techniques leads to a scalable solver for large multilevel Toeplitz systems, possibly with several right-hand sides. We show experimental results on matrices of size up to the order of one billion with nearly perfect scaling by using up to 1024 MPI processes. We also demonstrate an application of the solver in parameter estimation for analyzing large-scale climate data.

A structure preserving matrix factorization for solving general periodic pentadiagonal Toeplitz linear systems

In the current article, a new algorithm for solving linear systems with the periodic pentadiagonal Toeplitz (PPT) coefficient matrix is presented. This algorithm is based on a structure preserving matrix factorization and Sherman–Morrison–Woodbury formula. Some numerical examples are given in order to illustrate the effectiveness of the proposed algorithms. All of the experiments are performed on a computer with the aid of programs written in MATLAB.

On Structured Variants of Modified HSS Iteration Methods for Complex Toeplitz Linear Systems

The Modified Hermitian and skew-Hermitian splitting (MHSS) iteration method was presented and studied by Bai, Benzi and Chen (Computing, 87(2010), 93-111) for solving a class of complex symmetric linear systems. In this paper, using the properties of Toeplitz matrix, we propose a class of structured MHSS iteration methods for solving the complex Toeplitz linear system. Theoretical analysis shows that the structured MHSS iteration method is unconditionally convergent to the exact solution. When the MHSS iteration method is used directly to complex symmetric Toeplitz linear systems, the computational costs can be considerately reduced by use of Toeplitz structure. Finally, numerical experiments show that the structured MHSS iteration method and the structured MHSS preconditioner are efficient for solving the complex Toeplitz linear system.

AMG preconditioning for nonlinear degenerate parabolic equations on nonuniform grids with application to monument degradation

Motivated by the modelling of marble degradation by chemical pollutants, we consider the approximation by implicit finite differences schemes of nonlinear degenerate parabolic equations in which sharp boundary layers naturally occur. The latter suggests to consider various types of nonuniform griddings, when defining suitable approximation schemes. The resulting large nonlinear systems are treated by Newton methods, while the locally Toeplitz linear systems arising from the Jacobian have to be solved efficiently. To this end, we propose the use of AMG preconditioners and we study the related convergence issues, together with the associated spectral features. We present some numerical experiments supporting our theoretical results on the spectrum of the coefficient matrix of the linear systems, alongside others regarding the numerical simulations in the case of the specific model.

Parallelizing the Conjugate Gradient Algorithm for Multilevel Toeplitz Systems

Multilevel Toeplitz linear systems appear in a wide range of scientific and engineering applications. While several fast direct solvers exist for the basic 1-level Toeplitz matrices, in the multilevel case an iterative solver provides the most general and practical solution. This paper proposes several parallelization techniques that enable an efficient implementation of the conjugate gradient algorithm for solving multilevel Toeplitz systems on distributed-memory machines. The two major differences between this implementation and that for a general sparse linear solver are (1) a communication-efficient approach to handle data expansion and truncation and data transpose simultaneously; (2) the interleaving of matrix-vector multiplications and vector inner product calculations to reduce synchronization cost and latency. Similar ideas can be applied to the implementation of other iterative methods for Toeplitz systems that are not necessarily symmetric positive definite. Scaling results are shown to demonstrate the usefulness of the proposed techniques.

A multicore solution to Block–Toeplitz linear systems of equations

There exist algorithms, also called "fast" algorithms, which exploit the special structure of Toeplitz matrices so that, e.g., allow to solve a linear system of equations in O(n 2) flops. However, ...

Recursive self preconditioning method based on Schur complement for Toeplitz matrices

In this paper, we propose to solve the Toeplitz linear systems T n x = b by a recursive-based method. The method is based on repeatedly dividing the original problem into two subproblems that invol...

DCT- and DST-based splitting methods for Toeplitz systems

New splitting iterative methods for Toeplitz systems are proposed by means of recently developed matrix splittings based on discrete sine and cosine transforms due to Kailath and Olshevsky [Displacement structure approach to discrete-trigonometric transform-based preconditioners of G. Strang type and of T. Chan type, SIAM J. Matrix Anal. Appl. 26 (2005), pp. 706–734]. Theoretical analysis shows that new splitting iterative methods converge to the unique solution of a symmetric Toeplitz linear system. Moreover, an upper bound of the contraction factor of our new splitting iterations is derived. Numerical examples are reported to illustrate the effectiveness of new splitting iterative methods.

Smt: a Matlab toolbox for structured matrices

We introduce the smt toolbox for Matlab. It implements optimized storage and fast arithmetics for circulant and Toeplitz matrices, and is intended to be transparent to the user and easily extensible. It also provides a set of test matrices, computation of circulant preconditioners, and two fast algorithms for Toeplitz linear systems.