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:   options->fname[0]             = '\0';
 99:   options->includes_constraints = PETSC_TRUE;
100:   PetscOptionsBegin(comm, "", "PetscSectionView()/Load() in HDF5 Test Options", "DMPLEX");
101:   PetscOptionsString("-fname", "The output file", "ex5.c", options->fname, options->fname, sizeof(options->fname), NULL);
102:   PetscOptionsBool("-includes_constraints", "Flag for if global section is to include constrained DoFs or not", "ex5.c", options->includes_constraints, &options->includes_constraints, NULL);
103:   PetscOptionsEnd();
104:   return 0;
105: }

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

114:   PetscInitialize(&argc, &argv, NULL, help);
115:   ProcessOptions(PETSC_COMM_WORLD, &user);
116:   MPI_Comm_size(PETSC_COMM_WORLD, &size);
117:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);

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

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

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

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

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

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

241: /*TEST

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

255: TEST*/