3 ===============================
4 SMPI: Simulate MPI Applications
5 ===============================
9 <object id="TOC" data="graphical-toc.svg" width="100%" type="image/svg+xml"></object>
11 window.onload=function() { // Wait for the SVG to be loaded before changing it
12 var elem=document.querySelector("#TOC").contentDocument.getElementById("SMPIBox")
13 elem.style="opacity:0.93999999;fill:#ff0000;fill-opacity:0.1";
14 elem.style="opacity:0.93999999;fill:#ff0000;fill-opacity:0.1;stroke:#000000;stroke-width:0.35277778;stroke-linecap:round;stroke-linejoin:round;stroke-miterlimit:4;stroke-dasharray:none;stroke-dashoffset:0;stroke-opacity:1";
20 SMPI enables the study of MPI application by emulating them on top of
21 the SimGrid simulator. This is particularly interesting to study
22 existing MPI applications within the comfort of the simulator.
24 To get started with SMPI, you should head to :ref:`the SMPI tutorial
25 <usecase_smpi>`. You may also want to read the `SMPI reference
26 article <https://hal.inria.fr/hal-01415484>`_ or these `introductory
27 slides <http://simgrid.org/tutorials/simgrid-smpi-101.pdf>`_. If you
28 are new to MPI, you should first take our online `SMPI CourseWare
29 <https://simgrid.github.io/SMPI_CourseWare/>`_. It consists in several
30 projects that progressively introduce the MPI concepts. It proposes to
31 use SimGrid and SMPI to run the experiments, but the learning
32 objectives are centered on MPI itself.
34 Our goal is to enable the study of **unmodified MPI applications**.
35 Some constructs and features are still missing, but we can probably
36 add them on demand. If you already used MPI before, SMPI should sound
37 very familiar to you: Use smpicc instead of mpicc, and smpirun instead
38 of mpirun. The main difference is that smpirun takes a :ref:`simulated
39 platform <platform>` as an extra parameter.
41 For **further scalability**, you may modify your code to speed up your
42 studies or save memory space. Maximal **simulation accuracy**
43 requires some specific care from you.
51 In this mode, your application is actually executed. Every computation
52 occurs for real while every communication is simulated. In addition,
53 the executions are automatically benchmarked so that their timings can
54 be applied within the simulator.
56 SMPI can also go offline by replaying a trace. :ref:`Trace replay
57 <SMPI_offline>` is usually ways faster than online simulation (because
58 the computation are skipped), but it can only applied to applications
59 with constant execution and communication patterns (for the exact same
66 If your application is in C, then simply use ``smpicc`` as a
67 compiler just like you use mpicc with other MPI implementations. This
68 script still calls your default compiler (gcc, clang, ...) and adds
69 the right compilation flags along the way. If your application is in
70 C++, Fortran 77 or Fortran 90, use respectively ``smpicxx``,
71 ``smpiff`` or ``smpif90``.
77 Use the ``smpirun`` script as follows:
81 smpirun -hostfile my_hostfile.txt -platform my_platform.xml ./program -blah
83 - ``my_hostfile.txt`` is a classical MPI hostfile (that is, this file
84 lists the machines on which the processes must be dispatched, one
86 - ``my_platform.xml`` is a classical SimGrid platform file. Of course,
87 the hosts of the hostfile must exist in the provided platform.
88 - ``./program`` is the MPI program to simulate, that you compiled with ``smpicc``
89 - ``-blah`` is a command-line parameter passed to this program.
91 ``smpirun`` accepts other parameters, such as ``-np`` if you don't
92 want to use all the hosts defined in the hostfile, ``-map`` to display
93 on which host each rank gets mapped of ``-trace`` to activate the
94 tracing during the simulation. You can get the full list by running
97 ...............................
98 Debugging your Code within SMPI
99 ...............................
101 If you want to explore the automatic platform and deployment files
102 that are generated by ``smpirun``, add ``-keep-temps`` to the command
105 You can also run your simulation within valgrind or gdb using the
106 following commands. Once in GDB, each MPI ranks will be represented as
107 a regular thread, and you can explore the state of each of them as
110 .. code-block:: shell
112 smpirun -wrapper valgrind ...other args...
113 smpirun -wrapper "gdb --args" --cfg=contexts/factory:thread ...other args...
117 ................................
118 Simulating Collective Operations
119 ................................
121 MPI collective operations are crucial to the performance of MPI
122 applications and must be carefully optimized according to many
123 parameters. Every existing implementation provides several algorithms
124 for each collective operation, and selects by default the best suited
125 one, depending on the sizes sent, the number of nodes, the
126 communicator, or the communication library being used. These
127 decisions are based on empirical results and theoretical complexity
128 estimation, and are very different between MPI implementations. In
129 most cases, the users can also manually tune the algorithm used for
130 each collective operation.
132 SMPI can simulate the behavior of several MPI implementations:
133 OpenMPI, MPICH, `STAR-MPI <http://star-mpi.sourceforge.net/>`_, and
134 MVAPICH2. For that, it provides 115 collective algorithms and several
135 selector algorithms, that were collected directly in the source code
136 of the targeted MPI implementations.
138 You can switch the automatic selector through the
139 ``smpi/coll-selector`` configuration item. Possible values:
141 - **ompi:** default selection logic of OpenMPI (version 3.1.2)
142 - **mpich**: default selection logic of MPICH (version 3.3b)
143 - **mvapich2**: selection logic of MVAPICH2 (version 1.9) tuned
144 on the Stampede cluster
145 - **impi**: preliminary version of an Intel MPI selector (version
146 4.1.3, also tuned for the Stampede cluster). Due the closed source
147 nature of Intel MPI, some of the algorithms described in the
148 documentation are not available, and are replaced by mvapich ones.
149 - **default**: legacy algorithms used in the earlier days of
150 SimGrid. Do not use for serious perform performance studies.
152 .. todo:: default should not even exist.
158 You can also pick the algorithm used for each collective with the
159 corresponding configuration item. For example, to use the pairwise
160 alltoall algorithm, one should add ``--cfg=smpi/alltoall:pair`` to the
161 line. This will override the selector (if any) for this algorithm. It
162 means that the selected algorithm will be used
164 .. Warning:: Some collective may require specific conditions to be
165 executed correctly (for instance having a communicator with a power
166 of two number of nodes only), which are currently not enforced by
167 Simgrid. Some crashes can be expected while trying these algorithms
168 with unusual sizes/parameters
173 Most of these are best described in `STAR-MPI's white paper <https://doi.org/10.1145/1183401.1183431>`_.
175 - default: naive one, by default
176 - ompi: use openmpi selector for the alltoall operations
177 - mpich: use mpich selector for the alltoall operations
178 - mvapich2: use mvapich2 selector for the alltoall operations
179 - impi: use intel mpi selector for the alltoall operations
180 - automatic (experimental): use an automatic self-benchmarking algorithm
181 - bruck: Described by Bruck et.al. in <a href="http://ieeexplore.ieee.org/xpl/articleDetails.jsp?arnumber=642949">this paper</a>
182 - 2dmesh: organizes the nodes as a two dimensional mesh, and perform allgather
184 - 3dmesh: adds a third dimension to the previous algorithm
185 - rdb: recursive doubling: extends the mesh to a nth dimension, each one
187 - pair: pairwise exchange, only works for power of 2 procs, size-1 steps,
188 each process sends and receives from the same process at each step
189 - pair_light_barrier: same, with small barriers between steps to avoid
191 - pair_mpi_barrier: same, with MPI_Barrier used
192 - pair_one_barrier: only one barrier at the beginning
193 - ring: size-1 steps, at each step a process send to process (n+i)%size, and receives from (n-i)%size
194 - ring_light_barrier: same, with small barriers between some phases to avoid contention
195 - ring_mpi_barrier: same, with MPI_Barrier used
196 - ring_one_barrier: only one barrier at the beginning
197 - basic_linear: posts all receives and all sends,
198 starts the communications, and waits for all communication to finish
199 - mvapich2_scatter_dest: isend/irecv with scattered destinations, posting only a few messages at the same time
203 - default: naive one, by default
204 - ompi: use openmpi selector for the alltoallv operations
205 - mpich: use mpich selector for the alltoallv operations
206 - mvapich2: use mvapich2 selector for the alltoallv operations
207 - impi: use intel mpi selector for the alltoallv operations
208 - automatic (experimental): use an automatic self-benchmarking algorithm
209 - bruck: same as alltoall
210 - pair: same as alltoall
211 - pair_light_barrier: same as alltoall
212 - pair_mpi_barrier: same as alltoall
213 - pair_one_barrier: same as alltoall
214 - ring: same as alltoall
215 - ring_light_barrier: same as alltoall
216 - ring_mpi_barrier: same as alltoall
217 - ring_one_barrier: same as alltoall
218 - ompi_basic_linear: same as alltoall
223 - default: naive one, by default
224 - ompi: use openmpi selector for the gather operations
225 - mpich: use mpich selector for the gather operations
226 - mvapich2: use mvapich2 selector for the gather operations
227 - impi: use intel mpi selector for the gather operations
228 - automatic (experimental): use an automatic self-benchmarking algorithm which will iterate over all implemented versions and output the best
229 - ompi_basic_linear: basic linear algorithm from openmpi, each process sends to the root
230 - ompi_binomial: binomial tree algorithm
231 - ompi_linear_sync: same as basic linear, but with a synchronization at the
232 beginning and message cut into two segments.
233 - mvapich2_two_level: SMP-aware version from MVAPICH. Gather first intra-node (defaults to mpich's gather), and then exchange with only one process/node. Use mvapich2 selector to change these to tuned algorithms for Stampede cluster.
238 - default: naive one, by default
239 - ompi: use openmpi selector for the barrier operations
240 - mpich: use mpich selector for the barrier operations
241 - mvapich2: use mvapich2 selector for the barrier operations
242 - impi: use intel mpi selector for the barrier operations
243 - automatic (experimental): use an automatic self-benchmarking algorithm
244 - ompi_basic_linear: all processes send to root
245 - ompi_two_procs: special case for two processes
246 - ompi_bruck: nsteps = sqrt(size), at each step, exchange data with rank-2^k and rank+2^k
247 - ompi_recursivedoubling: recursive doubling algorithm
248 - ompi_tree: recursive doubling type algorithm, with tree structure
249 - ompi_doublering: double ring algorithm
250 - mvapich2_pair: pairwise algorithm
251 - mpich_smp: barrier intra-node, then inter-node
256 - default: naive one, by default
257 - ompi: use openmpi selector for the scatter operations
258 - mpich: use mpich selector for the scatter operations
259 - mvapich2: use mvapich2 selector for the scatter operations
260 - impi: use intel mpi selector for the scatter operations
261 - automatic (experimental): use an automatic self-benchmarking algorithm
262 - ompi_basic_linear: basic linear scatter
263 - ompi_binomial: binomial tree scatter
264 - mvapich2_two_level_direct: SMP aware algorithm, with an intra-node stage (default set to mpich selector), and then a basic linear inter node stage. Use mvapich2 selector to change these to tuned algorithms for Stampede cluster.
265 - mvapich2_two_level_binomial: SMP aware algorithm, with an intra-node stage (default set to mpich selector), and then a binomial phase. Use mvapich2 selector to change these to tuned algorithms for Stampede cluster.
270 - default: naive one, by default
271 - ompi: use openmpi selector for the reduce operations
272 - mpich: use mpich selector for the reduce operations
273 - mvapich2: use mvapich2 selector for the reduce operations
274 - impi: use intel mpi selector for the reduce operations
275 - automatic (experimental): use an automatic self-benchmarking algorithm
276 - arrival_pattern_aware: root exchanges with the first process to arrive
277 - binomial: uses a binomial tree
278 - flat_tree: uses a flat tree
279 - NTSL: Non-topology-specific pipelined linear-bcast function
280 0->1, 1->2 ,2->3, ....., ->last node: in a pipeline fashion, with segments
282 - scatter_gather: scatter then gather
283 - ompi_chain: openmpi reduce algorithms are built on the same basis, but the
284 topology is generated differently for each flavor
285 chain = chain with spacing of size/2, and segment size of 64KB
286 - ompi_pipeline: same with pipeline (chain with spacing of 1), segment size
287 depends on the communicator size and the message size
288 - ompi_binary: same with binary tree, segment size of 32KB
289 - ompi_in_order_binary: same with binary tree, enforcing order on the
291 - ompi_binomial: same with binomial algo (redundant with default binomial
293 - ompi_basic_linear: basic algorithm, each process sends to root
294 - mvapich2_knomial: k-nomial algorithm. Default factor is 4 (mvapich2 selector adapts it through tuning)
295 - mvapich2_two_level: SMP-aware reduce, with default set to mpich both for intra and inter communicators. Use mvapich2 selector to change these to tuned algorithms for Stampede cluster.
296 - rab: `Rabenseifner <https://fs.hlrs.de/projects/par/mpi//myreduce.html>`_'s reduce algorithm
301 - default: naive one, by default
302 - ompi: use openmpi selector for the allreduce operations
303 - mpich: use mpich selector for the allreduce operations
304 - mvapich2: use mvapich2 selector for the allreduce operations
305 - impi: use intel mpi selector for the allreduce operations
306 - automatic (experimental): use an automatic self-benchmarking algorithm
307 - lr: logical ring reduce-scatter then logical ring allgather
308 - rab1: variations of the <a href="https://fs.hlrs.de/projects/par/mpi//myreduce.html">Rabenseifner</a> algorithm: reduce_scatter then allgather
309 - rab2: variations of the <a href="https://fs.hlrs.de/projects/par/mpi//myreduce.html">Rabenseifner</a> algorithm: alltoall then allgather
310 - rab_rsag: variation of the <a href="https://fs.hlrs.de/projects/par/mpi//myreduce.html">Rabenseifner</a> algorithm: recursive doubling
311 reduce_scatter then recursive doubling allgather
312 - rdb: recursive doubling
313 - smp_binomial: binomial tree with smp: binomial intra
314 SMP reduce, inter reduce, inter broadcast then intra broadcast
315 - smp_binomial_pipeline: same with segment size = 4096 bytes
316 - smp_rdb: intra: binomial allreduce, inter: Recursive
317 doubling allreduce, intra: binomial broadcast
318 - smp_rsag: intra: binomial allreduce, inter: reduce-scatter,
319 inter:allgather, intra: binomial broadcast
320 - smp_rsag_lr: intra: binomial allreduce, inter: logical ring
321 reduce-scatter, logical ring inter:allgather, intra: binomial broadcast
322 - smp_rsag_rab: intra: binomial allreduce, inter: rab
323 reduce-scatter, rab inter:allgather, intra: binomial broadcast
324 - redbcast: reduce then broadcast, using default or tuned algorithms if specified
325 - ompi_ring_segmented: ring algorithm used by OpenMPI
326 - mvapich2_rs: rdb for small messages, reduce-scatter then allgather else
327 - mvapich2_two_level: SMP-aware algorithm, with mpich as intra algoritm, and rdb as inter (Change this behavior by using mvapich2 selector to use tuned values)
328 - rab: default `Rabenseifner <https://fs.hlrs.de/projects/par/mpi//myreduce.html>`_ implementation
333 - default: naive one, by default
334 - ompi: use openmpi selector for the reduce_scatter operations
335 - mpich: use mpich selector for the reduce_scatter operations
336 - mvapich2: use mvapich2 selector for the reduce_scatter operations
337 - impi: use intel mpi selector for the reduce_scatter operations
338 - automatic (experimental): use an automatic self-benchmarking algorithm
339 - ompi_basic_recursivehalving: recursive halving version from OpenMPI
340 - ompi_ring: ring version from OpenMPI
341 - mpich_pair: pairwise exchange version from MPICH
342 - mpich_rdb: recursive doubling version from MPICH
343 - mpich_noncomm: only works for power of 2 procs, recursive doubling for noncommutative ops
349 - default: naive one, by default
350 - ompi: use openmpi selector for the allgather operations
351 - mpich: use mpich selector for the allgather operations
352 - mvapich2: use mvapich2 selector for the allgather operations
353 - impi: use intel mpi selector for the allgather operations
354 - automatic (experimental): use an automatic self-benchmarking algorithm
355 - 2dmesh: see alltoall
356 - 3dmesh: see alltoall
357 - bruck: Described by Bruck et.al. in <a href="http://ieeexplore.ieee.org/xpl/articleDetails.jsp?arnumber=642949">
358 Efficient algorithms for all-to-all communications in multiport message-passing systems</a>
359 - GB: Gather - Broadcast (uses tuned version if specified)
360 - loosely_lr: Logical Ring with grouping by core (hardcoded, default
362 - NTSLR: Non Topology Specific Logical Ring
363 - NTSLR_NB: Non Topology Specific Logical Ring, Non Blocking operations
366 - rhv: only power of 2 number of processes
368 - SMP_NTS: gather to root of each SMP, then every root of each SMP node
369 post INTER-SMP Sendrecv, then do INTRA-SMP Bcast for each receiving message,
370 using logical ring algorithm (hardcoded, default processes/SMP: 8)
371 - smp_simple: gather to root of each SMP, then every root of each SMP node
372 post INTER-SMP Sendrecv, then do INTRA-SMP Bcast for each receiving message,
373 using simple algorithm (hardcoded, default processes/SMP: 8)
374 - spreading_simple: from node i, order of communications is i -> i + 1, i ->
375 i + 2, ..., i -> (i + p -1) % P
376 - ompi_neighborexchange: Neighbor Exchange algorithm for allgather.
377 Described by Chen et.al. in `Performance Evaluation of Allgather
378 Algorithms on Terascale Linux Cluster with Fast Ethernet <http://ieeexplore.ieee.org/xpl/articleDetails.jsp?tp=&arnumber=1592302>`_
379 - mvapich2_smp: SMP aware algorithm, performing intra-node gather, inter-node allgather with one process/node, and bcast intra-node
384 - default: naive one, by default
385 - ompi: use openmpi selector for the allgatherv operations
386 - mpich: use mpich selector for the allgatherv operations
387 - mvapich2: use mvapich2 selector for the allgatherv operations
388 - impi: use intel mpi selector for the allgatherv operations
389 - automatic (experimental): use an automatic self-benchmarking algorithm
390 - GB: Gatherv - Broadcast (uses tuned version if specified, but only for Bcast, gatherv is not tuned)
393 - ompi_neighborexchange: see allgather
394 - ompi_bruck: see allgather
395 - mpich_rdb: recursive doubling algorithm from MPICH
396 - mpich_ring: ring algorithm from MPICh - performs differently from the one from STAR-MPI
401 - default: naive one, by default
402 - ompi: use openmpi selector for the bcast operations
403 - mpich: use mpich selector for the bcast operations
404 - mvapich2: use mvapich2 selector for the bcast operations
405 - impi: use intel mpi selector for the bcast operations
406 - automatic (experimental): use an automatic self-benchmarking algorithm
407 - arrival_pattern_aware: root exchanges with the first process to arrive
408 - arrival_pattern_aware_wait: same with slight variation
409 - binomial_tree: binomial tree exchange
410 - flattree: flat tree exchange
411 - flattree_pipeline: flat tree exchange, message split into 8192 bytes pieces
412 - NTSB: Non-topology-specific pipelined binary tree with 8192 bytes pieces
413 - NTSL: Non-topology-specific pipelined linear with 8192 bytes pieces
414 - NTSL_Isend: Non-topology-specific pipelined linear with 8192 bytes pieces, asynchronous communications
415 - scatter_LR_allgather: scatter followed by logical ring allgather
416 - scatter_rdb_allgather: scatter followed by recursive doubling allgather
417 - arrival_scatter: arrival pattern aware scatter-allgather
418 - SMP_binary: binary tree algorithm with 8 cores/SMP
419 - SMP_binomial: binomial tree algorithm with 8 cores/SMP
420 - SMP_linear: linear algorithm with 8 cores/SMP
421 - ompi_split_bintree: binary tree algorithm from OpenMPI, with message split in 8192 bytes pieces
422 - ompi_pipeline: pipeline algorithm from OpenMPI, with message split in 128KB pieces
423 - mvapich2_inter_node: Inter node default mvapich worker
424 - mvapich2_intra_node: Intra node default mvapich worker
425 - mvapich2_knomial_intra_node: k-nomial intra node default mvapich worker. default factor is 4.
430 .. warning:: This is still very experimental.
432 An automatic version is available for each collective (or even as a selector). This specific
433 version will loop over all other implemented algorithm for this particular collective, and apply
434 them while benchmarking the time taken for each process. It will then output the quickest for
435 each process, and the global quickest. This is still unstable, and a few algorithms which need
436 specific number of nodes may crash.
441 To add a new algorithm, one should check in the src/smpi/colls folder
442 how other algorithms are coded. Using plain MPI code inside Simgrid
443 can't be done, so algorithms have to be changed to use smpi version of
444 the calls instead (MPI_Send will become smpi_mpi_send). Some functions
445 may have different signatures than their MPI counterpart, please check
446 the other algorithms or contact us using the `>SimGrid
447 developers mailing list <http://lists.gforge.inria.fr/mailman/listinfo/simgrid-devel>`_.
449 Example: adding a "pair" version of the Alltoall collective.
451 - Implement it in a file called alltoall-pair.c in the src/smpi/colls folder. This file should include colls_private.hpp.
453 - The name of the new algorithm function should be smpi_coll_tuned_alltoall_pair, with the same signature as MPI_Alltoall.
455 - Once the adaptation to SMPI code is done, add a reference to the file ("src/smpi/colls/alltoall-pair.c") in the SMPI_SRC part of the DefinePackages.cmake file inside buildtools/cmake, to allow the file to be built and distributed.
457 - To register the new version of the algorithm, simply add a line to the corresponding macro in src/smpi/colls/cools.h ( add a "COLL_APPLY(action, COLL_ALLTOALL_SIG, pair)" to the COLL_ALLTOALLS macro ). The algorithm should now be compiled and be selected when using --cfg=smpi/alltoall:pair at runtime.
459 - To add a test for the algorithm inside Simgrid's test suite, juste add the new algorithm name in the ALLTOALL_COLL list found inside teshsuite/smpi/CMakeLists.txt . When running ctest, a test for the new algorithm should be generated and executed. If it does not pass, please check your code or contact us.
461 - Please submit your patch for inclusion in SMPI, for example through a pull request on GitHub or directly per email.
464 Tracing of Internal Communications
465 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
467 By default, the collective operations are traced as a unique operation
468 because tracing all point-to-point communications composing them could
469 result in overloaded, hard to interpret traces. If you want to debug
470 and compare collective algorithms, you should set the
471 ``tracing/smpi/internals`` configuration item to 1 instead of 0.
473 Here are examples of two alltoall collective algorithms runs on 16 nodes,
474 the first one with a ring algorithm, the second with a pairwise one.
476 .. image:: /img/smpi_simgrid_alltoall_ring_16.png
479 Alltoall on 16 Nodes with the Ring Algorithm.
481 .. image:: /img/smpi_simgrid_alltoall_pair_16.png
484 Alltoall on 16 Nodes with the Pairwise Algorithm.
486 -------------------------
487 What can run within SMPI?
488 -------------------------
490 You can run unmodified MPI applications (both C/C++ and Fortran) within
491 SMPI, provided that you only use MPI calls that we implemented. Global
492 variables should be handled correctly on Linux systems.
498 Our coverage of the interface is very decent, but still incomplete;
499 Given the size of the MPI standard, we may well never manage to
500 implement absolutely all existing primitives. Currently, we have
501 almost no support for I/O primitives, but we still pass a very large
502 amount of the MPICH coverage tests.
504 The full list of not yet implemented functions is documented in the
505 file `include/smpi/smpi.h
506 <https://framagit.org/simgrid/simgrid/tree/master/include/smpi/smpi.h>`_
507 in your version of SimGrid, between two lines containing the ``FIXME``
508 marker. If you really miss a feature, please get in touch with us: we
509 can guide you though the SimGrid code to help you implementing it, and
510 we'd be glad to integrate your contribution to the main project.
512 .. _SMPI_what_globals:
514 .................................
515 Privatization of global variables
516 .................................
518 Concerning the globals, the problem comes from the fact that usually,
519 MPI processes run as real UNIX processes while they are all folded
520 into threads of a unique system process in SMPI. Global variables are
521 usually private to each MPI process while they become shared between
522 the processes in SMPI. The problem and some potential solutions are
523 discussed in this article: `Automatic Handling of Global Variables for
524 Multi-threaded MPI Programs
525 <http://charm.cs.illinois.edu/newPapers/11-23/paper.pdf>` (note that
526 this article does not deal with SMPI but with a competing solution
527 called AMPI that suffers of the same issue). This point used to be
528 problematic in SimGrid, but the problem should now be handled
529 automatically on Linux.
531 Older versions of SimGrid came with a script that automatically
532 privatized the globals through static analysis of the source code. But
533 our implementation was not robust enough to be used in production, so
534 it was removed at some point. Currently, SMPI comes with two
535 privatization mechanisms that you can :ref:`select at runtime
536 <cfg=smpi/privatization>`. The dlopen approach is used by
537 default as it is much faster and still very robust. The mmap approach
538 is an older approach that proves to be slower.
540 With the **mmap approach**, SMPI duplicates and dynamically switch the
541 ``.data`` and ``.bss`` segments of the ELF process when switching the
542 MPI ranks. This allows each ranks to have its own copy of the global
543 variables. No copy actually occures as this mechanism uses ``mmap()``
544 for efficiency. This mechanism is considered to be very robust on all
545 systems supporting ``mmap()`` (Linux and most BSDs). Its performance
546 is questionable since each context switch between MPI ranks induces
547 several syscalls to change the ``mmap`` that redirects the ``.data``
548 and ``.bss`` segments to the copies of the new rank. The code will
549 also be copied several times in memory, inducing a slight increase of
552 Another limitation is that SMPI only accounts for global variables
553 defined in the executable. If the processes use external global
554 variables from dynamic libraries, they won't be switched
555 correctly. The easiest way to solve this is to statically link against
556 the library with these globals. This way, each MPI rank will get its
557 own copy of these libraries. Of course you should never statically
558 link against the SimGrid library itself.
560 With the **dlopen approach**, SMPI loads several copies of the same
561 executable in memory as if it were a library, so that the global
562 variables get naturally dupplicated. It first requires the executable
563 to be compiled as a relocatable binary, which is less common for
564 programs than for libraries. But most distributions are now compiled
565 this way for security reason as it allows to randomize the address
566 space layout. It should thus be safe to compile most (any?) program
567 this way. The second trick is that the dynamic linker refuses to link
568 the exact same file several times, be it a library or a relocatable
569 executable. It makes perfectly sense in the general case, but we need
570 to circumvent this rule of thumb in our case. To that extend, the
571 binary is copied in a temporary file before being re-linked against.
572 ``dlmopen()`` cannot be used as it only allows 256 contextes, and as it
573 would also dupplicate simgrid itself.
575 This approach greatly speeds up the context switching, down to about
576 40 CPU cycles with our raw contextes, instead of requesting several
577 syscalls with the ``mmap()`` approach. Another advantage is that it
578 permits to run the SMPI contexts in parallel, which is obviously not
579 possible with the ``mmap()`` approach. It was tricky to implement, but
580 we are not aware of any flaws, so smpirun activates it by default.
582 In the future, it may be possible to further reduce the memory and
583 disk consumption. It seems that we could `punch holes
584 <https://lwn.net/Articles/415889/>`_ in the files before dl-loading
585 them to remove the code and constants, and mmap these area onto a
586 unique copy. If done correctly, this would reduce the disk- and
587 memory- usage to the bare minimum, and would also reduce the pressure
588 on the CPU instruction cache. See the `relevant bug
589 <https://github.com/simgrid/simgrid/issues/137>`_ on github for
590 implementation leads.\n
592 Also, currently, only the binary is copied and dlopen-ed for each MPI
593 rank. We could probably extend this to external dependencies, but for
594 now, any external dependencies must be statically linked into your
595 application. As usual, simgrid itself shall never be statically linked
596 in your app. You don't want to give a copy of SimGrid to each MPI rank:
597 that's ways too much for them to deal with.
599 .. todo: speak of smpi/privatize-libs here
601 ----------------------------------------------
602 Adapting your MPI code for further scalability
603 ----------------------------------------------
605 As detailed in the `reference article
606 <http://hal.inria.fr/hal-01415484>`_, you may want to adapt your code
607 to improve the simulation performance. But these tricks may seriously
608 hinder the result quality (or even prevent the app to run) if used
609 wrongly. We assume that if you want to simulate an HPC application,
610 you know what you are doing. Don't prove us wrong!
612 ..............................
613 Reducing your memory footprint
614 ..............................
616 If you get short on memory (the whole app is executed on a single node when
617 simulated), you should have a look at the SMPI_SHARED_MALLOC and
618 SMPI_SHARED_FREE macros. It allows to share memory areas between processes: The
619 purpose of these macro is that the same line malloc on each process will point
620 to the exact same memory area. So if you have a malloc of 2M and you have 16
621 processes, this macro will change your memory consumption from 2M*16 to 2M
622 only. Only one block for all processes.
624 If your program is ok with a block containing garbage value because all
625 processes write and read to the same place without any kind of coordination,
626 then this macro can dramatically shrink your memory consumption. For example,
627 that will be very beneficial to a matrix multiplication code, as all blocks will
628 be stored on the same area. Of course, the resulting computations will useless,
629 but you can still study the application behavior this way.
631 Naturally, this won't work if your code is data-dependent. For example, a Jacobi
632 iterative computation depends on the result computed by the code to detect
633 convergence conditions, so turning them into garbage by sharing the same memory
634 area between processes does not seem very wise. You cannot use the
635 SMPI_SHARED_MALLOC macro in this case, sorry.
637 This feature is demoed by the example file
638 `examples/smpi/NAS/dt.c <https://framagit.org/simgrid/simgrid/tree/master/examples/smpi/NAS/dt.c>`_
640 .........................
641 Toward Faster Simulations
642 .........................
644 If your application is too slow, try using SMPI_SAMPLE_LOCAL,
645 SMPI_SAMPLE_GLOBAL and friends to indicate which computation loops can
646 be sampled. Some of the loop iterations will be executed to measure
647 their duration, and this duration will be used for the subsequent
648 iterations. These samples are done per processor with
649 SMPI_SAMPLE_LOCAL, and shared between all processors with
650 SMPI_SAMPLE_GLOBAL. Of course, none of this will work if the execution
651 time of your loop iteration are not stable.
653 This feature is demoed by the example file
654 `examples/smpi/NAS/ep.c <https://framagit.org/simgrid/simgrid/tree/master/examples/smpi/NAS/ep.c>`_
656 .............................
657 Ensuring Accurate Simulations
658 .............................
660 Out of the box, SimGrid may give you fairly accurate results, but
661 there is a plenty of factors that could go wrong and make your results
662 inaccurate or even plainly wrong. Actually, you can only get accurate
663 results of a nicely built model, including both the system hardware
664 and your application. Such models are hard to pass over and reuse in
665 other settings, because elements that are not relevant to an
666 application (say, the latency of point-to-point communications,
667 collective operation implementation details or CPU-network
668 interaction) may be irrelevant to another application. The dream of
669 the perfect model, encompassing every aspects is only a chimera, as
670 the only perfect model of the reality is the reality. If you go for
671 simulation, then you have to ignore some irrelevant aspects of the
672 reality, but which aspects are irrelevant is actually
673 application-dependent...
675 The only way to assess whether your settings provide accurate results
676 is to double-check these results. If possible, you should first run
677 the same experiment in simulation and in real life, gathering as much
678 information as you can. Try to understand the discrepancies in the
679 results that you observe between both settings (visualization can be
680 precious for that). Then, try to modify your model (of the platform,
681 of the collective operations) to reduce the most preeminent differences.
683 If the discrepancies come from the computing time, try adapting the
684 ``smpi/host-speed``: reduce it if your simulation runs faster than in
685 reality. If the error come from the communication, then you need to
686 fiddle with your platform file.
688 Be inventive in your modeling. Don't be afraid if the names given by
689 SimGrid does not match the real names: we got very good results by
690 modeling multicore/GPU machines with a set of separate hosts
691 interconnected with very fast networks (but don't trust your model
692 because it has the right names in the right place either).
694 Finally, you may want to check `this article
695 <https://hal.inria.fr/hal-00907887>`_ on the classical pitfalls in
696 modeling distributed systems.
698 -------------------------
699 Troubleshooting with SMPI
700 -------------------------
702 .................................
703 ./configure refuses to use smpicc
704 .................................
706 If your ``./configure`` reports that the compiler is not
707 functional or that you are cross-compiling, try to define the
708 ``SMPI_PRETEND_CC`` environment variable before running the
711 .. code-block:: shell
713 SMPI_PRETEND_CC=1 ./configure # here come the configure parameters
716 Indeed, the programs compiled with ``smpicc`` cannot be executed
717 without ``smpirun`` (they are shared libraries and do weird things on
718 startup), while configure wants to test them directly. With
719 ``SMPI_PRETEND_CC`` smpicc does not compile as shared, and the SMPI
720 initialization stops and returns 0 before doing anything that would
721 fail without ``smpirun``.
725 Make sure that SMPI_PRETEND_CC is only set when calling ./configure,
726 not during the actual execution, or any program compiled with smpicc
727 will stop before starting.
729 ..............................................
730 ./configure does not pick smpicc as a compiler
731 ..............................................
733 In addition to the previous answers, some projects also need to be
734 explicitely told what compiler to use, as follows:
736 .. code-block:: shell
738 SMPI_PRETEND_CC=1 ./configure CC=smpicc # here come the other configure parameters
741 Maybe your configure is using another variable, such as ``cc`` (in
742 lower case) or similar. Just check the logs.
744 .....................................
745 error: unknown type name 'useconds_t'
746 .....................................
748 Try to add ``-D_GNU_SOURCE`` to your compilation line to get ride
751 The reason is that SMPI provides its own version of ``usleep(3)``
752 to override it and to block in the simulation world, not in the real
753 one. It needs the ``useconds_t`` type for that, which is declared
754 only if you declare ``_GNU_SOURCE`` before including
755 ``unistd.h``. If your project includes that header file before
756 SMPI, then you need to ensure that you pass the right configuration
757 defines as advised above.
763 -----------------------------
764 Trace Replay and Offline SMPI
765 -----------------------------
767 Although SMPI is often used for :ref:`online simulation
768 <SMPI_online>`, where the application is executed for real, you can
769 also go for offline simulation through trace replay.
771 SimGrid uses time-independent traces, in which each actor is given a
772 script of the actions to do sequentially. These trace files can
773 actually be captured with the online version of SMPI, as follows:
775 .. code-block:: shell
777 $ smpirun -trace-ti --cfg=tracing/filename:LU.A.32 -np 32 -platform ../cluster_backbone.xml bin/lu.A.32
779 The produced trace is composed of a file ``LU.A.32`` and a folder
780 ``LU.A.32_files``. The file names don't match with the MPI ranks, but
783 To replay this with SMPI, you need to first compile the provided
784 ``smpi_replay.cpp`` file, that comes from
785 `simgrid/examples/smpi/replay
786 <https://framagit.org/simgrid/simgrid/tree/master/examples/smpi/replay>`_.
788 .. code-block:: shell
790 $ smpicxx ../replay.cpp -O3 -o ../smpi_replay
792 Afterward, you can replay your trace in SMPI as follows:
794 $ smpirun -np 32 -platform ../cluster_torus.xml -ext smpi_replay ../smpi_replay LU.A.32
796 All the outputs are gone, as the application is not really simulated
797 here. Its trace is simply replayed. But if you visualize the live
798 simulation and the replay, you will see that the behavior is
799 unchanged. The simulation does not run much faster on this very
800 example, but this becomes very interesting when your application
801 is computationally hungry.