\clearpage
\section{Shifted power-GMRES method}

This section discusses the approach used by the authors of \cite{SHEN2022126799} to combine the shifted power method with the fast shifted \texttt{GMRES} method to create a hybrid algorithm for solving complex PageRank problems with multiple damping factors. The goal of this combination is to create an efficient and reliable algorithm for solving these types of problems. The details of this approach and how it was implemented are described in the cited paper.

\subsection{Restarted GMRES method}

The Restarted GMRES method (hereafter referred to as GMRES in short) is a non-symmetric Krylov subspace solver based on the Arnoldi decomposition procedure, that the authors sketch in algorithm \ref{alg:arnoldi}

\begin{algorithm}[H]
\caption{Arnoldi}
\label{alg:arnoldi}
\begin{algorithmic}[1]
    \Require $A, v_0, m$
    \Ensure $V_m, H_m, v_{m+1}, h_{m+1,m}, \beta, j$
    \State Compute $\beta = \lVert v_0 \rVert$
    \State $v_1 = v_0/ \beta$
    \For {$j = 1:m$} 
        \State Compute $w = Av_j$
        \For {$i = 1:j$} 
            \State Compute $h_{i,j} = v_i^T w$
            \State Compute $w = w - h_{i,j} v_i$
        \EndFor
        \State $h_{j+1,j} = \lVert w_i \rVert$
        \If {$h_{j+1,j} = 0$} 
            \State $m = j$, 
            \State $v_{m+1} = 0$
            \State \textbf{break}
        \Else
            \State $v_{j+1} = w / h_{j+1,j}$
        \EndIf
    \EndFor
\end{algorithmic}
\end{algorithm}

% \noindent Where $A \in \R^{n\times n}$ and $v_0 \in \R ^{n \times 1}$ is the initial vector. After $m$ iterations, the Arnoldi procedure produces the orthogonal basis $V_m = [v_1, \dots, v_m]$ and the upper Hessenberg matrix $H_m \in \R^{m\times m}$, and the residual vector $v_{m+1} \in \R^{n \times 1}$ and the residual norm $h_{m+1,m} \in \R$. Starting from $v_0 = b - Ax_0$ with an initial guess $x_0$, after running $m$ steps of the algorithm \ref{alg:arnoldi}, the \texttt{GMRES} method produces the approximate solution $\tilde x$ of the linear system $Ax = b$ that minimizes the residual norm $\lVert b - Ax \rVert$ in the Krylov subspace of dimension $m$. \vspace*{0.4cm}
\noindent The Arnoldi procedure, which is used as the basis for the GMRES method, involves iteratively constructing an orthogonal basis $V_m = [v_1, \dots, v_m]$ and an upper Hessenberg matrix $H_m \in \R^{m\times m}$ from an initial vector $v_0 \in \R ^{n \times 1}$ and a matrix $A \in \R^{n\times n}$. After $m$ iterations, it also produces a residual vector $v_{m+1} \in \R^{n \times 1}$ and a residual norm $h_{m+1,m} \in \R$. The GMRES method then uses these to approximate the solution $\tilde x$ of the linear system $Ax = b$ that minimizes the residual norm $\lVert b - Ax \rVert$ in the Krylov subspace of dimension $m$. To do this, it starts with an initial guess $x_0$ and a residual vector $v_0 = b$, and runs $m$ steps of the Arnoldi procedure outlined in algorithm \ref{alg:arnoldi}. \vspace*{0.4cm}

% \noindent We know that the accuracy of the approximate solution $\tilde x$ of \texttt{GMRES} depends heavily on the dimension $m$ of the search space. The authors in \cite{SHEN2022126799} propose to use the \texttt{GMRES} method as a preconditioner for the shifted power method presented in the previous section. The core idea of the method is to run standard GMRES on a seed system and to approximate the other solutions as by products. The theoretical basis is the shift-invariance property of the Krylov subspace that enables us to use only one Krylov subspace for all the shifted systems, provided that the residual vectors are collinear to one other. The algorithm proposed by the authors is presented in Algorithm \ref{alg:shifted_GMRES}.
% \noindent The accuracy of the approximate solution $\tilde x$ produced by the GMRES method depends significantly on the dimension $m$ of the search space used. In order to improve the efficiency of this method for solving the PageRank problem with multiple damping factors, the authors of \cite{SHEN2022126799} propose using a shifted-GMRES algorithm as an extension of the standard GMRES one, that can use the same Krylov subspace to solve the whole sequence of shifted linear systems simultaneously. The idea behind this approach is to run standard GMRES on a seed system and approximate the solutions for the other systems as byproducts. This is possible due to the shift-invariance property of the Krylov subspace, which allows the use of a single Krylov subspace for all the shifted systems as long as the residual vectors are collinear with one another. The authors present a specific algorithm for implementing this approach in Algorithm \ref{alg:shifted_GMRES}.

% \noindent The shifted GMRES method proposed by Frommer and Glassner is an extension of the standard GMRES algorithm that can use the same Krylov subspace to solve the whole  sequence of shifted linear systems
% \begin{equation*}
%     (A + \sigma_i I)x_i = b \qquad (\sigma_i \in \R, b_i \in \R^{n\times 1}, i = \dots s) 
% \end{equation*}
% simultaneously.  The core idea of the method is to run standard GMRES on a seed system and to approximate the other solutions as by products. The theoretical basis is the shift-invariance property of the Krylov subspace.  That enables us to use only one Krylov subspace to solve the whole sequence at once provided that the initial residual vectors of the shifted linear systems are collinear one another. However, after $m$ steps of Algorithm \ref{alg:arnoldi} the residuals vectors $r_m^i ~ (i = 2, 3, \dots, s)$  of the additional systems and the residual vector $r_m^1$ of the seed system will in general lose the collinearity, and consequently the shift-invariance property of the search space will not be maintained
% \begin{equation*}
%     \mathcal{K}_m(A + \sigma_i I, r_m^i) \neq \mathcal{K}_m(A + \sigma_1 I, r_m^1) \qquad i = 2, 3, \dots, s
% \end{equation*}
% This means: upon restarting the GMRES algorithm, the mutual collinearity of all the residual vectors needs to be enforced explicitly. More precisely, let us assume that the initial residual vectors are collinear, i.e. $r_0^i = \gamma_0^i r_0^1 ~~ (\gamma_0^i \in \R, ~ i = 2, 3, \dots, s)$  where $r_0^i$ is the initial residual vector corresponding
% to the $i$-th system

