Actual source code: mesh.c

  1: #ifdef PETSC_RCS_HEADER
  2: static char vcid[] = "$Id: mesh.c,v 1.19 2000/10/17 13:48:57 knepley Exp $";
  3: #endif

  5: /*
  6:      Defines the interface to the mesh functions
  7: */

 9:  #include src/mesh/meshimpl.h

 11: /* Logging support */
 12: int MESH_COOKIE;
 13: int MESH_Reform, MESH_IsDistorted, MESH_Partition, MESH_SetupBoundary, MESH_MoveMesh, MESH_CalcNodeVelocities;
 14: int MESH_CalcNodeAccelerations, MESH_CreateLocalCSR, MESH_CreateFullCSR, MESH_CreateDualCSR, MESH_LocatePoint;

 16: /*------------------------------------------------ Generic Operations ------------------------------------------------*/
 17: /*@
 18:   MeshSetUp - Set up any required internal data structures for a Mesh.

 20:   Collective on Mesh

 22:   Input Parameter:
 23: . mesh - The mesh

 25:   Level: beginner

 27: .keywords: Mesh, setup
 28: .seealso: MeshDestroy()
 29: @*/
 30: int MeshSetUp(Mesh mesh)
 31: {

 36:   if (mesh->setupcalled) return(0);
 37:   if (mesh->ops->setup) {
 38:     (*mesh->ops->setup)(mesh);
 39:   }
 40:   return(0);
 41: }

 43: /*
 44:   MeshSetTypeFromOptions - Sets the type of mesh generation from user options.

 46:   Collective on Mesh

 48:   Input Parameter:
 49: . mesh - The mesh

 51:   Level: intermediate

 53: .keywords: Mesh, set, options, database, type
 54: .seealso: MeshSetFromOptions(), MeshSetType()
 55: */
 56: static int MeshSetTypeFromOptions(Mesh mesh)
 57: {
 58:   PetscTruth opt;
 59:   char      *defaultType;
 60:   char       typeName[256];
 61:   int        dim;
 62:   int        ierr;

 65:   MeshGetDimension(mesh, &dim);
 66:   if (mesh->dim == -1) {
 67:     dim  = 2;
 68:     if (mesh->dict != PETSC_NULL) {
 69:       ParameterDictGetInteger(mesh->dict, "dim", &dim);
 70:     }
 71:     MeshSetDimension(mesh, dim);
 72:   }
 73:   if (mesh->type_name != PETSC_NULL) {
 74:     defaultType = mesh->type_name;
 75:   } else {
 76:     switch(dim)
 77:     {
 78:     case 2:
 79:       defaultType = MESH_TRIANGULAR_2D;
 80:       break;
 81:     default:
 82:       SETERRQ1(PETSC_ERR_SUP, "Mesh dimension %d is not supported", dim);
 83:     }
 84:   }

 86:   if (!MeshRegisterAllCalled) {
 87:     MeshRegisterAll(PETSC_NULL);
 88:   }
 89:   PetscOptionsList("-mesh_type", "Mesh generation method"," MeshSetType", MeshList, defaultType, typeName, 256, &opt);
 90: 
 91:   if (opt == PETSC_TRUE) {
 92:     MeshSetType(mesh, typeName);
 93:   } else {
 94:     MeshSetType(mesh, defaultType);
 95:   }
 96:   return(0);
 97: }

 99: /*@
100:   MeshSetFromOptions - Sets various Mesh parameters from user options.

102:   Collective on Mesh

104:   Input Parameter:
105: . mesh - The mesh

107:   Notes:  To see all options, run your program with the -help option, or consult the users manual.
108:           Must be called after MeshCreate() but before the Mesh is used.

110:   Level: intermediate

112: .keywords: Mesh, set, options, database
113: .seealso: MeshCreate(), MeshPrintHelp(), MeshSetOptionsPrefix()
114: @*/
115: int MeshSetFromOptions(Mesh mesh)
116: {
117:   Partition  part;
118:   PetscTruth opt;
119:   int        ierr;

123:   PetscOptionsBegin(mesh->comm, mesh->prefix, "Mesh options", "Mesh");

125:   /* Handle generic mesh options */
126:   PetscOptionsHasName(PETSC_NULL, "-help", &opt);
127:   if (opt == PETSC_TRUE) {
128:     MeshPrintHelp(mesh);
129:   }
130:   PetscOptionsGetInt(mesh->prefix, "-mesh_dim", &mesh->dim, &opt);
131:   PetscOptionsGetReal(mesh->prefix, "-mesh_max_aspect_ratio", &mesh->maxAspectRatio, &opt);

133:   /* Handle mesh type options */
134:   MeshSetTypeFromOptions(mesh);

136:   /* Handle specific mesh options */
137:   if (mesh->ops->setfromoptions != PETSC_NULL) {
138:     (*mesh->ops->setfromoptions)(mesh);
139:   }

141:   PetscOptionsEnd();

143:   /* Handle subobject options */
144:   MeshGetPartition(mesh, &part);
145:   PartitionSetFromOptions(part);

147:   MeshViewFromOptions(mesh, mesh->name);
148:   return(0);
149: }

