Numerical differentiation for two-dimensional scattered data on arbitrary domain base on Hermite extension with an implicit iteration process

1.   Introduction

A problem that occurs frequently in a variety of mathematical contexts, is to find the common invariant subspaces of a single matrix or set of matrices. In the case of a single endomorphism or matrix, it is relatively easy to find all the invariant subspaces by using the Jordan normal form. Also, some theoretical results are given only for the invariant subspaces of two matrices. However, when there are more than two matrices, the problem becomes much harder, and unexpected invariant subspaces may occur. No systematic method is known. In a recent article ^[1], we have provided a new algorithms to determine common invariant subspaces of a single matrix or of a set of matrices systematically.

In the present article we consider a more general version of this problem, that is, providing two algorithms for simultaneous block triangularization and block diagonalization of sets of matrices. One of the main steps in the first two proposed algorithms, consists of finding the common invariant subspaces of matrices using the new method proposed in the recent article ^[1]. It is worth mentioning that an efficient algorithm to explicitly compute a transfer matrix which realizes the simultaneous block diagonalization of unitary matrices whose decomposition in irreducible blocks (common invariant subspaces) is known from elsewhere is given in ^[2]. An application of simultaneous block-diagonalization of normal matrices in quantum theory is presented in ^[3].

In this article we shall be concerned with finite dimensions only. Of course the fact that a single complex matrix can always be put into triangular form follows readily from the Jordan normal form theorem ^[4]. For a set of matrices, Jacobson in ^[5] introduced the notion of a composition series for a collection of matrices. The idea of a composition series for a group is quite familiar. The Jordan-Hölder Theorem ^[4] states that any two composition series of the same group have the same length and the same composition factors (up to permutation). Jacobson in ^[5] characterized the simultaneous block triangularization of a set of matrices by the existence of a chain $\{0\} = V_0 \subset V_1 \subset... \subset V_t = \mathbb{C}^n$ of invariant subspaces with dimension $dim(V_i/V_{i-1}) = n_i$. Therefore, in the context of a collection of matrices $\Omega = \{A_i \}_{ i = 1}^N$, the idea is to locate a common invariant subspace $V$ of minimal dimension $d$ of a set of matrices $\Omega$. Assume $V$ is generated by the (linearly independent) set $\mathcal{B}_1 = \{u_1, u_2, ..., u_d \}$, and let $\mathcal{B} = \{u_1, u_2, ..., u_d, u_{d+1}, u_{d+2}, ..., u_n\}$ be a basis of $\mathbb{C}^n$ containing $\mathcal{B}_1$. Upon setting $S = (u_1, u_2, ..., u_d, u_{d+1}, u_{d+2}, ..., u_n)$, $S^{-1}A_i S$ has the block triangular form

for $i = 1, ..., n$. Thereafter, one may define a quotient of the ambient vector space, and each of the matrices in the given collection will pass to this quotient. As such, one defines

Then one may begin again the process of looking for a common invariant subspace of minimal dimension of a set of matrices $\{T_i \}_{ i = 1}^N$ and iterate the procedure. Since all spaces and matrices are of finite dimension, the procedure must terminate at some point. Again, any two such composition series will be isomorphic. When the various quotients and submatrices are lifted back to the original vector space, one obtains precisely the block-triangular form for the original set of matrices. It is important to find a composition series in the construction in order to make the set of matrices as "block-triangular as possible."

Dubi ^[6] gave an algorithmic approach to simultaneous triangularization of a set of matrices based on the idea of Jacobson in ^[5]. In the case of simultaneous triangularization, it can be understood as the existence of a chain $\{0\} = V_0 \subset V_1 \subset... \subset V_t = \mathbb{C}^n$ of invariant subspaces with dimension $dim(V_i) = i$. We generalize his study to cover simultaneous block triangularization of a set of matrices. The generalized algorithm depends on the novel algorithm for constructing invariant subspaces of a set of matrices given in the recent article ^[1].

