1 /* testIVL_alltoall.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 IVL_MPI_alltoall() method.
13
14 (1) each process creates a "send" IVL object
15 and fills the lists with random numbers
16 (3) the processes all-to-all gather the lists
17 to create a "recv" IVL object
18
19 created -- 98jul26, cca
20 -------------------------------------------------
21 */
22 {
23 char *buffer ;
24 double chksum, globalsum1, globalsum2, t1, t2 ;
25 Drand drand ;
26 int ilist, length, myid, msglvl, n, nlist,
27 nproc, rc, seed, size, tag ;
28 int *list, *vec ;
29 int stats[4], tstats[4] ;
30 IVL *recvIVL, *sendIVL ;
31 FILE *msgFile ;
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 != 5 ) {
43 fprintf(stdout,
44 "\n\n usage : %s msglvl msgFile n seed "
45 "\n msglvl -- message level"
46 "\n msgFile -- message file"
47 "\n n -- maximum size of any list"
48 "\n seed -- random number seed"
49 "\n", argv[0]) ;
50 return(0) ;
51 }
52 msglvl = atoi(argv[1]) ;
53 if ( strcmp(argv[2], "stdout") == 0 ) {
54 msgFile = stdout ;
55 } else {
56 length = strlen(argv[2]) + 1 + 4 ;
57 buffer = CVinit(length, '\0') ;
58 sprintf(buffer, "%s.%d", argv[2], myid) ;
59 if ( (msgFile = fopen(buffer, "w")) == NULL ) {
60 fprintf(stderr, "\n fatal error in %s"
61 "\n unable to open file %s\n",
62 argv[0], argv[2]) ;
63 return(-1) ;
64 }
65 CVfree(buffer) ;
66 }
67 n = atoi(argv[3]) ;
68 seed = atoi(argv[4]) ;
69 fprintf(msgFile,
70 "\n %s "
71 "\n msglvl -- %d"
72 "\n msgFile -- %s"
73 "\n n -- %d"
74 "\n seed -- %d"
75 "\n",
76 argv[0], msglvl, argv[2], n, seed) ;
77 fflush(msgFile) ;
78 /*
79 ----------------------------
80 set up the "send" IVL object
81 ----------------------------
82 */
83 MARKTIME(t1) ;
84 sendIVL = IVL_new() ;
85 nlist = nproc ;
86 IVL_init1(sendIVL, IVL_CHUNKED, nlist) ;
87 vec = IVinit(n, -1) ;
88 Drand_setDefaultFields(&drand) ;
89 Drand_setSeed(&drand, seed + myid) ;
90 for ( ilist = 0 ; ilist < nlist ; ilist++ ) {
91 Drand_setUniform(&drand, 0, n) ;
92 size = (int) Drand_value(&drand) ;
93 Drand_setUniform(&drand, 0, n*n) ;
94 Drand_fillIvector(&drand, size, vec) ;
95 IVqsortUp(size, vec) ;
96 IVL_setList(sendIVL, ilist, size, vec) ;
97 }
98 MARKTIME(t2) ;
99 fprintf(msgFile, "\n CPU %8.3f : initialize the IVL object",
100 t2 - t1) ;
101 fflush(msgFile) ;
102 if ( msglvl > 2 ) {
103 IVL_writeForHumanEye(sendIVL, msgFile) ;
104 } else {
105 IVL_writeStats(sendIVL, msgFile) ;
106 }
107 fflush(msgFile) ;
108 /*
109 --------------------------------------------
110 compute the local checksum of the ivl object
111 --------------------------------------------
112 */
113 for ( ilist = 0, chksum = 0.0 ; ilist < nlist ; ilist++ ) {
114 IVL_listAndSize(sendIVL, ilist, &size, &list) ;
115 chksum += 1 + ilist + size + IVsum(size, list) ;
116 }
117 fprintf(msgFile, "\n\n local partial chksum = %12.4e", chksum) ;
118 fflush(msgFile) ;
119 /*
120 -----------------------------
121 get the first global checksum
122 -----------------------------
123 */
124 rc = MPI_Allreduce((void *) &chksum, (void *) &globalsum1,
125 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD) ;
126 /*
127 --------------------------------
128 execute the all-to-all operation
129 --------------------------------
130 */
131 tag = 47 ;
132 IVzero(4, stats) ;
133 recvIVL = IVL_MPI_alltoall(sendIVL, NULL, stats,
134 msglvl, msgFile, tag, MPI_COMM_WORLD) ;
135 if ( msglvl > 0 ) {
136 fprintf(msgFile, "\n\n return from IVL_MPI_alltoall()") ;
137 fprintf(msgFile,
138 "\n local send stats : %10d messages with %10d bytes"
139 "\n local recv stats : %10d messages with %10d bytes",
140 stats[0], stats[2], stats[1], stats[3]) ;
141 fflush(msgFile) ;
142 }
143 MPI_Reduce((void *) stats, (void *) tstats, 4, MPI_INT,
144 MPI_SUM, 0, MPI_COMM_WORLD) ;
145 if ( myid == 0 ) {
146 fprintf(msgFile,
147 "\n total send stats : %10d messages with %10d bytes"
148 "\n total recv stats : %10d messages with %10d bytes",
149 tstats[0], tstats[2], tstats[1], tstats[3]) ;
150 fflush(msgFile) ;
151 }
152 if ( msglvl > 2 ) {
153 fprintf(msgFile, "\n\n recvIVL") ;
154 IVL_writeForHumanEye(recvIVL, msgFile) ;
155 fflush(msgFile) ;
156 }
157 /*
158 ----------------------------------
159 compute the second global checksum
160 ----------------------------------
161 */
162 for ( ilist = 0, chksum = 0.0 ; ilist < nlist ; ilist++ ) {
163 IVL_listAndSize(recvIVL, ilist, &size, &list) ;
164 chksum += 1 + ilist + size + IVsum(size, list) ;
165 }
166 rc = MPI_Allreduce((void *) &chksum, (void *) &globalsum2,
167 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD) ;
168 fprintf(msgFile,
169 "\n globalsum1 = %12.4e, globalsum2 = %12.4e, error = %12.4e",
170 globalsum1, globalsum2, fabs(globalsum1 - globalsum2)) ;
171 fflush(msgFile) ;
172 /*
173 ----------------
174 free the objects
175 ----------------
176 */
177 IVL_free(sendIVL) ;
178 IVL_free(recvIVL) ;
179 IVfree(vec) ;
180 /*
181 ------------------------
182 exit the MPI environment
183 ------------------------
184 */
185 MPI_Finalize() ;
186
187 fprintf(msgFile, "\n") ;
188 fclose(msgFile) ;
189
190 return(0) ; }
191
192 /*--------------------------------------------------------------------*/
193