-The code is organzed by what is named kernels, portions o code that are run on GPU devices. For step 3, there are two kernels, the
-first named \textit{save} is used to save vector $Z^{K-1}$ and the seconde one is named
-\textit{update} and is used to update the $Z^{K}$ vector. For step 4, a kernel
-tests the convergence of the method. In order to
-compute the function H, we have two possibilities: either to use
-the Jacobi mode, or the Gauss-Seidel mode of iterating which uses the
-most recent computed roots. It is well known that the Gauss-
-Seidel mode converges more quickly. So, we used the Gauss-Seidel mode of iteration. To
-parallelize the code, we created kernels and many functions to
-be executed on the GPU for all the operations dealing with the
+\begin{equation}
+\label{Log_H2}
+z^{k+1}_{i}=z_{i}^{k}-\exp \left(\ln \left(
+p(z_{i}^{k})\right)-\ln\left(p'(z^{k}_{i})\right)- \ln\left(1-Q(z^{k}_{i})\right)\right),
+\end{equation}
+
+where:
+
+\begin{eqnarray}
+\label{Log_H1}
+Q(z^{k}_{i})=\exp\left( \ln (p(z^{k}_{i}))-\ln(p'(z^{k}_{i}))+\ln \left(
+\sum_{i\neq j}^{n}\frac{1}{z^{k}_{i}-z^{k}_{j}}\right)\right) \nonumber \\
+i=1,...,n
+\end{eqnarray}
+
+
+%We propose to use the logarithm and the exponential of a complex in order to compute the power at a high exponent.
+Using the logarithm and the exponential operators, we can replace any multiplications and divisions with additions and subtractions. Consequently, computations manipulate lower absolute values~\cite{Karimall98}.
+
+%This problem was discussed earlier in~\cite{Karimall98} for the Durand-Kerner method. The authors
+%propose to use the logarithm and the exponential of a complex in order to compute the power at a high exponent. Using the logarithm and the exponential operators, we can replace any multiplications and divisions with additions and subtractions. Consequently, computations manipulate lower absolute values and the roots for large polynomial degrees can be looked for successfully~\cite{Karimall98}.
+
+\subsection{Ehrlich-Aberth parallel implementation on CUDA}
+%We introduced three paradigms of parallel programming.
+
+Our objective consists in implementing a root finding polynomial
+algorithm on multiple GPUs. To this end, it is primordial to know how
+to manage CUDA contexts of different GPUs. A direct method for
+controlling the various GPUs is to use as many threads or processes as
+GPU devices. We can choose the GPU index based on the identifier of
+OpenMP thread or the rank of the MPI process. Both approaches will be
+investigated.
+
+
+
+
+Like any parallel code, a GPU parallel implementation first requires
+to determine the sequential tasks and the parallelizable parts of the
+sequential version of the program/algorithm. In our case, all the
+operations that are easy to execute in parallel must be made by the
+GPU to accelerate the execution of the application, like the step 3
+and step 4. On the other hand, all the sequential operations and the
+operations that have data dependencies between threads or recursive
+computations must be executed by only one CUDA or CPU thread (step 1
+and step 2). Initially, we specify the organization of parallel
+threads, by specifying the dimension of the grid Dimgrid, the number
+of blocks per grid DimBlock and the number of threads per block.
+
+The code is organized kernels which are part of code that are run on
+GPU devices. For step 3, there are two kernels, the first named
+\textit{save} is used to save vector $Z^{K-1}$ and the second one is
+named \textit{update} and is used to update the $Z^{K}$ vector. For
+step 4, a kernel tests the convergence of the method. In order to
+compute the function H, we have two possibilities: either to use the
+Jacobi mode, or the Gauss-Seidel mode of iterating which uses the most
+recent computed roots. It is well known that the Gauss-Seidel mode
+converges more quickly. So, we use Gauss-Seidel iterations. To
+parallelize the code, we create kernels and many functions to be
+executed on the GPU for all the operations dealing with the