1 /* testScatterInpMtx.c */
2
3 #include "../spoolesMPI.h"
4 #include "../../Drand.h"
5 #include "../../timings.h"
6
7 /*--------------------------------------------------------------------*/
8
9 int
main(int argc,char * argv[])10 main ( int argc, char *argv[] )
11 /*
12 ---------------------------------------------------
13 (1) process root generates an InpMtx object
14 (2) using a random map, it distributes the object
15
16 created -- 98sep27, cca
17 ---------------------------------------------------
18 */
19 {
20 char *buffer ;
21 InpMtx *mtxA, *mtxB ;
22 double t1, t2 ;
23 double schecksums[3], pchecksums[3], tchecksums[3] ;
24 Drand *drand ;
25 int coordType, inputMode, iproc, length, myid,
26 msglvl, nent, neqns, nproc, rc, root, seed, tag, v ;
27 int ibuffer[3], stats[4], tstats[4] ;
28 int *map ;
29 IV *mapIV ;
30 FILE *msgFile ;
31 MPI_Status status ;
32 /*
33 ---------------------------------------------------------------
34 find out the identity of this process and the number of process
35 ---------------------------------------------------------------
36 */
37 MPI_Init(&argc, &argv) ;
38 MPI_Comm_rank(MPI_COMM_WORLD, &myid) ;
39 MPI_Comm_size(MPI_COMM_WORLD, &nproc) ;
40 fprintf(stdout, "\n process %d of %d, argc = %d", myid, nproc, argc) ;
41 fflush(stdout) ;
42 if ( argc != 9 ) {
43 fprintf(stdout,
44 "\n\n usage : %s msglvl msgFile neqns nent seed "
45 "\n coordType inputMode root"
46 "\n msglvl -- message level"
47 "\n msgFile -- message file"
48 "\n neqns -- number of equations"
49 "\n nent -- number of equations"
50 "\n seed -- random number seed"
51 "\n coordType -- coordinate type"
52 "\n 1 -- store by rows"
53 "\n 2 -- store by columns"
54 "\n 3 -- store by chevrons"
55 "\n inputMode -- input mode"
56 "\n 0 -- indices only"
57 "\n 1 -- indices and real entries"
58 "\n 2 -- indices and complex entries"
59 "\n root -- root processor that creates the global matrix"
60 "\n", argv[0]) ;
61 return(0) ;
62 }
63 msglvl = atoi(argv[1]) ;
64 if ( strcmp(argv[2], "stdout") == 0 ) {
65 msgFile = stdout ;
66 } else {
67 length = strlen(argv[2]) + 1 + 4 ;
68 buffer = CVinit(length, '\0') ;
69 sprintf(buffer, "%s.%d", argv[2], myid) ;
70 if ( (msgFile = fopen(buffer, "w")) == NULL ) {
71 fprintf(stderr, "\n fatal error in %s"
72 "\n unable to open file %s\n",
73 argv[0], argv[2]) ;
74 return(-1) ;
75 }
76 CVfree(buffer) ;
77 }
78 neqns = atoi(argv[3]) ;
79 nent = atoi(argv[4]) ;
80 seed = atoi(argv[5]) ;
81 coordType = atoi(argv[6]) ;
82 inputMode = atoi(argv[7]) ;
83 root = atoi(argv[8]) ;
84 fprintf(msgFile,
85 "\n %s "
86 "\n msglvl -- %d"
87 "\n msgFile -- %s"
88 "\n neqns -- %d"
89 "\n nent -- %d"
90 "\n seed -- %d"
91 "\n coordType -- %d"
92 "\n inputMode -- %d"
93 "\n root -- %d"
94 "\n",
95 argv[0], msglvl, argv[2], neqns, nent, seed,
96 coordType, inputMode, root) ;
97 fflush(msgFile) ;
98 /*
99 ----------------------------------
100 create the random number generator
101 ----------------------------------
102 */
103 if ( myid == root ) {
104 /*
105 ----------------------------------------------
106 processor root, generate a random input matrix
107 ----------------------------------------------
108 */
109 MARKTIME(t1) ;
110 mtxA = InpMtx_new() ;
111 InpMtx_init(mtxA, INPMTX_BY_ROWS, inputMode, nent, 0) ;
112 drand = Drand_new() ;
113 Drand_setSeed(drand, seed) ;
114 Drand_setUniform(drand, 0, neqns) ;
115 Drand_fillIvector(drand, nent, InpMtx_ivec1(mtxA)) ;
116 Drand_fillIvector(drand, nent, InpMtx_ivec2(mtxA)) ;
117 if ( inputMode == SPOOLES_REAL ) {
118 Drand_setUniform(drand, -1.0, 1.0) ;
119 Drand_fillDvector(drand, nent, InpMtx_dvec(mtxA)) ;
120 } else if ( inputMode == SPOOLES_COMPLEX ) {
121 Drand_setUniform(drand, -1.0, 1.0) ;
122 Drand_fillDvector(drand, 2*nent, InpMtx_dvec(mtxA)) ;
123 }
124 mtxA->nent = nent ;
125 InpMtx_sortAndCompress(mtxA) ;
126 InpMtx_changeCoordType(mtxA, coordType) ;
127 MARKTIME(t2) ;
128 fprintf(msgFile, "\n CPU %9.5f : generate random matrix", t2 - t1) ;
129 /*
130 -----------------------------------------------
131 change the coordinate type and the storage mode
132 -----------------------------------------------
133 */
134 InpMtx_changeCoordType(mtxA, coordType) ;
135 InpMtx_changeStorageMode(mtxA, INPMTX_BY_VECTORS) ;
136 mtxA->inputMode = inputMode ;
137 fprintf(msgFile, "\n\n random input matrix") ;
138 if ( msglvl > 2 ) {
139 InpMtx_writeForHumanEye(mtxA, msgFile) ;
140 } else {
141 InpMtx_writeStats(mtxA, msgFile) ;
142 }
143 fflush(msgFile) ;
144 /*
145 ------------------------------------------
146 compute the serial checksums of the matrix
147 ------------------------------------------
148 */
149 InpMtx_checksums(mtxA, schecksums) ;
150 fprintf(msgFile, "\n schecksums %12.4e %12.4e %12.4e",
151 schecksums[0], schecksums[1], schecksums[2]) ;
152 fflush(msgFile) ;
153 /*
154 ---------------------
155 generate a random map
156 ---------------------
157 */
158 MARKTIME(t1) ;
159 Drand_setSeed(drand, seed + 1) ;
160 Drand_setUniform(drand, 0.0, (double) neqns) ;
161 mapIV = IV_new() ;
162 IV_init(mapIV, neqns, NULL) ;
163 map = IV_entries(mapIV) ;
164 for ( v = 0 ; v < neqns ; v++ ) {
165 map[v] = ((int) Drand_value(drand)) % nproc ;
166 }
167 MARKTIME(t2) ;
168 fprintf(msgFile, "\n CPU %9.5f : random map set", t2 - t1) ;
169 if ( msglvl > 2 ) {
170 fprintf(msgFile, "\n\n map") ;
171 IV_writeForHumanEye(mapIV, msgFile) ;
172 fflush(msgFile) ;
173 }
174 } else {
175 mtxA = NULL ;
176 mapIV = NULL ;
177 }
178 /*
179 ------------------------------------------
180 now scatter the matrix with the random map
181 ------------------------------------------
182 */
183 stats[0] = stats[1] = stats[2] = stats[3] = 0 ;
184 tstats[0] = tstats[1] = tstats[2] = tstats[3] = 0 ;
185 MARKTIME(t1) ;
186 mtxB = InpMtx_MPI_splitFromGlobal(mtxA, NULL, mapIV, root, stats,
187 msglvl, msgFile, tag, MPI_COMM_WORLD) ;
188 MARKTIME(t2) ;
189 fprintf(msgFile, "\n CPU %9.5f : matrix split", t2 - t1) ;
190 fprintf(msgFile,
191 "\n send stats : %d messages with %d bytes"
192 "\n recv stats : %d messages with %d bytes",
193 stats[0], stats[2], stats[1], stats[3]) ;
194 MPI_Reduce((void *) stats, (void *) tstats, 4, MPI_INT,
195 MPI_SUM, root, MPI_COMM_WORLD) ;
196 if ( myid == root ) {
197 fprintf(msgFile,
198 "\n total send stats : %d messages with %d bytes"
199 "\n total recv stats : %d messages with %d bytes",
200 tstats[0], tstats[2], tstats[1], tstats[3]) ;
201 fflush(msgFile) ;
202 }
203 if ( mtxA != NULL ) {
204 InpMtx_free(mtxA) ;
205 }
206 mtxA = mtxB ;
207 fprintf(msgFile,
208 "\n\n after splitting the InpMtx object with the owners map") ;
209 if ( msglvl > 2 ) {
210 InpMtx_writeForHumanEye(mtxA, msgFile) ;
211 } else {
212 InpMtx_writeStats(mtxA, msgFile) ;
213 }
214 fflush(msgFile) ;
215 /*
216 -----------------------------------------------
217 compute the checksums of the distributed matrix
218 -----------------------------------------------
219 */
220 InpMtx_checksums(mtxA, pchecksums) ;
221 MPI_Reduce((void *) pchecksums, (void *) tchecksums, 3, MPI_DOUBLE,
222 MPI_SUM, root, MPI_COMM_WORLD) ;
223 if ( myid == root ) {
224 fprintf(msgFile,
225 "\n\n checksums for original matrix : %12.4e %12.4e %12.4e"
226 "\n checksums for distributed matrix : %12.4e %12.4e %12.4e"
227 "\n error in checksum : %12.4e %12.4e %12.4e",
228 schecksums[0], schecksums[1], schecksums[2],
229 tchecksums[0], tchecksums[1], tchecksums[2],
230 tchecksums[0] - schecksums[0],
231 tchecksums[1] - schecksums[1],
232 tchecksums[2] - schecksums[2]) ;
233 }
234 /*
235 ----------------
236 free the objects
237 ----------------
238 */
239 if ( myid == root ) {
240 IV_free(mapIV) ;
241 Drand_free(drand) ;
242 }
243 InpMtx_free(mtxA) ;
244
245 MPI_Finalize() ;
246
247 fprintf(msgFile, "\n") ;
248 fclose(msgFile) ;
249
250 return(0) ; }
251
252 /*--------------------------------------------------------------------*/
253