Actual source code: ex48.c
2: static char help[] = "Tests HDF5 attribute I/O.\n\n";
4: #include <petscviewerhdf5.h>
5: #include <petscvec.h>
7: static PetscInt n = 5; /* testing vector size */
8: static PetscBool verbose = PETSC_FALSE;
9: #define SLEN 128
11: /* sequence of unique absolute paths */
12: #define nap 9
13: static const char *apaths[nap] = {
14: /* 0 */
15: "/", "/g1", "/g1/g2", "/g1/nonExistingGroup1", "/g1/g3",
16: /* 5 */
17: "/g1/g3/g4", "/g1/nonExistingGroup2", "/g1/nonExistingGroup2/g5", "/g1/g6/g7"};
19: #define np 21
20: /* sequence of paths (absolute or relative); "<" encodes Pop */
21: static const char *paths[np] = {
22: /* 0 */
23: "/",
24: "/g1",
25: "/g1/g2",
26: "/g1/nonExistingGroup1",
27: "<",
28: /* 5 */
29: ".", /* /g1/g2 */
30: "<",
31: "<",
32: "g3", /* /g1/g3 */
33: "g4", /* /g1/g3/g4 */
34: /* 10 */
35: "<",
36: "<",
37: ".", /* /g1 */
38: "<",
39: "nonExistingGroup2", /* /g1/nonExistingG2 */
40: /* 15 */
41: "g5", /* /g1/nonExistingG2/g5 */
42: "<",
43: "<",
44: "g6/g7", /* /g1/g6/g7 */
45: "<",
46: /* 20 */
47: "<",
48: };
49: /* corresponding expected absolute paths - positions in abspath */
50: static const PetscInt paths2apaths[np] = {
51: /* 0 */
52: 0,
53: 1,
54: 2,
55: 3,
56: 2,
57: /* 5 */
58: 2,
59: 2,
60: 1,
61: 4,
62: 5,
63: /* 10 */
64: 4,
65: 1,
66: 1,
67: 1,
68: 6,
69: /* 15 */
70: 7,
71: 6,
72: 1,
73: 8,
74: 1,
75: /* 20 */
76: 0,
77: };
79: #define ns 4
80: /* for "" attribute will be stored to group, otherwise to given dataset */
81: static const char *datasets[ns] = {"", "x", "nonExistingVec", "y"};
83: /* beware this yields PETSC_FALSE for "" but group "" is interpreted as "/" */
84: static inline PetscErrorCode shouldExist(const char name[], PetscBool emptyExists, PetscBool *has)
85: {
86: size_t len = 0;
88: PetscFunctionBegin;
89: PetscCall(PetscStrlen(name, &len));
90: *has = emptyExists;
91: if (len) {
92: char *loc = NULL;
93: PetscCall(PetscStrstr(name, "nonExisting", &loc));
94: *has = PetscNot(loc);
95: }
96: PetscFunctionReturn(PETSC_SUCCESS);
97: }
99: static inline PetscErrorCode isPop(const char path[], PetscBool *has)
100: {
101: PetscFunctionBegin;
102: PetscCall(PetscStrcmp(path, "<", has));
103: PetscFunctionReturn(PETSC_SUCCESS);
104: }
106: static inline PetscErrorCode isDot(const char path[], PetscBool *has)
107: {
108: PetscFunctionBegin;
109: PetscCall(PetscStrcmp(path, ".", has));
110: PetscFunctionReturn(PETSC_SUCCESS);
111: }
113: static inline PetscErrorCode isRoot(const char path[], PetscBool *flg)
114: {
115: size_t len;
117: PetscFunctionBegin;
118: PetscCall(PetscStrlen(path, &len));
119: *flg = PetscNot(len);
120: if (!*flg) PetscCall(PetscStrcmp(path, "/", flg));
121: PetscFunctionReturn(PETSC_SUCCESS);
122: }
124: static inline PetscErrorCode compare(PetscDataType dt, void *ptr0, void *ptr1, PetscBool *flg)
125: {
126: PetscFunctionBegin;
127: switch (dt) {
128: case PETSC_INT:
129: *flg = (PetscBool)(*(PetscInt *)ptr0 == *(PetscInt *)ptr1);
130: if (verbose) {
131: if (*flg) {
132: PetscCall(PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT, *(PetscInt *)ptr0));
133: } else {
134: PetscCall(PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT " != %" PetscInt_FMT "\n", *(PetscInt *)ptr0, *(PetscInt *)ptr1));
135: }
136: }
137: break;
138: case PETSC_REAL:
139: *flg = (PetscBool)(*(PetscReal *)ptr0 == *(PetscReal *)ptr1);
140: if (verbose) {
141: if (*flg) {
142: PetscCall(PetscPrintf(PETSC_COMM_SELF, "%f", *(PetscReal *)ptr0));
143: } else {
144: PetscCall(PetscPrintf(PETSC_COMM_SELF, "%f != %f\n", *(PetscReal *)ptr0, *(PetscReal *)ptr1));
145: }
146: }
147: break;
148: case PETSC_BOOL:
149: *flg = (PetscBool)(*(PetscBool *)ptr0 == *(PetscBool *)ptr1);
150: if (verbose) {
151: if (*flg) {
152: PetscCall(PetscPrintf(PETSC_COMM_SELF, "%s", PetscBools[*(PetscBool *)ptr0]));
153: } else {
154: PetscCall(PetscPrintf(PETSC_COMM_SELF, "%s != %s\n", PetscBools[*(PetscBool *)ptr0], PetscBools[*(PetscBool *)ptr1]));
155: }
156: }
157: break;
158: case PETSC_STRING:
159: PetscCall(PetscStrcmp((const char *)ptr0, (const char *)ptr1, flg));
160: if (verbose) {
161: if (*flg) {
162: PetscCall(PetscPrintf(PETSC_COMM_SELF, "%s", (char *)ptr0));
163: } else {
164: PetscCall(PetscPrintf(PETSC_COMM_SELF, "%s != %s\n", (char *)ptr0, (char *)ptr1));
165: }
166: }
167: break;
168: default:
169: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "PetscDataType %s not handled here", PetscDataTypes[dt]);
170: }
171: PetscFunctionReturn(PETSC_SUCCESS);
172: }
174: static inline PetscErrorCode alterString(const char oldstr[], char str[])
175: {
176: size_t i, n;
178: PetscFunctionBegin;
179: PetscCall(PetscStrlen(oldstr, &n));
180: PetscCall(PetscStrncpy(str, oldstr, n + 1));
181: for (i = 0; i < n; i++) {
182: if (('A' <= str[i] && str[i] < 'Z') || ('a' <= str[i] && str[i] < 'z')) {
183: str[i]++;
184: break;
185: }
186: }
187: PetscFunctionReturn(PETSC_SUCCESS);
188: }
190: /* if name given, check dataset with this name exists under current group, otherwise just check current group exists */
191: /* flg: 0 doesn't exist, 1 group, 2 dataset */
192: static PetscErrorCode hasGroupOrDataset(PetscViewer viewer, const char path[], int *flg)
193: {
194: PetscBool has;
196: PetscFunctionBegin;
197: *flg = 0;
198: PetscCall(PetscViewerHDF5HasGroup(viewer, path, &has));
199: if (has) *flg = 1;
200: else {
201: PetscCall(PetscViewerHDF5HasDataset(viewer, path, &has));
202: if (has) *flg = 2;
203: }
204: PetscFunctionReturn(PETSC_SUCCESS);
205: }
207: #define nt 5 /* number of datatypes */
208: typedef struct _n_Capsule *Capsule;
209: struct _n_Capsule {
210: char names[nt][SLEN];
211: PetscDataType types[nt];
212: char typeNames[nt][SLEN];
213: size_t sizes[nt];
214: void *vals[nt];
215: PetscInt id, ntypes;
216: };
218: static PetscErrorCode CapsuleCreate(Capsule old, Capsule *newcapsule)
219: {
220: Capsule c;
221: PetscBool bool0 = PETSC_TRUE;
222: PetscInt int0 = -1;
223: PetscReal real0 = -1.1;
224: char str0[] = "Test String";
225: char nestr0[] = "NONEXISTING STRING"; /* this attribute shall be skipped for writing */
226: void *vals[nt] = {&bool0, &int0, &real0, str0, nestr0};
227: size_t sizes[nt] = {sizeof(bool0), sizeof(int0), sizeof(real0), sizeof(str0), sizeof(str0)};
228: PetscDataType types[nt] = {PETSC_BOOL, PETSC_INT, PETSC_REAL, PETSC_STRING, PETSC_STRING};
229: const char *tNames[nt] = {"bool", "int", "real", "str", "nonExisting"};
230: PetscInt t;
232: PetscFunctionBegin;
233: PetscCall(PetscNew(&c));
234: c->id = 0;
235: c->ntypes = nt;
236: if (old) {
237: /* alter values */
238: t = 0;
239: bool0 = PetscNot(*((PetscBool *)old->vals[t]));
240: t++;
241: int0 = *((PetscInt *)old->vals[t]) * -2;
242: t++;
243: real0 = *((PetscReal *)old->vals[t]) * -2.0;
244: t++;
245: PetscCall(alterString((const char *)old->vals[t], str0));
246: t++;
247: c->id = old->id + 1;
248: }
249: for (t = 0; t < nt; t++) {
250: c->sizes[t] = sizes[t];
251: c->types[t] = types[t];
252: PetscCall(PetscStrncpy(c->typeNames[t], tNames[t], sizeof(c->typeNames[t])));
253: PetscCall(PetscSNPrintf(c->names[t], SLEN, "attr_%" PetscInt_FMT "_%s", c->id, tNames[t]));
254: PetscCall(PetscMalloc(sizes[t], &c->vals[t]));
255: PetscCall(PetscMemcpy(c->vals[t], vals[t], sizes[t]));
256: }
257: *newcapsule = c;
258: PetscFunctionReturn(PETSC_SUCCESS);
259: }
260: #undef nt
262: static PetscErrorCode CapsuleWriteAttributes(Capsule c, PetscViewer v, const char parent[])
263: {
264: PetscInt t;
265: PetscBool flg = PETSC_FALSE;
267: PetscFunctionBegin;
268: for (t = 0; t < c->ntypes; t++) {
269: PetscCall(shouldExist(c->names[t], PETSC_FALSE, &flg));
270: if (!flg) continue;
271: PetscCall(PetscViewerHDF5WriteAttribute(v, parent, c->names[t], c->types[t], c->vals[t]));
272: }
273: PetscFunctionReturn(PETSC_SUCCESS);
274: }
276: static PetscErrorCode CapsuleReadAndCompareAttributes(Capsule c, PetscViewer v, const char parent[])
277: {
278: char *group;
279: int gd = 0;
280: PetscInt t;
281: PetscBool flg = PETSC_FALSE, hasAttr = PETSC_FALSE;
282: MPI_Comm comm;
284: PetscFunctionBegin;
285: PetscCall(PetscObjectGetComm((PetscObject)v, &comm));
286: PetscCall(PetscViewerHDF5GetGroup(v, NULL, &group));
287: PetscCall(hasGroupOrDataset(v, parent, &gd));
288: /* check correct existence of attributes */
289: for (t = 0; t < c->ntypes; t++) {
290: const char *attribute = c->names[t];
291: PetscCall(shouldExist(attribute, PETSC_FALSE, &flg));
292: PetscCall(PetscViewerHDF5HasAttribute(v, parent, attribute, &hasAttr));
293: if (verbose) {
294: PetscCall(PetscPrintf(comm, " %-24s = ", attribute));
295: if (!hasAttr) PetscCall(PetscPrintf(comm, "---"));
296: }
297: PetscCheck(gd || !hasAttr, comm, PETSC_ERR_PLIB, "Attribute %s/%s/%s exists while its parent %s/%s doesn't exist", group, parent, attribute, group, parent);
298: PetscCheck(flg == hasAttr, comm, PETSC_ERR_PLIB, "Attribute %s/%s should exist? %s Exists? %s", parent, attribute, PetscBools[flg], PetscBools[hasAttr]);
300: /* check loaded attributes are the same as original */
301: if (hasAttr) {
302: char buffer[SLEN];
303: char *str;
304: void *ptr0;
305: /* check the stored data is the same as original */
306: //TODO datatype should better be output arg, not input
307: //TODO string attributes should probably have a separate function since the handling is different;
308: //TODO or maybe it should just accept string buffer rather than pointer to string
309: if (c->types[t] == PETSC_STRING) {
310: PetscCall(PetscViewerHDF5ReadAttribute(v, parent, attribute, c->types[t], NULL, &str));
311: ptr0 = str;
312: } else {
313: PetscCall(PetscViewerHDF5ReadAttribute(v, parent, attribute, c->types[t], NULL, &buffer));
314: ptr0 = &buffer;
315: }
316: PetscCall(compare(c->types[t], ptr0, c->vals[t], &flg));
317: PetscCheck(flg, comm, PETSC_ERR_PLIB, "Value of attribute %s/%s/%s is not equal to the original value", group, parent, attribute);
318: if (verbose) PetscCall(PetscPrintf(comm, " (=)"));
319: if (c->types[t] == PETSC_STRING) PetscCall(PetscFree(str));
320: }
321: if (verbose && gd) PetscCall(PetscPrintf(comm, "\n"));
322: }
323: PetscCall(PetscFree(group));
324: PetscFunctionReturn(PETSC_SUCCESS);
325: }
327: static PetscErrorCode CapsuleDestroy(Capsule *c)
328: {
329: PetscInt t;
331: PetscFunctionBegin;
332: if (!*c) PetscFunctionReturn(PETSC_SUCCESS);
333: for (t = 0; t < (*c)->ntypes; t++) PetscCall(PetscFree((*c)->vals[t]));
334: PetscCall(PetscFree(*c));
335: PetscFunctionReturn(PETSC_SUCCESS);
336: }
338: static PetscErrorCode testGroupsDatasets(PetscViewer viewer)
339: {
340: char buf[PETSC_MAX_PATH_LEN];
341: Vec vecs[nap][ns];
342: PetscInt p, s;
343: PetscBool flg = PETSC_FALSE, flg1 = PETSC_FALSE, flg2 = PETSC_FALSE;
344: PetscRandom rand;
345: const char *filename;
346: MPI_Comm comm;
348: PetscFunctionBegin;
349: PetscCall(PetscObjectGetComm((PetscObject)viewer, &comm));
350: PetscCall(PetscViewerFileGetName(viewer, &filename));
351: if (verbose) PetscCall(PetscPrintf(comm, "# TEST testGroupsDatasets\n"));
352: /* store random vectors */
353: PetscCall(PetscRandomCreate(comm, &rand));
354: PetscCall(PetscRandomSetInterval(rand, 0.0, 10.0));
355: PetscCall(PetscRandomSetFromOptions(rand));
356: PetscCall(PetscMemzero(vecs, nap * ns * sizeof(Vec)));
358: /* test dataset writing */
359: if (verbose) PetscCall(PetscPrintf(comm, "## WRITE PHASE\n"));
360: for (p = 0; p < np; p++) {
361: PetscCall(isPop(paths[p], &flg));
362: PetscCall(isDot(paths[p], &flg1));
363: PetscCall(shouldExist(apaths[paths2apaths[p]], PETSC_FALSE, &flg2));
364: if (flg) {
365: PetscCall(PetscViewerHDF5PopGroup(viewer));
366: } else {
367: PetscCall(PetscViewerHDF5PushGroup(viewer, paths[p]));
368: }
369: if (verbose) PetscCall(PetscPrintf(comm, "%-32s => %4s => %-32s should exist? %s\n", paths[p], flg ? "pop" : "push", apaths[paths2apaths[p]], PetscBools[flg2]));
370: if (flg || flg1 || !flg2) continue;
372: for (s = 0; s < ns; s++) {
373: Vec v;
375: PetscCall(shouldExist(datasets[s], PETSC_FALSE, &flg));
376: if (!flg) continue;
378: PetscCall(VecCreate(comm, &v));
379: PetscCall(PetscObjectSetName((PetscObject)v, datasets[s]));
380: PetscCall(VecSetSizes(v, n, PETSC_DECIDE));
381: PetscCall(VecSetFromOptions(v));
382: PetscCall(VecSetRandom(v, rand));
383: if (verbose) {
384: PetscReal min, max;
385: PetscCall(VecMin(v, NULL, &min));
386: PetscCall(VecMax(v, NULL, &max));
387: PetscCall(PetscPrintf(comm, " Create dataset %s/%s, keep in memory in vecs[%" PetscInt_FMT "][%" PetscInt_FMT "], min %.3e max %.3e\n", apaths[paths2apaths[p]], datasets[s], paths2apaths[p], s, min, max));
388: }
390: PetscCall(VecView(v, viewer));
391: vecs[paths2apaths[p]][s] = v;
392: }
393: }
394: PetscCall(PetscViewerFlush(viewer));
395: PetscCall(PetscRandomDestroy(&rand));
397: if (verbose) PetscCall(PetscPrintf(comm, "\n## READ PHASE\n"));
398: /* check correct existence of groups in file */
399: for (p = 0; p < np; p++) {
400: char *group;
401: const char *expected = apaths[paths2apaths[p]];
403: /* check Push/Pop is correct */
404: PetscCall(isPop(paths[p], &flg));
405: if (flg) {
406: PetscCall(PetscViewerHDF5PopGroup(viewer));
407: } else {
408: PetscCall(PetscViewerHDF5PushGroup(viewer, paths[p]));
409: }
410: PetscCall(PetscViewerHDF5GetGroup(viewer, NULL, &group));
411: PetscCall(PetscViewerHDF5HasGroup(viewer, NULL, &flg1));
412: if (verbose) PetscCall(PetscPrintf(comm, "%-32s => %4s => %-32s exists? %s\n", paths[p], flg ? "pop" : "push", group, PetscBools[flg1]));
413: PetscCall(PetscStrcmp(group, expected, &flg2));
414: PetscCheck(flg2, comm, PETSC_ERR_PLIB, "Current group %s not equal to expected %s", group, expected);
415: PetscCall(shouldExist(group, PETSC_TRUE, &flg2));
416: PetscCheck(flg1 == flg2, comm, PETSC_ERR_PLIB, "Group %s should exist? %s Exists in %s? %s", group, PetscBools[flg2], filename, PetscBools[flg1]);
417: PetscCall(PetscFree(group));
418: }
420: /* check existence of datasets; compare loaded vectors with original ones */
421: for (p = 0; p < np; p++) {
422: char *group;
424: /* check Push/Pop is correct */
425: PetscCall(isPop(paths[p], &flg));
426: if (flg) {
427: PetscCall(PetscViewerHDF5PopGroup(viewer));
428: } else {
429: PetscCall(PetscViewerHDF5PushGroup(viewer, paths[p]));
430: }
431: PetscCall(PetscViewerHDF5GetGroup(viewer, NULL, &group));
432: PetscCall(PetscViewerHDF5HasGroup(viewer, NULL, &flg));
433: if (verbose) PetscCall(PetscPrintf(comm, "Has %s group? %s\n", group, PetscBools[flg]));
434: for (s = 0; s < ns; s++) {
435: const char *name = datasets[s];
436: char *fullname = buf;
438: /* check correct existence of datasets in file */
439: PetscCall(PetscSNPrintf(fullname, sizeof(buf), "%s/%s", group, name));
440: PetscCall(shouldExist(name, PETSC_FALSE, &flg1));
441: flg1 = (PetscBool)(flg && flg1); /* both group and dataset need to exist */
442: PetscCall(PetscViewerHDF5HasDataset(viewer, name, &flg2));
443: if (verbose) PetscCall(PetscPrintf(comm, " %s dataset? %s", fullname, PetscBools[flg2]));
444: PetscCheck(flg2 == flg1, comm, PETSC_ERR_PLIB, "Dataset %s should exist? %s Exists in %s? %s", fullname, PetscBools[flg1], filename, PetscBools[flg2]);
446: if (flg2) {
447: Vec v;
448: /* check loaded Vec is the same as original */
449: PetscCall(VecCreate(comm, &v));
450: PetscCall(PetscObjectSetName((PetscObject)v, name));
451: PetscCall(VecLoad(v, viewer));
452: PetscCall(VecEqual(v, vecs[paths2apaths[p]][s], &flg1));
453: PetscCheck(flg1, comm, PETSC_ERR_PLIB, "Dataset %s in %s is not equal to the original Vec", fullname, filename);
454: if (verbose) PetscCall(PetscPrintf(comm, " (=)"));
455: PetscCall(VecDestroy(&v));
456: }
457: if (verbose) PetscCall(PetscPrintf(comm, "\n"));
458: }
459: PetscCall(PetscFree(group));
460: }
461: PetscCall(PetscViewerFlush(viewer));
462: for (p = 0; p < nap; p++)
463: for (s = 0; s < ns; s++) PetscCall(VecDestroy(&vecs[p][s]));
464: if (verbose) PetscCall(PetscPrintf(comm, "# END testGroupsDatasets\n\n"));
465: PetscFunctionReturn(PETSC_SUCCESS);
466: }
468: static inline PetscErrorCode formPath(PetscBool relativize, const char path[], const char dataset[], char buf[], size_t bufsize)
469: {
470: PetscBool isroot = PETSC_FALSE;
472: PetscFunctionBegin;
473: PetscCall(isRoot(path, &isroot));
474: if (relativize) {
475: if (isroot) {
476: PetscCall(PetscStrncpy(buf, dataset, bufsize));
477: } else {
478: /* skip initial '/' in paths[p] if prefix given */
479: PetscCall(PetscSNPrintf(buf, bufsize, "%s/%s", path + 1, dataset));
480: }
481: } else {
482: PetscCall(PetscSNPrintf(buf, bufsize, "%s/%s", isroot ? "" : path, dataset));
483: }
484: PetscFunctionReturn(PETSC_SUCCESS);
485: }
487: /* test attribute writing, existence checking and reading, use absolute paths */
488: static PetscErrorCode testAttributesAbsolutePath(PetscViewer viewer, const char prefix[])
489: {
490: char buf[PETSC_MAX_PATH_LEN];
491: Capsule capsules[nap][ns], c = NULL, old = NULL;
492: PetscInt p, s;
493: PetscBool flg = PETSC_FALSE, flg1 = PETSC_FALSE;
494: MPI_Comm comm;
496: PetscFunctionBegin;
497: PetscCall(PetscObjectGetComm((PetscObject)viewer, &comm));
498: if (verbose) {
499: if (prefix) {
500: PetscCall(PetscPrintf(comm, "# TEST testAttributesAbsolutePath, prefix=\"%s\"\n", prefix));
501: } else {
502: PetscCall(PetscPrintf(comm, "# TEST testAttributesAbsolutePath\n"));
503: }
504: PetscCall(PetscPrintf(comm, "## WRITE PHASE\n"));
505: }
506: PetscCall(PetscMemzero(capsules, nap * ns * sizeof(Capsule)));
508: /* test attribute writing */
509: if (prefix) PetscCall(PetscViewerHDF5PushGroup(viewer, prefix));
510: for (p = 0; p < np; p++)
511: for (s = 0; s < ns; s++) {
512: /* we test only absolute paths here */
513: PetscCall(PetscViewerHDF5PathIsRelative(paths[p], PETSC_FALSE, &flg));
514: if (flg) continue;
515: {
516: char *group;
518: PetscCall(PetscViewerHDF5GetGroup(viewer, NULL, &group));
519: PetscCall(PetscStrcmp(group, prefix, &flg));
520: PetscCheck(flg, comm, PETSC_ERR_PLIB, "prefix %s not equal to pushed group %s", prefix, group);
521: PetscCall(PetscFree(group));
522: }
523: PetscCall(formPath((PetscBool) !!prefix, paths[p], datasets[s], buf, sizeof(buf)));
524: PetscCall(shouldExist(buf, PETSC_TRUE, &flg));
525: if (!flg) continue;
527: if (verbose) {
528: if (prefix) {
529: PetscCall(PetscPrintf(comm, "Write attributes to %s/%s\n", prefix, buf));
530: } else {
531: PetscCall(PetscPrintf(comm, "Write attributes to %s\n", buf));
532: }
533: }
535: PetscCall(CapsuleCreate(old, &c));
536: PetscCall(CapsuleWriteAttributes(c, viewer, buf));
537: PetscCheck(!capsules[paths2apaths[p]][s], comm, PETSC_ERR_PLIB, "capsules[%" PetscInt_FMT "][%" PetscInt_FMT "] gets overwritten for %s", paths2apaths[p], s, buf);
538: capsules[paths2apaths[p]][s] = c;
539: old = c;
540: }
541: if (prefix) PetscCall(PetscViewerHDF5PopGroup(viewer));
542: PetscCall(PetscViewerFlush(viewer));
544: if (verbose) PetscCall(PetscPrintf(comm, "\n## READ PHASE\n"));
545: if (prefix) PetscCall(PetscViewerHDF5PushGroup(viewer, prefix));
546: for (p = 0; p < np; p++)
547: for (s = 0; s < ns; s++) {
548: /* we test only absolute paths here */
549: PetscCall(PetscViewerHDF5PathIsRelative(paths[p], PETSC_FALSE, &flg));
550: if (flg) continue;
552: /* check existence of given group/dataset */
553: PetscCall(formPath((PetscBool) !!prefix, paths[p], datasets[s], buf, sizeof(buf)));
554: PetscCall(shouldExist(buf, PETSC_TRUE, &flg));
555: if (verbose) {
556: if (prefix) {
557: PetscCall(PetscPrintf(comm, "Has %s/%s? %s\n", prefix, buf, PetscBools[flg]));
558: } else {
559: PetscCall(PetscPrintf(comm, "Has %s? %s\n", buf, PetscBools[flg]));
560: }
561: }
563: /* check attribute capsule has been created for given path */
564: c = capsules[paths2apaths[p]][s];
565: flg1 = (PetscBool) !!c;
566: PetscCheck(flg == flg1, comm, PETSC_ERR_PLIB, "Capsule should exist for %s? %s Exists? %s", buf, PetscBools[flg], PetscBools[flg1]);
567: if (!flg) continue;
569: /* check correct existence and fidelity of attributes in file */
570: PetscCall(CapsuleReadAndCompareAttributes(c, viewer, buf));
571: }
572: if (prefix) PetscCall(PetscViewerHDF5PopGroup(viewer));
573: PetscCall(PetscViewerFlush(viewer));
574: for (p = 0; p < nap; p++)
575: for (s = 0; s < ns; s++) PetscCall(CapsuleDestroy(&capsules[p][s]));
576: if (verbose) PetscCall(PetscPrintf(comm, "# END testAttributesAbsolutePath\n\n"));
577: PetscFunctionReturn(PETSC_SUCCESS);
578: }
580: /* test attribute writing, existence checking and reading, use group push/pop */
581: static PetscErrorCode testAttributesPushedPath(PetscViewer viewer)
582: {
583: Capsule capsules[nap][ns], c = NULL, old = NULL;
584: PetscInt p, s;
585: int gd;
586: PetscBool flg = PETSC_FALSE, flg1 = PETSC_FALSE;
587: MPI_Comm comm;
589: PetscFunctionBegin;
590: PetscCall(PetscObjectGetComm((PetscObject)viewer, &comm));
591: if (verbose) {
592: PetscCall(PetscPrintf(comm, "# TEST testAttributesPushedPath\n"));
593: PetscCall(PetscPrintf(comm, "## WRITE PHASE\n"));
594: }
595: PetscCall(PetscMemzero(capsules, nap * ns * sizeof(Capsule)));
597: /* test attribute writing */
598: for (p = 0; p < np; p++) {
599: PetscCall(isPop(paths[p], &flg));
600: PetscCall(isDot(paths[p], &flg1));
601: if (flg) {
602: PetscCall(PetscViewerHDF5PopGroup(viewer));
603: } else {
604: PetscCall(PetscViewerHDF5PushGroup(viewer, paths[p]));
605: }
606: /* < and . have been already visited => skip */
607: if (flg || flg1) continue;
609: /* assume here that groups and datasets are already in the file */
610: for (s = 0; s < ns; s++) {
611: PetscCall(hasGroupOrDataset(viewer, datasets[s], &gd));
612: if (!gd) continue;
613: if (verbose) PetscCall(PetscPrintf(comm, "Write attributes to %s/%s\n", apaths[paths2apaths[p]], datasets[s]));
614: PetscCall(CapsuleCreate(old, &c));
615: PetscCall(CapsuleWriteAttributes(c, viewer, datasets[s]));
616: PetscCheck(!capsules[paths2apaths[p]][s], comm, PETSC_ERR_PLIB, "capsules[%" PetscInt_FMT "][%" PetscInt_FMT "] gets overwritten for %s/%s", paths2apaths[p], s, paths[p], datasets[s]);
617: capsules[paths2apaths[p]][s] = c;
618: old = c;
619: }
620: }
621: PetscCall(PetscViewerFlush(viewer));
623: if (verbose) PetscCall(PetscPrintf(comm, "\n## READ PHASE\n"));
624: for (p = 0; p < np; p++) {
625: char *group;
627: PetscCall(isPop(paths[p], &flg1));
628: if (flg1) {
629: PetscCall(PetscViewerHDF5PopGroup(viewer));
630: } else {
631: PetscCall(PetscViewerHDF5PushGroup(viewer, paths[p]));
632: }
633: PetscCall(PetscViewerHDF5GetGroup(viewer, NULL, &group));
634: for (s = 0; s < ns; s++) {
635: PetscCall(hasGroupOrDataset(viewer, datasets[s], &gd));
636: if (verbose) PetscCall(PetscPrintf(comm, "%s/%s %s\n", group, datasets[s], gd ? (gd == 1 ? "is group" : "is dataset") : "does not exist"));
638: /* check attribute capsule has been created for given path */
639: c = capsules[paths2apaths[p]][s];
640: flg = (PetscBool) !!gd;
641: flg1 = (PetscBool) !!c;
642: PetscCheck(flg == flg1, comm, PETSC_ERR_PLIB, "Capsule should exist for %s/%s? %s Exists? %s", group, datasets[s], PetscBools[flg], PetscBools[flg1]);
643: if (!flg) continue;
645: /* check correct existence of attributes in file */
646: PetscCall(CapsuleReadAndCompareAttributes(c, viewer, datasets[s]));
647: }
648: PetscCall(PetscFree(group));
649: }
650: PetscCall(PetscViewerFlush(viewer));
651: for (p = 0; p < nap; p++)
652: for (s = 0; s < ns; s++) PetscCall(CapsuleDestroy(&capsules[p][s]));
653: if (verbose) PetscCall(PetscPrintf(comm, "# END testAttributesPushedPath\n\n"));
654: PetscFunctionReturn(PETSC_SUCCESS);
655: }
657: /* test attribute writing, existence checking and reading, use group push/pop */
658: static PetscErrorCode testObjectAttributes(PetscViewer viewer)
659: {
660: Capsule capsules[nap][ns], c = NULL, old = NULL;
661: PetscInt p, s;
662: PetscBool flg = PETSC_FALSE, flg1 = PETSC_FALSE;
663: MPI_Comm comm;
665: PetscFunctionBegin;
666: PetscCall(PetscObjectGetComm((PetscObject)viewer, &comm));
667: if (verbose) {
668: PetscCall(PetscPrintf(comm, "# TEST testObjectAttributes\n"));
669: PetscCall(PetscPrintf(comm, "## WRITE PHASE\n"));
670: }
671: PetscCall(PetscMemzero(capsules, nap * ns * sizeof(Capsule)));
673: /* test attribute writing */
674: for (p = 0; p < np; p++) {
675: PetscCall(isPop(paths[p], &flg));
676: PetscCall(isDot(paths[p], &flg1));
677: if (flg) {
678: PetscCall(PetscViewerHDF5PopGroup(viewer));
679: } else {
680: PetscCall(PetscViewerHDF5PushGroup(viewer, paths[p]));
681: }
682: /* < and . have been already visited => skip */
683: if (flg || flg1) continue;
685: /* assume here that groups and datasets are already in the file */
686: for (s = 0; s < ns; s++) {
687: Vec v;
688: size_t len;
689: const char *name = datasets[s];
691: PetscCall(PetscStrlen(name, &len));
692: if (!len) continue;
693: PetscCall(VecCreate(comm, &v));
694: PetscCall(PetscObjectSetName((PetscObject)v, name));
695: PetscCall(PetscViewerHDF5HasObject(viewer, (PetscObject)v, &flg));
696: if (flg) {
697: if (verbose) PetscCall(PetscPrintf(comm, "Write attributes to %s/%s\n", apaths[paths2apaths[p]], name));
698: PetscCall(CapsuleCreate(old, &c));
699: PetscCall(CapsuleWriteAttributes(c, viewer, name));
700: PetscCheck(!capsules[paths2apaths[p]][s], comm, PETSC_ERR_PLIB, "capsules[%" PetscInt_FMT "][%" PetscInt_FMT "] gets overwritten for %s/%s", paths2apaths[p], s, paths[p], name);
701: capsules[paths2apaths[p]][s] = c;
702: old = c;
703: }
704: PetscCall(VecDestroy(&v));
705: }
706: }
707: PetscCall(PetscViewerFlush(viewer));
709: if (verbose) PetscCall(PetscPrintf(comm, "\n## READ PHASE\n"));
710: for (p = 0; p < np; p++) {
711: char *group;
713: PetscCall(isPop(paths[p], &flg));
714: if (flg) {
715: PetscCall(PetscViewerHDF5PopGroup(viewer));
716: } else {
717: PetscCall(PetscViewerHDF5PushGroup(viewer, paths[p]));
718: }
719: PetscCall(PetscViewerHDF5GetGroup(viewer, NULL, &group));
720: for (s = 0; s < ns; s++) {
721: Vec v;
722: size_t len;
723: const char *name = datasets[s];
725: PetscCall(PetscStrlen(name, &len));
726: if (!len) continue;
727: PetscCall(VecCreate(comm, &v));
728: PetscCall(PetscObjectSetName((PetscObject)v, name));
729: PetscCall(PetscViewerHDF5HasObject(viewer, (PetscObject)v, &flg));
730: if (verbose) PetscCall(PetscPrintf(comm, "Is %s/%s dataset? %s\n", group, name, PetscBools[flg]));
732: /* check attribute capsule has been created for given path */
733: c = capsules[paths2apaths[p]][s];
734: flg1 = (PetscBool) !!c;
735: PetscCheck(flg == flg1, comm, PETSC_ERR_PLIB, "Capsule should exist for %s/%s? %s Exists? %s", group, name, PetscBools[flg], PetscBools[flg1]);
737: /* check correct existence of attributes in file */
738: if (flg) PetscCall(CapsuleReadAndCompareAttributes(c, viewer, name));
739: PetscCall(VecDestroy(&v));
740: }
741: PetscCall(PetscFree(group));
742: }
743: PetscCall(PetscViewerFlush(viewer));
744: for (p = 0; p < nap; p++)
745: for (s = 0; s < ns; s++) PetscCall(CapsuleDestroy(&capsules[p][s]));
746: if (verbose) PetscCall(PetscPrintf(comm, "# END testObjectAttributes\n\n"));
747: PetscFunctionReturn(PETSC_SUCCESS);
748: }
750: static PetscErrorCode testAttributesDefaultValue(PetscViewer viewer)
751: {
752: #define nv 4
753: PetscBool bools[nv];
754: PetscInt ints[nv];
755: PetscReal reals[nv];
756: char *strings[nv];
757: PetscBool flg;
758: PetscInt i;
759: MPI_Comm comm;
761: PetscFunctionBegin;
762: PetscCall(PetscObjectGetComm((PetscObject)viewer, &comm));
763: if (verbose) PetscCall(PetscPrintf(comm, "# TEST testAttributesDefaultValue\n"));
765: PetscCall(PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_bool", PETSC_BOOL, NULL, &bools[0]));
766: bools[1] = PetscNot(bools[0]);
767: PetscCall(PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_bool", PETSC_BOOL, &bools[1], &bools[2]));
768: PetscCall(PetscViewerHDF5ReadAttribute(viewer, "/", "attr_nonExisting_bool", PETSC_BOOL, &bools[1], &bools[3]));
769: PetscCheck(bools[2] == bools[0], comm, PETSC_ERR_PLIB, "%s = bools[2] != bools[0] = %s", PetscBools[bools[2]], PetscBools[bools[0]]);
770: PetscCheck(bools[3] == bools[1], comm, PETSC_ERR_PLIB, "%s = bools[3] != bools[1] = %s", PetscBools[bools[3]], PetscBools[bools[1]]);
772: PetscCall(PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_int", PETSC_INT, NULL, &ints[0]));
773: ints[1] = ints[0] * -333;
774: PetscCall(PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_int", PETSC_INT, &ints[1], &ints[2]));
775: PetscCall(PetscViewerHDF5ReadAttribute(viewer, "/", "attr_nonExisting_int", PETSC_INT, &ints[1], &ints[3]));
776: PetscCheck(ints[2] == ints[0], comm, PETSC_ERR_PLIB, "%" PetscInt_FMT " = ints[2] != ints[0] = %" PetscInt_FMT, ints[2], ints[0]);
777: PetscCheck(ints[3] == ints[1], comm, PETSC_ERR_PLIB, "%" PetscInt_FMT " = ints[3] != ints[1] = %" PetscInt_FMT, ints[3], ints[1]);
778: if (verbose) PetscCall(PetscIntView(nv, ints, PETSC_VIEWER_STDOUT_WORLD));
780: PetscCall(PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_real", PETSC_REAL, NULL, &reals[0]));
781: reals[1] = reals[0] * -11.1;
782: PetscCall(PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_real", PETSC_REAL, &reals[1], &reals[2]));
783: PetscCall(PetscViewerHDF5ReadAttribute(viewer, "/", "attr_nonExisting_real", PETSC_REAL, &reals[1], &reals[3]));
784: PetscCheck(reals[2] == reals[0], comm, PETSC_ERR_PLIB, "%f = reals[2] != reals[0] = %f", reals[2], reals[0]);
785: PetscCheck(reals[3] == reals[1], comm, PETSC_ERR_PLIB, "%f = reals[3] != reals[1] = %f", reals[3], reals[1]);
786: if (verbose) PetscCall(PetscRealView(nv, reals, PETSC_VIEWER_STDOUT_WORLD));
788: PetscCall(PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_str", PETSC_STRING, NULL, &strings[0]));
789: PetscCall(PetscStrallocpy(strings[0], &strings[1]));
790: PetscCall(alterString(strings[0], strings[1]));
791: PetscCall(PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_str", PETSC_STRING, &strings[1], &strings[2]));
792: PetscCall(PetscViewerHDF5ReadAttribute(viewer, "/", "attr_nonExisting_str", PETSC_STRING, &strings[1], &strings[3]));
793: PetscCall(PetscStrcmp(strings[2], strings[0], &flg));
794: PetscCheck(flg, comm, PETSC_ERR_PLIB, "%s = strings[2] != strings[0] = %s", strings[2], strings[0]);
795: PetscCall(PetscStrcmp(strings[3], strings[1], &flg));
796: PetscCheck(flg, comm, PETSC_ERR_PLIB, "%s = strings[3] != strings[1] = %s", strings[3], strings[1]);
797: for (i = 0; i < nv; i++) PetscCall(PetscFree(strings[i]));
799: PetscCall(PetscViewerFlush(viewer));
800: if (verbose) PetscCall(PetscPrintf(comm, "# END testAttributesDefaultValue\n"));
801: #undef nv
802: PetscFunctionReturn(PETSC_SUCCESS);
803: }
805: int main(int argc, char **argv)
806: {
807: static char filename[PETSC_MAX_PATH_LEN] = "ex48.h5";
808: PetscMPIInt rank;
809: MPI_Comm comm;
810: PetscViewer viewer;
812: PetscFunctionBegin;
813: PetscFunctionBeginUser;
814: PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
815: comm = PETSC_COMM_WORLD;
816: PetscCallMPI(MPI_Comm_rank(comm, &rank));
817: PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL));
818: PetscCall(PetscOptionsGetBool(NULL, NULL, "-verbose", &verbose, NULL));
819: PetscCall(PetscOptionsGetString(NULL, NULL, "-filename", filename, sizeof(filename), NULL));
820: if (verbose) PetscCall(PetscPrintf(comm, "np ns " PetscStringize(np) " " PetscStringize(ns) "\n"));
822: PetscCall(PetscViewerHDF5Open(comm, filename, FILE_MODE_WRITE, &viewer));
823: PetscCall(testGroupsDatasets(viewer));
824: PetscCall(testAttributesAbsolutePath(viewer, "/"));
825: PetscCall(testAttributesAbsolutePath(viewer, "/prefix"));
826: PetscCall(PetscViewerDestroy(&viewer));
828: /* test reopening in update mode */
829: PetscCall(PetscViewerHDF5Open(comm, filename, FILE_MODE_UPDATE, &viewer));
830: PetscCall(testAttributesPushedPath(viewer));
831: PetscCall(testObjectAttributes(viewer));
832: PetscCall(testAttributesDefaultValue(viewer));
833: PetscCall(PetscViewerDestroy(&viewer));
834: PetscCall(PetscFinalize());
835: return 0;
836: }
838: /*TEST
840: build:
841: requires: hdf5
843: test:
844: nsize: {{1 4}}
846: TEST*/