Specht ^[7] (see also ^[8]) proved that if the associative algebra $\mathcal{L}$ generated by a set of matrices $\Omega$ over $\mathbb{C}$ satisfies $\mathcal{L} = \mathcal{L}^{*}$, then $\Omega$ admits simultaneous block triangularization if and only if it admits simultaneous block diagonalization, in both cases via a unitary matrix. Following a result of Specht, we prove that a set of matrices $\Omega$ admits simultaneous block diagonalization if and only if the set $\Gamma = \Omega \cup \Omega^{*}$ admits simultaneous block triangularization. Finally, an algorithmic approach to simultaneous block diagonalization of a set of matrices based on this fact is proposed.

The latter part of this paper presents an alternate approach for simultaneous block diagonalization of a set of $n \times n$ matrices $\{A_s\}_{ s = 1}^N$ by an invertible matrix that does not require finding the common invariant subspaces. Maehara et al. ^[9] introduced an algorithm for simultaneous block diagonalization of a set of matrices by a unitary matrix based on the existence of a Hermitian commuting matrix. Here, we extend their algorithm to simultaneous block diagonalization of a set of matrices by an invertible matrix based on the existence of a commuting matrix which is not necessarily Hermitian. For example, consider the set of matrices $\Omega = \{A_i \}_{i = 1}^2$ where

The only Hermitian matrix commuting with the set $\Omega$ is the identity matrix. Therefore, we cannot apply the proposed algorithm given in ^[9]. However, one can verify that the following non Hermitian matrix $C$ commutes with all the matrices $\{A_i \}_{ i = 1}^2$

The matrix $C$ has distinct eigenvalues $\lambda_1 = 0, \lambda_2 = 1$ with algebraic multiplicities $n_1 = 2, n_2 = 1$, respectively. Moreover, the matrix $C$ is not diagonalizable. Therefore, we cannot construct the eigenvalue decomposition for the matrix $C$. However, one can decompose the matrix $C$ by its generalized eigen vectors as follows:

where

Initially, it is noted that the matrices $\{A_i \}_{ i = 1}^2$ can be decomposed into two diagonal blocks by the constructed invertible matrix $S$ where

Then, a new algorithm is developed for simultaneous block diagonalization by an invertible matrix based on the generalized eigenvectors of a commuting matrix. Moreover, a new characterization is presented by proving that the existence of a commuting matrix that possesses at least two distinct eigenvalues is the necessary and sufficient condition to guarantee the simultaneous block diagonalization by an invertible matrix.

An outline of the paper is as follows. In Section 2 we review several definitions pertaining to block-triangular and block-diagonal matrices and state several elementary consequences that follow from them. In Section 3, following a result of Specht ^[7] (see also ^[8]), we provide conditions for putting a set of matrices into block-diagonal form simultaneously. Furthermore, we apply the theoretical results to provide two algorithms that enable a collection of matrices to be put into block-triangular form or block-diagonal form simultaneously by a unitary matrix based on the existence of invariant subspaces. In Section 4, a new characterization is presented by proving that the existence of a commuting matrix that possesses at least two distinct eigenvalues is the necessary and sufficient condition to guarantee the simultaneous block diagonalization by an invertible matrix. Furthermore, we apply the theoretical results to provide an algorithm that enables a collection of matrices to be put into block-diagonal form simultaneously by an invertible matrix based on the existence of a commuting matrix. Sections 3 and 4 also provide concrete examples using the symbolic manipulation system Maple.

2.   Preliminaries

Let $\Omega$ be a set of $n \times n$ matrices over an algebraically closed field $\mathcal{F}$, and let $\mathcal{L}$ denote the algebra generated by $\Omega$ over $\mathcal{F}$. Similarly, let $\Omega^{*}$ be the set of the conjugate transpose of each matrix in $\Omega$ and $\mathcal{L}^{*}$ denote the algebra generated by $\Omega^{*}$ over $\mathcal{F}$.

