The Inverses of Block Toeplitz Matrices

Abstract

We study the inverses of block Toeplitz matrices based on the analysis of the block cyclic displacement. New formulas for the inverses of block Toeplitz matrices are proposed. We show that the inverses of block Toeplitz matrices can be decomposed as a sum of products of block circulant matrices. In the scalar case, the inverse formulas are proved to be numerically forward stable, if the Toeplitz matrix is nonsingular and well conditioned.

1. Introduction

Let 𝑇𝑇=⁰𝑇₋₁⋯𝑇_1−𝑛𝑇₁𝑇₀⋱𝑇_2−𝑛𝑇⋮⋱⋱⋮_𝑛−1𝑇_𝑛−2⋯𝑇₀ (1) be an 𝑛×𝑛 block Toeplitz matrix with blocks of size 𝑚×𝑚. We use the shorthand 𝑇=(𝑇_𝑖−𝑗)^𝑛_𝑖,𝑗=1 for a block Toeplitz matrix. The block Toeplitz systems 𝑇𝑥=𝑏 arise in a variety of applications in mathematics, scientific computing, and engineering, for instance, image restoration problems in image processing, numerical differential equations and integral equations, time series analysis, and control theory [1–3]. If we want to solve more than one block Toeplitz linear system with the same coefficient matrix, then we usually solve four or so special block linear systems in order to determine the block Toeplitz inverse formula that expresses 𝑇⁻¹ as the sum of products of block upper and block lower Toeplitz matrices. For example, Van Barel and Bultheel [4] gave an inverse formula for a block Toeplitz matrix and then derived a weakly stable algorithm to solve a block Toeplitz system of linear equations. The special structure of block Toeplitz matrices has resulted in some closed formulas for their inverses.

In the scalar case, Gohberg and Semencul [5] have shown that if the (1,1)st entry of the inverse of a Toeplitz matrix 𝑇 is nonzero, then the first and the last columns of the inverse of the Toeplitz matrix are sufficient to reconstruct 𝑇⁻¹. In [6], an inverse formula can be obtained by the solutions of two equations (the so-called fundamental equations), where each right-hand side of them is a shifted column of the Toeplitz matrix. Later, Ben-Artzi and Shalom [7], Labahn and Shalom [8], Huckle [9], Ng et al. [10], and Heinig [11] have studied the Toeplitz matrix inverse formulas when the (1,1)st entry of the inverse of a Toeplitz matrix is zero. In [12], Cabay and Meleshko presented an efficient algorithm (NPADE) for numerically computing Padé approximants in a weakly stable fashion. As an application of NPADE, it has been shown that it can be used to compute stably, in a weak sense, the inverse of a Hankel or Toeplitz matrix.

When 𝑚>1, additional problems are encountered in obtaining the inverse formula of a block Toeplitz matrix. A well-known formula of Gohberg and Heinig can construct 𝑇⁻¹, provided that the first and last columns together with the first and last rows of the inverse are known [13]. In [14], a set of new formulas for the inverse of a block Hankel or block Toeplitz matrix is given by Labahn et al. The formulas are expressed in terms of certain matrix Padé forms, which approximate a matrix power series associated with the block Hankel matrix. We refer the reader to [15, 16] for the computation of Padé-Hermite and simultaneous Padé systems in detail. In [7], Ben-Artzi and Shalom have proved that each inverse of a Toeplitz matrix can be constructed via three of its columns, and thus, a parametrization of the set of inverses of Toeplitz matrices is obtained. Then they generalized these results to block Toeplitz matrices; see [17]. In [18], Gemignani has shown that the representation of 𝑇⁻¹ relies upon a strong structure-preserving property of the Schur complements of the nonsingular leading principal submatrices of a certain generalized Bezoutian of matrix polynomials.

In this paper, we focus our attention to the inverses of block Toeplitz matrices with the help of the block cyclic displacement. In [19], Ammar and Gader have shown that the inverse of a Toeplitz matrix can be represented as sums of products of lower triangular Toeplitz matrices and circulant matrices. The derivation of their results is based on the idea of cyclic displacement structure. In [20], Gohberg and Olshevsky also obtained new formulas for representation of matrices and their inverses in the form of sums of products of factor circulant, which are based on the analysis of the factor cyclic displacement of matrices. The results in applications to Toeplitz matrices generalized the Gohberg-Semencul, Ben-Artzi-Shalom, and Heinig-Rost formulas. Motivated by a number of related results on Toepltiz inverse formulas, we study the representation of the inverses of block Toeplitz matrices. Since block Toeplitz matrices have similar displacement structure as Toeplitz matrices, all results about Toeplitz matrices extend quite naturally to block Toeplitz matrices.

At first, for purposes of presentation, we adopt some notations that will be used throughout the paper. We denote the 𝑚×𝑚 identity matrix by 𝐼_𝑚. 𝐸_𝑖(𝑖=0,1,…,𝑛−1) means a matrix of dimension 𝑚𝑛×𝑚 in which the (𝑖+1)-th block is the identity matrix 𝐼_𝑚, and the remaining blocks are all zeros. 𝐽 is called an 𝑛×𝑛 block reverse identity matrix, if 𝐽 has the identity matrices on the block antidiagonal and zero elsewhere. Let 𝐴=(𝐴_𝑖𝑗) be a block matrix or a block vector with block size of 𝑚×𝑚. 𝐴^′=(𝐴_𝑗𝑖) is said to be the block transpose of 𝐴. By 𝐶(𝑟,𝜑) we denote the block 𝜑-circulant with the first column 𝑟=[𝑟₀,…,𝑟_𝑛−1]^′ with a block size of 𝑚×𝑚; that is, the matrix of form is defined by 𝑟₀𝜑𝑟_𝑛−1⋯⋯𝜑𝑟₁𝑟₁𝑟₀𝜑𝑟_𝑛−1⋱⋮⋮𝑟₁⋱⋱⋮⋮⋱⋱⋱𝜑𝑟_𝑛−1𝑟_𝑛−1⋯⋯𝑟₁𝑟₀. (2)

