Actual source code: ex1.c

  1: static char help[] = "Demonstration of creating and viewing DMFields objects.\n\n";

  3: #include <petscdmfield.h>
  4: #include <petscdmplex.h>
  5: #include <petscdmda.h>

  7: static PetscErrorCode ViewResults(PetscViewer viewer, PetscInt N, PetscInt dim, PetscScalar *B, PetscScalar *D, PetscScalar *H, PetscReal *rB, PetscReal *rD, PetscReal *rH)
  8: {
  9:   PetscFunctionBegin;
 10:   PetscCall(PetscViewerASCIIPrintf(viewer, "B:\n"));
 11:   PetscCall(PetscScalarView(N, B, viewer));
 12:   PetscCall(PetscViewerASCIIPrintf(viewer, "D:\n"));
 13:   PetscCall(PetscScalarView(N * dim, D, viewer));
 14:   PetscCall(PetscViewerASCIIPrintf(viewer, "H:\n"));
 15:   PetscCall(PetscScalarView(N * dim * dim, H, viewer));

 17:   PetscCall(PetscViewerASCIIPrintf(viewer, "rB:\n"));
 18:   PetscCall(PetscRealView(N, rB, viewer));
 19:   PetscCall(PetscViewerASCIIPrintf(viewer, "rD:\n"));
 20:   PetscCall(PetscRealView(N * dim, rD, viewer));
 21:   PetscCall(PetscViewerASCIIPrintf(viewer, "rH:\n"));
 22:   PetscCall(PetscRealView(N * dim * dim, rH, viewer));
 23:   PetscFunctionReturn(PETSC_SUCCESS);
 24: }

 26: static PetscErrorCode TestEvaluate(DMField field, PetscInt n, PetscRandom rand)
 27: {
 28:   DM           dm;
 29:   PetscInt     dim, i, nc;
 30:   PetscScalar *B, *D, *H;
 31:   PetscReal   *rB, *rD, *rH;
 32:   Vec          points;
 33:   PetscScalar *array;
 34:   PetscViewer  viewer;
 35:   MPI_Comm     comm;

 37:   PetscFunctionBegin;
 38:   comm = PetscObjectComm((PetscObject)field);
 39:   PetscCall(DMFieldGetNumComponents(field, &nc));
 40:   PetscCall(DMFieldGetDM(field, &dm));
 41:   PetscCall(DMGetDimension(dm, &dim));
 42:   PetscCall(VecCreateFromOptions(PetscObjectComm((PetscObject)field), NULL, 1, n * dim, PETSC_DETERMINE, &points));
 43:   PetscCall(VecSetBlockSize(points, dim));
 44:   PetscCall(VecGetArray(points, &array));
 45:   for (i = 0; i < n * dim; i++) PetscCall(PetscRandomGetValue(rand, &array[i]));
 46:   PetscCall(VecRestoreArray(points, &array));
 47:   PetscCall(PetscMalloc6(n * nc, &B, n * nc, &rB, n * nc * dim, &D, n * nc * dim, &rD, n * nc * dim * dim, &H, n * nc * dim * dim, &rH));
 48:   PetscCall(DMFieldEvaluate(field, points, PETSC_SCALAR, B, D, H));
 49:   PetscCall(DMFieldEvaluate(field, points, PETSC_REAL, rB, rD, rH));
 50:   viewer = PETSC_VIEWER_STDOUT_(comm);

 52:   PetscCall(PetscObjectSetName((PetscObject)points, "Test Points"));
 53:   PetscCall(VecView(points, viewer));
 54:   PetscCall(ViewResults(viewer, n * nc, dim, B, D, H, rB, rD, rH));

 56:   PetscCall(PetscFree6(B, rB, D, rD, H, rH));
 57:   PetscCall(VecDestroy(&points));
 58:   PetscFunctionReturn(PETSC_SUCCESS);
 59: }

 61: static PetscErrorCode TestEvaluateFE(DMField field, PetscInt n, PetscInt cStart, PetscInt cEnd, PetscQuadrature quad, PetscRandom rand)
 62: {
 63:   DM           dm;
 64:   PetscInt     dim, i, nc, nq;
 65:   PetscInt     N;
 66:   PetscScalar *B, *D, *H;
 67:   PetscReal   *rB, *rD, *rH;
 68:   PetscInt    *cells;
 69:   IS           cellIS;
 70:   PetscViewer  viewer;
 71:   MPI_Comm     comm;

 73:   PetscFunctionBegin;
 74:   comm = PetscObjectComm((PetscObject)field);
 75:   PetscCall(DMFieldGetNumComponents(field, &nc));
 76:   PetscCall(DMFieldGetDM(field, &dm));
 77:   PetscCall(DMGetDimension(dm, &dim));
 78:   PetscCall(PetscRandomSetInterval(rand, (PetscScalar)cStart, (PetscScalar)cEnd));
 79:   PetscCall(PetscMalloc1(n, &cells));
 80:   for (i = 0; i < n; i++) {
 81:     PetscReal rc;

 83:     PetscCall(PetscRandomGetValueReal(rand, &rc));
 84:     cells[i] = PetscFloorReal(rc);
 85:   }
 86:   PetscCall(ISCreateGeneral(comm, n, cells, PETSC_OWN_POINTER, &cellIS));
 87:   PetscCall(PetscObjectSetName((PetscObject)cellIS, "FE Test Cells"));
 88:   PetscCall(PetscQuadratureGetData(quad, NULL, NULL, &nq, NULL, NULL));
 89:   N = n * nq * nc;
 90:   PetscCall(PetscMalloc6(N, &B, N, &rB, N * dim, &D, N * dim, &rD, N * dim * dim, &H, N * dim * dim, &rH));
 91:   PetscCall(DMFieldEvaluateFE(field, cellIS, quad, PETSC_SCALAR, B, D, H));
 92:   PetscCall(DMFieldEvaluateFE(field, cellIS, quad, PETSC_REAL, rB, rD, rH));
 93:   viewer = PETSC_VIEWER_STDOUT_(comm);

 95:   PetscCall(PetscObjectSetName((PetscObject)quad, "Test quadrature"));
 96:   PetscCall(PetscQuadratureView(quad, viewer));
 97:   PetscCall(ISView(cellIS, viewer));
 98:   PetscCall(ViewResults(viewer, N, dim, B, D, H, rB, rD, rH));

100:   PetscCall(PetscFree6(B, rB, D, rD, H, rH));
101:   PetscCall(ISDestroy(&cellIS));
102:   PetscFunctionReturn(PETSC_SUCCESS);
103: }

