]> AND Private Git Repository - hpcc2014.git/blob - hpcc.tex
Logo AND Algorithmique Numérique Distribuée

Private GIT Repository
1e90946f706b73b3808b1bce4dbd397185e8ecf9
[hpcc2014.git] / hpcc.tex
1 \documentclass[conference]{IEEEtran}
2
3 \usepackage[T1]{fontenc}
4 \usepackage[utf8]{inputenc}
5 \usepackage{amsfonts,amssymb}
6 \usepackage{amsmath}
7 %\usepackage{algorithm}
8 \usepackage{algpseudocode}
9 %\usepackage{amsthm}
10 \usepackage{graphicx}
11 \usepackage[american]{babel}
12 % Extension pour les liens intra-documents (tagged PDF)
13 % et l'affichage correct des UR (commande \url{http://example.com})
14 %\usepackage{hyperref}
15
16 \usepackage{url}
17 \DeclareUrlCommand\email{\urlstyle{same}}
18
19 \usepackage[autolanguage,np]{numprint}
20 \AtBeginDocument{%
21   \renewcommand*\npunitcommand[1]{\text{#1}}
22   \npthousandthpartsep{}}
23
24 \usepackage{xspace}
25 \usepackage[textsize=footnotesize]{todonotes}
26 \newcommand{\AG}[2][inline]{%
27   \todo[color=green!50,#1]{\sffamily\textbf{AG:} #2}\xspace}
28 \newcommand{\DL}[2][inline]{%
29   \todo[color=yellow!50,#1]{\sffamily\textbf{DL:} #2}\xspace}
30 \newcommand{\LZK}[2][inline]{%
31   \todo[color=blue!10,#1]{\sffamily\textbf{LZK:} #2}\xspace}
32 \newcommand{\RC}[2][inline]{%
33   \todo[color=red!10,#1]{\sffamily\textbf{RC:} #2}\xspace}
34 \newcommand{\CER}[2][inline]{%
35   \todo[color=pink!10,#1]{\sffamily\textbf{CER:} #2}\xspace}
36
37 \algnewcommand\algorithmicinput{\textbf{Input:}}
38 \algnewcommand\Input{\item[\algorithmicinput]}
39
40 \algnewcommand\algorithmicoutput{\textbf{Output:}}
41 \algnewcommand\Output{\item[\algorithmicoutput]}
42
43 \newcommand{\MI}{\mathit{MaxIter}}
44 \newcommand{\Time}[1]{\mathit{Time}_\mathit{#1}}
45
46 \begin{document}
47
48 \title{Simulation of Asynchronous Iterative Algorithms Using SimGrid}
49
50 \author{%
51   \IEEEauthorblockN{%
52     Charles Emile Ramamonjisoa\IEEEauthorrefmark{1},
53     Lilia Ziane Khodja\IEEEauthorrefmark{2},
54     David Laiymani\IEEEauthorrefmark{1},
55     Arnaud Giersch\IEEEauthorrefmark{1} and
56     Raphaël Couturier\IEEEauthorrefmark{1}
57   }
58   \IEEEauthorblockA{\IEEEauthorrefmark{1}%
59     Femto-ST Institute -- DISC Department\\
60     Université de Franche-Comté,
61     IUT de Belfort-Montbéliard\\
62     19 avenue du Maréchal Juin, BP 527, 90016 Belfort cedex, France\\
63     Email: \email{{charles.ramamonjisoa,david.laiymani,arnaud.giersch,raphael.couturier}@univ-fcomte.fr}
64   }
65   \IEEEauthorblockA{\IEEEauthorrefmark{2}%
66     Inria Bordeaux Sud-Ouest\\
67     200 avenue de la Vieille Tour, 33405 Talence cedex, France \\
68     Email: \email{lilia.ziane@inria.fr}
69   }
70 }
71
72 \maketitle
73
74 \begin{abstract}
75
76 Synchronous  iterative  algorithms  are  often less  scalable  than  asynchronous
77 iterative  ones.  Performing  large  scale experiments  with  different kind  of
78 network parameters is not easy  because with supercomputers such parameters are
79 fixed. So, one  solution consists in using simulations first  in order to analyze
80 what parameters  could influence or not  the behavior of an  algorithm. In this
81 paper, we show  that it is interesting to use SimGrid  to simulate the behavior
82 of asynchronous  iterative algorithms. For that,  we compare the  behavior of a
83 synchronous  GMRES  algorithm  with  an  asynchronous  multisplitting  one  with
84 simulations  which let us easily choose  some parameters.   Both  codes  are real  MPI
85 codes and simulations allow us to see when the asynchronous multisplitting algorithm can be more
86 efficient than the GMRES one to solve a 3D Poisson problem.
87
88
89 % no keywords for IEEE conferences
90 % Keywords: Algorithm distributed iterative asynchronous simulation SimGrid
91 \end{abstract}
92
93 \section{Introduction}
94
95 Parallel computing and high performance computing (HPC) are becoming  more and more imperative to solve various
96 problems raised by  researchers on various scientific disciplines but also by industrialists in  the field. Indeed, the
97 increasing complexity of these requested  applications combined with a continuous increase of their sizes lead to  write
98 distributed and parallel algorithms requiring significant hardware  resources (grid computing, clusters, broadband
99 network, etc.) but also a non-negligible CPU execution time. We consider in this paper a class of highly efficient
100 parallel algorithms called \emph{iterative algorithms} executed in a distributed environment. As their name
101 suggests, these algorithms solve a given problem by successive iterations ($X_{n +1} = f(X_{n})$) from an initial value
102 $X_{0}$ to find an approximate value $X^*$ of the solution with a very low residual error. Several well-known methods
103 demonstrate the convergence of these algorithms~\cite{BT89,Bahi07}.
104
105 Parallelization of such algorithms generally involves the division of the problem
106 into  several  \emph{blocks}  that  will  be  solved  in  parallel  on  multiple
107 processing units. The latter will communicate each intermediate results before a
108 new  iteration starts  and until  the  approximate solution  is reached.   These
109 parallel computations can be performed either in a \emph{synchronous} mode, where a
110 new iteration  begins only  when all nodes  communications are completed,  or in an
111 \emph{asynchronous}  mode where  processors can  continue independently  with no
112 synchronization points~\cite{bcvc06:ij}. In this case, local computations do not
113 need to  wait for  required data. Processors  can then perform  their iterations
114 with the  data present at that time.  Even if the number of required iterations 
115 before  the convergence  is generally  greater  than in  the synchronous  case,
116 asynchronous  iterative algorithms  can significantly  reduce  overall execution
117 times by  suppressing idle  times due to  synchronizations especially in  a grid
118 computing context (see~\cite{Bahi07} for more details).
119
120 Parallel applications  based on a synchronous or  asynchronous iteration model
121 may have different configuration and deployment requirements.  Quantifying their
122 resource  allocation  policies and  application  scheduling  algorithms in  grid
123 computing environments under varying load,  CPU power and network speeds are very
124 costly,       very        labor       intensive       and        very       time
125 consuming~\cite{Calheiros:2011:CTM:1951445.1951450}.   The case  of asynchronous
126 iterative algorithms  is even more problematic  since they are  very sensitive to
127 the  execution environment  context.  For instance,  variations  in the  network
128 bandwidth (intra and  inter-clusters), in the number and the  power of nodes, in
129 the number  of clusters\dots{} can lead  to very different  number of iterations
130 and so  to very  different execution times.   Then, it  appears that the  use of
131 simulation tools to explore various  platform scenarios and to run large numbers
132 of  experiments quickly  can  be very  promising.  
133
134 Thus, using a simulation environment to execute parallel iterative algorithms can prove to be very interesting to reduce the highly cost  of access to computing resources: (1) for
135 the  applications  development life  cycle  and in  code  debugging  (2) and  in
136 production  to get  results  in a  reasonable  execution time  with a  simulated
137 infrastructure not  accessible with physical  resources. Indeed, to  find optimal  configurations
138 giving  the best  results  with  a lowest  residual  error and  in  the best 
139 execution time is very challenging  for large scale distributed  iterative asynchronous  algorithms
140
141
142 To our knowledge,  there is no existing work on the  large-scale simulation of a
143 real asynchronous  iterative application.  {\bf The contribution  of the present
144   paper can be  summarized in two main points}.  First we  give a first approach
145 of the simulation  of asynchronous iterative algorithms using  a simulation tool
146 (i.e.    the   SimGrid   toolkit~\cite{SimGrid}).    Second,  we   confirm   the
147 efficiency  of the  asynchronous  multisplitting algorithm  by comparing  its
148 performances   with  the   synchronous  GMRES   (Generalized   Minimal  Residual) method
149 \cite{ref1}.  Both  these codes can  be used to  solve large linear  systems. In
150 this  paper, we  focus  on  a 3D  Poisson  problem.  We  show  that, with  minor
151 modifications of the initial MPI code,  the SimGrid toolkit allows us to perform
152 a  test campaign  of  a  real asynchronous  iterative  application on  different
153 computing architectures.
154 % The  simulated results  we
155 %obtained are  in line with real  results exposed in  ??\AG[]{ref?}. 
156 SimGrid  has  allowed us  to  launch the  application  from  a modest  computing
157 infrastructure  by simulating  different distributed  architectures  composed by
158 clusters  nodes interconnected by  variable speed  networks.  Parameters  of the
159 network  platforms  are   the  bandwidth  and  the  latency   of  inter  cluster
160 network. Parameters on the cluster's architecture are the number of machines and
161 the  computation power  of a  machine.  Simulations show  that the  asynchronous
162 multisplitting algorithm  can solve the  3D Poisson problem  approximately twice
163 faster than GMRES with two distant clusters. In this way, we present an original solution to optimize the use of a simulation 
164 tool to run efficiently an  asynchronous iterative parallel algorithm in a grid architecture
165
166
167
168 This article is structured as follows: after this introduction, the next section
169 will  give a  brief  description  of the iterative  asynchronous  model.  Then,  the
170 simulation framework  SimGrid is presented  with the settings to  create various
171 distributed architectures.  Then, the  multisplitting method is presented, it is
172 based  on GMRES to  solve each  block obtained  from the  splitting. This  code is
173 written with MPI  primitives and its adaptation to  SimGrid with SMPI (Simulated
174 MPI) is  detailed in the next  section. At last, the  simulation results carried
175 out will be presented before some concluding remarks and future works.
176
177  
178 \section{Motivations and scientific context}
179
180 As described in the introduction, parallel iterative methods are now widely used
181 in many scientific domains.  They can be classified in three main classes
182 depending on how iterations and communications are managed (for more details
183 readers can refer to~\cite{bcvc06:ij}). In the synchronous iterations model,
184 data are exchanged at the end of each iteration. All the processors must begin
185 the same iteration at the same time and important idle times on processors are
186 generated.  It is possible to use asynchronous communications, in this case, the
187 model can be compared to the previous one except that data required on another
188 processor are sent asynchronously i.e.  without stopping current computations.
189 This technique allows communications to be partially overlapped by computations
190 but unfortunately, the overlapping is only partial and important idle times
191 remain.  It is clear that, in a grid computing context, where the number of
192 computational nodes is large, heterogeneous and widely distributed, the idle
193 times generated by synchronizations are very penalizing. One way to overcome
194 this problem is to use the asynchronous iterations model.  Here, local
195 computations do not need to wait for required data. Processors can then perform
196 their iterations with the data present at that time.  Figure~\ref{fig:aiac}
197 illustrates this model where the gray blocks represent the computation phases.
198 With this algorithmic model, the number of iterations required before the
199 convergence is generally greater than for the two former classes.  But, and as
200 detailed in~\cite{bcvc06:ij}, asynchronous iterative algorithms can
201 significantly reduce overall execution times by suppressing idle times due to
202 synchronizations especially in a grid computing context.
203
204 \begin{figure}[!t]
205   \centering
206     \includegraphics[width=8cm]{AIAC.pdf}
207   \caption{The asynchronous iterations model}
208   \label{fig:aiac}
209 \end{figure}
210
211
212 %% It is very challenging to develop efficient applications for large scale,
213 %% heterogeneous and distributed platforms such as computing grids. Researchers and
214 %% engineers have to develop techniques for maximizing application performance of
215 %% these multi-cluster platforms, by redesigning the applications and/or by using
216 %% novel algorithms that can account for the composite and heterogeneous nature of
217 %% the platform. Unfortunately, the deployment of such applications on these very
218 %% large scale systems is very costly, labor intensive and time consuming. In this
219 %% context, it appears that the use of simulation tools to explore various platform
220 %% scenarios at will and to run enormous numbers of experiments quickly can be very
221 %% promising. Several works\dots{}
222
223 %% \AG{Several works\dots{} what?\\
224 %  Le paragraphe suivant se trouve déjà dans l'intro ?}
225 In the context of asynchronous algorithms, the number of iterations to reach the
226 convergence depends on  the delay of the messages. With  synchronous iterations, the
227 number of  iterations is exactly  the same than  in the sequential mode  (if the
228 parallelization process does  not change the algorithm). So  the difficulty with
229 asynchronous iterative algorithms comes from the fact that it is necessary to run the algorithm
230 with real data. Indeed, from one execution to the other the order of messages will
231 change and the  number of iterations to reach the  convergence will also change.
232 According  to all  the parameters  of the  platform (number  of nodes,  power of
233 nodes,  inter  and  intra clusters  bandwidth  and  latency, etc.) and  of  the
234 algorithm  (number   of  splittings  with  the   multisplitting  algorithm),  the
235 multisplitting code  will obtain the solution  more or less  quickly. Of course,
236 the GMRES method also depends on the same parameters. As it is difficult to have
237 access to  many clusters,  grids or supercomputers  with many  different network
238 parameters,  it  is  interesting  to  be  able  to  simulate  the  behavior  of
239 asynchronous iterative algorithms before being able to run real experiments.
240
241
242
243
244
245
246 \section{SimGrid}
247
248 SimGrid~\cite{SimGrid,casanova+legrand+quinson.2008.simgrid,casanova+giersch+legrand+al.2014.versatile}
249 is a simulation framework to study the behavior of large-scale distributed
250 systems.  As its name suggests, it emanates from the grid computing community,
251 but is nowadays used to study grids, clouds, HPC or peer-to-peer systems.  The
252 early versions of SimGrid date back from 1999, but it is still actively
253 developed and distributed as an open source software.  Today, it is one of the
254 major generic tools in the field of simulation for large-scale distributed
255 systems.
256
257 SimGrid provides several programming interfaces: MSG to simulate Concurrent
258 Sequential Processes, SimDAG to simulate DAGs of (parallel) tasks, and SMPI to
259 run real applications written in MPI~\cite{MPI}.  Apart from the native C
260 interface, SimGrid provides bindings for the C++, Java, Lua and Ruby programming
261 languages.  SMPI is the interface that has been used for the work described in
262 this paper.  The SMPI interface implements about \np[\%]{80} of the MPI 2.0
263 standard~\cite{bedaride+degomme+genaud+al.2013.toward}, and supports
264 applications written in C or Fortran, with little or no modifications.
265
266 Within SimGrid, the execution of a distributed application is simulated by a
267 single process.  The application code is really executed, but some operations,
268 like communications, are intercepted, and their running time is computed
269 according to the characteristics of the simulated execution platform.  The
270 description of this target platform is given as an input for the execution, by
271 means of an XML file.  It describes the properties of the platform, such as
272 the computing nodes with their computing power, the interconnection links with
273 their bandwidth and latency, and the routing strategy.  The scheduling of the
274 simulated processes, as well as the simulated running time of the application
275 are computed according to these properties.
276
277 To compute the durations of the operations in the simulated world, and to take
278 into account resource sharing (e.g. bandwidth sharing between competing
279 communications), SimGrid uses a fluid model.  This allows users to run relatively fast
280 simulations, while still keeping accurate
281 results~\cite{bedaride+degomme+genaud+al.2013.toward,
282   velho+schnorr+casanova+al.2013.validity}.  Moreover, depending on the
283 simulated application, SimGrid/SMPI allows to skip long lasting computations and
284 to only take their duration into account.  When the real computations cannot be
285 skipped, but the results are unimportant for the simulation results, it is
286 also possible to share dynamically allocated data structures between
287 several simulated processes, and thus to reduce the whole memory consumption.
288 These two techniques can help to run simulations on a very large scale.
289
290 The validity of simulations with SimGrid has been asserted by several studies.
291 See, for example, \cite{velho+schnorr+casanova+al.2013.validity} and articles
292 referenced therein for the validity of the network models.  Comparisons between
293 real execution of MPI applications on the one hand, and their simulation with
294 SMPI on the other hand, are presented in~\cite{guermouche+renard.2010.first,
295   clauss+stillwell+genaud+al.2011.single,
296   bedaride+degomme+genaud+al.2013.toward}.  All these works conclude that
297 SimGrid is able to simulate pretty accurately the real behavior of the
298 applications.
299
300
301 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
302 \section{Simulation of the multisplitting method}
303
304 \subsection{The multisplitting method}
305 %Décrire le problème (algo) traité ainsi que le processus d'adaptation à SimGrid.
306 Let $Ax=b$ be a large sparse system of $n$ linear equations in $\mathbb{R}$, where $A$ is a sparse square and nonsingular matrix, $x$ is the solution vector and $b$ is the right-hand side vector. We use a multisplitting method based on the block Jacobi splitting to solve this linear system on a large scale platform composed of $L$ clusters of processors~\cite{o1985multi}. In this case, we apply a row-by-row splitting without overlapping  
307 \begin{equation*}
308   \left(\begin{array}{ccc}
309       A_{11} & \cdots & A_{1L} \\
310       \vdots & \ddots & \vdots\\
311       A_{L1} & \cdots & A_{LL}
312     \end{array} \right)
313   \times
314   \left(\begin{array}{c}
315       X_1 \\
316       \vdots\\
317       X_L
318     \end{array} \right)
319   =
320   \left(\begin{array}{c}
321       B_1 \\
322       \vdots\\
323       B_L
324     \end{array} \right)
325 \end{equation*}
326 in such a way that successive rows of matrix $A$ and both vectors $x$ and $b$
327 are assigned to one cluster, where for all $\ell,m\in\{1,\ldots,L\}$, $A_{\ell
328   m}$ is a rectangular block of $A$ of size $n_\ell\times n_m$, $X_\ell$ and
329 $B_\ell$ are sub-vectors of $x$ and $b$, respectively, of size $n_\ell$ each,
330 and $\sum_{\ell} n_\ell=\sum_{m} n_m=n$.
331
332 The multisplitting method proceeds by iteration to solve in parallel the linear system on $L$ clusters of processors, in such a way each sub-system
333 \begin{equation}
334   \label{eq:4.1}
335   \left\{
336     \begin{array}{l}
337       A_{\ell\ell}X_\ell = Y_\ell \text{, such that}\\
338       Y_\ell = B_\ell - \displaystyle\sum_{\substack{m=1\\ m\neq \ell}}^{L}A_{\ell m}X_m
339     \end{array}
340   \right.
341 \end{equation}
342 is solved independently by a cluster and communications are required to update
343 the right-hand side sub-vector $Y_\ell$, such that the sub-vectors $X_m$
344 represent the data dependencies between the clusters. As each sub-system
345 (\ref{eq:4.1}) is solved in parallel by a cluster of processors, our
346 multisplitting method uses an iterative method as an inner solver which is
347 easier to parallelize and more scalable than a direct method. In this work, we
348 use the parallel algorithm of GMRES method~\cite{ref1} which is one of the most
349 used iterative method by many researchers.
350
351 \begin{figure}[!t]
352   %%% IEEE instructions forbid to use an algorithm environment here, use figure
353   %%% instead
354 \begin{algorithmic}[1]
355 \Input $A_\ell$ (sparse sub-matrix), $B_\ell$ (right-hand side sub-vector)
356 \Output $X_\ell$ (solution sub-vector)\medskip
357
358 \State Load $A_\ell$, $B_\ell$
359 \State Set the initial guess $x^0$
360 \For {$k=0,1,2,\ldots$ until the global convergence}
361 \State Restart outer iteration with $x^0=x^k$
362 \State Inner iteration: \Call{InnerSolver}{$x^0$, $k+1$}
363 \State\label{algo:01:send} Send shared elements of $X_\ell^{k+1}$ to neighboring clusters
364 \State\label{algo:01:recv} Receive shared elements in $\{X_m^{k+1}\}_{m\neq \ell}$
365 \EndFor
366
367 \Statex
368
369 \Function {InnerSolver}{$x^0$, $k$}
370 \State Compute local right-hand side $Y_\ell$:
371        \begin{equation*}
372          Y_\ell = B_\ell - \sum\nolimits^L_{\substack{m=1\\ m\neq \ell}}A_{\ell m}X_m^0
373        \end{equation*}
374 \State Solving sub-system $A_{\ell\ell}X_\ell^k=Y_\ell$ with the parallel GMRES method
375 \State \Return $X_\ell^k$
376 \EndFunction
377 \end{algorithmic}
378 \caption{A multisplitting solver with GMRES method}
379 \label{algo:01}
380 \end{figure}
381
382 The algorithm in Figure~\ref{algo:01} shows the main key points of the
383 multisplitting method to solve a large sparse linear system. This algorithm is
384 based on an outer-inner iteration method where the parallel synchronous GMRES
385 method is used to solve the inner iteration. It is executed in parallel by each
386 cluster of processors. For all $\ell,m\in\{1,\ldots,L\}$, the matrices and
387 vectors with the subscript $\ell$ represent the local data for cluster $\ell$,
388 while $\{A_{\ell m}\}_{m\neq \ell}$ are off-diagonal matrices of sparse matrix
389 $A$ and $\{X_m\}_{m\neq \ell}$ contain vector elements of solution $x$ shared
390 with neighboring clusters. At every outer iteration $k$, asynchronous
391 communications are performed between processors of the local cluster and those
392 of distant clusters (lines~\ref{algo:01:send} and~\ref{algo:01:recv} in
393 Figure~\ref{algo:01}). The shared vector elements of the solution $x$ are
394 exchanged by message passing using MPI non-blocking communication routines.
395
396 \begin{figure}[!t]
397 \centering
398   \includegraphics[width=60mm,keepaspectratio]{clustering}
399 \caption{Example of three distant clusters of processors.}
400 \label{fig:4.1}
401 \end{figure}
402
403 The global convergence of the asynchronous multisplitting solver is detected
404 when the clusters of processors have all converged locally. We implemented the
405 global convergence detection process as follows. On each cluster a master
406 processor is designated (for example the processor with rank 1) and masters of
407 all clusters are interconnected by a virtual unidirectional ring network (see
408 Figure~\ref{fig:4.1}). During the resolution, a Boolean token circulates around
409 the virtual ring from a master processor to another until the global convergence
410 is achieved. So, starting from the cluster with rank 1, each master processor $\ell$
411 sets the token to \textit{True} if the local convergence is achieved or to
412 \textit{False} otherwise, and sends it to master processor $\ell+1$. Finally, the
413 global convergence is detected when the master of cluster 1 receives from the
414 master of cluster $L$ a token set to \textit{True}. In this case, the master of
415 cluster 1 broadcasts a stop message to the masters of other clusters. In this work,
416 the local convergence on each cluster $\ell$ is detected when the following
417 condition is satisfied
418 \begin{equation*}
419   (k=\MI) \text{ or } (\|X_\ell^k - X_\ell^{k+1}\|_{\infty}\leq\epsilon)
420 \end{equation*}
421 where $\MI$ is the maximum number of outer iterations and $\epsilon$ is the
422 tolerance threshold of the error computed between two successive local solution
423 $X_\ell^k$ and $X_\ell^{k+1}$. It should be noted that with asynchronous iterative algorithms, we cannot use a classical norm (which would require to synchronize all processors), such as $\|X_\ell^k - X_\ell^{k+1}\|_{2}$ for example. Nevertheless, in our experiments, we check that the final result is correct, for this we compute the precision with $max_i | A*x-b |_i$.
424
425
426
427 In this paper, we solve the 3D Poisson problem whose mathematical model is 
428 \begin{equation}
429 \left\{
430 \begin{array}{l}
431 \nabla^2 u = f \text{~in~} \Omega \\
432 u =0 \text{~on~} \Gamma =\partial\Omega
433 \end{array}
434 \right.
435 \label{eq:02}
436 \end{equation}
437 where $\nabla^2$ is the Laplace operator, $f$ and $u$ are real-valued functions, and $\Omega=[0,1]^3$. The spatial discretization with a finite differences scheme reduces problem~(\ref{eq:02}) to a system of sparse linear equations. Our multisplitting method solves the 3D Poisson problem using a seven point stencil whose general expression could be written as
438 \begin{equation}
439 \begin{array}{l}
440 u(x-1,y,z) + u(x,y-1,z) + u(x,y,z-1)\\+u(x+1,y,z)+u(x,y+1,z)+u(x,y,z+1) \\ -6u(x,y,z)=h^2f(x,y,z),
441 %u(x,y,z)= & \frac{1}{6}\times [u(x-1,y,z) + u(x+1,y,z) + \\
442  %         & u(x,y-1,z) + u(x,y+1,z) + \\
443   %        & u(x,y,z-1) + u(x,y,z+1) - \\ & h^2f(x,y,z)],
444 \end{array}
445 \label{eq:03}
446 \end{equation} 
447 where $h$ is the distance between two adjacent elements in the spatial discretization scheme and the iteration matrix $A$ of size $N_x\times N_y\times N_z$ of the discretized linear system is sparse, symmetric and positive definite. 
448
449 The parallel solving of the 3D Poisson problem with our multisplitting method requires a data partitioning of the problem between clusters and between processors within a cluster. We have chosen the 3D partitioning instead of the row-by-row partitioning one in order to reduce the data exchanges at sub-domain boundaries. Figure~\ref{fig:4.2} shows an example of the data partitioning of the 3D Poisson problem between two clusters of processors, where each sub-problem is assigned to a processor. In this context, a processor has at most six neighbors within a cluster or in distant clusters with which it shares data at sub-domain boundaries. 
450
451 \begin{figure}[!t]
452 \centering
453   \includegraphics[width=80mm,keepaspectratio]{partition}
454 \caption{Example of the 3D data partitioning between two clusters of processors.}
455 \label{fig:4.2}
456 \end{figure}
457
458
459 \subsection{Simulation of the multisplitting method using SimGrid and SMPI}
460
461
462
463 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
464 We did not encounter major blocking problems when adapting the multisplitting algorithm previously described to a simulation environment like SimGrid unless some code 
465 debugging. Indeed, apart from the review of the program sequence for asynchronous exchanges between processors within a cluster or between clusters, the algorithm was executed successfully with SMPI and provided identical outputs as those obtained with direct execution under MPI. For the synchronous GMRES method, the execution of the program raised no particular issue but in the asynchronous multisplitting method, the review of the sequence of \texttt{MPI\_Isend, MPI\_Irecv} and \texttt{MPI\_Waitall} instructions
466 and with the addition of the primitive \texttt{MPI\_Test} was needed to avoid a memory fault due to an infinite loop resulting from the non-convergence of the algorithm.
467 %\CER{On voulait en fait montrer la simplicité de l'adaptation de l'algo a SimGrid. Les problèmes rencontrés décrits dans ce paragraphe concerne surtout le mode async}\LZK{OK. J'aurais préféré avoir un peu plus de détails sur l'adaptation de la version async} 
468 %\CER{Le problème majeur sur l'adaptation MPI vers SMPI pour la partie asynchrone de l'algorithme a été le plantage en SMPI de Waitall après un Isend et Irecv. J'avais proposé un workaround en utilisant un MPI\_wait séparé pour chaque échange a la place d'un waitall unique pour TOUTES les échanges, une instruction qui semble bien fonctionner en MPI. Ce workaround aussi fonctionne bien. Mais après, tu as modifié le programme avec l'ajout d'un MPI\_Test, au niveau de la routine de détection de la convergence et du coup, l'échange global avec waitall a aussi fonctionné.}
469 Note here that the use of SMPI functions optimizer for memory footprint and CPU usage is not recommended knowing that one wants to get real results by simulation.
470 As mentioned, upon this adaptation, the algorithm is executed as in real life in the simulated environment after the following minor changes. The scope of all declared 
471 global variables have been moved to local subroutines. Indeed, global variables generate side effects arising from the concurrent access of 
472 shared memory used by threads simulating each computing unit in the SimGrid architecture. 
473 %Second, some compilation errors on MPI\_Waitall and MPI\_Finalize primitives have been fixed with the latest version of SimGrid.
474 %\AG{compilation or run-time error?}
475 In total, the initial MPI program running on the simulation environment SMPI gave after a very simple adaptation the same results as those obtained in a real 
476 environment. We have successfully executed the code for the synchronous GMRES algorithm compared with our asynchronous multisplitting algorithm after few modifications. 
477
478
479
480 \section{Simulation results}
481
482 When the \textit{real} application runs in the simulation environment and produces the expected results, varying the input
483 parameters and the program arguments allows us to compare outputs from the code execution. We have noticed from this
484 study that the results depend on the following parameters:  
485 \begin{itemize}
486 \item At the network level, we found that the most critical values are the
487   bandwidth and the network latency.
488 \item Host processor power (GFlops) can also influence the results.
489 \item Finally, when submitting job batches for execution, the arguments values
490   passed to the program like the maximum number of iterations or the precision are critical. They allow us to ensure not only the convergence of the
491   algorithm but also to get the main objective in getting an execution time with the asynchronous multisplitting  less than with synchronous GMRES. 
492   \end{itemize}
493
494 The ratio between the simulated execution time of synchronous GMRES algorithm
495 compared to the asynchronous multisplitting algorithm ($t_\text{GMRES} / t_\text{Multisplitting}$) is defined as the \emph{relative gain}. So,
496 our objective running the algorithm in SimGrid is to obtain a relative gain greater than 1.
497 A priori, obtaining a relative gain greater than 1 would be difficult in a local
498 area network configuration where the synchronous GMRES method will take advantage on the
499 rapid exchange of information on such high-speed links. Thus, the methodology
500 adopted was to launch the application on a clustered network. In this
501 configuration, degrading the inter-cluster network performance will penalize the
502 synchronous mode allowing to get a relative gain greater than 1.  This action
503 simulates the case of distant clusters linked with long distance networks as in grid computing context.
504
505
506
507 Both codes were simulated on a two clusters based network with 50 hosts each, totalling 100 hosts. Various combinations of the above
508 factors have provided the results shown in Table~\ref{tab.cluster.2x50}. The problem size of the 3D Poisson problem  ranges from $N=N_x = N_y = N_z = \text{62}$ to 150 elements (that is from
509 $\text{62}^\text{3} = \text{\np{238328}}$ to $\text{150}^\text{3} =
510 \text{\np{3375000}}$ entries). With the asynchronous multisplitting algorithm the simulated execution time is on average 2.5 times faster than with the synchronous GMRES one. 
511 %\AG{Expliquer comment lire les tableaux.}
512 %\CER{J'ai reformulé la phrase par la lecture du tableau. Plus de détails seront lus dans la partie Interprétations et commentaires}
513 % use the same column width for the following three tables
514 \newlength{\mytablew}\settowidth{\mytablew}{\footnotesize\np{E-11}}
515 \newenvironment{mytable}[1]{% #1: number of columns for data
516   \renewcommand{\arraystretch}{1.3}%
517   \begin{tabular}{|>{\bfseries}r%
518                   |*{#1}{>{\centering\arraybackslash}p{\mytablew}|}}}{%
519     \end{tabular}}
520
521 \begin{table}[!t]
522   \centering
523   \caption{Relative gain  of the multisplitting algorithm compared  to GMRES for
524     different configurations with 2 clusters, each one composed of 50 nodes. Latency = $20$ms}
525   \label{tab.cluster.2x50}
526
527   \begin{mytable}{5}
528     \hline
529     bandwidth (Mbit/s)
530     & 5         & 5         & 5         & 5         & 5         \\
531     \hline
532   %  latency (ms)
533    % & 20      &  20      & 20      & 20      & 20      \\
534     %\hline
535     power (GFlops)
536     & 1         & 1         & 1         & 1.5       & 1.5       \\
537     \hline
538     size $(N)$
539     & $62^3$        & $62^3$        & $62^3$        & $100^3$       & $100^3$       \\
540     \hline
541     Precision
542     & \np{E-5}  & \np{E-8}  & \np{E-9}  & \np{E-11} & \np{E-11} \\
543     \hline
544     \hline
545     Relative gain
546     & 2.52      & 2.55      & 2.52      & 2.57      & 2.54      \\
547     \hline
548   \end{mytable}
549
550   \bigskip
551
552   \begin{mytable}{5}
553     \hline
554     bandwidth (Mbit/s)
555     & 50        & 50        & 50        & 50        & 50 \\ %       & 10        & 10 \\
556     \hline
557     %latency (ms)
558     %& 20      & 20      & 20      & 20      & 20 \\ %      & 0.03      & 0.01 \\
559     %\hline
560     Power (GFlops)
561     & 1.5       & 1.5       & 1.5       & 1.5       & 1.5 \\ %      & 1         & 1.5 \\
562     \hline
563     size $(N)$
564     & $110^3$       & $120^3$       & $130^3$       & $140^3$       & $150^3$  \\ %     & 171       & 171 \\
565     \hline
566     Precision
567     & \np{E-11} & \np{E-11} & \np{E-11} & \np{E-11} & \np{E-11} \\ % & \np{E-5}  & \np{E-5} \\
568     \hline
569     \hline
570     Relative gain
571     & 2.53      & 2.51     & 2.58     & 2.55     & 2.54   \\ %  & 1.59      & 1.29 \\
572     \hline
573   \end{mytable}
574 \end{table}
575   
576 %\RC{Du coup la latence est toujours la même, pourquoi la mettre dans la table?}
577
578 %Then we have changed the network configuration using three clusters containing
579 %respectively 33, 33 and 34 hosts, or again by on hundred hosts for all the
580 %clusters. In the same way as above, a judicious choice of key parameters has
581 %permitted to get the results in Table~\ref{tab.cluster.3x33} which shows the
582 %relative gains greater than 1 with a matrix size from 62 to 100 elements.
583
584 %\CER{En accord avec RC, on a pour le moment enlevé les tableaux 2 et 3 sachant que les résultats obtenus sont limites. De même, on a enlevé aussi les deux dernières colonnes du tableau I en attendant une meilleure performance et une meilleure precision}
585 %\begin{table}[!t]
586 %  \centering
587 %  \caption{3 clusters, each with 33 nodes}
588 %  \label{tab.cluster.3x33}
589 %
590 %  \begin{mytable}{6}
591 %    \hline
592 %    bandwidth 
593 %    & 10       & 5        & 4        & 3        & 2        & 6 \\
594 %    \hline
595 %    latency
596 %    & 0.01     & 0.02     & 0.02     & 0.02     & 0.02     & 0.02 \\
597 %    \hline
598 %    power
599 %    & 1        & 1        & 1        & 1        & 1        & 1 \\
600 %    \hline
601 %    size
602 %    & 62       & 100      & 100      & 100      & 100      & 171 \\
603 %    \hline
604 %    Prec/Eprec
605 %    & \np{E-5} & \np{E-5} & \np{E-5} & \np{E-5} & \np{E-5} & \np{E-5} \\
606 %    \hline
607 %    \hline
608 %    Relative gain
609 %    & 1.003    & 1.01     & 1.08     & 1.19     & 1.28     & 1.01 \\
610 %    \hline
611 %  \end{mytable}
612 %\end{table}
613
614 %In a final step, results of an execution attempt to scale up the three clustered
615 %configuration but increasing by two hundreds hosts has been recorded in
616 %Table~\ref{tab.cluster.3x67}.
617
618 %\begin{table}[!t]
619 %  \centering
620 %  \caption{3 clusters, each with 66 nodes}
621 %  \label{tab.cluster.3x67}
622 %
623 %  \begin{mytable}{1}
624 %    \hline
625 %    bandwidth  & 1 \\
626 %    \hline
627 %    latency    & 0.02 \\
628 %    \hline
629 %    power      & 1 \\
630 %    \hline
631 %    size       & 62 \\
632 %    \hline
633 %    Prec/Eprec & \np{E-5} \\
634 %    \hline
635 %    \hline
636 %    Relative gain    & 1.11 \\
637 %    \hline
638 %  \end{mytable}
639 %\end{table}
640
641 Note that the program was run with the following parameters:
642
643 \paragraph*{SMPI parameters}
644
645 \begin{itemize}
646 \item HOSTFILE: Text file containing the list of the processors units name. Here 100 hosts;
647 \item PLATFORM: XML file description of the platform architecture with the
648   following characteristics:
649   % two clusters (cluster1 and cluster2) with the following characteristics:
650   \begin{itemize}
651   \item 2 clusters of 50 hosts each;
652   \item Processor unit power: \np[GFlops]{1} or \np[GFlops]{1.5};
653   \item Intra-cluster network bandwidth: \np[Gbit/s]{1.25} and latency: \np[$\mu$s]{50};
654   \item Inter-cluster network bandwidth: \np[Mbit/s]{5} or \np[Mbit/s]{50} and latency: \np[ms]{20};
655   \end{itemize}
656 \end{itemize}
657
658
659 \paragraph*{Arguments of the program}
660
661 \begin{itemize}
662 \item Description of the cluster architecture matching the format <Number of
663   clusters> <Number of hosts in cluster1> <Number of hosts in cluster2>;
664 \item Maximum numbers of outer and inner iterations;
665 \item Outer and inner precisions on the residual error;
666 \item Matrix size $N_x$, $N_y$ and $N_z$;
667 \item Matrix diagonal value: $6$ (see Equation~(\ref{eq:03}));
668 \item Matrix off-diagonal values: $-1$;
669 \item Communication mode: asynchronous.
670 \end{itemize}
671
672 \paragraph*{Interpretations and comments}
673
674 After analyzing the outputs, generally, for the two clusters including one hundred hosts configuration (Tables~\ref{tab.cluster.2x50}), some combinations of parameters affecting
675 the results, have given a relative gain of more than 2.5, showing the effectiveness of the
676 asynchronous multisplitting  compared to GMRES with two distant clusters.
677
678 With these settings, Table~\ref{tab.cluster.2x50} shows
679 that after setting the bandwidth of the  inter cluster network to  \np[Mbit/s]{5}, the latency to $20$ millisecond and the processor power
680 to one GFlops, an efficiency of about \np[\%]{40} is
681 obtained in asynchronous mode for a matrix size of $62^3$ elements. It is noticed that the result remains
682 stable even if the residual error precision varies from \np{E-5} to \np{E-9}. By
683 increasing the matrix size up to $100^3$ elements, it was necessary to increase the
684 CPU power by \np[\%]{50} to \np[GFlops]{1.5} to get the algorithm convergence and the same order of asynchronous mode efficiency.  Maintaining  a relative gain of $2.5$ and such processor power but increasing network throughput inter cluster up to \np[Mbit/s]{50},  is obtained with
685 high external precision of \np{E-11} for a matrix size from $110^3$ to $150^3$ side
686 elements. 
687
688 %For the 3 clusters architecture including a total of 100 hosts,
689 %Table~\ref{tab.cluster.3x33} shows that it was difficult to have a combination
690 %which gives a relative gain of asynchronous mode more than 1.2. Indeed, for a
691 %matrix size of 62 elements, equality between the performance of the two modes
692 %(synchronous and asynchronous) is achieved with an inter cluster of
693 %\np[Mbit/s]{10} and a latency of \np[ms]{E-1}. To challenge an efficiency greater than 1.2 with a matrix %size of 100 points, it was necessary to degrade the
694 %inter cluster network bandwidth from 5 to \np[Mbit/s]{2}.
695 %\AG{Conclusion, on prend une plateforme pourrie pour avoir un bon ratio sync/async ???
696   %Quelle est la perte de perfs en faisant ça ?}
697
698 %A last attempt was made for a configuration of three clusters but more powerful
699 %with 200 nodes in total. The convergence with a relative gain around 1.1 was
700 %obtained with a bandwidth of \np[Mbit/s]{1} as shown in
701 %Table~\ref{tab.cluster.3x67}.
702
703 %\RC{Est ce qu'on sait expliquer pourquoi il y a une telle différence entre les résultats avec 2 et 3 clusters... Avec 3 clusters, ils sont pas très bons... Je me demande s'il ne faut pas les enlever...}
704 %\RC{En fait je pense avoir la réponse à ma remarque... On voit avec les 2 clusters que le gain est d'autant plus grand qu'on choisit une bonne précision. Donc, plusieurs solutions, lancer rapidement un long test pour confirmer ca, ou enlever des tests... ou on ne change rien :-)}
705 %\LZK{Ma question est: le bandwidth et latency sont ceux inter-clusters ou pour les deux inter et intra cluster??}
706 %\CER{Définitivement, les paramètres réseaux variables ici se rapportent au réseau INTER cluster.}
707 \section{Conclusion}
708 The simulation of the execution of parallel asynchronous iterative algorithms on large scale  clusters has been presented. 
709 In this work, we show that SimGrid is an efficient simulation tool that has enabled us to 
710 reach the following two objectives: 
711
712 \begin{enumerate}
713 \item  To have  a flexible  configurable execution  platform that  allows  us to
714   simulate algorithms for  which execution of all parts of
715   the  code is  necessary. Using  simulations before  real executions  is  a nice
716   solution to detect potential scalability problems.
717
718 \item To test the combination of the cluster and network specifications permitting to execute an asynchronous algorithm faster than a synchronous one.
719 \end{enumerate}
720 Our results have shown that with two distant clusters, the asynchronous multisplitting method is faster by \np[\%]{40} compared to the synchronous GMRES method
721 which is not negligible for solving complex practical problems with ever increasing size.
722
723 Several studies have already addressed the performance execution time of 
724 this class of algorithm. The work presented in this paper has 
725 demonstrated an original solution to optimize the use of a simulation 
726 tool to run efficiently an iterative parallel algorithm in asynchronous 
727 mode in a grid architecture. 
728
729 In future works, we plan to extend our experimentations to larger scale platforms by increasing the number of computing cores and the number of clusters. 
730 We will also have to increase the size of the input problem which will require the use of a more powerful simulation platform. At last, we expect to compare our simulation results to real execution results on real architectures in order to better experimentally validate our study. Finally, we also plan to study other problems with the multisplitting method and other asynchronous iterative methods.
731
732 \section*{Acknowledgment}
733
734 This work is partially funded by the Labex ACTION program (contract ANR-11-LABX-01-01).
735 %\todo[inline]{The authors would like to thank\dots{}}
736
737 % trigger a \newpage just before the given reference
738 % number - used to balance the columns on the last page
739 % adjust value as needed - may need to be readjusted if
740 % the document is modified later
741 \bibliographystyle{IEEEtran}
742 \bibliography{IEEEabrv,hpccBib}
743
744
745
746 \end{document}
747
748 %%% Local Variables:
749 %%% mode: latex
750 %%% TeX-master: t
751 %%% fill-column: 80
752 %%% ispell-local-dictionary: "american"
753 %%% End:
754
755 % LocalWords:  Ramamonjisoa Laiymani Arnaud Giersch Ziane Khodja Raphaël Femto
756 % LocalWords:  Université Franche Comté IUT Montbéliard Maréchal Juin Inria Sud
757 % LocalWords:  Ouest Vieille Talence cedex scalability experimentations HPC MPI
758 % LocalWords:  Parallelization AIAC GMRES multi SMPI SISC SIAC SimDAG DAGs Lua
759 % LocalWords:  Fortran GFlops priori Mbit de du fcomte multisplitting scalable
760 % LocalWords:  SimGrid Belfort parallelize Labex ANR LABX IEEEabrv hpccBib Gbit
761 % LocalWords:  intra durations nonsingular Waitall discretization discretized
762 % LocalWords:  InnerSolver Isend Irecv parallelization