In this work, new formulas for the inverses of block Toeplitz matrices are proposed. As is well known, any block Toeplitz matrix 𝑇 has a property of block persymmetry, that is, 𝐽𝑇^′𝐽=𝑇, where 𝐽 is the block reverse identity matrix. Unfortunately, the matrix inverse does not hold this property. This can explain why the Gohberg-Heĭnig formula [13] needs the first and last columns together with the first and last rows of the inverse for constructing the inverse. For a block Toeplitz matrix 𝑇, let 𝑇_𝑐 and 𝑇_𝑠 be its block circulant and block skew-circulant parts, respectively. The representation 𝑇=𝑇_𝑐⋅𝐼+𝑇_𝑠⋅𝐼 shows that there will exist block circulant matrices 𝐴_𝑖(𝑖=1,2) and 𝐵_𝑖(𝑖=1,2) such that 𝑇⁻¹=𝐴₁𝐵₁+𝐴₂𝐵₂. In this paper, by solving four linear equations with the same coefficient matrix 𝑇, we decompose the inverses of block Toeplitz matrices as a sum of products of block circulant matrices. From the point of view of the computation, the formulas involving block circulants instead of block upper or lower Toeplitz matrices are more attractive. For example, in the scalar case, matrix-vector product by a circulant matrix is roughly twice as fast as matrix-vector product by a triangular Toeplitz matrix of the same size. In fact, the efficient multiplication of a Toeplitz matrix and a vector is achieved by embedding the Toeplitz matrix in a circulant matrix of twice the size. The derivation of the formulas we present is based on the factor cyclic displacement of block matrices. Later on, we give some corollaries about the representations of the inverse matrices. The results generalize the Gohberg-Olshevsky, Heinig-Rost, and Ben-Artzi-Shalom formulas. In the scalar case, we consider the stability of the Toeplitz matrix inverse formulas. We show that if the Toeplitz matrix is nonsingular and well conditioned, the inverse formulas are numerically forward stable.

We summarize briefly the content of our paper. Section 2 reviews some fundamental concepts such as the factor cyclic displacement structure of a block matrix and every square block matrix can be written as a sum of products of block circulant matrices. In Section 3, we decompose the inverses of block Toeplitz matrices as a sum of products of block circulant matrices by solving four equations with the same coefficient matrix 𝑇. Some corollaries are also obtained in this section. We discuss the stability of the inverse formulas of scalar Toeplitz matrices in Section 4. Some remarks are given in Section 5.

2. The Factor Block Cyclic Displacement

The inverse formulas derived in the next section depend on some results of the factor block cyclic displacement of a block matrix. The focus of this section lies upon the definitions and the properties of the factor block cyclic displacement of a block matrix. We start with the definition of the factor block cyclic displacement.

Let 𝑍_𝜑 be the block 𝜑-cyclic lower shift matrix defined by 0⋯⋯0𝜑𝐼𝐼0⋱⋱00𝐼⋱⋱⋮⋮⋱⋱⋱⋮0⋯0𝐼0. (3) The block 𝜑-cyclic displacement of an 𝑛×𝑛 block matrix 𝐴 with blocks size of 𝑚×𝑚 is defined as ∇_𝜑(𝐴)=𝐴−𝑍_𝜑𝐴𝑍^′_1/𝜑,𝜑≠0. (4) The number 𝛼=rank∇₁(𝐴) is referred to as block 𝜑-cyclic displacement rank of the block matrix 𝐴. Since all block Toeplitz matrices have a special structure, it is not difficult to obtain that ∇_𝜑(𝑇)=𝑂𝑇₋₁−𝜑𝑇_𝑛−1⋯𝑇_1−𝑛−𝜑𝑇₁𝑇₁−1𝜑𝑇_1−𝑛𝑇𝑂⋯𝑂⋮⋮⋮_𝑛−1−1𝜑𝑇₋₁𝑂⋯𝑂=𝐸₀⋅−𝑇_𝑛𝑇₋₁−𝜑𝑇_𝑛−1⋮𝑇_1−𝑛−𝜑𝑇₁^′+𝑇_𝑛𝑇₁−1𝜑𝑇_1−𝑛⋮𝑇_𝑛−1−1𝜑𝑇₋₁⋅𝐸^′₀, (5) where 𝑇_𝑛 may be an arbitrary 𝑚×𝑚 matrix. Therefore, the 𝜑-cyclic displacement rank of a block Toeplitz matrix is at most 2𝑚.

By straightforward computations, it is easy to see that the equality ∇_𝜑(𝐴)=0 holds if and only if 𝐴 is a block 𝜑-circulant matrix. For any invertible block matrix 𝐴 there exists a simple interrelation between the block 𝜑-cyclic displacement of the inverse matrix 𝐴⁻¹ and the block 𝜑-cyclic displacement of 𝐴. This connection is given by ∇_𝜑𝐴⁻¹=𝐴⁻¹∇_𝜑(𝐴)𝑍_𝜑𝐴⁻¹𝑍^′_1/𝜑. (6) From the previous equation, we know that the block 𝜑-cyclic displacement rank of the inverse matrix 𝐴⁻¹ is the same as that of 𝐴. Moreover, the previous equation also allows us to make use of the block 𝜑-cyclic displacement technique for the calculation of the inverse matrix. We will study the inverses of block Toeplitz matrices with the help of the block cyclic displacement. We need the following result.

