Logo AND Algorithmique Numérique Distribuée

Public GIT Repository
d995984fad8223a6bd94d7532b2721294cdf3a4a
[simgrid.git] / examples / smpi / trace / trace.c
1 /* Copyright (c) 2009-2010, 2012-2015. The SimGrid Team.
2  * All rights reserved.                                                     */
3
4 /* This program is free software; you can redistribute it and/or modify it
5  * under the terms of the license (GNU LGPL) which comes with this package. */
6
7 #include <stdio.h>
8 #include "mpi.h"
9 #include "simgrid/instr.h"
10
11 #define DATATOSENT 100000
12
13 int main(int argc, char *argv[])
14 {
15   int rank;
16   int numprocs;
17
18   MPI_Init(&argc, &argv);
19   MPI_Comm_size(MPI_COMM_WORLD, &numprocs);
20   MPI_Comm_rank(MPI_COMM_WORLD, &rank);
21
22   /** Tests:
23   * A: 0(isend -> wait) with 1(recv)
24   * B: 0(send) with 1(irecv -> wait)
25   * C: 0(N * isend -> N * wait) with 1(N * recv)
26   * D: 0(N * isend -> N * waitany) with 1(N * recv)
27   * E: 0(N*send) with 1(N*irecv, N*wait)
28   * F: 0(N*send) with 1(N*irecv, N*waitany)
29   * G: 0(N* isend -> waitall) with 1(N*recv)
30   * H: 0(N*send) with 1(N*irecv, waitall)
31   * I: 0(2*N*send, 2*N*Irecv, Waitall) with
32   *    1(N*irecv, waitall, N*isend, N*waitany) with
33   *    2(N*irecv, N*waitany, N*isend, waitall)
34   * J: 0(N*isend, N*test, N*wait) with (N*irecv, N*test, N*wait)
35   s*/
36   int N = 13;
37
38   int tag = 12345;
39 /////////////////////////////////////////
40 ////////////////// RANK 0
41 ///////////////////////////////////
42   if (rank == 0) {
43     MPI_Request request;
44     MPI_Status status;
45     MPI_Request req[2 * N];
46     MPI_Status sta[2 * N];
47     int *r = (int *) malloc(sizeof(int) * DATATOSENT);
48
49     /** Test A */
50     TRACE_smpi_set_category("A");
51     MPI_Isend(r, DATATOSENT, MPI_INT, 1, tag, MPI_COMM_WORLD, &request);
52     MPI_Wait(&request, &status);
53     MPI_Barrier(MPI_COMM_WORLD);
54
55     /** Test B */
56     TRACE_smpi_set_category("B");
57     MPI_Send(r, DATATOSENT, MPI_INT, 1, tag, MPI_COMM_WORLD);
58     MPI_Barrier(MPI_COMM_WORLD);
59
60     /** Test C */
61     TRACE_smpi_set_category("C");
62     for (int i = 0; i < N; i++) {
63       MPI_Isend(r, DATATOSENT, MPI_INT, 1, tag, MPI_COMM_WORLD, &req[i]);
64     }
65     for (int i = 0; i < N; i++) {
66       MPI_Wait(&req[i], &sta[i]);
67     }
68     MPI_Barrier(MPI_COMM_WORLD);
69
70     TRACE_smpi_set_category("D");
71     for (int i = 0; i < N; i++) {
72       MPI_Isend(r, DATATOSENT, MPI_INT, 1, tag, MPI_COMM_WORLD, &req[i]);
73     }
74     for (int i = 0; i < N; i++) {
75       int completed;
76       MPI_Waitany(N, req, &completed, sta);
77     }
78     MPI_Barrier(MPI_COMM_WORLD);
79
80     TRACE_smpi_set_category("E");
81     for (int i = 0; i < N; i++) {
82       MPI_Send(r, DATATOSENT, MPI_INT, 1, tag, MPI_COMM_WORLD);
83     }
84     MPI_Barrier(MPI_COMM_WORLD);
85
86     TRACE_smpi_set_category("F");
87     for (int i = 0; i < N; i++) {
88       MPI_Send(r, DATATOSENT, MPI_INT, 1, tag, MPI_COMM_WORLD);
89     }
90     MPI_Barrier(MPI_COMM_WORLD);
91
92     TRACE_smpi_set_category("G");
93     for (int i = 0; i < N; i++) {
94       MPI_Isend(r, DATATOSENT, MPI_INT, 1, tag, MPI_COMM_WORLD, &req[i]);
95     }
96     MPI_Waitall(N, req, sta);
97     MPI_Barrier(MPI_COMM_WORLD);
98
99     TRACE_smpi_set_category("H");
100     for (int i = 0; i < N; i++) {
101       MPI_Send(r, DATATOSENT, MPI_INT, 1, tag, MPI_COMM_WORLD);
102     }
103     MPI_Barrier(MPI_COMM_WORLD);
104
105     TRACE_smpi_set_category("I");
106     for (int i = 0; i < 2 * N; i++) {
107       if (i < N) {
108         MPI_Send(r, DATATOSENT, MPI_INT, 2, tag, MPI_COMM_WORLD);
109       } else {
110         MPI_Send(r, DATATOSENT, MPI_INT, 1, tag, MPI_COMM_WORLD);
111       }
112     }
113     MPI_Barrier(MPI_COMM_WORLD);
114     for (int i = 0; i < 2 * N; i++) {
115       if (i < N) {
116         MPI_Irecv(r, DATATOSENT, MPI_INT, 1, tag, MPI_COMM_WORLD, &req[i]);
117       } else {
118         MPI_Irecv(r, DATATOSENT, MPI_INT, 2, tag, MPI_COMM_WORLD, &req[i]);
119       }
120     }
121     MPI_Waitall(2 * N, req, sta);
122     MPI_Barrier(MPI_COMM_WORLD);
123
124     TRACE_smpi_set_category("J");
125     for (int i = 0; i < N; i++) {
126       MPI_Isend(r, DATATOSENT, MPI_INT, 1, tag, MPI_COMM_WORLD, &req[i]);
127     }
128     for (int i = 0; i < N; i++) {
129       int flag;
130       MPI_Test(&req[i], &flag, &sta[i]);
131     }
132     for (int i = 0; i < N; i++) {
133       MPI_Wait(&req[i], &sta[i]);
134     }
135     free(r);
136 /////////////////////////////////////////
137 ////////////////// RANK 1
138 ///////////////////////////////////
139   } else if (rank == 1) {
140     MPI_Request request;
141     MPI_Status status;
142     MPI_Request req[N];
143     MPI_Status sta[N];
144     int *r = (int *) malloc(sizeof(int) * DATATOSENT);
145
146     TRACE_smpi_set_category("A");
147     MPI_Recv(r, DATATOSENT, MPI_INT, 0, tag, MPI_COMM_WORLD, &status);
148     MPI_Barrier(MPI_COMM_WORLD);
149
150     TRACE_smpi_set_category("B");
151     MPI_Irecv(r, DATATOSENT, MPI_INT, 0, tag, MPI_COMM_WORLD, &request);
152     MPI_Wait(&request, &status);
153     MPI_Barrier(MPI_COMM_WORLD);
154
155     TRACE_smpi_set_category("C");
156     for (int i = 0; i < N; i++) {
157       MPI_Recv(r, DATATOSENT, MPI_INT, 0, tag, MPI_COMM_WORLD, &sta[i]);
158     }
159     MPI_Barrier(MPI_COMM_WORLD);
160
161     TRACE_smpi_set_category("D");
162     for (int i = 0; i < N; i++) {
163       MPI_Recv(r, DATATOSENT, MPI_INT, 0, tag, MPI_COMM_WORLD, &sta[i]);
164     }
165     MPI_Barrier(MPI_COMM_WORLD);
166
167     TRACE_smpi_set_category("E");
168     for (int i = 0; i < N; i++) {
169       MPI_Irecv(r, DATATOSENT, MPI_INT, 0, tag, MPI_COMM_WORLD, &req[i]);
170     }
171     for (int i = 0; i < N; i++) {
172       MPI_Wait(&req[i], &sta[i]);
173     }
174     MPI_Barrier(MPI_COMM_WORLD);
175
176     TRACE_smpi_set_category("F");
177     for (int i = 0; i < N; i++) {
178       MPI_Irecv(r, DATATOSENT, MPI_INT, 0, tag, MPI_COMM_WORLD, &req[i]);
179     }
180     for (int i = 0; i < N; i++) {
181       int completed;
182       MPI_Waitany(N, req, &completed, sta);
183     }
184     MPI_Barrier(MPI_COMM_WORLD);
185
186     TRACE_smpi_set_category("G");
187     for (int i = 0; i < N; i++) {
188       MPI_Recv(r, DATATOSENT, MPI_INT, 0, tag, MPI_COMM_WORLD, &sta[i]);
189     }
190     MPI_Barrier(MPI_COMM_WORLD);
191
192     TRACE_smpi_set_category("H");
193     for (int i = 0; i < N; i++) {
194       MPI_Irecv(r, DATATOSENT, MPI_INT, 0, tag, MPI_COMM_WORLD, &req[i]);
195     }
196     MPI_Waitall(N, req, sta);
197     MPI_Barrier(MPI_COMM_WORLD);
198
199     TRACE_smpi_set_category("I");
200     for (int i = 0; i < N; i++) {
201       MPI_Irecv(r, DATATOSENT, MPI_INT, 0, tag, MPI_COMM_WORLD, &req[i]);
202     }
203     MPI_Waitall(N, req, sta);
204
205     MPI_Barrier(MPI_COMM_WORLD);
206     for (int i = 0; i < N; i++) {
207       MPI_Isend(r, DATATOSENT, MPI_INT, 0, tag, MPI_COMM_WORLD, &req[i]);
208     }
209     MPI_Waitall(N, req, sta);
210     MPI_Barrier(MPI_COMM_WORLD);
211
212     TRACE_smpi_set_category("J");
213     for (int i = 0; i < N; i++) {
214       MPI_Irecv(r, DATATOSENT, MPI_INT, 0, tag, MPI_COMM_WORLD, &req[i]);
215     }
216     for (int i = 0; i < N; i++) {
217       int flag;
218       MPI_Test(&req[i], &flag, &sta[i]);
219     }
220     for (int i = 0; i < N; i++) {
221       MPI_Wait(&req[i], &sta[i]);
222     }
223     free(r);
224 /////////////////////////////////////////
225 ////////////////// RANK 2
226 ///////////////////////////////////
227   } else if (rank == 2) {
228     MPI_Request req[N];
229     MPI_Status sta[N];
230     int *r = (int *) malloc(sizeof(int) * DATATOSENT);
231
232     MPI_Barrier(MPI_COMM_WORLD);
233     MPI_Barrier(MPI_COMM_WORLD);
234     MPI_Barrier(MPI_COMM_WORLD);
235     MPI_Barrier(MPI_COMM_WORLD);
236     MPI_Barrier(MPI_COMM_WORLD);
237     MPI_Barrier(MPI_COMM_WORLD);
238     MPI_Barrier(MPI_COMM_WORLD);
239     MPI_Barrier(MPI_COMM_WORLD);
240     TRACE_smpi_set_category("I");
241     for (int i = 0; i < N; i++) {
242       MPI_Irecv(r, DATATOSENT, MPI_INT, 0, tag, MPI_COMM_WORLD, &req[i]);
243     }
244     for (int i = 0; i < N; i++) {
245       int completed;
246       MPI_Waitany(N, req, &completed, sta);
247     }
248     MPI_Barrier(MPI_COMM_WORLD);
249
250     for (int i = 0; i < N; i++) {
251       MPI_Send(r, DATATOSENT, MPI_INT, 0, tag, MPI_COMM_WORLD);
252     }
253     MPI_Barrier(MPI_COMM_WORLD);
254     free(r);
255   }
256   MPI_Barrier(MPI_COMM_WORLD);
257   MPI_Finalize();
258   return 0;
259 }