This source file includes following definitions.
- main
1
2
3
4
5
6 #include "mpi.h"
7 #include <stdio.h>
8 #include <string.h>
9 #include <stdlib.h>
10
11
12
13
14
15
16
17
18
19
20
21 int main(int argc, char **argv)
22 {
23 MPI_Datatype newtype;
24 int i, ndims, array_of_gsizes[3], array_of_distribs[3];
25 int order, nprocs, len, *buf, mynod;
26 MPI_Count bufcount;
27 int array_of_dargs[3], array_of_psizes[3];
28 MPI_File fh;
29 MPI_Status status;
30 double stim, write_tim, new_write_tim, write_bw;
31 double read_tim, new_read_tim, read_bw;
32 char *filename;
33
34 MPI_Init(&argc,&argv);
35 MPI_Comm_rank(MPI_COMM_WORLD, &mynod);
36 MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
37
38
39
40 if (!mynod) {
41 i = 1;
42 while ((i < argc) && strcmp("-fname", *argv)) {
43 i++;
44 argv++;
45 }
46 if (i >= argc) {
47 fprintf(stderr, "\n*# Usage: coll_perf -fname filename\n\n");
48 MPI_Abort(MPI_COMM_WORLD, 1);
49 }
50 argv++;
51 len = strlen(*argv);
52 filename = (char *) malloc(len+1);
53 strcpy(filename, *argv);
54 MPI_Bcast(&len, 1, MPI_INT, 0, MPI_COMM_WORLD);
55 MPI_Bcast(filename, len+1, MPI_CHAR, 0, MPI_COMM_WORLD);
56 }
57 else {
58 MPI_Bcast(&len, 1, MPI_INT, 0, MPI_COMM_WORLD);
59 filename = (char *) malloc(len+1);
60 MPI_Bcast(filename, len+1, MPI_CHAR, 0, MPI_COMM_WORLD);
61 }
62
63
64 ndims = 3;
65 order = MPI_ORDER_C;
66
67 array_of_gsizes[0] = 128*17;
68 array_of_gsizes[1] = 128*9;
69 array_of_gsizes[2] = 128*11;
70
71 array_of_distribs[0] = MPI_DISTRIBUTE_BLOCK;
72 array_of_distribs[1] = MPI_DISTRIBUTE_BLOCK;
73 array_of_distribs[2] = MPI_DISTRIBUTE_BLOCK;
74
75 array_of_dargs[0] = MPI_DISTRIBUTE_DFLT_DARG;
76 array_of_dargs[1] = MPI_DISTRIBUTE_DFLT_DARG;
77 array_of_dargs[2] = MPI_DISTRIBUTE_DFLT_DARG;
78
79 for (i=0; i<ndims; i++) array_of_psizes[i] = 0;
80 MPI_Dims_create(nprocs, ndims, array_of_psizes);
81
82 MPI_Type_create_darray(nprocs, mynod, ndims, array_of_gsizes,
83 array_of_distribs, array_of_dargs,
84 array_of_psizes, order, MPI_INT, &newtype);
85 MPI_Type_commit(&newtype);
86
87 MPI_Type_size_x(newtype, &bufcount);
88 bufcount = bufcount/sizeof(int);
89 buf = (int *) malloc(bufcount * sizeof(int));
90
91
92
93
94 MPI_File_open(MPI_COMM_WORLD, filename, MPI_MODE_CREATE | MPI_MODE_RDWR,
95 MPI_INFO_NULL, &fh);
96 MPI_File_set_view(fh, 0, MPI_INT, newtype, "native", MPI_INFO_NULL);
97 MPI_File_write_all(fh, buf, bufcount, MPI_INT, &status);
98 MPI_File_seek(fh, 0, MPI_SEEK_SET);
99 MPI_File_read_all(fh, buf, bufcount, MPI_INT, &status);
100 MPI_File_close(&fh);
101
102 MPI_Barrier(MPI_COMM_WORLD);
103
104
105 MPI_File_open(MPI_COMM_WORLD, filename, MPI_MODE_CREATE | MPI_MODE_RDWR,
106 MPI_INFO_NULL, &fh);
107 MPI_File_set_view(fh, 0, MPI_INT, newtype, "native", MPI_INFO_NULL);
108
109 MPI_Barrier(MPI_COMM_WORLD);
110 stim = MPI_Wtime();
111 MPI_File_write_all(fh, buf, bufcount, MPI_INT, &status);
112 write_tim = MPI_Wtime() - stim;
113 MPI_File_close(&fh);
114
115 MPI_Allreduce(&write_tim, &new_write_tim, 1, MPI_DOUBLE, MPI_MAX,
116 MPI_COMM_WORLD);
117
118 if (mynod == 0) {
119 write_bw = (array_of_gsizes[0]*array_of_gsizes[1]*array_of_gsizes[2]*sizeof(int))/(new_write_tim*1024.0*1024.0);
120 fprintf(stderr, "Global array size %d x %d x %d integers\n", array_of_gsizes[0], array_of_gsizes[1], array_of_gsizes[2]);
121 fprintf(stderr, "Collective write time = %f sec, Collective write bandwidth = %f Mbytes/sec\n", new_write_tim, write_bw);
122 }
123
124 MPI_Barrier(MPI_COMM_WORLD);
125
126
127 MPI_File_open(MPI_COMM_WORLD, filename, MPI_MODE_CREATE | MPI_MODE_RDWR,
128 MPI_INFO_NULL, &fh);
129 MPI_File_set_view(fh, 0, MPI_INT, newtype, "native", MPI_INFO_NULL);
130
131 MPI_Barrier(MPI_COMM_WORLD);
132 stim = MPI_Wtime();
133 MPI_File_read_all(fh, buf, bufcount, MPI_INT, &status);
134 read_tim = MPI_Wtime() - stim;
135 MPI_File_close(&fh);
136
137 MPI_Allreduce(&read_tim, &new_read_tim, 1, MPI_DOUBLE, MPI_MAX,
138 MPI_COMM_WORLD);
139
140 if (mynod == 0) {
141 read_bw = (array_of_gsizes[0]*array_of_gsizes[1]*array_of_gsizes[2]*sizeof(int))/(new_read_tim*1024.0*1024.0);
142 fprintf(stderr, "Collective read time = %f sec, Collective read bandwidth = %f Mbytes/sec\n", new_read_tim, read_bw);
143 }
144
145 MPI_Type_free(&newtype);
146 free(buf);
147 free(filename);
148
149 MPI_Finalize();
150 return 0;
151 }