Theorem 1. If the block 𝜑-cyclic displacement of 𝐴 is given as the sum 𝐴−𝑍_𝜑𝐴𝑍^′_1/𝜑=^𝛼_𝑚=1𝑔_𝑚ℎ^′_𝑚, (7) where 𝑔_𝑖 and ℎ_𝑖(𝑖=1,…,𝑚) are 𝑛×1 block vectors with blocks of size 𝑚×𝑚, then +𝜑𝐴=Clr,𝜑𝜑−𝜓^𝛼_𝑚=1C𝑔_𝑚,𝜓⋅Cℎ_𝑚,1𝜑^′ (8) or +𝜓𝐴=C𝑙𝑐,𝜑𝜓−𝜑^𝛼_𝑚=1C𝑔_𝑚,𝜑⋅Cℎ_𝑚,1𝜓^′, (9) where C(𝑙𝑟,𝜑) is a block 𝜑-circulant matrix with the same last block row as that of 𝐴 and C(𝑙𝑐,𝜑) is the block 𝜑-circulant matrix with the same last block column as that of 𝐴.

Proof. From computing the block 𝜑-cyclic displacement of matrix 𝐴 defined by (4), it then follows that ∇_𝜑(𝐴)=∇_𝜑𝐶+𝜑lr,𝜑∇𝜑−𝜓^𝜑𝛼_𝑚=1𝐶𝑔_𝑚,𝜓⋅𝐶ℎ_𝑚,1𝜑^′=𝜑𝜑−𝜓^𝛼_𝑚=1∇_𝜑𝐶𝑔_𝑚,𝜓⋅𝐶ℎ_𝑚,1𝜑^′. (10) Note that the 𝜑-cyclic displacement of any block matrix 𝐴 is a linear operator and two block 𝜑-circulant matrices can commute; we obtain ∇_𝜑(𝐴)=𝜑𝜑−𝜓^𝛼_𝑚=1∇_𝜑𝐶𝑔_𝑚,𝜓⋅𝐶ℎ_𝑚,1𝜑^′. (11) By computing the block 𝜑-cyclic displacement of block 𝜓-cyclic matrix 𝐶(𝑔_𝑚,𝜓) for all 𝑚, we have ∇_𝜑(𝐴)=𝜑𝜑−𝜓^𝛼_𝑚=1∇_𝜑𝐶𝑔_𝑚,𝜓⋅𝐶ℎ_𝑚,1𝜑^′=^𝛼_𝑚=1𝑔_𝑚⋅𝐸^′₀⋅𝐶ℎ_𝑚,1𝜑^′−^𝛼_𝑚=1𝐸₀⋅𝑔^′_𝑚⋅𝐶ℎ_𝑚,1𝜑^′, (12) where 𝑔^′_𝑚 for 𝑚=1,…,𝛼 are the first row of the block 𝜑-circulant matrices 𝐶(𝑔_𝑚,𝜑), respectively. Since 𝐸^′₀⋅𝐶(ℎ_𝑚,1/𝜑)^′=ℎ^′_𝑚, it is easy to observe that the sum of the first 𝛼 terms in the previous equality is equal to ∑^𝛼_𝑚=1𝑔_𝑚ℎ^′_𝑚.
In the following, we need to show that the sum of the last 𝛼 terms in the previous equality is equal to the zero matrix. Using the fact that 𝑍^𝑛_𝜑=𝜑𝐼 and 𝑍_𝜑⁻¹=𝑍^′_1/𝜑, we check ^𝑛−1_𝑘=0𝑍^𝑘_𝜑𝑇⁻¹−𝑍_𝜑𝑇⁻¹𝑍_𝜑⁻¹𝑍_𝜑^−𝑘=0. (13) Thus we obtain ^𝑛−1_𝑘=0𝑍^𝑘_𝜑𝑇⁻¹−𝑍_𝜑𝑇⁻¹𝑍_𝜑⁻¹𝑍_𝜑^−𝑘=^𝑛−1𝛼_𝑘=0𝑚=1𝑍^𝑘_𝜑𝑔_𝑚ℎ^′_𝑚𝑍_𝜑^−𝑘=^𝛼_𝑚=1𝐶𝑔_𝑚,𝜑⋅𝐶ℎ_𝑚,1𝜑^′=0. (14) By comparing ∑^𝛼_𝑚=1𝐸₀⋅𝑔^′_𝑚⋅𝐶(ℎ_𝑚,1/𝜑)^′ with ∑^𝛼_𝑚=1𝐶(𝑔_𝑚,𝜑)⋅𝐶(ℎ_𝑚,1/𝜑)^′, we obtain ∑^𝛼_𝑚=1𝐸₀⋅𝑔^′_𝑚⋅𝐶(ℎ_𝑚,1/𝜑)^′=0.
Finally, we discuss why 𝐶(lr,𝜑) is the block 𝜑-circulant with the same last block rows as the matrix 𝐴. We know that the last block row of the matrix 𝐶(𝑔_𝑚,𝜓) is independent of the number 𝜓. In view of ∑^𝛼_𝑚=1𝐶(𝑔_𝑚,𝜑)⋅𝐶(ℎ_𝑚,1/𝜑)=0, the last block rows of the matrix 𝐴 and 𝐶(lr,𝜑) are the same.
Equation (8) is now completely proved. A similar argument shows that (9) can be obtained from (7).

