1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
2 * Copyright by The HDF Group. *
3 * Copyright by the Board of Trustees of the University of Illinois. *
4 * All rights reserved. *
5 * *
6 * This file is part of HDF5. The full HDF5 copyright notice, including *
7 * terms governing use, modification, and redistribution, is contained in *
8 * the COPYING file, which can be found at the root of the source code *
9 * distribution tree, or in https://support.hdfgroup.org/ftp/HDF5/releases. *
10 * If you do not have access to either file, you may request a copy from *
11 * help@hdfgroup.org. *
12 * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
13
14 #include "H5private.h"
15 #include "h5diff.h"
16 #include "ph5diff.h"
17 #include "h5diff_common.h"
18 #include "h5tools.h"
19 #include "h5tools_utils.h"
20
21 /* Name of tool */
22 #define PROGRAMNAME "h5diff"
23
24 static void ph5diff_worker(int );
25
26 /*-------------------------------------------------------------------------
27 * Function: main
28 *
29 * Purpose: ph5diff main program
30 *
31 * Return: An exit status of 0 means no differences were found, 1 means some
32 * differences were found.
33 *
34 * Programmer: Pedro Vicente, pvn@ncsa.uiuc.edu
35 *
36 * Date: May 9, 2003
37 *
38 * Comments:
39 *
40 * This function drives the diff process and will do a serial or parallel diff depending
41 * on the value of the global variable g_Parallel (default is 0), set to 1 when the program
42 * is run as "ph5diff"
43 *-------------------------------------------------------------------------
44 */
45
main(int argc,const char * argv[])46 int main(int argc, const char *argv[])
47 {
48 int nID = 0;
49 const char *fname1 = NULL;
50 const char *fname2 = NULL;
51 const char *objname1 = NULL;
52 const char *objname2 = NULL;
53 diff_opt_t opts;
54
55 h5tools_setprogname(PROGRAMNAME);
56 h5tools_setstatus(EXIT_SUCCESS);
57
58 /* Initialize h5tools lib */
59 h5tools_init();
60
61 outBuffOffset = 0;
62 g_Parallel = 1;
63
64 MPI_Init(&argc, (char***) &argv);
65
66 MPI_Comm_rank(MPI_COMM_WORLD, &nID);
67 MPI_Comm_size(MPI_COMM_WORLD, &g_nTasks);
68
69 if(g_nTasks == 1)
70 {
71 HDprintf("Only 1 task available...doing serial diff\n");
72
73 g_Parallel = 0;
74
75 parse_command_line(argc, argv, &fname1, &fname2, &objname1, &objname2, &opts);
76
77 h5diff(fname1, fname2, objname1, objname2, &opts);
78
79 print_info(&opts);
80 }
81 /* Parallel h5diff */
82 else {
83
84 /* Have the manager process the command-line */
85 if(nID == 0)
86 {
87 parse_command_line(argc, argv, &fname1, &fname2, &objname1, &objname2, &opts);
88
89 h5diff(fname1, fname2, objname1, objname2, &opts);
90
91 MPI_Barrier(MPI_COMM_WORLD);
92
93 print_info(&opts);
94 print_manager_output();
95 }
96 /* All other tasks become workers and wait for assignments. */
97 else {
98 ph5diff_worker(nID);
99
100 MPI_Barrier(MPI_COMM_WORLD);
101 } /* end else */
102
103 } /* end else */
104
105 MPI_Finalize();
106
107 return 0;
108 }
109
110 /*-------------------------------------------------------------------------
111 * Function: ph5diff_worker
112 *
113 * Purpose: worker process of ph5diff
114 *
115 * Return: none
116 *
117 * Programmer: Leon Arber
118 * Date: January 2005
119 *
120 * Comments:
121 *
122 * Modifications:
123 *
124 *-------------------------------------------------------------------------
125 */
126 static void
ph5diff_worker(int nID)127 ph5diff_worker(int nID)
128 {
129 hid_t file1_id = -1, file2_id = -1;
130
131 while(1)
132 {
133 MPI_Status Status;
134
135 MPI_Probe(0, MPI_ANY_TAG, MPI_COMM_WORLD, &Status);
136
137 /* Check for filenames */
138 if(Status.MPI_TAG == MPI_TAG_PARALLEL)
139 {
140 char filenames[2][MAX_FILENAME];
141
142 /* Retrieve filenames */
143 MPI_Recv(filenames, MAX_FILENAME*2, MPI_CHAR, 0, MPI_ANY_TAG, MPI_COMM_WORLD, &Status);
144
145 /* disable error reporting */
146 H5E_BEGIN_TRY
147 {
148 /* Open the files */
149 if ((file1_id = H5Fopen (filenames[0], H5F_ACC_RDONLY, H5P_DEFAULT)) < 0)
150 {
151 HDprintf("h5diff Task [%d]: <%s>: unable to open file\n", nID, filenames[0]);
152 MPI_Abort(MPI_COMM_WORLD, 0);
153 }
154 if ((file2_id = H5Fopen (filenames[1], H5F_ACC_RDONLY, H5P_DEFAULT)) < 0)
155 {
156 HDprintf("h5diff Task [%d]: <%s>: unable to open file\n", nID, filenames[1]);
157 MPI_Abort(MPI_COMM_WORLD, 0);
158 }
159 /* enable error reporting */
160 }
161 H5E_END_TRY;
162 }
163 /* Check for work */
164 else if(Status.MPI_TAG == MPI_TAG_ARGS)
165 {
166 struct diff_mpi_args args;
167 struct diffs_found diffs;
168 int i;
169
170 /* Make certain we've received the filenames and opened the files already */
171 if(file1_id < 0 || file2_id < 0)
172 {
173 HDprintf("ph5diff_worker: ERROR: work received before/without filenames\n");
174 break;
175 }
176
177 /* Recv parameters for diff from manager task */
178 MPI_Recv(&args, sizeof(args), MPI_BYTE, 0, MPI_TAG_ARGS, MPI_COMM_WORLD, &Status);
179
180 /* Do the diff */
181 diffs.nfound = diff(file1_id, args.name1, file2_id, args.name2, &(args.opts), &(args.argdata));
182 diffs.not_cmp = args.opts.not_cmp;
183
184 /* If print buffer has something in it, request print token.*/
185 if(outBuffOffset>0)
186 {
187 MPI_Send(NULL, 0, MPI_BYTE, 0, MPI_TAG_TOK_REQUEST, MPI_COMM_WORLD);
188
189 /* Wait for print token. */
190 MPI_Recv(NULL, 0, MPI_BYTE, 0, MPI_TAG_PRINT_TOK, MPI_COMM_WORLD, &Status);
191
192 /* When get token, send all of our output to the manager task and then return the token */
193 for(i=0; i<outBuffOffset; i+=PRINT_DATA_MAX_SIZE)
194 MPI_Send(outBuff+i, PRINT_DATA_MAX_SIZE, MPI_BYTE, 0, MPI_TAG_PRINT_DATA, MPI_COMM_WORLD);
195
196 /* An overflow file exists, so we send it's output to the manager too and then delete it */
197 if(overflow_file)
198 {
199 char out_data[PRINT_DATA_MAX_SIZE];
200 int tmp;
201
202 HDmemset(out_data, 0, PRINT_DATA_MAX_SIZE);
203 i=0;
204
205 rewind(overflow_file);
206 while((tmp = getc(overflow_file)) >= 0)
207 {
208 *(out_data + i++) = (char)tmp;
209 if(i==PRINT_DATA_MAX_SIZE)
210 {
211 MPI_Send(out_data, PRINT_DATA_MAX_SIZE, MPI_BYTE, 0, MPI_TAG_PRINT_DATA, MPI_COMM_WORLD);
212 i=0;
213 HDmemset(out_data, 0, PRINT_DATA_MAX_SIZE);
214 }
215 }
216
217 if(i>0)
218 MPI_Send(out_data, PRINT_DATA_MAX_SIZE, MPI_BYTE, 0, MPI_TAG_PRINT_DATA, MPI_COMM_WORLD);
219
220 fclose(overflow_file);
221 overflow_file = NULL;
222 }
223
224 HDfflush(stdout);
225 HDmemset(outBuff, 0, OUTBUFF_SIZE);
226 outBuffOffset = 0;
227
228 MPI_Send(&diffs, sizeof(diffs), MPI_BYTE, 0, MPI_TAG_TOK_RETURN, MPI_COMM_WORLD);
229 }
230 else
231 MPI_Send(&diffs, sizeof(diffs), MPI_BYTE, 0, MPI_TAG_DONE, MPI_COMM_WORLD);
232 }
233 /* Check for leaving */
234 else if(Status.MPI_TAG == MPI_TAG_END)
235 {
236 MPI_Recv(NULL, 0, MPI_BYTE, 0, MPI_TAG_END, MPI_COMM_WORLD, &Status);
237 break;
238 }
239 else
240 {
241 HDprintf("ph5diff_worker: ERROR: invalid tag (%d) received\n", Status.MPI_TAG);
242 break;
243 }
244
245 }
246
247 return;
248 }
249
250 /*-------------------------------------------------------------------------
251 * Function: print_manager_output
252 *
253 * Purpose: special function that prints any output accumulated by the
254 * manager task.
255 *
256 * Return: none
257 *
258 * Programmer: Leon Arber
259 *
260 * Date: Feb 7, 2005
261 *
262 *-------------------------------------------------------------------------
263 */
print_manager_output(void)264 void print_manager_output(void)
265 {
266 /* If there was something we buffered, let's print it now */
267 if( (outBuffOffset>0) && g_Parallel)
268 {
269 HDprintf("%s", outBuff);
270
271 if(overflow_file)
272 {
273 int tmp;
274 rewind(overflow_file);
275 while((tmp = HDgetc(overflow_file)) >= 0)
276 HDputchar(tmp);
277 fclose(overflow_file);
278 overflow_file = NULL;
279 }
280
281 HDfflush(stdout);
282 HDmemset(outBuff, 0, OUTBUFF_SIZE);
283 outBuffOffset = 0;
284 }
285 else if( (outBuffOffset>0) && !g_Parallel)
286 {
287 HDfprintf(stderr, "h5diff error: outBuffOffset>0, but we're not in parallel!\n");
288 }
289 }
290
291 /*-------------------------------------------------------------------------
292 * Function: h5diff_exit
293 *
294 * Purpose: dismiss phdiff worker processes and exit
295 *
296 * Return: none
297 *
298 * Programmer: Albert Cheng
299 * Date: Feb 6, 2005
300 *
301 * Comments:
302 *
303 * Modifications:
304 *
305 *-------------------------------------------------------------------------
306 */
h5diff_exit(int status)307 void h5diff_exit(int status)
308 {
309 /* if in parallel mode, dismiss workers, close down MPI, then exit */
310 if(g_Parallel) {
311 if(g_nTasks > 1) {
312 phdiff_dismiss_workers();
313 MPI_Barrier(MPI_COMM_WORLD);
314 }
315 MPI_Finalize();
316 status = EXIT_SUCCESS; /* Reset exit status, since some mpiexec commands generate output on failure status */
317 }
318
319 h5tools_close();
320
321 /* Always exit(0), since MPI implementations do weird stuff when they
322 * receive a non-zero exit value. - QAK
323 */
324 HDexit(0);
325 }
326
327