151: /*@
152:   MeshViewFromOptions - This function visualizes the mesh based upon user options.

154:   Collective on Mesh

156:   Input Parameter:
157: . mesh - The mesh

159:   Level: intermediate

161: .keywords: Mesh, view, options, database
162: .seealso: MeshSetFromOptions(), MeshView()
163: @*/
164: int MeshViewFromOptions(Mesh mesh, char *title)
165: {
166:   PetscViewer viewer;
167:   PetscDraw   draw;
168:   PetscTruth  opt;
169:   char       *titleStr;
170:   char        typeName[1024];
171:   char        fileName[1024];
172:   int         len;
173:   int         ierr;

176:   PetscOptionsHasName(mesh->prefix, "-mesh_view", &opt);
177:   if (opt == PETSC_TRUE) {
178:     PetscOptionsGetString(mesh->prefix, "-mesh_view", typeName, 1024, &opt);
179:     PetscStrlen(typeName, &len);
180:     if (len > 0) {
181:       PetscViewerCreate(mesh->comm, &viewer);
182:       PetscViewerSetType(viewer, typeName);
183:       PetscOptionsGetString(mesh->prefix, "-mesh_view_file", fileName, 1024, &opt);
184:       if (opt == PETSC_TRUE) {
185:         PetscViewerSetFilename(viewer, fileName);
186:       } else {
187:         PetscViewerSetFilename(viewer, mesh->name);
188:       }
189:       MeshView(mesh, viewer);
190:       PetscViewerFlush(viewer);
191:       PetscViewerDestroy(viewer);
192:     } else {
193:       MeshView(mesh, PETSC_NULL);
194:     }
195:   }
196:   PetscOptionsHasName(mesh->prefix, "-mesh_view_draw", &opt);
197:   if (opt == PETSC_TRUE) {
198:     PetscViewerDrawOpen(mesh->comm, 0, 0, 0, 0, 300, 300, &viewer);
199:     PetscViewerDrawGetDraw(viewer, 0, &draw);
200:     if (title != PETSC_NULL) {
201:       titleStr = title;
202:     } else {
203:       PetscObjectName((PetscObject) mesh);                                                         CHKERRQ(ierr) ;
204:       titleStr = mesh->name;
205:     }
206:     PetscDrawSetTitle(draw, titleStr);
207:     MeshView(mesh, viewer);
208:     PetscViewerFlush(viewer);
209:     PetscDrawPause(draw);
210:     PetscViewerDestroy(viewer);
211:   }
212:   return(0);
213: }

215: /*@
216:   MeshView - Views a mesh object.

218:   Collective on Mesh

220:   Input Parameters:
221: + mesh   - The mesh
222: - viewer - The viewer with which to view the mesh

224:   Level: beginner

226: .keywords: mesh, view
227: .seealso: MeshDestroy(), PetscViewerDrawOpen()
228: @*/
229: int MeshView(Mesh mesh, PetscViewer viewer)
230: {

235:   if (viewer == PETSC_NULL) {
236:     viewer = PETSC_VIEWER_STDOUT_SELF;
237:   } else {
239:   }
240:   (*mesh->ops->view)(mesh, viewer);
241:   return(0);
242: }

244: /*
245:   MeshPreCopy_Private - Executes hook for generic pre-copying actions.

247:   Collective on Mesh

249:   Input Parameter:
250: . mesh - The mesh

252:   Level: developer

254: .keywords: Mesh, copy
255: .seealso: MeshPosCopy_Private(), MeshCopy(), MeshDuplicate(), MeshRefine(), MeshReform()
256: */
257: int MeshPreCopy_Private(Mesh mesh) {
258:   int (*precopy)(Mesh);
259:   int   ierr;

263:   /* Specific operations */
264:   PetscObjectQueryFunction((PetscObject) mesh, "PreCopy", (void (**)(void)) &precopy);
265:   if (precopy != PETSC_NULL) {
266:     (*precopy)(mesh);
267:   }
268:   return(0);
269: }

