Actual source code: ex25.c

  1: static char help[]     = "Test DMCreateLocalVector_Plex, DMPlexGetCellFields and DMPlexRestoreCellFields work properly for 0 fields/cells/DS dimension\n\n";
  2: static char FILENAME[] = "ex25.c";

  4: #include <petscdmplex.h>
  5: #include <petscds.h>
  6: #include <petscsnes.h>

  8: typedef struct {
  9:   PetscInt test;
 10: } AppCtx;

 12: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
 13: {
 14:   options->test = 0;
 15:   PetscOptionsBegin(comm, "", "Zero-sized DMPlexGetCellFields Test Options", "DMPLEX");
 16:   PetscOptionsBoundedInt("-test", "Test to run", FILENAME, options->test, &options->test, NULL, 0);
 17:   PetscOptionsEnd();
 18:   return 0;
 19: }

 21: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *options, DM *dm)
 22: {
 23:   DMCreate(comm, dm);
 24:   DMSetType(*dm, DMPLEX);
 25:   DMSetFromOptions(*dm);
 26:   DMViewFromOptions(*dm, NULL, "-dm_view");
 27:   return 0;
 28: }

 30: /* no discretization is given so DMGetNumFields yields 0 */
 31: static PetscErrorCode test0(DM dm, AppCtx *options)
 32: {
 33:   Vec locX;

 35:   DMGetLocalVector(dm, &locX);
 36:   DMRestoreLocalVector(dm, &locX);
 37:   return 0;
 38: }

 40: /* no discretization is given so DMGetNumFields and PetscDSGetTotalDimension yield 0 */
 41: static PetscErrorCode test1(DM dm, AppCtx *options)
 42: {
 43:   IS           cells;
 44:   Vec          locX, locX_t, locA;
 45:   PetscScalar *u, *u_t, *a;

 47:   ISCreateStride(PETSC_COMM_SELF, 0, 0, 1, &cells);
 48:   DMGetLocalVector(dm, &locX);
 49:   DMGetLocalVector(dm, &locX_t);
 50:   DMGetLocalVector(dm, &locA);
 51:   DMPlexGetCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a);
 52:   DMPlexRestoreCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a);
 53:   DMRestoreLocalVector(dm, &locX);
 54:   DMRestoreLocalVector(dm, &locX_t);
 55:   DMRestoreLocalVector(dm, &locA);
 56:   ISDestroy(&cells);
 57:   return 0;
 58: }

 60: /* no discretization is given so DMGetNumFields and PetscDSGetTotalDimension yield 0 */
 61: static PetscErrorCode test2(DM dm, AppCtx *options)
 62: {
 63:   IS           cells;
 64:   Vec          locX, locX_t, locA;
 65:   PetscScalar *u, *u_t, *a;
 66:   PetscMPIInt  rank;

 68:   MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank);
 69:   ISCreateStride(PETSC_COMM_SELF, rank ? 0 : 1, 0, 1, &cells);
 70:   DMGetLocalVector(dm, &locX);
 71:   DMGetLocalVector(dm, &locX_t);
 72:   DMGetLocalVector(dm, &locA);
 73:   DMPlexGetCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a);
 74:   DMPlexRestoreCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a);
 75:   DMRestoreLocalVector(dm, &locX);
 76:   DMRestoreLocalVector(dm, &locX_t);
 77:   DMRestoreLocalVector(dm, &locA);
 78:   ISDestroy(&cells);
 79:   return 0;
 80: }

 82: static PetscErrorCode test3(DM dm, AppCtx *options)
 83: {
 84:   PetscDS   ds;
 85:   PetscFE   fe;
 86:   PetscInt  dim;
 87:   PetscBool simplex;

 89:   DMGetDimension(dm, &dim);
 90:   DMPlexIsSimplex(dm, &simplex);
 91:   DMGetDS(dm, &ds);
 92:   PetscFECreateDefault(PetscObjectComm((PetscObject)dm), dim, 1, simplex, NULL, -1, &fe);
 93:   PetscDSSetDiscretization(ds, 0, (PetscObject)fe);
 94:   PetscFEDestroy(&fe);
 95:   test1(dm, options);
 96:   return 0;
 97: }

 99: static PetscErrorCode test4(DM dm, AppCtx *options)
100: {
101:   PetscFE   fe;
102:   PetscInt  dim;
103:   PetscBool simplex;

105:   DMGetDimension(dm, &dim);
106:   DMPlexIsSimplex(dm, &simplex);
107:   PetscFECreateDefault(PetscObjectComm((PetscObject)dm), dim, 1, simplex, NULL, -1, &fe);
108:   DMSetField(dm, 0, NULL, (PetscObject)fe);
109:   PetscFEDestroy(&fe);
110:   DMCreateDS(dm);
111:   test2(dm, options);
112:   return 0;
113: }

