Actual source code: ex5.c

  1: static char help[] = "Tests PetscSectionView()/Load() with HDF5.\n\n";

  3: #include <petscdmshell.h>
  4: #include <petscdmplex.h>
  5: #include <petscsection.h>
  6: #include <petscsf.h>
  7: #include <petsclayouthdf5.h>

  9: /* Save/Load abstract sections

 11: =====================
 12:  Save on 2 processes
 13: =====================

 15: section:
 16:                          0   1   2   3
 17:   rank 0: Dof (Field 0)  2   3   5   7
 18:           Dof (Field 1)  1   0   0   0

 20:                          0   1   2
 21:   rank 1: Dof (Field 0)  7   5  11 <- DoF 7 is constrained
 22:           Dof (Field 1)  0   0   2

 24: sf:
 25:   [0] 3 <- (1, 0)
 26:   [1] 1 <- (0, 2)

 28: global section (includesConstraints = PETSC_FALSE):
 29:                          0   1   2   3
 30:   rank 0: Dof (Field 0)  2   3   5  -8
 31:           Off (Field 0)  0   3   6  -12
 32:           Dof (Field 1)  1   0   0  -1
 33:           Off (Field 1)  2   6  11  -19

 35:                          0   1   2
 36:   rank 1: Dof (Field 0)  7  -6  11
 37:           Off (Field 0) 11  -7  18
 38:           Dof (Field 1)  0  -1   2
 39:           Off (Field 1) 18 -12  28

 41: global section (includesConstraints = PETSC_TRUE):
 42:                          0   1   2   3
 43:   rank 0: Dof (Field 0)  2   3   5  -8
 44:           Off (Field 0)  0   3   6  -12
 45:           Dof (Field 1)  1   0   0  -1
 46:           Off (Field 1)  2   6  11  -19

 48:                          0   1   2
 49:   rank 1: Dof (Field 0)  7  -6  11
 50:           Off (Field 0) 11  -7  18
 51:           Dof (Field 1)  0  -1   2
 52:           Off (Field 1) 18 -12  29

 54: =====================
 55:  Load on 3 Processes
 56: =====================

 58: (Set chartSize = 4, 0, 1 for rank 0, 1, 2, respectively)

 60: global section (includesConstraints = PETSC_FALSE):

 62:   rank 0: Dof (Field 0)  2   3   5   7
 63:           Off (Field 0)  0   3   6  11
 64:           Dof (Field 1)  1   0   0   0
 65:           Off (Field 1)  2   6  11  18

 67:   rank 1: Dof (Field 0)
 68:           Dof (Field 1)

 70:   rank 2: Dof (Field 0) 11
 71:           Off (Field 0) 18
 72:           Dof (Field 1)  2
 73:           Off (Field 1) 28

 75: global section (includesConstraints = PETSC_TRUE):

 77:   rank 0: Dof (Field 0)  2   3   5   7
 78:           Off (Field 0)  0   3   6  11
 79:           Dof (Field 1)  1   0   0   0
 80:           Off (Field 1)  2   6  11  18

 82:   rank 1: Dof (Field 0)
 83:           Dof (Field 1)

 85:   rank 2: Dof (Field 0) 11
 86:           Off (Field 0) 18
 87:           Dof (Field 1)  2
 88:           Off (Field 1) 29
 89: */

 91: typedef struct {
 92:   char      fname[PETSC_MAX_PATH_LEN]; /* Output mesh filename */
 93:   PetscBool includes_constraints;      /* Flag for if global section is to include constrained DoFs or not */
 94: } AppCtx;

 96: PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
 97: {
 98:   PetscFunctionBegin;
 99:   options->fname[0]             = '\0';
100:   options->includes_constraints = PETSC_TRUE;
101:   PetscOptionsBegin(comm, "", "PetscSectionView()/Load() in HDF5 Test Options", "DMPLEX");
102:   PetscCall(PetscOptionsString("-fname", "The output file", "ex5.c", options->fname, options->fname, sizeof(options->fname), NULL));
103:   PetscCall(PetscOptionsBool("-includes_constraints", "Flag for if global section is to include constrained DoFs or not", "ex5.c", options->includes_constraints, &options->includes_constraints, NULL));
104:   PetscOptionsEnd();
105:   PetscFunctionReturn(PETSC_SUCCESS);
106: }

108: int main(int argc, char **argv)
109: {
110:   MPI_Comm    comm;
111:   PetscMPIInt size, rank, mycolor;
112:   AppCtx      user;

114:   PetscFunctionBeginUser;
115:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
116:   PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user));
117:   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
118:   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
119:   PetscCheck(size >= 3, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "Example only works with three or more processes");