3. Inverses of Block Toeplitz Matrices

In this section, we focus our attention to new formulas on the inverses of block Toeplitz matrices. It is easy to find that any block Toeplitz matrix has a property of block persymmetry, that is, 𝐽𝑇^′𝐽=𝑇, where 𝐽 is the block reverse identity matrix. Unfortunately, the matrix inverse does not have this property. This can explain why the Gohberg-Heĭnig formula [13] needs the first and last block columns together with the first and last block rows of the inverse to construct it. Unlike the Gohberg-Semencul formula [5], it only needs the first and the last columns of the inverse if the (1,1)st entry of the inverse is nonzero.

Let 𝑇𝑥₀,𝑥₁,…,𝑥_𝑛−1^′=𝐸₀,𝑥₀,𝑥₁,…,𝑥_𝑛−1𝑇=𝐸^′_𝑛−1,𝑇𝑢₀,𝑢₁,…,𝑢_𝑛−1^′=𝑇_𝑛,𝑇₁−1𝜑𝑇_1−𝑛,…,𝑇_𝑛−1−1𝜑𝑇₋₁^′,𝑢₀,𝑢₁,…,𝑢_𝑛−1𝑇=𝑇_𝑛−1−1𝜑𝑇₋₁,…,𝑇₁−1𝜑𝑇_1−𝑛,𝑇_𝑛. (15)(16)(17)(18)

For convenience, we will write 𝑥,𝑥,𝑢 and 𝑢 as 𝑥=𝑥₀,𝑥₁,…,𝑥_𝑛−1^′,𝑥=𝑥₀,𝑥₁,…,𝑥_𝑛−1^′,𝑢=𝑢₀,𝑢₁,…,𝑢_𝑛−1^′,𝑢=𝑢₀,𝑢₁,…,𝑢_𝑛−1^′. (19) From (5), (6), (15), and (17), by simple computations we obtain the following equality: ∇_𝜑𝑇⁻¹=−𝑇⁻¹∇_𝜑(𝑇)𝑍_𝜑𝑇⁻¹𝑍^′_1/𝜑𝑥=−⁰𝑥₁⋮𝑥_𝑛−1⋅−𝑇_𝑛,𝑇₋₁−𝜑𝑇_𝑛−1,…,𝑇_1−𝑛−𝜑𝑇₁⋅𝑍_𝜑𝑇⁻¹𝑍^′_1/𝜑−𝑢₀𝑢₁⋮𝑢_𝑛−1⋅𝐼_𝑚,0,…,0⋅𝑍_𝜑𝑇⁻¹𝑍^′_1/𝜑. (20) Substituting (16) and (18) into the previous equality and computing, we obtain ∇_𝜑𝑇⁻¹𝑥=𝜑⁰𝑥₁⋮𝑥_𝑛−1⋅𝑇_𝑛−1−1𝜑𝑇₋₁,…,𝑇₁−1𝜑𝑇_1−𝑛,𝑇_𝑛𝑇⁻¹𝑍^′_1/𝜑𝑢−𝜑⁰𝑢₁⋮𝑢_𝑛−1⋅0,…,0,𝐼_𝑚⋅𝑇⁻¹𝑍^′_1/𝜑=𝑥₀𝑥₁⋮𝑥_𝑛−1𝑢_𝑛−1,𝜑𝑢₀,…,𝜑𝑢_𝑛−2−𝑢₀𝑢₁⋮𝑢_𝑛−1𝑥_𝑛−1,𝜑𝑥₀,…,𝜑𝑥_𝑛−2. (21)

Theorem 2. Let 𝑇 be a block Toeplitz matrix given by (1). If (15)–(18) are solvable for some 𝜑≠0 and an arbitrary matrix 𝑇_𝑛 of order 𝑚, then(a)𝑇is invertible;(b)𝑇⁻¹ can be denoted as 𝑇⁻¹=𝜑C𝐽𝜑−𝜓𝑥,𝜓⋅C𝑢,𝜑−C𝑢−𝜑−𝜓𝜑𝐸₀,𝜓⋅C𝐽𝑥,𝜑, (22)where 𝜓(≠𝜑) is an arbitrary number.

Proof. (a) Assume that there exists nonzero vector 𝑝 of dimension 𝑚𝑛 such that 𝑇𝑝=0. Then by computing (16) and (18), we have (−𝑇_𝑛,𝑇₋₁−𝜑𝑇_𝑛−1,…,𝑇_1−𝑛−𝜑𝑇₁)𝑍_𝜑𝑝=0 and (𝐼_𝑚,0,…,0)𝑍_𝜑𝑝=0. These imply ∇_𝜑(𝑇)𝑍_𝜑𝑝=0. Since ∇_𝜑(𝑇)=𝑇−𝑍_𝜑𝑇𝑍^′_1/𝜑, we obtain 𝑇𝑍_𝜑𝑝=0. A similar argument shows that 𝑇𝑍^𝑘_𝜑𝑝=0(𝑘=2,3,…,𝑛−1). Premultiplying the equalities 𝑇𝑍^𝑘_𝜑𝑝=0(𝑘=0,1,…,𝑛−1) by (𝑥₀,𝑥₁,…,𝑥_𝑛−1), we can conclude that all elements of the vector are equal to zero. Therefore, the invertibility of 𝑇 is now proved.
(b) According to Theorem 1, we get 𝑇⁻¹+𝜑=𝐶lr,𝜑𝐶𝜑−𝜓𝑥,𝜓⋅𝐶𝜑𝑍_1/𝜑1𝑢,𝜑^′+𝐶𝑢,𝜓⋅𝐶𝜑𝑍_1/𝜑1𝑥,𝜑^′. (23)𝐶(lr,𝜑) is the block 𝜑-circulant matrix with the same last block row as that of 𝑇⁻¹. It is not difficult to see that 𝐶(lr,𝜑)=𝐶(𝜑𝑍_1/𝜑𝑥,1/𝜑)^′=𝐶(𝐽𝑥,𝜑). After computations and arrangement, we obtain the formula (22) on the inverses of block Toeplitz matrices.