271: /*
272:   MeshPostCopy_Private - Executes hook for generic post-copying actions.

274:   Collective on Mesh

276:   Input Parameters:
277: + mesh    - The mesh
278: - newMesh - The new mesh

280:   Level: developer

282: .keywords: Mesh, copy
283: .seealso: MeshPosCopy_Private(), MeshCopy(), MeshDuplicate(), MeshRefine(), MeshReform()
284: */
285: int MeshPostCopy_Private(Mesh mesh, Mesh newMesh) {
286:   int      (*postcopy)(Mesh, Mesh);
287:   void      *ctx;
288:   PetscTruth isMoving;
289:   int        highlightElement;
290:   int        ierr;

295:   /* Generic operations */
296:   MeshGetUserContext(mesh, &ctx);
297:   MeshSetUserContext(newMesh, ctx);
298:   MeshGetHighlightElement(mesh, &highlightElement);
299:   MeshSetHighlightElement(newMesh, highlightElement);
300:   MeshGetMovement(mesh, &isMoving);
301:   MeshSetMovement(newMesh, isMoving);
302:   /* Specific operations */
303:   PetscObjectQueryFunction((PetscObject) mesh, "PostCopy", (void (**)(void)) &postcopy);
304:   if (postcopy != PETSC_NULL) {
305:     (*postcopy)(mesh, newMesh);
306:   }
307:   return(0);
308: }

310: /*@
311:   MeshCopy - Copies all the data from one mesh to another.

313:   Collective on Mesh

315:   Input Parameter:
316: . mesh    - The mesh

318:   Output Parameter:
319: . newmesh - The identical new mesh

321:   Level: beginner

323: .keywords: Mesh, copy
324: .seealso: MeshDuplicate()
325: @*/
326: int MeshCopy(Mesh mesh, Mesh newmesh)
327: {

334:   if (mesh->ops->copy == PETSC_NULL) SETERRQ(PETSC_ERR_SUP, " ");
335:   MeshPreCopy_Private(mesh);
336:   (*mesh->ops->copy)(mesh, newmesh);
337:   MeshPostCopy_Private(mesh, newmesh);
338:   return(0);
339: }

341: /*@
342:   MeshDuplicate - Duplicates the Mesh.

344:   Collective on Mesh

346:   Input Parameter:
347: . mesh    - The mesh

349:   Output Parameter:
350: . newmesh - The duplicate mesh

352:   Level: beginner

354: .keywords: mesh, duplicate, copy
355: .seealso: MeshCopy()
356: @*/
357: int MeshDuplicate(Mesh mesh, Mesh *newmesh)
358: {
359:   int   ierr;

364:   if (mesh->ops->duplicate == PETSC_NULL) SETERRQ(PETSC_ERR_SUP, " ");
365:   MeshPreCopy_Private(mesh);
366:   (*mesh->ops->duplicate)(mesh, newmesh);
367:   MeshPostCopy_Private(mesh, *newmesh);
368:   return(0);
369: }

371: /*@
372:   MeshDestroy - Destroys a mesh object.

374:   Collective on Mesh

376:   Input Parameter:
377: . mesh - The mesh

379:   Level: beginner

381: .keywords: mesh, destroy
382: .seealso MeshCreate()
383: @*/
384: int MeshDestroy(Mesh mesh)
385: {

390:   if (--mesh->refct > 0) return(0);
391:   (*mesh->ops->destroy)(mesh);
392:   if (mesh->holes != PETSC_NULL) {
393:     PetscFree(mesh->holes);
394:   }
395:   PetscLogObjectDestroy(mesh);
396:   PetscHeaderDestroy(mesh);
397:   return(0);
398: }

400: /*@
401:   MeshSetUserContext - Sets the optional user-defined context for a mesh object.

403:   Collective on Mesh

405:   Input Parameters:
406: + mesh   - The mesh
407: - usrCtx - The optional user context

409:   Level: intermediate

411: .keywords: mesh, set, context
412: .seealso: MeshGetUserContext(), GridSetMeshContext()
413: @*/
414: int MeshSetUserContext(Mesh mesh, void *usrCtx)
415: {
418:   mesh->usr = usrCtx;
419:   return(0);
420: }