115: static PetscErrorCode test5(DM dm, AppCtx *options)
116: {
117:   IS           cells;
118:   Vec          locX, locX_t, locA;
119:   PetscScalar *u, *u_t, *a;

121:   locX_t = NULL;
122:   locA   = NULL;
123:   ISCreateStride(PETSC_COMM_SELF, 0, 0, 1, &cells);
124:   DMGetLocalVector(dm, &locX);
125:   DMPlexGetCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a);
126:   DMPlexRestoreCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a);
127:   DMRestoreLocalVector(dm, &locX);
128:   ISDestroy(&cells);
129:   return 0;
130: }

132: static PetscErrorCode test6(DM dm, AppCtx *options)
133: {
134:   IS           cells;
135:   Vec          locX, locX_t, locA;
136:   PetscScalar *u, *u_t, *a;
137:   PetscMPIInt  rank;

139:   MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank);
140:   locX_t = NULL;
141:   locA   = NULL;
142:   ISCreateStride(PETSC_COMM_SELF, rank ? 0 : 1, 0, 1, &cells);
143:   DMGetLocalVector(dm, &locX);
144:   DMPlexGetCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a);
145:   DMPlexRestoreCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a);
146:   DMRestoreLocalVector(dm, &locX);
147:   ISDestroy(&cells);
148:   return 0;
149: }

151: static PetscErrorCode test7(DM dm, AppCtx *options)
152: {
153:   PetscFE   fe;
154:   PetscInt  dim;
155:   PetscBool simplex;

157:   DMGetDimension(dm, &dim);
158:   DMPlexIsSimplex(dm, &simplex);
159:   PetscFECreateDefault(PetscObjectComm((PetscObject)dm), dim, 1, simplex, NULL, -1, &fe);
160:   DMSetField(dm, 0, NULL, (PetscObject)fe);
161:   PetscFEDestroy(&fe);
162:   DMCreateDS(dm);
163:   test5(dm, options);
164:   return 0;
165: }

167: static PetscErrorCode test8(DM dm, AppCtx *options)
168: {
169:   PetscFE   fe;
170:   PetscInt  dim;
171:   PetscBool simplex;

173:   DMGetDimension(dm, &dim);
174:   DMPlexIsSimplex(dm, &simplex);
175:   PetscFECreateDefault(PetscObjectComm((PetscObject)dm), dim, 1, simplex, NULL, -1, &fe);
176:   DMSetField(dm, 0, NULL, (PetscObject)fe);
177:   PetscFEDestroy(&fe);
178:   DMCreateDS(dm);
179:   test6(dm, options);
180:   return 0;
181: }

183: int main(int argc, char **argv)
184: {
185:   MPI_Comm comm;
186:   DM       dm;
187:   AppCtx   options;

190:   PetscInitialize(&argc, &argv, NULL, help);
191:   comm = PETSC_COMM_WORLD;
192:   ProcessOptions(comm, &options);
193:   CreateMesh(comm, &options, &dm);

195:   switch (options.test) {
196:   case 0:
197:     test0(dm, &options);
198:     break;
199:   case 1:
200:     test1(dm, &options);
201:     break;
202:   case 2:
203:     test2(dm, &options);
204:     break;
205:   case 3:
206:     test3(dm, &options);
207:     break;
208:   case 4:
209:     test4(dm, &options);
210:     break;
211:   case 5:
212:     test5(dm, &options);
213:     break;
214:   case 6:
215:     test6(dm, &options);
216:     break;
217:   case 7:
218:     test7(dm, &options);
219:     break;
220:   case 8:
221:     test8(dm, &options);
222:     break;
223:   default:
224:     SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No such test: %" PetscInt_FMT, options.test);
225:   }

227:   DMDestroy(&dm);
228:   PetscFinalize();
229:   return 0;
230: }

232: /*TEST

234:   testset:
235:     args: -dm_plex_dim 3 -dm_plex_interpolate 0

237:     test:
238:       suffix: 0
239:       requires: ctetgen
240:       args: -test 0
241:     test:
242:       suffix: 1
243:       requires: ctetgen
244:       args: -test 1
245:     test:
246:       suffix: 2
247:       requires: ctetgen
248:       args: -test 2
249:     test:
250:       suffix: 3
251:       requires: ctetgen
252:       args: -test 3
253:     test:
254:       suffix: 4
255:       requires: ctetgen
256:       args: -test 4
257:     test:
258:       suffix: 5
259:       requires: ctetgen
260:       args: -test 5
261:     test:
262:       suffix: 6
263:       requires: ctetgen
264:       args: -test 6
265:     test:
266:       suffix: 7
267:       requires: ctetgen
268:       args: -test 7
269:     test:
270:       suffix: 8
271:       requires: ctetgen
272:       args: -test 8
273:     test:
274:       suffix: 9
275:       requires: ctetgen
276:       nsize: 2
277:       args: -test 1
278:     test:
279:       suffix: 10
280:       requires: ctetgen
281:       nsize: 2
282:       args: -test 2
283:     test:
284:       suffix: 11
285:       requires: ctetgen
286:       nsize: 2
287:       args: -test 3
288:     test:
289:       suffix: 12
290:       requires: ctetgen
291:       nsize: 2
292:       args: -test 4

294: TEST*/