Partially doubly symmetric solutions of general Sylvester matrix equations

Abstract

The study of linear matrix equations is extremely important in many scientific fields such as control systems and stability analysis. In this work, we aim to design the Hestenes-Stiefel (HS) version of biconjugate residual (Bi-CR) algorithm for computing the (least Frobenius norm) partially doubly symmetric solution $X$ of the general Sylvester matrix equations $\sum_{j = 1}^{f} A_{ij} X B_{ij} = M_{i}$ for $i = 1, 2, . . ., g$ . We show that the proposed algorithm converges in a finite number of iterations. Finally, numerical results compare the proposed algorithm to alternative algorithms.

Keywords

Hestenes-Stiefel bi-conjugate residual (Bi-CR) algorithm general Sylvester matrix equations partially doubly symmetric solution CGNR algorithm CGNE algorithm

Introduction

It is well known that linear matrix equations play crucial roles in engineering and applied mathematics (Ke and Ma, 2018; Simoncini, 2016). For example, the Stein matrix equation

X + AX A^{T} = B,

(1.1)

and the Lyapunov matrix equation

X A^{T} + A^{T} X = B,

(1.2)

occur in the study of discrete-time dynamical systems and the analysis of continuous-time dynamical systems (DeAlba and Johnson, 1995). The generalized Lyapunov matrix equation

A_{0} X + {XA}_{0}^{T} + A_{1} {XA}_{1}^{T} = - Q,

(1.3)

is used to check the stability of an Itô stochastic system (Zhang et al., 2018)

dx (t) = A_{0} x (t) dt + A_{1} x (t) d ω (t) .

(1.4)

The generalized Sylvester matrix equation

A_{1} {XB}_{1}^{T} + A_{2} {XB}_{2}^{T} + . . . + A_{8} {XB}_{8}^{T} = C,

(1.5)

appears in the application of the Sinc-Galerkin method to the Poisson problem in the curvilinear coordinate domain (Yamamoto, 2006). In the last 15 years, growing attention is seen in the literature to the presentation of efficient methods for solving various linear matrix equations (Hajarian, 2017; Larin, 2009a; Liao et al., 2005; Wang and He, 2014; Xie et al., 2014; Yuan and Liao, 2011). Peng et al. (2007) studied the necessary and sufficient conditions for the solvability of matrix equation

A^{H} XB = C,

over the reflexive or anti-reflexive matrices. In Wang and Zhang (2008); Wang et al. (2009); and Wang and Li (2009), the least-norm and reflexive re-nonnegative definite solutions of the quaternion matrix equations were studied. Zhou et al. (2010) investigated the gradient-based iterative (GI) method to obtain the solution X of the general class of linear matrix equations

A_{1} {XB}_{1}^{T} + A_{2} {XB}_{2}^{T} + . . . + A_{p} {XB}_{p}^{T} = C .

(1.6)

Based on the Bass relation, Larin (2009b) proposed algorithms for solving the Lyapunov and Sylvester matrix equations. In Xiong and Qin (2011), the common Re-nonnegative definite and Re-positive definite solutions of the matrix equations

AX = C and XB = D,

(1.7)

were given. Xie et al. (2010) introduced a gradient based and a least squares based iterative algorithms to solve the generalized Sylvester-transpose matrix equation

AXB + C X^{T} D = F .

(1.8)

In Huang and Ma (2014), a modification of the conjugate gradient (CG) method was investigated to find the solutions of the generalized coupled Sylvester-transpose matrix equations

A_{1} X B_{1} + C_{1} Y^{T} D_{1} = F_{1}, A_{2} Y B_{2} + C_{2} X^{T} D_{2} = F_{2} .

(1.9)

Li et al. (2010) proposed the conjugate gradient least squares method (CGLS) for minimizing

| | AX A^{T} - B | |,