422: /*@
423:   MeshGetUserContext - Gets the optional user-defined context for a mesh object.

425:   Not collective

427:   Input Parameters:
428: . mesh - The mesh

430:   Output Parameters:
431: . usrCtx   - The optional user context

433:   Level: intermediate

435: .keywords: mesh, set, context
436: .seealso: MeshSetUserContext(), GridGetMeshContext()
437: @*/
438: int MeshGetUserContext(Mesh mesh, void **usrCtx)
439: {
443:   *usrCtx = mesh->usr;
444:   return(0);
445: }

447: /*@C
448:   MeshSetOptionsPrefix - Sets the prefix used for searching for all Mesh options in the database.

450:   Not collective

452:   Input Parameters:
453: + mesh   - The mesh
454: - prefix - The prefix to prepend to all option names

456:   Notes:
457:   A hyphen (-) must NOT be given at the beginning of the prefix name.
458:   The first character of all runtime options is AUTOMATICALLY the hyphen.

460:   Level: intermediate

462: .keywords: Mesh, set, options, prefix, database
463: .seealso: MeshGetOptionsPrefix(), MeshAppendOptionsPrefix(), MeshSetFromOptions()
464: @*/
465: int MeshSetOptionsPrefix(Mesh mesh, char *prefix)
466: {

471:   PetscObjectSetOptionsPrefix((PetscObject) mesh, prefix);
472:   return(0);
473: }

475: /*@C
476:   MeshAppendOptionsPrefix - Appends to the prefix used for searching for all Mesh options in the database.

478:   Not collective

480:   Input Parameters:
481: + mesh   - The mesh
482: - prefix - The prefix to prepend to all option names

484:   Notes:
485:   A hyphen (-) must NOT be given at the beginning of the prefix name.
486:   The first character of all runtime options is AUTOMATICALLY the hyphen.

488:   Level: intermediate

490: .keywords: Mesh, append, options, prefix, database
491: .seealso: MeshGetOptionsPrefix(), MeshSetOptionsPrefix(), MeshSetFromOptions()
492: @*/
493: int MeshAppendOptionsPrefix(Mesh mesh, char *prefix)
494: {
496: 
499:   PetscObjectAppendOptionsPrefix((PetscObject) mesh, prefix);
500:   return(0);
501: }

503: /*@
504:   MeshGetOptionsPrefix - Returns the prefix used for searching for all Mesh options in the database.

506:   Input Parameter:
507: . mesh   - The mesh

509:   Output Parameter:
510: . prefix - A pointer to the prefix string used

512:   Level: intermediate

514: .keywords: Mesh, get, options, prefix, database
515: .seealso: MeshSetOptionsPrefix(), MeshSetOptionsPrefix(), MeshSetFromOptions()
516: @*/
517: int MeshGetOptionsPrefix(Mesh mesh, char **prefix)
518: {

523:   PetscObjectGetOptionsPrefix((PetscObject) mesh, prefix);
524:   return(0);
525: }

527: /*@
528:   MeshPrintHelp - Prints all options for the Mesh.

530:   Input Parameter:
531: . mesh - The mesh

533:   Options Database Keys:
534: $  -help, -h

536:   Level: intermediate

538: .keywords: Mesh, help
539: .seealso: MeshSetFromOptions()
540: @*/
541: int MeshPrintHelp(Mesh mesh)
542: {
543:   char p[64];
544:   int  ierr;


549:   PetscStrcpy(p, "-");
550:   if (mesh->prefix != PETSC_NULL) {
551:     PetscStrcat(p, mesh->prefix);
552:   }

554:   (*PetscHelpPrintf)(mesh->comm, "Mesh options ------------------------------------------------n");
555:   (*PetscHelpPrintf)(mesh->comm,"   %smesh_type <typename> : Sets the mesh typen", p);
556:   (*PetscHelpPrintf)(mesh->comm,"   %smesh_dim  <num>      : Sets the dimension of the meshn", p);
557:   return(0);
558: }

560: /*--------------------------------------------- Mesh-Specific Operations ---------------------------------------------*/

