1 static char help[] = "Tests for cell geometry\n\n"; 2 3 #include <petscdmplex.h> 4 5 typedef enum {RUN_REFERENCE, RUN_FILE, RUN_DISPLAY} RunType; 6 7 typedef struct { 8 DM dm; 9 RunType runType; /* Type of mesh to use */ 10 char filename[PETSC_MAX_PATH_LEN]; /* Import mesh from file */ 11 PetscBool interpolate; /* Interpolate the mesh */ 12 PetscBool transform; /* Use random coordinate transformations */ 13 /* Data for input meshes */ 14 PetscReal *v0, *J, *invJ, *detJ; /* FEM data */ 15 PetscReal *centroid, *normal, *vol; /* FVM data */ 16 } AppCtx; 17 18 static PetscErrorCode ReadMesh(MPI_Comm comm, const char *filename, AppCtx *user, DM *dm) 19 { 20 PetscErrorCode ierr; 21 22 PetscFunctionBegin; 23 ierr = DMPlexCreateFromFile(comm, filename, user->interpolate, dm);CHKERRQ(ierr); 24 ierr = PetscObjectSetName((PetscObject) *dm, "Input Mesh");CHKERRQ(ierr); 25 ierr = DMViewFromOptions(*dm, NULL, "-dm_view");CHKERRQ(ierr); 26 PetscFunctionReturn(0); 27 } 28 29 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 30 { 31 const char *runTypes[3] = {"reference", "file", "display"}; 32 PetscInt run; 33 PetscErrorCode ierr; 34 35 PetscFunctionBeginUser; 36 options->runType = RUN_REFERENCE; 37 options->filename[0] = '\0'; 38 options->interpolate = PETSC_FALSE; 39 options->transform = PETSC_FALSE; 40 41 ierr = PetscOptionsBegin(comm, "", "Geometry Test Options", "DMPLEX");CHKERRQ(ierr); 42 run = options->runType; 43 ierr = PetscOptionsEList("-run_type", "The run type", "ex8.c", runTypes, 3, runTypes[options->runType], &run, NULL);CHKERRQ(ierr); 44 options->runType = (RunType) run; 45 ierr = PetscOptionsString("-filename", "The mesh file", "ex8.c", options->filename, options->filename, PETSC_MAX_PATH_LEN, NULL);CHKERRQ(ierr); 46 ierr = PetscOptionsBool("-interpolate", "Interpolate the mesh", "ex8.c", options->interpolate, &options->interpolate, NULL);CHKERRQ(ierr); 47 ierr = PetscOptionsBool("-transform", "Use random transforms", "ex8.c", options->transform, &options->transform, NULL);CHKERRQ(ierr); 48 49 if (options->runType == RUN_FILE) { 50 PetscInt dim, cStart, cEnd, numCells, n; 51 PetscBool flag; 52 53 ierr = ReadMesh(PETSC_COMM_WORLD, options->filename, options, &options->dm);CHKERRQ(ierr); 54 ierr = DMGetDimension(options->dm, &dim);CHKERRQ(ierr); 55 ierr = DMPlexGetHeightStratum(options->dm, 0, &cStart, &cEnd);CHKERRQ(ierr); 56 numCells = cEnd-cStart; 57 ierr = PetscMalloc7(numCells*dim,&options->v0,numCells*dim*dim,&options->J,numCells*dim*dim,&options->invJ,numCells,&options->detJ, 58 numCells*dim,&options->centroid,numCells*dim,&options->normal,numCells,&options->vol);CHKERRQ(ierr); 59 n = numCells*dim; 60 ierr = PetscOptionsRealArray("-v0", "Input v0 for each cell", "ex8.c", options->v0, &n, &flag);CHKERRQ(ierr); 61 if (flag && n != numCells*dim) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid size of v0 %D should be %D", n, numCells*dim); 62 n = numCells*dim*dim; 63 ierr = PetscOptionsRealArray("-J", "Input Jacobian for each cell", "ex8.c", options->J, &n, &flag);CHKERRQ(ierr); 64 if (flag && n != numCells*dim*dim) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid size of J %D should be %D", n, numCells*dim*dim); 65 n = numCells*dim*dim; 66 ierr = PetscOptionsRealArray("-invJ", "Input inverse Jacobian for each cell", "ex8.c", options->invJ, &n, &flag);CHKERRQ(ierr); 67 if (flag && n != numCells*dim*dim) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid size of invJ %D should be %D", n, numCells*dim*dim); 68 n = numCells; 69 ierr = PetscOptionsRealArray("-detJ", "Input Jacobian determinant for each cell", "ex8.c", options->detJ, &n, &flag);CHKERRQ(ierr); 70 if (flag && n != numCells) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid size of detJ %D should be %D", n, numCells); 71 n = numCells*dim; 72 ierr = PetscOptionsRealArray("-centroid", "Input centroid for each cell", "ex8.c", options->centroid, &n, &flag);CHKERRQ(ierr); 73 if (flag && n != numCells*dim) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid size of centroid %D should be %D", n, numCells*dim); 74 n = numCells*dim; 75 ierr = PetscOptionsRealArray("-normal", "Input normal for each cell", "ex8.c", options->normal, &n, &flag);CHKERRQ(ierr); 76 if (flag && n != numCells*dim) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid size of normal %D should be %D", n, numCells*dim); 77 n = numCells; 78 ierr = PetscOptionsRealArray("-vol", "Input volume for each cell", "ex8.c", options->vol, &n, &flag);CHKERRQ(ierr); 79 if (flag && n != numCells) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid size of vol %D should be %D", n, numCells); 80 } else if (options->runType == RUN_DISPLAY) { 81 ierr = ReadMesh(PETSC_COMM_WORLD, options->filename, options, &options->dm);CHKERRQ(ierr); 82 } 83 ierr = PetscOptionsEnd(); 84 85 if (options->transform) {ierr = PetscPrintf(comm, "Using random transforms\n");CHKERRQ(ierr);} 86 PetscFunctionReturn(0); 87 } 88 89 static PetscErrorCode ChangeCoordinates(DM dm, PetscInt spaceDim, PetscScalar vertexCoords[]) 90 { 91 PetscSection coordSection; 92 Vec coordinates; 93 PetscScalar *coords; 94 PetscInt vStart, vEnd, v, d, coordSize; 95 PetscErrorCode ierr; 96 97 PetscFunctionBegin; 98 ierr = DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd);CHKERRQ(ierr); 99 ierr = DMGetCoordinateSection(dm, &coordSection);CHKERRQ(ierr); 100 ierr = PetscSectionSetNumFields(coordSection, 1);CHKERRQ(ierr); 101 ierr = PetscSectionSetFieldComponents(coordSection, 0, spaceDim);CHKERRQ(ierr); 102 ierr = PetscSectionSetChart(coordSection, vStart, vEnd);CHKERRQ(ierr); 103 for (v = vStart; v < vEnd; ++v) { 104 ierr = PetscSectionSetDof(coordSection, v, spaceDim);CHKERRQ(ierr); 105 ierr = PetscSectionSetFieldDof(coordSection, v, 0, spaceDim);CHKERRQ(ierr); 106 } 107 ierr = PetscSectionSetUp(coordSection);CHKERRQ(ierr); 108 ierr = PetscSectionGetStorageSize(coordSection, &coordSize);CHKERRQ(ierr); 109 ierr = VecCreate(PETSC_COMM_SELF, &coordinates);CHKERRQ(ierr); 110 ierr = PetscObjectSetName((PetscObject) coordinates, "coordinates");CHKERRQ(ierr); 111 ierr = VecSetSizes(coordinates, coordSize, PETSC_DETERMINE);CHKERRQ(ierr); 112 ierr = VecSetFromOptions(coordinates);CHKERRQ(ierr); 113 ierr = VecGetArray(coordinates, &coords);CHKERRQ(ierr); 114 for (v = vStart; v < vEnd; ++v) { 115 PetscInt off; 116 117 ierr = PetscSectionGetOffset(coordSection, v, &off);CHKERRQ(ierr); 118 for (d = 0; d < spaceDim; ++d) { 119 coords[off+d] = vertexCoords[(v-vStart)*spaceDim+d]; 120 } 121 } 122 ierr = VecRestoreArray(coordinates, &coords);CHKERRQ(ierr); 123 ierr = DMSetCoordinateDim(dm, spaceDim);CHKERRQ(ierr); 124 ierr = DMSetCoordinatesLocal(dm, coordinates);CHKERRQ(ierr); 125 ierr = VecDestroy(&coordinates);CHKERRQ(ierr); 126 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 127 PetscFunctionReturn(0); 128 } 129 130 #define RelativeError(a,b) PetscAbs(a-b)/(1.0+PetscMax(PetscAbs(a),PetscAbs(b))) 131 132 static PetscErrorCode CheckFEMGeometry(DM dm, PetscInt cell, PetscInt spaceDim, PetscReal v0Ex[], PetscReal JEx[], PetscReal invJEx[], PetscReal detJEx) 133 { 134 PetscReal v0[3], J[9], invJ[9], detJ; 135 PetscInt d, i, j; 136 PetscErrorCode ierr; 137 138 PetscFunctionBegin; 139 ierr = DMPlexComputeCellGeometryFEM(dm, cell, NULL, v0, J, invJ, &detJ);CHKERRQ(ierr); 140 for (d = 0; d < spaceDim; ++d) { 141 if (v0[d] != v0Ex[d]) { 142 switch (spaceDim) { 143 case 2: SETERRQ4(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Invalid v0 (%g, %g) != (%g, %g)", (double)v0[0], (double)v0[1], (double)v0Ex[0], (double)v0Ex[1]);break; 144 case 3: SETERRQ6(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Invalid v0 (%g, %g, %g) != (%g, %g, %g)", (double)v0[0], (double)v0[1], (double)v0[2], (double)v0Ex[0], (double)v0Ex[1], (double)v0Ex[2]);break; 145 default: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Invalid space dimension %D", spaceDim); 146 } 147 } 148 } 149 for (i = 0; i < spaceDim; ++i) { 150 for (j = 0; j < spaceDim; ++j) { 151 if (RelativeError(J[i*spaceDim+j],JEx[i*spaceDim+j]) > 10*PETSC_SMALL) SETERRQ4(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Invalid J[%D,%D]: %g != %g", i, j, (double)J[i*spaceDim+j], (double)JEx[i*spaceDim+j]); 152 if (RelativeError(invJ[i*spaceDim+j],invJEx[i*spaceDim+j]) > 10*PETSC_SMALL) SETERRQ4(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Invalid invJ[%D,%D]: %g != %g", i, j, (double)invJ[i*spaceDim+j], (double)invJEx[i*spaceDim+j]); 153 } 154 } 155 if (RelativeError(detJ,detJEx) > 10*PETSC_SMALL) SETERRQ3(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Invalid |J| = %g != %g diff %g", (double)detJ, (double)detJEx,(double)(detJ - detJEx)); 156 PetscFunctionReturn(0); 157 } 158 159 static PetscErrorCode CheckFVMGeometry(DM dm, PetscInt cell, PetscInt spaceDim, PetscReal centroidEx[], PetscReal normalEx[], PetscReal volEx) 160 { 161 PetscReal centroid[3], normal[3], vol; 162 PetscInt d; 163 PetscErrorCode ierr; 164 165 PetscFunctionBegin; 166 ierr = DMPlexComputeCellGeometryFVM(dm, cell, &vol, centroid, normal);CHKERRQ(ierr); 167 for (d = 0; d < spaceDim; ++d) { 168 if (RelativeError(centroid[d],centroidEx[d]) > 10*PETSC_SMALL) SETERRQ5(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Cell %D, Invalid centroid[%D]: %g != %g diff %g", cell, d, (double)centroid[d], (double)centroidEx[d],(double)(centroid[d]-centroidEx[d])); 169 if (RelativeError(normal[d],normalEx[d]) > 10*PETSC_SMALL) SETERRQ4(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Cell %D, Invalid normal[%D]: %g != %g", cell, d, (double)normal[d], (double)normalEx[d]); 170 } 171 if (RelativeError(volEx,vol) > 10*PETSC_SMALL) SETERRQ4(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Cell %D, Invalid volume = %g != %g diff %g", cell, (double)vol, (double)volEx,(double)(vol - volEx)); 172 PetscFunctionReturn(0); 173 } 174 175 static PetscErrorCode CheckCell(DM dm, PetscInt cell, PetscBool transform, PetscReal v0Ex[], PetscReal JEx[], PetscReal invJEx[], PetscReal detJEx, PetscReal centroidEx[], PetscReal normalEx[], PetscReal volEx, PetscReal faceCentroidEx[], PetscReal faceNormalEx[], PetscReal faceVolEx[]) 176 { 177 const PetscInt *cone; 178 PetscInt coneSize, c; 179 PetscInt dim, depth, cdim; 180 PetscErrorCode ierr; 181 182 PetscFunctionBegin; 183 ierr = DMPlexGetDepth(dm, &depth);CHKERRQ(ierr); 184 ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr); 185 ierr = DMGetCoordinateDim(dm, &cdim);CHKERRQ(ierr); 186 if (v0Ex) { 187 ierr = CheckFEMGeometry(dm, cell, cdim, v0Ex, JEx, invJEx, detJEx);CHKERRQ(ierr); 188 } 189 if (dim == depth && centroidEx) { 190 ierr = CheckFVMGeometry(dm, cell, cdim, centroidEx, normalEx, volEx);CHKERRQ(ierr); 191 if (faceCentroidEx) { 192 ierr = DMPlexGetConeSize(dm, cell, &coneSize);CHKERRQ(ierr); 193 ierr = DMPlexGetCone(dm, cell, &cone);CHKERRQ(ierr); 194 for (c = 0; c < coneSize; ++c) { 195 ierr = CheckFVMGeometry(dm, cone[c], dim, &faceCentroidEx[c*dim], &faceNormalEx[c*dim], faceVolEx[c]);CHKERRQ(ierr); 196 } 197 } 198 } 199 if (transform) { 200 Vec coordinates; 201 PetscSection coordSection; 202 PetscScalar *coords = NULL, *origCoords, *newCoords; 203 PetscReal *v0ExT, *JExT, *invJExT, detJExT=0, *centroidExT, *normalExT, volExT=0; 204 PetscReal *faceCentroidExT, *faceNormalExT, faceVolExT; 205 PetscRandom r, ang, ang2; 206 PetscInt coordSize, numCorners, t; 207 208 ierr = DMGetCoordinatesLocal(dm, &coordinates);CHKERRQ(ierr); 209 ierr = DMGetCoordinateSection(dm, &coordSection);CHKERRQ(ierr); 210 ierr = DMPlexVecGetClosure(dm, coordSection, coordinates, cell, &coordSize, &coords);CHKERRQ(ierr); 211 ierr = PetscMalloc2(coordSize, &origCoords, coordSize, &newCoords);CHKERRQ(ierr); 212 ierr = PetscMalloc5(cdim, &v0ExT, cdim*cdim, &JExT, cdim*cdim, &invJExT, cdim, ¢roidExT, cdim, &normalExT);CHKERRQ(ierr); 213 ierr = PetscMalloc2(cdim, &faceCentroidExT, cdim, &faceNormalExT);CHKERRQ(ierr); 214 for (c = 0; c < coordSize; ++c) origCoords[c] = coords[c]; 215 ierr = DMPlexVecRestoreClosure(dm, coordSection, coordinates, cell, &coordSize, &coords);CHKERRQ(ierr); 216 numCorners = coordSize/cdim; 217 218 ierr = PetscRandomCreate(PETSC_COMM_SELF, &r);CHKERRQ(ierr); 219 ierr = PetscRandomSetFromOptions(r);CHKERRQ(ierr); 220 ierr = PetscRandomSetInterval(r, 0.0, 10.0);CHKERRQ(ierr); 221 ierr = PetscRandomCreate(PETSC_COMM_SELF, &ang);CHKERRQ(ierr); 222 ierr = PetscRandomSetFromOptions(ang);CHKERRQ(ierr); 223 ierr = PetscRandomSetInterval(ang, 0.0, 2*PETSC_PI);CHKERRQ(ierr); 224 ierr = PetscRandomCreate(PETSC_COMM_SELF, &ang2);CHKERRQ(ierr); 225 ierr = PetscRandomSetFromOptions(ang2);CHKERRQ(ierr); 226 ierr = PetscRandomSetInterval(ang2, 0.0, PETSC_PI);CHKERRQ(ierr); 227 for (t = 0; t < 1; ++t) { 228 PetscScalar trans[3]; 229 PetscReal R[9], rot[3], rotM[9]; 230 PetscReal scale, phi, theta, psi = 0.0, norm; 231 PetscInt d, e, f, p; 232 233 for (c = 0; c < coordSize; ++c) newCoords[c] = origCoords[c]; 234 ierr = PetscRandomGetValueReal(r, &scale);CHKERRQ(ierr); 235 ierr = PetscRandomGetValueReal(ang, &phi);CHKERRQ(ierr); 236 ierr = PetscRandomGetValueReal(ang2, &theta);CHKERRQ(ierr); 237 for (d = 0; d < cdim; ++d) { 238 ierr = PetscRandomGetValue(r, &trans[d]);CHKERRQ(ierr); 239 } 240 switch (cdim) { 241 case 2: 242 R[0] = PetscCosReal(phi); R[1] = -PetscSinReal(phi); 243 R[2] = PetscSinReal(phi); R[3] = PetscCosReal(phi); 244 break; 245 case 3: 246 { 247 const PetscReal ct = PetscCosReal(theta), st = PetscSinReal(theta); 248 const PetscReal cp = PetscCosReal(phi), sp = PetscSinReal(phi); 249 const PetscReal cs = PetscCosReal(psi), ss = PetscSinReal(psi); 250 R[0] = ct*cs; R[1] = sp*st*cs - cp*ss; R[2] = sp*ss + cp*st*cs; 251 R[3] = ct*ss; R[4] = cp*cs + sp*st*ss; R[5] = cp*st*ss - sp*cs; 252 R[6] = -st; R[7] = sp*ct; R[8] = cp*ct; 253 break; 254 } 255 default: SETERRQ1(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "Invalid coordinate dimension %D", cdim); 256 } 257 if (v0Ex) { 258 detJExT = detJEx; 259 for (d = 0; d < cdim; ++d) { 260 v0ExT[d] = v0Ex[d]; 261 for (e = 0; e < cdim; ++e) { 262 JExT[d*cdim+e] = JEx[d*cdim+e]; 263 invJExT[d*cdim+e] = invJEx[d*cdim+e]; 264 } 265 } 266 for (d = 0; d < cdim; ++d) { 267 v0ExT[d] *= scale; 268 v0ExT[d] += PetscRealPart(trans[d]); 269 /* Only scale dimensions in the manifold */ 270 for (e = 0; e < dim; ++e) { 271 JExT[d*cdim+e] *= scale; 272 invJExT[d*cdim+e] /= scale; 273 } 274 if (d < dim) detJExT *= scale; 275 } 276 /* Do scaling and translation before rotation, so that we can leave out the normal dimension for lower dimensional manifolds */ 277 for (d = 0; d < cdim; ++d) { 278 for (e = 0, rot[d] = 0.0; e < cdim; ++e) { 279 rot[d] += R[d*cdim+e] * v0ExT[e]; 280 } 281 } 282 for (d = 0; d < cdim; ++d) v0ExT[d] = rot[d]; 283 for (d = 0; d < cdim; ++d) { 284 for (e = 0; e < cdim; ++e) { 285 for (f = 0, rotM[d*cdim+e] = 0.0; f < cdim; ++f) { 286 rotM[d*cdim+e] += R[d*cdim+f] * JExT[f*cdim+e]; 287 } 288 } 289 } 290 for (d = 0; d < cdim; ++d) { 291 for (e = 0; e < cdim; ++e) { 292 JExT[d*cdim+e] = rotM[d*cdim+e]; 293 } 294 } 295 for (d = 0; d < cdim; ++d) { 296 for (e = 0; e < cdim; ++e) { 297 for (f = 0, rotM[d*cdim+e] = 0.0; f < cdim; ++f) { 298 rotM[d*cdim+e] += invJExT[d*cdim+f] * R[e*cdim+f]; 299 } 300 } 301 } 302 for (d = 0; d < cdim; ++d) { 303 for (e = 0; e < cdim; ++e) { 304 invJExT[d*cdim+e] = rotM[d*cdim+e]; 305 } 306 } 307 } 308 if (centroidEx) { 309 volExT = volEx; 310 for (d = 0; d < cdim; ++d) { 311 centroidExT[d] = centroidEx[d]; 312 normalExT[d] = normalEx[d]; 313 } 314 for (d = 0; d < cdim; ++d) { 315 centroidExT[d] *= scale; 316 centroidExT[d] += PetscRealPart(trans[d]); 317 normalExT[d] /= scale; 318 /* Only scale dimensions in the manifold */ 319 if (d < dim) volExT *= scale; 320 } 321 /* Do scaling and translation before rotation, so that we can leave out the normal dimension for lower dimensional manifolds */ 322 for (d = 0; d < cdim; ++d) { 323 for (e = 0, rot[d] = 0.0; e < cdim; ++e) { 324 rot[d] += R[d*cdim+e] * centroidExT[e]; 325 } 326 } 327 for (d = 0; d < cdim; ++d) centroidExT[d] = rot[d]; 328 for (d = 0; d < cdim; ++d) { 329 for (e = 0, rot[d] = 0.0; e < cdim; ++e) { 330 rot[d] += R[d*cdim+e] * normalExT[e]; 331 } 332 } 333 for (d = 0; d < cdim; ++d) normalExT[d] = rot[d]; 334 for (d = 0, norm = 0.0; d < cdim; ++d) norm += PetscSqr(normalExT[d]); 335 norm = PetscSqrtReal(norm); 336 for (d = 0; d < cdim; ++d) normalExT[d] /= norm; 337 } 338 for (d = 0; d < cdim; ++d) { 339 for (p = 0; p < numCorners; ++p) { 340 newCoords[p*cdim+d] *= scale; 341 newCoords[p*cdim+d] += trans[d]; 342 } 343 } 344 for (p = 0; p < numCorners; ++p) { 345 for (d = 0; d < cdim; ++d) { 346 for (e = 0, rot[d] = 0.0; e < cdim; ++e) { 347 rot[d] += R[d*cdim+e] * PetscRealPart(newCoords[p*cdim+e]); 348 } 349 } 350 for (d = 0; d < cdim; ++d) newCoords[p*cdim+d] = rot[d]; 351 } 352 353 ierr = ChangeCoordinates(dm, cdim, newCoords);CHKERRQ(ierr); 354 if (v0Ex) { 355 ierr = CheckFEMGeometry(dm, 0, cdim, v0ExT, JExT, invJExT, detJExT);CHKERRQ(ierr); 356 } 357 if (dim == depth && centroidEx) { 358 ierr = CheckFVMGeometry(dm, cell, cdim, centroidExT, normalExT, volExT);CHKERRQ(ierr); 359 if (faceCentroidEx) { 360 ierr = DMPlexGetConeSize(dm, cell, &coneSize);CHKERRQ(ierr); 361 ierr = DMPlexGetCone(dm, cell, &cone);CHKERRQ(ierr); 362 for (c = 0; c < coneSize; ++c) { 363 PetscInt off = c*cdim; 364 365 faceVolExT = faceVolEx[c]; 366 for (d = 0; d < cdim; ++d) { 367 faceCentroidExT[d] = faceCentroidEx[off+d]; 368 faceNormalExT[d] = faceNormalEx[off+d]; 369 } 370 for (d = 0; d < cdim; ++d) { 371 faceCentroidExT[d] *= scale; 372 faceCentroidExT[d] += PetscRealPart(trans[d]); 373 faceNormalExT[d] /= scale; 374 /* Only scale dimensions in the manifold */ 375 if (d < dim-1) { 376 faceVolExT *= scale; 377 } 378 } 379 for (d = 0; d < cdim; ++d) { 380 for (e = 0, rot[d] = 0.0; e < cdim; ++e) { 381 rot[d] += R[d*cdim+e] * faceCentroidExT[e]; 382 } 383 } 384 for (d = 0; d < cdim; ++d) faceCentroidExT[d] = rot[d]; 385 for (d = 0; d < cdim; ++d) { 386 for (e = 0, rot[d] = 0.0; e < cdim; ++e) { 387 rot[d] += R[d*cdim+e] * faceNormalExT[e]; 388 } 389 } 390 for (d = 0; d < cdim; ++d) faceNormalExT[d] = rot[d]; 391 for (d = 0, norm = 0.0; d < cdim; ++d) norm += PetscSqr(faceNormalExT[d]); 392 norm = PetscSqrtReal(norm); 393 for (d = 0; d < cdim; ++d) faceNormalExT[d] /= norm; 394 395 ierr = CheckFVMGeometry(dm, cone[c], cdim, faceCentroidExT, faceNormalExT, faceVolExT);CHKERRQ(ierr); 396 } 397 } 398 } 399 } 400 ierr = PetscRandomDestroy(&r);CHKERRQ(ierr); 401 ierr = PetscRandomDestroy(&ang);CHKERRQ(ierr); 402 ierr = PetscRandomDestroy(&ang2);CHKERRQ(ierr); 403 ierr = PetscFree2(origCoords, newCoords);CHKERRQ(ierr); 404 ierr = PetscFree5(v0ExT, JExT, invJExT, centroidExT, normalExT);CHKERRQ(ierr); 405 ierr = PetscFree2(faceCentroidExT, faceNormalExT);CHKERRQ(ierr); 406 } 407 PetscFunctionReturn(0); 408 } 409 410 static PetscErrorCode TestTriangle(MPI_Comm comm, PetscBool interpolate, PetscBool transform) 411 { 412 DM dm; 413 PetscInt dim; 414 PetscErrorCode ierr; 415 416 PetscFunctionBegin; 417 /* Create reference triangle */ 418 dim = 2; 419 ierr = DMCreate(comm, &dm);CHKERRQ(ierr); 420 ierr = PetscObjectSetName((PetscObject) dm, "triangle");CHKERRQ(ierr); 421 ierr = DMSetType(dm, DMPLEX);CHKERRQ(ierr); 422 ierr = DMSetDimension(dm, dim);CHKERRQ(ierr); 423 { 424 PetscInt numPoints[2] = {3, 1}; 425 PetscInt coneSize[4] = {3, 0, 0, 0}; 426 PetscInt cones[3] = {1, 2, 3}; 427 PetscInt coneOrientations[3] = {0, 0, 0}; 428 PetscScalar vertexCoords[6] = {-1.0, -1.0, 1.0, -1.0, -1.0, 1.0}; 429 430 ierr = DMPlexCreateFromDAG(dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 431 if (interpolate) { 432 DM idm; 433 434 ierr = DMPlexInterpolate(dm, &idm);CHKERRQ(ierr); 435 ierr = DMDestroy(&dm);CHKERRQ(ierr); 436 dm = idm; 437 } 438 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 439 } 440 /* Check reference geometry: determinant is scaled by reference volume (2.0) */ 441 { 442 PetscReal v0Ex[2] = {-1.0, -1.0}; 443 PetscReal JEx[4] = {1.0, 0.0, 0.0, 1.0}; 444 PetscReal invJEx[4] = {1.0, 0.0, 0.0, 1.0}; 445 PetscReal detJEx = 1.0; 446 PetscReal centroidEx[2] = {-((PetscReal)1.)/((PetscReal)3.), -((PetscReal)1.)/((PetscReal)3.)}; 447 PetscReal normalEx[2] = {0.0, 0.0}; 448 PetscReal volEx = 2.0; 449 450 ierr = CheckCell(dm, 0, transform, v0Ex, JEx, invJEx, detJEx, centroidEx, normalEx, volEx, NULL, NULL, NULL);CHKERRQ(ierr); 451 } 452 /* Move to 3D: Check reference geometry: determinant is scaled by reference volume (2.0) */ 453 dim = 3; 454 { 455 PetscScalar vertexCoords[9] = {-1.0, -1.0, 0.0, 1.0, -1.0, 0.0, -1.0, 1.0, 0.0}; 456 PetscReal v0Ex[3] = {-1.0, -1.0, 0.0}; 457 PetscReal JEx[9] = {1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0}; 458 PetscReal invJEx[9] = {1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0}; 459 PetscReal detJEx = 1.0; 460 PetscReal centroidEx[3] = {-((PetscReal)1.)/((PetscReal)3.), -((PetscReal)1.)/((PetscReal)3.), 0.0}; 461 PetscReal normalEx[3] = {0.0, 0.0, 1.0}; 462 PetscReal volEx = 2.0; 463 464 ierr = ChangeCoordinates(dm, dim, vertexCoords);CHKERRQ(ierr); 465 ierr = CheckCell(dm, 0, transform, v0Ex, JEx, invJEx, detJEx, centroidEx, normalEx, volEx, NULL, NULL, NULL);CHKERRQ(ierr); 466 } 467 /* Cleanup */ 468 ierr = DMDestroy(&dm);CHKERRQ(ierr); 469 PetscFunctionReturn(0); 470 } 471 472 static PetscErrorCode TestQuadrilateral(MPI_Comm comm, PetscBool interpolate, PetscBool transform) 473 { 474 DM dm; 475 PetscInt dim; 476 PetscErrorCode ierr; 477 478 PetscFunctionBegin; 479 /* Create reference quadrilateral */ 480 dim = 2; 481 ierr = DMCreate(comm, &dm);CHKERRQ(ierr); 482 ierr = PetscObjectSetName((PetscObject) dm, "quadrilateral");CHKERRQ(ierr); 483 ierr = DMSetType(dm, DMPLEX);CHKERRQ(ierr); 484 ierr = DMSetDimension(dm, dim);CHKERRQ(ierr); 485 { 486 PetscInt numPoints[2] = {4, 1}; 487 PetscInt coneSize[5] = {4, 0, 0, 0, 0}; 488 PetscInt cones[4] = {1, 2, 3, 4}; 489 PetscInt coneOrientations[4] = {0, 0, 0, 0}; 490 PetscScalar vertexCoords[8] = {-1.0, -1.0, 1.0, -1.0, 1.0, 1.0, -1.0, 1.0}; 491 492 ierr = DMPlexCreateFromDAG(dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 493 if (interpolate) { 494 DM idm; 495 496 ierr = DMPlexInterpolate(dm, &idm);CHKERRQ(ierr); 497 ierr = DMDestroy(&dm);CHKERRQ(ierr); 498 dm = idm; 499 } 500 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 501 } 502 /* Check reference geometry: determinant is scaled by reference volume (2.0) */ 503 { 504 PetscReal v0Ex[2] = {-1.0, -1.0}; 505 PetscReal JEx[4] = {1.0, 0.0, 0.0, 1.0}; 506 PetscReal invJEx[4] = {1.0, 0.0, 0.0, 1.0}; 507 PetscReal detJEx = 1.0; 508 PetscReal centroidEx[2] = {0.0, 0.0}; 509 PetscReal normalEx[2] = {0.0, 0.0}; 510 PetscReal volEx = 4.0; 511 512 ierr = CheckCell(dm, 0, transform, v0Ex, JEx, invJEx, detJEx, centroidEx, normalEx, volEx, NULL, NULL, NULL);CHKERRQ(ierr); 513 } 514 /* Move to 3D: Check reference geometry: determinant is scaled by reference volume (4.0) */ 515 dim = 3; 516 { 517 PetscScalar vertexCoords[12] = {-1.0, -1.0, 0.0, 1.0, -1.0, 0.0, 1.0, 1.0, 0.0, -1.0, 1.0, 0.0}; 518 PetscReal v0Ex[3] = {-1.0, -1.0, 0.0}; 519 PetscReal JEx[9] = {1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0}; 520 PetscReal invJEx[9] = {1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0}; 521 PetscReal detJEx = 1.0; 522 PetscReal centroidEx[3] = {0.0, 0.0, 0.0}; 523 PetscReal normalEx[3] = {0.0, 0.0, 1.0}; 524 PetscReal volEx = 4.0; 525 526 ierr = ChangeCoordinates(dm, dim, vertexCoords);CHKERRQ(ierr); 527 ierr = CheckCell(dm, 0, transform, v0Ex, JEx, invJEx, detJEx, centroidEx, normalEx, volEx, NULL, NULL, NULL);CHKERRQ(ierr); 528 } 529 /* Cleanup */ 530 ierr = DMDestroy(&dm);CHKERRQ(ierr); 531 PetscFunctionReturn(0); 532 } 533 534 static PetscErrorCode TestTetrahedron(MPI_Comm comm, PetscBool interpolate, PetscBool transform) 535 { 536 DM dm; 537 PetscInt dim; 538 PetscErrorCode ierr; 539 540 PetscFunctionBegin; 541 /* Create reference tetrahedron */ 542 dim = 3; 543 ierr = DMCreate(comm, &dm);CHKERRQ(ierr); 544 ierr = PetscObjectSetName((PetscObject) dm, "tetrahedron");CHKERRQ(ierr); 545 ierr = DMSetType(dm, DMPLEX);CHKERRQ(ierr); 546 ierr = DMSetDimension(dm, dim);CHKERRQ(ierr); 547 { 548 PetscInt numPoints[2] = {4, 1}; 549 PetscInt coneSize[5] = {4, 0, 0, 0, 0}; 550 PetscInt cones[4] = {1, 2, 3, 4}; 551 PetscInt coneOrientations[4] = {0, 0, 0, 0}; 552 PetscScalar vertexCoords[12] = {-1.0, -1.0, -1.0, -1.0, 1.0, -1.0, 1.0, -1.0, -1.0, -1.0, -1.0, 1.0}; 553 554 ierr = DMPlexCreateFromDAG(dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 555 if (interpolate) { 556 DM idm; 557 558 ierr = DMPlexInterpolate(dm, &idm);CHKERRQ(ierr); 559 ierr = DMDestroy(&dm);CHKERRQ(ierr); 560 dm = idm; 561 } 562 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 563 } 564 /* Check reference geometry: determinant is scaled by reference volume (4/3) */ 565 { 566 PetscReal v0Ex[3] = {-1.0, -1.0, -1.0}; 567 PetscReal JEx[9] = {1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0}; 568 PetscReal invJEx[9] = {1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0}; 569 PetscReal detJEx = 1.0; 570 PetscReal centroidEx[3] = {-0.5, -0.5, -0.5}; 571 PetscReal normalEx[3] = {0.0, 0.0, 0.0}; 572 PetscReal volEx = (PetscReal)4.0/(PetscReal)3.0; 573 574 ierr = CheckCell(dm, 0, transform, v0Ex, JEx, invJEx, detJEx, centroidEx, normalEx, volEx, NULL, NULL, NULL);CHKERRQ(ierr); 575 } 576 /* Cleanup */ 577 ierr = DMDestroy(&dm);CHKERRQ(ierr); 578 PetscFunctionReturn(0); 579 } 580 581 static PetscErrorCode TestHexahedron(MPI_Comm comm, PetscBool interpolate, PetscBool transform) 582 { 583 DM dm; 584 PetscInt dim; 585 PetscErrorCode ierr; 586 587 PetscFunctionBegin; 588 /* Create reference hexahedron */ 589 dim = 3; 590 ierr = DMCreate(comm, &dm);CHKERRQ(ierr); 591 ierr = PetscObjectSetName((PetscObject) dm, "hexahedron");CHKERRQ(ierr); 592 ierr = DMSetType(dm, DMPLEX);CHKERRQ(ierr); 593 ierr = DMSetDimension(dm, dim);CHKERRQ(ierr); 594 { 595 PetscInt numPoints[2] = {8, 1}; 596 PetscInt coneSize[9] = {8, 0, 0, 0, 0, 0, 0, 0, 0}; 597 PetscInt cones[8] = {1, 2, 3, 4, 5, 6, 7, 8}; 598 PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0}; 599 PetscScalar vertexCoords[24] = {-1.0, -1.0, -1.0, -1.0, 1.0, -1.0, 1.0, 1.0, -1.0, 1.0, -1.0, -1.0, 600 -1.0, -1.0, 1.0, 1.0, -1.0, 1.0, 1.0, 1.0, 1.0, -1.0, 1.0, 1.0}; 601 602 ierr = DMPlexCreateFromDAG(dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 603 if (interpolate) { 604 DM idm; 605 606 ierr = DMPlexInterpolate(dm, &idm);CHKERRQ(ierr); 607 ierr = DMDestroy(&dm);CHKERRQ(ierr); 608 dm = idm; 609 } 610 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 611 } 612 /* Check reference geometry: determinant is scaled by reference volume 8.0 */ 613 { 614 PetscReal v0Ex[3] = {-1.0, -1.0, -1.0}; 615 PetscReal JEx[9] = {1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0}; 616 PetscReal invJEx[9] = {1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0}; 617 PetscReal detJEx = 1.0; 618 PetscReal centroidEx[3] = {0.0, 0.0, 0.0}; 619 PetscReal normalEx[3] = {0.0, 0.0, 0.0}; 620 PetscReal volEx = 8.0; 621 622 ierr = CheckCell(dm, 0, transform, v0Ex, JEx, invJEx, detJEx, centroidEx, normalEx, volEx, NULL, NULL, NULL);CHKERRQ(ierr); 623 } 624 /* Cleanup */ 625 ierr = DMDestroy(&dm);CHKERRQ(ierr); 626 PetscFunctionReturn(0); 627 } 628 629 /* This wedge is a tensor product cell, rather than a normal wedge */ 630 static PetscErrorCode TestWedge(MPI_Comm comm, PetscBool interpolate, PetscBool transform) 631 { 632 DM dm; 633 PetscInt dim; 634 PetscErrorCode ierr; 635 636 PetscFunctionBegin; 637 /* Create reference wedge */ 638 dim = 3; 639 ierr = DMCreate(comm, &dm);CHKERRQ(ierr); 640 ierr = PetscObjectSetName((PetscObject) dm, "Triangular Prism");CHKERRQ(ierr); 641 ierr = DMSetType(dm, DMPLEX);CHKERRQ(ierr); 642 ierr = DMSetDimension(dm, dim);CHKERRQ(ierr); 643 { 644 PetscInt numPoints[2] = {6, 1}; 645 PetscInt coneSize[7] = {6, 0, 0, 0, 0, 0, 0}; 646 PetscInt cones[6] = {1, 2, 3, 4, 5, 6}; 647 PetscInt coneOrientations[6] = {0, 0, 0, 0, 0, 0}; 648 PetscScalar vertexCoords[18] = {-1.0, -1.0, -1.0, 1.0, -1.0, -1.0, -1.0, 1.0, -1.0, 649 -1.0, -1.0, 1.0, 1.0, -1.0, 1.0, -1.0, 1.0, 1.0}; 650 651 ierr = DMPlexCreateFromDAG(dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 652 ierr = DMPlexSetHybridBounds(dm, 0, PETSC_DETERMINE, PETSC_DETERMINE, PETSC_DETERMINE);CHKERRQ(ierr); 653 if (interpolate) { 654 DM idm; 655 656 ierr = DMPlexInterpolate(dm, &idm);CHKERRQ(ierr); 657 ierr = DMDestroy(&dm);CHKERRQ(ierr); 658 dm = idm; 659 } 660 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 661 } 662 /* Check reference geometry: determinant is scaled by reference volume 4.0 */ 663 { 664 #if 0 665 /* FEM geometry is not functional for wedges */ 666 PetscReal v0Ex[3] = {-1.0, -1.0, -1.0}; 667 PetscReal JEx[9] = {1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0}; 668 PetscReal invJEx[9] = {1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0}; 669 PetscReal detJEx = 1.0; 670 #endif 671 672 if (interpolate) { 673 PetscReal centroidEx[3] = {-((PetscReal)1.)/((PetscReal)3.), -((PetscReal)1.)/((PetscReal)3.), 0.0}; 674 PetscReal normalEx[3] = {0.0, 0.0, 0.0}; 675 PetscReal volEx = 4.0; 676 PetscReal faceVolEx[5] = {2.0, 2.0, 4.0, PETSC_SQRT2*4.0, 4.0}; 677 PetscReal faceNormalEx[15] = {0.0, 0.0, 1.0, 0.0, 0.0, 1.0, 0.0, -1.0, 0.0, PETSC_SQRT2/2.0, PETSC_SQRT2/2.0, 0.0, -1.0, 0.0, 0.0}; 678 PetscReal faceCentroidEx[15] = {-((PetscReal)1.)/((PetscReal)3.), -((PetscReal)1.)/((PetscReal)3.), -1.0, 679 -((PetscReal)1.)/((PetscReal)3.), -((PetscReal)1.)/((PetscReal)3.), 1.0, 680 0.0, -1.0, 0.0, 0.0, 0.0, 0.0, -1.0, 0.0, 0.0}; 681 682 ierr = CheckCell(dm, 0, transform, NULL, NULL, NULL, 0.0, centroidEx, normalEx, volEx, faceCentroidEx, faceNormalEx, faceVolEx);CHKERRQ(ierr); 683 } 684 } 685 /* Cleanup */ 686 ierr = DMDestroy(&dm);CHKERRQ(ierr); 687 PetscFunctionReturn(0); 688 } 689 690 int main(int argc, char **argv) 691 { 692 AppCtx user; 693 PetscErrorCode ierr; 694 695 ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr; 696 ierr = ProcessOptions(PETSC_COMM_WORLD, &user);CHKERRQ(ierr); 697 if (user.runType == RUN_REFERENCE) { 698 ierr = TestTriangle(PETSC_COMM_SELF, user.interpolate, user.transform);CHKERRQ(ierr); 699 ierr = TestQuadrilateral(PETSC_COMM_SELF, user.interpolate, user.transform);CHKERRQ(ierr); 700 ierr = TestTetrahedron(PETSC_COMM_SELF, user.interpolate, user.transform);CHKERRQ(ierr); 701 ierr = TestHexahedron(PETSC_COMM_SELF, user.interpolate, user.transform);CHKERRQ(ierr); 702 ierr = TestWedge(PETSC_COMM_SELF, user.interpolate, user.transform);CHKERRQ(ierr); 703 } else if (user.runType == RUN_FILE) { 704 PetscInt dim, cStart, cEnd, c; 705 706 ierr = DMGetDimension(user.dm, &dim);CHKERRQ(ierr); 707 ierr = DMPlexGetHeightStratum(user.dm, 0, &cStart, &cEnd);CHKERRQ(ierr); 708 for (c = 0; c < cEnd-cStart; ++c) { 709 ierr = CheckCell(user.dm, c+cStart, PETSC_FALSE, &user.v0[c*dim], &user.J[c*dim*dim], &user.invJ[c*dim*dim], user.detJ[c], &user.centroid[c*dim], &user.normal[c*dim], user.vol[c], NULL, NULL, NULL);CHKERRQ(ierr); 710 } 711 ierr = PetscFree7(user.v0,user.J,user.invJ,user.detJ,user.centroid,user.normal,user.vol);CHKERRQ(ierr); 712 ierr = DMDestroy(&user.dm);CHKERRQ(ierr); 713 } else if (user.runType == RUN_DISPLAY) { 714 DM gdm, dmCell; 715 Vec cellgeom, facegeom; 716 const PetscScalar *cgeom; 717 PetscInt dim, d, cStart, cEnd, cEndInterior, c; 718 719 ierr = DMGetCoordinateDim(user.dm, &dim);CHKERRQ(ierr); 720 ierr = DMPlexConstructGhostCells(user.dm, NULL, NULL, &gdm);CHKERRQ(ierr); 721 if (gdm) { 722 ierr = DMDestroy(&user.dm);CHKERRQ(ierr); 723 user.dm = gdm; 724 } 725 ierr = DMPlexComputeGeometryFVM(user.dm, &cellgeom, &facegeom);CHKERRQ(ierr); 726 ierr = DMPlexGetHeightStratum(user.dm, 0, &cStart, &cEnd);CHKERRQ(ierr); 727 ierr = DMPlexGetGhostCellStratum(user.dm, &cEndInterior, NULL);CHKERRQ(ierr); 728 if (cEndInterior >= 0) cEnd = cEndInterior; 729 ierr = VecGetDM(cellgeom, &dmCell);CHKERRQ(ierr); 730 ierr = VecGetArrayRead(cellgeom, &cgeom);CHKERRQ(ierr); 731 for (c = 0; c < cEnd-cStart; ++c) { 732 PetscFVCellGeom *cg; 733 734 ierr = DMPlexPointLocalRead(dmCell, c, cgeom, &cg);CHKERRQ(ierr); 735 ierr = PetscPrintf(PETSC_COMM_SELF, "Cell %4D: Centroid (", c);CHKERRQ(ierr); 736 for (d = 0; d < dim; ++d) { 737 if (d > 0) {ierr = PetscPrintf(PETSC_COMM_SELF, ", ");CHKERRQ(ierr);} 738 ierr = PetscPrintf(PETSC_COMM_SELF, "%12.2g", cg->centroid[d]);CHKERRQ(ierr); 739 } 740 ierr = PetscPrintf(PETSC_COMM_SELF, ") Vol %12.2g\n", cg->volume);CHKERRQ(ierr); 741 } 742 ierr = VecRestoreArrayRead(cellgeom, &cgeom);CHKERRQ(ierr); 743 ierr = VecDestroy(&cellgeom);CHKERRQ(ierr); 744 ierr = VecDestroy(&facegeom);CHKERRQ(ierr); 745 ierr = DMDestroy(&user.dm);CHKERRQ(ierr); 746 } 747 ierr = PetscFinalize(); 748 return ierr; 749 } 750 751 /*TEST 752 753 test: 754 suffix: 0 755 args: -dm_view ascii::ascii_info_detail 756 test: 757 suffix: 1 758 args: -interpolate -dm_view ascii::ascii_info_detail 759 test: 760 suffix: 2 761 args: -transform 762 test: 763 suffix: 3 764 args: -interpolate -transform 765 test: 766 suffix: 4 767 requires: exodusii 768 args: -run_type file -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/simpleblock-100.exo -dm_view ascii::ascii_info_detail -v0 -1.5,-0.5,0.5,-0.5,-0.5,0.5,0.5,-0.5,0.5 -J 0.0,0.0,0.5,0.0,0.5,0.0,-0.5,0.0,0.0,0.0,0.0,0.5,0.0,0.5,0.0,-0.5,0.0,0.0,0.0,0.0,0.5,0.0,0.5,0.0,-0.5,0.0,0.0 -invJ 0.0,0.0,-2.0,0.0,2.0,0.0,2.0,0.0,0.0,0.0,0.0,-2.0,0.0,2.0,0.0,2.0,0.0,0.0,0.0,0.0,-2.0,0.0,2.0,0.0,2.0,0.0,0.0 -detJ 0.125,0.125,0.125 769 test: 770 suffix: 5 771 requires: exodusii 772 args: -interpolate -run_type file -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/simpleblock-100.exo -dm_view ascii::ascii_info_detail -v0 -1.5,-0.5,0.5,-0.5,-0.5,0.5,0.5,-0.5,0.5 -J 0.0,0.0,0.5,0.0,0.5,0.0,-0.5,0.0,0.0,0.0,0.0,0.5,0.0,0.5,0.0,-0.5,0.0,0.0,0.0,0.0,0.5,0.0,0.5,0.0,-0.5,0.0,0.0 -invJ 0.0,0.0,-2.0,0.0,2.0,0.0,2.0,0.0,0.0,0.0,0.0,-2.0,0.0,2.0,0.0,2.0,0.0,0.0,0.0,0.0,-2.0,0.0,2.0,0.0,2.0,0.0,0.0 -detJ 0.125,0.125,0.125 -centroid -1.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0 -normal 0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 -vol 1.0,1.0,1.0 773 774 TEST*/ 775