xref: /petsc/src/dm/field/tutorials/ex1.c (revision c4762a1b19cd2af06abeed90e8f9d34fb975dd94) 
1 static char help[] = "Demonstration of creating and viewing DMFields objects.\n\n";
2 
3 #include <petscdmfield.h>
4 #include <petscdmplex.h>
5 #include <petscdmda.h>
6 
7 static PetscErrorCode ViewResults(PetscViewer viewer, PetscInt N, PetscInt dim, PetscScalar *B, PetscScalar *D, PetscScalar *H, PetscReal *rB, PetscReal *rD, PetscReal *rH)
8 {
9   PetscErrorCode ierr;
10 
11   PetscFunctionBegin;
12   ierr = PetscViewerASCIIPrintf(viewer,"B:\n");CHKERRQ(ierr);
13   ierr = PetscScalarView(N,B,viewer);CHKERRQ(ierr);
14   ierr = PetscViewerASCIIPrintf(viewer,"D:\n");CHKERRQ(ierr);
15   ierr = PetscScalarView(N*dim,D,viewer);CHKERRQ(ierr);
16   ierr = PetscViewerASCIIPrintf(viewer,"H:\n");CHKERRQ(ierr);
17   ierr = PetscScalarView(N*dim*dim,H,viewer);CHKERRQ(ierr);
18 
19   ierr = PetscViewerASCIIPrintf(viewer,"rB:\n");CHKERRQ(ierr);
20   ierr = PetscRealView(N,rB,viewer);CHKERRQ(ierr);
21   ierr = PetscViewerASCIIPrintf(viewer,"rD:\n");CHKERRQ(ierr);
22   ierr = PetscRealView(N*dim,rD,viewer);CHKERRQ(ierr);
23   ierr = PetscViewerASCIIPrintf(viewer,"rH:\n");CHKERRQ(ierr);
24   ierr = PetscRealView(N*dim*dim,rH,viewer);CHKERRQ(ierr);
25   PetscFunctionReturn(0);
26 }
27 
28 static PetscErrorCode TestEvaluate(DMField field, PetscInt n, PetscRandom rand)
29 {
30   DM             dm;
31   PetscInt       dim, i, nc;
32   PetscScalar    *B, *D, *H;
33   PetscReal      *rB, *rD, *rH;
34   Vec            points;
35   PetscScalar    *array;
36   PetscViewer    viewer;
37   MPI_Comm       comm;
38   PetscErrorCode ierr;
39 
40   PetscFunctionBegin;
41   comm = PetscObjectComm((PetscObject)field);
42   ierr = DMFieldGetNumComponents(field,&nc);CHKERRQ(ierr);
43   ierr = DMFieldGetDM(field,&dm);CHKERRQ(ierr);
44   ierr = DMGetDimension(dm,&dim);CHKERRQ(ierr);
45   ierr = VecCreateMPI(PetscObjectComm((PetscObject)field),n * dim,PETSC_DETERMINE,&points);CHKERRQ(ierr);
46   ierr = VecSetBlockSize(points,dim);CHKERRQ(ierr);
47   ierr = VecGetArray(points,&array);CHKERRQ(ierr);
48   for (i = 0; i < n * dim; i++) {ierr = PetscRandomGetValue(rand,&array[i]);CHKERRQ(ierr);}
49   ierr = VecRestoreArray(points,&array);CHKERRQ(ierr);
50   ierr = 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);CHKERRQ(ierr);
51   ierr = DMFieldEvaluate(field,points,PETSC_SCALAR,B,D,H);CHKERRQ(ierr);
52   ierr = DMFieldEvaluate(field,points,PETSC_REAL,rB,rD,rH);CHKERRQ(ierr);
53   viewer = PETSC_VIEWER_STDOUT_(comm);
54 
55   ierr = PetscObjectSetName((PetscObject)points,"Test Points");CHKERRQ(ierr);
56   ierr = VecView(points,viewer);CHKERRQ(ierr);
57   ierr = ViewResults(viewer,n*nc,dim,B,D,H,rB,rD,rH);CHKERRQ(ierr);
58 
59   ierr = PetscFree6(B,rB,D,rD,H,rH);CHKERRQ(ierr);
60   ierr = VecDestroy(&points);CHKERRQ(ierr);
61   PetscFunctionReturn(0);
62 }
63 
64 static PetscErrorCode TestEvaluateFE(DMField field, PetscInt n, PetscInt cStart, PetscInt cEnd, PetscQuadrature quad, PetscRandom rand)
65 {
66   DM             dm;
67   PetscInt       dim, i, nc, nq;
68   PetscInt       N;
69   PetscScalar    *B, *D, *H;
70   PetscReal      *rB, *rD, *rH;
71   PetscInt       *cells;
72   IS             cellIS;
73   PetscViewer    viewer;
74   MPI_Comm       comm;
75   PetscErrorCode ierr;
76 
77   PetscFunctionBegin;
78   comm = PetscObjectComm((PetscObject)field);
79   ierr = DMFieldGetNumComponents(field,&nc);CHKERRQ(ierr);
80   ierr = DMFieldGetDM(field,&dm);CHKERRQ(ierr);
81   ierr = DMGetDimension(dm,&dim);CHKERRQ(ierr);
82   ierr = PetscRandomSetInterval(rand,(PetscScalar) cStart, (PetscScalar) cEnd);CHKERRQ(ierr);
83   ierr = PetscMalloc1(n,&cells);CHKERRQ(ierr);
84   for (i = 0; i < n; i++) {
85     PetscReal rc;
86 
87     ierr = PetscRandomGetValueReal(rand,&rc);CHKERRQ(ierr);
88     cells[i] = PetscFloorReal(rc);
89   }
90   ierr = ISCreateGeneral(comm,n,cells,PETSC_OWN_POINTER,&cellIS);CHKERRQ(ierr);
91   ierr = PetscObjectSetName((PetscObject)cellIS,"FE Test Cells");CHKERRQ(ierr);
92   ierr = PetscQuadratureGetData(quad,NULL,NULL,&nq,NULL,NULL);CHKERRQ(ierr);
93   N    = n * nq * nc;
94   ierr = PetscMalloc6(N,&B,N,&rB,N*dim,&D,N*dim,&rD,N*dim*dim,&H,N*dim*dim,&rH);CHKERRQ(ierr);
95   ierr = DMFieldEvaluateFE(field,cellIS,quad,PETSC_SCALAR,B,D,H);CHKERRQ(ierr);
96   ierr = DMFieldEvaluateFE(field,cellIS,quad,PETSC_REAL,rB,rD,rH);CHKERRQ(ierr);
97   viewer = PETSC_VIEWER_STDOUT_(comm);
98 
99   ierr = PetscObjectSetName((PetscObject)quad,"Test quadrature");CHKERRQ(ierr);
100   ierr = PetscQuadratureView(quad,viewer);CHKERRQ(ierr);
101   ierr = ISView(cellIS,viewer);CHKERRQ(ierr);
102   ierr = ViewResults(viewer,N,dim,B,D,H,rB,rD,rH);CHKERRQ(ierr);
103 
104   ierr = PetscFree6(B,rB,D,rD,H,rH);CHKERRQ(ierr);
105   ierr = ISDestroy(&cellIS);CHKERRQ(ierr);
106   PetscFunctionReturn(0);
107 }
108 
109 static PetscErrorCode TestEvaluateFV(DMField field, PetscInt n, PetscInt cStart, PetscInt cEnd, PetscRandom rand)
110 {
111   DM             dm;
112   PetscInt       dim, i, nc;
113   PetscInt       N;
114   PetscScalar    *B, *D, *H;
115   PetscReal      *rB, *rD, *rH;
116   PetscInt       *cells;
117   IS             cellIS;
118   PetscViewer    viewer;
119   MPI_Comm       comm;
120   PetscErrorCode ierr;
121 
122   PetscFunctionBegin;
123   comm = PetscObjectComm((PetscObject)field);
124   ierr = DMFieldGetNumComponents(field,&nc);CHKERRQ(ierr);
125   ierr = DMFieldGetDM(field,&dm);CHKERRQ(ierr);
126   ierr = DMGetDimension(dm,&dim);CHKERRQ(ierr);
127   ierr = PetscRandomSetInterval(rand,(PetscScalar) cStart, (PetscScalar) cEnd);CHKERRQ(ierr);
128   ierr = PetscMalloc1(n,&cells);CHKERRQ(ierr);
129   for (i = 0; i < n; i++) {
130     PetscReal rc;
131 
132     ierr = PetscRandomGetValueReal(rand,&rc);CHKERRQ(ierr);
133     cells[i] = PetscFloorReal(rc);
134   }
135   ierr = ISCreateGeneral(comm,n,cells,PETSC_OWN_POINTER,&cellIS);CHKERRQ(ierr);
136   ierr = PetscObjectSetName((PetscObject)cellIS,"FV Test Cells");CHKERRQ(ierr);
137   N    = n * nc;
138   ierr = PetscMalloc6(N,&B,N,&rB,N*dim,&D,N*dim,&rD,N*dim*dim,&H,N*dim*dim,&rH);CHKERRQ(ierr);
139   ierr = DMFieldEvaluateFV(field,cellIS,PETSC_SCALAR,B,D,H);CHKERRQ(ierr);
140   ierr = DMFieldEvaluateFV(field,cellIS,PETSC_REAL,rB,rD,rH);CHKERRQ(ierr);
141   viewer = PETSC_VIEWER_STDOUT_(comm);
142 
143   ierr = ISView(cellIS,viewer);CHKERRQ(ierr);
144   ierr = ViewResults(viewer,N,dim,B,D,H,rB,rD,rH);CHKERRQ(ierr);
145 
146   ierr = PetscFree6(B,rB,D,rD,H,rH);CHKERRQ(ierr);
147   ierr = ISDestroy(&cellIS);CHKERRQ(ierr);
148   PetscFunctionReturn(0);
149 }
150 
151 static PetscErrorCode radiusSquared(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar u[], void *ctx)
152 {
153   PetscInt       i;
154   PetscReal      r2 = 0.;
155 
156   PetscFunctionBegin;
157   for (i = 0; i < dim; i++) {r2 += PetscSqr(x[i]);}
158   for (i = 0; i < Nf; i++) {
159     u[i] = (i + 1) * r2;
160   }
161   PetscFunctionReturn(0);
162 }
163 
164 static PetscErrorCode TestShellEvaluate(DMField field, Vec points, PetscDataType type, void *B, void *D, void *H)
165 {
166   Vec                ctxVec = NULL;
167   const PetscScalar *mult;
168   PetscInt           dim;
169   const PetscScalar *x;
170   PetscInt           Nc, n, i, j, k, l;
171   PetscErrorCode     ierr;
172 
173   PetscFunctionBegin;
174   ierr = DMFieldGetNumComponents(field, &Nc);CHKERRQ(ierr);
175   ierr = DMFieldShellGetContext(field, (void *) &ctxVec);CHKERRQ(ierr);
176   ierr = VecGetBlockSize(points, &dim);CHKERRQ(ierr);
177   ierr = VecGetLocalSize(points, &n);CHKERRQ(ierr);
178   n /= Nc;
179   ierr = VecGetArrayRead(ctxVec, &mult);CHKERRQ(ierr);
180   ierr = VecGetArrayRead(points, &x);CHKERRQ(ierr);
181   for (i = 0; i < n; i++) {
182     PetscReal r2 = 0.;
183 
184     for (j = 0; j < dim; j++) {r2 += PetscSqr(PetscRealPart(x[i * dim + j]));}
185     for (j = 0; j < Nc; j++) {
186       PetscReal m = PetscRealPart(mult[j]);
187       if (B) {
188         if (type == PETSC_SCALAR) {
189           ((PetscScalar *)B)[i * Nc + j] = m * r2;
190         } else {
191           ((PetscReal *)B)[i * Nc + j] = m * r2;
192         }
193       }
194       if (D) {
195         if (type == PETSC_SCALAR) {
196           for (k = 0; k < dim; k++) ((PetscScalar *)D)[(i * Nc + j) * dim + k] = 2. * m * x[i * dim + k];
197         } else {
198           for (k = 0; k < dim; k++) ((PetscReal   *)D)[(i * Nc + j) * dim + k] = 2. * m * PetscRealPart(x[i * dim + k]);
199         }
200       }
201       if (H) {
202         if (type == PETSC_SCALAR) {
203           for (k = 0; k < dim; k++) for (l = 0; l < dim; l++) ((PetscScalar *)H)[((i * Nc + j) * dim + k) * dim + l] = (k == l) ? 2. * m : 0.;
204         } else {
205           for (k = 0; k < dim; k++) for (l = 0; l < dim; l++) ((PetscReal   *)H)[((i * Nc + j) * dim + k) * dim + l] = (k == l) ? 2. * m : 0.;
206         }
207       }
208     }
209   }
210   ierr = VecRestoreArrayRead(points, &x);CHKERRQ(ierr);
211   ierr = VecRestoreArrayRead(ctxVec, &mult);CHKERRQ(ierr);
212   PetscFunctionReturn(0);
213 }
214 
215 static PetscErrorCode TestShellDestroy(DMField field)
216 {
217   Vec                ctxVec = NULL;
218   PetscErrorCode     ierr;
219 
220   PetscFunctionBegin;
221   ierr = DMFieldShellGetContext(field, (void *) &ctxVec);CHKERRQ(ierr);
222   ierr = VecDestroy(&ctxVec);CHKERRQ(ierr);
223   PetscFunctionReturn(0);
224 }
225 
226 int main(int argc, char **argv)
227 {
228   DM              dm = NULL;
229   MPI_Comm        comm;
230   char            type[256] = DMPLEX;
231   PetscBool       isda, isplex;
232   PetscInt        dim = 2;
233   DMField         field = NULL;
234   PetscInt        nc = 1;
235   PetscInt        cStart = -1, cEnd = -1;
236   PetscRandom     rand;
237   PetscQuadrature quad = NULL;
238   PetscInt        pointsPerEdge = 2;
239   PetscInt        numPoint = 0, numFE = 0, numFV = 0;
240   PetscBool       testShell = PETSC_FALSE;
241   PetscErrorCode  ierr;
242 
243   ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr;
244   comm = PETSC_COMM_WORLD;
245   ierr = PetscOptionsBegin(comm, "", "DMField Tutorial Options", "DM");CHKERRQ(ierr);
246   ierr = PetscOptionsFList("-dm_type","DM implementation on which to define field","ex1.c",DMList,type,type,256,NULL);CHKERRQ(ierr);
247   ierr = PetscOptionsRangeInt("-dim","DM intrinsic dimension", "ex1.c", dim, &dim, NULL,1,3);CHKERRQ(ierr);
248   ierr = PetscOptionsBoundedInt("-num_components","Number of components in field", "ex1.c", nc, &nc, NULL,1);CHKERRQ(ierr);
249   ierr = PetscOptionsBoundedInt("-num_quad_points","Number of quadrature points per dimension", "ex1.c", pointsPerEdge, &pointsPerEdge, NULL,1);CHKERRQ(ierr);
250   ierr = PetscOptionsBoundedInt("-num_point_tests", "Number of test points for DMFieldEvaluate()", "ex1.c", numPoint, &numPoint, NULL,0);CHKERRQ(ierr);
251   ierr = PetscOptionsBoundedInt("-num_fe_tests", "Number of test cells for DMFieldEvaluateFE()", "ex1.c", numFE, &numFE, NULL,0);CHKERRQ(ierr);
252   ierr = PetscOptionsBoundedInt("-num_fv_tests", "Number of test cells for DMFieldEvaluateFV()", "ex1.c", numFV, &numFV, NULL,0);CHKERRQ(ierr);
253   ierr = PetscOptionsBool("-test_shell", "Test the DMFIELDSHELL implementation of DMField", "ex1.c", testShell, &testShell, NULL);CHKERRQ(ierr);
254   ierr = PetscOptionsEnd();CHKERRQ(ierr);
255 
256   if (dim > 3) SETERRQ1(comm,PETSC_ERR_ARG_OUTOFRANGE,"This examples works for dim <= 3, not %D",dim);
257   ierr = PetscStrncmp(type,DMPLEX,256,&isplex);CHKERRQ(ierr);
258   ierr = PetscStrncmp(type,DMDA,256,&isda);CHKERRQ(ierr);
259 
260   ierr = PetscRandomCreate(PETSC_COMM_SELF,&rand);CHKERRQ(ierr);
261   ierr = PetscRandomSetFromOptions(rand);CHKERRQ(ierr);
262   if (isplex) {
263     PetscBool simplex = PETSC_TRUE;
264     PetscInt  overlap = 0;
265     Vec       fieldvec;
266     PetscBool interpolate = PETSC_TRUE;
267     PetscInt       inttypes[3];
268     DMBoundaryType types[3] = {DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE};
269     PetscInt       cells[3] = {3,3,3};
270     PetscBool      flags[3] = {PETSC_FALSE, PETSC_FALSE, PETSC_FALSE};
271     PetscInt       n = 3, i;
272     PetscFE   fe;
273 
274     ierr = PetscOptionsBegin(comm, "", "DMField DMPlex Options", "DM");CHKERRQ(ierr);
275     ierr = PetscOptionsBool("-simplex","Create a simplicial DMPlex","ex1.c",simplex,&simplex,NULL);CHKERRQ(ierr);
276     ierr = PetscOptionsBoundedInt("-overlap","DMPlex parallel overlap","ex1.c",overlap,&overlap,NULL,0);CHKERRQ(ierr);
277     ierr = PetscOptionsBool("-interpolate","Interpolate the DMPlex","ex1.c",interpolate,&interpolate,NULL);CHKERRQ(ierr);
278     ierr = PetscOptionsEList("-boundary_x","type of boundary in x direction","ex1.c",DMBoundaryTypes,5,DMBoundaryTypes[types[0]],&inttypes[0],&flags[0]);CHKERRQ(ierr);
279     ierr = PetscOptionsEList("-boundary_y","type of boundary in y direction","ex1.c",DMBoundaryTypes,5,DMBoundaryTypes[types[1]],&inttypes[1],&flags[1]);CHKERRQ(ierr);
280     ierr = PetscOptionsEList("-boundary_z","type of boundary in z direction","ex1.c",DMBoundaryTypes,5,DMBoundaryTypes[types[2]],&inttypes[2],&flags[2]);CHKERRQ(ierr);
281     ierr = PetscOptionsIntArray("-cells","Cells per dimension","ex1.c",cells,&n,NULL);CHKERRQ(ierr);
282     ierr = PetscOptionsEnd();CHKERRQ(ierr);
283     for (i = 0; i < 3; i++) {
284       if (flags[i]) {types[i] = (DMBoundaryType) inttypes[i];}
285     }
286     ierr = DMPlexCreateBoxMesh(comm,dim,simplex,cells,NULL,NULL,types,interpolate,&dm);CHKERRQ(ierr);
287     if (simplex) {
288       ierr = PetscDTStroudConicalQuadrature(dim, 1, pointsPerEdge, -1.0, 1.0, &quad);CHKERRQ(ierr);
289     } else {
290       ierr = PetscDTGaussTensorQuadrature(dim, 1, pointsPerEdge, -1.0, 1.0, &quad);CHKERRQ(ierr);
291     }
292     {
293       PetscPartitioner part;
294       DM               dmDist;
295 
296       ierr = DMPlexGetPartitioner(dm,&part);CHKERRQ(ierr);
297       ierr = PetscPartitionerSetType(part,PETSCPARTITIONERSIMPLE);CHKERRQ(ierr);
298       ierr = DMPlexDistribute(dm,overlap,NULL,&dmDist);CHKERRQ(ierr);
299       if (dmDist) {
300         ierr = DMDestroy(&dm);CHKERRQ(ierr);
301         dm = dmDist;
302       }
303     }
304     ierr = DMPlexGetHeightStratum(dm,0,&cStart,&cEnd);CHKERRQ(ierr);
305     if (testShell) {
306       Vec ctxVec;
307       PetscInt i;
308       PetscScalar *array;
309 
310       ierr = VecCreateSeq(PETSC_COMM_SELF, nc, &ctxVec);CHKERRQ(ierr);
311       ierr = VecSetUp(ctxVec);CHKERRQ(ierr);
312       ierr = VecGetArray(ctxVec,&array);CHKERRQ(ierr);
313       for (i = 0; i < nc; i++) array[i] = i + 1.;
314       ierr = VecRestoreArray(ctxVec,&array);CHKERRQ(ierr);
315       ierr = DMFieldCreateShell(dm, nc, DMFIELD_VERTEX, (void *) ctxVec, &field);CHKERRQ(ierr);
316       ierr = DMFieldShellSetEvaluate(field, TestShellEvaluate);CHKERRQ(ierr);
317       ierr = DMFieldShellSetDestroy(field, TestShellDestroy);CHKERRQ(ierr);
318     } else {
319       ierr = PetscFECreateDefault(comm,dim,nc,simplex,NULL,PETSC_DEFAULT,&fe);CHKERRQ(ierr);
320       ierr = PetscFESetName(fe,"MyPetscFE");CHKERRQ(ierr);
321       ierr = DMSetField(dm,0,NULL,(PetscObject)fe);CHKERRQ(ierr);
322       ierr = DMCreateDS(dm);CHKERRQ(ierr);
323       ierr = PetscFEDestroy(&fe);CHKERRQ(ierr);
324       ierr = DMCreateLocalVector(dm,&fieldvec);CHKERRQ(ierr);
325       {
326         PetscErrorCode (*func[1]) (PetscInt,PetscReal,const PetscReal [],PetscInt, PetscScalar *,void *);
327         void            *ctxs[1];
328 
329         func[0] = radiusSquared;
330         ctxs[0] = NULL;
331 
332         ierr = DMProjectFunctionLocal(dm,0.0,func,ctxs,INSERT_ALL_VALUES,fieldvec);CHKERRQ(ierr);
333       }
334       ierr = DMFieldCreateDS(dm,0,fieldvec,&field);CHKERRQ(ierr);
335       ierr = VecDestroy(&fieldvec);CHKERRQ(ierr);
336     }
337   } else if (isda) {
338     PetscInt       i;
339     PetscScalar    *cv;
340 
341     switch (dim) {
342     case 1:
343       ierr = DMDACreate1d(comm, DM_BOUNDARY_NONE, 3, 1, 1, NULL, &dm);CHKERRQ(ierr);
344       break;
345     case 2:
346       ierr = DMDACreate2d(comm, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, 3, 3, PETSC_DETERMINE, PETSC_DETERMINE, 1, 1, NULL, NULL, &dm);CHKERRQ(ierr);
347       break;
348     default:
349       ierr = 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);CHKERRQ(ierr);
350       break;
351     }
352     ierr = DMSetUp(dm);CHKERRQ(ierr);
353     ierr = DMDAGetHeightStratum(dm,0,&cStart,&cEnd);CHKERRQ(ierr);
354     ierr = PetscMalloc1(nc * (1 << dim),&cv);CHKERRQ(ierr);
355     for (i = 0; i < nc * (1 << dim); i++) {
356       PetscReal rv;
357 
358       ierr = PetscRandomGetValueReal(rand,&rv);CHKERRQ(ierr);
359       cv[i] = rv;
360     }
361     ierr = DMFieldCreateDA(dm,nc,cv,&field);CHKERRQ(ierr);
362     ierr = PetscFree(cv);CHKERRQ(ierr);
363     ierr = PetscDTGaussTensorQuadrature(dim, 1, pointsPerEdge, -1.0, 1.0, &quad);CHKERRQ(ierr);
364   } else SETERRQ1(comm,PETSC_ERR_SUP,"This test does not run for DM type %s",type);
365 
366   ierr = PetscObjectSetName((PetscObject)dm,"mesh");CHKERRQ(ierr);
367   ierr = DMViewFromOptions(dm,NULL,"-dm_view");CHKERRQ(ierr);
368   ierr = DMDestroy(&dm);CHKERRQ(ierr);
369   ierr = PetscObjectSetName((PetscObject)field,"field");CHKERRQ(ierr);
370   ierr = PetscObjectViewFromOptions((PetscObject)field,NULL,"-dmfield_view");CHKERRQ(ierr);
371   if (numPoint) {ierr = TestEvaluate(field,numPoint,rand);CHKERRQ(ierr);}
372   if (numFE) {ierr = TestEvaluateFE(field,numFE,cStart,cEnd,quad,rand);CHKERRQ(ierr);}
373   if (numFV) {ierr = TestEvaluateFV(field,numFV,cStart,cEnd,rand);CHKERRQ(ierr);}
374   ierr = PetscQuadratureDestroy(&quad);CHKERRQ(ierr);
375   ierr = DMFieldDestroy(&field);CHKERRQ(ierr);
376   ierr = PetscRandomDestroy(&rand);CHKERRQ(ierr);
377   ierr = PetscFinalize();
378   return ierr;
379 }
380 
381 /*TEST
382 
383   test:
384     suffix: da
385     requires: !complex
386     args: -dm_type da -dim 2 -num_components 2 -num_point_tests 2 -num_fe_tests 2 -num_fv_tests 2 -dmfield_view
387 
388   test:
389     suffix: da_1
390     requires: !complex
391     args: -dm_type da -dim 1  -num_fe_tests 2
392 
393   test:
394     suffix: da_2
395     requires: !complex
396     args: -dm_type da -dim 2  -num_fe_tests 2
397 
398   test:
399     suffix: da_3
400     requires: !complex
401     args: -dm_type da -dim 3  -num_fe_tests 2
402 
403   test:
404     suffix: ds
405     requires: !complex triangle
406     args: -dm_type plex -dim 2 -num_components 2 -num_point_tests 2 -num_fe_tests 2 -num_fv_tests 2 -dmfield_view -petscspace_degree 2 -num_quad_points 1
407 
408   test:
409     suffix: ds_simplex_0
410     requires: !complex triangle
411     args: -dm_type plex -dim 2  -num_fe_tests 2  -petscspace_degree 0
412 
413   test:
414     suffix: ds_simplex_1
415     requires: !complex triangle
416     args: -dm_type plex -dim 2  -num_fe_tests 2  -petscspace_degree 1
417 
418   test:
419     suffix: ds_simplex_2
420     requires: !complex triangle
421     args: -dm_type plex -dim 2  -num_fe_tests 2  -petscspace_degree 2
422 
423   test:
424     suffix: ds_tensor_2_0
425     requires: !complex
426     args: -dm_type plex -dim 2  -num_fe_tests 2  -petscspace_poly_tensor 1 -petscspace_degree 0 -simplex 0
427 
428   test:
429     suffix: ds_tensor_2_1
430     requires: !complex
431     args: -dm_type plex -dim 2  -num_fe_tests 2  -petscspace_poly_tensor 1 -petscspace_degree 1 -simplex 0
432 
433   test:
434     suffix: ds_tensor_2_2
435     requires: !complex
436     args: -dm_type plex -dim 2  -num_fe_tests 2  -petscspace_poly_tensor 1 -petscspace_degree 2 -simplex 0
437 
438   test:
439     suffix: ds_tensor_3_0
440     requires: !complex
441     args: -dm_type plex -dim 3  -num_fe_tests 2  -petscspace_poly_tensor 1 -petscspace_degree 0 -simplex 0
442 
443   test:
444     suffix: ds_tensor_3_1
445     requires: !complex
446     args: -dm_type plex -dim 3  -num_fe_tests 2  -petscspace_poly_tensor 1 -petscspace_degree 1 -simplex 0
447 
448   test:
449     suffix: ds_tensor_3_2
450     requires: !complex
451     args: -dm_type plex -dim 3  -num_fe_tests 2  -petscspace_poly_tensor 1 -petscspace_degree 2 -simplex 0
452 
453   test:
454     suffix: shell
455     requires: !complex triangle
456     args: -dm_type plex -dim 2 -num_components 2 -num_point_tests 2 -num_fe_tests 2 -num_fv_tests 2 -dmfield_view -num_quad_points 1 -test_shell
457 
458 TEST*/
459