On Convergents Infinite Products and Some Generalized Inverses of Matrix Sequences

Abstract

The definition of convergence of an infinite product of scalars is extended to the infinite usual and Kronecker products of matrices. The new definitions are less restricted invertibly convergence. Whereas the invertibly convergence is based on the invertible of matrices; in this study, we assume that matrices are not invertible. Some sufficient conditions for these kinds of convergence are studied. Further, some matrix sequences which are convergent to the Moore-Penrose inverses 𝐴⁺ and outer inverses 𝐴⁽²⁾_𝑇,𝑆 as a general case are also studied. The results are derived here by considering the related well-known methods, namely, Euler-Knopp, Newton-Raphson, and Tikhonov methods. Finally, we provide some examples for computing both generalized inverses 𝐴⁽²⁾_𝑇,𝑆 and 𝐴⁺ numerically for any arbitrary matrix A_𝑚,𝑛 of large dimension by using MATLAB and comparing the results between some of different methods.

1. Introduction and Preliminaries

A scalar infinite product 𝑝=∏^∞_𝑚=1𝑏_𝑚 of complex numbers is said to converge if 𝑏_𝑚 is nonzero for 𝑚 sufficiently large, say 𝑚≥𝑁, and 𝑞=lim_𝑚→∞∏^∞_𝑚=1𝑏_𝑚 exists and is nonzero. If this is so then 𝑝 is defined to be 𝑝=𝑞=∏^𝑁−1_𝑚=1𝑏_𝑚. With this definition, a convergent infinite product vanishes if and only if one of its factors vanishes.

Let {𝐵_𝑚} be a sequence of 𝑘×𝑘 matrices, then ^𝑠_𝑚=𝑟𝐵_𝑚=𝐵_𝑠𝐵_𝑠−1⋯𝐵_𝑟if𝑟≤𝑠𝐼if𝑟>𝑠. (1.1) In [1], Daubechies and Lagarias defined the converges of an infinite product of matrices without the adverb “invertibly” as follows.(i)An infinite product ∏^∞_𝑚=1𝐵_𝑚 of 𝑘×𝑘 matrices is said to be right converges if lim_𝑚→∞𝐵₁𝐵₂⋯𝐵_𝑚 exists, in which case ^∞_𝑚=1𝐵_𝑚=lim_𝑚→∞𝐵₁𝐵₂⋯𝐵_𝑚. (1.2)(ii)An infinite product ∏^∞_𝑚=1𝐵_𝑚 of 𝑘×𝑘 matrices is said to be left converges if lim_𝑚→∞𝐵_𝑚⋯𝐵₂𝐵₁ exists, in which case ^∞_𝑚=1𝐵_𝑚=lim_𝑚→∞𝐵_𝑚⋯𝐵₂𝐵₁. (1.3)

The idea of invertibly convergence of sequence of matrices was introduced by Trench [2, 3] as follows. An infinite product ∏^∞_𝑚=1𝐵_𝑚 of 𝑘×𝑘 matrices is said to be invertibly converged if there is an integer 𝑁 such that 𝐵_𝑚 is invertible for 𝑚≥𝑁, and 𝑄=lim^𝑛_{𝑛→∞𝑚=𝑁}𝐵_𝑚, (1.4) exists and is invertible. In this case, ^∞_𝑚=1𝐵_𝑚=𝑄^𝑁−1_𝑚=1𝐵_𝑚. (1.5) Let us recall some concepts that will be used below. Before starting, throughout we consider matrices over the field of complex numbers ℂ or real numbers ℝ. The set of 𝑚-by-𝑛 complex matrices is denoted by 𝑀_𝑚,𝑛(ℂ)=ℂ^𝑚×𝑛. For simplicity, we write 𝑀_𝑚,𝑛 instead of 𝑀_𝑚,𝑛(ℂ) and when 𝑚=𝑛, we write 𝑀_𝑛 instead of 𝑀_𝑛,𝑛. The notations 𝐴^𝑇, 𝐴^∗, 𝐴⁺, 𝐴⁽²⁾_𝑇,𝑆, rank(𝐴), rang(𝐴), null(𝐴), 𝜌(𝐴), ‖𝐴‖_𝑠, ‖𝐴‖_𝑝, and 𝜎(𝐴) stand, respectively, for the transpose, conjugate transpose, Moore-Penrose inverse, outer inverse, rank, range, null space, spectral radius, spectrum norm, 𝑝-norm, and the set of all eigenvalues of matrix 𝐴.

The Moore-Penrose and outer inverses of an arbitrary matrix (including singular and rectangular) are very useful in various applications in control system analysis, statistics, singular differential and difference equations, Markov chains, iterative methods, least-square problem, perturbation theory, neural networks problem, and many other subjects were found in the literature (see, e.g., [4–14]).

It is well known that Moore-Penrose inverse (MPI) of a matrix 𝐴∈𝑀_𝑚,𝑛 is defined to be the unique solution of the following four matrix equations (see, e.g. [4, 11, 14–20]): 𝐴𝑋𝐴=𝐴,𝑋𝐴𝑋=𝑋,(𝐴𝑋)^∗=𝐴𝑋,(𝑋𝐴)^∗=𝑋𝐴 (1.6) and is often denoted by 𝑋=𝐴⁺∈𝑀_𝑛,𝑚. In particular, when 𝐴 is a square and nonsingular matrix, then 𝐴⁺ reduce to 𝐴⁻¹.

