Logo AND Algorithmique Numérique Distribuée

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