1 /* testGraph_Bcast.c */
2
3 #include "../spoolesMPI.h"
4 #include "../../Drand.h"
5 #include "../../timings.h"
6
7 /*--------------------------------------------------------------------*/
8 int
main(int argc,char * argv[])9 main ( int argc, char *argv[] )
10 /*
11 --------------------------------------------------------------------
12 this program tests the Graph_MPI_Bcast() method
13
14 (1) process root generates a random Graph object
15 and computes its checksum
16 (2) process root broadcasts the Graph object to the other processors
17 (3) each process computes the checksum of its Graph object
18 (4) the checksums are compared on root
19
20 created -- 98sep10, cca
21 --------------------------------------------------------------------
22 */
23 {
24 char *buffer ;
25 double chksum, t1, t2 ;
26 double *sums ;
27 Drand drand ;
28 int iproc, length, loc, msglvl, myid, nitem, nproc,
29 nvtx, root, seed, size, type, v ;
30 int *list ;
31 FILE *msgFile ;
32 Graph *graph ;
33 /*
34 ---------------------------------------------------------------
35 find out the identity of this process and the number of process
36 ---------------------------------------------------------------
37 */
38 MPI_Init(&argc, &argv) ;
39 MPI_Comm_rank(MPI_COMM_WORLD, &myid) ;
40 MPI_Comm_size(MPI_COMM_WORLD, &nproc) ;
41 fprintf(stdout, "\n process %d of %d, argc = %d", myid, nproc, argc) ;
42 fflush(stdout) ;
43 if ( argc != 8 ) {
44 fprintf(stdout,
45 "\n\n usage : %s msglvl msgFile type nvtx nitem root seed "
46 "\n msglvl -- message level"
47 "\n msgFile -- message file"
48 "\n type -- type of graph"
49 "\n nvtx -- # of vertices"
50 "\n nitem -- # of items used to generate graph"
51 "\n root -- root processor for broadcast"
52 "\n seed -- random number seed"
53 "\n", argv[0]) ;
54 return(0) ;
55 }
56 msglvl = atoi(argv[1]) ;
57 if ( strcmp(argv[2], "stdout") == 0 ) {
58 msgFile = stdout ;
59 } else {
60 length = strlen(argv[2]) + 1 + 4 ;
61 buffer = CVinit(length, '\0') ;
62 sprintf(buffer, "%s.%d", argv[2], myid) ;
63 if ( (msgFile = fopen(buffer, "w")) == NULL ) {
64 fprintf(stderr, "\n fatal error in %s"
65 "\n unable to open file %s\n",
66 argv[0], argv[2]) ;
67 return(-1) ;
68 }
69 CVfree(buffer) ;
70 }
71 type = atoi(argv[3]) ;
72 nvtx = atoi(argv[4]) ;
73 nitem = atoi(argv[5]) ;
74 root = atoi(argv[6]) ;
75 seed = atoi(argv[7]) ;
76 fprintf(msgFile,
77 "\n %s "
78 "\n msglvl -- %d"
79 "\n msgFile -- %s"
80 "\n type -- %d"
81 "\n nvtx -- %d"
82 "\n nitem -- %d"
83 "\n root -- %d"
84 "\n seed -- %d"
85 "\n",
86 argv[0], msglvl, argv[2], type, nvtx, nitem, root, seed) ;
87 fflush(msgFile) ;
88 /*
89 -----------------------
90 set up the Graph object
91 -----------------------
92 */
93 MARKTIME(t1) ;
94 graph = Graph_new() ;
95 if ( myid == root ) {
96 InpMtx *inpmtx ;
97 int nedges, totewght, totvwght, v ;
98 int *adj, *vwghts ;
99 IVL *adjIVL, *ewghtIVL ;
100 /*
101 -----------------------
102 generate a random graph
103 -----------------------
104 */
105 inpmtx = InpMtx_new() ;
106 InpMtx_init(inpmtx, INPMTX_BY_ROWS, INPMTX_INDICES_ONLY, nitem, 0) ;
107 Drand_setDefaultFields(&drand) ;
108 Drand_setSeed(&drand, seed) ;
109 Drand_setUniform(&drand, 0, nvtx) ;
110 Drand_fillIvector(&drand, nitem, InpMtx_ivec1(inpmtx)) ;
111 Drand_fillIvector(&drand, nitem, InpMtx_ivec2(inpmtx)) ;
112 InpMtx_setNent(inpmtx, nitem) ;
113 InpMtx_changeStorageMode(inpmtx, INPMTX_BY_VECTORS) ;
114 if ( msglvl > 2 ) {
115 fprintf(msgFile, "\n\n inpmtx mtx filled with raw entries") ;
116 InpMtx_writeForHumanEye(inpmtx, msgFile) ;
117 fflush(msgFile) ;
118 }
119 adjIVL = InpMtx_fullAdjacency(inpmtx) ;
120 if ( msglvl > 2 ) {
121 fprintf(msgFile, "\n\n full adjacency structure") ;
122 IVL_writeForHumanEye(adjIVL, msgFile) ;
123 fflush(msgFile) ;
124 }
125 nedges = adjIVL->tsize ;
126 if ( type == 1 || type == 3 ) {
127 Drand_setUniform(&drand, 1, 10) ;
128 vwghts = IVinit(nvtx, 0) ;
129 Drand_fillIvector(&drand, nvtx, vwghts) ;
130 totvwght = IVsum(nvtx, vwghts) ;
131 if ( msglvl > 2 ) {
132 fprintf(msgFile, "\n\n vertex weights") ;
133 IVfprintf(msgFile, nvtx, vwghts) ;
134 fflush(msgFile) ;
135 }
136 } else {
137 vwghts = NULL ;
138 totvwght = nvtx ;
139 }
140 if ( msglvl > 2 ) {
141 fprintf(msgFile, "\n\n totvwght %d", totvwght) ;
142 fflush(msgFile) ;
143 }
144 if ( type == 2 || type == 3 ) {
145 ewghtIVL = IVL_new() ;
146 IVL_init1(ewghtIVL, IVL_CHUNKED, nvtx) ;
147 Drand_setUniform(&drand, 1, 100) ;
148 totewght = 0 ;
149 for ( v = 0 ; v < nvtx ; v++ ) {
150 IVL_listAndSize(adjIVL, v, &size, &adj) ;
151 IVL_setList(ewghtIVL, v, size, NULL) ;
152 IVL_listAndSize(ewghtIVL, v, &size, &adj) ;
153 Drand_fillIvector(&drand, size, adj) ;
154 totewght += IVsum(size, adj) ;
155 }
156 if ( msglvl > 2 ) {
157 fprintf(msgFile, "\n\n ewghtIVL") ;
158 IVL_writeForHumanEye(ewghtIVL, msgFile) ;
159 fflush(msgFile) ;
160 }
161 } else {
162 ewghtIVL = NULL ;
163 totewght = nedges ;
164 }
165 if ( msglvl > 2 ) {
166 fprintf(msgFile, "\n\n totewght %d", totewght) ;
167 fflush(msgFile) ;
168 }
169 Graph_init2(graph, type, nvtx, 0, nedges, totvwght, totewght,
170 adjIVL, vwghts, ewghtIVL) ;
171 InpMtx_free(inpmtx) ;
172 }
173 MARKTIME(t2) ;
174 fprintf(msgFile,
175 "\n CPU %8.3f : initialize the Graph object", t2 - t1) ;
176 fflush(msgFile) ;
177 if ( msglvl > 2 ) {
178 Graph_writeForHumanEye(graph, msgFile) ;
179 } else {
180 Graph_writeStats(graph, msgFile) ;
181 }
182 fflush(msgFile) ;
183 if ( myid == root ) {
184 /*
185 ----------------------------------------
186 compute the checksum of the Graph object
187 ----------------------------------------
188 */
189 chksum = graph->type + graph->nvtx + graph->nvbnd
190 + graph->nedges + graph->totvwght + graph->totewght ;
191 for ( v = 0 ; v < nvtx ; v++ ) {
192 IVL_listAndSize(graph->adjIVL, v, &size, &list) ;
193 chksum += 1 + v + size + IVsum(size, list) ;
194 }
195 if ( graph->vwghts != NULL ) {
196 chksum += IVsum(nvtx, graph->vwghts) ;
197 }
198 if ( graph->ewghtIVL != NULL ) {
199 for ( v = 0 ; v < nvtx ; v++ ) {
200 IVL_listAndSize(graph->ewghtIVL, v, &size, &list) ;
201 chksum += 1 + v + size + IVsum(size, list) ;
202 }
203 }
204 fprintf(msgFile, "\n\n local chksum = %12.4e", chksum) ;
205 fflush(msgFile) ;
206 }
207 /*
208 --------------------------
209 broadcast the Graph object
210 --------------------------
211 */
212 MARKTIME(t1) ;
213 graph = Graph_MPI_Bcast(graph, root, msglvl, msgFile, MPI_COMM_WORLD) ;
214 MARKTIME(t2) ;
215 fprintf(msgFile, "\n CPU %8.3f : broadcast the Graph object", t2 - t1) ;
216 if ( msglvl > 2 ) {
217 Graph_writeForHumanEye(graph, msgFile) ;
218 } else {
219 Graph_writeStats(graph, msgFile) ;
220 }
221 /*
222 ----------------------------------------
223 compute the checksum of the Graph object
224 ----------------------------------------
225 */
226 chksum = graph->type + graph->nvtx + graph->nvbnd
227 + graph->nedges + graph->totvwght + graph->totewght ;
228 for ( v = 0 ; v < nvtx ; v++ ) {
229 IVL_listAndSize(graph->adjIVL, v, &size, &list) ;
230 chksum += 1 + v + size + IVsum(size, list) ;
231 }
232 if ( graph->vwghts != NULL ) {
233 chksum += IVsum(nvtx, graph->vwghts) ;
234 }
235 if ( graph->ewghtIVL != NULL ) {
236 for ( v = 0 ; v < nvtx ; v++ ) {
237 IVL_listAndSize(graph->ewghtIVL, v, &size, &list) ;
238 chksum += 1 + v + size + IVsum(size, list) ;
239 }
240 }
241 fprintf(msgFile, "\n\n local chksum = %12.4e", chksum) ;
242 fflush(msgFile) ;
243 /*
244 ---------------------------------------
245 gather the checksums from the processes
246 ---------------------------------------
247 */
248 sums = DVinit(nproc, 0.0) ;
249 MPI_Gather((void *) &chksum, 1, MPI_DOUBLE,
250 (void *) sums, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD) ;
251 if ( myid == 0 ) {
252 fprintf(msgFile, "\n\n sums") ;
253 DVfprintf(msgFile, nproc, sums) ;
254 for ( iproc = 0 ; iproc < nproc ; iproc++ ) {
255 sums[iproc] -= chksum ;
256 }
257 fprintf(msgFile, "\n\n errors") ;
258 DVfprintf(msgFile, nproc, sums) ;
259 fprintf(msgFile, "\n\n maxerror = %12.4e", DVmax(nproc, sums, &loc));
260 }
261 /*
262 ----------------
263 free the objects
264 ----------------
265 */
266 DVfree(sums) ;
267 Graph_free(graph) ;
268 /*
269 ------------------------
270 exit the MPI environment
271 ------------------------
272 */
273 MPI_Finalize() ;
274
275 fprintf(msgFile, "\n") ;
276 fclose(msgFile) ;
277
278 return(0) ; }
279
280 /*--------------------------------------------------------------------*/
281