Remark 3. In fact, the inverse formula (22) is a generalization of the formula of Gohberg and Olshevsky [20], where they consider scalar Toeplitz matrices. Since the matrix inverse does not hold the property of persymmetry, in the block case, for constructing 𝑇⁻¹ we need to solve four linear equations with the same coefficient matrix 𝑇. Labahn et al. also gave this conclusion in terms of certain matrix Padé forms [14].

Remark 4. For the inverse formula (22), we take an interest in two particular cases. When 𝜑=1 and 𝜓=0, the inverse formula (22) can be written as 𝑇⁻¹=𝑥₀𝑥0⋯0₁𝑥₀𝑥⋱⋮⋮⋱⋱0_𝑛−1⋮𝑥₁𝑥₀×𝑢_𝑛−1𝑢₀⋯𝑢_𝑛−2𝑢_𝑛−2𝑢_𝑛−1𝑢⋱⋮⋮⋱⋱₀𝑢₀⋮𝑢_𝑛−2𝑢_𝑛−1−𝑢₀−𝐼_𝑚𝑢0⋯0₁𝑢₀−𝐼_𝑚𝑢⋱⋮⋮⋱⋱0_𝑛−1⋮𝑢₁𝑢₀−𝐼_𝑚×𝑥_𝑛−1𝑥₀⋯𝑥_𝑛−2𝑥_𝑛−2𝑥_𝑛−1𝑥⋱⋮⋮⋱⋱₀𝑥₀⋮𝑥_𝑛−2𝑥_𝑛−1. (24) So we obtain that the inverse of a block Toeplitz matrix can be represented as a sum of products of block lower triangular Toeplitz matrices and block circulant matrices. This is a generalization of the formula of Ammar and Gader [19]. And when 𝜑=1 and 𝜓=−1, the inverse formula (22) can be written as a sum of products of block skew-circulant matrices and block circulant matrices. 𝑇⁻¹=12𝑥₀−𝑥_𝑛−1⋯−𝑥₁𝑥₁𝑥₀⋱⋮⋮⋱⋱−𝑥_𝑛−1𝑥_𝑛−1⋮𝑥₁𝑥₀×𝑢_𝑛−1𝑢₀⋯𝑢_𝑛−2𝑢_𝑛−2𝑢_𝑛−1𝑢⋱⋮⋮⋱⋱₀𝑢₀⋮𝑢_𝑛−2𝑢_𝑛−1−𝑢₀−2𝐼_𝑚−𝑢_𝑛−1⋯−𝑢₁𝑢₁𝑢₀−2𝐼_𝑚⋱⋮⋮⋱⋱−𝑢_𝑛−1𝑢_𝑛−1⋮𝑢₁𝑢₀−2𝐼_𝑚×𝑥_𝑛−1𝑥₀⋯𝑥_𝑛−2𝑥_𝑛−2𝑥_𝑛−1𝑥⋱⋮⋮⋱⋱₀𝑥₀⋮𝑥_𝑛−2𝑥_𝑛−1. (25) This is a generalization of the formula of Ammar and Gader [21], where they considered the Hermitian Toeplitz matrix.

Let 𝑦=(𝑦₀,𝑦₁,…,𝑦_𝑛−1)^′, 𝑧=(𝑧₀,𝑧₁,…,𝑧_𝑛−1)^′, 𝑦^′=(𝑦₀,𝑦₁,…,𝑦_𝑛−1) and 𝑧^′=(𝑧₀,𝑧₁,…,𝑧_𝑛−1) be the solutions of the equations 𝑇𝑦=𝐸_𝑘, 𝑇𝑧=𝐸_𝑘+1, 𝑦^′𝑇=𝐸^′_{𝑛−𝑘−1}, and 𝑧^′𝑇=𝐸^′_{𝑛−𝑘−2} for some integer 𝑘(0≤𝑘≤𝑛−1), respectively. In addition, 𝑦_𝑛−1 and 𝑦₀ are both nonsingular, 𝑦₀⁻¹∑^𝑛−1_𝑖=1𝑦_𝑖𝑇_−𝑛+𝑖=∑^𝑛−1_𝑖=1𝑇_−𝑖𝑦_𝑖−1𝑦⁻¹_𝑛−1 for 0≤𝑘≤𝑛−2 and (1/𝜑)𝑦₀⁻¹+𝑦₀⁻¹∑^𝑛−1_𝑖=1𝑦_𝑖𝑇_−𝑛+𝑖=(1/𝜑)𝑦⁻¹_𝑛−1+∑^𝑛−1_𝑖=1𝑇_−𝑖𝑦_𝑖−1𝑦⁻¹_𝑛−1 for 𝑘=𝑛−1.