Definition 2.1. An $n \times n$ matrix $A$ is given the notation $BT(n_1, ..., n_t)$ provided $A$ is block upper triangular with $t$ square blocks on the diagonal, of sizes $n_1, ..., n_t$, where $t \geq 2$ and $n_1+... +n_t = n$. That is, a block upper triangular matrix $A$ has the form

where ${\bf{A}}_{i, j}$ is a square matrix for all $i = 1, ..., t$ and $j = i, ..., t$.

Definition 2.2. A set of $n \times n$ matrices $\Omega$ is $BT(n_1, ..., n_t)$ if all of the matrices in $\Omega$ are $BT(n_1, ..., n_t)$.

Remark 2.3. A set of $n \times n$ matrices $\Omega$ admits a simultaneous triangularization if it is $BT(n_1, ..., n_t)$ with $n_i = 1$ for $i = 1, ..., t$.

Remark 2.4. A set of $n \times n$ matrices $\Omega$ is $BT(n_1, ..., n_t)$ if and only if the algebra $\mathcal{L}$ generated by $\Omega$ is $BT(n_1, ..., n_t)$.

Proposition 2.5. ^[7] (see also ^[8]) Let $\Omega$ be a nonempty set of complex $n \times n$ matrices. Then, there is a nonsingular matrix $S$ such that $S \Omega S^{-1}$ is $BT(n_1, ..., n_t)$ if and only if there is a unitary matrix $U$ such that $U \Omega U^{*}$ is $BT(n_1, ..., n_t)$.

Theorem 2.6. [5,Chapter Ⅳ] Let $\Omega$ be a nonempty set of complex $n \times n$ matrices. Then, there is a unitary matrix $U$ such that $U \Omega U^{*}$ is $BT(n_1, ..., n_t)$ if and only if the set $\Omega$ has a chain $\{0\} = V_0 \subset V_1 \subset... \subset V_t = \mathbb{C}^n$ of invariant subspaces with dimension $dim(V_i/V_{i-1}) = n_i$.

Definition 2.7. An $n \times n$ matrix $A$ is given the notation $BD(n_1, ..., n_t)$ provided $A$ is block diagonal with $t$ square blocks on the diagonal, of sizes $n_1, ..., n_t$, where $t \geq 2$, $n_1+... +n_t = n$, and the blocks off the diagonal are the zero matrices. That is, a block diagonal matrix $A$ has the form

where ${\bf{A}}_k$ is a square matrix for all $k = 1, ..., t$. In other words, matrix ${\bf{A}}$ is the direct sum of ${\bf{A}}_1, ..., {\bf{A}}_t$. It can also be indicated as ${\bf{A}}_{\text{1}} \oplus {\bf{A}}_{\text{2}} \oplus... \oplus {\bf{A}}_{\text{t}}$.

Definition 2.8. A set of $n \times n$ matrices $\Omega$ is $BD(n_1, ..., n_t)$ if all of the matrices in $\Omega$ are $BD(n_1, ..., n_t)$.

Remark 2.9. A set of $n \times n$ matrices $\Omega$ admits a simultaneous diagonalization if it is $BD(n_1, ..., n_t)$ with $n_i = 1$ for $i = 1, ..., t$.

Remark 2.10. A set of $n \times n$ matrices $\Omega$ is $BD(n_1, ..., n_t)$ if and only if the algebra $\mathcal{L}$ generated by $\Omega$ is $BD(n_1, ..., n_t)$.

Proposition 2.11. ^[7] (see also ^[8]) Let $\Omega$ be a nonempty set of complex $n \times n$ matrices and let $\mathcal{L}$ be the algebra generated by $\Omega$ over $\mathbb{C}$. Suppose $\mathcal{L} = \mathcal{L}^{*}$. Then, there is a nonsingular matrix $S$ such that $S \mathcal{L} S^{-1}$ is $BT(n_1, ..., n_t)$ if and only if there is a unitary matrix $U$ such that $U \mathcal{L} U^{*}$ is $BD(n_1, ..., n_t)$.