where matrix $X$ is a symmetric centrosymmetric with a specified central submatrix. Dehghan and Hajarian (2011, 2010a,b, 2012) introduced CG algorithms for computing the reflexive and anti-reflexive solutions of the generalized coupled Sylvester matrix equations. Based on the hierarchical identification principle, gradient based iterative algorithms were proposed for finding the solutions of several (coupled) linear matrix equations (Ding and Chen, 2005b, 2006, 2005a; Ding and Zhang, 2014. In Hajarian (2016a), the generalized conjugate direction (CD) algorithm was derived for computing the symmetric solution group of the general coupled matrix equations. Hajarian (2013b,a, 2014a,b, 2015) developed the modern Krylov space methods for finding solutions of linear matrix equations.

In this paper, we develop the Hestenes-Stiefel (HS) version of biconjugate residual (Bi-CR) algorithm to propose a matrix algorithm for solving the general Sylvester matrix equations

\sum_{j = 1}^{f} A_{ij} X B_{ij} = M_{i}, i = 1, 2, . . ., g .

(1.10)

over the partially doubly symmetric matrix $X$ . The HS version of Bi-CR algorithm was originally established in Vespucci and Broyden (2001) to solve the nonsymmetric linear system of equations $Ax = b$ . It is worth pointing out that the general Sylvester matrix equations (1.10) contain the Lyapunov and Sylvester matrix equations as special cases.

The layout of the paper is as follows. In Section 2, we cover the basic concepts, notation and definitions needed throughout the paper. In Section 3, first we briefly introduce the HS version of Bi-CR algorithm for the linear system $Ax = b$ and then we extend this algorithm for solving the general Sylvester matrix equations (1.10) over the partially doubly symmetric matrix $X$ . Also the convergence of the proposed algorithm is discussed. The results of numerical examples are presented in Section 4, illustrating the performance of the algorithm. Section 5 ends this paper with brief conclusions.

Preliminary

First, let us introduce the following notation Hajarian (2018d, 2019). Given a matrix $A \in R^{m \times n}$ , let $tr (A)$ , $| | A | |$ and $R (A)$ , respectively, denote the trace, the Frobenius norm and the column space of $A$ . The vectorization of the matrix $A \in R^{m \times n}$ , denoted by $vec (A)$ , is defined by $vec (A) = (a_{1}^{T}, a_{2}^{T}, . . ., a_{n}^{T})^{T}$ where $a_{k}$ is the $k$ -th column of $A$ . Let $A \otimes B$ be the Kronecker product of matrices $A$ and $B$ . Denote by $I_{n} = (e_{1}, e_{2}, . . ., e_{n})$ and $J_{n} = (e_{n}, e_{n - 1}, . . ., e_{1})$ , respectively, the $n \times n$ unit matrix and $n \times n$ flip matrix with ones on the secondary diagonal and zeros elsewhere. Let $S_{p, n} = {α = (α_{1}, α_{2}, . . ., α_{p}, n + 1 - α_{p}, . . ., n + 1 - α_{2}, n + 1 - α_{1}) : 1 \leq α_{1} < α_{2} < . . . < α_{p} < \frac{n + 1}{2}}$ denote strictly increasing sequences of $2 p$ elements from $1, 2, . . ., n$ . For $α = (α_{1}, α_{2}, . . ., α_{p}, n + 1 - α_{p}, . . ., n + 1 - α_{2}, n + 1 - α_{1}) \in S_{p, n}$ , we consider $E_{α} = (e_{α_{1}}, e_{α_{2}}, . . ., e_{α_{p}}, e_{n + 1 - α_{p}}, . . ., e_{n + 1 - α_{2}}, e_{n + 1 - α_{1}}) \in R^{n \times 2 p}$ . We use $A [α | β]$ to represent the $p \times q$ submatrix of $A$ determined by rows indexed by $α$ and columns indexed by $β = (β_{1}, β_{2}, . . ., β_{q}, n + 1 - β_{q}, . . ., n + 1 - β_{2}, n + 1 - β_{1}) \in S_{q, n}$ . We use $A [\bar{α} | \bar{β}]$ to denote the $(m - 2 p) \times (n - 2 q)$ submatrix of $A$ determined by deleting rows indexed by $α$ and columns indexed by $β$ . A real $n \times n$ matrix $A = (a_{i, j})$ is said to be doubly symmetric if $a_{i, j} = a_{j, i}$ and $a_{i, j} = a_{n - j + 1, n - i + 1}$ for $1 \leq i, j \leq n$ . Let $BS R^{n \times n}$ be the set of $n \times n$ doubly symmetric matrices. It can be clearly seen that a matrix $X \in BS R^{n \times n}$ iff $X = X^{T} = J_{n} X J_{n}$ .

Definition 1: Given $α = (α_{1}, α_{2}, . . ., α_{p}, n + 1 - α_{p}, . . ., n + 1 - α_{2}, n + 1 - α_{1}) \in S_{p, n}$ and $X_{p} \in R^{2 p \times 2 p}$ . Let ${\hat{X}}_{p}$ be the matrix satisfying ${\hat{X}}_{p} [α | α] = X_{p}$ and zeros elsewhere. The set of $n \times n$ partially doubly symmetric matrices, also called the partially bisymmetric matrices, is defined as:

T_{α} = {Y \in R^{n \times n} | Y [α | α] = X_{p}, (Y - {\hat{X}}_{p}) \in BS R^{n \times n}} .

Doubly symmetric matrices containing the symmetric Toeplitz matrices, centrosymmetric Hankel matrices and symmetric circulant matrices as special types play a fundamental role in many important application fields like the pattern recognition, signal reconstruction, communication theory, differential equations, magic squares and harmonic differential quadrature (Andrew, 1998; Cantoni and Butler, 1976; Datta and Morgera, 1989; Good, 1970; Hanna and Mansoori, 2003; Mattingly, 2000). The present paper is concerned with the following problem related to the general Sylvester matrix equations (1.10) over the partially doubly symmetric matrices:

Problem 1: Given $A_{ij} \in R^{s_{i} \times n}$ , $B_{ij} \in R^{n \times t_{i}}$ , $M_{i} \in R^{s_{i} \times t_{i}}$ , for $j = 1, 2, . . ., f$ and $i = 1, 2, . . ., g$ , $α = (α_{1}, α_{2}, . . ., α_{p}, n + 1 - α_{p}, . . ., n + 1 - α_{2}, n + 1 - α_{1}) \in S_{p, n}$ and $X_{p} \in R^{2 p \times 2 p}$ . Find the partially doubly symmetric matrix $X \in T_{α}$ such that (1.10) holds.

Since doubly symmetric matrices are complicated and we do not generate them directly, Problem 1 is very difficult to solve in practice. In order to overcome this complication, first we divide the set of $n \times n$ partially doubly symmetric matrices into two sets with simple structures and second by using these sets we transform Problem 1 into a simpler problem and then solve it. Setting $S_{α} = U_{α} + V_{α}$ where

\begin{matrix} U_{α} = {Y \in R^{n \times n} | Y [α | α] = Xp, Y [\bar{α} | \bar{α}] \\ = 0, Y [α | \bar{α}] = 0, Y [\bar{α} | α] = 0}, \end{matrix}

and

V_{α} = {Y \in BS R^{n \times n} | Y [α | α] = 0},

gives us $Y \in S_{α} \Rightarrow Y \in T_{α}$ . Based on the above discussion, instead of Problem 1 we study the following problem.

Problem 2: Given $A_{ij} \in R^{s_{i} \times n}$ , $B_{ij} \in R^{n \times t_{i}}$ , $M_{i} \in R^{s_{i} \times t_{i}}$ , for $j = 1, 2, . . ., f$ and $i = 1, 2, . . ., g$ , $α = (α_{1}, α_{2}, . . ., α_{p}, n + 1 - α_{p}, . . ., n + 1 - α_{2}, n + 1 - α_{1}) \in S_{p, n}$ and $X_{p} \in R^{2 p \times 2 p}$ . Find $Z \in V_{α}$ such that

$\sum_{j = 1}^{f} A_{ij} Z B_{ij} = N_{i}, i = 1, 2, . . ., g .$

where $N_{i} = M_{i} - \sum_{j = 1}^{f} A_{ij} {\hat{X}}_{p} B_{ij}$ for $i = 1, 2, . . ., g$ .

Remark 1: It is straightforward to verify that $Z^{*}$ is the solution of Problem 2 iff $Z^{*} + {\hat{X}}_{p}$ is the solution of Problem 1.

Algorithm and its convergence analysis

The aim of this section is to generalize the HS version of Bi-CR algorithm developed in Vespucci and Broyden (2001) for solving Problem 2. Vespucci and Broyden (2001) established various versions of Bi-CR algorithm that approximately solve the nonsymmetric linear system of equations $Ax = b$ . The convergence theories of the proposed algorithms were not studied in Vespucci and Broyden (2001) and only numerical results were presented to show the efficiency and performance of the algorithms. Recently, Hajarian (2016b, 2018a,b,c,d; 2019) generalized the Bi-CR algorithms for solving Sylvester matrix equations and constrained quadratic inverse eigenvalue problems. One of the algorithms proposed in Vespucci and Broyden (2001) is the HS version of Bi-CR algorithm, which can be summarized as follows.

HS version of Bi-CR algorithm

Initial values: $x (1)$ and $s (1)$ arbitrary, $r (1) = Ax (1) - b$ ,

$u (1) = s (1)$ , $v (1) = r (1)$ , $w (1) = Au (1)$ , and $y (1) = A^{T} v (1)$ .

Recursions

α (k) = \frac{r {(k)}^{T} As (k)}{w {(k)}^{T} w (k)}, x (k + 1) = x (k) - α (k) u (k),

r (k + 1) = r (k) - α (k) w (k),

β (k) = \frac{r {(k)}^{T} As (k)}{y {(k)}^{T} y (k)}, s (k + 1) = s (k) - β (k) y (k),

γ (k) = \frac{r {(k + 1)}^{T} As (k + 1)}{r {(k)}^{T} As (k)}, u (k + 1) = s (k + 1) + γ (k) u (k),

v (k + 1) = r (k + 1) + γ (k) v (k),

w (k + 1) = Au (k + 1),

y (k + 1) = A^{T} v (k + 1) .

In order to develop the HS version of Bi-CR algorithm for solving Problem 2, we first need a linear system that is equivalent to Problem 2. In view of $V_{α}$ , it can be readily seen that the solvability of Problem 2 is equivalent to the solvability of the following general Sylvester matrix equations

{\begin{matrix} \sum_{j = 1}^{f} A_{ij} Z B_{ij} = N_{i}, i = 1, 2, . . ., g, \\ \sum_{j = 1}^{f} B_{ij}^{T} {ZA}_{ij}^{T} = N_{i}^{T}, i = 1, 2, . . ., g, \\ \sum_{j = 1}^{f} A_{ij} J_{n} Z J_{n} B_{ij} = N_{i}, i = 1, 2, . . ., g, \\ \sum_{j = 1}^{f} B_{ij}^{T} J_{n} Z J_{n} A_{ij}^{T} = N_{i}^{T}, i = 1, 2, . . ., g, \\ E_{α}^{T} Z E_{α} = 0, \\ E_{α}^{T} J_{n} Z J_{n} E_{α} = 0 . \end{matrix}

(3.1)

By means of Kronecker product and vectorization, we rewrite the general Sylvester matrix equations (3.1) in an equivalent form

\underset{A}{\underset{︸}{(\begin{matrix} \sum_{j = 1}^{f} B_{1 j}^{T} \otimes A_{1 j} \\ ⋮ ⋮ ⋮ \\ \sum_{j = 1}^{f} B_{gj}^{T} \otimes A_{gj} \\ \sum_{j = 1}^{f} A_{1 j} \otimes B_{1 j}^{T} \\ ⋮ ⋮ ⋮ \\ \sum_{j = 1}^{f} A_{gj} \otimes B_{gj}^{T} \\ \sum_{j = 1}^{f} B_{1 j}^{T} J_{n} \otimes A_{1 j} J_{n} \\ ⋮ ⋮ ⋮ \\ \sum_{j = 1}^{f} B_{gj}^{T} J_{n} \otimes A_{gj} J_{n} \\ \sum_{j = 1}^{f} A_{1 j} J_{n} \otimes B_{1 j}^{T} J_{n} \\ ⋮ ⋮ ⋮ \\ \sum_{j = 1}^{f} A_{gj} J_{n} \otimes B_{gj}^{T} J_{n} \\ E_{α}^{T} \otimes E_{α}^{T} \\ E_{α}^{T} J_{n} \otimes E_{α}^{T} J_{n} \end{matrix})}} \underset{x}{\underset{︸}{vec (Z)}} = \underset{b}{\underset{︸}{(\begin{matrix} vec (N_{1}) \\ ⋮ \\ vec (N_{g}) \\ vec (N_{1}^{T}) \\ ⋮ \\ vec (N_{g}^{T}) \\ vec (N_{1}) \\ ⋮ \\ vec (N_{g}) \\ vec (N_{1}^{T}) \\ ⋮ \\ vec (N_{g}^{T}) \\ 0 \\ 0 \end{matrix})}} .

(3.2)

This leads to the following proposition.

Proposition 1: Problem 2 has a unique solution if and only if $rank ((A, b)) = rank (A)$ and $A$ has a full column rank.

It must be pointed out that in this step we can apply the above HS version of Bi-CR algorithm and some other iterative methods (Chronopoulos, 1991, 1994; Chronopoulos and Kincaid, 2001; Chronopoulos and Kucherov, 2010; Saad, 2003) for solving the linear system (3.2) but the size of the coefficient matrix of this system is very large in comparison with the size of the coefficient matrices of (1.10). Therefore we announce the matrix form of HS version of Bi-CR algorithm to solve (1.10). For this purpose, we substitute the parameters $A$ and $b$ of linear system (3.2) into the HS version of Bi-CR algorithm. With the above discussion as background, we are now ready to propose the extended HS version of Bi-CR algorithm for finding the solution of Problem 2 by the following.

Algorithm 1: Given the arbitrary doubly symmetric matrices $Z (1) \in BS R^{n \times n}$ and $S (1) \in BS R^{n \times n}$ ;

Compute

R_{i} (1) = \sum_{j = 1}^{f} A_{ij} Z (1) B_{ij} - N_{i}, i = 1, 2, . . ., g,

\hat{R} (1) = E_{α}^{T} Z (1) E_{α},

U (1) = S (1), V_{i} (1) = R_{i} (1), i = 1, 2, . . ., g, \hat{V} (1) = \hat{R} (1),

W_{i} (1) = \sum_{j = 1}^{f} A_{ij} U (1) B_{ij}, i = 1, 2, . . ., g,

\hat{W} (1) = E_{α}^{T} U (1) E_{α},

\begin{matrix} Y (1) = \frac{1}{4} [\sum_{i = 1}^{g} \sum_{j = 1}^{f} (A_{ij}^{T} V_{i} (1) B_{ij}^{T} + (A_{ij}^{T} V_{i} (1) B_{ij}^{T})^{T} \\ + J_{n} (A_{ij}^{T} V_{i} (1) B_{ij}^{T} + (A_{ij}^{T} V_{i} (1) B_{ij}^{T})^{T}) J_{n}) \end{matrix}

\begin{matrix} + E_{α} \hat{V} (1) E_{α}^{T} + (E_{α} \hat{V} (1) E_{α}^{T})^{T} \\ + J_{n} (E_{α} \hat{V} (1) E_{α}^{T} + (E_{α} \hat{V} (1) E_{α}^{T})^{T}) J_{n}]; \end{matrix}

For $k = 1, 2, . . .,$ until convergence, do

α (k) = \frac{\sum_{i = 1}^{g} tr (R_{i} {(k)}^{T} (\sum_{j = 1}^{f} A_{ij} S (k) B_{ij})) + tr (\hat{R} {(k)}^{T} (E_{α}^{T} S (k) E_{α}))}{\sum_{i = 1}^{g} | | W_{i} (k {) | |}^{2} + | | \hat{W} (k {) | |}^{2}},

Z (k + 1) = Z (k) - α (k) U (k),

R_{i} (k + 1) = R_{i} (k) - α (k) W_{i} (k), i = 1, 2, . . ., g, \hat{R} (k + 1) = \hat{R} (k) - α (k) \hat{W} (k),

β (k) = \frac{\sum_{i = 1}^{g} tr (R_{i} {(k)}^{T} (\sum_{j = 1}^{f} A_{ij} S (k) B_{ij})) + tr (\hat{R} {(k)}^{T} (E_{α}^{T} S (k) E_{α}))}{| | Y (k {) | |}^{2}},

S (k + 1) = S (k) - β (k) Y (k),

\begin{matrix} γ (k) = \\ \frac{\sum_{i = 1}^{g} tr (R_{i} {(k + 1)}^{T} (\sum_{j = 1}^{f} A_{ij} S (k + 1) B_{ij})) + tr (\hat{R} {(k + 1)}^{T} (E_{α}^{T} S (k + 1) E_{α}))}{\sum_{i = 1}^{g} tr (R_{i} {(k)}^{T} (\sum_{j = 1}^{f} A_{ij} S (k) B_{ij})) + tr (\hat{R} {(k)}^{T} (E_{α}^{T} S (k) E_{α}))} \end{matrix}

U (k + 1) = S (k + 1) + γ (k) U (k),

\begin{matrix} V_{i} (k + 1) = R_{i} (k + 1) + γ (k) V_{i} (k), \\ i = 1, 2, . . ., g, \hat{V} (k + 1) = \hat{R} (k + 1) + γ (k) \hat{V} (k), \end{matrix}

W_{i} (k + 1) = \sum_{j = 1}^{f} A_{ij} U (k + 1) B_{ij}, i = 1, 2, . . ., g,

\hat{W} (k + 1) = E_{α}^{T} U (k + 1) E_{α},

\begin{matrix} Y (k + 1) = \frac{1}{4} [\sum_{i = 1}^{g} \sum_{j = 1}^{f} (A_{ij}^{T} V_{i} (k + 1) B_{ij}^{T} + (A_{ij}^{T} V_{i} (k + 1) B_{ij}^{T})^{T} \\ + J_{n} (A_{ij}^{T} V_{i} (k + 1) B_{ij}^{T} + (A_{ij}^{T} V_{i} (k + 1) B_{ij}^{T})^{T}) J_{n}) \end{matrix}

\begin{matrix} + E_{α} \hat{V} (k + 1) E_{α}^{T} + (E_{α} \hat{V} (k + 1) E_{α}^{T})^{T} \\ + J_{n} (E_{α} \hat{V} (k + 1) E_{α}^{T} + (E_{α} \hat{V} (k + 1) E_{α}^{T})^{T}) J_{n}] . \end{matrix}

Stopping criterion: As a stopping criterion, we choose to satisfy $\sqrt{\sum_{i = 1}^{g} | | R_{i} (k + {1) | |}^{2} + | | \hat{R} (k + {1) | |}^{2}} = 0$ . So our convergence criterion for Algorithm 1 becomes

\sqrt{\sum_{i = 1}^{g} | | R_{i} (k + {1) | |}^{2} + | | \hat{R} (k + {1) | |}^{2}} \leq tol,

where $tol$ is a chosen fixed threshold.

Remark 2: From Algorithm 1, we can observe that

γ (k) = 1 - \frac{\sum_{i = 1}^{g} tr (R_{i} {(k)}^{T} (\sum_{j = 1}^{f} (A_{ij} Y (k) B_{ij}))) + tr (\hat{R} {(k)}^{T} (E_{α}^{T} Y (k) E_{α}))}{| | Y (k {) | |}^{2}}

- \frac{\sum_{i = 1}^{g} tr (W_{i} {(k)}^{T} (\sum_{j = 1}^{f} (A_{ij} Y (k) B_{ij}))) + tr (\hat{W} {(k)}^{T} (E_{α}^{T} Y (k) E_{α}))}{\sum_{i = 1}^{g} | | W_{i} (k {) | |}^{2} + | | \hat{W} (k {) | |}^{2}}

+ [\frac{\sum_{i = 1}^{g} tr (R_{i} {(k)}^{T} (\sum_{j = 1}^{f} (A_{ij} S (k) B_{ij}))) + tr (\hat{R} {(k)}^{T} (E_{α}^{T} S (k) E_{α}))}{\sum_{i = 1}^{g} | | W_{i} (k {) | |}^{2} + | | \hat{W} (k {) | |}^{2}}

\times \frac{\sum_{i = 1}^{g} tr (W_{i} {(k)}^{T} (\sum_{j = 1}^{f} (A_{ij} Y (k) B_{ij}))) + tr (\hat{W} {(k)}^{T} (E_{α}^{T} Y (k) E_{α}))}{| | Y (k {) | |}^{2}}],

W_{i} (k + 1) = \sum_{j = 1}^{f} A_{ij} S (k + 1) B_{ij} + γ (k) W_{i} (k), i = 1, 2, . . ., g,

\hat{W} (k + 1) = E_{α}^{T} S (k + 1) E_{α} + γ (k) \hat{W} (k),

\begin{matrix} Y (k + 1) = \frac{1}{4} [\sum_{i = 1}^{g} \sum_{j = 1}^{f} (A_{ij}^{T} R_{i} (k + 1) B_{ij}^{T} \\ + (A_{ij}^{T} R_{i} (k + 1) B_{ij}^{T})^{T} + J_{n} (A_{ij}^{T} R_{i} (k + 1) B_{ij}^{T} \\ + (A_{ij}^{T} R_{i} (k + 1) B_{ij}^{T})^{T}) J_{n}) \end{matrix}

\begin{matrix} + E_{α} \hat{R} (k + 1) E_{α}^{T} + (E_{α} \hat{R} (k + 1) E_{α}^{T})^{T} + J_{n} (E_{α} \hat{R} (k + 1) E_{α}^{T} \\ + (E_{α} \hat{R} (k + 1) E_{α}^{T})^{T}) J_{n}] + γ (k) Y (k) . \end{matrix}

Now we are going to analyze the convergence of Algorithm 1. First, we deal with some properties of the sequences generated by Algorithm 1, which will be needed to prove the convergence result of this algorithm.

Lemma 1: If there exists a positive integer number $r$ such that $\sum_{i = 1}^{g} | | W_{i} (k) | |^{2} + | | \hat{W} (k) | |^{2} \neq 0$ and $\sum_{i = 1}^{g} | | R_{i} (k) | |^{2} + | | \hat{R} (k) | |^{2} \neq 0$ for all $k = 1, 2, . . ., r$ then

\begin{matrix} \sum_{i = 1}^{g} tr (R_{i} (v)^{T} W_{i} (u)) + tr (\hat{R} (v)^{T} \hat{W} (u)) = 0, \\ for u, v = 1, 2, . . ., r, u < v, \end{matrix}

(3.3)

tr (S (v)^{T} Y (u)) = 0, for u, v = 1, 2, . . ., r, u < v,

(3.4)

tr (Y (v)^{T} Y (u)) = 0, for u, v = 1, 2, . . ., r, u \neq v,

(3.5)

\begin{matrix} \sum_{i = 1}^{g} tr (W_{i} (v)^{T} W_{i} (u)) + tr (\hat{W} (v)^{T} \hat{W} (u)) = 0, \\ for u, v = 1, 2, . . ., r, u \neq v . \end{matrix}

(3.6)

The proof of Lemma 1 is long and can be found in Appendix A.

Our next theorem shows that Algorithm 1 can compute the solution of Problem 2 within a finite number of iterations in the absence of round-off errors.

Theorem 1: Let Problem 2 be consistent. Algorithm 1 obtains the solution of Problem 2 within a finite number of iterations in the absence of round-off errors.

Proof: In view of (3.3) and (3.6), it is not difficult to see that Algorithm 1 converges to the solution of Problem 2 within a finite number of iterations in the absence of round-off errors.

In what follows, we show that the least Frobenius norm solution of Problem 2 can be obtained by Algorithm 1 when the appropriate initial matrices are chosen.

Theorem 2: Let Problem 2 be consistent. If we take the initial matrix

\begin{matrix} Z (1) = \sum_{i = 1}^{g} \sum_{j = 1}^{f} (A_{ij}^{T} K_{i} (1) B_{ij}^{T} + (A_{ij}^{T} K_{i} (1) B_{ij}^{T})^{T} \\ + J_{n} (A_{ij}^{T} K_{i} (1) B_{ij}^{T} + (A_{ij}^{T} K_{i} (1) B_{ij}^{T})^{T}) J_{n}) \\ + E_{α} \hat{K} (1) E_{α}^{T} + (E_{α} \hat{K} (1) E_{α}^{T})^{T} + J_{n} (E_{α} \hat{K} (1) E_{α}^{T} \\ + (E_{α} \hat{K} (1) E_{α}^{T})^{T}) J_{n}, \end{matrix}

(3.7)

and matrix

\begin{matrix} S (1) = \sum_{i = 1}^{g} \sum_{j = 1}^{f} (A_{ij}^{T} L_{i} (1) B_{ij}^{T} + (A_{ij}^{T} L_{i} (1) B_{ij}^{T})^{T} \\ + J_{n} (A_{ij}^{T} L_{i} (1) B_{ij}^{T} + (A_{ij}^{T} L_{i} (1) B_{ij}^{T})^{T}) J_{n}) \\ + E_{α} \hat{L} (1) E_{α}^{T} + (E_{α} \hat{L} (1) E_{α}^{T})^{T} + J_{n} (E_{α} \hat{L} (1) E_{α}^{T} \\ + (E_{α} \hat{L} (1) E_{α}^{T})^{T}) J_{n}, \end{matrix}

(3.8)

where $K_{i} (1), L_{i} (1) \in R^{s_{i} \times t_{i}}$ and $\hat{K} (1), \hat{L} (1) \in R^{n \times n}$ for $i = 1, 2, . . ., g$ are arbitrary, then the solution $Z^{*}$ generated by Algorithm 1 is the least Frobenius norm solution of Problem 2.

Proof: By Algorithm 1, (3.7) and (3.8), it can be shown that there exist the matrices $K_{i} (k) \in R^{s_{i} \times t_{i}}$ and $\hat{K} (k) \in R^{n \times n}$ for $i = 1, 2, . . ., g$ such that

\begin{matrix} Z (k) = \sum_{i = 1}^{g} \sum_{j = 1}^{f} (A_{ij}^{T} K_{i} (k) B_{ij}^{T} + (A_{ij}^{T} K_{i} (k) B_{ij}^{T})^{T} \\ + J_{n} (A_{ij}^{T} K_{i} (k) B_{ij}^{T} + (A_{ij}^{T} K_{i} (k) B_{ij}^{T})^{T}) J_{n}) \\ + E_{α} \hat{K} (k) E_{α}^{T} + (E_{α} \hat{K} (k) E_{α}^{T})^{T} + J_{n} (E_{α} \hat{K} (k) E_{α}^{T} \\ + (E_{α} \hat{K} (k) E_{α}^{T})^{T}) J_{n} . \end{matrix}

(3.9)

Now we have using (3.9) and (3.2) that

vec (Z (k)) \in R (\begin{matrix} {(\begin{matrix} \sum_{j = 1}^{f} B_{1 j}^{T} \otimes A_{1 j} \\ ⋮ ⋮ ⋮ \\ \sum_{j = 1}^{f} B_{gj}^{T} \otimes A_{gj} \\ \sum_{j = 1}^{f} A_{1 j} \otimes B_{1 j}^{T} \\ ⋮ ⋮ ⋮ \\ \sum_{j = 1}^{f} A_{gj} \otimes B_{gj}^{T} \\ \sum_{j = 1}^{f} B_{1 j}^{T} J_{n} \otimes A_{1 j} J_{n} \\ ⋮ ⋮ ⋮ \\ \sum_{j = 1}^{f} B_{gj}^{T} J_{n} \otimes A_{gj} J_{n} \\ \sum_{j = 1}^{f} A_{1 j} J_{n} \otimes B_{1 j}^{T} J_{n} \\ ⋮ ⋮ ⋮ \\ \sum_{j = 1}^{f} A_{gj} J_{n} \otimes B_{gj}^{T} J_{n} \\ E_{α}^{T} \otimes E_{α}^{T} \\ E_{α}^{T} J_{n} \otimes E_{α}^{T} J_{n} \end{matrix})}^{T} \end{matrix}) \in R (A^{T}) .

This implies that the generated solution by Algorithm 1 is the least Frobenius norm solution of Problem 2.

Remark 3: To solve Problem 2, we can similarity obtain the generalization of other iterative methods which were proposed for computing the solution of nonsymmetric linear systems $Ax = b$ .

Here, we generalize the conjugate gradient normal equation residual (CGNR) method and the conjugate gradient normal equation error (CGNE) method for solving Problem 2. In recent years, Ramadan and Bayoumi (2018); Xie et al. (2014); and Dehghan and Hajarian (2011, 2010a,b) extended the CGNR and CGNE methods to solve several Sylvester matrix equations. As is well known, the CGNR and CGNE methods are obtained by applying the CG algorithm for solving either linear systems

A^{T} Ax = A^{T} b,

A A^{T} u = b,

where $x = A^{T} u$ Saad (2003). The CGNR and CGNE methods can be extended to solve Problem 2 as follows.

Extended form of CGNR method

Choose the arbitrary doubly symmetric matrix $Z (1) \in BS R^{n \times n}$ .

Compute

R_{i} (1) = N_{i} - \sum_{j = 1}^{f} A_{ij} Z (1) B_{ij}, i = 1, 2, . . ., g,

\hat{R} (1) = - E_{α}^{T} Z (1) E_{α},

\begin{matrix} Y (1) = \frac{1}{4} [\sum_{i = 1}^{g} \sum_{j = 1}^{f} (A_{ij}^{T} R_{i} (1) B_{ij}^{T} + (A_{ij}^{T} R_{i} (1) B_{ij}^{T})^{T} + J_{n} (A_{ij}^{T} R_{i} (1) B_{ij}^{T} \\ + (A_{ij}^{T} R_{i} (1) B_{ij}^{T})^{T}) J_{n}) \end{matrix}

+ E_{α} \hat{R} (1) E_{α}^{T} + (E_{α} \hat{R} (1) E_{α}^{T})^{T} + J_{n} (E_{α} \hat{R} (1) E_{α}^{T} + (E_{α} \hat{R} (1) E_{α}^{T})^{T}) J_{n}];

Set $P (1) = Y (1)$ .

For $k = 1, 2, . . .$ , until convergence, do

W_{i} (k) = \sum_{j = 1}^{f} A_{ij} P (k) B_{ij}, i = 1, 2, . . ., g, \hat{W} (k) = E_{α}^{T} P (k) E_{α},

α (k) = \frac{| | Y (k {) | |}^{2}}{\sum_{i = 1}^{g} | | W_{i} (k {) | |}^{2} + | | \hat{W} (k {) | |}^{2}},

Z (k + 1) = Z (k) + α (k) P (k),

\begin{matrix} R_{i} (k + 1) = R_{i} (k) - α (k) W_{i} (k), \\ i = 1, 2, . . ., g, \hat{R} (k + 1) = \hat{R} (k) - α (k) \hat{W} (k), \end{matrix}

\begin{matrix} Y (k + 1) = \frac{1}{4} [\sum_{i = 1}^{g} \sum_{j = 1}^{f} (A_{ij}^{T} R_{i} (k + 1) B_{ij}^{T} + (A_{ij}^{T} R_{i} (k + 1) B_{ij}^{T})^{T} \\ + J_{n} (A_{ij}^{T} R_{i} (k + 1) B_{ij}^{T} + (A_{ij}^{T} R_{i} (k + 1) B_{ij}^{T})^{T}) J_{n}) \end{matrix}

\begin{matrix} + E_{α} \hat{R} (k + 1) E_{α}^{T} + (E_{α} \hat{R} (k + 1) E_{α}^{T})^{T} \\ + J_{n} (E_{α} \hat{R} (k + 1) E_{α}^{T} + (E_{α} \hat{R} (k + 1) E_{α}^{T})^{T}) J_{n}]; \end{matrix}

β (k) = \frac{| | Y (k + {1) | |}^{2}}{| | Y (k {) | |}^{2}},

P (k + 1) = Y (k + 1) + β (k) P (k) .

Extended form of CGNE method

Choose the arbitrary doubly symmetric matrix $Z (1) \in BS R^{n \times n}$ .

Compute

R_{i} (1) = N_{i} - \sum_{j = 1}^{f} A_{ij} Z (1) B_{ij}, i = 1, 2, . . ., g,

\hat{R} (1) = - E_{α}^{T} Z (1) E_{α},

\begin{matrix} P (1) = \frac{1}{4} [\sum_{i = 1}^{g} \sum_{j = 1}^{f} (A_{ij}^{T} R_{i} (1) B_{ij}^{T} + (A_{ij}^{T} R_{i} (1) B_{ij}^{T})^{T} \\ + J_{n} (A_{ij}^{T} R_{i} (1) B_{ij}^{T} + (A_{ij}^{T} R_{i} (1) B_{ij}^{T})^{T}) J_{n}) \\ + E_{α} \hat{R} (1) E_{α}^{T} + (E_{α} \hat{R} (1) E_{α}^{T})^{T} + J_{n} (E_{α} \hat{R} (1) E_{α}^{T} \\ + (E_{α} \hat{R} (1) E_{α}^{T})^{T}) J_{n}]; \end{matrix}

For $k = 1, 2, . . .$ , until convergence, do

α (k) = \frac{\sum_{i = 1}^{g} | | R_{i} (k {) | |}^{2} + | | \hat{R} (k {) | |}^{2}}{| | P (k) | |},

W_{i} (k) = \sum_{j = 1}^{f} A_{ij} P (k) B_{ij}, i = 1, 2, . . ., g, \hat{W} (k) = E_{α}^{T} P (k) E_{α},

Z (k + 1) = Z (k) + α (k) P (k),

\begin{matrix} R_{i} (k + 1) = R_{i} (k) - α (k) W_{i} (k), i = 1, 2, . . ., g, \\ \hat{R} (k + 1) = \hat{R} (k) - α (k) \hat{W} (k), \end{matrix}

\begin{matrix} Y (k + 1) = \frac{1}{4} [\sum_{i = 1}^{g} \sum_{j = 1}^{f} (A_{ij}^{T} R_{i} (k + 1) B_{ij}^{T} + (A_{ij}^{T} R_{i} (k + 1) B_{ij}^{T})^{T} \\ + J_{n} (A_{ij}^{T} R_{i} (k + 1) B_{ij}^{T} + (A_{ij}^{T} R_{i} (k + 1) B_{ij}^{T})^{T}) J_{n}) \end{matrix}

\begin{matrix} + E_{α} \hat{R} (k + 1) E_{α}^{T} + (E_{α} \hat{R} (k + 1) E_{α}^{T})^{T} \\ + J_{n} (E_{α} \hat{R} (k + 1) E_{α}^{T} + (E_{α} \hat{R} (k + 1) E_{α}^{T})^{T}) J_{n}], \end{matrix}

β (k) = \frac{\sum_{i = 1}^{g} | | R_{i} (k + {1) | |}^{2} + | | \hat{R} (k + {1) | |}^{2}}{\sum_{i = 1}^{g} | | R_{i} (k {) | |}^{2} + | | \hat{R} (k {) | |}^{2}},

P (k + 1) = Y (k + 1) + β (k) P (k) .

Illustrative examples

In this section, the computational results are presented to support our theoretical discussion developed in the previous section. In what follows, we compare the performance of Algorithm 1 with the extended form of CGNR and CGNE algorithms for solving Problem 1. All numerical results were done with MATLAB.

Example 1: The first example involves the Sylvester matrix equation $AXB + CXD = M$ with the following parameters

A = (\begin{matrix} - 2 & 4000 & 10 & 12 & 2 & - 12 & - 15 & - 24 & - 18 & - 900 & - 3 & 6 \\ 16 & 16 & 18 & 4 & 6 & - 2 & - 3 & - 3 & - 6 & - 6 & - 9 & - 9 \\ 12 & 12 & 4 & 6 & - 18 & - 18 & - 27 & - 24 & - 21 & - 15 & - 6 & - 6 \\ 6 & 10 & 12 & 16 & 8 & - 18 & 15 & 15 & - 27 & - 24 & - 12 & - 6 \\ 6 & 6 & 12 & 12 & 4 & - 14 & - 12 & - 12 & - 18 & - 27 & - 6 & - 9 \\ 14 & 16 & 10 & 18 & 12 & 6 & - 3 & - 6 & - 15 & - 24 & - 27 & 18 \\ 2 & - 4000 & - 10 & - 12 & - 2 & 12 & - 2 & 4000 & 10 & 12 & 2 & - 12 \\ - 16 & - 16 & - 18 & - 4 & - 6 & 2 & 16 & 16 & 18 & 4 & 6 & - 2 \\ - 12 & - 12 & - 4 & - 6 & 18 & 18 & 12 & 12 & 4 & 6 & - 18 & - 18 \\ - 6 & - 10 & - 12 & - 16 & - 8 & 18 & 6 & 10 & 12 & 16 & 8 & - 18 \\ - 6 & - 6 & - 12 & - 12 & - 4 & 14 & 6 & 6 & 12 & 12 & 4 & - 14 \\ - 14 & - 16 & - 10 & - 18 & - 12 & - 6 & 14 & 16 & 10 & 18 & 12 & 6 \end{matrix}),

B = (\begin{matrix} 15 & 24 & 18 & 900 & 3 & - 6 & - 2 & 4000 & 10 & 12 & 2 & - 12 \\ 3 & 3 & 6 & 6 & 9 & 9 & 16 & 16 & 18 & 4 & 6 & - 2 \\ 27 & 24 & 21 & 15 & 6 & 6 & 12 & 12 & 4 & 6 & - 18 & - 18 \\ - 15 & - 15 & 27 & 24 & 12 & 6 & 6 & 10 & 12 & 16 & 8 & - 18 \\ 12 & 12 & 18 & 27 & 6 & 9 & 6 & 6 & 12 & 12 & 4 & - 14 \\ 3 & 6 & 15 & 24 & 27 & - 18 & 14 & 16 & 10 & 18 & 12 & 6 \\ 15 & 24 & 18 & 900 & 3 & - 6 & 2 & - 4000 & - 10 & - 12 & - 2 & 12 \\ 3 & 3 & 6 & 6 & 9 & 9 & - 16 & - 16 & - 18 & - 4 & - 6 & 2 \\ 27 & 24 & 21 & 15 & 6 & 6 & - 12 & - 12 & - 4 & - 6 & 18 & 18 \\ - 15 & - 15 & 27 & 24 & 12 & 6 & - 6 & - 10 & - 12 & - 16 & - 8 & 18 \\ 12 & 12 & 18 & 27 & 6 & 9 & - 6 & - 6 & - 12 & - 12 & - 4 & 14 \\ 3 & 6 & 15 & 24 & 27 & - 18 & - 14 & - 16 & - 10 & - 18 & - 12 & - 6 \end{matrix}),

C = (\begin{matrix} - 2 & 4000 & 10 & 12 & 2 & - 12 & - 17 & 3976 & - 8 & - 888 & - 1 & - 6 \\ 16 & 16 & 18 & 4 & 6 & - 2 & 13 & 13 & 12 & - 2 & - 3 & - 11 \\ 12 & 12 & 4 & 6 & - 18 & - 18 & - 15 & - 12 & - 17 & - 9 & - 24 & - 24 \\ 6 & 10 & 12 & 16 & 8 & - 18 & 21 & 25 & - 15 & - 8 & - 4 & - 24 \\ 6 & 6 & 12 & 12 & 4 & - 14 & - 6 & - 6 & - 6 & - 15 & - 2 & - 23 \\ 14 & 16 & 10 & 18 & 12 & 6 & 11 & 10 & - 5 & - 6 & - 15 & 24 \\ 17 & - 3976 & 8 & 888 & 1 & 6 & - 2 & 4000 & 10 & 12 & 2 & - 12 \\ - 13 & - 13 & - 12 & 2 & 3 & 11 & 16 & 16 & 18 & 4 & 6 & - 2 \\ 15 & 12 & 17 & 9 & 24 & 24 & 12 & 12 & 4 & 6 & - 18 & - 18 \\ - 21 & - 25 & 15 & 8 & 4 & 24 & 6 & 10 & 12 & 16 & 8 & - 18 \\ 6 & 6 & 6 & 15 & 2 & 23 & 6 & 6 & 12 & 12 & 4 & - 14 \\ - 11 & - 10 & 5 & 6 & 15 & - 24 & 14 & 16 & 10 & 18 & 12 & 6 \end{matrix}),

D = (\begin{matrix} - 2 & 4000 & 10 & 12 & 2 & - 12 & - 2 & 4000 & 10 & 12 & 2 & - 12 \\ 16 & 16 & 18 & 4 & 6 & - 2 & 16 & 16 & 18 & 4 & 6 & - 2 \\ 12 & 12 & 4 & 6 & - 18 & - 18 & 12 & 12 & 4 & 6 & - 18 & - 18 \\ 6 & 10 & 12 & 16 & 8 & - 18 & 6 & 10 & 12 & 16 & 8 & - 18 \\ 6 & 6 & 12 & 12 & 4 & - 14 & 6 & 6 & 12 & 12 & 4 & - 14 \\ 14 & 16 & 10 & 18 & 12 & 6 & 14 & 16 & 10 & 18 & 12 & 6 \\ 2 & - 4000 & - 10 & - 12 & - 2 & 12 & - 2 & 4000 & 10 & 12 & 2 & - 12 \\ - 16 & - 16 & - 18 & - 4 & - 6 & 2 & 16 & 16 & 18 & 4 & 6 & - 2 \\ - 12 & - 12 & - 4 & - 6 & 18 & 18 & 12 & 12 & 4 & 6 & - 18 & - 18 \\ - 6 & - 10 & - 12 & - 16 & - 8 & 18 & 6 & 10 & 12 & 16 & 8 & - 18 \\ - 6 & - 6 & - 12 & - 12 & - 4 & 14 & 6 & 6 & 12 & 12 & 4 & - 14 \\ - 14 & - 16 & - 10 & - 18 & - 12 & - 6 & 14 & 16 & 10 & 18 & 12 & 6 \end{matrix}),

M = (\begin{matrix} 779812 & 710496 & 7636206 & - 78639900 & 8263420 & - 3811200 & 7476128 & - 784932024 & 5080232 & 14801268 & 6296148 & - 3515132 \\ 2402 & - 19842 & 55629 & - 483006 & 35696 & - 6654 & 40868 & - 2684240 & 54844 & 111268 & 39344 & - 78564 \\ - 53376 & - 16464 & - 112932 & - 1400688 & - 63384 & - 26280 & - 77240 & - 12196112 & - 107960 & - 35120 & - 920 & 57760 \\ 2328 & 73836 & 44586 & - 646572 & 32292 & - 34692 & - 3600 & - 1898832 & 6840 & 96048 & 24648 & - 63992 \\ - 29930 & 38154 & - 51111 & - 164838 & - 37376 & - 19242 & - 36548 & 37368 & - 30556 & 6096 & - 2692 & - 11152 \\ - 4720 & - 164028 & 1254 & - 1112604 & 29180 & - 4212 & 70000 & - 2900704 & 70280 & 79832 & 30640 & - 46512 \\ 260 & - 19620 & 330 & 9270 & 350 & - 150 & 659012 & 27571496 & 693908 & 523410 & - 353654 & - 1013298 \\ 84 & 60300 & 594 & 16542 & 390 & - 702 & 25668 & 2560416 & 33612 & 21366 & 2238 & - 25326 \\ 492 & 12852 & 990 & 27702 & 870 & - 774 & 59724 & - 303072 & 74340 & 99378 & 28266 & - 80538 \\ - 192 & 156096 & 792 & 21942 & 330 & - 1278 & 33180 & 5015904 & 43644 & 558 & - 19290 & - 39726 \\ 232 & 32472 & 660 & 18432 & 520 & - 624 & 56800 & 1207240 & 66832 & 66228 & 8156 & - 68772 \\ 64 & 20160 & 264 & 7362 & 190 & - 282 & 27460 & 7455664 & 52132 & 53178 & - 8494 & - 79290 \end{matrix}) .

We can easily see that $κ (B^{T} \otimes A + D^{T} \otimes C) = | | B^{T} \otimes A + D^{T} \otimes C | | | | (B^{T} \otimes A + D^{T} \otimes C)^{†} | | = 5.7691 \times 10^{6}$ . This implies that the above Sylvester matrix equation is ill-conditioned. As far as we know, solving any ill-conditioned problem is a challenge in the areas of applied sciences. Let $α = {1, 12}$ and

X_{p} = (\begin{matrix} 1 - 3 \\ 5 - 7 \end{matrix}) .

Hence, we have

{\hat{X}}_{p} = (\begin{matrix} 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & - 3 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 5 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & - 7 \end{matrix}), E_{α} = (\begin{matrix} 1 & 0 \\ 0 & 0 \\ 0 & 0 \\ 0 & 0 \\ 0 & 0 \\ 0 & 0 \\ 0 & 0 \\ 0 & 0 \\ 0 & 0 \\ 0 & 0 \\ 0 & 0 \\ 0 & 1 \end{matrix}),

and

N = M - A {\hat{X}}_{p} B - C {\hat{X}}_{p} D

= (\begin{matrix} 779516 & 582184 & 7635370 & - 78616812 & 8261892 & - 3810624 & 7477312 & - 784946680 & 5081032 & 14802732 & 6297148 & - 3514580 \\ 1684 & - 176912 & 54652 & - 509736 & 35588 & - 6456 & 41200 & - 2956016 & 54304 & 110704 & 39376 & - 77664 \\ - 55260 & - 450972 & - 115566 & - 1420416 & - 65076 & - 25776 & - 75224 & - 12698096 & - 107960 & - 34364 & 340 & 60028 \\ 168 & - 384996 & 41874 & - 671664 & 30840 & - 34512 & - 1560 & - 2448816 & 6720 & 96660 & 25884 & - 61580 \\ - 32164 & - 400800 & - 53658 & - 202740 & - 38212 & - 19368 & - 34880 & - 553064 & - 31056 & 6084 & - 1812 & - 8788 \\ - 992 & 376820 & 4046 & - 1017648 & 27744 & - 2160 & 68936 & - 1949376 & 72240 & 81932 & 30612 & - 49620 \\ - 752 & - 193072 & - 304 & - 42168 & 1800 & - 888 & 658584 & 27166776 & 692448 & 521388 & - 354396 & - 1012356 \\ - 816 & - 32800 & 296 & - 6600 & 1304 & - 1704 & 25640 & 2364272 & 33032 & 20616 & 2032 & - 24792 \\ - 1536 & - 289920 & 24 & - 62568 & 3816 & - 2664 & 58944 & - 1012368 & 71760 & 95796 & 26940 & - 78900 \\ - 3624 & - 292368 & - 1776 & - 65736 & 1440 & - 3096 & 34224 & 4188624 & 42024 & - 1116 & - 19164 & - 36972 \\ - 1552 & - 225936 & - 408 & - 49392 & 2336 & - 1968 & 56576 & 646664 & 65072 & 63900 & 7444 & - 67308 \\ 392 & 96688 & 832 & 22152 & 384 & - 552 & 27264 & 7595520 & 52392 & 53436 & - 8532 & - 79764 \end{matrix}) .

By the mentioned algorithms with the initial matrix $Z (1) = 0$ , we obtain the matrix sequence ${Z (k)}$ to solve Problems 1 and 2. The numerical results are shown in Figures 1 and 2 where

r (k) = \underset{10}{\log} \sqrt{{\underset{R (k)}{\underset{︸}{| | N - AZ (k) B - CZ (k) D | |}}}^{2} + {\underset{\hat{R} (k)}{\underset{︸}{| | E_{α}^{T} Z (k) E_{α} | |}}}^{2}},

and

r_{1} (k) = \underset{10}{\log} | | M - A (Z (k) + {\hat{X}}_{p}) B - C (Z (k) + {\hat{X}}_{p}) D | | .

The numerical results illustrate that Algorithm 1 is capable of solving the ill-conditioned matrix equation. Now we can obtain the solutions of Problems 2 and 1 as follows

\begin{array}{l} Z^{*} \\ = (\begin{matrix} \begin{matrix} - 0.0000 \end{matrix} & - 12.0000 & - 4.0000 & 4.0000 & 36.0000 & 1.0000 & 1.0000 & - 12.0000 & - 4.0000 & 4.0000 & 36.0000 & \begin{matrix} - 0.0000 \end{matrix} \\ - 12.0000 & - 20.0000 & 6.0000 & 10.0000 & 28.0000 & 36.0000 & - 12.0000 & - 20.0000 & 6.0000 & 10.0000 & 28.0000 & 36.0000 \\ - 4.0000 & 6.0000 & 4.0000 & 16.0000 & 10.0000 & 4.0000 & - 4.0000 & 6.0000 & 4.0000 & 16.0000 & 10.0000 & 4.0000 \\ 4.0000 & 10.0000 & 16.0000 & 4.0000 & 6.0000 & - 4.0000 & 4.0000 & 10.0000 & 16.0000 & 4.0000 & 6.0000 & - 4.0000 \\ 36.0000 & 28.0000 & 10.0000 & 6.0000 & - 20.0000 & - 12.0000 & 36.0000 & 28.0000 & 10.0000 & 6.0000 & - 20.0000 & - 12.0000 \\ 1.0000 & 36.0000 & 4.0000 & - 4.0000 & - 12.0000 & 1.0000 & 1.0000 & 36.0000 & 4.0000 & - 4.0000 & - 12.0000 & 1.0000 \\ 1.0000 & - 12.0000 & - 4.0000 & 4.0000 & 36.0000 & 1.0000 & 1.0000 & - 12.0000 & - 4.0000 & 4.0000 & 36.0000 & 1.0000 \\ - 12.0000 & - 20.0000 & 6.0000 & 10.0000 & 28.0000 & 36.0000 & - 12.0000 & - 20.0000 & 6.0000 & 10.0000 & 28.0000 & 36.0000 \\ - 4.0000 & 6.0000 & 4.0000 & 16.0000 & 10.0000 & 4.0000 & - 4.0000 & 6.0000 & 4.0000 & 16.0000 & 10.0000 & 4.0000 \\ 4.0000 & 10.0000 & 16.0000 & 4.0000 & 6.0000 & - 4.0000 & 4.0000 & 10.0000 & 16.0000 & 4.0000 & 6.0000 & - 4.0000 \\ 36.0000 & 28.0000 & 10.0000 & 6.0000 & - 20.0000 & - 12.0000 & 36.0000 & 28.0000 & 10.0000 & 6.0000 & - 20.0000 & - 12.0000 \\ \begin{matrix} - 0.0000 \end{matrix} & 36.0000 & 4.0000 & - 4.0000 & - 12.0000 & 1.0000 & 1.0000 & 36.0000 & 4.0000 & - 4.0000 & - 12.0000 & \begin{matrix} - 0.0000 \end{matrix} \end{matrix}) \in V_{α}, \end{array}

and

X^{*} = Z^{*} + {\hat{X}}_{p}

= (\begin{matrix} - 0.0000 & - 12.0000 & - 4.0000 & 4.0000 & 36.0000 & 1.0000 & 1.0000 & - 12.0000 & - 4.0000 & 4.0000 & 36.0000 & - 3.0000 \\ - 12.0000 & - 20.0000 & 6.0000 & 10.0000 & 28.0000 & 36.0000 & - 12.0000 & - 20.0000 & 6.0000 & 10.0000 & 28.0000 & 36.0000 \\ - 4.0000 & 6.0000 & 4.0000 & 16.0000 & 10.0000 & 4.0000 & - 4.0000 & 6.0000 & 4.0000 & 16.0000 & 10.0000 & 4.0000 \\ 4.0000 & 10.0000 & 16.0000 & 4.0000 & 6.0000 & - 4.0000 & 4.0000 & 10.0000 & 16.0000 & 4.0000 & 6.0000 & - 4.0000 \\ 36.0000 & 28.0000 & 10.0000 & 6.0000 & - 20.0000 & - 12.0000 & 36.0000 & 28.0000 & 10.0000 & 6.0000 & - 20.0000 & - 12.0000 \\ 1.0000 & 36.0000 & 4.0000 & - 4.0000 & - 12.0000 & 1.0000 & 1.0000 & 36.0000 & 4.0000 & - 4.0000 & - 12.0000 & 1.0000 \\ 1.0000 & - 12.0000 & - 4.0000 & 4.0000 & 36.0000 & 1.0000 & 1.0000 & - 12.0000 & - 4.0000 & 4.0000 & 36.0000 & 1.0000 \\ - 12.0000 & - 20.0000 & 6.0000 & 10.0000 & 28.0000 & 36.0000 & - 12.0000 & - 20.0000 & 6.0000 & 10.0000 & 28.0000 & 36.0000 \\ - 4.0000 & 6.0000 & 4.0000 & 16.0000 & 10.0000 & 4.0000 & - 4.0000 & 6.0000 & 4.0000 & 16.0000 & 10.0000 & 4.0000 \\ 4.0000 & 10.0000 & 16.0000 & 4.0000 & 6.0000 & - 4.0000 & 4.0000 & 10.0000 & 16.0000 & 4.0000 & 6.0000 & - 4.0000 \\ 36.0000 & 28.0000 & 10.0000 & 6.0000 & - 20.0000 & - 12.0000 & 36.0000 & 28.0000 & 10.0000 & 6.0000 & - 20.0000 & - 12.0000 \\ 5.0000 & 36.0000 & 4.0000 & - 4.0000 & - 12.0000 & 1.0000 & 1.0000 & 36.0000 & 4.0000 & - 4.0000 & - 12.0000 & - 7.0000 \end{matrix}) \in T_{α} .

where

| | M - A X^{*} B - C X^{*} D | | = 8.9834 \times 10^{- 07} .

The obtained results show that the speed of convergence for Algorithm 1 is high. Also, Figures 1 and 2 show that our computations are stable.

Figure 1.

Comparing the performance of Algorithm 1, CGNR algorithm and CGNE algorithm for Example 1.

Figure 2.

Comparing the performance of Algorithm 1, CGNR algorithm and CGNE algorithm for Example 1.

Example 2: As another example, consider the pair of matrix equations $A_{1} X B_{1} = M_{1}$ and $A_{2} X B_{2} = M_{2}$ where

A_{1} = (\begin{matrix} - 3000 & 600 & 15 & 18 & 3 & - 18 & - 10 & - 16 & - 12 & - 600 & - 2 & 4 \\ 24 & 24 & 27 & 6 & 9 & - 3 & - 2 & - 2 & - 4 & - 4 & - 6 & - 6 \\ 18 & 18 & 6 & 9 & - 27 & - 27 & - 18 & - 16 & - 14 & - 10 & - 4 & - 4 \\ 9 & 15 & 18 & 24 & 12 & - 27 & 10 & 10 & - 18 & - 16 & - 8 & - 4 \\ 9 & 9 & 18 & 18 & 6 & - 21 & - 8 & - 8 & - 12 & - 18 & - 4 & - 6 \\ 21 & 24 & 15 & 27 & 18 & 9 & - 2 & - 4 & - 10 & - 16 & - 18 & 12 \\ 2980 & - 632 & - 39 & - 1218 & - 7 & 26 & - 3000 & 600 & 15 & 18 & 3 & - 18 \\ - 28 & - 28 & - 35 & - 14 & - 21 & - 9 & 24 & 24 & 27 & 6 & 9 & - 3 \\ - 54 & - 50 & - 34 & - 29 & 19 & 19 & 18 & 18 & 6 & 9 & - 27 & - 27 \\ 11 & 5 & - 54 & - 56 & - 28 & 19 & 9 & 15 & 18 & 24 & 12 & - 27 \\ - 25 & - 25 & - 42 & - 54 & - 14 & 9 & 9 & 9 & 18 & 18 & 6 & - 21 \\ - 25 & - 32 & - 35 & - 59 & - 54 & 15 & 21 & 24 & 15 & 27 & 18 & 9 \end{matrix}),

B_{1} = (\begin{matrix} 10 & 16 & 12 & 600 & 2 & - 4 & - 2980 & 632 & 39 & 1218 & 7 & - 26 \\ 2 & 2 & 4 & 4 & 6 & 6 & 28 & 28 & 35 & 14 & 21 & 9 \\ 18 & 16 & 14 & 10 & 4 & 4 & 54 & 50 & 34 & 29 & - 19 & - 19 \\ - 10 & - 10 & 18 & 16 & 8 & 4 & - 11 & - 5 & 54 & 56 & 28 & - 19 \\ 8 & 8 & 12 & 18 & 4 & 6 & 25 & 25 & 42 & 54 & 14 & - 9 \\ 2 & 4 & 10 & 16 & 18 & - 12 & 25 & 32 & 35 & 59 & 54 & - 15 \\ 10 & 16 & 12 & 600 & 2 & - 4 & 2980 & - 632 & - 39 & - 1218 & - 7 & 26 \\ 2 & 2 & 4 & 4 & 6 & 6 & - 28 & - 28 & - 35 & - 14 & - 21 & - 9 \\ 18 & 16 & 14 & 10 & 4 & 4 & - 54 & - 50 & - 34 & - 29 & 19 & 19 \\ - 10 & - 10 & 18 & 16 & 8 & 4 & 11 & 5 & - 54 & - 56 & - 28 & 19 \\ 8 & 8 & 12 & 18 & 4 & 6 & - 25 & - 25 & - 42 & - 54 & - 14 & 9 \\ 2 & 4 & 10 & 16 & 18 & - 12 & - 25 & - 32 & - 35 & - 59 & - 54 & 15 \end{matrix}),

A_{2} = (\begin{matrix} 2980 & - 632 & - 39 & - 1218 & - 7 & 26 & 2990 & - 616 & - 27 & - 618 & - 5 & 22 \\ - 28 & - 28 & - 35 & - 14 & - 21 & - 9 & - 26 & - 26 & - 31 & - 10 & - 15 & - 3 \\ - 54 & - 50 & - 34 & - 29 & 19 & 19 & - 36 & - 34 & - 20 & - 19 & 23 & 23 \\ 11 & 5 & - 54 & - 56 & - 28 & 19 & 1 & - 5 & - 36 & - 40 & - 20 & 23 \\ - 25 & - 25 & - 42 & - 54 & - 14 & 9 & - 17 & - 17 & - 30 & - 36 & - 10 & 15 \\ - 25 & - 32 & - 35 & - 59 & - 54 & 15 & - 23 & - 28 & - 25 & - 43 & - 36 & 3 \\ - 2990 & 616 & 27 & 618 & 5 & - 22 & - 2980 & 632 & 39 & 1218 & 7 & - 26 \\ 26 & 26 & 31 & 10 & 15 & 3 & 28 & 28 & 35 & 14 & 21 & 9 \\ 36 & 34 & 20 & 19 & - 23 & - 23 & 54 & 50 & 34 & 29 & - 19 & - 19 \\ - 1 & 5 & 36 & 40 & 20 & - 23 & - 11 & - 5 & 54 & 56 & 28 & - 19 \\ 17 & 17 & 30 & 36 & 10 & - 15 & 25 & 25 & 42 & 54 & 14 & - 9 \\ 23 & 28 & 25 & 43 & 36 & - 3 & 25 & 32 & 35 & 59 & 54 & - 15 \end{matrix}),

B_{2} = (\begin{matrix} - 2980 & 632 & 39 & 1218 & 7 & - 26 & - 2980 & 632 & 39 & 1218 & 7 & - 26 \\ 28 & 28 & 35 & 14 & 21 & 9 & 28 & 28 & 35 & 14 & 21 & 9 \\ 54 & 50 & 34 & 29 & - 19 & - 19 & 54 & 50 & 34 & 29 & - 19 & - 19 \\ - 11 & - 5 & 54 & 56 & 28 & - 19 & - 11 & - 5 & 54 & 56 & 28 & - 19 \\ 25 & 25 & 42 & 54 & 14 & - 9 & 25 & 25 & 42 & 54 & 14 & - 9 \\ 25 & 32 & 35 & 59 & 54 & - 15 & 25 & 32 & 35 & 59 & 54 & - 15 \\ 2980 & - 632 & - 39 & - 1218 & - 7 & 26 & - 2980 & 632 & 39 & 1218 & 7 & - 26 \\ - 28 & - 28 & - 35 & - 14 & - 21 & - 9 & 28 & 28 & 35 & 14 & 21 & 9 \\ - 54 & - 50 & - 34 & - 29 & 19 & 19 & 54 & 50 & 34 & 29 & - 19 & - 19 \\ 11 & 5 & - 54 & - 56 & - 28 & 19 & - 11 & - 5 & 54 & 56 & 28 & - 19 \\ - 25 & - 25 & - 42 & - 54 & - 14 & 9 & 25 & 25 & 42 & 54 & 14 & - 9 \\ - 25 & - 32 & - 35 & - 59 & - 54 & 15 & 25 & 32 & 35 & 59 & 54 & - 15 \end{matrix}),

M_{1} = (\begin{matrix} - 9530538 & - 10994604 & - 12626164 & 33313238 & - 19878206 & 7955066 & 9721987 & - 1962707 & - 18655 & - 3707878 & 127731 & 7199 \\ 84244 & 93972 & 138080 & - 815332 & 214076 & - 101416 & - 1097968 & 240068 & 23698 & 473704 & 12942 & - 17540 \\ - 49352 & - 108396 & - 56760 & - 6150288 & 46040 & - 25860 & 23862226 & - 4879438 & - 113576 & - 9275180 & 214552 & 70794 \\ 41882 & 16776 & 50412 & - 3634014 & 154142 & - 85758 & 10674983 & - 2477947 & - 353081 & - 4669358 & - 259851 & 124111 \\ - 22182 & - 62616 & - 38084 & - 3877190 & 12062 & - 10646 & 9959783 & - 2073091 & - 81065 & - 3905654 & 52293 & 26419 \\ 113514 & 169464 & 210236 & 3237578 & 247998 & - 103030 & - 10379957 & 2238033 & 171507 & 4270266 & 53585 & - 82649 \\ 3927694 & 5310912 & 8843068 & 52035490 & 10151698 & - 3918014 & - 130473693 & 52382365 & 28912973 & 104386178 & 29378431 & - 11992025 \\ - 69172 & - 106488 & - 162976 & - 2058304 & - 184116 & 79712 & 10598604 & - 2554704 & - 458586 & - 4889736 & - 390390 & 194964 \\ - 290812 & - 404496 & - 489040 & - 7103632 & - 523180 & 231320 & - 17090388 & 3423392 & 7234 & 6505984 & - 243058 & - 22672 \\ 39290 & 47460 & 59052 & - 1494810 & 190254 & - 145662 & - 2391079 & 452139 & - 33201 & 799338 & - 83459 & 36767 \\ - 121378 & - 167556 & - 224028 & - 2859366 & - 232502 & 73614 & - 899833 & 115909 & - 71107 & 185726 & - 99601 & 39945 \\ - 12604 & 16632 & - 64528 & 3743804 & - 139764 & 60884 & 7690806 & - 1808922 & - 295458 & - 3541416 & - 246090 & 162690 \end{matrix}),

M_{2} = (\begin{matrix} 107088997 & - 47116789 & - 28308593 & - 94452698 & - 29003699 & 11774505 & - 351837631 & 175581831 & 78761889 & 288855726 & 88854427 & - 20481627 \\ 12536504 & - 2326204 & 186214 & - 4429736 & 405210 & - 49436 & 28280964 & - 7375656 & - 1609366 & - 13697788 & - 1377518 & 538460 \\ - 15127826 & 3690190 & 698060 & 7042676 & 605940 & - 285094 & - 56936070 & 9949698 & - 903044 & 19980304 & - 1658272 & 158458 \\ - 1298795 & 368591 & 103315 & 633586 & 105337 & - 11517 & - 42989911 & 9527523 & 962469 & 18277314 & 934951 & - 324261 \\ 107875 & 234729 & 264189 & 424926 & 314111 & - 87143 & - 29477417 & 5206173 & - 485821 & 10433702 & - 684687 & 94481 \\ 9917638 & - 1687914 & 294198 & - 3296184 & 495038 & - 58502 & 10246474 & - 3706590 & - 1490734 & - 6584092 & - 1348866 & 538106 \\ - 109436885 & 47476949 & 28195501 & 95130562 & 28828975 & - 11731713 & 353969519 & - 176203239 & - 78966429 & - 290032470 & - 89059223 & 20555379 \\ - 7437364 & 1210016 & - 287594 & 2304520 & - 450726 & 94156 & - 33398304 & 8468268 & 1685546 & 15779852 & 1382122 & - 572020 \\ 18011656 & - 4530500 & - 975022 & - 8652712 & - 891234 & 394844 & 53737740 & - 9489240 & 716998 & - 19096268 & 1394006 & - 98972 \\ 6392083 & - 1322407 & - 36647 & - 2449034 & 48667 & - 5655 & 38104223 & - 8323995 & - 722481 & - 15983370 & - 732059 & 230601 \\ 4131415 & - 1247595 & - 436995 & - 2353890 & - 457105 & 156985 & 25103827 & - 4357119 & 462515 & - 8815906 & 583833 & - 90679 \\ - 11069306 & 1951758 & - 263622 & 3759528 & - 479638 & 63010 & - 9121406 & 3415362 & 1416254 & 6062732 & 1314090 & - 532366 \end{matrix}) .

Noting the condition number

κ (\begin{matrix} (\begin{matrix} B_{1}^{T} \otimes A_{1} \\ B_{2}^{T} \otimes A_{2} \end{matrix}) \end{matrix}) = 1.1657 \times 10^{5},

the pair of matrix equations is ill-conditioned. Suppose that $α = {5, 6, 7, 8}$ and

X_{p} = (\begin{matrix} 1 & 2 & 3 & 4 \\ 5 & 6 & 7 & 8 \\ 9 & 10 & 11 & 12 \\ 13 & 14 & 15 & 16 \end{matrix}) .

This leads to

{\hat{X}}_{p} = (\begin{matrix} 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 1 & 2 & 3 & 4 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 5 & 6 & 7 & 8 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 9 & 10 & 11 & 12 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 13 & 14 & 15 & 16 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \end{matrix}), E_{α} = (\begin{matrix} 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \end{matrix}),

N_{1} = M_{1} - A_{1} {\hat{X}}_{p} B_{1} =

(\begin{matrix} - 9520920 & - 10981332 & - 12609648 & 33609216 & - 19865016 & 7953444 & 11119698 & - 2248818 & - 23568 & - 4237872 & 142188 & 4914 \\ 85288 & 95388 & 139888 & - 785888 & 215496 & - 101612 & - 959670 & 212550 & 24144 & 422592 & 14988 & - 17910 \\ - 35184 & - 88740 & - 32496 & - 5702808 & 65520 & - 28164 & 25978578 & - 5316042 & - 125004 & - 10083360 & 233808 & 67950 \\ 39608 & 13656 & 46496 & - 3702040 & 151032 & - 85360 & 10354284 & - 2412876 & - 352632 & - 4548720 & - 263616 & 124740 \\ - 15264 & - 53064 & - 26208 & - 3663672 & 21552 & - 11808 & 10968264 & - 2279712 & - 84828 & - 4288368 & 62580 & 24804 \\ 113080 & 168792 & 209536 & 3216256 & 247368 & - 103016 & - 10482864 & 2261448 & 174636 & 4313232 & 54348 & - 82908 \\ 4414056 & 5983224 & 9677520 & 67130484 & 10819248 & - 3999084 & - 59154258 & 37745538 & 28617816 & 77279376 & 30086748 & - 12101754 \\ - 79840 & - 121104 & - 181360 & - 2375164 & - 198696 & 81596 & 9105462 & - 2252430 & - 457308 & - 4329216 & - 408456 & 198018 \\ - 304104 & - 422904 & - 511824 & - 7519884 & - 541440 & 233508 & - 19058082 & 3828306 & 16656 & 7255680 & - 261756 & - 19842 \\ 31324 & 36540 & 45328 & - 1732064 & 179364 & - 144260 & - 3509280 & 678708 & - 32010 & 1219440 & - 96834 & 39018 \\ - 126780 & - 174972 & - 233328 & - 3021432 & - 239892 & 74556 & - 1664004 & 271104 & - 69870 & 473424 & - 108462 & 41418 \\ - 24256 & 720 & - 84640 & 3403400 & - 155664 & 62984 & 6088356 & - 1486260 & - 296136 & - 2942784 & - 266832 & 166284 \end{matrix}),

and

N_{2} = M_{2} - A_{2} {\hat{X}}_{p} B_{2} =

(\begin{matrix} 35854854 & - 32495526 & - 28011504 & - 67375056 & - 29709660 & 11883750 & - 282625638 & 158638566 & 75367440 & 257196060 & 86714544 & - 19605882 \\ 14935662 & - 2818830 & 175992 & - 5342016 & 428844 & - 53082 & 25949706 & - 6805050 & - 1495152 & - 12631644 & - 1305528 & 509022 \\ - 12981690 & 3254922 & 695280 & 6235680 & 631284 & - 289338 & - 59015406 & 10461342 & - 797616 & 20938428 & - 1591488 & 130446 \\ - 1255104 & 358116 & 101358 & 614448 & 104598 & - 11310 & - 43034052 & 9537576 & 963642 & 18295500 & 935598 & - 324330 \\ 1367700 & - 23424 & 259446 & - 53232 & 326934 & - 89154 & - 30700992 & 5505924 & - 425526 & 10993980 & - 646650 & 78858 \\ 12099828 & - 2138004 & 282528 & - 4129344 & 514968 & - 61452 & 8123784 & - 3188040 & - 1388064 & - 5616012 & - 1284276 & 511956 \\ - 39316602 & 33083082 & 27901608 & 68474256 & 29522940 & - 11839026 & 285838386 & - 159525186 & - 75625800 & - 258868860 & - 86953344 & 19693710 \\ - 9813330 & 1696938 & - 278580 & 3206400 & - 474864 & 97938 & - 31090638 & 7902966 & 1571868 & 14723244 & 1310412 & - 542574 \\ 14673402 & - 3852378 & - 969432 & - 7395648 & - 929820 & 401250 & 56973294 & - 10284846 & 553656 & - 20585724 & 1290600 & - 55734 \\ 7281720 & - 1503588 & - 38682 & - 2784816 & 58590 & - 7278 & 37241436 & - 8112072 & - 679230 & - 15586812 & - 704706 & 219222 \\ 2409900 & - 895920 & - 431790 & - 1702320 & - 475470 & 159930 & 26774592 & - 4766964 & 379470 & - 9582396 & 531378 & - 68994 \\ - 13587720 & 2469384 & - 252288 & 4718016 & - 504048 & 66744 & - 6673692 & 2816508 & 1296672 & 4944012 & 1238748 & - 501624 \end{matrix}) .

In Figures 3 and 4, we present the numerical results of the presented algorithms with $Z (1) = 0$ where

r (k) = \underset{10}{\log} \sqrt{{\underset{R_{1} (k)}{\underset{︸}{| | N_{1} - A_{1} Z (k) B_{1} | |}}}^{2} + {\underset{R_{2} (k)}{\underset{︸}{| | N_{2} - A_{2} Z (k) B_{2} | |}}}^{2} + {\underset{\hat{R} (k)}{\underset{︸}{| | E_{α}^{T} Z (k) E_{α} | |}}}^{2}},

and

r_{1} (k) = \underset{10}{\log} \sqrt{| | M_{1} - A_{1} (Z (k) + {\hat{X}}_{p}) B_{1} {| |}^{2} + | | M_{2} - A_{2} (Z (k) + {\hat{X}}_{p}) B_{2} {| |}^{2}} .

Here, the solutions of Problems 2 and 1 can be computed as follows

\begin{array}{l} Z^{*} = \\ (\begin{matrix} 0.0000 & 36.0000 & - 0.0000 & - 48.0000 & 108.0000 & 192.0000 & - 0.0000 & 36.0000 & 0.0000 & - 48.0000 & 108.0000 & 192.0000 \\ 36.0000 & 6.0000 & 0.0000 & 0.0000 & 12.0000 & 108.0000 & 36.0000 & 6.0000 & 0.0000 & 0.0000 & 12.0000 & 108.0000 \\ - 0.0000 & 0.0000 & 7.0000 & 7.0000 & 0.0000 & - 48.0000 & - 0.0000 & 0.0000 & 7.0000 & 7.0000 & 0.0000 & - 48.0000 \\ - 48.0000 & 0.0000 & 7.0000 & 7.0000 & 0.0000 & - 0.0000 & - 48.0000 & - 0.0000 & 7.0000 & 7.0000 & 0.0000 & 0.0000 \\ 108.0000 & 12.0000 & 0.0000 & 0.0000 & \begin{matrix} - 0.0000 \end{matrix} & \begin{matrix} 0.0000 \end{matrix} & \begin{matrix} - 0.0000 \end{matrix} & \begin{matrix} - 0.0000 \end{matrix} & - 0.0000 & 0.0000 & 6.0000 & 36.0000 \\ 192.0000 & 108.0000 & - 48.0000 & - 0.0000 & \begin{matrix} 0.0000 \end{matrix} & \begin{matrix} - 0.0000 \end{matrix} & \begin{matrix} - 0.0000 \end{matrix} & \begin{matrix} - 0.0000 \end{matrix} & - 48.0000 & - 0.0000 & 36.0000 & - 0.0000 \\ - 0.0000 & 36.0000 & - 0.0000 & - 48.0000 & \begin{matrix} - 0.0000 \end{matrix} & \begin{matrix} - 0.0000 \end{matrix} & \begin{matrix} - 0.0000 \end{matrix} & \begin{matrix} 0.0000 \end{matrix} & - 0.0000 & - 48.0000 & 108.0000 & 192.0000 \\ 36.0000 & 6.0000 & 0.0000 & - 0.0000 & \begin{matrix} - 0.0000 \end{matrix} & \begin{matrix} - 0.0000 \end{matrix} & \begin{matrix} 0.0000 \end{matrix} & \begin{matrix} - 0.0000 \end{matrix} & 0.0000 & 0.0000 & 12.0000 & 108.0000 \\ 0.0000 & 0.0000 & 7.0000 & 7.0000 & - 0.0000 & - 48.0000 & - 0.0000 & 0.0000 & 7.0000 & 7.0000 & 0.0000 & - 48.0000 \\ - 48.0000 & 0.0000 & 7.0000 & 7.0000 & 0.0000 & - 0.0000 & - 48.0000 & 0.0000 & 7.0000 & 7.0000 & 0.0000 & - 0.0000 \\ 108.0000 & 12.0000 & 0.0000 & 0.0000 & 6.0000 & 36.0000 & 108.0000 & 12.0000 & 0.0000 & 0.0000 & 6.0000 & 36.0000 \\ 192.0000 & 108.0000 & - 48.0000 & 0.0000 & 36.0000 & - 0.0000 & 192.0000 & 108.0000 & - 48.0000 & - 0.0000 & 36.0000 & 0.0000 \end{matrix}) \in V_{α}, \end{array}

and

X^{*} = Z^{*} + {\hat{X}}_{p} =

(\begin{matrix} 0.0000 & 36.0000 & - 0.0000 & - 48.0000 & 108.0000 & 192.0000 & - 0.0000 & 36.0000 & 0.0000 & - 48.0000 & 108.0000 & 192.0000 \\ 36.0000 & 6.0000 & 0.0000 & 0.0000 & 12.0000 & 108.0000 & 36.0000 & 6.0000 & 0.0000 & 0.0000 & 12.0000 & 108.0000 \\ - 0.0000 & 0.0000 & 7.0000 & 7.0000 & 0.0000 & - 48.0000 & - 0.0000 & 0.0000 & 7.0000 & 7.0000 & 0.0000 & - 48.0000 \\ - 48.0000 & 0.0000 & 7.0000 & 7.0000 & 0.0000 & - 0.0000 & - 48.0000 & - 0.0000 & 7.0000 & 7.0000 & 0.0000 & 0.0000 \\ 108.0000 & 12.0000 & 0.0000 & 0.0000 & 1.0000 & 2.0000 & 3.0000 & 4.0000 & - 0.0000 & 0.0000 & 6.0000 & 36.0000 \\ 192.0000 & 108.0000 & - 48.0000 & - 0.0000 & 5.0000 & 6.0000 & 7.0000 & 8.0000 & - 48.0000 & - 0.0000 & 36.0000 & - 0.0000 \\ - 0.0000 & 36.0000 & - 0.0000 & - 48.0000 & 9.0000 & 10.0000 & 11.0000 & 12.0000 & - 0.0000 & - 48.0000 & 108.0000 & 192.0000 \\ 36.0000 & 6.0000 & 0.0000 & - 0.0000 & 13.0000 & 14.0000 & 15.0000 & 16.0000 & 0.0000 & 0.0000 & 12.0000 & 108.0000 \\ 0.0000 & 0.0000 & 7.0000 & 7.0000 & - 0.0000 & - 48.0000 & - 0.0000 & 0.0000 & 7.0000 & 7.0000 & 0.0000 & - 48.0000 \\ - 48.0000 & 0.0000 & 7.0000 & 7.0000 & 0.0000 & - 0.0000 & - 48.0000 & 0.0000 & 7.0000 & 7.0000 & 0.0000 & - 0.0000 \\ 108.0000 & 12.0000 & 0.0000 & 0.0000 & 6.0000 & 36.0000 & 108.0000 & 12.0000 & 0.0000 & 0.0000 & 6.0000 & 36.0000 \\ 192.0000 & 108.0000 & - 48.0000 & 0.0000 & 36.0000 & - 0.0000 & 192.0000 & 108.0000 & - 48.0000 & - 0.0000 & 36.0000 & 0.0000 \end{matrix}) \in T_{α},

where

\sqrt{| | M_{1} - A_{1} X^{*} B_{1} {| |}^{2} + | | M_{2} - A_{2} X^{*} B_{2} {| |}^{2}} = 4.4315 \times 10^{- 07} .

Figures 3 and 4 display the computational superiority of Algorithm 1 because of its good accuracy.

Figure 3.

Comparing the performance of Algorithm 1, CGNR algorithm and CGNE algorithm for Example 2.

Figure 4.

Comparing the performance of Algorithm 1, CGNR algorithm and CGNE algorithm for Example 2.

Conclusions

In this contribution, we developed the HS version of Bi-CR algorithm for finding the least Frobenius norm partially doubly symmetric solution of the general Sylvester matrix equations (1.10) within a finite number of iterations in the absence of round-off errors. To show the efficiency of Algorithm 1, we gave two numerical examples. Overall, the numerical examples illustrated that the results obtained by Algorithm 1 have good accuracy and high speed.

Footnotes

Appendix A. Proof of Lemma 1

This lemma can be proved via induction on $v$ and $u$ . We begin by showing (3.3)–(3.6) for $v = 2$ and $u = 1$ . Observe that

\begin{matrix} \sum_{i = 1}^{g} tr (R_{i} (2)^{T} W_{i} (1)) + tr (\hat{R} (2)^{T} \hat{W} (1)) \\ = \sum_{i = 1}^{g} tr ((R_{i} (1) - α (1) W_{i} (1))^{T} W_{i} (1)) \end{matrix}

+ tr ((\hat{R} (1) - α (1) \hat{W} (1))^{T} \hat{W} (1)) = 0,

tr (S (2)^{T} Y (1)) = tr ((S (1) - β (1) Y (1))^{T} Y (1)) = 0,

\begin{matrix} tr (Y (2)^{T} Y (1)) = | | Y (1) | |^{2} \\ - α (1) tr ((\frac{1}{4} [\sum_{i = 1}^{g} \sum_{j = 1}^{f} (A_{ij}^{T} R_{i} (1) B_{ij}^{T} + (A_{ij}^{T} R_{i} (1) B_{ij}^{T})^{T} \\ + J_{n} (A_{ij}^{T} R_{i} (1) B_{ij}^{T} + (A_{ij}^{T} R_{i} (1) B_{ij}^{T})^{T}) J_{n}) \end{matrix}

\begin{matrix} + E_{α} \hat{R} (1) E_{α}^{T} + (E_{α} \hat{R} (1) E_{α}^{T})^{T} + J_{n} (E_{α} \hat{R} (1) E_{α}^{T} \\ + (E_{α} \hat{R} (1) E_{α}^{T})^{T}) J_{n} {])}^{T} Y (1)) \end{matrix}

\begin{matrix} + | | Y (1) | |^{2} \\ [1 - \frac{\sum_{i = 1}^{g} tr (R_{i} {(1)}^{T} (\sum_{j = 1}^{f} (A_{ij} Y (1) B_{ij}))) + tr (\hat{R} {(1)}^{T} (E_{α}^{T} Y (1) E_{α}))}{| | Y {(1) | |}^{2}} \end{matrix}

- \frac{\sum_{i = 1}^{g} tr (W_{i} {(1)}^{T} (\sum_{j = 1}^{f} (A_{ij} Y (1) B_{ij}))) + tr (\hat{W} {(1)}^{T} (E_{α}^{T} Y (1) E_{α}))}{\sum_{i = 1}^{g} | | W_{i} {(1) | |}^{2} + | | \hat{W} {(1) | |}^{2}}

+ [\frac{\sum_{i = 1}^{g} tr (R_{i} {(1)}^{T} (\sum_{j = 1}^{f} (A_{ij} S (1) B_{ij}))) + tr (\hat{R} {(1)}^{T} (E_{α}^{T} S (1) E_{α}))}{\sum_{i = 1}^{g} | | W_{i} {(1) | |}^{2} + | | \hat{W} {(1) | |}^{2}}

\times \frac{\sum_{i = 1}^{g} tr (W_{i} {(1)}^{T} (\sum_{j = 1}^{f} (A_{ij} Y (1) B_{ij}))) + tr (\hat{W} {(1)}^{T} (E_{α}^{T} Y (1) E_{α}))}{| | Y {(1) | |}^{2}}]] = 0,

\sum_{i = 1}^{g} tr (W_{i} (2)^{T} W_{i} (1)) + tr (\hat{W} (2)^{T} \hat{W} (1))

= \sum_{i = 1}^{g} | | W_{i} (1) | |^{2} + | | \hat{W} (1) | |^{2} - β (1) \sum_{i = 1}^{g} tr ((\sum_{j = 1}^{f} A_{ij} Y (1) B_{ij}))^{T} W_{i} (1))

- β (1) tr ((E_{α}^{T} Y (1) E_{α}))^{T} \hat{W} (1))

