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: PetscStrlen(name, &len);
89: *has = emptyExists;
90: if (len) {
91: char *loc = NULL;
92: PetscStrstr(name, "nonExisting", &loc);
93: *has = PetscNot(loc);
94: }
95: return 0;
96: }
98: static inline PetscErrorCode isPop(const char path[], PetscBool *has)
99: {
100: PetscStrcmp(path, "<", has);
101: return 0;
102: }
104: static inline PetscErrorCode isDot(const char path[], PetscBool *has)
105: {
106: PetscStrcmp(path, ".", has);
107: return 0;
108: }
110: static inline PetscErrorCode isRoot(const char path[], PetscBool *flg)
111: {
112: size_t len;
114: PetscStrlen(path, &len);
115: *flg = PetscNot(len);
116: if (!*flg) PetscStrcmp(path, "/", flg);
117: return 0;
118: }
120: static inline PetscErrorCode compare(PetscDataType dt, void *ptr0, void *ptr1, PetscBool *flg)
121: {
122: switch (dt) {
123: case PETSC_INT:
124: *flg = (PetscBool)(*(PetscInt *)ptr0 == *(PetscInt *)ptr1);
125: if (verbose) {
126: if (*flg) {
127: PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT, *(PetscInt *)ptr0);
128: } else {
129: PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT " != %" PetscInt_FMT "\n", *(PetscInt *)ptr0, *(PetscInt *)ptr1);
130: }
131: }
132: break;
133: case PETSC_REAL:
134: *flg = (PetscBool)(*(PetscReal *)ptr0 == *(PetscReal *)ptr1);
135: if (verbose) {
136: if (*flg) {
137: PetscPrintf(PETSC_COMM_SELF, "%f", *(PetscReal *)ptr0);
138: } else {
139: PetscPrintf(PETSC_COMM_SELF, "%f != %f\n", *(PetscReal *)ptr0, *(PetscReal *)ptr1);
140: }
141: }
142: break;
143: case PETSC_BOOL:
144: *flg = (PetscBool)(*(PetscBool *)ptr0 == *(PetscBool *)ptr1);
145: if (verbose) {
146: if (*flg) {
147: PetscPrintf(PETSC_COMM_SELF, "%s", PetscBools[*(PetscBool *)ptr0]);
148: } else {
149: PetscPrintf(PETSC_COMM_SELF, "%s != %s\n", PetscBools[*(PetscBool *)ptr0], PetscBools[*(PetscBool *)ptr1]);
150: }
151: }
152: break;
153: case PETSC_STRING:
154: PetscStrcmp((const char *)ptr0, (const char *)ptr1, flg);
155: if (verbose) {
156: if (*flg) {
157: PetscPrintf(PETSC_COMM_SELF, "%s", (char *)ptr0);
158: } else {
159: PetscPrintf(PETSC_COMM_SELF, "%s != %s\n", (char *)ptr0, (char *)ptr1);
160: }
161: }
162: break;
163: default:
164: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "PetscDataType %s not handled here", PetscDataTypes[dt]);
165: }
166: return 0;
167: }
169: static inline PetscErrorCode alterString(const char oldstr[], char str[])
170: {
171: size_t i, n;
173: PetscStrcpy(str, oldstr);
174: PetscStrlen(oldstr, &n);
175: for (i = 0; i < n; i++) {
176: if (('A' <= str[i] && str[i] < 'Z') || ('a' <= str[i] && str[i] < 'z')) {
177: str[i]++;
178: break;
179: }
180: }
181: return 0;
182: }
184: /* if name given, check dataset with this name exists under current group, otherwise just check current group exists */
185: /* flg: 0 doesn't exist, 1 group, 2 dataset */
186: static PetscErrorCode hasGroupOrDataset(PetscViewer viewer, const char path[], int *flg)
187: {
188: PetscBool has;
190: *flg = 0;
191: PetscViewerHDF5HasGroup(viewer, path, &has);
192: if (has) *flg = 1;
193: else {
194: PetscViewerHDF5HasDataset(viewer, path, &has);
195: if (has) *flg = 2;
196: }
197: return 0;
198: }
200: #define nt 5 /* number of datatypes */
201: typedef struct _n_Capsule *Capsule;
202: struct _n_Capsule {
203: char names[nt][SLEN];
204: PetscDataType types[nt];
205: char typeNames[nt][SLEN];
206: size_t sizes[nt];
207: void *vals[nt];
208: PetscInt id, ntypes;
209: };
211: static PetscErrorCode CapsuleCreate(Capsule old, Capsule *newcapsule)
212: {
213: Capsule c;
214: PetscBool bool0 = PETSC_TRUE;
215: PetscInt int0 = -1;
216: PetscReal real0 = -1.1;
217: char str0[] = "Test String";
218: char nestr0[] = "NONEXISTING STRING"; /* this attribute shall be skipped for writing */
219: void *vals[nt] = {&bool0, &int0, &real0, str0, nestr0};
220: size_t sizes[nt] = {sizeof(bool0), sizeof(int0), sizeof(real0), sizeof(str0), sizeof(str0)};
221: PetscDataType types[nt] = {PETSC_BOOL, PETSC_INT, PETSC_REAL, PETSC_STRING, PETSC_STRING};
222: const char *tNames[nt] = {"bool", "int", "real", "str", "nonExisting"};
223: PetscInt t;
225: PetscNew(&c);
226: c->id = 0;
227: c->ntypes = nt;
228: if (old) {
229: /* alter values */
230: t = 0;
231: bool0 = PetscNot(*((PetscBool *)old->vals[t]));
232: t++;
233: int0 = *((PetscInt *)old->vals[t]) * -2;
234: t++;
235: real0 = *((PetscReal *)old->vals[t]) * -2.0;
236: t++;
237: alterString((const char *)old->vals[t], str0);
238: t++;
239: c->id = old->id + 1;
240: }
241: for (t = 0; t < nt; t++) {
242: c->sizes[t] = sizes[t];
243: c->types[t] = types[t];
244: PetscStrcpy(c->typeNames[t], tNames[t]);
245: PetscSNPrintf(c->names[t], SLEN, "attr_%" PetscInt_FMT "_%s", c->id, tNames[t]);
246: PetscMalloc(sizes[t], &c->vals[t]);
247: PetscMemcpy(c->vals[t], vals[t], sizes[t]);
248: }
249: *newcapsule = c;
250: return 0;
251: }
252: #undef nt
254: static PetscErrorCode CapsuleWriteAttributes(Capsule c, PetscViewer v, const char parent[])
255: {
256: PetscInt t;
257: PetscBool flg = PETSC_FALSE;
259: for (t = 0; t < c->ntypes; t++) {
260: shouldExist(c->names[t], PETSC_FALSE, &flg);
261: if (!flg) continue;
262: PetscViewerHDF5WriteAttribute(v, parent, c->names[t], c->types[t], c->vals[t]);
263: }
264: return 0;
265: }
267: static PetscErrorCode CapsuleReadAndCompareAttributes(Capsule c, PetscViewer v, const char parent[])
268: {
269: char *group;
270: int gd = 0;
271: PetscInt t;
272: PetscBool flg = PETSC_FALSE, hasAttr = PETSC_FALSE;
273: MPI_Comm comm;
275: PetscObjectGetComm((PetscObject)v, &comm);
276: PetscViewerHDF5GetGroup(v, NULL, &group);
277: hasGroupOrDataset(v, parent, &gd);
278: /* check correct existence of attributes */
279: for (t = 0; t < c->ntypes; t++) {
280: const char *attribute = c->names[t];
281: shouldExist(attribute, PETSC_FALSE, &flg);
282: PetscViewerHDF5HasAttribute(v, parent, attribute, &hasAttr);
283: if (verbose) {
284: PetscPrintf(comm, " %-24s = ", attribute);
285: if (!hasAttr) PetscPrintf(comm, "---");
286: }
290: /* check loaded attributes are the same as original */
291: if (hasAttr) {
292: char buffer[SLEN];
293: char *str;
294: void *ptr0;
295: /* check the stored data is the same as original */
296: //TODO datatype should better be output arg, not input
297: //TODO string attributes should probably have a separate function since the handling is different;
298: //TODO or maybe it should just accept string buffer rather than pointer to string
299: if (c->types[t] == PETSC_STRING) {
300: PetscViewerHDF5ReadAttribute(v, parent, attribute, c->types[t], NULL, &str);
301: ptr0 = str;
302: } else {
303: PetscViewerHDF5ReadAttribute(v, parent, attribute, c->types[t], NULL, &buffer);
304: ptr0 = &buffer;
305: }
306: compare(c->types[t], ptr0, c->vals[t], &flg);
308: if (verbose) PetscPrintf(comm, " (=)");
309: if (c->types[t] == PETSC_STRING) PetscFree(str);
310: }
311: if (verbose && gd) PetscPrintf(comm, "\n");
312: }
313: PetscFree(group);
314: return 0;
315: }
317: static PetscErrorCode CapsuleDestroy(Capsule *c)
318: {
319: PetscInt t;
321: if (!*c) return 0;
322: for (t = 0; t < (*c)->ntypes; t++) PetscFree((*c)->vals[t]);
323: PetscFree(*c);
324: return 0;
325: }
327: static PetscErrorCode testGroupsDatasets(PetscViewer viewer)
328: {
329: char buf[PETSC_MAX_PATH_LEN];
330: Vec vecs[nap][ns];
331: PetscInt p, s;
332: PetscBool flg = PETSC_FALSE, flg1 = PETSC_FALSE, flg2 = PETSC_FALSE;
333: PetscRandom rand;
334: const char *filename;
335: MPI_Comm comm;
337: PetscObjectGetComm((PetscObject)viewer, &comm);
338: PetscViewerFileGetName(viewer, &filename);
339: if (verbose) PetscPrintf(comm, "# TEST testGroupsDatasets\n");
340: /* store random vectors */
341: PetscRandomCreate(comm, &rand);
342: PetscRandomSetInterval(rand, 0.0, 10.0);
343: PetscRandomSetFromOptions(rand);
344: PetscMemzero(vecs, nap * ns * sizeof(Vec));
346: /* test dataset writing */
347: if (verbose) PetscPrintf(comm, "## WRITE PHASE\n");
348: for (p = 0; p < np; p++) {
349: isPop(paths[p], &flg);
350: isDot(paths[p], &flg1);
351: shouldExist(apaths[paths2apaths[p]], PETSC_FALSE, &flg2);
352: if (flg) {
353: PetscViewerHDF5PopGroup(viewer);
354: } else {
355: PetscViewerHDF5PushGroup(viewer, paths[p]);
356: }
357: if (verbose) PetscPrintf(comm, "%-32s => %4s => %-32s should exist? %s\n", paths[p], flg ? "pop" : "push", apaths[paths2apaths[p]], PetscBools[flg2]);
358: if (flg || flg1 || !flg2) continue;
360: for (s = 0; s < ns; s++) {
361: Vec v;
363: shouldExist(datasets[s], PETSC_FALSE, &flg);
364: if (!flg) continue;
366: VecCreate(comm, &v);
367: PetscObjectSetName((PetscObject)v, datasets[s]);
368: VecSetSizes(v, n, PETSC_DECIDE);
369: VecSetFromOptions(v);
370: VecSetRandom(v, rand);
371: if (verbose) {
372: PetscReal min, max;
373: VecMin(v, NULL, &min);
374: VecMax(v, NULL, &max);
375: 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);
376: }
378: VecView(v, viewer);
379: vecs[paths2apaths[p]][s] = v;
380: }
381: }
382: PetscViewerFlush(viewer);
383: PetscRandomDestroy(&rand);
385: if (verbose) PetscPrintf(comm, "\n## READ PHASE\n");
386: /* check correct existence of groups in file */
387: for (p = 0; p < np; p++) {
388: char *group;
389: const char *expected = apaths[paths2apaths[p]];
391: /* check Push/Pop is correct */
392: isPop(paths[p], &flg);
393: if (flg) {
394: PetscViewerHDF5PopGroup(viewer);
395: } else {
396: PetscViewerHDF5PushGroup(viewer, paths[p]);
397: }
398: PetscViewerHDF5GetGroup(viewer, NULL, &group);
399: PetscViewerHDF5HasGroup(viewer, NULL, &flg1);
400: if (verbose) PetscPrintf(comm, "%-32s => %4s => %-32s exists? %s\n", paths[p], flg ? "pop" : "push", group, PetscBools[flg1]);
401: PetscStrcmp(group, expected, &flg2);
403: shouldExist(group, PETSC_TRUE, &flg2);
405: PetscFree(group);
406: }
408: /* check existence of datasets; compare loaded vectors with original ones */
409: for (p = 0; p < np; p++) {
410: char *group;
412: /* check Push/Pop is correct */
413: isPop(paths[p], &flg);
414: if (flg) {
415: PetscViewerHDF5PopGroup(viewer);
416: } else {
417: PetscViewerHDF5PushGroup(viewer, paths[p]);
418: }
419: PetscViewerHDF5GetGroup(viewer, NULL, &group);
420: PetscViewerHDF5HasGroup(viewer, NULL, &flg);
421: if (verbose) PetscPrintf(comm, "Has %s group? %s\n", group, PetscBools[flg]);
422: for (s = 0; s < ns; s++) {
423: const char *name = datasets[s];
424: char *fullname = buf;
426: /* check correct existence of datasets in file */
427: PetscSNPrintf(fullname, sizeof(buf), "%s/%s", group, name);
428: shouldExist(name, PETSC_FALSE, &flg1);
429: flg1 = (PetscBool)(flg && flg1); /* both group and dataset need to exist */
430: PetscViewerHDF5HasDataset(viewer, name, &flg2);
431: if (verbose) PetscPrintf(comm, " %s dataset? %s", fullname, PetscBools[flg2]);
434: if (flg2) {
435: Vec v;
436: /* check loaded Vec is the same as original */
437: VecCreate(comm, &v);
438: PetscObjectSetName((PetscObject)v, name);
439: VecLoad(v, viewer);
440: VecEqual(v, vecs[paths2apaths[p]][s], &flg1);
442: if (verbose) PetscPrintf(comm, " (=)");
443: VecDestroy(&v);
444: }
445: if (verbose) PetscPrintf(comm, "\n");
446: }
447: PetscFree(group);
448: }
449: PetscViewerFlush(viewer);
450: for (p = 0; p < nap; p++)
451: for (s = 0; s < ns; s++) VecDestroy(&vecs[p][s]);
452: if (verbose) PetscPrintf(comm, "# END testGroupsDatasets\n\n");
453: return 0;
454: }
456: static inline PetscErrorCode formPath(PetscBool relativize, const char path[], const char dataset[], char buf[], size_t bufsize)
457: {
458: PetscBool isroot = PETSC_FALSE;
460: isRoot(path, &isroot);
461: if (relativize) {
462: if (isroot) {
463: PetscStrncpy(buf, dataset, bufsize);
464: } else {
465: /* skip initial '/' in paths[p] if prefix given */
466: PetscSNPrintf(buf, bufsize, "%s/%s", path + 1, dataset);
467: }
468: } else {
469: PetscSNPrintf(buf, bufsize, "%s/%s", isroot ? "" : path, dataset);
470: }
471: return 0;
472: }
474: /* test attribute writing, existence checking and reading, use absolute paths */
475: static PetscErrorCode testAttributesAbsolutePath(PetscViewer viewer, const char prefix[])
476: {
477: char buf[PETSC_MAX_PATH_LEN];
478: Capsule capsules[nap][ns], c = NULL, old = NULL;
479: PetscInt p, s;
480: PetscBool flg = PETSC_FALSE, flg1 = PETSC_FALSE;
481: MPI_Comm comm;
483: PetscObjectGetComm((PetscObject)viewer, &comm);
484: if (verbose) {
485: if (prefix) {
486: PetscPrintf(comm, "# TEST testAttributesAbsolutePath, prefix=\"%s\"\n", prefix);
487: } else {
488: PetscPrintf(comm, "# TEST testAttributesAbsolutePath\n");
489: }
490: PetscPrintf(comm, "## WRITE PHASE\n");
491: }
492: PetscMemzero(capsules, nap * ns * sizeof(Capsule));
494: /* test attribute writing */
495: if (prefix) PetscViewerHDF5PushGroup(viewer, prefix);
496: for (p = 0; p < np; p++)
497: for (s = 0; s < ns; s++) {
498: /* we test only absolute paths here */
499: PetscViewerHDF5PathIsRelative(paths[p], PETSC_FALSE, &flg);
500: if (flg) continue;
501: {
502: char *group;
504: PetscViewerHDF5GetGroup(viewer, NULL, &group);
505: PetscStrcmp(group, prefix, &flg);
507: PetscFree(group);
508: }
509: formPath((PetscBool) !!prefix, paths[p], datasets[s], buf, sizeof(buf));
510: shouldExist(buf, PETSC_TRUE, &flg);
511: if (!flg) continue;
513: if (verbose) {
514: if (prefix) {
515: PetscPrintf(comm, "Write attributes to %s/%s\n", prefix, buf);
516: } else {
517: PetscPrintf(comm, "Write attributes to %s\n", buf);
518: }
519: }
521: CapsuleCreate(old, &c);
522: CapsuleWriteAttributes(c, viewer, buf);
524: capsules[paths2apaths[p]][s] = c;
525: old = c;
526: }
527: if (prefix) PetscViewerHDF5PopGroup(viewer);
528: PetscViewerFlush(viewer);
530: if (verbose) PetscPrintf(comm, "\n## READ PHASE\n");
531: if (prefix) PetscViewerHDF5PushGroup(viewer, prefix);
532: for (p = 0; p < np; p++)
533: for (s = 0; s < ns; s++) {
534: /* we test only absolute paths here */
535: PetscViewerHDF5PathIsRelative(paths[p], PETSC_FALSE, &flg);
536: if (flg) continue;
538: /* check existence of given group/dataset */
539: formPath((PetscBool) !!prefix, paths[p], datasets[s], buf, sizeof(buf));
540: shouldExist(buf, PETSC_TRUE, &flg);
541: if (verbose) {
542: if (prefix) {
543: PetscPrintf(comm, "Has %s/%s? %s\n", prefix, buf, PetscBools[flg]);
544: } else {
545: PetscPrintf(comm, "Has %s? %s\n", buf, PetscBools[flg]);
546: }
547: }
549: /* check attribute capsule has been created for given path */
550: c = capsules[paths2apaths[p]][s];
551: flg1 = (PetscBool) !!c;
553: if (!flg) continue;
555: /* check correct existence and fidelity of attributes in file */
556: CapsuleReadAndCompareAttributes(c, viewer, buf);
557: }
558: if (prefix) PetscViewerHDF5PopGroup(viewer);
559: PetscViewerFlush(viewer);
560: for (p = 0; p < nap; p++)
561: for (s = 0; s < ns; s++) CapsuleDestroy(&capsules[p][s]);
562: if (verbose) PetscPrintf(comm, "# END testAttributesAbsolutePath\n\n");
563: return 0;
564: }
566: /* test attribute writing, existence checking and reading, use group push/pop */
567: static PetscErrorCode testAttributesPushedPath(PetscViewer viewer)
568: {
569: Capsule capsules[nap][ns], c = NULL, old = NULL;
570: PetscInt p, s;
571: int gd;
572: PetscBool flg = PETSC_FALSE, flg1 = PETSC_FALSE;
573: MPI_Comm comm;
575: PetscObjectGetComm((PetscObject)viewer, &comm);
576: if (verbose) {
577: PetscPrintf(comm, "# TEST testAttributesPushedPath\n");
578: PetscPrintf(comm, "## WRITE PHASE\n");
579: }
580: PetscMemzero(capsules, nap * ns * sizeof(Capsule));
582: /* test attribute writing */
583: for (p = 0; p < np; p++) {
584: isPop(paths[p], &flg);
585: isDot(paths[p], &flg1);
586: if (flg) {
587: PetscViewerHDF5PopGroup(viewer);
588: } else {
589: PetscViewerHDF5PushGroup(viewer, paths[p]);
590: }
591: /* < and . have been already visited => skip */
592: if (flg || flg1) continue;
594: /* assume here that groups and datasets are already in the file */
595: for (s = 0; s < ns; s++) {
596: hasGroupOrDataset(viewer, datasets[s], &gd);
597: if (!gd) continue;
598: if (verbose) PetscPrintf(comm, "Write attributes to %s/%s\n", apaths[paths2apaths[p]], datasets[s]);
599: CapsuleCreate(old, &c);
600: CapsuleWriteAttributes(c, viewer, datasets[s]);
602: capsules[paths2apaths[p]][s] = c;
603: old = c;
604: }
605: }
606: PetscViewerFlush(viewer);
608: if (verbose) PetscPrintf(comm, "\n## READ PHASE\n");
609: for (p = 0; p < np; p++) {
610: char *group;
612: isPop(paths[p], &flg1);
613: if (flg1) {
614: PetscViewerHDF5PopGroup(viewer);
615: } else {
616: PetscViewerHDF5PushGroup(viewer, paths[p]);
617: }
618: PetscViewerHDF5GetGroup(viewer, NULL, &group);
619: for (s = 0; s < ns; s++) {
620: hasGroupOrDataset(viewer, datasets[s], &gd);
621: if (verbose) PetscPrintf(comm, "%s/%s %s\n", group, datasets[s], gd ? (gd == 1 ? "is group" : "is dataset") : "does not exist");
623: /* check attribute capsule has been created for given path */
624: c = capsules[paths2apaths[p]][s];
625: flg = (PetscBool) !!gd;
626: flg1 = (PetscBool) !!c;
628: if (!flg) continue;
630: /* check correct existence of attributes in file */
631: CapsuleReadAndCompareAttributes(c, viewer, datasets[s]);
632: }
633: PetscFree(group);
634: }
635: PetscViewerFlush(viewer);
636: for (p = 0; p < nap; p++)
637: for (s = 0; s < ns; s++) CapsuleDestroy(&capsules[p][s]);
638: if (verbose) PetscPrintf(comm, "# END testAttributesPushedPath\n\n");
639: return 0;
640: }
642: /* test attribute writing, existence checking and reading, use group push/pop */
643: static PetscErrorCode testObjectAttributes(PetscViewer viewer)
644: {
645: Capsule capsules[nap][ns], c = NULL, old = NULL;
646: PetscInt p, s;
647: PetscBool flg = PETSC_FALSE, flg1 = PETSC_FALSE;
648: MPI_Comm comm;
650: PetscObjectGetComm((PetscObject)viewer, &comm);
651: if (verbose) {
652: PetscPrintf(comm, "# TEST testObjectAttributes\n");
653: PetscPrintf(comm, "## WRITE PHASE\n");
654: }
655: PetscMemzero(capsules, nap * ns * sizeof(Capsule));
657: /* test attribute writing */
658: for (p = 0; p < np; p++) {
659: isPop(paths[p], &flg);
660: isDot(paths[p], &flg1);
661: if (flg) {
662: PetscViewerHDF5PopGroup(viewer);
663: } else {
664: PetscViewerHDF5PushGroup(viewer, paths[p]);
665: }
666: /* < and . have been already visited => skip */
667: if (flg || flg1) continue;
669: /* assume here that groups and datasets are already in the file */
670: for (s = 0; s < ns; s++) {
671: Vec v;
672: size_t len;
673: const char *name = datasets[s];
675: PetscStrlen(name, &len);
676: if (!len) continue;
677: VecCreate(comm, &v);
678: PetscObjectSetName((PetscObject)v, name);
679: PetscViewerHDF5HasObject(viewer, (PetscObject)v, &flg);
680: if (flg) {
681: if (verbose) PetscPrintf(comm, "Write attributes to %s/%s\n", apaths[paths2apaths[p]], name);
682: CapsuleCreate(old, &c);
683: CapsuleWriteAttributes(c, viewer, name);
685: capsules[paths2apaths[p]][s] = c;
686: old = c;
687: }
688: VecDestroy(&v);
689: }
690: }
691: PetscViewerFlush(viewer);
693: if (verbose) PetscPrintf(comm, "\n## READ PHASE\n");
694: for (p = 0; p < np; p++) {
695: char *group;
697: isPop(paths[p], &flg);
698: if (flg) {
699: PetscViewerHDF5PopGroup(viewer);
700: } else {
701: PetscViewerHDF5PushGroup(viewer, paths[p]);
702: }
703: PetscViewerHDF5GetGroup(viewer, NULL, &group);
704: for (s = 0; s < ns; s++) {
705: Vec v;
706: size_t len;
707: const char *name = datasets[s];
709: PetscStrlen(name, &len);
710: if (!len) continue;
711: VecCreate(comm, &v);
712: PetscObjectSetName((PetscObject)v, name);
713: PetscViewerHDF5HasObject(viewer, (PetscObject)v, &flg);
714: if (verbose) PetscPrintf(comm, "Is %s/%s dataset? %s\n", group, name, PetscBools[flg]);
716: /* check attribute capsule has been created for given path */
717: c = capsules[paths2apaths[p]][s];
718: flg1 = (PetscBool) !!c;
721: /* check correct existence of attributes in file */
722: if (flg) CapsuleReadAndCompareAttributes(c, viewer, name);
723: VecDestroy(&v);
724: }
725: PetscFree(group);
726: }
727: PetscViewerFlush(viewer);
728: for (p = 0; p < nap; p++)
729: for (s = 0; s < ns; s++) CapsuleDestroy(&capsules[p][s]);
730: if (verbose) PetscPrintf(comm, "# END testObjectAttributes\n\n");
731: return 0;
732: }
734: static PetscErrorCode testAttributesDefaultValue(PetscViewer viewer)
735: {
736: #define nv 4
737: PetscBool bools[nv];
738: PetscInt ints[nv];
739: PetscReal reals[nv];
740: char *strings[nv];
741: PetscBool flg;
742: PetscInt i;
743: MPI_Comm comm;
745: PetscObjectGetComm((PetscObject)viewer, &comm);
746: if (verbose) PetscPrintf(comm, "# TEST testAttributesDefaultValue\n");
748: PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_bool", PETSC_BOOL, NULL, &bools[0]);
749: bools[1] = PetscNot(bools[0]);
750: PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_bool", PETSC_BOOL, &bools[1], &bools[2]);
751: PetscViewerHDF5ReadAttribute(viewer, "/", "attr_nonExisting_bool", PETSC_BOOL, &bools[1], &bools[3]);
755: PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_int", PETSC_INT, NULL, &ints[0]);
756: ints[1] = ints[0] * -333;
757: PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_int", PETSC_INT, &ints[1], &ints[2]);
758: PetscViewerHDF5ReadAttribute(viewer, "/", "attr_nonExisting_int", PETSC_INT, &ints[1], &ints[3]);
761: if (verbose) PetscIntView(nv, ints, PETSC_VIEWER_STDOUT_WORLD);
763: PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_real", PETSC_REAL, NULL, &reals[0]);
764: reals[1] = reals[0] * -11.1;
765: PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_real", PETSC_REAL, &reals[1], &reals[2]);
766: PetscViewerHDF5ReadAttribute(viewer, "/", "attr_nonExisting_real", PETSC_REAL, &reals[1], &reals[3]);
769: if (verbose) PetscRealView(nv, reals, PETSC_VIEWER_STDOUT_WORLD);
771: PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_str", PETSC_STRING, NULL, &strings[0]);
772: PetscStrallocpy(strings[0], &strings[1]);
773: alterString(strings[0], strings[1]);
774: PetscViewerHDF5ReadAttribute(viewer, "/", "attr_0_str", PETSC_STRING, &strings[1], &strings[2]);
775: PetscViewerHDF5ReadAttribute(viewer, "/", "attr_nonExisting_str", PETSC_STRING, &strings[1], &strings[3]);
776: PetscStrcmp(strings[2], strings[0], &flg);
778: PetscStrcmp(strings[3], strings[1], &flg);
780: for (i = 0; i < nv; i++) PetscFree(strings[i]);
782: PetscViewerFlush(viewer);
783: if (verbose) PetscPrintf(comm, "# END testAttributesDefaultValue\n");
784: #undef nv
785: return 0;
786: }
788: int main(int argc, char **argv)
789: {
790: static char filename[PETSC_MAX_PATH_LEN] = "ex48.h5";
791: PetscMPIInt rank;
792: MPI_Comm comm;
793: PetscViewer viewer;
796: PetscInitialize(&argc, &argv, (char *)0, help);
797: comm = PETSC_COMM_WORLD;
798: MPI_Comm_rank(comm, &rank);
799: PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL);
800: PetscOptionsGetBool(NULL, NULL, "-verbose", &verbose, NULL);
801: PetscOptionsGetString(NULL, NULL, "-filename", filename, sizeof(filename), NULL);
802: if (verbose) PetscPrintf(comm, "np ns " PetscStringize(np) " " PetscStringize(ns) "\n");
804: PetscViewerHDF5Open(comm, filename, FILE_MODE_WRITE, &viewer);
805: testGroupsDatasets(viewer);
806: testAttributesAbsolutePath(viewer, "/");
807: testAttributesAbsolutePath(viewer, "/prefix");
808: PetscViewerDestroy(&viewer);
810: /* test reopening in update mode */
811: PetscViewerHDF5Open(comm, filename, FILE_MODE_UPDATE, &viewer);
812: testAttributesPushedPath(viewer);
813: testObjectAttributes(viewer);
814: testAttributesDefaultValue(viewer);
815: PetscViewerDestroy(&viewer);
816: PetscFinalize();
817: return 0;
818: }
820: /*TEST
822: build:
823: requires: hdf5
825: test:
826: nsize: {{1 4}}
828: TEST*/