3.   Algorithms for simultaneous block triangularization and block diagonalization of a set of matrices based on the invariant subspaces

Dubi ^[6] gave an algorithmic approach to simultaneous triangularization of a set of $n \times n$ matrices. In this section, we will generalize his study to cover simultaneous block triangularization and simultaneous block diagonalization of a set of $n \times n$ matrices. The generalized algorithms depend on the novel algorithm for constructing invariant subspaces of a set of matrices given in the recent article ^[1] and Theorem 3.3.

Lemma 3.1. Let $\Omega$ be a nonempty set of complex $n \times n$ matrices, $\Omega^{*}$ be the set of the conjugate transpose of each matrix in $\Omega$ and $\mathcal{L}$ be the algebra generated by $\Gamma = \Omega \cup \Omega^{*}$. Then, $\mathcal{L} = \mathcal{L}^{*}$.

Proof. Let $A$ be a matrix in $\mathcal{L}$. Then, $A = P(B_1, ..., B_m)$ for some multivariate noncommutative polynomial $P(x_1, ..., x_m)$ and matrices $\{B_i\}_{i = 1}^m\in \Gamma$. Therefore, $A^{*} = P^*(B_1, ..., B_m) = Q(B_1^*, ..., B_m^*)$ for some multivariate noncommutative polynomial $Q(x_1, ..., x_m)$ where the matrices $\{B_i^*\}_{i = 1}^m\in \Gamma^* = \Gamma$. Hence, the matrix $A^* \in \mathcal{L}$

Lemma 3.2. Let $\Omega$ be a nonempty set of complex $n \times n$ matrices and $\Omega^{*}$ be the set of the conjugate transpose of each matrix in $\Omega$, and $\Gamma = \Omega \cup \Omega^{*}$. Then, there is a unitary matrix $U$ such that $U \Gamma U^{*}$ is $BD(n_1, ..., n_t)$ if and only if there is a unitary matrix $U$ such that $U \Omega U^{*}$ is $BD(n_1, ..., n_t)$.

Proof. Assume that there exists a unitary matrix $U$ such that $U \Omega U^{*}$ is $BD(n_1, ..., n_t)$. Then, $(U \Omega U^{*})^{*} = U \Omega^{*} U^{*}$ is $BD(n_1, ..., n_t)$. Hence, $U \Gamma U^{*}$ is $BD(n_1, ..., n_t)$.

Theorem 3.3. Let $\Omega$ be a nonempty set of complex $n \times n$ matrices and $\Omega^{*}$ be the set of the conjugate transpose of each matrix in $\Omega$, and $\Gamma = \Omega \cup \Omega^{*}$. Then, there is a unitary matrix $U$ such that $U \Omega U^{*}$ is $BD(n_1, ..., n_t)$ if and only if there is a unitary matrix $U$ such that $U \Gamma U^{*}$ is $BT(n_1, ..., n_t)$.

Proof. Let $\mathcal{L}$ be the algebra generated by $\Gamma$. Then, $\mathcal{L} = \mathcal{L}^{*}$ using Lemma 3.1. Now, by applying Proposition 2.11 and Lemma 3.2, the following statements are equivalent :

There is a unitary matrix $U$ such that $U \Gamma U^{*}$ is $BT(n_1, ..., n_t)$.

$\iff$ There is a unitary matrix $U$ such that $U \mathcal{L} U^{*}$ is $BT(n_1, ..., n_t)$.

$\iff$ There is a unitary matrix $U$ such that $U \mathcal{L} U^{*}$ is $BD(n_1, ..., n_t)$.

$\iff$ There is a unitary matrix $U$ such that $U \Gamma U^{*}$ is $BD(n_1, ..., n_t)$.

$\iff$ There is a unitary matrix $U$ such that $U \Omega U^{*}$ is $BD(n_1, ..., n_t)$.

3.1.   Algorithm $A$: Simultaneous block triangularization of a set of $n \times n$ matrices $\{A_i \}_{ i=1}^N$.