105: static PetscErrorCode TestEvaluateFV(DMField field, PetscInt n, PetscInt cStart, PetscInt cEnd, PetscRandom rand)
106: {
107:   DM           dm;
108:   PetscInt     dim, i, nc;
109:   PetscInt     N;
110:   PetscScalar *B, *D, *H;
111:   PetscReal   *rB, *rD, *rH;
112:   PetscInt    *cells;
113:   IS           cellIS;
114:   PetscViewer  viewer;
115:   MPI_Comm     comm;

117:   PetscFunctionBegin;
118:   comm = PetscObjectComm((PetscObject)field);
119:   PetscCall(DMFieldGetNumComponents(field, &nc));
120:   PetscCall(DMFieldGetDM(field, &dm));
121:   PetscCall(DMGetDimension(dm, &dim));
122:   PetscCall(PetscRandomSetInterval(rand, (PetscScalar)cStart, (PetscScalar)cEnd));
123:   PetscCall(PetscMalloc1(n, &cells));
124:   for (i = 0; i < n; i++) {
125:     PetscReal rc;

127:     PetscCall(PetscRandomGetValueReal(rand, &rc));
128:     cells[i] = PetscFloorReal(rc);
129:   }
130:   PetscCall(ISCreateGeneral(comm, n, cells, PETSC_OWN_POINTER, &cellIS));
131:   PetscCall(PetscObjectSetName((PetscObject)cellIS, "FV Test Cells"));
132:   N = n * nc;
133:   PetscCall(PetscMalloc6(N, &B, N, &rB, N * dim, &D, N * dim, &rD, N * dim * dim, &H, N * dim * dim, &rH));
134:   PetscCall(DMFieldEvaluateFV(field, cellIS, PETSC_SCALAR, B, D, H));
135:   PetscCall(DMFieldEvaluateFV(field, cellIS, PETSC_REAL, rB, rD, rH));
136:   viewer = PETSC_VIEWER_STDOUT_(comm);

138:   PetscCall(ISView(cellIS, viewer));
139:   PetscCall(ViewResults(viewer, N, dim, B, D, H, rB, rD, rH));

141:   PetscCall(PetscFree6(B, rB, D, rD, H, rH));
142:   PetscCall(ISDestroy(&cellIS));
143:   PetscFunctionReturn(PETSC_SUCCESS);
144: }

