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 /*
15 * Programmer: Leon Arber <larber@uiuc.edu>
16 * Sept. 28, 2006.
17 *
18 * Purpose: This is the second half of a two-part test that makes sure
19 * that a file can be read after a parallel application crashes as long
20 * as the file was flushed first. We simulate a crash by
21 * calling _exit(0) since this doesn't flush HDF5 caches but
22 * still exits with success.
23 */
24
25 #include "h5test.h"
26
27 const char *FILENAME[] = {
28 "flush",
29 "noflush",
30 NULL
31 };
32
33 static double the_data[100][100];
34
35
36 /*-------------------------------------------------------------------------
37 * Function: check_file
38 *
39 * Purpose: Part 2 of a two-part H5Fflush() test.
40 *
41 * Return: Success: 0
42 *
43 * Failure: 1
44 *
45 * Programmer: Leon Arber
46 * Sept. 26, 2006.
47 *
48 *-------------------------------------------------------------------------
49 */
50 static int
check_file(char * name,hid_t fapl)51 check_file(char* name, hid_t fapl)
52 {
53 hid_t file, space, dset, groups, grp, plist;
54 hsize_t ds_size[2];
55 double error;
56 hsize_t i, j;
57
58 plist = H5Pcreate(H5P_DATASET_XFER);
59 H5Pset_dxpl_mpio(plist, H5FD_MPIO_COLLECTIVE);
60 if((file = H5Fopen(name, H5F_ACC_RDONLY, fapl)) < 0) goto error;
61
62 /* Open the dataset */
63 if((dset = H5Dopen2(file, "dset", H5P_DEFAULT)) < 0) goto error;
64 if((space = H5Dget_space(dset)) < 0) goto error;
65 if(H5Sget_simple_extent_dims(space, ds_size, NULL) < 0) goto error;
66 assert(100==ds_size[0] && 100==ds_size[1]);
67
68 /* Read some data */
69 if (H5Dread(dset, H5T_NATIVE_DOUBLE, space, space, plist,
70 the_data) < 0) goto error;
71 for (i=0; i<ds_size[0]; i++) {
72 for (j=0; j<ds_size[1]; j++) {
73 /*
74 * The extra cast in the following statement is a bug workaround
75 * for the Win32 version 5.0 compiler.
76 * 1998-11-06 ptl
77 */
78 error = fabs(the_data[i][j]-(double)(hssize_t)i/((hssize_t)j+1));
79 if (error>0.0001) {
80 H5_FAILED();
81 printf(" dset[%lu][%lu] = %g\n",
82 (unsigned long)i, (unsigned long)j, the_data[i][j]);
83 printf(" should be %g\n",
84 (double)(hssize_t)i/(hssize_t)(j+1));
85 goto error;
86 }
87 }
88 }
89
90 /* Open some groups */
91 if((groups = H5Gopen2(file, "some_groups", H5P_DEFAULT)) < 0) goto error;
92 for(i = 0; i < 100; i++) {
93 sprintf(name, "grp%02u", (unsigned)i);
94 if((grp = H5Gopen2(groups, name, H5P_DEFAULT)) < 0) goto error;
95 if(H5Gclose(grp) < 0) goto error;
96 }
97
98 if(H5Gclose(groups) < 0) goto error;
99 if(H5Dclose(dset) < 0) goto error;
100 if(H5Fclose(file) < 0) goto error;
101 if(H5Pclose(plist) < 0) goto error;
102 if(H5Sclose(space) < 0) goto error;
103
104 return 0;
105
106 error:
107 H5E_BEGIN_TRY {
108 H5Pclose(plist);
109 H5Gclose(groups);
110 H5Dclose(dset);
111 H5Fclose(file);
112 H5Sclose(space);
113 } H5E_END_TRY;
114 return 1;
115 }
116
117 /*-------------------------------------------------------------------------
118 * Function: main
119 *
120 * Purpose: Part 2 of a two-part H5Fflush() test.
121 *
122 * Return: Success: 0
123 *
124 * Failure: 1
125 *
126 * Programmer: Robb Matzke
127 * Friday, October 23, 1998
128 *
129 * Modifications:
130 * Leon Arber
131 * Sept. 26, 2006, expand to check for case where the was file not flushed.
132 *
133 *-------------------------------------------------------------------------
134 */
135 int
main(int argc,char * argv[])136 main(int argc, char* argv[])
137 {
138 hid_t fapl1, fapl2;
139 H5E_auto2_t func;
140
141 char name[1024];
142 const char *envval = NULL;
143
144 int mpi_size, mpi_rank;
145 MPI_Comm comm = MPI_COMM_WORLD;
146 MPI_Info info = MPI_INFO_NULL;
147
148 MPI_Init(&argc, &argv);
149 MPI_Comm_size(comm, &mpi_size);
150 MPI_Comm_rank(comm, &mpi_rank);
151
152 fapl1 = H5Pcreate(H5P_FILE_ACCESS);
153 H5Pset_fapl_mpio(fapl1, comm, info);
154
155 fapl2 = H5Pcreate(H5P_FILE_ACCESS);
156 H5Pset_fapl_mpio(fapl2, comm, info);
157
158
159 if(mpi_rank == 0)
160 TESTING("H5Fflush (part2 with flush)");
161
162 /* Don't run this test using the core or split file drivers */
163 envval = HDgetenv("HDF5_DRIVER");
164 if (envval == NULL)
165 envval = "nomatch";
166 if (HDstrcmp(envval, "core") && HDstrcmp(envval, "split")) {
167 /* Check the case where the file was flushed */
168 h5_fixname(FILENAME[0], fapl1, name, sizeof name);
169 if(check_file(name, fapl1))
170 {
171 H5_FAILED()
172 goto error;
173 }
174 else if(mpi_rank == 0)
175 {
176 PASSED()
177 }
178
179 /* Check the case where the file was not flushed. This should give an error
180 * so we turn off the error stack temporarily */
181 if(mpi_rank == 0)
182 TESTING("H5Fflush (part2 without flush)");
183 H5Eget_auto2(H5E_DEFAULT,&func,NULL);
184 H5Eset_auto2(H5E_DEFAULT, NULL, NULL);
185
186 h5_fixname(FILENAME[1], fapl2, name, sizeof name);
187 if(check_file(name, fapl2))
188 {
189 if(mpi_rank == 0)
190 {
191 PASSED()
192 }
193 }
194 else
195 {
196 H5_FAILED()
197 goto error;
198 }
199 H5Eset_auto2(H5E_DEFAULT, func, NULL);
200
201
202 h5_clean_files(&FILENAME[0], fapl1);
203 h5_clean_files(&FILENAME[1], fapl2);
204 }
205 else
206 {
207 SKIPPED();
208 puts(" Test not compatible with current Virtual File Driver");
209 }
210
211 MPI_Finalize();
212 return 0;
213
214 error:
215 return 1;
216 }
217
218
219
220
221