562: /*
563:   MeshCheckBoundary_Triangular_2D - Checks a boundary for validity.

565:   Collective on Mesh

567:   Input Parameter:
568: . mesh        - The mesh

570:   Input Parameters from bdCtx:
571: + numBD       - The number of closed boundaries in the geometry, or different markers
572: . numVertices - The umber of boundary points
573: . vertices    - The (x,y) coordinates of the boundary points
574: . markers     - The boundary markers for nodes, 0 indicates an interior point, each boundary must have a different marker
575: . numSegments - The number of boundary segments
576: . segments    - The endpoints of boundary segments or PETSC_NULL
577: . segMarkers  - The boundary markers for each segment
578: . numHoles    - The number of holes
579: - holes       - The (x,y) coordinates of holes or PETSC_NULL

581:   Level: advanced

583: .keywords: mesh, boundary
584: .seealso: MeshSetBoundary(), MeshSetReformBoundary()
585: */
586: static int MeshCheckBoundary_Triangular_2D(Mesh mesh, MeshBoundary2D *bdCtx) {
587:   int rank;

592:   MPI_Comm_rank(mesh->comm, &rank);
593:   if (rank == 0) {
594:     if (bdCtx->numVertices < 3) {
595:       SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE, "Insufficient number of boundary points %d", bdCtx->numVertices);
596:     }
597:     if (bdCtx->numSegments < 3) {
598:       SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE, "Invalid number of boundary segments %d", bdCtx->numSegments);
599:     }
600:     if (bdCtx->numHoles    < 0) {
601:       SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE, "Invalid number of holes %d", bdCtx->numHoles);
602:     }
603:     if (bdCtx->numVertices > 0) {
606:     }
607:     if (bdCtx->numSegments > 0) {
610:     }
611:     if (bdCtx->numHoles    > 0) {
613:     }
614:   }
615:   return(0);
616: }

618: /*@
619:   MeshCheckBoundary - Checks the mesh boundary for validity.

621:   Collective on Mesh

623:   Input Parameters:
624: + mesh  - The mesh
625: - bdCtx - The MeshBoundary

627:   Level: intermediate

629: .keywords: mesh, partition
630: .seealso: MeshSetBoundary(), MeshSetReformBoundary()
631: @*/
632: int MeshCheckBoundary(Mesh mesh, MeshBoundary2D *bdCtx) {

637:   /* This needs to go in MeshBoundary later */
638:   MeshCheckBoundary_Triangular_2D(mesh, bdCtx);
639:   return(0);
640: }

642: /*@
643:   MeshSetBoundary - Store a boundary to use for mesh generation.

645:   Not collective

647:   Input Parameters:
648: + mesh  - The mesh
649: - bdCtx - The MeshBoundary

651:   Level: intermediate

653: .keywords: mesh, boundary, refinement
654: .seealso: MeshReform()
655: @*/
656: int MeshSetBoundary(Mesh mesh, MeshBoundary2D *bdCtx) {

662:   MeshCheckBoundary(mesh, bdCtx);
663:   mesh->bdCtx = bdCtx;
664:   return(0);
665: }

667: /*@
668:   MeshSetReformBoundary - Store an alternate boundary to use for reforming.

670:   Not collective

672:   Input Parameters:
673: + mesh  - The mesh
674: - bdCtx - The MeshBoundary

676:   Level: intermediate

678: .keywords: mesh, boundary, refinement
679: .seealso: MeshReform()
680: @*/
681: int MeshSetReformBoundary(Mesh mesh, MeshBoundary2D *bdCtx) {

687:   MeshCheckBoundary(mesh, bdCtx);
688:   mesh->bdCtxNew = bdCtx;
689:   return(0);
690: }

692: /*@
693:   MeshPartition - Partitions the mesh among the processors

695:   Collective on Mesh

697:   Input Parameter:
698: . mesh - the mesh

700:   Level: beginner

702: .keywords: mesh, partition
703: .seealso: MeshGetBoundaryStart()
704: @*/
705: int MeshPartition(Mesh mesh)
706: {

711:   PetscLogEventBegin(MESH_Partition, mesh, 0, 0, 0);
712:   (*mesh->ops->partition)(mesh);
713:   PetscLogEventEnd(MESH_Partition, mesh, 0, 0, 0);
714:   MeshUpdateBoundingBox(mesh);
715:   return(0);
716: }
717: 
718: /*@
719:   MeshCoarsen - Coarsen a mesh based on area constraints.

721:   Collective on Mesh

723:   Input Parameters:
724: + mesh    - The initial mesh
725: - area    - A function which gives an area constraint when evaluated inside an element

727:   Output Parameter:
728: . newmesh - The coarse mesh

730:   Note:
731:   If PETSC_NULL is used for the 'area' argument, then the mesh consisting only of vertices
732:   is returned.

734:   Level: beginner

736: .keywords: mesh, coarsening
737: .seealso: MeshRefine(), MeshDestroy()
738: @*/
739: int MeshCoarsen(Mesh mesh, PointFunction area, Mesh *newmesh)
740: {

746:   if (mesh->ops->coarsen == PETSC_NULL) SETERRQ(PETSC_ERR_SUP, " ");
747:   MeshPreCopy_Private(mesh);
748:   (*mesh->ops->coarsen)(mesh, area, newmesh);
749:   MeshPostCopy_Private(mesh, *newmesh);
750:   return(0);
751: }
752: 
753: /*@
754:   MeshRefine - Refine a mesh based on area constraints.

756:   Input Parameters:
757: + mesh    - The initial mesh
758: - area    - A function which gives an area constraint when evaluated inside an element

760:   Output Parameter:
761: . newmesh - The refined mesh

763:   Level: beginner

765: .keywords: mesh, refinement
766: .seealso: MeshCoarsen(), MeshDestroy()
767: @*/
768: int MeshRefine(Mesh mesh, PointFunction area, Mesh *newmesh)
769: {

775:   if (mesh->ops->refine == PETSC_NULL) SETERRQ(PETSC_ERR_SUP, " ");
776:   MeshPreCopy_Private(mesh);
777:   (*mesh->ops->refine)(mesh, area, newmesh);
778:   MeshPostCopy_Private(mesh, *newmesh);
779:   return(0);
780: }