(1) Input: the set $\Omega = \{A_i \}_{ i = 1}^N$.

(2) Set $k = 0, \mathcal{B} = \phi, s = n, T_i = A_i, S_2 = I$.

(3) Search for a $d$-dimensional invariant subspace $V = \langle v_1, v_2, ..., v_d \rangle$ of a set of matrices $\{T_i \}_{ i = 1}^N$ starting from $d = 1$ up to $d = s-1$. If one does not exist and $k = 0$, abort and print "no simultaneous block triangularization". Else, if one does not exist and $k\ne 0$, go to step (8). Else, go to next step.

(4) Set $V_{k+1} = (S_2 v_1\; S_2 v_2\; ...\; S_2 v_d), \mathcal{B} = \mathcal{B} \cup \{S_2 v_1, S_2 v_2, ..., S_2 v_d\}, S_1 = (V_1\; V_2\; ...\; V_{k+1})$.

(5) Find a basis $\{u_1, u_2, ..., u_l \}$ for the orthogonal complement of $\mathcal{B}$.

(6) Set $S_2 = (u_1\; u_2\; ...\; u_l), S = (S_1\; S_2)$, and

$T_i = \left({\begin{array}{*{20}{c}} \textbf{0}_{(s-d)\times d} & \textbf{I}_{s-d} \end{array}} \right)S^{-1}A_i S \left({\begin{array}{*{20}{c}} \textbf{0}_{d \times (s-d)} \\ \textbf{I}_{s-d} \end{array}} \right)$.

(7) Set $k = k+1, s = s-d$, and return to step (3).

(8) Compute the QR decomposition of the invertible matrix $S$, by means of the Gram–Schmidt process, to convert it to a unitary matrix $Q$.

(9) Output: a unitary matrix $U$ as the conjugate transpose of the resulting matrix $Q$.

Remark 3.4. If one uses any non-orthogonal complement in step 5 of Algorithm $A$, then the matrix $S$ is invertible such that $S^{-1} \Omega S$ is $BT(n_1, ..., n_t)$. However, in such a case, one cannot guarantee that $U \Omega U^{*}$ is $BT(n_1, ..., n_t)$.

Example 3.5. The set of matrices $\Omega = \{A_i \}_{i = 1}^2$ admits simultaneous block triangularization where

Applying Algorithm $A$ to the set $\Omega$ can be summarized as follows:

● Input: $\Omega$.

● Initiation step:

We have $k = 0, \mathcal{B} = \phi, s = 6, T_1 = A_1, T_2 = A_2, S_2 = I$.

● In the first iteration:

We found two-dimensional invariant subspace $V = \langle e_1, e_4 \rangle$ of a set of matrices $\{T_i \}_{ i = 1}^2$. Therefore, $\mathcal{B} = \{e_1, e_4\}, S_1 = (e_1, e_4), S_2 = (e_2, e_3, e_5, e_6)$,

$k = 1$, and $s = 4$.

● In the second iteration: We found two-dimensional invariant subspace $V = \langle e_2, e_3 \rangle$ of a set of matrices $\{T_i \}_{ i = 1}^2$. Therefore, $\mathcal{B} = \{e_1, e_4, e_3, e_5\}, S_1 = (e_1, e_4, e_3, e_5), S_2 = (e_2, e_6)$,

$k = 2$, and $s = 2$.

● In the third iteration: There is no one-dimensional invariant subspace of a set of matrices $\{T_i \}_{ i = 1}^2$. Therefore, $S = (e_1\; e_4\; e_3\; e_5\; e_2\; e_6)$, and the corresponding unitary matrix is

such that the set $U \Omega U^{*} = \{U A_i U^{*}\}_{i = 1}^2$ is $BT(2, 2, 2)$ where

3.2.   Algorithm $B$: Simultaneous block diagonalization of a set of $n \times n$ matrices $\{A_i \}_{ i=1}^N$.