146: static PetscErrorCode radiusSquared(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar u[], void *ctx)
147: {
148:   PetscInt  i;
149:   PetscReal r2 = 0.;

151:   PetscFunctionBegin;
152:   for (i = 0; i < dim; i++) r2 += PetscSqr(x[i]);
153:   for (i = 0; i < Nf; i++) u[i] = (i + 1) * r2;
154:   PetscFunctionReturn(PETSC_SUCCESS);
155: }

157: static PetscErrorCode TestShellEvaluate(DMField field, Vec points, PetscDataType type, void *B, void *D, void *H)
158: {
159:   Vec                ctxVec = NULL;
160:   const PetscScalar *mult;
161:   PetscInt           dim;
162:   const PetscScalar *x;
163:   PetscInt           Nc, n, i, j, k, l;

165:   PetscFunctionBegin;
166:   PetscCall(DMFieldGetNumComponents(field, &Nc));
167:   PetscCall(DMFieldShellGetContext(field, &ctxVec));
168:   PetscCall(VecGetBlockSize(points, &dim));
169:   PetscCall(VecGetLocalSize(points, &n));
170:   n /= Nc;
171:   PetscCall(VecGetArrayRead(ctxVec, &mult));
172:   PetscCall(VecGetArrayRead(points, &x));
173:   for (i = 0; i < n; i++) {
174:     PetscReal r2 = 0.;

176:     for (j = 0; j < dim; j++) r2 += PetscSqr(PetscRealPart(x[i * dim + j]));
177:     for (j = 0; j < Nc; j++) {
178:       PetscReal m = PetscRealPart(mult[j]);
179:       if (B) {
180:         if (type == PETSC_SCALAR) {
181:           ((PetscScalar *)B)[i * Nc + j] = m * r2;
182:         } else {
183:           ((PetscReal *)B)[i * Nc + j] = m * r2;
184:         }
185:       }
186:       if (D) {
187:         if (type == PETSC_SCALAR) {
188:           for (k = 0; k < dim; k++) ((PetscScalar *)D)[(i * Nc + j) * dim + k] = 2. * m * x[i * dim + k];
189:         } else {
190:           for (k = 0; k < dim; k++) ((PetscReal *)D)[(i * Nc + j) * dim + k] = 2. * m * PetscRealPart(x[i * dim + k]);
191:         }
192:       }
193:       if (H) {
194:         if (type == PETSC_SCALAR) {
195:           for (k = 0; k < dim; k++)
196:             for (l = 0; l < dim; l++) ((PetscScalar *)H)[((i * Nc + j) * dim + k) * dim + l] = (k == l) ? 2. * m : 0.;
197:         } else {
198:           for (k = 0; k < dim; k++)
199:             for (l = 0; l < dim; l++) ((PetscReal *)H)[((i * Nc + j) * dim + k) * dim + l] = (k == l) ? 2. * m : 0.;
200:         }
201:       }
202:     }
203:   }
204:   PetscCall(VecRestoreArrayRead(points, &x));
205:   PetscCall(VecRestoreArrayRead(ctxVec, &mult));
206:   PetscFunctionReturn(PETSC_SUCCESS);
207: }

209: static PetscErrorCode TestShellDestroy(DMField field)
210: {
211:   Vec ctxVec = NULL;

213:   PetscFunctionBegin;
214:   PetscCall(DMFieldShellGetContext(field, &ctxVec));
215:   PetscCall(VecDestroy(&ctxVec));
216:   PetscFunctionReturn(PETSC_SUCCESS);
217: }

