1 #include <petsc/private/dmpleximpl.h> /*I "petscdmplex.h" I*/
2 #include <../src/sys/classes/viewer/impls/vtk/vtkvimpl.h>
3
DMPlexVTKGetCellType_Internal(DM dm,PetscInt dim,PetscInt corners,PetscInt * cellType)4 PetscErrorCode DMPlexVTKGetCellType_Internal(DM dm, PetscInt dim, PetscInt corners, PetscInt *cellType)
5 {
6 PetscFunctionBegin;
7 *cellType = -1;
8 switch (dim) {
9 case 0:
10 switch (corners) {
11 case 1:
12 *cellType = 1; /* VTK_VERTEX */
13 break;
14 default:
15 break;
16 }
17 break;
18 case 1:
19 switch (corners) {
20 case 2:
21 *cellType = 3; /* VTK_LINE */
22 break;
23 case 3:
24 *cellType = 21; /* VTK_QUADRATIC_EDGE */
25 break;
26 default:
27 break;
28 }
29 break;
30 case 2:
31 switch (corners) {
32 case 3:
33 *cellType = 5; /* VTK_TRIANGLE */
34 break;
35 case 4:
36 *cellType = 9; /* VTK_QUAD */
37 break;
38 case 6:
39 *cellType = 22; /* VTK_QUADRATIC_TRIANGLE */
40 break;
41 case 9:
42 *cellType = 23; /* VTK_QUADRATIC_QUAD */
43 break;
44 default:
45 break;
46 }
47 break;
48 case 3:
49 switch (corners) {
50 case 4:
51 *cellType = 10; /* VTK_TETRA */
52 break;
53 case 5:
54 *cellType = 14; /* VTK_PYRAMID */
55 break;
56 case 6:
57 *cellType = 13; /* VTK_WEDGE */
58 break;
59 case 8:
60 *cellType = 12; /* VTK_HEXAHEDRON */
61 break;
62 case 10:
63 *cellType = 24; /* VTK_QUADRATIC_TETRA */
64 break;
65 case 27:
66 *cellType = 29; /* VTK_QUADRATIC_HEXAHEDRON */
67 break;
68 default:
69 break;
70 }
71 }
72 PetscFunctionReturn(PETSC_SUCCESS);
73 }
74
75 /*@
76 DMPlexVTKWriteAll - Write a file containing all the fields that have been provided to the viewer
77
78 Collective
79
80 Input Parameters:
81 + odm - The `DMPLEX` specifying the mesh, passed as a `PetscObject`
82 - viewer - viewer of type `PETSCVIEWERVTK`
83
84 Level: developer
85
86 Note:
87 This function is a callback used by the VTK viewer to actually write the file.
88 The reason for this odd model is that the VTK file format does not provide any way to write one field at a time.
89 Instead, metadata for the entire file needs to be available up-front before you can start writing the file.
90
91 .seealso: [](ch_unstructured), `DM`, `PETSCVIEWEREXODUSII`, `DMPLEX`, `PETSCVIEWERVTK`
92 @*/
DMPlexVTKWriteAll(PetscObject odm,PetscViewer viewer)93 PetscErrorCode DMPlexVTKWriteAll(PetscObject odm, PetscViewer viewer)
94 {
95 DM dm = (DM)odm;
96
97 PetscFunctionBegin;
98 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
99 PetscValidHeaderSpecificType(viewer, PETSC_VIEWER_CLASSID, 2, PETSCVIEWERVTK);
100 switch (viewer->format) {
101 case PETSC_VIEWER_DEFAULT:
102 case PETSC_VIEWER_VTK_VTU:
103 PetscCall(DMPlexVTKWriteAll_VTU(dm, viewer));
104 break;
105 default:
106 SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "No support for format '%s'", PetscViewerFormats[viewer->format]);
107 }
108 PetscFunctionReturn(PETSC_SUCCESS);
109 }
110