-followings steps. First roots are shared with OpenMP. Each thread
-starts by copying all the previous roots inside its GPU. Then each GPU
-will compute an iteration of the EA method on its own roots. For that
-all the other roots are used. At the end of an iteration, the updated
-roots are copied from the GPU to the CPU. The convergence is checked
-on the new roots. Finally each CPU will update its own roots in the
-shared memory arrays containing all the roots.
-
-%In principle a grid is set by two parameter DimGrid, the number of block per grid, DimBloc: the number of threads per block. The following schema shows the architecture of (CUDA,OpenMP).
-
-%\begin{figure}[htbp]
-%\centering
- % \includegraphics[angle=-90,width=0.5\textwidth]{OpenMP-CUDA}
-%\caption{The OpenMP-CUDA architecture}
-%\label{fig:03}
-%\end{figure}
-%Each thread OpenMP compute the kernels on GPUs,than after each iteration they copy out the data from GPU memory to CPU shared memory. The kernels are re-runs is up to the roots converge sufficiently. Here are below the corresponding algorithm:
-
-%% \RC{Surement à virer ou réécrire pour etre compris sans algo}
-%% $num\_gpus$ OpenMP threads are created using
-%% \verb=omp_set_num_threads();=function (step $3$, Algorithm
-%% \ref{alg2-cuda-openmp}), the shared memory is created using
-%% \verb=#pragma omp parallel shared()= OpenMP function (line $5$,
-%% Algorithm\ref{alg2-cuda-openmp}), then each OpenMP thread allocates
-%% memory and copies initial data from CPU memory to GPU global memory,
-%% executes the kernels on GPU, but computes only his portion of roots
-%% indicated with variable \textit{index} initialized in (line 5,
-%% Algorithm \ref{alg2-cuda-openmp}), used as input data in the
-%% $kernel\_update$ (line 10, Algorithm \ref{alg2-cuda-openmp}). After
-%% each iteration, all OpenMP threads synchronize using
-%% \verb=#pragma omp barrier;= to gather all the correct values of
-%% $\Delta z$, thus allowing the computation the maximum stop condition
-%% on vector $\Delta z$ (line 12, Algorithm
-%% \ref{alg2-cuda-openmp}). Finally, threads copy the results from GPU
-%% memories to CPU memory. The OpenMP threads execute kernels until the
-%% roots sufficiently converge.
-
-
-%% \begin{enumerate}
-%% \begin{algorithm}[htpb]
-%% \label{alg2-cuda-openmp}
-%% %\LinesNumbered
-%% \caption{CUDA-OpenMP Algorithm to find roots with the Ehrlich-Aberth method}
-
-%% \KwIn{$Z^{0}$ (Initial root's vector), $\varepsilon$ (Error tolerance
-%% threshold), P (Polynomial to solve), Pu (Derivative of P), $n$ (Polynomial degree), $\Delta z$ ( Vector of errors for stop condition), $num_gpus$ (number of OpenMP threads/ Number of GPUs), $Size$ (number of roots)}
-
-%% \KwOut {$Z$ ( Root's vector), $ZPrec$ (Previous root's vector)}
-
-%% \BlankLine
-
-%% \item Initialization of P\;
-%% \item Initialization of Pu\;
-%% \item Initialization of the solution vector $Z^{0}$\;
-%% \verb=omp_set_num_threads(num_gpus);=
-%% \verb=#pragma omp parallel shared(Z,$\Delta$ z,P);=
-%% \verb=cudaGetDevice(gpu_id);=
-%% \item Allocate and copy initial data from CPU memory to the GPU global memories\;
-%% \item index= $Size/num\_gpus$\;
-%% \item k=0\;
-%% \While {$error > \epsilon$}{
-%% \item Let $\Delta z=0$\;
-%% \item $ kernel\_save(ZPrec,Z)$\;
-%% \item k=k+1\;
-%% \item $ kernel\_update(Z,P,Pu,index)$\;
-%% \item $kernel\_testConverge(\Delta z[gpu\_id],Z,ZPrec)$\;
-%% %\verb=#pragma omp barrier;=
-%% \item error= Max($\Delta z$)\;
-%% }
-
-%% \item Copy results from GPU memories to CPU memory\;
-%% \end{algorithm}
-%% \end{enumerate}
-%% ~\\
-%% \RC{C'est encore pire ici, on ne voit pas les comm CPU <-> GPU }
-
-
-\subsection{Multi-GPU : an MPI-CUDA approach}
-%\begin{figure}[htbp]
-%\centering
- % \includegraphics[angle=-90,width=0.2\textwidth]{MPI-CUDA}
-%\caption{The MPI-CUDA architecture }
-%\label{fig:03}
-%\end{figure}
-Our parallel implementation of EA to find root of polynomials using a CUDA-MPI approach is a data parallel approach. It splits input data of the polynomial to solve among MPI processes. In Algorithm \ref{alg2-cuda-mpi}, input data are the polynomial to solve $P$, the solution vector $Z$, the previous solution vector $ZPrev$, and the value of errors of stop condition $\Delta z$. Let $p$ denote the number of MPI processes on and $n$ the degree of the polynomial to be solved. The algorithm performs a simple data partitioning by creating $p$ portions, of at most $\lceil n/p \rceil$ roots to find per MPI process, for each $Z$ and $ZPrec$. Consequently, each MPI process of rank $k$ will have its own solution vector $Z_{k}$ and $ZPrec$, the error related to the stop condition $\Delta z_{k}$, enabling each MPI process to compute $\lceil n/p \rceil$ roots.
-
-Since a GPU works only on data already allocated in its memory, all local input data, $Z_{k}$, $ZPrec$ and $\Delta z_{k}$, must be transferred from CPU memories to the corresponding GPU memories. Afterwards, the same EA algorithm (Algorithm \ref{alg1-cuda}) is run by all processes but on different polynomial subset of roots $ p(x)_{k}=\sum_{i=1}^{n} a_{i}x^{i}, k=1,...,p$. Each MPI process executes the loop \verb=(While(...)...do)= containing the CUDA kernels but each MPI process computes only its own portion of the roots according to the rule ``''owner computes``''. The local range of roots is indicated with the \textit{index} variable initialized at (line 5, Algorithm \ref{alg2-cuda-mpi}), and passed as an input variable to $kernel\_update$ (line 10, Algorithm \ref{alg2-cuda-mpi}). After each iteration, MPI processes synchronize (\verb=MPI_Allreduce= function) by a reduction on $\Delta z_{k}$ in order to compute the maximum error related to the stop condition. Finally, processes copy the values of new computed roots from GPU memories to CPU memories, then communicate their results to other processes with \verb=MPI_Alltoall= broadcast. If the stop condition is not verified ($error > \epsilon$) then processes stay withing the loop \verb= while(...)...do= until all the roots sufficiently converge.
-
-%% \begin{enumerate}
-%% \begin{algorithm}[htpb]
-%% \label{alg2-cuda-mpi}
-%% %\LinesNumbered
-%% \caption{CUDA-MPI Algorithm to find roots with the Ehrlich-Aberth method}
-
-%% \KwIn{$Z^{0}$ (Initial root's vector), $\varepsilon$ (Error tolerance
-%% threshold), P (Polynomial to solve), Pu (Derivative of P), $n$ (Polynomial degrees), $\Delta z$ ( error of stop condition), $num_gpus$ (number of MPI processes/ number of GPUs), Size (number of roots)}
-
-%% \KwOut {$Z$ (Solution root's vector), $ZPrec$ (Previous solution root's vector)}
-
-%% \BlankLine
-%% \item Initialization of P\;
-%% \item Initialization of Pu\;
-%% \item Initialization of the solution vector $Z^{0}$\;
-%% \item Allocate and copy initial data from CPU memories to GPU global memories\;
-%% \item $index= Size/num_gpus$\;
-%% \item k=0\;
-%% \While {$error > \epsilon$}{
-%% \item Let $\Delta z=0$\;
-%% \item $kernel\_save(ZPrec,Z)$\;
-%% \item k=k+1\;
-%% \item $kernel\_update(Z,P,Pu,index)$\;
-%% \item $kernel\_testConverge(\Delta z,Z,ZPrec)$\;
-%% \item ComputeMaxError($\Delta z$,error)\;
-%% \item Copy results from GPU memories to CPU memories\;
-%% \item Send $Z[id]$ to all processes\;
-%% \item Receive $Z[j]$ from every other process j\;
-%% }
-%% \end{algorithm}
-%% \end{enumerate}
-%% ~\\
-
-%% \RC{ENCORE ENCORE PIRE}
+followings steps. First roots are shared with OpenMP and the
+computation of the local size for each GPU is performed (lines 5-7 in
+Algo\ref{alg2-cuda-openmp}). Each thread starts by copying all the
+previous roots inside its GPU (line 9). Then each GPU will copy the
+previous roots (line 10) and it will compute an iteration of the EA
+method on its own roots (line 11). For that all the other roots are
+used. The convergence is checked on the new roots (line 12). At the end
+of an iteration, the updated roots are copied from the GPU to the
+CPU (line 14) by direcly updating its own roots in the shared memory
+arrays containing all the roots.
+
+
+
+\begin{algorithm}[htpb]
+\LinesNumbered
+\SetAlgoNoLine
+\caption{Finding roots of polynomials with the Ehrlich-Aberth method on multiple GPUs using OpenMP}
+\KwIn{$n$ (polynomial's degree), $\epsilon$ (tolerance threshold), $ngpu$ (number of GPUs)}
+\KwOut{$Z$ (solution vector of roots)}
+Initialize the polynomial $P$ and its derivative $P'$\;
+Set the initial values of vector $Z$\;
+Start of a parallel part with OpenMP ($Z$, $\Delta Z$, $\Delta
+Z_{max}$, $P$, $P'$ are shared variables)\;
+$id_{gpu}$ = cudaGetDevice()\;
+$n_{loc}$ = $n/ngpu$ (local size)\;
+%$idx$ = $id_{gpu}\times n_{loc}$ (local offset)\;
+Copy $P$, $P'$ from CPU to GPU\;
+\While{\emph{not convergence}}{
+ Copy $Z$ from CPU to GPU\;
+ $Z^{prev}$ = KernelSave($Z,n$)\;
+ $Z_{loc}$ = KernelUpdate($P,P',Z^{prev},n_{loc}$)\;
+ $\Delta Z_{loc}$ = KernelComputeError($Z_{loc},Z^{prev}_{loc},n_{loc}$)\;
+ $\Delta Z_{max}[id_{gpu}]$ = CudaMaxFunction($\Delta Z_{loc},n_{loc}$)\;
+ Copy $Z_{loc}$ from GPU to $Z$ in CPU\;
+ $max$ = MaxFunction($\Delta Z_{max},ngpu$)\;
+ TestConvergence($max,\epsilon$)\;
+}
+\label{alg2-cuda-openmp}
+\LZK{J'ai modifié l'algo. Le $P$ est mis shared. Qu'en est-il pour
+ $P'$?}\RC{Je l'ai rajouté. Bon sinon le n\_loc ne remplace pas
+ vraiment un offset et une taille mais bon... et là il y a 4 lignes
+ pour la convergence, c'est bcp ... Zloc, Zmax, max et
+ testconvergence. On pourrait faire mieux}
+\end{algorithm}
+
+
+
+
+
+\subsection{a MPI-CUDA approach}
+
+Our parallel implementation of EA to find root of polynomials using a
+CUDA-MPI approach follows a similar computing approach to the one used
+in CUDA-OpenMP. Each process is responsible to compute its own part of
+roots using all the roots computed by other processors at the previous
+iteration. The difference between both approaches lies in the way
+processes communicate and exchange data. With MPI processors need to
+send and receive data explicitely. So in
+Algorithm~\ref{alg2-cuda-mpi}, after the initialization all the
+processors have the same $Z$ vector. Then they need to compute the
+parameters used by the $MPI\_AlltoAll$ routines (line 4). In practise,
+each processor needs to compute its offset and its local size. Then
+processors need to allocate memory on their GPU (line 5). At the
+beginning of each iteration, a processor starts by transfering the
+whole vector Z from the CPU to the GPU (line 7). Then only the local
+part of $Z^{prev}$ is saved (line 8). After that, a processor is able
+to compute its own roots (line 9). Next, the local error can be
+computed (ligne 10) and the global error (line 11). Then the local
+roots are transfered from the GPU memory to the CPU memory (line 12)
+before being exchanged between all processors (linge 13) in order to
+give to all processors the last version of the roots. If the
+convergence is not statisfied, an new iteration is executed.
+
+
+
+\begin{algorithm}[htpb]
+\label{alg2-cuda-mpi}
+\LinesNumbered
+\caption{CUDA-MPI Algorithm to find roots with the Ehrlich-Aberth method}
+
+\KwIn{$Z^{0}$ (Initial root's vector), $\varepsilon$ (Error tolerance
+ threshold), P (Polynomial to solve), P' (Derivative of P), $n$ (Polynomial degrees), $\Delta z$ ( error of stop condition)}
+
+\KwOut {$Z$ (Solution root's vector)}
+
+\BlankLine
+Initialization of P\;
+Initialization of Pu\;
+Initialization of the solution vector $Z^{0}$\;
+Computation of the parameters for the $MPI\_AlltoAll$\;
+Allocate memory to GPU\;
+\While {$error > \epsilon$}{
+copy Z from CPU to GPU\;
+$Z^{Prev}_{loc}=kernel\_save(Z_{loc})$\;
+$Z_{loc}=kernel\_update(Z,P,P')$\;
+$\Delta z=kernel\_testConv(Z_{loc},Z^{prev}_{loc})$\;
+$error=MPI\_Reduce(\Delta z)$\;
+Copy $Z_{loc}$ from GPU to CPU\;
+$Z=MPI\_AlltoAll(Z_{loc})$\;
+}
+\RC{A uniformiser avec les autres algos, mais les grandes lignes sont là}
+\end{algorithm}
+