(1) Input: the set $\Omega = \{A_i \}_{ i = 1}^N$.

(2) Construct the set $\Gamma = \Omega \cup \Omega^{*}$.

(3) Find a unitary matrix $U$ such that $U \Gamma U^{*}$ is $BT(n_1, ..., n_t)$ using Algorithm $A$.

(4) Output: a unitary matrix $U$.

Remark 3.6. Algorithm $B$ provides the finest block-diagonalization. Moreover, the number of the blocks equals the number the of the invariant subspaces, and the size of each block is $n_i \times n_i$, where $n_i$ is the dimension of the invariant subspace.

Example 3.7. The set of matrices $\Omega = \{A_i \}_{i = 1}^2$ admits simultaneous block diagonalization where

Applying Algorithm $B$ to the set $\Omega$ can be summarized as follows:

● Input: $\Gamma = \Omega \cup \Omega^{*}$.

● Initiation step:

We have $k = 0, \mathcal{B} = \phi, s = 7, T_1 = A_1, T_2 = A_2, T_3 = A_2^T, S_2 = I$.

● In the first iteration:

We found one-dimensional invariant subspace $V = \langle e_5 \rangle$ of a set of matrices $\{T_i \}_{ i = 1}^3$. Therefore, $\mathcal{B} = \{e_5\}, S_1 = (e_5), S_2 = (e_1, e_2, e_3, e_4, e_6, e_7)$,

$k = 1$, and $s = 6$.

● In the second iteration: We found two-dimensional invariant subspace $V = \langle e_4, e_5 \rangle$ of a set of matrices $\{T_i \}_{ i = 1}^3$. Therefore, $\mathcal{B} = \{e_5, e_4, e_6\}, S_1 = (e_5\; e_4\; e_6), S_2 = (e_1, e_2, e_3, e_7)$,

$k = 2$, and $s = 4$.

● In the third iteration: We found two-dimensional invariant subspace $V = \langle e_2, e_3 \rangle$ of a set of matrices $\{T_i \}_{ i = 1}^3$. Therefore, $\mathcal{B} = \{e_5, e_4, e_6, e_2, e_3\}, S_1 = (e_5\; e_4\; e_6\; e_2\; e_3), S_2 = (e_1, e_7)$,

$k = 3$, and $s = 2$.

● In the fourth iteration: There is no one-dimensional invariant subspace of a set of matrices $\{T_i \}_{ i = 1}^3$. Therefore, $S = (e_5\; e_4\; e_6\; e_2\; e_3\; e_1\; e_7)$, and the corresponding unitary matrix is

such that the set $U \Omega U^{*} = \{U A_i U^{*}\}_{i = 1}^2$ is $BD(1, 2, 2, 2)$ where

Example 3.8. The set of matrices $\Omega = \{A_i \}_{i = 1}^2$ admits simultaneous block diagonalization where

Similarly, applying Algorithm $B$ to the set $\Omega$ provides the matrix $S = (e_6\; e_5\; e_7\; e_1\; e_3\; e_2\; e_4)$. Therefore, the corresponding unitary matrix is

such that the set $U \Omega U^{*} = \{U A_i U^{*}\}_{i = 1}^2$ is $BD(2, 2, 3)$ where

Example 3.9. The set of matrices $\Omega = \{A_i \}_{i = 1}^3$ admits simultaneous block diagonalization where

Similarly, applying Algorithm $B$ to the set $\Omega$ provides the matrix $S = (e_1+e_5\; e_9\; e_3\; e_6\; e_8\; -e_7\; e_1-e_5, e_2\; e_4)$. Therefore, the corresponding unitary matrix is

such that the set $U \Omega U^{*} = \{U A_i U^{*}\}_{i = 1}^3$ is $BD(1, 1, 2, 2, 3)$ where

4.   Algorithm for simultaneous block diagonalization of a set of matrices based on a commuting matrix

