Commit e-Offprint provided by Springer.

[GMRES_For_Journal.git] / GMRES_Journal.tex

diff --git a/GMRES_Journal.tex b/GMRES_Journal.tex
index a4a026adef83ca556aee81f8b7641b40dc7f5795..3add2567e3fcbc522d7296230a2f6da9121bca8b 100644 (file)
--- a/GMRES_Journal.tex
+++ b/GMRES_Journal.tex
@@ -1,5 +1,5 @@
 \documentclass[11pt]{article}
 \documentclass[11pt]{article}

-%\documentclass{acmconf}
+\usepackage{multicol}

 
 \usepackage[paper=a4paper,dvips,top=1.5cm,left=1.5cm,right=1.5cm,foot=1cm,bottom=1.5cm]{geometry}
 \usepackage{times}
 
 \usepackage[paper=a4paper,dvips,top=1.5cm,left=1.5cm,right=1.5cm,foot=1cm,bottom=1.5cm]{geometry}
 \usepackage{times}


@@ -18,28 +18,38 @@
 \usepackage{url}
 \usepackage{mdwlist}
 \usepackage{multirow}
 \usepackage{url}
 \usepackage{mdwlist}
 \usepackage{multirow}

-\usepackage{color}

 
 \date{}
 
 \title{Parallel sparse linear solver with GMRES method using minimization techniques of communications for GPU clusters}
 
 \author{
 
 \date{}
 
 \title{Parallel sparse linear solver with GMRES method using minimization techniques of communications for GPU clusters}
 
 \author{

-\textsc{Jacques M. Bahi}
+\textsc{Lilia Ziane Khodja}

 \qquad
 \textsc{Rapha\"el Couturier}\thanks{Contact author}
 \qquad
 \qquad
 \textsc{Rapha\"el Couturier}\thanks{Contact author}
 \qquad

-\textsc{Lilia Ziane Khodja}
+\textsc{Arnaud Giersch}
+\qquad
+\textsc{Jacques M. Bahi}

 \mbox{}\\ %
 \\
 FEMTO-ST Institute, University of Franche-Comte\\
 IUT Belfort-Montb\'eliard\\
 \mbox{}\\ %
 \\
 FEMTO-ST Institute, University of Franche-Comte\\
 IUT Belfort-Montb\'eliard\\

-Rue Engel Gros, BP 527, 90016 Belfort, \underline{France}\\
+19 Av. du Maréchal Juin, BP 527, 90016 Belfort, France\\

 \mbox{}\\ %
 \normalsize
 \mbox{}\\ %
 \normalsize

-\{\texttt{jacques.bahi},~\texttt{raphael.couturier},~\texttt{lilia.ziane\_khoja}\}\texttt{@univ-fcomte.fr}
+\{\texttt{lilia.ziane\_khoja},~\texttt{raphael.couturier},~\texttt{arnaud.giersch},~\texttt{jacques.bahi}\}\texttt{@univ-fcomte.fr}

 }
 
 }
 

+\newcommand{\BW}{\mathit{bw}}
+\newcommand{\Iter}{\mathit{iter}}
+\newcommand{\Max}{\mathit{max}}
+\newcommand{\Offset}{\mathit{offset}}
+\newcommand{\Prec}{\mathit{prec}}
+\newcommand{\Ratio}{\mathit{Ratio}}
+\newcommand{\Time}[1]{\mathit{Time}_\mathit{#1}}
+\newcommand{\Wavg}{W_{\mathit{avg}}}
+

 \begin{document}
 
 \maketitle
 \begin{document}
 
 \maketitle


@@ -47,7 +57,7 @@ Rue Engel Gros, BP 527, 90016 Belfort, \underline{France}\\
 \begin{abstract}
 In this paper, we aim at exploiting the power computing of a GPU cluster for solving large sparse
 linear systems. We implement the parallel algorithm of the GMRES iterative method using the CUDA
 \begin{abstract}
 In this paper, we aim at exploiting the power computing of a GPU cluster for solving large sparse
 linear systems. We implement the parallel algorithm of the GMRES iterative method using the CUDA

-programming language and the MPI parallel environment. The experiments shows that a GPU cluster
+programming language and the MPI parallel environment. The experiments show that a GPU cluster

 is more efficient than a CPU cluster. In order to optimize the performances, we use a compressed
 storage format for the sparse vectors and the hypergraph partitioning. These solutions improve
 the spatial and temporal localization of the shared data between the computing nodes of the GPU
 is more efficient than a CPU cluster. In order to optimize the performances, we use a compressed
 storage format for the sparse vectors and the hypergraph partitioning. These solutions improve
 the spatial and temporal localization of the shared data between the computing nodes of the GPU


@@ -125,9 +135,8 @@ This method firstly uses the Reverse Cuthill-McKee reordering to reduce the tota
 volume. In addition, the performances of the parallel FEM algorithm are improved by overlapping
 the communication with computation. 
 
 volume. In addition, the performances of the parallel FEM algorithm are improved by overlapping
 the communication with computation. 
 

-%%% MODIF %%%
-\textcolor{red}{ \bf Our main contribution in this work is to show the difficulties to implement the GMRES method for solving sparse linear systems on a cluster of GPUs. First, we show the main key points of the parallel GMRES algorithm on a GPU cluster. Then, we discuss the improvements of the algorithm which are mainly performed on the sparse matrix-vector multiplication when the matrix is distributed on several GPUs. In fact, on a cluster of GPUs the influence of the communications is greater than on clusters of CPUs due to the CPU/GPU communications between two GPUs that are not on the same machines. We propose to perform a hypergraph partitioning on the problem to be solved, then we reorder the matrix columns according to the partitioning scheme, and we use a compressed format for storing the vectors in such a way to minimize the communication overheads between two GPUs.}
-%%% END %%%
+ Our main contribution in this work is to show the difficulties of implementing the GMRES method to solve sparse linear systems on a cluster of GPUs. First, we show the main key points of the parallel GMRES algorithm on a GPU cluster. Then, we discuss the improvements of the algorithm which are mainly performed on the sparse matrix-vector multiplication when the matrix is distributed on several GPUs. In fact, on a cluster of GPUs the influence of the communications is greater than on clusters of CPUs due to the CPU/GPU communications between two GPUs that are not on the same machines. We propose to perform a hypergraph partitioning on the problem to be solved, then we reorder the matrix columns according to the partitioning scheme, and we use a compressed format to store the vectors in order to minimize the communication overheads between two GPUs.
+

 
 %%--------------------%%
 %%      SECTION 3     %%
 
 %%--------------------%%
 %%      SECTION 3     %%


@@ -142,7 +151,7 @@ read-only \emph{constant} and \emph{texture} caches per multiprocessor and a rea
 \emph{global memory} shared by all its multiprocessors. The new architectures (Fermi, Kepler,
 etc) have also L1 and L2 caches to improve the accesses to the global memory.
 
 \emph{global memory} shared by all its multiprocessors. The new architectures (Fermi, Kepler,
 etc) have also L1 and L2 caches to improve the accesses to the global memory.
 

-NVIDIA has released the CUDA platform (Compute Unified Device Architecture)~\cite{Nvi10}
+NVIDIA has released the CUDA platform (Compute Unified Device Architecture)~\cite{ref19}

 which provides a high level GPGPU-based programming language (General-Purpose computing
 on GPUs), allowing to program GPUs for general purpose computations. In CUDA programming
 environment, all data-parallel and compute intensive portions of an application running
 which provides a high level GPGPU-based programming language (General-Purpose computing
 on GPUs), allowing to program GPUs for general purpose computations. In CUDA programming
 environment, all data-parallel and compute intensive portions of an application running


@@ -178,7 +187,7 @@ the transfer of data between the GPU and its host.
 %%--------------------%%
 \section{{GMRES} method}
 \label{sec:04}
 %%--------------------%%
 \section{{GMRES} method}
 \label{sec:04}

-%%% MODIF %%%
+

 The generalized minimal residual method (GMRES) is an iterative method designed by Saad and Schultz in 1986~\cite{Saa86}. It is a generalization of the minimal residual method (MNRES)~\cite{Pai75} to deal with asymmetric and non Hermitian problems and indefinite symmetric problems. 
 
 Let us consider the following sparse linear system of $n$ equations:
 The generalized minimal residual method (GMRES) is an iterative method designed by Saad and Schultz in 1986~\cite{Saa86}. It is a generalization of the minimal residual method (MNRES)~\cite{Pai75} to deal with asymmetric and non Hermitian problems and indefinite symmetric problems. 
 
 Let us consider the following sparse linear system of $n$ equations:


@@ -200,20 +209,23 @@ such that the Petrov-Galerkin condition is satisfied:
 Algorithm~\ref{alg:01} illustrates the main key points of the GMRES method with restarts. The linear system to be solved in this algorithm is left-preconditioned where $M$ is the preconditioning matrix. The Arnoldi process~\cite{Arn51} is used (from line~$7$ to line~$17$ of algorithm~\ref{alg:01}) to construct an orthonormal basis $V_m$ and a Hessenberg matrix $\bar{H}_m$ of order $(m+1)\times m$ such that $m\ll n$. Then, the least-squares problem is solved (line~$18$) to find the vector $y\in\mathbb{R}^m$ which minimizes the residual. Finally, the solution $x_m$ is computed in the Krylov sub-space spanned by $V_m$ (line~$19$). In practice, the GMRES algorithm stops when the Euclidean norm of the residual is small enough and/or the maximum number of iterations is reached.
 %%% END %%%
 
 Algorithm~\ref{alg:01} illustrates the main key points of the GMRES method with restarts. The linear system to be solved in this algorithm is left-preconditioned where $M$ is the preconditioning matrix. The Arnoldi process~\cite{Arn51} is used (from line~$7$ to line~$17$ of algorithm~\ref{alg:01}) to construct an orthonormal basis $V_m$ and a Hessenberg matrix $\bar{H}_m$ of order $(m+1)\times m$ such that $m\ll n$. Then, the least-squares problem is solved (line~$18$) to find the vector $y\in\mathbb{R}^m$ which minimizes the residual. Finally, the solution $x_m$ is computed in the Krylov sub-space spanned by $V_m$ (line~$19$). In practice, the GMRES algorithm stops when the Euclidean norm of the residual is small enough and/or the maximum number of iterations is reached.
 %%% END %%%
 

-\begin{algorithm}[!h]
+\begin{algorithm}[!ht]
+  \newcommand{\Convergence}{\mathit{convergence}}
+  \newcommand{\False}{\mathit{false}}
+  \newcommand{\True}{\mathit{true}}

   %\SetAlgoLined
   %\SetAlgoLined

-  \Entree{$A$ (matrix), $b$ (vector), $M$ (preconditioning matrix),
-$x_{0}$ (initial guess), $\varepsilon$ (tolerance threshold), $max$ (maximum number of iterations),
+  \KwIn{$A$ (matrix), $b$ (vector), $M$ (preconditioning matrix),
+$x_{0}$ (initial guess), $\varepsilon$ (tolerance threshold), $\Max$ (maximum number of iterations),

 $m$ (number of iterations of the Arnoldi process)}
 $m$ (number of iterations of the Arnoldi process)}

-  \Sortie{$x$ (solution vector)}
+  \KwOut{$x$ (solution vector)}

   \BlankLine
   $r_{0} \leftarrow M^{-1}(b - Ax_{0})$\;
   $\beta \leftarrow \|r_{0}\|_{2}$\;
   $\alpha \leftarrow \|M^{-1}b\|_{2}$\;
   \BlankLine
   $r_{0} \leftarrow M^{-1}(b - Ax_{0})$\;
   $\beta \leftarrow \|r_{0}\|_{2}$\;
   $\alpha \leftarrow \|M^{-1}b\|_{2}$\;

-  $convergence \leftarrow false$\;
+  $\Convergence \leftarrow \False$\;

   $k \leftarrow 1$\;
   \BlankLine
   $k \leftarrow 1$\;
   \BlankLine

-  \While{$(\neg convergence)$}{
+  \While{$(\neg \Convergence)$}{

     $v_{1} \leftarrow r_{0} / \beta$\;
     \For{$j=1$ {\bf to} $m$}{
       $w_{j} \leftarrow M^{-1}Av_{j}$\;
     $v_{1} \leftarrow r_{0} / \beta$\;
     \For{$j=1$ {\bf to} $m$}{
       $w_{j} \leftarrow M^{-1}Av_{j}$\;


@@ -226,14 +238,14 @@ $m$ (number of iterations of the Arnoldi process)}
     }
     \BlankLine
     Put $V_{m}=\{v_{j}\}_{1\leq j \leq m}$ and $\bar{H}_{m}=(h_{i,j})$ Hessenberg matrix of order $(m+1)\times m$\;
     }
     \BlankLine
     Put $V_{m}=\{v_{j}\}_{1\leq j \leq m}$ and $\bar{H}_{m}=(h_{i,j})$ Hessenberg matrix of order $(m+1)\times m$\;

-    Solve the least-squares problem of size $m$: $\underset{y\in\mathbb{R}^{m}}{min}\|\beta e_{1}-\bar{H}_{m}y\|_{2}$\;
+    Solve the least-squares problem of size $m$: $\displaystyle\min_{y\in\mathbb{R}^{m}} \|\beta e_{1}-\bar{H}_{m}y\|_{2}$\;

     \BlankLine
     $x_{m} \leftarrow x_{0} + V_{m}y$\;
     $r_{m} \leftarrow M^{-1}(b-Ax_{m})$\;
     $\beta \leftarrow \|r_{m}\|_{2}$\;
     \BlankLine
     \BlankLine
     $x_{m} \leftarrow x_{0} + V_{m}y$\;
     $r_{m} \leftarrow M^{-1}(b-Ax_{m})$\;
     $\beta \leftarrow \|r_{m}\|_{2}$\;
     \BlankLine

-    \eIf{$(\frac{\beta}{\alpha}<\varepsilon)$ {\bf or} $(k\geq max)$}{
-      $convergence \leftarrow true$\;
+    \eIf{$(\frac{\beta}{\alpha}<\varepsilon)$ {\bf or} $(k\geq \Max)$}{
+      $\Convergence \leftarrow \True$\;

     }{
       $x_{0} \leftarrow x_{m}$\;
       $r_{0} \leftarrow r_{m}$\;
     }{
       $x_{0} \leftarrow x_{m}$\;
       $r_{0} \leftarrow r_{m}$\;


@@ -302,7 +314,7 @@ The most important operation in the GMRES method is the sparse matrix-vector mul
 It is quite expensive for large size matrices in terms of execution time and memory space. In
 addition, it performs irregular memory accesses to read the nonzero values of the sparse matrix,
 implying non coalescent accesses to the GPU global memory which slow down the performances of
 It is quite expensive for large size matrices in terms of execution time and memory space. In
 addition, it performs irregular memory accesses to read the nonzero values of the sparse matrix,
 implying non coalescent accesses to the GPU global memory which slow down the performances of

-the GPUs. So we use the HYB kernel developed and optimized by Nvidia~\cite{CUSP} which gives on
+the GPUs. So we use the HYB kernel developed and optimized by NVIDIA~\cite{CUSP} which gives on

 average the best performance in sparse matrix-vector multiplications on GPUs~\cite{Bel09}. The
 HYB (Hybrid) storage format is the combination of two sparse storage formats: Ellpack format
 (ELL) and Coordinate format (COO). It stores a typical number of nonzero values per row in ELL
 average the best performance in sparse matrix-vector multiplications on GPUs~\cite{Bel09}. The
 HYB (Hybrid) storage format is the combination of two sparse storage formats: Ellpack format
 (ELL) and Coordinate format (COO). It stores a typical number of nonzero values per row in ELL


@@ -318,10 +330,10 @@ requires the vector elements of its neighboring nodes corresponding to the colum
 which its local sub-matrix has nonzero values. Consequently, each computing node manages a global
 vector composed of a local vector of size $\frac{n}{p}$ and a shared vector of size $S$:
 \begin{equation}
 which its local sub-matrix has nonzero values. Consequently, each computing node manages a global
 vector composed of a local vector of size $\frac{n}{p}$ and a shared vector of size $S$:
 \begin{equation}

-  S = bw - \frac{n}{p},
+  S = \BW - \frac{n}{p},

 \label{eq:11}
 \end{equation}
 \label{eq:11}
 \end{equation}

-where $\frac{n}{p}$ is the size of the local vector and $bw$ is the bandwidth of the local sparse
+where $\frac{n}{p}$ is the size of the local vector and $\BW$ is the bandwidth of the local sparse

 sub-matrix which represents the number of columns between the minimum and the maximum column indices
 (see Figure~\ref{fig:01}). In order to improve memory accesses, we use the texture memory to
 cache elements of the global vector.
 sub-matrix which represents the number of columns between the minimum and the maximum column indices
 (see Figure~\ref{fig:01}). In order to improve memory accesses, we use the texture memory to
 cache elements of the global vector.


@@ -345,20 +357,20 @@ is composed of $240$ cores running at $1.3$ GHz and has $4$~GB of global memory
 of $102$~GB/s. The GPU is connected to the CPU via a PCI-Express 16x Gen2.0 with a throughput of $8$~GB/s.
 Figure~\ref{fig:02} shows the general scheme of the GPU cluster.
 
 of $102$~GB/s. The GPU is connected to the CPU via a PCI-Express 16x Gen2.0 with a throughput of $8$~GB/s.
 Figure~\ref{fig:02} shows the general scheme of the GPU cluster.
 

-\begin{figure}[!h]
+\begin{figure}[!ht]

 \centering
   \includegraphics[width=80mm,keepaspectratio]{Figures/clusterGPU}
 \caption{A cluster composed of six machines, each equipped with two Tesla C1060 GPUs}
 \label{fig:02}
 \end{figure}
 
 \centering
   \includegraphics[width=80mm,keepaspectratio]{Figures/clusterGPU}
 \caption{A cluster composed of six machines, each equipped with two Tesla C1060 GPUs}
 \label{fig:02}
 \end{figure}
 

-Linux cluster version 2.6.18 OS is installed on the six machines. The C programming language is used for
+Scientific Linux 5.10, with Linux version 2.6.18, is installed on the six machines. The C programming language is used for

 coding the GMRES algorithm on both the CPU and the GPU versions. CUDA version 4.0~\cite{ref19} is used for programming 
 the GPUs, using CUBLAS library~\cite{ref37} to deal with the functions operating on vectors. Finally, MPI routines
 of OpenMPI 1.3.3 are used to carry out the communication between the CPU cores.
 
 The experiments are done on linear systems associated to sparse matrices chosen from the Davis collection of the
 coding the GMRES algorithm on both the CPU and the GPU versions. CUDA version 4.0~\cite{ref19} is used for programming 
 the GPUs, using CUBLAS library~\cite{ref37} to deal with the functions operating on vectors. Finally, MPI routines
 of OpenMPI 1.3.3 are used to carry out the communication between the CPU cores.
 
 The experiments are done on linear systems associated to sparse matrices chosen from the Davis collection of the

-university of Florida~\cite{Dav97}. They are matrices arising in real-world applications. Table~\ref{tab:01} shows
+University of Florida~\cite{Dav97}. They are matrices arising in real-world applications. Table~\ref{tab:01} shows

 the main characteristics of these sparse matrices and Figure~\ref{fig:03} shows their sparse structures. For
 each matrix, we give the number of rows (column~$3$ in Table~\ref{tab:01}), the number of nonzero values (column~$4$)
 and the bandwidth (column~$5$).
 the main characteristics of these sparse matrices and Figure~\ref{fig:03} shows their sparse structures. For
 each matrix, we give the number of rows (column~$3$ in Table~\ref{tab:01}), the number of nonzero values (column~$4$)
 and the bandwidth (column~$5$).


@@ -395,21 +407,19 @@ Matrix type                 & Name              & \# Rows   & \# Nonzeros  & Ban
 
 All the experiments are performed on double-precision data. The parameters of the parallel
 GMRES algorithm are as follows: the tolerance threshold $\varepsilon=10^{-12}$, the maximum
 
 All the experiments are performed on double-precision data. The parameters of the parallel
 GMRES algorithm are as follows: the tolerance threshold $\varepsilon=10^{-12}$, the maximum

-number of iterations $max=500$, the Arnoldi process is limited to $m=16$ iterations, the elements 
+number of iterations $\Max=500$, the Arnoldi process is limited to $m=16$ iterations, the elements 

 of the guess solution $x_0$ is initialized to $0$ and those of the right-hand side vector are
 of the guess solution $x_0$ is initialized to $0$ and those of the right-hand side vector are

-initialized to $1$. For simplicity sake, we chose the matrix preconditioning $M$ as the
+initialized to $1$. For simplicity's sake, we chose the matrix preconditioning $M$ as the

 main diagonal of the sparse matrix $A$. Indeed, it allows us to easily compute the required inverse
 matrix $M^{-1}$ and it provides relatively good preconditioning in most cases. Finally, we set 
 the size of a thread-block in GPUs to $512$ threads. 
 main diagonal of the sparse matrix $A$. Indeed, it allows us to easily compute the required inverse
 matrix $M^{-1}$ and it provides relatively good preconditioning in most cases. Finally, we set 
 the size of a thread-block in GPUs to $512$ threads. 

-%%% MODIF %%%
-\textcolor{red}{\bf It would be noted that the same optimizations done on the GPU version of the parallel GMRES algorithm are performed on the CPU version.}
-%%% END %%%
+ It should be noted that the same optimizations are performed on the CPU version and on the GPU version of the parallel GMRES algorithm.

 
 

-\begin{table}[!h]
+\begin{table}[!ht]

 \begin{center}
 \begin{tabular}{|c|c|c|c|c|c|c|} 
 \hline
 \begin{center}
 \begin{tabular}{|c|c|c|c|c|c|c|} 
 \hline

-Matrix            & $Time_{cpu}$  & $Time_{gpu}$  & $\tau$  & \# $iter$ & $prec$ & $\Delta$   \\ \hline \hline
+Matrix            & $\Time{cpu}$  & $\Time{gpu}$  & $\tau$  & \# $\Iter$ & $\Prec$ & $\Delta$   \\ \hline \hline

 2cubes\_sphere    & 0.234 s     & 0.124 s      & 1.88  & 21     & 2.10e-14      & 3.47e-18 \\
 ecology2          & 0.076 s     & 0.035 s      & 2.15  & 21     & 4.30e-13      & 4.38e-15 \\
 finan512          & 0.073 s     & 0.052 s      & 1.40  & 17     & 3.21e-12      & 5.00e-16 \\
 2cubes\_sphere    & 0.234 s     & 0.124 s      & 1.88  & 21     & 2.10e-14      & 3.47e-18 \\
 ecology2          & 0.076 s     & 0.035 s      & 2.15  & 21     & 4.30e-13      & 4.38e-15 \\
 finan512          & 0.073 s     & 0.052 s      & 1.40  & 17     & 3.21e-12      & 5.00e-16 \\


@@ -431,24 +441,24 @@ torso3            & 4.242 s     & 2.030 s      & 2.09  & 175    & 2.69e-10
 In Table~\ref{tab:02}, we give the performances of the parallel GMRES algorithm for solving the linear
 systems associated to the sparse matrices shown in Table~\ref{tab:01}. The second and third columns show
 the execution times in seconds obtained on a cluster of 24 CPU cores and on a cluster of 12 GPUs, respectively.
 In Table~\ref{tab:02}, we give the performances of the parallel GMRES algorithm for solving the linear
 systems associated to the sparse matrices shown in Table~\ref{tab:01}. The second and third columns show
 the execution times in seconds obtained on a cluster of 24 CPU cores and on a cluster of 12 GPUs, respectively.

-The fourth column shows the ratio $\tau$ between the CPU time $Time_{cpu}$ and the GPU time $Time_{gpu}$ that
+The fourth column shows the ratio $\tau$ between the CPU time $\Time{cpu}$ and the GPU time $\Time{gpu}$ that

 is computed as follows:
 \begin{equation}
 is computed as follows:
 \begin{equation}

-  \tau = \frac{Time_{cpu}}{Time_{gpu}}.
+  \tau = \frac{\Time{cpu}}{\Time{gpu}}.

 \end{equation}
 From these ratios, we can notice that the use of many GPUs is not interesting to solve small sparse linear
 systems. Solving these sparse linear systems on a cluster of 12 GPUs is as fast as on a cluster of 24 CPU
 cores. Indeed, the small sizes of the sparse matrices do not allow to maximize the utilization of the GPU
 cores of the cluster. The fifth, sixth and seventh columns show, respectively, the number of iterations performed
 \end{equation}
 From these ratios, we can notice that the use of many GPUs is not interesting to solve small sparse linear
 systems. Solving these sparse linear systems on a cluster of 12 GPUs is as fast as on a cluster of 24 CPU
 cores. Indeed, the small sizes of the sparse matrices do not allow to maximize the utilization of the GPU
 cores of the cluster. The fifth, sixth and seventh columns show, respectively, the number of iterations performed

-by the parallel GMRES algorithm to converge, the precision of the solution, $prec$, computed on the GPU
+by the parallel GMRES algorithm to converge, the precision of the solution, $\Prec$, computed on the GPU

 cluster and the difference, $\Delta$, between the solutions computed on the GPU and the GPU clusters. The
 last two parameters are used to validate the results obtained on the GPU cluster and they are computed as
 follows:
 \begin{equation}
 cluster and the difference, $\Delta$, between the solutions computed on the GPU and the GPU clusters. The
 last two parameters are used to validate the results obtained on the GPU cluster and they are computed as
 follows:
 \begin{equation}

-\begin{array}{c}
-  prec = \|M^{-1}(b-Ax^{gpu})\|_{\infty}, \\
-  \Delta = \|x^{cpu}-x^{gpu}\|_{\infty},
-\end{array}
+  \begin{aligned}
+    \Prec &= \|M^{-1}(b-Ax^{gpu})\|_{\infty}, \\
+    \Delta &= \|x^{cpu}-x^{gpu}\|_{\infty},
+  \end{aligned}

 \end{equation}
 where $x^{cpu}$ and $x^{gpu}$ are the solutions computed, respectively, on the CPU cluster and on the GPU cluster.
 We can see that the precision of the solutions computed on the GPU cluster are sufficient, they are about $10^{-10}$,
 \end{equation}
 where $x^{cpu}$ and $x^{gpu}$ are the solutions computed, respectively, on the CPU cluster and on the GPU cluster.
 We can see that the precision of the solutions computed on the GPU cluster are sufficient, they are about $10^{-10}$,


@@ -460,18 +470,17 @@ developed in C programming language a generator of large sparse matrices having
 in  most numerical problems. This generator uses the sparse matrices of the Davis collection as the initial
 matrices to build the large band matrices. It is executed in parallel by all the MPI processes of the cluster
 so that each process constructs its own sub-matrix as a rectangular block of the global sparse matrix. Each process
 in  most numerical problems. This generator uses the sparse matrices of the Davis collection as the initial
 matrices to build the large band matrices. It is executed in parallel by all the MPI processes of the cluster
 so that each process constructs its own sub-matrix as a rectangular block of the global sparse matrix. Each process

-$i$ computes the size $n_i$ and the offset $offset_i$ of its sub-matrix in the global sparse matrix according to the
+$i$ computes the size $n_i$ and the offset $\Offset_i$ of its sub-matrix in the global sparse matrix according to the

 size $n$ of the linear system to be solved and the number of the GPU computing nodes $p$ as follows:
 \begin{equation}
   n_i = \frac{n}{p},
 \end{equation} 
 \begin{equation}
 size $n$ of the linear system to be solved and the number of the GPU computing nodes $p$ as follows:
 \begin{equation}
   n_i = \frac{n}{p},
 \end{equation} 
 \begin{equation}

-  offset_i = \left\{
-  \begin{array}{l}
-  0\mbox{~if~}i=0,\\
-  offset_{i-1}+n_{i-1}\mbox{~otherwise.}
-  \end{array}
-  \right.
+  \Offset_i =
+  \begin{cases}
+    0 & \text{if $i=0$,}\\
+    \Offset_{i-1}+n_{i-1} & \text{otherwise.}
+  \end{cases}

 \end{equation} 
 So each process $i$ performs several copies of the same initial matrix chosen from the Davis collection and it
 puts all these copies on the main diagonal of the global matrix in order to construct a band matrix. Moreover,
 \end{equation} 
 So each process $i$ performs several copies of the same initial matrix chosen from the Davis collection and it
 puts all these copies on the main diagonal of the global matrix in order to construct a band matrix. Moreover,


@@ -495,7 +504,7 @@ that for solving large sparse matrices the GPU cluster is more efficient (about
 cluster. The computing power of the GPUs allows to accelerate the computation of the functions operating
 on large vectors of the parallel GMRES algorithm.
 
 cluster. The computing power of the GPUs allows to accelerate the computation of the functions operating
 on large vectors of the parallel GMRES algorithm.
 

-\begin{table}[!h]
+\begin{table}[!ht]

 \begin{center}
 \begin{tabular}{|c|c|r|r|} 
 \hline
 \begin{center}
 \begin{tabular}{|c|c|r|r|} 
 \hline


@@ -519,11 +528,11 @@ Matrix type                 & Name              & \# nonzeros & Bandwidth \\ \hl
 \end{table}
 
 
 \end{table}
 
 

-\begin{table}[!h]
+\begin{table}[!ht]

 \begin{center}
 \begin{tabular}{|c|c|c|c|c|c|c|} 
 \hline
 \begin{center}
 \begin{tabular}{|c|c|c|c|c|c|c|} 
 \hline

-Matrix            & $Time_{cpu}$ & $Time_{gpu}$ & $\tau$ & \# $iter$& $prec$   & $\Delta$   \\ \hline \hline
+Matrix            & $\Time{cpu}$ & $\Time{gpu}$ & $\tau$ & \# $\Iter$& $\Prec$   & $\Delta$   \\ \hline \hline

 2cubes\_sphere    & 3.683 s     & 0.870 s      & 4.23   & 21       & 2.11e-14 & 8.67e-18 \\
 ecology2          & 2.570 s     & 0.424 s      & 6.06   & 21       & 4.88e-13 & 2.08e-14 \\
 finan512          & 2.727 s     & 0.533 s      & 5.11   & 17       & 3.22e-12 & 8.82e-14 \\
 2cubes\_sphere    & 3.683 s     & 0.870 s      & 4.23   & 21       & 2.11e-14 & 8.67e-18 \\
 ecology2          & 2.570 s     & 0.424 s      & 6.06   & 21       & 4.88e-13 & 2.08e-14 \\
 finan512          & 2.727 s     & 0.533 s      & 5.11   & 17       & 3.22e-12 & 8.82e-14 \\


@@ -570,13 +579,20 @@ and 9) and two elements from  node 3 (elements 13 and 14). Therefore to reduce t
 of the unused vector elements, the GPU computing nodes must exchange between them only the vector elements necessary
 to perform their local sparse matrix-vector multiplications.
 
 of the unused vector elements, the GPU computing nodes must exchange between them only the vector elements necessary
 to perform their local sparse matrix-vector multiplications.
 

-\begin{figure}
+\begin{figure}[!ht]

 \centering
   \includegraphics[width=120mm,keepaspectratio]{Figures/compress}
 \caption{Example of data exchanges between node 1 and its neighbors 0, 2 and 3.}
 \label{fig:05}
 \end{figure}
 
 \centering
   \includegraphics[width=120mm,keepaspectratio]{Figures/compress}
 \caption{Example of data exchanges between node 1 and its neighbors 0, 2 and 3.}
 \label{fig:05}
 \end{figure}
 

+\begin{figure}[!ht]
+\centering
+  \includegraphics[width=100mm,keepaspectratio]{Figures/reorder}
+\caption{Reordering of the columns of a local sparse matrix.}
+\label{fig:06}
+\end{figure}
+

 We propose to use a compressed storage format of the sparse global vector. In Figure~\ref{fig:05}, we show an
 example of the data exchanges between node 1 and its neighbors to construct the compressed global vector.
 First, the neighboring nodes 0, 2 and 3 determine the vector elements needed by node 1 and, then, they send
 We propose to use a compressed storage format of the sparse global vector. In Figure~\ref{fig:05}, we show an
 example of the data exchanges between node 1 and its neighbors to construct the compressed global vector.
 First, the neighboring nodes 0, 2 and 3 determine the vector elements needed by node 1 and, then, they send


@@ -588,13 +604,6 @@ For performance purposes, the computation of the shared data to send to the neig
 by the GPU as a kernel. In addition, we use the MPI point-to-point communication routines: a blocking send routine
 \verb+MPI_Send()+ and a nonblocking receive routine \verb+MPI_Irecv()+.
 
 by the GPU as a kernel. In addition, we use the MPI point-to-point communication routines: a blocking send routine
 \verb+MPI_Send()+ and a nonblocking receive routine \verb+MPI_Irecv()+.
 

-\begin{figure}
-\centering
-  \includegraphics[width=100mm,keepaspectratio]{Figures/reorder}
-\caption{Reordering of the columns of a local sparse matrix.}
-\label{fig:06}
-\end{figure}
-

 Table~\ref{tab:05} shows the performances of the parallel GMRES algorithm using the compressed storage format
 of the sparse global vector. The results are obtained from solving large linear systems associated to the sparse
 band matrices presented in Table~\ref{tab:03}. We can see from Table~\ref{tab:05} that the execution times 
 Table~\ref{tab:05} shows the performances of the parallel GMRES algorithm using the compressed storage format
 of the sparse global vector. The results are obtained from solving large linear systems associated to the sparse
 band matrices presented in Table~\ref{tab:03}. We can see from Table~\ref{tab:05} that the execution times 


@@ -604,11 +613,11 @@ and those spent on the cluster of 12 GPUs have increased. Indeed, the reordering
 the use of a compressed storage format for the sparse vectors minimize the communication overheads between the
 GPU computing nodes.
 
 the use of a compressed storage format for the sparse vectors minimize the communication overheads between the
 GPU computing nodes.
 

-\begin{table}[!h]
+\begin{table}[!ht]

 \begin{center}
 \begin{tabular}{|c|c|c|c|c|c|c|} 
 \hline
 \begin{center}
 \begin{tabular}{|c|c|c|c|c|c|c|} 
 \hline

-Matrix            & $Time_{cpu}$ & $Time_{gpu}$ & $\tau$& \# $iter$& $prec$   & $\Delta$   \\ \hline \hline
+Matrix            & $\Time{cpu}$ & $\Time{gpu}$ & $\tau$& \# $\Iter$& $\Prec$   & $\Delta$   \\ \hline \hline

 2cubes\_sphere    & 3.597 s     & 0.514 s      & 6.99  & 21       & 2.11e-14 & 8.67e-18 \\
 ecology2          & 2.549 s     & 0.288 s      & 8.83  & 21       & 4.88e-13 & 2.08e-14 \\
 finan512          & 2.660 s     & 0.377 s      & 7.05  & 17       & 3.22e-12 & 8.82e-14 \\
 2cubes\_sphere    & 3.597 s     & 0.514 s      & 6.99  & 21       & 2.11e-14 & 8.67e-18 \\
 ecology2          & 2.549 s     & 0.288 s      & 8.83  & 21       & 4.88e-13 & 2.08e-14 \\
 finan512          & 2.660 s     & 0.377 s      & 7.05  & 17       & 3.22e-12 & 8.82e-14 \\


@@ -629,7 +638,6 @@ a cluster of 12 GPUs.}
 \end{center}
 \end{table}
 
 \end{center}
 \end{table}
 

-

 \subsection{Hypergraph partitioning}
 \label{sec:06.02}
 In this section, we use another structure of the sparse matrices. We are interested in sparse matrices
 \subsection{Hypergraph partitioning}
 \label{sec:06.02}
 In this section, we use another structure of the sparse matrices. We are interested in sparse matrices


@@ -641,15 +649,14 @@ off-diagonals on the left and right of the main diagonal. Table~\ref{tab:06} sho
 of sparse matrices of size 25 million of rows and generated from those of the Davis collection. We can
 see in the fourth column that the bandwidths of these matrices are as large as their sizes.
 
 of sparse matrices of size 25 million of rows and generated from those of the Davis collection. We can
 see in the fourth column that the bandwidths of these matrices are as large as their sizes.
 

-\begin{figure}
+\begin{figure}[!ht]

 \centering
   \includegraphics[width=100mm,keepaspectratio]{Figures/generation_1}
 \caption{Example of the generation of a large sparse matrix having five bands by four computing nodes.}
 \label{fig:07}
 \end{figure}
 
 \centering
   \includegraphics[width=100mm,keepaspectratio]{Figures/generation_1}
 \caption{Example of the generation of a large sparse matrix having five bands by four computing nodes.}
 \label{fig:07}
 \end{figure}
 

-
-\begin{table}[!h]
+\begin{table}[!ht]

 \begin{center}
 \begin{tabular}{|c|c|r|r|} 
 \hline
 \begin{center}
 \begin{tabular}{|c|c|r|r|} 
 \hline


@@ -681,11 +688,11 @@ In fact, the naive partitioning row-by-row or column-by-column of this type of s
 a GPU node to many neighboring nodes and produces a significant number of data dependencies between
 the different GPU nodes.   
 
 a GPU node to many neighboring nodes and produces a significant number of data dependencies between
 the different GPU nodes.   
 

-\begin{table}[!h]
+\begin{table}[!ht]

 \begin{center}
 \begin{tabular}{|c|c|c|c|c|c|c|} 
 \hline
 \begin{center}
 \begin{tabular}{|c|c|c|c|c|c|c|} 
 \hline

-Matrix            & $Time_{cpu}$ & $Time_{gpu}$ & $\tau$ & \# $iter$& $prec$ & $\Delta$   \\ \hline \hline
+Matrix            & $\Time{cpu}$ & $\Time{gpu}$ & $\tau$ & \# $\Iter$& $\Prec$ & $\Delta$   \\ \hline \hline

 2cubes\_sphere    & 15.963 s    & 7.250 s      & 2.20  & 58     & 6.23e-16 & 3.25e-19 \\
 ecology2          & 3.549 s     & 2.176 s      & 1.63  & 21     & 4.78e-15 & 1.06e-15 \\
 finan512          & 3.862 s     & 1.934 s      & 1.99  & 17     & 3.21e-14 & 8.43e-17 \\
 2cubes\_sphere    & 15.963 s    & 7.250 s      & 2.20  & 58     & 6.23e-16 & 3.25e-19 \\
 ecology2          & 3.549 s     & 2.176 s      & 1.63  & 21     & 4.78e-15 & 1.06e-15 \\
 finan512          & 3.862 s     & 1.934 s      & 1.99  & 17     & 3.21e-14 & 8.43e-17 \\


@@ -756,9 +763,9 @@ where $\mathcal{E}_C$ is the set of the cut hyperedges issued from the partition
 and $\lambda_j$ are, respectively, the cost and the connectivity of the cut hyperedge $e_j$. In addition,
 the hypergraph partitioning is constrained to maintain the load balance between the $K$ parts:
 \begin{equation}
 and $\lambda_j$ are, respectively, the cost and the connectivity of the cut hyperedge $e_j$. In addition,
 the hypergraph partitioning is constrained to maintain the load balance between the $K$ parts:
 \begin{equation}

-W_k\leq (1+\epsilon)W_{avg}\mbox{,~}(1\leq k\leq K)\mbox{~and~}(0<\epsilon<1),
+W_k\leq (1+\epsilon)\Wavg\mbox{,~}(1\leq k\leq K)\mbox{~and~}(0<\epsilon<1),

 \end{equation}
 \end{equation}

-where $W_k$ is the sum of the vertex weights in the subset $\mathcal{V}_k$, $W_{avg}$ is the average part's
+where $W_k$ is the sum of the vertex weights in the subset $\mathcal{V}_k$, $\Wavg$ is the average part's

 weight and $\epsilon$ is the maximum allowed imbalanced ratio.
 
 The hypergraph partitioning is a NP-complete problem but software tools using heuristics are developed, for
 weight and $\epsilon$ is the maximum allowed imbalanced ratio.
 
 The hypergraph partitioning is a NP-complete problem but software tools using heuristics are developed, for


@@ -774,18 +781,18 @@ compressed storage format of the sparse vectors and the parallel hypergraph part
 Zoltan tool~\cite{ref20,ref21}. The parameters of the hypergraph partitioning are initialized as follows:
 \begin{itemize}
 \item The weight $w_i$ of each vertex $v_i$ is set to the number of the nonzero values on the matrix row $i$,
 Zoltan tool~\cite{ref20,ref21}. The parameters of the hypergraph partitioning are initialized as follows:
 \begin{itemize}
 \item The weight $w_i$ of each vertex $v_i$ is set to the number of the nonzero values on the matrix row $i$,

-\item For simplicity sake, the cost $c_j$ of each hyperedge $e_j$ is set to 1,
+\item For simplicity's sake, the cost $c_j$ of each hyperedge $e_j$ is set to 1,

 \item The maximum imbalanced ratio $\epsilon$ is limited to 10\%.  
 \end{itemize} 
 We can notice from Table~\ref{tab:08} that the execution times on the cluster of 12 GPUs are significantly
 improved compared to those presented in Table~\ref{tab:07}. The hypergraph partitioning applied on the large
 sparse matrices having large bandwidths have improved the execution times on the GPU cluster by about 65\%.
 
 \item The maximum imbalanced ratio $\epsilon$ is limited to 10\%.  
 \end{itemize} 
 We can notice from Table~\ref{tab:08} that the execution times on the cluster of 12 GPUs are significantly
 improved compared to those presented in Table~\ref{tab:07}. The hypergraph partitioning applied on the large
 sparse matrices having large bandwidths have improved the execution times on the GPU cluster by about 65\%.
 

-\begin{table}[!h]
+\begin{table}[!ht]

 \begin{center}
 \begin{tabular}{|c|c|c|c|c|c|c|} 
 \hline
 \begin{center}
 \begin{tabular}{|c|c|c|c|c|c|c|} 
 \hline

-Matrix            & $Time_{cpu}$ & $Time_{gpu}$ & $\tau$  & \# iter & $prec$     & $\Delta$   \\ \hline \hline
+Matrix            & $\Time{cpu}$ & $\Time{gpu}$ & $\tau$  & \# $\Iter$ & $\Prec$     & $\Delta$   \\ \hline \hline

 2cubes\_sphere    & 16.430 s    & 2.840 s      & 5.78 & 58     & 6.23e-16 & 3.25e-19 \\
 ecology2          & 3.152 s     & 0.367 s      & 8.59 & 21     & 4.78e-15 & 1.06e-15 \\
 finan512          & 3.672 s     & 0.723 s      & 5.08 & 17     & 3.21e-14 & 8.43e-17 \\
 2cubes\_sphere    & 16.430 s    & 2.840 s      & 5.78 & 58     & 6.23e-16 & 3.25e-19 \\
 ecology2          & 3.152 s     & 0.367 s      & 8.59 & 21     & 4.78e-15 & 1.06e-15 \\
 finan512          & 3.672 s     & 0.723 s      & 5.08 & 17     & 3.21e-14 & 8.43e-17 \\


@@ -806,24 +813,14 @@ having five bands on a cluster of 24 CPU cores vs. a cluster of 12 GPUs.}
 \end{center}
 \end{table}
 
 \end{center}
 \end{table}
 

-Table~\ref{tab:09} shows in the second, third and fourth columns the total communication volume on a cluster of 12 GPUs by using row-by-row partitioning or hypergraph partitioning and compressed format. The total communication volume defines the total number of the vector elements exchanged between the 12 GPUs. From these columns we can see that the two heuristics, compressed format for the vectors and the hypergraph partitioning, minimize the number the vector elements to be exchanged over the GPU cluster. The compressed format applied  
-
+Table~\ref{tab:09} shows in the second, third and fourth columns the total communication volume on a cluster of 12 GPUs by using row-by-row partitioning or hypergraph partitioning and compressed format. The total communication volume defines the total number of the vector elements exchanged between the 12 GPUs. From these columns we can see that the two heuristics, compressed format for the vectors and the hypergraph partitioning, minimize the number of vector elements to be exchanged over the GPU cluster. The compressed format allows the GPUs to exchange the needed vector elements without any communication overheads. The hypergraph partitioning allows to split the large sparse matrices so as to minimize  data dependencies between the GPU computing nodes. However, we can notice in the fifth column that the hypergraph partitioning takes longer than the computation times. As we have mentioned before, the hypergraph partitioning method is less efficient in terms of memory consumption and partitioning time than its graph counterpart. So for the applications which often use the same sparse matrices, we can perform the hypergraph partitioning only once and, then, we save the traces in files to be reused several times. Therefore, this allows us to avoid the partitioning of the sparse matrices at each resolution of the linear systems.

 
 

-
-This table shows that the hypergraph partitioning allows to split the large sparse matrices so as to minimize  data
-dependencies between the GPU computing nodes. However, we can notice in the fourth column that the hypergraph
-partitioning takes longer than the computation times. As we have mentioned before, the hypergraph partitioning
-method is less efficient in terms of memory consumption and partitioning time than its graph counterpart.
-So for the applications which often use the same sparse matrices, we can perform the hypergraph partitioning
-only once and, then, we save the traces in files to be reused several times. Therefore, this allows us to
-avoid the partitioning of the sparse matrices at each resolution of the linear systems. 
-
-\begin{table}[!h]
+\begin{table}

 \begin{center}
 \begin{tabular}{|c|c|c|c|c|} 
 \hline
 \begin{center}
 \begin{tabular}{|c|c|c|c|c|} 
 \hline

-\multirow{3}{*}{Matrix} & Total comm. vol. & Total comm. vol. & Total comm. vol.              & Time of hypergraph \\
-                        & using row-by row & using compressed & using hypergraph partitioning & partitioning \\
+\multirow{3}{*}{Matrix} & Total comm. vol. & Total comm. vol. & Total comm. vol. using        & Time of hypergraph \\
+                        & using row-by row & using compressed & hypergraph partitioning       & partitioning \\

                         & partitioning     & format           & and compressed format         & in minutes \\ \hline \hline
 2cubes\_sphere          & 182 061 791      & 25 360 543       & 240 679                       & 68.98         \\
 ecology2                & 181 267 000      & 26 044 002       & 73 021                        & 4.92          \\
                         & partitioning     & format           & and compressed format         & in minutes \\ \hline \hline
 2cubes\_sphere          & 182 061 791      & 25 360 543       & 240 679                       & 68.98         \\
 ecology2                & 181 267 000      & 26 044 002       & 73 021                        & 4.92          \\


@@ -843,16 +840,112 @@ torso3                  & 183 863 292      & 25 682 514       & 613 250
 \end{center}
 \end{table}
 
 \end{center}
 \end{table}
 

-%%% MODIF %%%
-\textcolor{red}{\bf In order to show the influence of the communications on a GPU cluster
-In tables, we compute the ratios of the computation time over the communication time to show the influence of the communications on a GPU cluster compared to a CPU cluster}

 
 

-\textcolor{red}{\bf Finally, the parallel solving of a linear system can be easy to optimize when the associated matrix is regular. This is unfortunately not the case of many real-world applications. When the matrix has an irregular structure, the amount of communication between processors is not the same. Another important parameter is the size of the matrix bandwidth which has a huge influence on the amount of communications. In this work, we have generated different kinds of matrices in order to analyze different difficulties. With as a large bandwidth as possible involving communications between all processors, which is the most difficult situation, we proposed to use two heuristics. Unfortunately, there is no fast method that optimizes the communication in any situation. For systems of non linear equations, there are different algorithms but most of them consist in linearizing the system of equations. In this case, a linear system needs to be solved. The big interest is that the matrix is the same at each step of the non linear system solving, so the partitioning method which is a time consuming step is performed once only.
-}

 
 

-\textcolor{red}{\bf 
-Another very important issue is that the communications have a greater influence on a cluster of GPUs than on a cluster of CPUs. There are two reasons for this. The first one comes from the fact that with a cluster of GPUs, the CPU/GPU data transfers slow down communications between two GPUs that are not on the same machines. The second one is due to the fact that with GPUs the ratio of the computation time over the communication time decreases since the computation time is reduced. So the impact of the communications between GPUs might be a very important issue that can limit the scalability of a parallel algorithm.}
-%%% END %%%
+
+
+
+
+
+
+
+
+
+
+
+
+Hereafter, we show the influence of the communications on a GPU cluster compared to a CPU cluster. In Tables~\ref{tab:10},~\ref{tab:11} and~\ref{tab:12}, we compute the ratios between the computation time over the communication time of three versions of the parallel GMRES algorithm to solve sparse linear systems associated to matrices of Table~\ref{tab:06}. These tables show that the hypergraph partitioning and the compressed format of the vectors increase the ratios either on the GPU cluster or on the CPU cluster. That means that the two optimization techniques allow the minimization of the total communication volume between the computing nodes. However, we can notice that the ratios obtained on the GPU cluster are lower than those obtained on the CPU cluster. Indeed, GPUs compute faster than CPUs but with GPUs there are more communications due to CPU/GPU communications, so communications are more time-consuming while the computation time remains unchanged. Furthermore, we can notice that the GPU computation times on Tables~\ref{tab:11} and~\ref{tab:12} are about 10\% lower than those on Table~\ref{tab:10}. Indeed, the compression of the vectors and the reordering of matrix columns allow to perform coalesced accesses to the GPU memory and thus accelerate the sparse matrix-vector multiplication.  
+
+\begin{table}[p]
+\begin{center}
+\begin{tabular}{|c||c|c|c||c|c|c|} 
+\hline
+\multirow{2}{*}{Matrix} & \multicolumn{3}{c||}{GPU version} & \multicolumn{3}{c|}{CPU version}  \\ \cline{2-7}
+                  & $\Time{comput}$ & $\Time{comm}$ & $\Ratio$ & $\Time{comput}$ & $\Time{comm}$ & $\Ratio$ \\ \hline \hline
+2cubes\_sphere    & 37.067 s       & 1434.512 s   & {\bf 0.026}   & 312.061 s      & 3453.931 s   & {\bf 0.090}\\
+ecology2          & 4.116 s        & 501.327 s    & {\bf 0.008}   & 60.776 s       & 1216.607 s   & {\bf 0.050}\\
+finan512          & 7.170 s        & 386.742 s    & {\bf 0.019}   & 72.464 s       & 932.538 s    & {\bf 0.078}\\
+G3\_circuit       & 4.797 s        & 537.343 s    & {\bf 0.009}   & 66.011 s       & 1407.378 s   & {\bf 0.047}\\
+shallow\_water2   & 3.620 s        & 411.208 s    & {\bf 0.009}   & 51.294 s       & 973.446 s    & {\bf 0.053}\\ 
+thermal2          & 6.902 s        & 511.618 s    & {\bf 0.013}   & 77.255 s       & 1281.979 s   & {\bf 0.060}\\ \hline \hline
+cage13            & 12.837 s       & 625.175 s    & {\bf 0.021}   & 139.178 s      & 1518.349 s   & {\bf 0.092}\\
+crashbasis        & 48.532 s       & 3195.183 s   & {\bf 0.015}   & 623.686 s      & 7741.777 s   & {\bf 0.081}\\
+FEM\_3D\_thermal2 & 37.211 s       & 1584.650 s   & {\bf 0.023}   & 370.297 s      & 3810.255 s   & {\bf 0.097}\\
+language          & 22.912 s       & 2242.897 s   & {\bf 0.010}   & 286.682 s      & 5348.733 s   & {\bf 0.054}\\
+poli\_large       & 13.618 s       & 1722.304 s   & {\bf 0.008}   & 190.302 s      & 4059.642 s   & {\bf 0.047}\\
+torso3            & 74.194 s       & 4454.936 s   & {\bf 0.017}   & 897.440 s      & 10800.787 s  & {\bf 0.083}\\ \hline       
+\end{tabular}
+\caption{Ratios of the computation time over the communication time obtained from the parallel GMRES algorithm using row-by-row partitioning on 12 GPUs and 24 CPUs.}
+\label{tab:10}
+\end{center}
+\end{table}
+
+
+\begin{table}[p]
+\begin{center}
+\begin{tabular}{|c||c|c|c||c|c|c|} 
+\hline
+\multirow{2}{*}{Matrix} & \multicolumn{3}{c||}{GPU version} & \multicolumn{3}{c|}{CPU version}  \\ \cline{2-7}
+                  & $\Time{comput}$ & $\Time{comm}$ & $\Ratio$ & $\Time{comput}$ & $\Time{comm}$ & $\Ratio$ \\ \hline \hline
+2cubes\_sphere    & 27.386 s       & 154.861 s   & {\bf 0.177}   & 342.255 s      & 42.100 s   & {\bf 8.130}\\
+ecology2          & 3.822 s        & 53.131 s    & {\bf 0.072}   & 69.956 s       & 15.019 s   & {\bf 4.658}\\
+finan512          & 6.366 s        & 41.155 s    & {\bf 0.155}   & 79.592 s       & 8.604 s    & {\bf 9.251}\\
+G3\_circuit       & 4.543 s        & 63.132 s    & {\bf 0.072}   & 76.540 s       & 27.371 s   & {\bf 2.796}\\
+shallow\_water2   & 3.282 s        & 43.080 s    & {\bf 0.076}   & 58.348 s       & 8.088 s    & {\bf 7.214}\\ 
+thermal2          & 5.986 s        & 57.100 s    & {\bf 0.105}   & 87.682 s       & 28.544 s   & {\bf 3.072}\\ \hline \hline
+cage13            & 10.227 s       & 70.388 s    & {\bf 0.145}   & 152.718 s      & 30.785 s   & {\bf 4.961}\\
+crashbasis        & 41.527 s       & 369.071 s   & {\bf 0.113}   & 701.040 s      & 158.916 s  & {\bf 4.411}\\
+FEM\_3D\_thermal2 & 28.691 s       & 167.140 s   & {\bf 0.172}   & 403.510 s      & 50.935 s   & {\bf 7.922}\\
+language          & 22.408 s       & 242.589 s   & {\bf 0.092}   & 333.119 s      & 64.409 s   & {\bf 5.172}\\
+poli\_large       & 13.710 s       & 179.208 s   & {\bf 0.077}   & 215.934 s      & 30.903 s   & {\bf 6.987}\\
+torso3            & 58.455 s       & 480.315 s   & {\bf 0.122}   & 993.609 s      & 152.173 s  & {\bf 6.529}\\ \hline       
+\end{tabular}
+\caption{Ratios of the computation time over the communication time obtained from the parallel GMRES algorithm using row-by-row partitioning and compressed format for vectors on 12 GPUs and 24 CPUs.}
+\label{tab:11}
+\end{center}
+\end{table}
+
+
+\begin{table}[p]
+\begin{center}
+\begin{tabular}{|c||c|c|c||c|c|c|} 
+\hline
+\multirow{2}{*}{Matrix} & \multicolumn{3}{c||}{GPU version} & \multicolumn{3}{c|}{CPU version}  \\ \cline{2-7}
+                  & $\Time{comput}$ & $\Time{comm}$ & $\Ratio$ & $\Time{comput}$ & $\Time{comm}$ & $\Ratio$ \\ \hline \hline
+2cubes\_sphere    & 28.440 s       & 7.768 s      & {\bf 3.661}   & 327.109 s      & 63.788 s   & {\bf 5.128}\\
+ecology2          & 3.652 s        & 0.757 s      & {\bf 4.823}   & 63.632 s       & 13.520 s   & {\bf 4.707}\\
+finan512          & 7.579 s        & 4.569 s      & {\bf 1.659}   & 74.120 s       & 22.505 s   & {\bf 3.294}\\
+G3\_circuit       & 4.876 s        & 8.745 s      & {\bf 0.558}   & 72.280 s       & 28.395 s   & {\bf 2.546}\\
+shallow\_water2   & 3.146 s        & 0.606 s      & {\bf 5.191}   & 52.903 s       & 11.177 s   & {\bf 4.733}\\ 
+thermal2          & 6.473 s        & 4.325 s      & {\bf 1.497}   & 81.171 s       & 20.907 s   & {\bf 3.882}\\ \hline \hline
+cage13            & 11.676 s       & 7.723 s      & {\bf 1.512}   & 145.755 s      & 46.547 s   & {\bf 3.131}\\
+crashbasis        & 42.799 s       & 29.399 s     & {\bf 1.456}   & 650.386 s      & 203.918 s  & {\bf 3.189}\\
+FEM\_3D\_thermal2 & 29.875 s       & 8.915 s      & {\bf 3.351}   & 382.887 s      & 93.252 s   & {\bf 4.106}\\
+language          & 20.991 s       & 11.197 s     & {\bf 1.875}   & 310.679 s      & 82.480 s   & {\bf 3.767}\\
+poli\_large       & 13.817 s       & 102.760 s    & {\bf 0.134}   & 197.508 s      & 151.672 s  & {\bf 1.302}\\
+torso3            & 57.469 s       & 16.828 s     & {\bf 3.415}   & 926.588 s      & 242.721 s  & {\bf 3.817}\\ \hline       
+\end{tabular}
+\caption{Ratios of the computation time over the communication time obtained from the parallel GMRES algorithm using hypergraph partitioning and compressed format for vectors on 12 GPUs and 24 CPUs.}
+\label{tab:12}
+\end{center}
+\end{table}
+
+\begin{figure}
+\centering
+  \includegraphics[width=120mm,keepaspectratio]{Figures/weak}
+\caption{Weak scaling of the parallel GMRES algorithm on a GPU cluster.}
+\label{fig:09}
+\end{figure}
+
+ Figure~\ref{fig:09} presents the weak scaling of four versions of the parallel GMRES algorithm on a GPU cluster. We fixed the size of a sub-matrix to 5 million of rows per GPU computing node. We used matrices having five bands generated from the symmetric matrix thermal2. This figure shows that the parallel GMRES algorithm, in its naive version or using either the compression format for vectors or the hypergraph partitioning, is not scalable on a GPU cluster due to the large amount of communications between GPUs. In contrast, we can see that the algorithm using both optimization techniques is fairly scalable. That means that in this version the cost of communications is relatively constant regardless the number of computing nodes in the cluster.
+
+\bigskip
+
+ Finally, as far as we are concerned, the parallel solving of a linear system can be easy to optimize when the associated matrix is regular. This is unfortunately not the case for many real-world applications. When the matrix has an irregular structure, the amount of communications between processors is not the same. Another important parameter is the size of the matrix bandwidth which has a huge influence on the amount of communications. In this work, we have generated different kinds of matrices in order to analyze several difficulties. With a bandwidth as large as possible, involving communications between all processors, which is the most difficult situation, we proposed to use two heuristics. Unfortunately, there is no fast method that optimizes the communications in any situation. For systems of non linear equations, there are different algorithms but most of them consist in linearizing the system of equations. In this case, a linear system needs to be solved. The big interest is that the matrix is the same at each step of the non linear system solving, so the partitioning method which is a time consuming step is performed only once.
+
+
+ 
+Another very important issue, which might be ignored by too many people, is that the communications have a greater influence on a cluster of GPUs than on a cluster of CPUs. There are two reasons for that. The first one comes from the fact that with a cluster of GPUs, the CPU/GPU data transfers slow down communications between two GPUs that are not on the same machines. The second one is due to the fact that with GPUs the ratio of the computation time over the communication time decreases since the computation time is reduced. So the impact of the communications between GPUs might be a very important issue that can limit the scalability of parallel algorithms.

 
 %%--------------------%%
 %%      SECTION 7     %%
 
 %%--------------------%%
 %%      SECTION 7     %%


@@ -886,7 +979,7 @@ and to improve the performances of the parallel GMRES algorithm as the multispli
 The recent GPU hardware and software architectures provide the GPU-Direct system which allows
 two GPUs, placed in the same machine or in two remote machines, to exchange data without using
 CPUs. This improves the data transfers between GPUs. Finally, it would be interesting to implement
 The recent GPU hardware and software architectures provide the GPU-Direct system which allows
 two GPUs, placed in the same machine or in two remote machines, to exchange data without using
 CPUs. This improves the data transfers between GPUs. Finally, it would be interesting to implement

-other iterative methods on GPU clusters for solving large sparse linear or nonlinear systems.
+other iterative methods on GPU clusters for solving large sparse linear or non linear systems.

 
 \paragraph{Acknowledgments}
 This paper is based upon work supported by the R\'egion de Franche-Comt\'e.
 
 \paragraph{Acknowledgments}
 This paper is based upon work supported by the R\'egion de Franche-Comt\'e.