• PDF
• Cite
• Share
Article Contents  Article Contents

# Extended Krylov subspace methods for solving Sylvester and Stein tensor equations

• * Corresponding author: Smahane El-Halouy
• This paper deals with Sylvester and Stein tensor equations with low rank right hand sides. It proposes extended Krylov-like methods for solving Sylvester and Stein tensor equations. The expressions of residual norms are presented. To show the performance of the proposed approaches, some numerical experiments are given.

Mathematics Subject Classification: 2020 MSC.

 Citation: • • Figure 1.  Residual norm of the global and the extended global Arnoldi algorithms for example $3$ with the first set of matrices and $n = 225$

Figure 2.  Residual norm of the global and the extended global Arnoldi algorithms for example $3$ with the second set of matrices and $n = 225$

Figure 3.  Residual norm of the extended block and global Arnoldi algorithms for example $4$ with $n = 300$

 Algorithm 1 Extended Block Arnoldi (EBA) 1: Input: $n\times n$ matrix $A$, $n\times r$ matrix $B$ and an integer $m$. 2: Compute the QR decomposition of $\left[ B,A^{-1}B\right]=U_1\Lambda, \mathbb{U}_0=\left[ \;\right]$ for $j=1,\ldots ,m$ 3: Set $U_j^{(1)}$: first $r$ columns of $U_j$, $U_j^{(2)}$: second $r$ columns of $U_j$ 4: $\mathbb{U}_j=\left[ \mathbb{U}_{j-1},U_j\right]$; $W=\left[ AU_j^{(1)},A^{-1}U_j^{(2)}\right]$ 5: Orthogonalize $W$ w.r. to $\mathbb{U}_j$ to get $U_{j+1}$ i.e., 6: for $i=1,\ldots ,j$ $H_{i,j}=U_i^TW$; $W=W-U_iH_{i,j}$ 7: Compute the QR decomposition of $W=U_{j+1}H_{j+1,j}$
 Algorithm 2 Extended Global Arnoldi (EGA) 1: Input: $n\times n$ matrix $A$, $n\times r$ matrix $B$ and an integer $m$. 2: Compute the global QR decomposition of $\left[ B,A^{-1}B\right]=V_1(\Omega\otimes I_r)$; $\mathbb{V}_0=\left[ \;\right]$ for $j=1,\ldots ,m$ 3: Set $V_j^{(1)}$: first $r$ columns of $V_j$, $V_j^{(2)}$: second $r$ columns of $V_j$ 4: $\mathbb{V}_j=\left[ \mathbb{V}_{j-1},V_j\right]$; $U=\left[ AV_j^{(1)},A^{-1}V_j^{(2)}\right]$ 5: Orthogonalize $U$ w.r. to $\mathbb{V}_j$ to get $V_{j+1}$ i.e., 6: for $i=1,\ldots ,j$ $H_{i,j}=V_i^T\lozenge U$; $U=U-V_iH_{i,j}$ 7:     Compute the QR decomposition of $U=V_{j+1}(H_{j+1,j}\otimes I_r)$
 Algorithm 3 Sylvester Stein EBA Process 1: Coefficient matrices $A^{(i)}$, $i=1,\ldots ,N$, and the right hand side in low rank representation $\mathcal{B}=\left[B^{(1)},\ldots ,B^{(N)} \right]$. 2: Output: Approximate solutions, $\mathcal{X}_m$, to the equations (1) / (2). 3: Choose a tolerance $\epsilon>0$, a maximum number of iterations $itermax_i$, a step-size parameter $k$, set for $i=1,\ldots ,N$, $m_i=k$. 4: For $i=1,\ldots ,N$ do   For $m_i=k$, construct the orthonormal basis $\mathbb{V}_{m_i}$ and the restriction matrices $\mathbb{T}_{m_i}$ by EBA algorithm (1). 5: Compute $\mathcal{Y}_m$ / $\mathcal{Z}_m$ the solution of the low dimensional equation (12) / (13). 6: Compute the residual norms $r_m=\|\mathcal{R}_m\|$ as in theorem (4.1)/ $q_m=\|\mathcal{Q}_m\|$ as in theorem (4.2). 7: If $r_m\geq\epsilon$ / $q_m\geq\epsilon$, set for $i=1,\ldots ,N$, $m_i=mi+k$ and go to step 2. 8: The approximate solutions are given by $\mathcal{X}_m=\mathcal{Y}_m\times_1\mathbb{V}_{m_1}\ldots\times_N\mathbb{V}_{m_N}$ / $\mathcal{X}_m=\mathcal{Z}_m\times_1\mathbb{V}_{m_1}\ldots\times_N\mathbb{V}_{m_N}$.
 Algorithm 4 Sylvester / Stein EGA Process 1: Input: Coefficient matrices $A^{(i)}$, $i=1,\ldots ,N$, and the right hand side in low rank representation $\mathcal{B}=\left[B^{(1)},\ldots ,B^{(N)} \right]$. 2: Output: Approximate solutions, $\mathcal{X}_m$, to the equations (1) / (2). 3: Choose a tolerance $\epsilon>0$, a maximum number of iterations $itermax_i$, a step-size parameter $k$, set for $i=1,\ldots ,N$, $m_i=k$. 4: For $i=1,\ldots ,N$ do   For $m_i=k$, construct the orthonormal basis $\mathbb{W}_{m_i}$ and the restriction matrices $T_{m_i}$ by EGA algorithm (2). 5: Compute $\mathcal{Y}_m$ / $\mathcal{Z}_m$ the solutions of the low dimensional equations (16) / (19). 6: Compute the upper bound for the residual norms $r_m=\left( {\sum_{i=1}^N\|\mathcal{Y}_m\times_iT_{m_i+1,m_i}E_{m_i}^T\|^2}\right)^{1/2}$ as in theorem (4.4)/ $q_m=\left( \beta_1^2+\sum_{i=1}^N{\beta_2^{(i)}}^2+\sum_{i=1}^N{\beta_3^{(i)}}^2\right)^{1/2}$ as in theorem (4.5). 7: If $r_m\geq\epsilon$ / $q_m\geq\epsilon$, set for $i=1,\ldots ,N$, $m_i=mi+k$ and go to step 2. 8: The approximate solutions are given by $\mathcal{X}_m=(\mathcal{Y}_m\otimes\mathcal{I}_R)\times_1\mathbb{W}_{m_1}\ldots\times_N\mathbb{W}_{m_N}$ / $\mathcal{X}_m=(\mathcal{Z}_m\otimes\mathcal{I}_R)\times_1\mathbb{W}_{m_1}\ldots\times_N\mathbb{W}_{m_N}$.