\noindent The modified GMRES algorithm represents an extension of the classical GMRES method, allowing for the solution of a series of shifted linear systems using the same Krylov subspace.
\begin{equation*}
    (A + \sigma_i I)x_i = b \qquad (\sigma_i \in \R, b_i \in \R^{n\times 1}, i = \dots s) 
\end{equation*}
The central concept behind the shifted GMRES method is the utilization of the standard GMRES algorithm on a seed system, with the solutions to the other systems being approximated as byproducts. This approach is made possible through the shift-invariance property of the Krylov subspace, which allows for the utilization of a single Krylov subspace to solve the entire sequence of systems provided that the initial residual vectors of the shifted linear systems are collinear with one another. However, after a certain number of steps of the Arnoldi algorithm, the residual vectors of the additional systems and the residual vector of the seed system will typically lose collinearity, resulting in the loss of the shift-invariance property of the search space.
\begin{equation*}
    \mathcal{K}_m(A + \sigma_i I, r_m^i) \neq \mathcal{K}_m(A + \sigma_1 I, r_m^1) \qquad i = 2, 3, \dots, s
\end{equation*}
This implies that the mutual collinearity of all the residual vectors must be enforced explicitly upon restarting the GMRES algorithm. Specifically, if we assume that the initial residual vectors are collinear, meaning $r_0^i = \gamma_0^i r_0^1$ for $i = 2, 3, \dots, s$ where $r_0^i$ is the initial residual vector corresponding to the $i$-th system and $\gamma_0^i \in \R$, then the shift-invariance property of the Krylov subspace can be maintained.

\begin{algorithm}[H]
\caption{Shifted GMRES}
\label{alg:shifted_GMRES}
\begin{algorithmic}[1]
    \Require $\tilde P, v, m, \alpha_i, maxit, x_0^i ~~ (1 \leq i \leq s)$
    \Ensure $x^i, res_i ~~(1 \leq i \leq s), mv$
    \State Set $_0^i = \frac{1-\alpha_i}{\alpha_i} v - \Big(\frac{1}{\alpha_i} I - \tilde P \Big) x_0^i$, iter = 1
    \State Set $res_i = \alpha_i \lVert v \rVert ~~ (1 \leq i \leq s)$
    \State Set mv = 0
    \While {$\max (res_i) \geq \tau~~ \&\& ~~ iter \leq maxit$}
        \State Find $k$ that satisfies $res_k = \max (res_i)$
        \State Compute $\gamma^i = \frac{res_i \alpha_k}{res_k \alpha_i}$ for all $i \neq k$
        \State Run Arnoldi by $ [V_m, \bar H_m^k, v_{m+1}, \bar h_{m+1,m}, \beta, j] = Arnoldi(\frac{1}{\alpha_k}I - \tilde P, r_0^k, m)$
        \State Set $mv = mv + j$
        \State Compute $y_k$, the minimizer of $\lVert \beta e_1 - \bar H_m^k y_k \rVert_2$
        \State Compute $x^k = x_0^k + V_m y_k$
        \State Compute $res_k = \alpha_k \lVert \beta e_1 - \bar H_m^k y^k \rVert$
        \For {i = 1, 2, \dots , k-1, k+1, \dots , s}
            \If {$res_i \geq \tau$} 
                \State Set $\bar H_m^i = \bar H_m^k + \Big( \frac{1-\alpha_i}{\alpha_i} - \frac{1-\alpha_k}{\alpha_k} \Big) I_m$
                \State Solve $y_i$ and $\gamma_i$ from $\begin{bmatrix} \bar H_m^i & z \end{bmatrix} \begin{bmatrix} y^i \\ \gamma^i \end{bmatrix} = \gamma^i \beta e_1$
                \State Set $x^i = x_0^i + V_m y^i$
                \State Set $res_i = \frac{\alpha_i}{\alpha_k} \gamma_k^i res_k$
            \EndIf
        \EndFor
        \State Set $iter = iter + 1$
        \State Set $x_0^i = x^i$
    \EndWhile
\end{algorithmic}
\end{algorithm}

\noindent Where $z = \beta e_1 - H_m^1 y_m^1$. In line 15, by solving this small size system, we can obtain the vector $y_m^i$ and scalar $\gamma_m^i$ that ensures the collinearity of the shifted results.

\paragraph{Problems:} There have been significant issues with the implementation of the shifted GMRES algorithm. A key component of the algorithm is the use of the seed choosing strategy outlined in \cite{SHEN2022126799}. However, during testing, it was observed that after the second iteration, the value of $k$ remained constant and the $res$ vector did not change, resulting in a stall situation where the program runs without updating values until the maximum number of iterations is reached. The cause of this problem is currently under investigation. Although a notebook containing the code for the algorithm has been included in the GitHub repository for completeness, it is not functional at this time. It is believed that the issue may be related to a misunderstanding of the algorithm as presented in the pseudo-code, but this has not yet been determined. As a result, there will be no test results for this algorithm in the following section.