For 𝑘<𝑛−1, we have 𝑇_𝑛=𝑇₀+(1/𝜑)∑^𝑛−1_𝑖=1𝑇_−𝑖𝑦_𝑖−1𝑦⁻¹_𝑛−1=𝑇₀+(1/𝜑)𝑦₀⁻¹∑^𝑛−1_𝑖=1𝑘_𝑖𝑇_−𝑛+𝑖. Then we obtain 1𝑇𝑢=𝑇𝜑𝑍_𝜑𝑦−𝑧𝑦⁻¹_𝑛−1=1𝜑𝑇𝜑𝑦_𝑛−10⋮0𝑦⁻¹_𝑛−1+1𝜑𝑇0𝑦₀⋮𝑦_𝑛−2𝑦⁻¹_𝑛−1−1𝜑𝑇𝑧𝑦⁻¹_𝑛−1=𝑇₀𝑇₁⋮𝑇_𝑛−1+1𝜑𝜑𝑇_𝑛−𝑇₀𝑦_𝑛−1−𝑇_−𝑛+1𝑦_𝑛−1⋮−𝑇₋₁𝑦_𝑛−1×𝑦⁻¹_𝑛−1+1𝜑𝐸_𝑘+1𝑦⁻¹_𝑛−1−1𝜑𝐸_𝑘+1𝑦⁻¹_𝑛−1=𝑇_𝑛𝑇₁−1𝜑𝑇_−𝑛+1⋮𝑇_𝑛−1−1𝜑𝑇₋₁,𝑢^′1𝑇=𝜑𝑦₀⁻¹𝑦^′𝑍_𝜑−𝑧^′𝑇=1𝜑𝑦₀⁻¹0,…,0,𝜑𝑦₀𝑇+1𝜑𝑦₀⁻¹𝑦₁,…,𝑦_𝑛−1−1,0𝜑𝑦⁻¹𝑧^′𝑇=𝑇_𝑛−1,…,𝑇₁,𝑇₀+1𝜑𝑦₀⁻¹−𝑦₀𝑇₋₁,…,−𝑦₀𝑇_−𝑛+1,𝜑𝑦₀𝑇_𝑛−𝑇₀+1𝜑𝑦₀⁻¹𝐸^′_{𝑛−𝑘−1}−1𝜑𝑦₀⁻¹𝐸^′_{𝑛−𝑘−1}=𝑇_𝑛−1−1𝜑𝑇₋₁,…,𝑇₁−1𝜑𝑇_−𝑛+1,𝑇_𝑛. (26) For 𝑘=𝑛−1, we have 𝑇_𝑛−𝑇₀=1𝜑𝑦⁻¹_𝑛−1+1𝜑∑^𝑛−1_𝑘=1𝑇_−𝑖𝑦_𝑖−1𝑦⁻¹_𝑛−1=1𝜑𝑦₀⁻¹+(1/𝜑)𝑦₀⁻¹∑^𝑛−1_𝑘=1𝑦_𝑘𝑇_−𝑛+𝑘. A similar argument shows that 𝑇𝑇𝑢=^𝑛𝑇₁−1𝜑𝑇_−𝑛+1⋮𝑇_𝑛−1−1𝜑𝑇₋₁,𝑢^′𝑇=𝑇_𝑛−1−1𝜑𝑇₋₁,…,𝑇₁−1𝜑𝑇_−𝑛+1,𝑇_𝑛. (27)

As a result, we obtain that 𝑢=(1/𝜑)(𝑍_𝜑𝑦−𝑧)𝑦⁻¹_𝑛−1 and 𝑢^′=(1/𝜑)𝑦₀⁻¹(𝑦^′𝑍_𝜑−𝑧^′) are the solutions of (17) and (18), respectively. This result together with Theorem 2 yields the following corollary.

Corollary 5. Let 𝑇 be the block Toeplitz matrix given by (1). If (15) and (16), 𝑇𝑦=𝐸_𝑘, 𝑇𝑧=𝐸_𝑘+1, 𝑦^′𝑇=𝐸^′_{𝑛−𝑘−1}, and 𝑧^′𝑇=𝐸^′_{𝑛−𝑘−2} are solvable for some integer 𝑘(0≤𝑘≤n−1), in addition, 𝑦_𝑛−1 and 𝑦₀ are both nonsingular, 𝑦₀⁻¹∑^𝑛−1_𝑖=1𝑦_𝑖𝑇_−𝑛+𝑖=∑^𝑛−1_𝑖=1𝑇_−𝑖𝑦_𝑖−1𝑦⁻¹_𝑛−1 for 0≤𝑘≤𝑛−2 and (1/𝜑)𝑦₀⁻¹+𝑦₀⁻¹∑^𝑛−1_𝑖=1𝑦_𝑖𝑇_−𝑛+𝑖=(1/𝜑)𝑦⁻¹_𝑛−1+∑^𝑛−1_𝑖=1𝑇_−𝑖𝑦_𝑖−1𝑦⁻¹_𝑛−1 for 𝑘=𝑛−1, then 𝑇 is invertible and 𝑇⁻¹=1𝐶𝐽𝑦𝜑−𝜓𝑥,𝜓⋅𝐶₀⁻¹𝑍^′_𝜑𝑧𝑍𝑦−,𝜑−𝐶_𝜓𝑦−𝑧𝑦⁻¹_𝑛−1,𝜓⋅𝐶𝐽𝑥,𝜑. (28)

Corollary 5 shows how to construct the inverse of a block Toeplitz matrix, provided that the first block column and last block row, along with two successive block columns and block rows of the inverse, are given. More precisely, the choice of 𝑘=𝑛−1 gives that 𝐸_𝑘+1=𝐸₀ and 𝐸_{𝑛−𝑘−2}=𝐸_𝑛−1. This means that 𝑧=𝑥 and 𝑧=𝑥. As a consequence, we have the following result.