Table 1.  Example $1$

 Algorithm $R$ $n$ Cycles $\mathcal{R}_m$ CPU Times (Second) Algorithm 3 $3$ $500$ $5$ $9.34\cdot 10^{-7}$ $16.1$ $5$ $300$ $4$ $4.03\cdot 10^{-7}$ $31.38$ Algorithm 4 $3$ $500$ $9$ $6.69\cdot 10^{-7}$ $4.57$ $5$ $300$ $8$ $4.34\cdot 10^{-7}$ $3.76$

Table 2.  Example 2

 Algorithm $R$ $n$ Cycles $\|\mathcal{R}_m\|$ $\|\mathcal{X}^{\ast}-\mathcal{X}_m\|$ CPU Times (Second) Algorithm 3 (case $1$) $3$ $300$ $2$ $4.08\cdot 10^{-9}$ $3.85\cdot 10^{-6}$ $1.58$ (case $2$) $3$ $300$ $2$ $1.26\cdot 10^{-8}$ $2.82\cdot 10^{-7}$ $1.53$ Algorithm 4 (case $1$) $3$ $300$ $6$ $5.82\cdot 10^{-8}$ $4.21\cdot 10^{-6}$ $1.46$ (case $2$) $3$ $300$ $7$ $5.39\cdot 10^{-8}$ $4.17\cdot 10^{-8}$ $2.22$ Algorithm 4 $5$ $500$ $10$ $4.54\cdot 10^{-7}$ $-$ $26.01$

Table 3.  Example 3, $N = 3$ and $R = 5$

 Matrices Algorithm $\|\mathcal{R}_m\|$ CPU Times (Second) fdm Algorithm $2$  $1.45\cdot 10^{-6}$ $56.0$ Algorithm 4 $3.6\cdot 10^{-9}$ $0.74$ poisson Algorithm $2$  $2.82\cdot 10^{-6}$ $14.5$ Algorithm 4 $2.24\cdot 10^{-9}$ $0.67$

Table 4.  Example $4$

 Algorithm $R$ $n$ Cycles $\mathcal{R}_m$ CPU Times (Second) Algorithm 3 $3$ $300$ $6$ $1.58\cdot 10^{-6}$ $27.49$ Algorithm 4 $3$ $300$ $4$ $1.56\cdot 10^{-7}$ $0.83$ $5$ $300$ $4$ $8.77\cdot 10^{-7}$ $0.85$
• Figures(3)

Tables(8)

## Article Metrics  DownLoad:  Full-Size Img  PowerPoint