For 𝑥=𝐴⁺𝑏, 𝑥^′∈ℂ^𝑛⧵{𝑥} arbitrary, it holds, (see, e.g., [14, 18]), ‖𝑏−𝐴𝑥‖²₂=(𝑏−𝐴𝑥)^∗(𝑏−𝐴𝑥)≤𝑏−𝐴𝑥^′2₂, (1.7) and ‖𝑏−𝐴𝑥‖₂=‖𝑏−𝐴𝑥^′‖₂, then ‖𝑥‖²₂=𝑥^∗𝑥<𝑥^′2₂. (1.8) Thus, 𝑥=𝐴⁺𝑏 is the unique minimum least-squares solution of the following linear squares problem, (see, e.g., [14, 21, 22]), ‖𝑏−𝐴𝑥‖₂=min_{𝑧∈𝐶𝑛}‖𝑏−𝐴𝑧‖₂. (1.9) It is well known also that the singular value decomposition of any rectangular matrix 𝐴∈𝑀_𝑚,𝑛 with rank(𝐴)=𝑟≠0 is given by 𝐴=𝑈𝑉𝐷000^∗∶𝑈^∗𝑈=𝐼_𝑚,𝑉^∗𝑉=𝐼_𝑛, (1.10) where 𝐷=diag(𝜇₁,𝜇₂,…,𝜇_𝑟)∈𝑀_𝑟 is a diagonal matrix with diagonal entries 𝛿_𝑖(𝑖=1,2,…,𝑟), and 𝜇₁≥𝜇₂≥⋯≥𝜇_𝑟>0 are the singular values of 𝐴, that is, 𝜇²_𝑖(𝑖=1,2,…,𝑟) are the nonzero eigenvalues of 𝐴^∗𝐴. This decomposition is extremely useful to represent the MPI of 𝐴∈𝑀_𝑚,𝑛 by [20, 23] 𝐴⁺𝐷=𝑉⁻¹0𝑈00^∗, (1.11) where 𝐷⁻¹=diag(𝜇₁⁻¹,𝜇₂⁻¹,…,𝜇_𝑟⁻¹)∈𝑀_𝑟 is a diagonal matrix with diagonal entries 𝜇_𝑖⁻¹(𝑖=1,…,𝑟).

Furthermore, the spectral norm of 𝐴 is defined by ‖𝐴‖_𝑠=max_{1≤𝑖≤𝑟}𝜇_𝑖=𝜇₁;𝐴⁺_𝑠=1𝜇_𝑟, (1.12) where 𝜇₁and 𝜇_𝑟 are, respectively, the largest and smallest singular value of 𝐴.

Generally speaking, the outer inverse 𝐴⁽²⁾_𝑇,𝑆 of a matrix 𝐴∈𝑀_𝑚,𝑛, which is a unique matrix 𝑋∈𝑀_𝑛,𝑚 satisfying the following equations (see, e.g., [20, 24–27]): 𝐴𝑋𝐴=𝐴,rang(𝐴)=𝑇,null(𝐴)=𝑆, (1.13) where 𝑇 is a subspace of ℂ^𝑛 of 𝑠≤𝑟, and 𝑆 is a subspace of ℂ^𝑚 of dimension 𝑚−𝑠.

As we see in [13, 20, 24–29], it is well-known fact that several important generalized inverses, such as the Moore-Penrose inverse 𝐴⁺, the weighted Moore-Penrose inverse 𝐴⁺_𝑀,𝑁, the Drazin inverse 𝐴^𝐷, and so forth, are all the generalized inverse 𝐴⁽²⁾_𝑇,𝑆, which is having the prescribed range 𝑇 and null space 𝑆 of outer inverse of 𝐴. In this case, the Moore-Penrose inverse 𝐴⁺ can be represented in outer inverse form as follows [27]: 𝐴⁺=𝐴⁽²⁾_{rang(𝐴∗),null(𝐴∗)}. (1.14) Also, the representation and characterization for the outer generalized inverse 𝐴⁽²⁾_𝑇,𝑆 have been considered by many authors (see, e.g., [15, 16, 20, 27, 30, 31])

Finally, given two matrices 𝐴=[𝑎_𝑖𝑗]∈𝑀_𝑚,𝑛, and 𝐵=[𝑏_𝑘𝑙]∈𝑀_𝑝,𝑞, then the Kronecker product of 𝐴 and 𝐵 is defined by (see, e.g., [5, 7, 32–35]) 𝐴⊗𝐵=𝑎_𝑖𝑗𝐵∈𝑀_{𝑚𝑝,𝑛𝑞}. (1.15) Furthermore, the Kronecker product enjoys the following well-known and important properties:(i)The Kronecker product is associative and distributive with respect to matrix addition.(ii)If 𝐴∈𝑀_𝑚,𝑛, 𝐵∈𝑀_𝑝,𝑞, 𝐶∈𝑀_𝑛,𝑟 and 𝐷∈𝑀_𝑞,𝑠, then (𝐴⊗𝐵)(𝐶⊗𝐷)=𝐴𝐶⊗𝐵𝐷 (1.16)(iii)If 𝐴∈𝑀_𝑚, 𝐵∈𝑀_𝑝 are positive definite matrices, then for any real number 𝑟, we have: (𝐴⊗𝐵)^𝑟=𝐴^𝑟⊗𝐵^𝑟. (1.17)(iv)If 𝐴∈𝑀_𝑚,𝑛, 𝐵∈𝑀_𝑝,𝑞, then (𝐴⊗𝐵)⁺=𝐴⁺⊗𝐵⁺. (1.18)

2. Convergent Moore-Penrose Inverse of Matrices

First, the need to compute 𝐴⁺ by using sequences method. The key to such results below is the following two Lemmas, due to Wei [23] and Wei and Wu [17], respectively.