Corollary 6. Let 𝑇 be the block Toeplitz matrix given by (1). If (15) and (16), 𝑇𝑦=𝐸_𝑛−1, and 𝑦^′𝑇=𝐸_0′ are solvable, in addition, 𝑦_𝑛−1 and 𝑦₀ are both nonsingular, 𝑦₀⁻¹∑^𝑛−1_𝑘=1𝑦_𝑘𝑇_−𝑛+𝑘=∑^𝑛−1_𝑘=1𝑇_−𝑖𝑦_𝑖−1𝑦⁻¹_𝑛−1, then 𝑇 is invertible and 𝑇⁻¹=1𝐶𝐽𝑦𝜑−𝜓𝑥,𝜓⋅𝐶₀⁻¹𝑍^′_𝜑𝑥𝑍𝑦−,𝜑−𝐶_𝜓𝑦−𝑥𝑦⁻¹_𝑛−1,𝜓⋅𝐶𝐽𝑥,𝜑. (29)

Corollary 6 is a generalization of the well-known formula of Gohberg and Heĭnig. In [13], the authors have shown that the first and last block columns together with the first and last block rows of the inverse are sufficient to construct the inverse 𝑇⁻¹. In addition, we decompose the inverses of block Toeplitz matrices as a sum of products of block circulant matrices.

4. Stability Analysis

In the scalar case, the inverse formula (22) in Section 3 is just the well-known formula of Gohberg and Olshevsky; see Theorem 3.1 in [20]. Hence, solving two linear systems are sufficient to reconstruct the inverse of the Toeplitz matrix 𝑇⁻¹. Let 𝑥=(𝑥₀,𝑥₁,…,𝑥_𝑛−1)^′ and 𝑢=(𝑢₀,𝑢₁,…,𝑢_𝑛−1)^′ be the solutions of (15) and (17), respectively. In the scalar case, formula (22) can be written as 𝑇⁻¹=𝜑𝐶𝜑−𝜓𝑥,𝜓⋅𝐶𝑢,𝜑−𝐶𝑢−𝜑−𝜓𝜑𝐸₀,𝜓⋅𝐶𝑥,𝜑, (30) where 𝜑≠0 and 𝜓(≠𝜑) is an arbitrary number. For convenience, we rewrite 𝑇⁻¹ as 𝑇⁻¹≐𝜑𝐶𝜑−𝜓¹¹𝐶₁₂−𝐶₂₁𝐶₂₂. (31) In this section, we want to show that the evaluation of some scalar inverse formulas in Section 3 is forward stable. An algorithm is called forward stable if for all well-conditioned problems, the computed solution 𝑥 is close to the true solution 𝑥 in the sense that the relative error ‖𝑥−𝑥‖₂/‖𝑥‖₂ is small. In the matrix computation, roundoff errors occur. Let 𝐴,𝐵∈𝐶^𝑛,𝑛 and 𝛼∈𝐶; if we neglect the 𝑂(𝜀²) terms, then for any floating-point arithmetic with machine precision 𝜀, there are fl()‖𝐸‖𝛼𝐴=𝛼𝐴+𝐸,_𝐹≤𝜀|𝛼|‖𝐴‖_𝐹≤𝜀√𝑛|𝛼|‖𝐴‖₂,fl()‖𝐸‖𝐴+𝐵=𝐴+𝐵+𝐸,_𝐹≤𝜀‖𝐴+𝐵‖_𝐹≤𝜀√𝑛‖𝐴+𝐵‖₂,fl(𝐴𝐵)=𝐴𝐵+𝐸,‖𝐸‖_𝐹≤𝜀𝑛‖𝐴‖_𝐹‖𝐵‖_𝐹. (32) See [22]. According to the floating-point arithmetic, we have the following result.

Theorem 7. Let 𝑇 be a nonsingular Toeplitz matrix given by (1) with 𝑚=1 and let it be well conditioned then the inversion formula (30) is forward stable for 𝜑≠0 and a given number 𝜓(≠𝜑).