782: /*@
783:   MeshResetNodes - Using the vertex and edge structure, move auxilliary nodes back to
784:   their default positions.

786:   Input Parameters:
787: + mesh    - The mesh
788: - resetBd - The flag which indicates whether boundaries should also be reset

790:   Note:
791:   This function was developed to allow midnodes to be moved back onto the edges
792:   with which they were associated. Midnodes on edges between nodes on the same
793:   boundary are only reset if resetBd == PETSC_TRUE (to allow for curved boundaries).

795:   Level: advanced

797: .keywords: mesh, reset, node
798: .seealso: MeshReform(), MeshSetBounday()
799: @*/
800: int MeshResetNodes(Mesh mesh, PetscTruth resetBd)
801: {

806:   if (!mesh->ops->resetnodes) {
807:     SETERRQ(PETSC_ERR_SUP, "Not supported by this Mesh object");
808:   }
809:   (*mesh->ops->resetnodes)(mesh, resetBd);
810:   return(0);
811: }

813: /*@
814:   MeshSaveMesh - This function saves the mesh coordinates.

816:   Collective on Mesh

818:   Input Parameter:
819: . mesh - The mesh

821:   Level: advanced

823:   Note:
824:   This function is meant to be used in conjuction with MeshMoveMesh(),
825:   and MeshRestoreMesh() so that the mesh may be moved, checked for
826:   distortion, and if necessary moved back.

828: .keywords: mesh, movement
829: .seealso: MeshMoveMesh(), MeshRestoreMesh()
830: @*/
831: int MeshSaveMesh(Mesh mesh) {

836:   (*mesh->ops->savemesh)(mesh);
837:   return(0);
838: }

840: /*@
841:   MeshRestoreMesh - This function restores the mesh coordinates which were
842:   previously saved.

844:   Collective on Mesh

846:   Input Parameter:
847: . mesh - The mesh

849:   Level: advanced

851:   Note:
852:   This function is meant to be used in conjuction with MeshMoveMesh(),
853:   and MeshRestoreMesh() so that the mesh may be moved, checked for
854:   distortion, and if necessary moved back.

856: .keywords: mesh, movement
857: .seealso: MeshMoveMesh(), MeshSaveMesh()
858: @*/
859: int MeshRestoreMesh(Mesh mesh) {

864:   (*mesh->ops->restoremesh)(mesh);
865:   return(0);
866: }

868: /*@
869:   MeshIsDistorted - This function checks for more than the maximum
870:   level of distortion in the mesh. It also returns an error when
871:   encountering elements with negative area.

873:   Collective on Mesh

875:   Input Parameter:
876: . mesh - The mesh

878:   Output Parameter:
879: . flag - Signals a distorted mesh

881:   Level: intermediate

883: .keywords: mesh, distortion, movement
884: .seealso: MeshMoveMesh()
885: @*/
886: int MeshIsDistorted(Mesh mesh, PetscTruth *flag) {

892:   PetscLogEventBegin(MESH_IsDistorted, mesh, 0, 0, 0);
893:   /* Do not use CHKERRQ since we want the return value at the top level */
894:   (*mesh->ops->isdistorted)(mesh, flag);
895:   PetscLogEventEnd(MESH_IsDistorted, mesh, 0, 0, 0);
896:   PetscFunctionReturn(ierr);
897: }