121:   /* Save */
122:   mycolor = (PetscMPIInt)(rank >= 2);
123:   PetscCallMPI(MPI_Comm_split(PETSC_COMM_WORLD, mycolor, rank, &comm));
124:   if (mycolor == 0) {
125:     PetscSection section, gsection;
126:     PetscSF      sf;
127:     PetscInt     nroots = -1, nleaves = -1, *ilocal;
128:     PetscSFNode *iremote;
129:     PetscViewer  viewer;

131:     /* Create section */
132:     PetscCall(PetscSectionCreate(comm, &section));
133:     PetscCall(PetscSectionSetNumFields(section, 2));
134:     switch (rank) {
135:     case 0:
136:       PetscCall(PetscSectionSetChart(section, 0, 4));
137:       PetscCall(PetscSectionSetDof(section, 0, 3));
138:       PetscCall(PetscSectionSetDof(section, 1, 3));
139:       PetscCall(PetscSectionSetDof(section, 2, 5));
140:       PetscCall(PetscSectionSetDof(section, 3, 7));
141:       PetscCall(PetscSectionSetFieldDof(section, 0, 0, 2));
142:       PetscCall(PetscSectionSetFieldDof(section, 1, 0, 3));
143:       PetscCall(PetscSectionSetFieldDof(section, 2, 0, 5));
144:       PetscCall(PetscSectionSetFieldDof(section, 3, 0, 7));
145:       PetscCall(PetscSectionSetFieldDof(section, 0, 1, 1));
146:       break;
147:     case 1:
148:       PetscCall(PetscSectionSetChart(section, 0, 3));
149:       PetscCall(PetscSectionSetDof(section, 0, 7));
150:       PetscCall(PetscSectionSetDof(section, 1, 5));
151:       PetscCall(PetscSectionSetDof(section, 2, 13));
152:       PetscCall(PetscSectionSetConstraintDof(section, 2, 1));
153:       PetscCall(PetscSectionSetFieldDof(section, 0, 0, 7));
154:       PetscCall(PetscSectionSetFieldDof(section, 1, 0, 5));
155:       PetscCall(PetscSectionSetFieldDof(section, 2, 0, 11));
156:       PetscCall(PetscSectionSetFieldDof(section, 2, 1, 2));
157:       PetscCall(PetscSectionSetFieldConstraintDof(section, 2, 0, 1));
158:       break;
159:     }
160:     PetscCall(PetscSectionSetUp(section));
161:     if (rank == 1) {
162:       const PetscInt indices[]  = {7};
163:       const PetscInt indices0[] = {7};

165:       PetscCall(PetscSectionSetConstraintIndices(section, 2, indices));
166:       PetscCall(PetscSectionSetFieldConstraintIndices(section, 2, 0, indices0));
167:     }
168:     /* Create sf */
169:     switch (rank) {
170:     case 0:
171:       nroots  = 4;
172:       nleaves = 1;
173:       PetscCall(PetscMalloc1(nleaves, &ilocal));
174:       PetscCall(PetscMalloc1(nleaves, &iremote));
175:       ilocal[0]        = 3;
176:       iremote[0].rank  = 1;
177:       iremote[0].index = 0;
178:       break;
179:     case 1:
180:       nroots  = 3;
181:       nleaves = 1;
182:       PetscCall(PetscMalloc1(nleaves, &ilocal));
183:       PetscCall(PetscMalloc1(nleaves, &iremote));
184:       ilocal[0]        = 1;
185:       iremote[0].rank  = 0;
186:       iremote[0].index = 2;
187:       break;
188:     }
189:     PetscCall(PetscSFCreate(comm, &sf));
190:     PetscCall(PetscSFSetGraph(sf, nroots, nleaves, ilocal, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER));
191:     /* Create global section*/
192:     PetscCall(PetscSectionCreateGlobalSection(section, sf, PETSC_TRUE, user.includes_constraints, PETSC_FALSE, &gsection));
193:     PetscCall(PetscSFDestroy(&sf));
194:     /* View */
195:     PetscCall(PetscViewerHDF5Open(comm, user.fname, FILE_MODE_WRITE, &viewer));
196:     PetscCall(PetscSectionView(gsection, viewer));
197:     PetscCall(PetscViewerDestroy(&viewer));
198:     PetscCall(PetscObjectSetName((PetscObject)section, "Save: local section"));
199:     PetscCall(PetscSectionView(section, PETSC_VIEWER_STDOUT_(comm)));
200:     PetscCall(PetscObjectSetName((PetscObject)gsection, "Save: global section"));
201:     PetscCall(PetscSectionView(gsection, PETSC_VIEWER_STDOUT_(comm)));
202:     PetscCall(PetscSectionDestroy(&gsection));
203:     PetscCall(PetscSectionDestroy(&section));
204:   }
205:   PetscCallMPI(MPI_Comm_free(&comm));

207:   /* Load */
208:   mycolor = (PetscMPIInt)(rank >= 3);
209:   PetscCallMPI(MPI_Comm_split(PETSC_COMM_WORLD, mycolor, rank, &comm));
210:   if (mycolor == 0) {
211:     PetscSection section;
212:     PetscInt     chartSize = -1;
213:     PetscViewer  viewer;

215:     PetscCall(PetscSectionCreate(comm, &section));
216:     switch (rank) {
217:     case 0:
218:       chartSize = 4;
219:       break;
220:     case 1:
221:       chartSize = 0;
222:       break;
223:     case 2:
224:       chartSize = 1;
225:       break;
226:     }
227:     PetscCall(PetscSectionSetChart(section, 0, chartSize));
228:     PetscCall(PetscViewerHDF5Open(comm, user.fname, FILE_MODE_READ, &viewer));
229:     PetscCall(PetscSectionLoad(section, viewer));
230:     PetscCall(PetscViewerDestroy(&viewer));
231:     PetscCall(PetscObjectSetName((PetscObject)section, "Load: section"));
232:     PetscCall(PetscSectionView(section, PETSC_VIEWER_STDOUT_(comm)));
233:     PetscCall(PetscSectionDestroy(&section));
234:   }
235:   PetscCallMPI(MPI_Comm_free(&comm));

237:   /* Finalize */
238:   PetscCall(PetscFinalize());
239:   return 0;
240: }

242: /*TEST

244:   build:
245:     requires: hdf5
246:     requires: !complex
247:   testset:
248:     nsize: 4
249:     test:
250:       suffix: 0
251:       args: -fname ex5_dump.h5 -includes_constraints 0
252:     test:
253:       suffix: 1
254:       args: -fname ex5_dump.h5 -includes_constraints 1

256: TEST*/