\begin{matrix} + \sum_{i = 1}^{g} | | W_{i} (1) | |^{2} \\ [1 - \frac{\sum_{i = 1}^{g} tr (R_{i} {(1)}^{T} (\sum_{j = 1}^{f} (A_{ij} Y (1) B_{ij}))) + tr (\hat{R} {(1)}^{T} (E_{α}^{T} Y (1) E_{α}))}{| | Y {(1) | |}^{2}} \end{matrix}

- \frac{\sum_{i = 1}^{g} tr (W_{i} {(1)}^{T} (\sum_{j = 1}^{f} (A_{ij} Y (1) B_{ij}))) + tr (\hat{W} {(1)}^{T} (E_{α}^{T} Y (1) E_{α}))}{\sum_{i = 1}^{g} | | W_{i} {(1) | |}^{2} + | | \hat{W} {(1) | |}^{2}}

+ [\frac{\sum_{i = 1}^{g} tr (R_{i} {(1)}^{T} (\sum_{j = 1}^{f} (A_{ij} S (1) B_{ij}))) + tr (\hat{R} {(1)}^{T} (E_{α}^{T} S (1) E_{α}))}{\sum_{i = 1}^{g} | | W_{i} {(1) | |}^{2} + | | \hat{W} {(1) | |}^{2}}

\times \frac{\sum_{i = 1}^{g} tr (W_{i} {(1)}^{T} (\sum_{j = 1}^{f} (A_{ij} Y (1) B_{ij}))) + tr (\hat{W} {(1)}^{T} (E_{α}^{T} Y (1) E_{α}))}{| | Y {(1) | |}^{2}}]] = 0 .