219: int main(int argc, char **argv)
220: {
221:   DM              dm = NULL;
222:   MPI_Comm        comm;
223:   char            type[256] = DMPLEX;
224:   PetscBool       isda, isplex;
225:   PetscInt        dim    = 2;
226:   DMField         field  = NULL;
227:   PetscInt        nc     = 1;
228:   PetscInt        cStart = -1, cEnd = -1;
229:   PetscRandom     rand;
230:   PetscQuadrature quad          = NULL;
231:   PetscInt        pointsPerEdge = 2;
232:   PetscInt        numPoint = 0, numFE = 0, numFV = 0;
233:   PetscBool       testShell = PETSC_FALSE;

235:   PetscFunctionBeginUser;
236:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
237:   comm = PETSC_COMM_WORLD;
238:   PetscOptionsBegin(comm, "", "DMField Tutorial Options", "DM");
239:   PetscCall(PetscOptionsFList("-dm_type", "DM implementation on which to define field", "ex1.c", DMList, type, type, 256, NULL));
240:   PetscCall(PetscOptionsRangeInt("-dim", "DM intrinsic dimension", "ex1.c", dim, &dim, NULL, 1, 3));
241:   PetscCall(PetscOptionsBoundedInt("-num_components", "Number of components in field", "ex1.c", nc, &nc, NULL, 1));
242:   PetscCall(PetscOptionsBoundedInt("-num_quad_points", "Number of quadrature points per dimension", "ex1.c", pointsPerEdge, &pointsPerEdge, NULL, 1));
243:   PetscCall(PetscOptionsBoundedInt("-num_point_tests", "Number of test points for DMFieldEvaluate()", "ex1.c", numPoint, &numPoint, NULL, 0));
244:   PetscCall(PetscOptionsBoundedInt("-num_fe_tests", "Number of test cells for DMFieldEvaluateFE()", "ex1.c", numFE, &numFE, NULL, 0));
245:   PetscCall(PetscOptionsBoundedInt("-num_fv_tests", "Number of test cells for DMFieldEvaluateFV()", "ex1.c", numFV, &numFV, NULL, 0));
246:   PetscCall(PetscOptionsBool("-test_shell", "Test the DMFIELDSHELL implementation of DMField", "ex1.c", testShell, &testShell, NULL));
247:   PetscOptionsEnd();

249:   PetscCheck(dim <= 3, comm, PETSC_ERR_ARG_OUTOFRANGE, "This examples works for dim <= 3, not %" PetscInt_FMT, dim);
250:   PetscCall(PetscStrncmp(type, DMPLEX, 256, &isplex));
251:   PetscCall(PetscStrncmp(type, DMDA, 256, &isda));

253:   PetscCall(PetscRandomCreate(PETSC_COMM_SELF, &rand));
254:   PetscCall(PetscRandomSetFromOptions(rand));
255:   if (isplex) {
256:     PetscInt  overlap = 0;
257:     Vec       fieldvec;
258:     PetscInt  cells[3] = {3, 3, 3};
259:     PetscBool simplex;
260:     PetscFE   fe;

262:     PetscOptionsBegin(comm, "", "DMField DMPlex Options", "DM");
263:     PetscCall(PetscOptionsBoundedInt("-overlap", "DMPlex parallel overlap", "ex1.c", overlap, &overlap, NULL, 0));
264:     PetscOptionsEnd();
265:     if (0) {
266:       PetscCall(DMPlexCreateBoxMesh(comm, 2, PETSC_TRUE, cells, NULL, NULL, NULL, PETSC_TRUE, &dm));
267:     } else {
268:       PetscCall(DMCreate(comm, &dm));
269:       PetscCall(DMSetType(dm, DMPLEX));
270:       PetscCall(DMSetFromOptions(dm));
271:       CHKMEMQ;
272:     }
273:     PetscCall(DMGetDimension(dm, &dim));
274:     PetscCall(DMPlexIsSimplex(dm, &simplex));
275:     if (simplex) PetscCall(PetscDTStroudConicalQuadrature(dim, 1, pointsPerEdge, -1.0, 1.0, &quad));
276:     else PetscCall(PetscDTGaussTensorQuadrature(dim, 1, pointsPerEdge, -1.0, 1.0, &quad));
277:     PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd));
278:     if (testShell) {
279:       Vec          ctxVec;
280:       PetscInt     i;
281:       PetscScalar *array;

283:       PetscCall(VecCreateSeq(PETSC_COMM_SELF, nc, &ctxVec));
284:       PetscCall(VecSetUp(ctxVec));
285:       PetscCall(VecGetArray(ctxVec, &array));
286:       for (i = 0; i < nc; i++) array[i] = i + 1.;
287:       PetscCall(VecRestoreArray(ctxVec, &array));
288:       PetscCall(DMFieldCreateShell(dm, nc, DMFIELD_VERTEX, (void *)ctxVec, &field));
289:       PetscCall(DMFieldShellSetEvaluate(field, TestShellEvaluate));
290:       PetscCall(DMFieldShellSetDestroy(field, TestShellDestroy));
291:     } else {
292:       PetscCall(PetscFECreateDefault(PETSC_COMM_SELF, dim, nc, simplex, NULL, PETSC_DEFAULT, &fe));
293:       PetscCall(PetscFESetName(fe, "MyPetscFE"));
294:       PetscCall(DMSetField(dm, 0, NULL, (PetscObject)fe));
295:       PetscCall(PetscFEDestroy(&fe));
296:       PetscCall(DMCreateDS(dm));
297:       PetscCall(DMCreateLocalVector(dm, &fieldvec));
298:       {
299:         PetscErrorCode (*func[1])(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar *, void *);
300:         void *ctxs[1];

302:         func[0] = radiusSquared;
303:         ctxs[0] = NULL;

305:         PetscCall(DMProjectFunctionLocal(dm, 0.0, func, ctxs, INSERT_ALL_VALUES, fieldvec));
306:       }
307:       PetscCall(DMFieldCreateDS(dm, 0, fieldvec, &field));
308:       PetscCall(VecDestroy(&fieldvec));
309:     }
310:   } else if (isda) {
311:     PetscInt     i;
312:     PetscScalar *cv;

314:     switch (dim) {
315:     case 1:
316:       PetscCall(DMDACreate1d(comm, DM_BOUNDARY_NONE, 3, 1, 1, NULL, &dm));
317:       break;
318:     case 2:
319:       PetscCall(DMDACreate2d(comm, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, 3, 3, PETSC_DETERMINE, PETSC_DETERMINE, 1, 1, NULL, NULL, &dm));
320:       break;
321:     default:
322:       PetscCall(DMDACreate3d(comm, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, 3, 3, 3, PETSC_DETERMINE, PETSC_DETERMINE, PETSC_DETERMINE, 1, 1, NULL, NULL, NULL, &dm));
323:       break;
324:     }
325:     PetscCall(DMSetUp(dm));
326:     PetscCall(DMDAGetHeightStratum(dm, 0, &cStart, &cEnd));
327:     PetscCall(PetscMalloc1(nc * (1 << dim), &cv));
328:     for (i = 0; i < nc * (1 << dim); i++) {
329:       PetscReal rv;

331:       PetscCall(PetscRandomGetValueReal(rand, &rv));
332:       cv[i] = rv;
333:     }
334:     PetscCall(DMFieldCreateDA(dm, nc, cv, &field));
335:     PetscCall(PetscFree(cv));
336:     PetscCall(PetscDTGaussTensorQuadrature(dim, 1, pointsPerEdge, -1.0, 1.0, &quad));
337:   } else SETERRQ(comm, PETSC_ERR_SUP, "This test does not run for DM type %s", type);