Lemma 2.1. Let 𝐴∈𝑀_𝑚,𝑛 ba a matrix. Then 𝐴⁺=𝐴⁻¹𝐴^∗, (2.1) where 𝐴=𝐴^∗𝐴⇃_rang(𝐴^∗₎ is the restriction of 𝐴^∗𝐴 on rang(𝐴^∗).

Lemma 2.2. Let 𝐴∈𝑀_𝑚,𝑛 with rank(𝐴)=𝑟 and 𝐴=𝐴^∗𝐴⇃_rang(𝐴^∗₎. Suppose Ω is an open set such that 𝜎(𝐴)⊂Ω⊂(0,∞). Let {𝑆_𝑛(𝑥)} be a family of continuous real valued function on Ω with lim_𝑛→∞𝑆_𝑛(𝑥)=1/𝑥 uniformly on𝜎(𝐴). Then 𝐴⁺=lim_𝑛→∞𝑆_𝑛𝐴𝐴^∗. (2.2) Furthermore, 𝑆_𝑛𝐴𝐴^∗−𝐴⁺₂≤sup_𝑥∈𝜎(𝐴₎𝑆_𝑛(𝑥)𝑥−1𝐴⁺₂. (2.3) Moreover, for each 𝜆∈𝜎(𝐴), we have 𝜇²_𝑟≤𝜆≤𝜇²₁. (2.4)

It is well known that the inverse of an invertible operator can be calculated by interpolating the function 1/𝑥, in a similar manner we will approximate the Moore-Penrose inverse by interpolating function 1/𝑥 and using Lemmas 2.1 and 2.2.

One way to produce a family of functions {𝑆_𝑛(𝑥)} which is suitable for use in the Lemma 2.2 is to employ the well known Euler-Knopp method. A series ∑^∞_𝑛=0𝑎_𝑛 is said to be Euler-Knopp summable with parameter 𝛼>0 to the value 𝑎 if the sequence defined by 𝑆_𝑛=𝛼^𝑛𝑘_𝑘=0𝑗=0𝑘𝑗(1−𝛼)^𝑘−𝑗𝛼^𝑗𝛼_𝑗 (2.5) converges to 𝑎. If 𝑎_𝑘=(1−𝑥) for 𝑘=0,1,2,…, then we obtain as the Euler-Knopp transform of the series ∑^∞_𝑘=0(1−𝑥)^𝑘, the sequence given by 𝑆_𝑛(𝑥)=𝛼^𝑛_𝑘=0(1−𝛼𝑥)^𝑘. (2.6) Clearly lim_𝑛→∞𝑆_𝑛(𝑥)=1/𝑥 uniformly on any compact subset of the set 𝐸_𝛼={||}=2𝑥∶1−𝛼𝑥<1𝑥∶0<𝑥<𝛼. (2.7) Another way to produce a family functions {𝑆_𝑛(𝑥)} which is suitable also for use in the Lemma 2.2 is to employ the well-known Newton-Raphson method. This can be done by generating a sequence 𝑦_𝑛, where 𝑦_𝑛+1=𝑦_𝑛−𝑠𝑦_𝑛𝑠^′𝑦_𝑛=𝑦_𝑛2−𝑥𝑦_𝑛, (2.8) for suitable 𝑦₀. Suppose that for 𝛼>0 we define a sequence of functions {𝑆_𝑛(𝑥)} by 𝑆₀(𝑥)=𝛼;𝑆_𝑛+1(𝑥)=𝑆_𝑛(𝑥)2−𝑥𝑆_𝑛(𝑥). (2.9) In fact, 𝑥𝑆_𝑛+1(𝑥)−1=−𝑥𝑆_𝑛(𝑥)−1². (2.10) Iterating on this equality, it follows that if 𝑥 is confined to a compact subset of 𝐸_𝛼={𝑥∶0<𝑥<2/𝛼}. Then there is a constant 𝛽 (defining on this compact set) with 0<𝛽<1 and 𝑥𝑆_𝑛(𝑥)−1=|𝛼𝑥−1|^2𝑛≤𝛽^2𝑛⟶0(𝑛⟶∞). (2.11) According to the variational definition, 𝐴⁺𝑏 is the vector 𝑥∈ℂ^𝑛 which minimizes the functional ‖𝐴𝑥−𝑏‖₂ and also has the smallest 2 norm among all such minimizing vectors. The idea of Tikhonov’s regularization [36, 37] of order zero is to approximately minimize both the functional ‖𝐴𝑥−𝑏‖₂ and the norm ‖𝑥‖₂ by minimizing the functional 𝑔∶ℂ^𝑛→ℝ defined by 𝑔(𝑥)=‖𝐴𝑥−𝑏‖²₂+𝑡‖𝑥‖²₂, (2.12) where 𝑡>0. The minimum of this functional will occur at the unique stationary point 𝑢 of 𝑔, that is, the vector 𝑢 which satisfies ∇𝑔(𝑢)=0. The gradient of 𝑔 is given by ∇𝑔(𝑥)=2𝐴^∗𝐴𝑥−𝐴^∗𝑏+2𝑡𝑥, (2.13) and hence the unique minimizer 𝑢_𝑡 satisfies 𝑢_𝑡=𝐴^∗𝐴+𝑡𝐼⁻¹𝐴^∗𝑏. (2.14) On intuitive grounds, it seems reasonable to expect that lim_𝑡→0+𝑢_𝑡=𝐴^∗𝐴⁻¹𝐴^∗𝑏=𝐴⁺𝑏. (2.15) Therefore, if we define a sequence of functions {𝑆_𝑛(𝑥)} by using Euler-Knopp method, Newton-Raphson method and the idea of Tikhonov’s regularization that mentioned above, then we get the following nice Theorem.

Theorem 2.3. Let 𝐴∈𝑀_𝑚,𝑛 with rank(𝐴)=𝑟 and 0<𝛼<2𝜇₁⁻². Then(i)the sequence {𝐴_𝑛} defined by 𝐴₀=𝛼𝐴^∗;𝐴_𝑛+1=1−𝛼𝐴^∗𝐴𝐴_𝑛+𝛼𝐴^∗ (2.16) converges to 𝐴⁺. Furthermore, the error estimate is given by 𝐴_𝑛−𝐴⁺₂≤𝛽^𝑛+1𝐴⁺₂, (2.17) where 0<𝛽<1.(ii)The sequence {𝐴_𝑛} defined by 𝐴₀=𝛼𝐴^∗;𝐴_𝑛+1=𝐴_𝑛2𝐼−𝐴𝐴_𝑛 (2.18) converges to 𝐴⁺. Furthermore, the error estimate is given by 𝐴_𝑛−𝐴⁺₂≤𝛽^2𝑛𝐴⁺₂. (2.19) where 0<𝛽<1.(iii)for 𝑡>0, 𝐴⁺=lim_𝑡→0+𝑡𝐼+𝐴^∗𝐴⁻¹𝐴^∗. (2.20)Thus, the error estimate is given by 𝑡𝐼+𝐴^∗𝐴⁻¹𝐴^∗−𝐴⁺₂≤𝑡𝜇²_𝑟𝐴+𝑡⁺². (2.21)

Proof. (i) It follows from 𝜎(𝐴)⊆[𝜇²_𝑟,𝜇²₁] that 𝜎(𝐴)⊂(0,𝜇²₁], and hence we apply Lemma 2.2 if we choose the parameter 𝛼 is such a way that (0,𝜇²₁]⊆𝐸_𝛼, where 𝐸_𝛼 is defined by (2.7). We may choose 𝛼 such that 0<𝛼<2𝜇₁⁻². If we use the sequence defined by 𝑆₀(𝑥)=𝛼;𝑆_𝑛+1(𝑥)=(1−𝛼𝑥)𝑆_𝑛(𝑥)+𝛼, (2.22) it is easy to verify that lim_𝑛→∞𝑆_𝑛(𝑥)=1𝑥, (2.23) uniformly on any compact subset of 𝐸_𝛼. Hence, if 0<𝛼<2𝜇₁⁻², then, applying Lemma 2.2, we get lim_𝑛→∞𝑆_𝑛𝐴𝐴^∗=𝐴⁺. (2.24) But it is easy to see from (2.22) that 𝑆_𝑛(𝐴)𝐴^∗=𝐴_𝑛, where 𝐴_𝑛 is given by (2.16). This is surely the case if 0<𝛼<2𝜇₁⁻², then, for such 𝛼, we have the representation 𝐴⁺=𝛼^𝑛_𝑗=01−𝛼𝐴^∗𝐴^𝑗𝐴^∗. (2.25) Note that if we set 𝐴_𝑛=𝛼^𝑛_𝑗=01−𝛼𝐴^∗𝐴^𝑗𝐴^∗, (2.26) then we get (2.16).
To derive an error estimate for the Euler-Knopp method, suppose that 0<𝛼<2𝜇₁⁻². If the sequence 𝑆_𝑛(𝑥) is defined as in (2.22), then 𝑆_𝑛+1(𝑥)𝑥−1=(1−𝛼𝑥)𝑆_𝑛(𝑥)𝑥−1. (2.27) Therefore, since 𝑆₀=𝛼, 𝑆_𝑛(𝑥)𝑥−1=|1−𝛼𝑥|^𝑛+1. (2.28) By 𝜇²_𝑟≤𝑥≤𝜇²₁ for 𝑥∈𝜎(𝐴) and 0<𝛼<2𝜇₁⁻², it follows that |1−𝛼𝑥|<𝛽, (2.29) where 𝛽 is given by 𝛽=max1−𝛼𝜇²₁,1−𝛼𝜇²_𝑟. (2.30) Clearly, 2>𝛼𝜇²₁≥𝛼𝜇²_𝑟>0, (2.31) and therefore 0<𝛽<1. From Lemma 2.2, we establish (2.17).
(ii) Using the Newton-Raphson iterations in (2.8)–(2.11) in conjunction with Lemma 2.2, we see that the sequence of {𝑆_𝑛(𝐴)} defined by 𝑆₀(𝑥)=𝛼𝐼;𝑆_𝑛+1𝐴=𝑆_𝑛𝐴2𝐼−𝐴^∗𝐴𝑆_𝑛𝐴 (2.32) has the property that lim_𝑛→∞𝑆_𝑛(𝐴)𝐴^∗=𝐴⁺ uniformly in ℂ^𝑚×𝑛. If we set 𝐴_𝑛=𝑆_𝑛(𝐴)𝐴^∗, then we get (2.18).
If 𝑥∈𝜎(𝐴) and 0<𝛼<2𝜇₁⁻², then we see that |1−𝛼𝑥|≤𝛽, where 𝛽 is given by (2.30). It follows as in (2.11) and hence from Lemma 2.2, then we get the error bound as in (2.19).
(iii) If we set 𝑆_𝑡(𝑥)=(𝑡+𝑥)⁻¹∶𝑡>0 (2.33) in Lemma 2.2 and using the idea of Tikhonov’s regularization as in (2.12)–(2.15), then it is easy to get (2.20) and (2.21).

Huang and Zhang [28] presented the following sequence which is convergent to 𝐴⁺: 𝐴₀=𝛼₀𝐴^∗;𝐴_𝑛+1=𝛼_𝑛+12𝐴_𝑛−𝐴_𝑛𝐴𝐴_𝑛,𝑛=1,2,…, (2.34) where 𝛼_𝑛+1∈[1,2] is called an acceleration parameter and chosen so as to minimize that which is bound on the maximum distance of any nonzero singular value of 𝐴_𝑛+1𝐴 from 1. They chose 𝐴₀ according to the first term of sequence (2.34) with 𝛼₀=2/(𝜇₁+𝜇_𝑟), and let 𝑝₀=𝛼₀𝜇₁, then the acceleration parameters 𝛼_𝑛+1 and 𝑝_𝑛+1 have the following sequences: 𝛼_𝑛+1=21+2−𝑝_𝑛𝑝_𝑛,𝑝_𝑛+1=𝛼_𝑛+12−𝑝_𝑛𝑝_𝑛. (2.35) We point out that the iteration (2.18) is a special case of an acceleration iteration (2.34). Further, we note that the above methods and the first-order iterative methods used by Ben-Israel and Greville [4] for computing 𝐴⁺ are a set of instructions for generating a sequence {𝐴_𝑛∶𝑛=1,2,3,…} converging to 𝐴⁺.

Similarly, Liu et al. [25] introduced some necessary and sufficient conditions for iterative convergence to the generalized inverse 𝐴⁽²⁾_𝑇,𝑆 and its existence and estimated the error bounds of the iterative methods for approximating 𝐴⁽²⁾_𝑇,𝑆 by defining the sequence {𝐴_𝑛} in the following way: 𝐴_𝑛=𝐴_𝑛−1+𝛽𝑋𝐼−𝐴𝐴_𝑛−1,𝑛=1,2,…, (2.36) where 𝛽≠0 with 𝑋≠𝑋𝐴𝐴₀. Then the iteration (2.36) converges if and only 𝜌(𝐼−𝛽𝐴𝑋)<1, equivalently, 𝜌(𝐼−𝛽𝐴𝑋)<1. In which case, if rang(𝑋)=𝑇, null(𝑋)=𝑆, and rang(𝐴₀)⊂𝑇. Then, 𝐴⁽²⁾_𝑇,𝑆 exists and {𝐴_𝑛} converges to 𝐴⁽²⁾_𝑇,𝑆. Furthermore, the error estimate is given by 𝐴⁽²⁾_𝑇,𝑆−𝐴_𝑛≤𝛽𝑞^𝑛‖𝑋1−𝑞‖𝐼−𝐴𝐴₀=𝑅𝛽,𝑛, (2.37) where 𝑞=min{‖𝐼−𝛽𝑋𝐴‖,‖𝐼−𝛽𝐴𝑋‖}.

What is the best value 𝛽_opt such that 𝜌(𝐼−𝛽𝐴𝑋) minimize in order to achieve good convergence? Unfortunately, it may be very difficult and still require further studies. If 𝜎(𝐴𝑋) is a subset of ℝ and 𝜆_min=min{𝜆∶𝜆∈𝜎(𝐴𝑋)}>0, analogous to [38, Example  4.1], we can have 𝛽_opt=2𝜆_min+𝜌(𝐴𝑋). (2.38)

3. Convergent Infinite Products of Matrices

Trench [3, Definition  1] defined invertibly convergence of an infinite products matrices ∏^∞_𝑚=1𝐵_𝑚 as in the invertible of matrices 𝐵_𝑚 for all 𝑚>𝑁 (where 𝑁 is an integer number). Here, we define the less restricted definitions of convergence of an infinite products ∏^∞_𝑚=1𝐵_𝑚 and ∏^∞_𝑚=1⊗𝐵_𝑚 for 𝑘×𝑘 complex matrices such that ^𝑠_𝑚=𝑟⊗𝐵_𝑚=𝐵_𝑠⊗𝐵_𝑠−1⊗⋯⊗𝐵_𝑟if𝑟≤𝑠𝐼if𝑟>𝑠 (3.1) as follows.

Definition 3.1. Let {𝐵_𝑚} be a sequence of 𝑘×𝑘 matrices. Then An infinite product ∏^∞_𝑚=1𝐵_𝑚 is said to be convergent if there is an integer 𝑁 such that 𝐵_𝑚≠0 (may 𝐵_𝑚 is invertible or not) for 𝑚≥𝑁 and 𝑄=lim^𝑛_{𝑛→∞𝑚=𝑁}𝐵_𝑚 (3.2) exists and is nonzero. In this case, we define ^∞_𝑚=1𝐵_𝑚=𝑄^𝑁−1_𝑚=1𝐵_𝑚. (3.3) Similarly, an infinite product ∏^∞_𝑚=1⊗𝐵_𝑚 converges if there is an integer 𝑁 such that 𝐵_𝑚≠0 for 𝑚≥𝑁, and 𝑅=lim^𝑛_{𝑛→∞𝑚=𝑁}⊗𝐵_𝑚 (3.4) exists and is nonzero. In this case, we define ^∞_𝑚=1⊗𝐵_𝑚=𝑅⊗^𝑁−1_𝑚=1⊗𝐵_𝑚. (3.5) In the above Definition 3.1, the matrix 𝑅 may be singular even if 𝐵_𝑚 is nonsingular for all 𝑚≥1, but 𝑅 may be singular if 𝐵_𝑚 is singular for some 𝑚≥1. However, this definition does not require that 𝐵_𝑚 is invertible for large 𝑚.

Definition 3.2. Let {𝐵_𝑚} be a sequence of 𝑘×𝑘 matrices. Then an infinite product ∏^∞_𝑚=1⊗𝐵_𝑚 is said to be invertibly convergent if there is an integer 𝑁 such that 𝐵_𝑚 is invertible for 𝑚≥𝑁, and 𝑅=lim^𝑛_{𝑛→∞𝑚=𝑁}⊗𝐵_𝑚 (3.6) exists and invertible. In this case, we define ^∞_𝑚=1𝐵_𝑚=𝑅⊗^𝑁−1_𝑚=1⊗𝐵_𝑚. (3.7) The Definitions 3.1 and 3.2 have the following obvious consequence.

Theorem 3.3. Let {𝐵_𝑚} be a sequence of k×k matrices such that the infinite products ∏^∞_𝑚=1𝐵_𝑚 and ∏^∞_𝑚=1⊗𝐵_𝑚 are invertibly convergent. Then, both infinite products are convergent, but the converse, in general, is not true.

Theorem 3.4. Let {𝐵_𝑚} be a sequence of 𝑘×𝑘 matrices such that(i)if the infinite product ∏^∞_𝑚=1𝐵_𝑚 is convergent, then lim_𝑚→∞𝐵_𝑚=𝑄𝑄⁺, (3.8) where 𝑄=lim_𝑛→∞∏^𝑛_𝑚=𝑁𝐵_𝑚.(ii)If the infinite product ∏^∞_𝑚=1⊗𝐵_𝑚 is converge, then lim_𝑚→∞𝐵_𝑚=𝑅⊗𝑅⁺. (3.9)

Proof. (i) Suppose that ∏^∞_𝑚=1𝐵_𝑚 is convergent such that 𝐵_𝑚≠0 when 𝑚≥𝑁. Let 𝑄=∏^𝑛_𝑚=𝑁𝐵_𝑚. Then lim_𝑛→∞𝑄_𝑛=𝑄, where 𝑄≠0. Therefore, lim_𝑛→∞𝑄⁺_𝑛=𝑄⁺. Since 𝐵_𝑛=𝑄_𝑛𝑄⁺_𝑛−1, we then have lim_𝑛→∞𝐵_𝑛=lim_𝑛→∞𝑄_𝑛lim_𝑛→∞𝑄⁺_𝑛−1=𝑄𝑄⁺. (3.10) But if 𝐵_𝑚 is invertible when 𝑚≥𝑁, then 𝑄 is invertible and lim_𝑛→∞𝐵_𝑛=lim_𝑛→∞𝑄_𝑛lim_𝑛→∞𝑄⁻¹_𝑛−1=𝑄𝑄⁻¹=𝐼. (3.11) Similarly, it is easy to prove (ii).

If the infinite products ∏^∞_𝑚=1𝐵_𝑚 and ∏^∞_𝑚=1⊗𝐵_𝑚 are invertibly convergent in Theorem 3.4, then we get the following corollary.

Corollary 3.5. Let {𝐵_𝑚} be a sequence of 𝑘×𝑘 matrices such that(i)if the infinite product ∏^∞_𝑚=1𝐵_𝑚 is invertibly convergent, then lim_𝑚→∞𝐵_𝑚=𝐼 (3.12)(ii)if the infinite product ∏^∞_𝑚=1⊗𝐵_𝑚 is invertibly convergent, then lim_𝑚→∞𝐵_𝑚=𝑅⊗𝑅⁻¹, (3.13) where 𝑅=lim_𝑛→∞∏^𝑛_𝑚=𝑁⊗𝐵_𝑚.

The main reason for interest in these products above is to generate matrix sequences for solving such matrix problems such as singular linear systems and singular coupled matrix equations. For example, Cao [21] and Shi et al. [22] constructed the general stationary and nonstationary iterative process generated by {𝐵_𝑚}^∞_𝑚=1 for solving the singular linear system 𝐴𝑥=𝑏, and Leizarowitz [39] established conditions for weak ergodicity of products, existence of optimal strategies for controlled Markov chains, and growth properties of certain linear nonautonomous differential equations based on a sequence (an infinite product) of stochastic matrices {𝐵_𝑚}^∞_𝑚=1(∏^∞_𝑚=1𝐵_𝑚). Also, as discussed in [2], the motivation of Definition 3.1 stems from a question about linear systems of difference equations and coupled matrix equations, under what conditions on {𝐵_𝑚}^∞_𝑚=1 of, for instance, the system 𝑥_𝑚=𝐵_𝑚𝑥_𝑚−1,𝑚=1,2,…, approach a finite nonzero limit whenever 𝑥₀≠0? A system with property linear asymptotic equilibrium if and only if 𝐵_𝑚 is invertible for every 𝑚≥1 and ∏^∞_𝑚=1𝐵_𝑚 invertibly convergent, but a system with the so-called least-squares linear asymptotic equilibrium if 𝐵_𝑚≠0 for every 𝑚≥1 and ∏^∞_𝑚=1𝐵_𝑚 converges.

Because of Theorem 3.4, we consider only infinite product of the form ∏^∞_𝑚=1(𝐼+𝐴_𝑚), where lim_𝑚→∞𝐴_𝑚=0. We will write 𝑃_𝑛=^𝑛_𝑚=1𝐼+𝐴_𝑚;𝑃=^∞_𝑚=1𝐼+𝐴_𝑚. (3.14) The following Theorem provides the convergence and invertibly convergence of the infinite product ∏^∞_𝑚=1(𝐼+𝐴_𝑚), and the proof here is omitted.

Theorem 3.6. The infinite product ∏^∞_𝑚=1(𝐼+𝐴_𝑚) converges (invertibly converges) if ∑^∞_𝑚=1‖𝐴_𝑚‖_𝑝<∞, for some 𝑝-norm ‖⋅‖.

The following theorem relates convergence of an infinite product to the asymptotic behavior of least-square solutions of a related system of difference equations.

Theorem 3.7. The infinite product ∏^∞_𝑚=1(𝐼+𝐴_𝑚) converges if and only if for some integer 𝑁≥1 the matrix difference equation 𝑋_𝑛+1=𝐼+𝐴_𝑛𝑋_𝑛∶𝑛≥𝑁 (3.15) has a least-square solution {𝑋_𝑛}^∞_𝑛=𝑁 such that lim_𝑛→∞𝑋_𝑛=𝐼. In this case, 𝑃_𝑛=𝑋_𝑛+1𝑋⁺_𝑁^𝑁−1_𝑚=1𝐼+𝐴_𝑚;𝑃=𝑋⁺_𝑁^𝑁−1_𝑚=1𝐼+𝐴_𝑚. (3.16)

Proof. Suppose that ∏^∞_𝑚=1(𝐼+𝐴_𝑚) converges. Choose 𝑁 so that 𝐼+𝐴_𝑚≠0 for 𝑚≥𝑁, and let 𝑄=∏^𝑛_𝑚=𝑁(𝐼+𝐴_𝑚). Define 𝑋_𝑛=^𝑛−1_𝑚=𝑁𝐼+𝐴_𝑚𝑄⁺∶𝑛≥𝑁. (3.17) Then {𝑋_𝑛}^∞_𝑛=𝑁 is a solution of (3.15) such that lim_𝑛→∞𝑋_𝑛=𝐼.
Conversely, suppose that (3.15) has a least-square solution {𝑋_𝑛}^∞_𝑛=𝑁 such that lim_𝑛→∞𝑋_𝑛=𝐼. Then𝑋_𝑛=^𝑛−1_𝑚=𝑁𝐼+𝐴_𝑚𝑋_𝑁∶𝑛≥𝑁, (3.18) where 𝑋_𝑁≠0. Therefore, ^𝑛−1_𝑚=𝑁𝐼+𝐴_𝑚=𝑋_𝑛𝑋⁺_𝑁∶𝑛≥𝑁. (3.19) Letting 𝑛→∞ shows that ^∞_𝑚=𝑁𝐼+𝐴_𝑚=𝑋⁺_𝑁 (3.20) which implies that ∏^∞_𝑚=1(𝐼+𝐴_𝑚) converges.
From (3.14) and (3.18), we get𝑃_𝑛=^𝑛_𝑚=1𝐼+𝐴_𝑚=𝑋_𝑛+1𝑋⁺_𝑁^𝑁−1_𝑚=1𝐼+𝐴_𝑚, (3.21) which proves the first expression in (3.16). From (3.14) and (3.20), we get 𝑃=^∞_𝑚=1𝐼+𝐴_𝑚=𝑋⁺_𝑁^𝑁−1_𝑚=1𝐼+𝐴_𝑚. (3.22) which proves the second expression in (3.16).

Remark 3.8. If the infinite product ∏^∞_𝑚=1(𝐼+𝐴_𝑚) is invertibly convergent in Theorem 3.7, then 𝑃_𝑛=𝑋_𝑛+1𝑋_𝑁^−1𝑁−1_𝑚=1𝐼+𝐴_𝑚;𝑃=𝑋_𝑁^−1𝑁−1_𝑚=1𝐼+𝐴_𝑚. (3.23)

Theorem 3.7 indicates the connection between convergence of an infinite product of 𝑘×𝑘 matrices ∏^∞_𝑚=1(𝐼+𝐴_𝑚) and the asymptotic properties of least-squares solutions of matrix difference equation as defined in (3.15). We can say that (3.15) has a least-squares linear asymptotic equilibrium if every least-squares solution of (3.15) for which 𝑋_𝑁≠0 that approaches 𝐼 as 𝑛→∞.

For example, Ding and Chen [5] and generalized later by Kılıçman and Al-Zhour [8] studied the convergence of least-square solutions to the coupled Sylvester matrix equations 𝐴𝑋+𝑌𝐵=𝐶;𝐷𝑋+𝑌𝐸=𝐹, (3.24) where 𝐴,𝐷∈𝑀_𝑚, 𝐵,𝐸∈𝑀_𝑛, 𝐶,𝐹∈𝑀_𝑚,𝑛 are given constant matrices, and 𝑋,𝑌∈𝑀_𝑚,𝑛 are the unknown matrices to be solved. If the coupled Sylvester matrix equation determined by (3.24) has a unique solution 𝑋 and 𝑌, then the following iterative solution 𝑋_𝑛+1 and 𝑌_𝑛+1 given by [5, 8]: 𝑋_𝑛+1=𝑋_𝑛+𝜇𝐺⁺𝐶−𝐴𝑋_𝑛−𝑌_𝑛𝐵𝐹−𝐷𝑋_𝑛−𝑌_𝑛𝐸;𝑌_𝑛+1=𝑌_𝑛+𝜇𝐶−𝐴𝑋_𝑛−𝑌_𝑛𝐵𝐹−𝐷𝑋_𝑛−𝑌_𝑛𝐸𝐻⁺, (3.25) where 𝐺=^𝐴_𝐷 and 𝐻=[^𝐵𝐸] are full column and full row-rank matrices, respectively; 1𝜇=1𝑚+𝑛or𝜇=𝜆_max𝐺(𝐺^∗𝐺)⁻¹𝐺^∗+𝜆_max𝐻^∗(𝐻𝐻^∗)⁻¹𝐻^∗. (3.26) converges to 𝑋 and 𝑌 for any finite initial values 𝑋₀ and 𝑌₀.

The convergence factor 𝜇 in (3.26) may not be the best and may be conservative. In fact, there exists a best 𝜇 such that the fast convergence rate of 𝑋_𝑘 to 𝑋 and 𝑌_𝑘 to 𝑌 can be obtained as in numerical examples given by Cao [21] and Kılıçman and Al-Zhour [8]. How to find the connections between convergence of an infinite products of 𝑘×𝑘 matrices and least-square solutions of coupled Sylvester matrix equation in (3.24) requires further research.

4. Numerical Examples

Here, we give some numerical example for computing outer inverse 𝐴⁽²⁾_𝑇,𝑆 and Moore-Penrose inverse 𝐴⁽²⁾_𝑇,𝑆 by applying sequences methods which are studied and derived in Section 2. Our results are obtained in this Section by choosing Frobenius norm (‖⋅‖₂) and using MATLAB software.

Example 4.1. Consider the matrix 𝐴=20.40.40.40200002000020000∈ℂ^5×4. (4.1) Let 𝑇=ℂ⁴,𝑒=[0,0,0,0,1]^𝑇∈ℂ⁵,𝑆=span{𝑒}.
Take 𝑋=𝐴₀=0.4000000.4000000.4000000.40∈ℂ^4×5. (4.2)

Here 𝑅(𝛽,𝑛)=|𝛽|𝑞^𝑛(1−𝑞)⁻¹‖𝑋‖₂‖𝐼−𝐴𝐴₀‖₂. Clearly rang(𝑋)=𝑇, null(𝑋)=𝑆, and rang(𝐴₀)⊂𝑇. By computing, we have 𝐴⁽²⁾_𝑇,𝑆=0.5−0.1−0.1−0.1000.5000000.5000000.50∈ℂ^4×5. (4.3) In order to satisfy 𝑞=min{‖𝐼−𝛽𝑋𝐴‖₂,‖𝐼−𝛽𝐴𝑋‖₂}<1, we get that 𝛽 should satisfy the following 0.63474563020816<𝛽<1.79243883581125.

From the iteration (2.5) in [24, Theorem  2.2], Let 𝐴∈ℂ^𝑚×𝑟, and 𝑇 and 𝑆 be given subspaces of ℂ^𝑚×𝑟 such that there exists 𝐴⁽²⁾_𝑇,𝑆. Then the sequence (𝐴_𝑛)_𝑛 in ℂ^𝑚×𝑟 defined in the following way: 𝑅_𝑛=𝑃_{𝐴(𝑇),𝑆}−𝑃_{𝐴(𝑇),𝑆}𝐴𝐴_𝑛,𝐴_𝑛+1=𝐴₀𝑅_𝑛+𝐴_𝑛𝑛=0,1,2,… (4.4) converges to 𝐴⁽²⁾_𝑇,𝑆 if and only if rang(𝐴₀)⊂𝑇 and 𝜌(𝑅₀)<1 (where 𝑅(𝑛)=(‖𝑅₀‖^𝑛+1/(1−‖𝑅₀‖))‖𝐴₀‖ and 𝑅₀=𝑃_{𝐴(𝑇),𝑆}−𝑃_{𝐴(𝑇),𝑆}𝐴𝐴₀).

In this case, if ‖𝑅₀‖=𝑞<1, then 𝐴⁽²⁾_𝑇,𝑆−𝐴_𝑛≤𝑞^𝑛+1𝐴1−𝑞⁰. (4.5)

Thus we have Tables 1 and 2 respectively, where 𝑅=𝛽𝑞𝛽,𝑛^𝑛‖𝑋1−𝑞‖𝐼−𝐴𝐴₀(𝑛)=𝑅,𝑅^0𝑛+11−𝑅₀𝐴₀. (4.6) Table 1 illustrates that 𝛽=1.25 is best value such that ‖𝐴⁽²⁾_𝑇,𝑆−𝐴_𝑛‖ reaches 2.020636405220133×10⁻¹⁶ iterating the least number of steps, the reason for which is that such a 𝛽 is calculating by using (2.38). Thus, for an appropriate 𝛽, the iteration is better than the iteration (4.4) (cf. Tables 1 and 2). And with respect to the error bound, the iterations for almost all are also better. Let us take the error bound smaller than 10⁻¹⁶; for instance, the number of steps of iterations in Table 1 is smaller than that of the iterations in Table 2. But, in practice, we consider also the quantity ‖𝐴_𝑛−𝐴_𝑛−1‖ in order to cease iteration since there exist such cases as 𝛽=1.25. For example, for ‖𝐴_𝑛−𝐴_𝑛−1‖<𝜇‖𝐴_𝑛‖, where 𝜇 is the machine precision, the iteration for 𝛽=1.25 only needs 3 steps. Therefore, in general, the iteration of (2.36) is better than the iteration (4.4) for an appropriate 𝛽. Note that the iterations in both Tables 1 and 2 indicate a fast convergence for the quantity ‖𝐴⁽²⁾_𝑇,𝑆−𝐴_𝑛‖ more than the quantity 𝑅(𝛽,𝑛) in Table 1 and the quantity 𝑅(𝑛) in Table 2 since each of 𝑅(𝛽,𝑛) and 𝑅(𝑛) is an upper bound of the quantity ‖𝐴⁽²⁾_𝑇,𝑆−𝐴_𝑛‖, and, to find the best or least upper bound for the quantity, ‖𝐴⁽²⁾_𝑇,𝑆−𝐴_𝑛‖ requires further research.

Table 1 

Results for Example 4.1 using the iterations (2.36) and (2.37) for 𝛽=1.25.

Table 2 

Results for Example 4.1 using the iterations (4.4) and (4.5).

Example 4.2. Consider the matrix 𝐴=22102371471008−113−1−113−3−213−24981−2491−75−12−6651450−22∈ℂ^8×5. (4.7) Then by computing we have 𝐴⁺=2.1129808×10⁻²4.6153846×10⁻³−2.1073718×10⁻³7.6041667×10⁻³3.8060897×10⁻³9.3108974×10⁻³2.2115385×10⁻³2.0528846×10⁻³−2.0833333×10⁻³1.0016026×10⁻³−1.1097756×10⁻²2.7403846×10⁻²−3.8862179×10⁻³−2.7604167×10⁻²4.2067308×10⁻³−7.9166667×10⁻³−5.0000000×10⁻³3.3750000×10⁻²−5.4166667×10⁻³1.0416667×10⁻³5.5128205×10⁻³9.8076923×10⁻³−8.9743590×10⁻⁴−5.0000000×10⁻³3.2051282×10⁻³1.4318910×10⁻²−2.5961538×10⁻³−2.0136218×10⁻²1.2812500×10⁻³−6.2099359×10⁻³4.8958333×10⁻³−1.5000000×10⁻²1.5312500×10⁻²1.2395833×10⁻²2.4604166×10⁻³1.5064103×10⁻³7.4038462×10⁻³−1.6987179×10⁻³−5.0000000×10⁻³1.6025641×10^−3𝑇∈ℂ^5×8. (4.8)

Thus, (see Tables 3 and 4).

Table 3 

Results for Example 4.2 using an acceleration iteration (2.34).

Table 4 

Results for Example 4.2 using the iteration (2.18).

Example 4.3. We generate a random matrix, 𝐴∈ℂ^100×80 by using MATLAB, and then we obtain the results as in Tables 5 and 6.

Table 5 

Results for a random matrix 𝐴∈ℂ^100×80 using an acceleration iteration (2.34).

Table 6 

Results for a random matrix 𝐴∈ℂ^100×80 using the iteration (2.18).

Note that from Tables 3, 4, 5, and 6, it is clear that the quantities ‖𝑋𝐴𝑋−𝑋‖, and ‖𝐴𝑋𝐴−𝐴‖, ‖(𝐴𝑋)^∗−𝐴𝑋‖, ‖(𝑋𝐴)^∗−𝑋𝐴‖ are becoming smaller and smaller and goes to zero as 𝑛 increases in both iterations (2.34) and (2.18). We can also conclude that both iterations almost have same fast of convergence when the dimension of any arbitrary matrix 𝐴 is not so large, but the acceleration iteration (2.34) is better more than the iteration (2.18) when the dimension of any arbitrary matrix 𝐴 is so large with an appropriate acceleration parameter 𝛼_𝑛+1∈[1,2].

5. Concluding Remarks

In this paper, we have studied some matrix sequences convergence to the Moore-Penrose inverse 𝐴⁺ and outer inverse 𝐴⁽²⁾_𝑇,𝑆 of an arbitrary matrix 𝐴∈𝑀_𝑚,𝑛. The key to derive matrix sequences which are convergent to weighted Drazin and weighted Moore-Penrose inverses is the Lemma 2.2. Some sufficient conditions for infinite products ∏^∞_𝑚=1𝐵_𝑚 and ∏^∞_𝑚=1⊗𝐵_𝑚 of 𝑘×𝑘 matrices are also derived. In our opinion, it is worth establishing some connections between convergence of an infinite products of 𝑘×𝑘 matrices and least-square solutions of such linear singular systems as well as the singular coupled matrix equations.

6. Acknowledgments

The authors express their sincere thanks to the referee(s) for careful reading of the manuscript and several helpful suggestions. The authors also gratefully acknowledge that this research was partially supported by Ministry of Science, Technology and Innovations (MOSTI), Malaysia under the e-Science Grant 06-01-04-SF1050.