It follows from the above that the relations (3.3)–(3.6) hold for $v = 2$ and $u = 1$ . Let us now for $u < w < r$ assume that

(5.1)

\sum_{i = 1}^{g} tr (R_{i} (w)^{T} W_{i} (u)) + tr (\hat{R} (w)^{T} \hat{W} (u)) = 0, tr (S (w)^{T} Y (u)) = 0,

(5.2)

tr (Y {(w)}^{T} Y (u)) = 0, \sum_{i = 1}^{g} tr (W_{i} {(w)}^{T} W_{i} (u)) + tr (\hat{W} {(w)}^{T} \hat{W} (u)) = 0 .

These, together with Remark 2, imply

\sum_{i = 1}^{g} tr (R_{i} (w + 1)^{T} W_{i} (u)) = \sum_{i = 1}^{g} tr ((R_{i} (w) - α (w) W_{i} (w))^{T} W_{i} (u))

+ tr ((\hat{R} (w) - α (w) \hat{W} (w))^{T} \hat{W} (u)) = 0,

tr (S (w + 1)^{T} Y (u)) = tr ((S (w) - β (w) Y (w))^{T} Y (u)) = 0,

tr (Y (w + 1)^{T} Y (u))

\begin{matrix} = \frac{1}{β (u)} [\sum_{i = 1}^{g} tr (R_{i} (w + 1)^{T} (- W_{i} (u + 1) + γ (u) W_{i} (u) \\ + W_{i} (u) - γ (u - 1) W_{i} (u - 1))) \end{matrix}

+ tr (\hat{R} (w + 1)^{T} (- \hat{W} (u + 1) + γ (u) \hat{W} (u) + \hat{W} (u) - γ (u - 1) \hat{W} (u - 1))]

(5.3)

= - \frac{1}{β (u)} [\sum_{i = 1}^{g} tr (R_{i} (w + 1)^{T} (W_{i} (u + 1)) + tr (\hat{R} (w + 1)^{T} (\hat{W} (u + 1))],

\sum_{i = 1}^{g} tr (W_{i} (w + 1)^{T} W_{i} (u)) + tr (\hat{W} (w + 1)^{T} \hat{W} (u))

\begin{matrix} = \frac{1}{α (u)} tr ((S (w + 1))^{T} (Y (u) - γ (u - 1) \\ Y (u - 1) - Y (u + 1) + γ (u) Y (u)))) \end{matrix}

(5.4)

= - \frac{1}{α (u)} tr ((S (w + 1)^{T} Y (u + 1)) .

Letting $u = w$ , we get

\begin{matrix} \sum_{i = 1}^{g} tr (R_{i} (w + 1)^{T} W_{i} (w)) \\ = \sum_{i = 1}^{g} tr ((R_{i} (w) - α (w) W_{i} (w))^{T} W_{i} (w)) + tr ((\hat{R} (w) - α (w) \hat{W} (w))^{T} \hat{W} (w)) \end{matrix}

= \sum_{i = 1}^{g} tr (R_{i} (w)^{T} W_{i} (w)) + tr (\hat{R} (w)^{T} \hat{W} (w))

\begin{matrix} - \sum_{i = 1}^{g} tr (R_{i} (w)^{T} (W_{i} (w) - γ (w - 1) W_{i} (w - 1))) \\ - tr (\hat{R} (w)^{T} (\hat{W} (w) - γ (w - 1) \hat{W} (w - 1))) = 0, \end{matrix}

tr (S (w + 1)^{T} Y (w)) = tr ((S (w) - β (w) Y (w))^{T} Y (w))

= tr (S (w)^{T} Y (w)) - tr ((Y (w) - γ (w - 1) Y (w - 1)))^{T} S (w)) = 0,

tr (Y (w + 1)^{T} Y (w)) = | | Y (w) | |^{2}

+ | | Y (w) | |^{2} [1 - \frac{tr ((Y (w) - γ (w - 1) Y (w - {1))}^{T} Y (w))}{tr (Y {(w)}^{T} Y (w))}

- \frac{\sum_{i = 1}^{g} tr ((W_{i} (w) - γ (w - 1) W_{i} (w - {1))}^{T} W_{i} (w)) + tr ((\hat{W} (w) - γ (w - 1) \hat{W} (w - {1))}^{T} \hat{W} (w))}{\sum_{i = 1}^{g} tr (W_{i} {(w)}^{T} W_{i} (w)) + tr (\hat{W} {(w)}^{T} \hat{W} (w))}] = 0,

\sum_{i = 1}^{g} tr (W_{i} (w + 1)^{T} W_{i} (w)) + tr (\hat{W} (w + 1)^{T} \hat{W} (w))

= \sum_{i = 1}^{g} | | W_{i} (w) | |^{2} + | | \hat{W} (w) | |^{2} + [\sum_{i = 1}^{g} | | W_{i} (w) | |^{2} + | | \hat{W} (w) | |^{2}] [1 - \frac{tr ((Y (w) - γ (w - 1) Y (w - {1))}^{T} Y (w))}{tr (Y {(w)}^{T} Y (w))}

- \frac{\sum_{i = 1}^{g} tr ((W_{i} (w) - γ (w - 1) W_{i} (w - {1))}^{T} W_{i} (w)) + tr ((\hat{W} (w) - γ (w - 1) \hat{W} (w - {1))}^{T} \hat{W} (w))}{\sum_{i = 1}^{g} tr (W_{i} {(w)}^{T} W_{i} (w)) + tr (\hat{W} {(w)}^{T} \hat{W} (w))}] = 0 .

From $tr (Y (w)^{T} Y (u)) = 0$ , $\sum_{i = 1}^{g} tr (R_{i} (w + 1)^{T} (W_{i} (w)) + tr (\hat{R} (w + 1)^{T} (\hat{W} (w)) = 0$ and (5.3), it can be obtained that

tr (Y (w + 1)^{T} Y (u)) = 0 .

Also, it follows from $\sum_{i = 1}^{g} tr (W_{i} (w)^{T} W_{i} (u)) + tr (\hat{W} (w)^{T} \hat{W} (u)) = 0$ , $tr ((S (w + 1)^{T} Y (w)) = 0$ and (5.4) that

\sum_{i = 1}^{g} tr (W_{i} (w + 1)^{T} W_{i} (u)) + tr (\hat{W} (w + 1)^{T} \hat{W} (u)) = 0 .

By combining the above relations, the proof is now complete by induction.

Acknowledgements

The author would like to express his great gratitude to the associate editor and three anonymous referees for very useful comments and constructive suggestions that led to a significant improvement in the quality and presentation of this paper.

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) received no financial support for the research, authorship, and/or publication of this article.

ORCID iD

Masoud Hajarian

References

Andrew

(1998) Classroom note: Centrosymmetric matrices. SIAM Review 40(3): 697–698.

Cantoni

Butler

(1976) Properties of the eigenvectors of persymmetric matrices with applications to communication theory. IEEE Transactions on Communications 24(8): 804–809.

Chronopoulos

(1991) S-step iterative methods for (non)symmetric (in)definite linear systems. SIAM Journal on Numerical Analysis 28(6): 1776–1789.

Chronopoulos

(1994) On the squared unsymmetric Lanczos method. Journal of Computational and Applied Mathematics 54(1): 65–78.

Chronopoulos

Kincaid

(2001) On the Odir iterative method for non-symmetric indefinite linear systems. Numerical Linear Algebra with Applications 8(2): 71–82.

Chronopoulos

Kucherov

(2010) Block s-step Krylov iterative methods. Numerical Linear Algebra with Applications 17(1): 3–15.

Datta

Morgera

(1989) On the reducibility of centrosymmetric matrices– Applications in engineering problems. Circuits, Systems and Signal Processing 8(1): 71–96.

DeAlba

Johnson

(1995) Possible inertia combinations in the Stein and Lyapunov equations. Linear Algebra and its Applications 222: 227–240.

Dehghan

Hajarian

(2010a) An iterative method for solving the generalized coupled Sylvester matrix equations over generalized bisymmetric matrices. Applied Mathematical Modelling 34(3): 639–654.

10.

Dehghan

Hajarian

(2010b) The general coupled matrix equations over generalized bisymmetric matrices. Linear Algebra and its Applications 432(6): 1531–1552.

11.

Dehghan

Hajarian

(2011) Analysis of an iterative algorithm to solve the generalized coupled Sylvester matrix equations. Applied Mathematical Modelling 35(7): 3285–3300.

12.

Dehghan

Hajarian

(2012) On the generalized reflexive and anti-reflexive solutions to a system of matrix equations. Linear Algebra and its Applications 437(11): 2793–2812.

13.

Ding

Chen

(2005a) Gradient based iterative algorithms for solving a class of matrix equations. IEEE Transactions on Automatic Control 50(8): 1216–1221.

14.

Ding

Chen

(2005b) Iterative least squares solutions of coupled Sylvester matrix equations. Systems & Control Letters 54(2): 95–107.

15.

Ding

Chen

(2006) On iterative solutions of general coupled matrix equations. SIAM Journal on Control and Optimization 44(6): 2269–2284.

16.

Ding

Zhang

(2014) Gradient-based iterative algorithm for a class of the coupled matrix equations related to control systems. IET Control Theory & Applications 8(15): 1588–1595.

17.

Good

(1970) The inverse of a centrosymmetric matrix. Technometrics 12(4): 925–928.

18.

Hajarian

(2013a) Matrix form of the Bi-CGSTAB method for solving the coupled Sylvester matrix equations. IET Control Theory & Applications 7(14): 1828–1833.

19.

Hajarian

(2013b) Matrix iterative methods for solving the Sylvester-transpose and periodic Sylvester matrix equations. Journal of the Franklin Institute 350(10): 3328–3341.

20.

Hajarian

(2014a) Matrix form of the CGS method for solving general coupled matrix equations. Applied Mathematics Letters 34: 37–42.

21.

Hajarian

(2014b) Matrix GPBiCG algorithms for solving the general coupled matrix equations. IET Control Theory & Applications 9(1): 74–81.

22.

Hajarian

(2015) Developing BiCOR and CORS methods for coupled Sylvester-transpose and periodic Sylvester matrix equations. Applied Mathematical Modelling 39(19): 6073–6084.

23.

Hajarian

(2016a) Generalized conjugate direction algorithm for solving the general coupled matrix equations over symmetric matrices. Numerical Algorithms 73(3): 591–609.

24.

Hajarian

(2016b) Lanczos version of BCR algorithm for solving the generalised second-order Sylvester matrix equation EVF+GVH+BVC=DWE+M. IET Control Theory & Applications 11(2): 273–281.

25.

Hajarian

(2017) Convergence analysis of the MCGNR algorithm for the least squares solution group of discrete-time periodic coupled matrix equations. Transactions of the Institute of Measurement and Control 39(1): 29–42.

26.

Hajarian

(2018a) Computing symmetric solutions of general Sylvester matrix equations via Lanczos version of biconjugate residual algorithm. Computers & Mathematics with Applications 76(4): 686–700.

27.

Hajarian

(2018b) Convergence properties of BCR method for generalized Sylvester matrix equation over generalized reflexive and anti-reflexive matrices. Linear and Multilinear Algebra 66(10): 1975–1990.

28.

Hajarian

(2018c) Finding solutions for periodic discrete-time generalized coupled Sylvester matrix equations via the generalized BCR method. Transactions of the Institute of Measurement and Control 40(2): 647–656.

29.

Hajarian

(2018d) Solving constrained quadratic inverse eigenvalue problem via conjugate direction method. Computers & Mathematics with Applications 76(10): 2384–2401.

30.

Hajarian

(2019) An efficient algorithm based on Lanczos type of BCR to solve constrained quadratic inverse eigenvalue problems. Journal of Computational and Applied Mathematics 346: 418–431.

31.

Hanna

Mansoori

(2003) A centrosymmetric matrix based technique for the interpolation of a Hermitian signal. Numerical Linear Algebra with Applications 10(8): 701–720.

32.

Huang

(2014) The modified conjugate gradient methods for solving a class of generalized coupled Sylvester-transpose matrix equations. Computers & Mathematics with Applications 67(8): 1545–1558.

33.

(2018) The unified frame of alternating direction method of multipliers for three classes of matrix equations arising in control theory. Asian Journal of Control 20(1): 437–454.

34.

Larin

(2009a) On solution of Sylvester equation. Journal of Automation and Information Sciences 41: 1–7.

35.

Larin

(2009b) Solutions of matrix equations in problems of mechanics and control. International Applied Mechanics 45(8): 847–872.

36.

Zhang

(2010) The nearness problems for symmetric centrosymmetric with a special submatrix constraint. Numerical Algorithms 55(1): 39–57.

37.

Liao

Bai

Lei

(2005) Best approximate solution of matrix equation AXB+CYD=E. SIAM. SIAM Journal on Matrix Analysis and Applications 27(3): 675–688.

38.

Mattingly

(2000) Even order regular magic squares are singular. The American Mathematical Monthly 107(9): 777–782.

39.

Peng

Zhang

(2007) The reflexive and anti-reflexive solutions of the matrix equation AHXB = C. Journal of Computational and Applied Mathematics 200(2): 749–760.

40.

Ramadan

Bayoumi

(2018) A modified gradient-based algorithm for solving extended Sylvester-conjugate matrix equations. Asian Journal of Control 20(1): 228–235.

41.

Saad

(2003) Iterative Methods for Sparse Linear Systems (2nd ed.). Philadelpha: SIAM.

42.

Simoncini

(2016) Computational methods for linear matrix equations. SIAM Review 58(3): 377–441.

43.

Vespucci

Broyden

(2001) Implementation of different computational variations of biconjugate residual methods. Computers & Mathematics with Applications 42(8-9): 1239–1253.

44.

Wang

(2014) Systems of coupled generalized Sylvester matrix equations. Automatica 50(11): 2840–2844.

45.

Wang

(2009) Ranks and the least-norm of the general solution to a system of quaternion matrix equations. Linear Algebra and its Applications 430(5-6): 1626–1640.

46.

Wang

Zhang

(2008) The reflexive re-nonnegative definite solution to a quaternion matrix equation. Electronic Journal of Linear Algebra 17(1): 88–101.

47.

Wang

Woude

Chang

(2009) A system of real quaternion matrix equations with applications. Linear Algebra and its Applications 431(12): 2291–2303.

48.

Xie

Liu

Yang

(2010) Gradient based and least squares based iterative algorithms for matrix equations AXB + CXTD = F. Applied Mathematics and Computation 217(5): 2191–2199.

49.

Xie

Huang

(2014) Iterative method to solve the generalized coupled Sylvester-transpose linear matrix equations over reflexive or anti-reflexive matrix. Computers & Mathematics with Applications 67(11): 2071–2084.

50.

Xiong

Qin

(2011) The common Re-nnd and Re-pd solutions to the matrix equations AX = C and XB = D. Applied Mathematics and Computation 218(7): 3330–3337.

51.

Yamamoto

(2006) Toward the Sinc-Galerkin method for the Poisson problem in one type of curvilinear coordinate domain. Electronic Transactions on Numerical Analysis 23: 63–75.

52.

Yuan

Liao

(2011) Least squares solution of the quaternion matrix equation X−A ^ B = C with the least norm. Linear and Multilinear Algebra 59(9): 985–998.

53.

Zhang

A-G

Sun

H-J

(2018) An implicit iterative algorithm with a parameter for Itô Lyapunov matrix equations. International Journal of Systems Science 49(2): 425–434.

54.

Zhou

Lam

Duan

(2010) Gradient-based maximal convergence rate iterative method for solving linear matrix equations. International Journal of Computer Mathematics 87(3): 515–527.