339:   PetscCall(PetscObjectSetName((PetscObject)dm, "mesh"));
340:   PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
341:   PetscCall(DMDestroy(&dm));
342:   PetscCall(PetscObjectSetName((PetscObject)field, "field"));
343:   PetscCall(PetscObjectViewFromOptions((PetscObject)field, NULL, "-dmfield_view"));
344:   if (numPoint) PetscCall(TestEvaluate(field, numPoint, rand));
345:   if (numFE) PetscCall(TestEvaluateFE(field, numFE, cStart, cEnd, quad, rand));
346:   if (numFV) PetscCall(TestEvaluateFV(field, numFV, cStart, cEnd, rand));
347:   PetscCall(DMFieldDestroy(&field));
348:   PetscCall(PetscQuadratureDestroy(&quad));
349:   PetscCall(PetscRandomDestroy(&rand));
350:   PetscCall(PetscFinalize());
351:   return 0;
352: }

354: /*TEST

356:   test:
357:     suffix: da
358:     requires: !complex
359:     args: -dm_type da -dim 2 -num_components 2 -num_point_tests 2 -num_fe_tests 2 -num_fv_tests 2 -dmfield_view

361:   test:
362:     suffix: da_1
363:     requires: !complex
364:     args: -dm_type da -dim 1  -num_fe_tests 2

366:   test:
367:     suffix: da_2
368:     requires: !complex
369:     args: -dm_type da -dim 2  -num_fe_tests 2

371:   test:
372:     suffix: da_3
373:     requires: !complex
374:     args: -dm_type da -dim 3  -num_fe_tests 2

376:   test:
377:     suffix: ds
378:     requires: !complex triangle
379:     args: -dm_type plex -dm_plex_dim 2 -dm_plex_box_faces 3,3 -num_components 2 -num_point_tests 2 -num_fe_tests 2 -num_fv_tests 2 -dmfield_view -petscspace_degree 2 -num_quad_points 1

381:   test:
382:     suffix: ds_simplex_0
383:     requires: !complex triangle
384:     args: -dm_type plex -dm_plex_dim 2 -dm_plex_box_faces 3,3  -num_fe_tests 2  -petscspace_degree 0

386:   test:
387:     suffix: ds_simplex_1
388:     requires: !complex triangle
389:     args: -dm_type plex -dm_plex_dim 2 -dm_plex_box_faces 3,3  -num_fe_tests 2  -petscspace_degree 1

391:   test:
392:     suffix: ds_simplex_2
393:     requires: !complex triangle
394:     args: -dm_type plex -dm_plex_dim 2 -dm_plex_box_faces 3,3  -num_fe_tests 2  -petscspace_degree 2

396:   test:
397:     suffix: ds_tensor_2_0
398:     requires: !complex
399:     args: -dm_type plex -dm_plex_dim 2 -dm_plex_box_faces 3,3  -num_fe_tests 2  -petscspace_poly_tensor 1 -petscspace_degree 0 -dm_plex_simplex 0

401:   test:
402:     suffix: ds_tensor_2_1
403:     requires: !complex
404:     args: -dm_type plex -dm_plex_dim 2 -dm_plex_box_faces 3,3  -num_fe_tests 2  -petscspace_poly_tensor 1 -petscspace_degree 1 -dm_plex_simplex 0

406:   test:
407:     suffix: ds_tensor_2_2
408:     requires: !complex
409:     args: -dm_type plex -dm_plex_dim 2 -dm_plex_box_faces 3,3  -num_fe_tests 2  -petscspace_poly_tensor 1 -petscspace_degree 2 -dm_plex_simplex 0

411:   test:
412:     suffix: ds_tensor_3_0
413:     requires: !complex
414:     args: -dm_type plex -dm_plex_dim 3 -dm_plex_box_faces 3,3,3  -num_fe_tests 2  -petscspace_poly_tensor 1 -petscspace_degree 0 -dm_plex_simplex 0

416:   test:
417:     suffix: ds_tensor_3_1
418:     requires: !complex
419:     args: -dm_type plex -dm_plex_dim 3 -dm_plex_box_faces 3,3,3  -num_fe_tests 2  -petscspace_poly_tensor 1 -petscspace_degree 1 -dm_plex_simplex 0

421:   test:
422:     suffix: ds_tensor_3_2
423:     requires: !complex
424:     args: -dm_type plex -dm_plex_dim 3 -dm_plex_box_faces 3,3,3  -num_fe_tests 2  -petscspace_poly_tensor 1 -petscspace_degree 2 -dm_plex_simplex 0

426:   test:
427:     suffix: shell
428:     requires: !complex triangle
429:     args: -dm_coord_space 0 -dm_type plex -dm_plex_dim 2 -dm_plex_box_faces 3,3 -num_components 2 -num_point_tests 2 -num_fe_tests 2 -num_fv_tests 2 -dmfield_view -num_quad_points 1 -test_shell

431: TEST*/