This section focuses on an alternate approach for simultaneous block diagonalization of a set of $n \times n$ matrices $\{A_s\}_{ s = 1}^N$ by an invertible matrix that does not require finding the common invariant subspaces as Algorithm $B$ given in the previous section. Maehara et al. ^[9] introduced an algorithm for simultaneous block diagonalization of a set of matrices by a unitary matrix based on the eigenvalue decomposition of a Hermitian commuting matrix. Here, we extend their algorithm to be applicable for a non-Hermitian commuting matrix by considering its generalized eigen vectors. Moreover, a new characterization is presented by proving that the existence of a commuting matrix that possesses at least two distinct eigenvalues is the necessary and sufficient condition to guarantee the simultaneous block diagonalization by an invertible matrix.

Proposition 4.1. Let $V$ be a vector space, and let $T:V \rightarrow V$ be a linear operator. Let $\lambda_1, ..., \lambda_k$ be distinct eigenvalues of $T$. Then, each generalized eigenspace $G_{\lambda_i}(T)$ is $T$-invariant, and we have the direct sum decomposition

Lemma 4.2. Let $V$ be a vector space, and let $T:V \rightarrow V$, $L:V \rightarrow V$ be linear commuting operators. Let $\lambda_1, ..., \lambda_k$ be distinct eigenvalues of $T$. Then, each generalized eigenspace $G_{\lambda_i}(T)$ is $L$-invariant.

Proof. Let $V$ be a vector space and $\lambda_1, ..., \lambda_k$ be distinct eigenvalues of $T$ with the minimal polynomial $\mu(x) = (x-\lambda_1)^{n_1}(x-\lambda_2)^{n_2}...(x-\lambda_k)^{n_k}$. Then, we have the direct sum decomposition $V = G_{\lambda_1}(T)\oplus G_{\lambda_2}(T) \oplus...\oplus G_{\lambda_k}(T)$.

For each $i = 1, .., k$, let $x \in G_{\lambda_i}(T)$, and then $(T-\lambda_i I)^{n_i}x = 0$. Then, $(T-\lambda_i I)^{n_i}Lx = L(T-\lambda_i I)^{n_i}x = 0$. Hence, $L x \in G_{\lambda_i}(T)$.

Theorem 4.3. Let $\{A_s\}_{ s = 1}^N$ be a set of $n \times n$ matrices. Then, the set $\{A_s\}_{ s = 1}^N$ admits simultaneous block diagonalization by an invertible matrix $S$ if and only if the set $\{A_s\}_{ s = 1}^N$ commutes with a matrix $C$ that possesses two distinct eigenvalues.

Proof. $\Rightarrow$ Assume that the set $\{A_s\}_{ s = 1}^N$ admits simultaneous block diagonalization by the an invertible matrix $S$ such that

where the number of blocks $k\geq 2$, and the matrices $B_{s, 1}, B_{s, 2}, ..., B_{s, k}$ have sizes $n_1 \times n_1, n_2 \times n_2, ..., n_k \times n_k$, respectively, for all $s = 1, .., N$.

Now, define the matrix $C$ as

where $\lambda_1, \lambda_2, ..., \lambda_k$ are any distinct numbers.

Clearly, the matrix $C$ commutes with the set $\{A_s\}_{ s = 1}^N$. Moreover, it has the distinct eigenvalues $\lambda_1, \lambda_2, ..., \lambda_k$.

$\Leftarrow$ Assume that the set $\{A_s\}_{ s = 1}^N$ commutes with a matrix $C$ that posseses distinct eigenvalues $\lambda_1, \lambda_2, ..., \lambda_k$.

Using Proposition 4.1, one can use the generalized eigenspace $G_{\lambda_i}(C)$ of the matrix $C$ associated to these distinct eigenvalues to decompose the matrix $C$ as a direct sum of $k$ matrices. This can be achieved by restricting the matrix $C$ on the invariant subspaces $G_{\lambda_i}(C)$ as follows:

where

Using Lemma 4.2, one can restrict each matrix $A_s$ on the invariant subspaces $G_{\lambda_i}(C)$ to decompose the matrix $A_s$ as a direct sum of $k$ matrices as follows:

Remark 4.4. For a given set of $n \times n$ matrices $\{A_s\}_{ s = 1}^N$, if the set $\{A_s\}_{ s = 1}^N$ commutes only with the matrices having only one eigenvalue, then it does not admit a simultaneous block diagonalization by an invertible matrix.

Algorithm $C$:

(1) Input: the set $\Omega = \{A_s \}_{ s = 1}^N$.

(2) Construct the the following matrix:

(3) Compute the null space of the matrix $X$ and reshape the obtained vectors as $n \times n$ matrices. These matrices commute with all the matrices $\{A_s \}_{ s = 1}^N$.

(4) Choose a matrix $C$ from the obtained matrices that possesses two distinct eigenvalues.

(5) Find the distinct eigenvalues $\lambda_1, ..., \lambda_k$ of the matrix $C$ and the corresponding algebraic multiplicity $n_1, n_2, ..., n_k$.

(6) Find each generalized eigenspace $G_{\lambda_i}(C)$ of the matrix $C$ associated to the eigenvalue $\lambda_i$ by computing the null space of $(C-\lambda_i I)^{n_i}$.

(7) Construct the invertible matrix $S$ as

(8) Verify that

where the matrices $B_{s, 1}, B_{s, 2}, ..., B_{s, k}$ have sizes $n_1 \times n_1, n_2 \times n_2, ..., n_k \times n_k$, respectively, for all $s = 1, .., N$.

(9) Output: an invertible matrix $S$.

Remark 4.5. Algorithm $C$ provides the finest block-diagonalization if one chooses a matrix $C$ with maximum number of distinct eigenvalues. Moreover, the number of the blocks equals the number the of the distinct eigenvalues, and the size of each block is $n_i \times n_i$, where $n_i$ is the algebraic multiplicity of the eigenvalue $\lambda_i$.

Example 4.6. Consider the set of matrices $\Omega = \{A_i \}_{i = 1}^6$ where

The set $\Omega$ admits simultaneous block diagonalization by an invertible matrix. An invertible matrix can be obtained by applying algorithm $C$ to the set $\Omega$ as summarized below:

● A matrix $C$ that commutes with all the matrices $\{A_i \}_{ i = 1}^6$ can be obtained as

.

● The distinct eigenvalues of the matrix $C$ are $\lambda_1 = -1, \lambda_2 = 1$ with algebraic multiplicities $n_1 = 3, n_2 = 3$, respectively..

● The generalized eigenspaces of the matrix $C$ associated to the distinct eigenvalues are

● The invertible matrix $S = \big(G_{\lambda_1}(C), G_{\lambda_2}(C) \big)$ is

● The set $S^{-1} \Omega S = \{S^{-1} A_i S\}_{i = 1}^6$ contains block diagonal matrices where

5.   Conclusions

It is well known that a set of non-defective matrices can be simultaneously diagonalized if and only if the matrices commute. In the case of non-commuting matrices, the best that can be achieved is simultaneous block diagonalization. Both Algorithm B and the Maehara et al. ^[9] algorithm are applicable for simultaneous block diagonalization of a set of matrices by a unitary matrix. Algorithm C can be applied for block diagonalization by an invertible matrix when finding a unitary matrix is not possible. In case block diagonalization of a set of matrices is not possible by a unitary or an invertible matrix, then one may utilize block triangularization by Algorithm A. Algorithms A and B are based on the existence of invariant subspaces; however, Algorithm C is based on the existence of a commuting matrix which is not necessarily Hermitian, unlike the Maehara et al. algorithm.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgments

Ahmad Y. Al-Dweik and M. T. Mustafa would like to thank Qatar University for its support and excellent research facilities. R. Ghanam and G. Thompson are grateful to VCU Qatar and Qatar Foundation for their support.

Conflict of Interest

The authors declare that they have no conflicts of interest.

Appendix: Maple procedures