Proof. Assume that we have computed the solutions 𝑥, 𝑢 in the inverse formula (30) which are perturbed by the normwise relative errors bounded by 𝜀̂‖𝑥‖₂≤‖𝑥‖₂(1+𝜀)̂,‖𝑢‖₂≤‖𝑢‖₂(1+𝜀)̂. (33) Let 𝑁=max{1,𝜓} and 𝑀=max{1,𝜑}. Then, we have 𝐶_11𝐹≤√𝑛𝑁‖𝑥‖₂,𝐶_12𝐹≤√𝑛𝑀‖𝑢‖₂,𝐶_21𝐹≤𝑛𝜑−𝜓𝜑²+𝑛𝑁²‖𝑢‖²₂≤√𝑛𝜑−𝜓𝜑+√𝑛𝑁‖𝑢‖₂,𝐶_22𝐹≤√𝑛𝑀‖𝑥‖₂. (34) By using the perturbed solutions 𝑥, 𝑢, the inverse formula (30) can be expressed as 𝑇⁻¹𝜑=fl𝐶𝜑−𝜓¹¹𝐶₁₂+𝐶₂₁𝐶₂₂𝜑=fl𝐶𝜑−𝜓¹¹+Δ𝐶₁₁𝐶₁₂+Δ𝐶₁₂+𝐶₂₁+Δ𝐶₂₁𝐶₂₂+Δ𝐶₂₂=𝑇⁻¹+𝜑𝜑−𝜓Δ𝐶₁₁𝐶₁₂+𝐶₁₁Δ𝐶₁₂+Δ𝐶₂₁𝐶₂₂+𝐶₂₁Δ𝐶₂₂+𝐸+𝐹+𝐺. (35) Here, 𝐸 is the matrix containing the error which results from computing the matrix products, 𝐹 contains the error from subtracting the matrices, and 𝐺 represents the error of the multiplication by 𝜑/(𝜑−𝜓). For the error matrices Δ𝐶₁₁, Δ𝐶₁₂, Δ𝐶₂₁, and Δ𝐶₂₂, we have Δ𝐶_11𝐹≤𝜀̂𝐶_11𝐹≤𝜀√̂𝑛𝑁‖𝑥‖₂,Δ𝐶_12𝐹≤𝜀̂𝐶_12𝐹≤𝜀√̂𝑛𝑀‖𝑢‖₂,Δ𝐶_21𝐹≤𝜀̂𝐶_21𝐹≤𝜀√̂𝑛𝜑−𝜓𝜑+𝑁‖𝑢‖₂,Δ𝐶_22𝐹≤𝜀̂𝐶_22𝐹=𝜀√̂𝑛𝑀‖𝑥‖₂. (36) It follows that ‖𝐸‖₂≤‖𝐸‖_𝐹≤𝜀𝑛𝐶_11𝐹𝐶_12𝐹+𝐶_21𝐹𝐶_22𝐹≤𝜀𝑛²𝑀‖𝑥‖₂𝜑−𝜓𝜑+2𝑁‖𝑢‖₂,‖𝐹‖₂≤𝜀√𝑛𝑇⁻¹₂,‖𝐺‖₂≤‖𝐺‖_𝐹𝜑≤𝜀𝐶𝜑−𝜓^11𝐹𝐶_12𝐹+𝐶_21𝐹𝐶_22𝐹≤𝜀𝑛𝑀𝜑‖𝑥‖𝜑−𝜓²𝜑−𝜓𝜑+2𝑁‖𝑢‖₂. (37) After adding all these error bounds, we have 𝑇⁻¹−𝑇⁻¹₂≤𝜑2̂𝜑−𝜓𝜀𝑛+𝜀𝑛²‖𝑥‖+𝜀𝑛×𝑀₂𝜑−𝜓𝜑+2𝑁‖𝑢‖₂+𝜀√𝑛𝜑𝑇𝜑−𝜓⁻¹². (38) Note that 𝑥 and 𝑢 are the solutions of (15) and (17), respectively; then ‖𝑥‖₂≤‖𝑇⁻¹‖₂ and ‖𝑢‖₂≤‖𝑇⁻¹‖₂‖𝑓‖₂, where 𝑓=(𝑇_𝑛,𝑇₁−(1/𝜑)𝑇_1−𝑛,…,𝑇_𝑛−1−(1/𝜑)𝑇₋₁)^′. Thus, the relative error is 𝑇⁻¹−𝑇⁻¹₂𝑇⁻¹₂≤𝜑2̂𝜑−𝜓𝜀𝑛+𝜀𝑛²+𝜀𝑛×𝑀𝜑−𝜓𝜑+2𝑁𝑇⁻¹₂𝑓₂+𝜀√𝑛𝜑.𝜑−𝜓 (39) As 𝑇 is well conditioned, thus, ‖𝑇⁻¹‖₂ is finite. Obviously, ‖𝑓‖₂ is finite. Therefore, the formula (30) presented in this section is forward stable.

For example, there are two particular cases of the inverse formula (30). When 𝜑=1 and 𝜓=−1, the inverse formula 𝑇⁻¹=(1/2)(𝐶(𝑥,−1)⋅𝐶(𝑢,1)−𝐶(𝑢−2𝐸₀,−1)⋅𝐶(𝑥,1)), which denotes 𝑇⁻¹ as a sum of products of skew-circulant and circulant matrices, is numerically forward stable for a well-conditioned and nonsingular Toeplitz matrix 𝑇. In fact, this is the formula of Ammar and Gader in [21]. And when 𝜑=1 and 𝜓=0, this is the formula of Ammar and Gader in [19]. The inverse formula 𝑇⁻¹=(𝐶(𝑥,0)⋅𝐶(𝑢,1)−𝐶(𝑢−𝐸₀,0)⋅𝐶(𝑥,1)) in which skew-circulant matrices are replaced by upper triangular factors is also forward stable, if 𝑇 is nonsingular and well conditioned.

We also find that even when 𝜓=∞ and a given number 𝜑≠0, the inverse formula (30) is still numerically forward stable. A similar argument can show that the inverse formulas in Corollaries 5 and 6 for the scalar case are all numerically forward stable when 𝑇 is nonsingular and well conditioned.

5. Concluding Remarks

In this paper, with the help of the block cyclic displacement, the inverses of block Toeplitz matrices are discussed. New formulas for inversion of block Toeplitz matrices are obtained, which can be decomposed in a sum of products of block circulant matrices. In the scalar case, we analyze the stability of the inverse formulas. In addition, it should be remarked that any mention of matrix Padé forms for representing inverses of Toeplitz matrices could also be given in terms of solving block linear systems. Also the proofs of the inverses formulas in [14] could also be extended to prove the inverse formulas in this work.

Acknowledgments

This research is supported by NSFC (61170311), 973 Program (2013CB329404), and Sichuan Province Sci. & Tech. Research Project (2012GZX0080). The authors are grateful to the unknown referee for pertinent comments, especially for suggesting four references, which contributed substantially to this revised work.