Actual source code: sectionhdf5.c
1: #include <petsc/private/sectionimpl.h>
2: #include <petscsf.h>
3: #include <petscis.h>
4: #include <petscviewerhdf5.h>
5: #include <petsclayouthdf5.h>
7: #if defined(PETSC_HAVE_HDF5)
8: static PetscErrorCode PetscSectionView_HDF5_SingleField(PetscSection s, PetscViewer viewer)
9: {
10: MPI_Comm comm;
11: PetscInt pStart, pEnd, p, n;
12: PetscBool hasConstraints, includesConstraints;
13: IS dofIS, offIS, cdofIS, coffIS, cindIS;
14: PetscInt *dofs, *offs, *cdofs, *coffs, *cinds, dof, cdof, m, moff, i;
16: PetscObjectGetComm((PetscObject)s, &comm);
17: PetscSectionGetChart(s, &pStart, &pEnd);
18: hasConstraints = (s->bc) ? PETSC_TRUE : PETSC_FALSE;
19: MPIU_Allreduce(MPI_IN_PLACE, &hasConstraints, 1, MPIU_BOOL, MPI_LOR, comm);
20: for (p = pStart, n = 0, m = 0; p < pEnd; ++p) {
21: PetscSectionGetDof(s, p, &dof);
22: if (dof >= 0) {
23: if (hasConstraints) {
24: PetscSectionGetConstraintDof(s, p, &cdof);
25: m += cdof;
26: }
27: n++;
28: }
29: }
30: PetscMalloc1(n, &dofs);
31: PetscMalloc1(n, &offs);
32: if (hasConstraints) {
33: PetscMalloc1(n, &cdofs);
34: PetscMalloc1(n, &coffs);
35: PetscMalloc1(m, &cinds);
36: }
37: for (p = pStart, n = 0, m = 0; p < pEnd; ++p) {
38: PetscSectionGetDof(s, p, &dof);
39: if (dof >= 0) {
40: dofs[n] = dof;
41: PetscSectionGetOffset(s, p, &offs[n]);
42: if (hasConstraints) {
43: const PetscInt *cpinds;
45: PetscSectionGetConstraintDof(s, p, &cdof);
46: PetscSectionGetConstraintIndices(s, p, &cpinds);
47: cdofs[n] = cdof;
48: coffs[n] = m;
49: for (i = 0; i < cdof; ++i) cinds[m++] = cpinds[i];
50: }
51: n++;
52: }
53: }
54: if (hasConstraints) {
55: MPI_Scan(&m, &moff, 1, MPIU_INT, MPI_SUM, comm);
56: moff -= m;
57: for (p = 0; p < n; ++p) coffs[p] += moff;
58: }
59: PetscViewerHDF5WriteAttribute(viewer, NULL, "hasConstraints", PETSC_BOOL, (void *)&hasConstraints);
60: PetscSectionGetIncludesConstraints(s, &includesConstraints);
61: PetscViewerHDF5WriteAttribute(viewer, NULL, "includesConstraints", PETSC_BOOL, (void *)&includesConstraints);
62: ISCreateGeneral(comm, n, dofs, PETSC_OWN_POINTER, &dofIS);
63: PetscObjectSetName((PetscObject)dofIS, "atlasDof");
64: ISView(dofIS, viewer);
65: ISDestroy(&dofIS);
66: ISCreateGeneral(comm, n, offs, PETSC_OWN_POINTER, &offIS);
67: PetscObjectSetName((PetscObject)offIS, "atlasOff");
68: ISView(offIS, viewer);
69: ISDestroy(&offIS);
70: if (hasConstraints) {
71: PetscViewerHDF5PushGroup(viewer, "bc");
72: ISCreateGeneral(comm, n, cdofs, PETSC_OWN_POINTER, &cdofIS);
73: PetscObjectSetName((PetscObject)cdofIS, "atlasDof");
74: ISView(cdofIS, viewer);
75: ISDestroy(&cdofIS);
76: ISCreateGeneral(comm, n, coffs, PETSC_OWN_POINTER, &coffIS);
77: PetscObjectSetName((PetscObject)coffIS, "atlasOff");
78: ISView(coffIS, viewer);
79: ISDestroy(&coffIS);
80: PetscViewerHDF5PopGroup(viewer);
81: ISCreateGeneral(comm, m, cinds, PETSC_OWN_POINTER, &cindIS);
82: PetscObjectSetName((PetscObject)cindIS, "bcIndices");
83: ISView(cindIS, viewer);
84: ISDestroy(&cindIS);
85: }
86: return 0;
87: }
89: PetscErrorCode PetscSectionView_HDF5_Internal(PetscSection s, PetscViewer viewer)
90: {
91: PetscInt numFields, f;
93: PetscViewerHDF5PushGroup(viewer, "section");
94: PetscSectionGetNumFields(s, &numFields);
95: PetscViewerHDF5WriteAttribute(viewer, NULL, "numFields", PETSC_INT, (void *)&numFields);
96: PetscSectionView_HDF5_SingleField(s, viewer);
97: for (f = 0; f < numFields; ++f) {
98: char fname[PETSC_MAX_PATH_LEN];
99: const char *fieldName;
100: PetscInt fieldComponents, c;
102: PetscSNPrintf(fname, sizeof(fname), "field%" PetscInt_FMT, f);
103: PetscViewerHDF5PushGroup(viewer, fname);
104: PetscSectionGetFieldName(s, f, &fieldName);
105: PetscViewerHDF5WriteAttribute(viewer, NULL, "fieldName", PETSC_STRING, fieldName);
106: PetscSectionGetFieldComponents(s, f, &fieldComponents);
107: PetscViewerHDF5WriteAttribute(viewer, NULL, "fieldComponents", PETSC_INT, (void *)&fieldComponents);
108: for (c = 0; c < fieldComponents; ++c) {
109: char cname[PETSC_MAX_PATH_LEN];
110: const char *componentName;
112: PetscSNPrintf(cname, sizeof(cname), "component%" PetscInt_FMT, c);
113: PetscViewerHDF5PushGroup(viewer, cname);
114: PetscSectionGetComponentName(s, f, c, &componentName);
115: PetscViewerHDF5WriteAttribute(viewer, NULL, "componentName", PETSC_STRING, componentName);
116: PetscViewerHDF5PopGroup(viewer);
117: }
118: PetscSectionView_HDF5_SingleField(s->field[f], viewer);
119: PetscViewerHDF5PopGroup(viewer);
120: }
121: PetscViewerHDF5PopGroup(viewer);
122: return 0;
123: }
125: static PetscErrorCode PetscSectionLoad_HDF5_SingleField_SetConstraintIndices(PetscSection s, IS cindIS, IS coffIS)
126: {
127: MPI_Comm comm;
128: PetscInt pStart, pEnd, p, M, m, i, cdof;
129: const PetscInt *data;
130: PetscInt *cinds;
131: const PetscInt *coffs;
132: PetscInt *coffsets;
133: PetscSF sf;
134: PetscLayout layout;
136: PetscObjectGetComm((PetscObject)s, &comm);
137: PetscSectionGetChart(s, &pStart, &pEnd);
138: ISGetSize(cindIS, &M);
139: ISGetLocalSize(cindIS, &m);
140: PetscMalloc1(m, &coffsets);
141: ISGetIndices(coffIS, &coffs);
142: for (p = pStart, m = 0; p < pEnd; ++p) {
143: PetscSectionGetConstraintDof(s, p, &cdof);
144: for (i = 0; i < cdof; ++i) coffsets[m++] = coffs[p - pStart] + i;
145: }
146: ISRestoreIndices(coffIS, &coffs);
147: PetscSFCreate(comm, &sf);
148: PetscLayoutCreate(comm, &layout);
149: PetscLayoutSetSize(layout, M);
150: PetscLayoutSetLocalSize(layout, m);
151: PetscLayoutSetBlockSize(layout, 1);
152: PetscLayoutSetUp(layout);
153: PetscSFSetGraphLayout(sf, layout, m, NULL, PETSC_OWN_POINTER, coffsets);
154: PetscLayoutDestroy(&layout);
155: PetscFree(coffsets);
156: PetscMalloc1(m, &cinds);
157: ISGetIndices(cindIS, &data);
158: PetscSFBcastBegin(sf, MPIU_INT, data, cinds, MPI_REPLACE);
159: PetscSFBcastEnd(sf, MPIU_INT, data, cinds, MPI_REPLACE);
160: ISRestoreIndices(cindIS, &data);
161: PetscSFDestroy(&sf);
162: PetscSectionSetUpBC(s);
163: for (p = pStart, m = 0; p < pEnd; ++p) {
164: PetscSectionGetConstraintDof(s, p, &cdof);
165: PetscSectionSetConstraintIndices(s, p, &cinds[m]);
166: m += cdof;
167: }
168: PetscFree(cinds);
169: return 0;
170: }
172: static PetscErrorCode PetscSectionLoad_HDF5_SingleField(PetscSection s, PetscViewer viewer)
173: {
174: MPI_Comm comm;
175: PetscInt pStart, pEnd, p, N, n, M, m;
176: #if defined(PETSC_USE_DEBUG)
177: PetscInt N1, M1;
178: #endif
179: PetscBool hasConstraints, includesConstraints;
180: IS dofIS, offIS, cdofIS, coffIS, cindIS;
181: const PetscInt *dofs, *offs, *cdofs;
182: PetscLayout map;
184: PetscObjectGetComm((PetscObject)s, &comm);
185: PetscViewerHDF5ReadAttribute(viewer, NULL, "includesConstraints", PETSC_BOOL, NULL, (void *)&includesConstraints);
186: PetscSectionSetIncludesConstraints(s, includesConstraints);
187: PetscSectionGetChart(s, &pStart, &pEnd);
188: n = pEnd - pStart;
189: #if defined(PETSC_USE_DEBUG)
190: MPIU_Allreduce(&n, &N1, 1, MPIU_INT, MPI_SUM, comm);
191: #endif
192: ISCreate(comm, &dofIS);
193: PetscObjectSetName((PetscObject)dofIS, "atlasDof");
194: PetscViewerHDF5ReadSizes(viewer, "atlasDof", NULL, &N);
195: #if defined(PETSC_USE_DEBUG)
197: #endif
198: ISGetLayout(dofIS, &map);
199: PetscLayoutSetSize(map, N);
200: PetscLayoutSetLocalSize(map, n);
201: ISLoad(dofIS, viewer);
202: ISCreate(comm, &offIS);
203: PetscObjectSetName((PetscObject)offIS, "atlasOff");
204: PetscViewerHDF5ReadSizes(viewer, "atlasOff", NULL, &N);
205: #if defined(PETSC_USE_DEBUG)
207: #endif
208: ISGetLayout(offIS, &map);
209: PetscLayoutSetSize(map, N);
210: PetscLayoutSetLocalSize(map, n);
211: ISLoad(offIS, viewer);
212: ISGetIndices(dofIS, &dofs);
213: ISGetIndices(offIS, &offs);
214: for (p = pStart, n = 0; p < pEnd; ++p, ++n) {
215: PetscSectionSetDof(s, p, dofs[n]);
216: PetscSectionSetOffset(s, p, offs[n]);
217: }
218: ISRestoreIndices(dofIS, &dofs);
219: ISRestoreIndices(offIS, &offs);
220: ISDestroy(&dofIS);
221: ISDestroy(&offIS);
222: PetscViewerHDF5ReadAttribute(viewer, NULL, "hasConstraints", PETSC_BOOL, NULL, (void *)&hasConstraints);
223: if (hasConstraints) {
224: PetscViewerHDF5PushGroup(viewer, "bc");
225: ISCreate(comm, &cdofIS);
226: PetscObjectSetName((PetscObject)cdofIS, "atlasDof");
227: PetscViewerHDF5ReadSizes(viewer, "atlasDof", NULL, &N);
228: #if defined(PETSC_USE_DEBUG)
230: #endif
231: ISGetLayout(cdofIS, &map);
232: PetscLayoutSetSize(map, N);
233: PetscLayoutSetLocalSize(map, n);
234: ISLoad(cdofIS, viewer);
235: ISGetIndices(cdofIS, &cdofs);
236: for (p = pStart, n = 0; p < pEnd; ++p, ++n) PetscSectionSetConstraintDof(s, p, cdofs[n]);
237: ISRestoreIndices(cdofIS, &cdofs);
238: ISDestroy(&cdofIS);
239: ISCreate(comm, &coffIS);
240: PetscObjectSetName((PetscObject)coffIS, "atlasOff");
241: PetscViewerHDF5ReadSizes(viewer, "atlasOff", NULL, &N);
242: #if defined(PETSC_USE_DEBUG)
244: #endif
245: ISGetLayout(coffIS, &map);
246: PetscLayoutSetSize(map, N);
247: PetscLayoutSetLocalSize(map, n);
248: ISLoad(coffIS, viewer);
249: PetscViewerHDF5PopGroup(viewer);
250: ISCreate(comm, &cindIS);
251: PetscObjectSetName((PetscObject)cindIS, "bcIndices");
252: PetscViewerHDF5ReadSizes(viewer, "bcIndices", NULL, &M);
253: if (!s->bc) m = 0;
254: else PetscSectionGetStorageSize(s->bc, &m);
255: #if defined(PETSC_USE_DEBUG)
256: MPIU_Allreduce(&m, &M1, 1, MPIU_INT, MPI_SUM, comm);
258: #endif
259: ISGetLayout(cindIS, &map);
260: PetscLayoutSetSize(map, M);
261: PetscLayoutSetLocalSize(map, m);
262: ISLoad(cindIS, viewer);
263: PetscSectionLoad_HDF5_SingleField_SetConstraintIndices(s, cindIS, coffIS);
264: ISDestroy(&coffIS);
265: ISDestroy(&cindIS);
266: }
267: return 0;
268: }
270: PetscErrorCode PetscSectionLoad_HDF5_Internal(PetscSection s, PetscViewer viewer)
271: {
272: MPI_Comm comm;
273: PetscInt N, n, numFields, f;
275: PetscObjectGetComm((PetscObject)s, &comm);
276: PetscViewerHDF5PushGroup(viewer, "section");
277: PetscViewerHDF5ReadAttribute(viewer, NULL, "numFields", PETSC_INT, NULL, (void *)&numFields);
278: if (s->pStart < 0 && s->pEnd < 0) n = PETSC_DECIDE;
279: else {
282: n = s->pEnd;
283: }
284: if (numFields > 0) PetscSectionSetNumFields(s, numFields);
285: PetscViewerHDF5ReadSizes(viewer, "atlasDof", NULL, &N);
286: if (n == PETSC_DECIDE) PetscSplitOwnership(comm, &n, &N);
287: PetscSectionSetChart(s, 0, n);
288: PetscSectionLoad_HDF5_SingleField(s, viewer);
289: for (f = 0; f < numFields; ++f) {
290: char fname[PETSC_MAX_PATH_LEN];
291: char *fieldName;
292: PetscInt fieldComponents, c;
294: PetscSNPrintf(fname, sizeof(fname), "field%" PetscInt_FMT, f);
295: PetscViewerHDF5PushGroup(viewer, fname);
296: PetscViewerHDF5ReadAttribute(viewer, NULL, "fieldName", PETSC_STRING, NULL, &fieldName);
297: PetscSectionSetFieldName(s, f, fieldName);
298: PetscFree(fieldName);
299: PetscViewerHDF5ReadAttribute(viewer, NULL, "fieldComponents", PETSC_INT, NULL, (void *)&fieldComponents);
300: PetscSectionSetFieldComponents(s, f, fieldComponents);
301: for (c = 0; c < fieldComponents; ++c) {
302: char cname[PETSC_MAX_PATH_LEN];
303: char *componentName;
305: PetscSNPrintf(cname, sizeof(cname), "component%" PetscInt_FMT, c);
306: PetscViewerHDF5PushGroup(viewer, cname);
307: PetscViewerHDF5ReadAttribute(viewer, NULL, "componentName", PETSC_STRING, NULL, &componentName);
308: PetscSectionSetComponentName(s, f, c, componentName);
309: PetscFree(componentName);
310: PetscViewerHDF5PopGroup(viewer);
311: }
312: PetscSectionLoad_HDF5_SingleField(s->field[f], viewer);
313: PetscViewerHDF5PopGroup(viewer);
314: }
315: PetscViewerHDF5PopGroup(viewer);
316: return 0;
317: }
319: #endif