1\chapter{{\tt testWrapperMPI.c} --- A MPI Driver Program}
2\label{chapter:MPI_driver}
3
4\begin{verbatim}
5/*  testWrapperMPI.c  */
6
7#include "../BridgeMPI.h"
8
9/*--------------------------------------------------------------------*/
10int
11main ( int argc, char *argv[] ) {
12/*
13   -----------------------------------------------------------
14   purpose -- main driver program to solve a linear system
15      where the matrix and rhs are read in from files and
16      the solution is written to a file.
17      NOTE: MPI version
18
19   created -- 98sep25, cca and pjs
20   -----------------------------------------------------------
21*/
22BridgeMPI   *bridge ;
23char        *mtxFileName, *rhsFileName, *solFileName ;
24double      nfactorops ;
25FILE        *msgFile ;
26InpMtx      *mtxA ;
27int         error, msglvl, myid, neqns, nfent, nfind, nfront,
28            nproc, nrhs, nrow, nsolveops, permuteflag, rc, seed,
29            symmetryflag, type ;
30int         tstats[6] ;
31DenseMtx    *mtxX, *mtxY ;
32/*--------------------------------------------------------------------*/
33/*
34   ---------------------------------------------------------------
35   find out the identity of this process and the number of process
36   ---------------------------------------------------------------
37*/
38MPI_Init(&argc, &argv) ;
39MPI_Comm_rank(MPI_COMM_WORLD, &myid) ;
40MPI_Comm_size(MPI_COMM_WORLD, &nproc) ;
41/*--------------------------------------------------------------------*/
42/*
43    --------------------
44    get input parameters
45    --------------------
46*/
47if ( argc != 10 ) {
48   fprintf(stdout,
49           "\n\n usage : %s msglvl msgFile neqns type symmetryflag"
50           "\n         mtxFile rhsFile solFile seed"
51           "\n   msglvl  -- message level"
52           "\n      0 -- no output"
53           "\n      1 -- timings and statistics"
54           "\n      2 and greater -- lots of output"
55           "\n   msgFile -- message file"
56           "\n   neqns   -- # of equations"
57           "\n   type    -- type of entries"
58           "\n      1 -- real"
59           "\n      2 -- complex"
60           "\n   symmetryflag -- symmetry flag"
61           "\n      0 -- symmetric"
62           "\n      1 -- hermitian"
63           "\n      2 -- nonsymmetric"
64           "\n   mtxFile -- input file for A matrix InpMtx object"
65           "\n      must be *.inpmtxf or *.inpmtxb"
66           "\n   rhsFile -- input file for Y DenseMtx object"
67           "\n      must be *.densemtxf or *.densemtxb"
68           "\n   solFile -- output file for X DenseMtx object"
69           "\n      must be none, *.densemtxf or *.densemtxb"
70           "\n   seed -- random number seed"
71           "\n",
72        argv[0]) ;
73   return(0) ;
74}
75msglvl = atoi(argv[1]) ;
76if ( strcmp(argv[2], "stdout") == 0 ) {
77   msgFile = stdout ;
78} else {
79   int    length = strlen(argv[2]) + 1 + 4 ;
80   char   *buffer = CVinit(length, '\0') ;
81   sprintf(buffer, "%s.%d", argv[2], myid) ;
82   if ( (msgFile = fopen(buffer, "w")) == NULL ) {
83      fprintf(stderr, "\n fatal error in %s"
84              "\n unable to open file %s\n",
85              argv[0], argv[2]) ;
86      MPI_Finalize() ;
87      return(0) ;
88   }
89   CVfree(buffer) ;
90}
91neqns        = atoi(argv[3]) ;
92type         = atoi(argv[4]) ;
93symmetryflag = atoi(argv[5]) ;
94mtxFileName  = argv[6] ;
95rhsFileName  = argv[7] ;
96solFileName  = argv[8] ;
97seed         = atoi(argv[9]) ;
98fprintf(msgFile,
99        "\n\n %s input :"
100        "\n msglvl       = %d"
101        "\n msgFile      = %s"
102        "\n neqns        = %d"
103        "\n type         = %d"
104        "\n symmetryflag = %d"
105        "\n mtxFile      = %s"
106        "\n rhsFile      = %s"
107        "\n solFile      = %s"
108        "\n",
109        argv[0], msglvl, argv[2], neqns, type, symmetryflag,
110        mtxFileName, rhsFileName, solFileName) ;
111/*--------------------------------------------------------------------*/
112/*
113   -----------------------------------
114   processor zero reads in the matrix.
115   if an error is found,
116   all processors exit cleanly
117   -----------------------------------
118*/
119if ( myid != 0 ) {
120   mtxA = NULL ;
121} else {
122/*
123   ----------------------------------------------------
124   open the file, read in the matrix and close the file
125   ----------------------------------------------------
126*/
127   mtxA = InpMtx_new() ;
128   rc = InpMtx_readFromFile(mtxA, mtxFileName) ;
129   if ( rc != 1 ) {
130      fprintf(msgFile,
131              "\n fatal error reading mtxA from file %s, rc = %d",
132              mtxFileName, rc) ;
133      fflush(msgFile) ;
134   }
135}
136/*
137   ---------------------------------------------------------------
138   processor 0 broadcasts the error return to the other processors
139   ---------------------------------------------------------------
140*/
141MPI_Bcast((void *) &rc, 1, MPI_INT, 0, MPI_COMM_WORLD) ;
142if ( rc != 1 ) {
143   MPI_Finalize() ;
144   return(-1) ;
145}
146/*--------------------------------------------------------------------*/
147/*
148   ---------------------------------------------------
149   processor zero reads in the right hand side matrix.
150   if an error is found, all processors exit cleanly
151   ---------------------------------------------------
152*/
153if ( myid != 0 ) {
154   mtxY = NULL ;
155} else {
156/*
157   ----------------------------------
158   read in the right hand side matrix
159   ----------------------------------
160*/
161   mtxY = DenseMtx_new() ;
162   rc = DenseMtx_readFromFile(mtxY, rhsFileName) ;
163   if ( rc != 1 ) {
164      fprintf(msgFile,
165              "\n fatal error reading mtxY from file %s, rc = %d",
166              rhsFileName, rc) ;
167      fflush(msgFile) ;
168   } else {
169      DenseMtx_dimensions(mtxY, &nrow, &nrhs) ;
170   }
171}
172/*
173   ---------------------------------------------------------------
174   processor 0 broadcasts the error return to the other processors
175   ---------------------------------------------------------------
176*/
177MPI_Bcast((void *) &rc, 1, MPI_INT, 0, MPI_COMM_WORLD) ;
178if ( rc != 1 ) {
179   MPI_Finalize() ;
180   return(-1) ;
181}
182/*--------------------------------------------------------------------*/
183/*
184   ------------------------------------------
185   create and setup a BridgeMPI object
186   set the MPI, matrix and message parameters
187   ------------------------------------------
188*/
189bridge = BridgeMPI_new() ;
190BridgeMPI_setMPIparams(bridge, nproc, myid, MPI_COMM_WORLD) ;
191BridgeMPI_setMatrixParams(bridge, neqns, type, symmetryflag) ;
192BridgeMPI_setMessageInfo(bridge, msglvl, msgFile) ;
193/*
194   -----------------
195   setup the problem
196   -----------------
197*/
198rc = BridgeMPI_setup(bridge, mtxA) ;
199fprintf(msgFile,
200        "\n\n ----- SETUP -----\n"
201        "\n    CPU %8.3f : time to construct Graph"
202        "\n    CPU %8.3f : time to compress Graph"
203        "\n    CPU %8.3f : time to order Graph"
204        "\n    CPU %8.3f : time for symbolic factorization"
205        "\n    CPU %8.3f : time to broadcast front tree"
206        "\n    CPU %8.3f : time to broadcast symbolic factorization"
207        "\n CPU %8.3f : total setup time\n",
208        bridge->cpus[0], bridge->cpus[1], bridge->cpus[2],
209        bridge->cpus[3], bridge->cpus[4], bridge->cpus[5],
210        bridge->cpus[6]) ;
211rc = BridgeMPI_factorStats(bridge, type, symmetryflag, &nfront,
212                           &nfind, &nfent, &nsolveops, &nfactorops) ;
213if ( rc != 1 ) {
214   fprintf(stderr,
215           "\n error return %d from BridgeMPI_factorStats()", rc) ;
216   MPI_Finalize() ;
217   exit(-1) ;
218}
219fprintf(msgFile,
220        "\n\n factor matrix statistics"
221        "\n %d fronts, %d indices, %d entries"
222        "\n %d solve operations, %12.4e factor operations",
223        nfront, nfind, nfent, nsolveops, nfactorops) ;
224fflush(msgFile) ;
225/*--------------------------------------------------------------------*/
226/*
227   --------------------------------
228   setup the parallel factorization
229   --------------------------------
230*/
231rc = BridgeMPI_factorSetup(bridge, 0, 0.0) ;
232if ( rc != 1 ) {
233   fprintf(stderr,
234           "\n error return %d from BridgeMPI_factorSetup()", rc) ;
235   MPI_Finalize() ;
236   exit(-1) ;
237}
238fprintf(msgFile, "\n\n ----- PARALLEL FACTOR SETUP -----\n") ;
239fprintf(msgFile,
240        "\n    CPU %8.3f : time to setup parallel factorization",
241        bridge->cpus[7]) ;
242if ( msglvl > 0 ) {
243   fprintf(msgFile, "\n total factor operations = %.0f"
244           "\n upper bound on speedup due to load balance = %.2f",
245           DV_sum(bridge->cumopsDV),
246           DV_sum(bridge->cumopsDV)/DV_max(bridge->cumopsDV)) ;
247   fprintf(msgFile, "\n operations distributions over processors") ;
248   DV_writeForHumanEye(bridge->cumopsDV, msgFile) ;
249   fflush(msgFile) ;
250}
251/*--------------------------------------------------------------------*/
252/*
253   ------------------------------------------------------
254   set the factorization parameters and factor the matrix
255   ------------------------------------------------------
256*/
257permuteflag = 1 ;
258rc = BridgeMPI_factor(bridge, mtxA, permuteflag, &error) ;
259fprintf(msgFile, "\n\n ----- FACTORIZATION -----\n") ;
260if ( rc == 1 ) {
261   fprintf(msgFile, "\n\n factorization completed successfully\n") ;
262} else {
263   fprintf(msgFile, "\n"
264           "\n return code from factorization = %d\n"
265           "\n error code                     = %d\n",
266           rc, error) ;
267   MPI_Finalize() ;
268   exit(-1) ;
269}
270fprintf(msgFile,
271        "\n    CPU %8.3f : time to permute original matrix"
272        "\n    CPU %8.3f : time to distribute original matrix"
273        "\n    CPU %8.3f : time to initialize factor matrix"
274        "\n    CPU %8.3f : time to compute factorization"
275        "\n    CPU %8.3f : time to post-process factorization"
276        "\n CPU %8.3f : total factorization time\n",
277        bridge->cpus[8], bridge->cpus[9], bridge->cpus[10],
278        bridge->cpus[11], bridge->cpus[12], bridge->cpus[13]) ;
279IVzero(6, tstats) ;
280MPI_Reduce((void *) bridge->stats, (void *) tstats, 6, MPI_INT,
281           MPI_SUM, 0, bridge->comm) ;
282fprintf(msgFile,
283        "\n\n   factorization statistics"
284        "\n   %d pivots, %d pivot tests, %d delayed vertices"
285        "\n   %d entries in D, %d entries in L, %d entries in U",
286        tstats[0], tstats[1], tstats[2],
287        tstats[3], tstats[4], tstats[5]) ;
288fprintf(msgFile,
289        "\n\n   factorization: raw mflops %8.3f, overall mflops %8.3f",
290        1.e-6*nfactorops/bridge->cpus[11],
291        1.e-6*nfactorops/bridge->cpus[13]) ;
292fflush(msgFile) ;
293/*--------------------------------------------------------------------*/
294/*
295   ------------------------
296   setup the parallel solve
297   ------------------------
298*/
299rc = BridgeMPI_solveSetup(bridge) ;
300fprintf(msgFile, "\n\n ----- PARALLEL SOLVE SETUP -----\n"
301        "\n    CPU %8.3f : time to setup parallel solve",
302        bridge->cpus[14]) ;
303if ( rc != 1 ) {
304   fprintf(stderr,
305           "\n error return %d from BridgeMPI_solveSetup()", rc) ;
306   MPI_Finalize() ;
307   exit(-1) ;
308}
309/*--------------------------------------------------------------------*/
310/*
311   -----------------------------------------
312   processor 0 initializes a DenseMtx object
313   to hold the global solution matrix
314   -----------------------------------------
315*/
316if ( myid == 0 ) {
317   mtxX = DenseMtx_new() ;
318   DenseMtx_init(mtxX, type, 0, 0, neqns, nrhs, 1, neqns) ;
319   DenseMtx_zero(mtxX) ;
320} else {
321   mtxX = NULL ;
322}
323/*
324   ---------------------------------------------
325   the processors solve the system cooperatively
326   ---------------------------------------------
327*/
328permuteflag = 1 ;
329rc = BridgeMPI_solve(bridge, permuteflag, mtxX, mtxY) ;
330if ( rc == 1 ) {
331   fprintf(msgFile, "\n\n solve complete successfully\n") ;
332} else {
333   fprintf(msgFile, "\n" " return code from solve = %d\n", rc) ;
334}
335fprintf(msgFile, "\n\n ----- SOLVE -----\n"
336        "\n    CPU %8.3f : time to permute rhs into new ordering"
337        "\n    CPU %8.3f : time to distribute rhs "
338        "\n    CPU %8.3f : time to initialize solution matrix "
339        "\n    CPU %8.3f : time to solve linear system"
340        "\n    CPU %8.3f : time to gather solution "
341        "\n    CPU %8.3f : time to permute solution into old ordering"
342        "\n CPU %8.3f : total solve time"
343        "\n\n solve: raw mflops %8.3f, overall mflops %8.3f",
344        bridge->cpus[15], bridge->cpus[16], bridge->cpus[17],
345        bridge->cpus[18], bridge->cpus[19], bridge->cpus[20],
346        bridge->cpus[21],
347        1.e-6*nsolveops/bridge->cpus[18],
348        1.e-6*nsolveops/bridge->cpus[21]) ;
349fflush(msgFile) ;
350if ( myid == 0 ) {
351   if ( msglvl > 0 ) {
352      fprintf(msgFile, "\n\n solution matrix in original ordering") ;
353      DenseMtx_writeForHumanEye(mtxX, msgFile) ;
354      fflush(msgFile) ;
355   }
356}
357/*--------------------------------------------------------------------*/
358/*
359   ---------------------
360   free the working data
361   ---------------------
362*/
363if ( myid == 0 ) {
364   InpMtx_free(mtxA) ;
365   DenseMtx_free(mtxX) ;
366   DenseMtx_free(mtxY) ;
367}
368BridgeMPI_free(bridge) ;
369
370/*--------------------------------------------------------------------*/
371
372MPI_Finalize() ;
373
374return(1) ; }
375
376/*--------------------------------------------------------------------*/
377\end{verbatim}
378