Actual source code: ex47.c
1: static char help[] = "Test VTK structured grid (.vts) viewer support\n\n";
3: #include <petscdm.h>
4: #include <petscdmda.h>
6: /*
7: Write 3D DMDA vector with coordinates in VTK .vts format
9: */
10: PetscErrorCode test_3d(const char filename[])
11: {
12: MPI_Comm comm = MPI_COMM_WORLD;
13: const PetscInt M = 10, N = 15, P = 30, dof = 1, sw = 1;
14: const PetscScalar Lx = 1.0, Ly = 1.0, Lz = 1.0;
15: DM da;
16: Vec v;
17: PetscViewer view;
18: DMDALocalInfo info;
19: PetscScalar ***va;
20: PetscInt i, j, k;
22: DMDACreate3d(comm, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_STAR, M, N, P, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, dof, sw, NULL, NULL, NULL, &da);
23: DMSetFromOptions(da);
24: DMSetUp(da);
26: DMDASetUniformCoordinates(da, 0.0, Lx, 0.0, Ly, 0.0, Lz);
27: DMDAGetLocalInfo(da, &info);
28: DMCreateGlobalVector(da, &v);
29: DMDAVecGetArray(da, v, &va);
30: for (k = info.zs; k < info.zs + info.zm; k++) {
31: for (j = info.ys; j < info.ys + info.ym; j++) {
32: for (i = info.xs; i < info.xs + info.xm; i++) {
33: PetscScalar x = (Lx * i) / M;
34: PetscScalar y = (Ly * j) / N;
35: PetscScalar z = (Lz * k) / P;
36: va[k][j][i] = PetscPowScalarInt(x - 0.5 * Lx, 2) + PetscPowScalarInt(y - 0.5 * Ly, 2) + PetscPowScalarInt(z - 0.5 * Lz, 2);
37: }
38: }
39: }
40: DMDAVecRestoreArray(da, v, &va);
41: PetscViewerVTKOpen(comm, filename, FILE_MODE_WRITE, &view);
42: VecView(v, view);
43: PetscViewerDestroy(&view);
44: VecDestroy(&v);
45: DMDestroy(&da);
46: return 0;
47: }
49: /*
50: Write 2D DMDA vector with coordinates in VTK .vts format
52: */
53: PetscErrorCode test_2d(const char filename[])
54: {
55: MPI_Comm comm = MPI_COMM_WORLD;
56: const PetscInt M = 10, N = 20, dof = 1, sw = 1;
57: const PetscScalar Lx = 1.0, Ly = 1.0, Lz = 1.0;
58: DM da;
59: Vec v;
60: PetscViewer view;
61: DMDALocalInfo info;
62: PetscScalar **va;
63: PetscInt i, j;
65: DMDACreate2d(comm, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_STAR, M, N, PETSC_DECIDE, PETSC_DECIDE, dof, sw, NULL, NULL, &da);
66: DMSetFromOptions(da);
67: DMSetUp(da);
68: DMDASetUniformCoordinates(da, 0.0, Lx, 0.0, Ly, 0.0, Lz);
69: DMDAGetLocalInfo(da, &info);
70: DMCreateGlobalVector(da, &v);
71: DMDAVecGetArray(da, v, &va);
72: for (j = info.ys; j < info.ys + info.ym; j++) {
73: for (i = info.xs; i < info.xs + info.xm; i++) {
74: PetscScalar x = (Lx * i) / M;
75: PetscScalar y = (Ly * j) / N;
76: va[j][i] = PetscPowScalarInt(x - 0.5 * Lx, 2) + PetscPowScalarInt(y - 0.5 * Ly, 2);
77: }
78: }
79: DMDAVecRestoreArray(da, v, &va);
80: PetscViewerVTKOpen(comm, filename, FILE_MODE_WRITE, &view);
81: VecView(v, view);
82: PetscViewerDestroy(&view);
83: VecDestroy(&v);
84: DMDestroy(&da);
85: return 0;
86: }
88: /*
89: Write 2D DMDA vector without coordinates in VTK .vts format
91: */
92: PetscErrorCode test_2d_nocoord(const char filename[])
93: {
94: MPI_Comm comm = MPI_COMM_WORLD;
95: const PetscInt M = 10, N = 20, dof = 1, sw = 1;
96: const PetscScalar Lx = 1.0, Ly = 1.0;
97: DM da;
98: Vec v;
99: PetscViewer view;
100: DMDALocalInfo info;
101: PetscScalar **va;
102: PetscInt i, j;
104: DMDACreate2d(comm, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_STAR, M, N, PETSC_DECIDE, PETSC_DECIDE, dof, sw, NULL, NULL, &da);
105: DMSetFromOptions(da);
106: DMSetUp(da);
107: DMDAGetLocalInfo(da, &info);
108: DMCreateGlobalVector(da, &v);
109: DMDAVecGetArray(da, v, &va);
110: for (j = info.ys; j < info.ys + info.ym; j++) {
111: for (i = info.xs; i < info.xs + info.xm; i++) {
112: PetscScalar x = (Lx * i) / M;
113: PetscScalar y = (Ly * j) / N;
114: va[j][i] = PetscPowScalarInt(x - 0.5 * Lx, 2) + PetscPowScalarInt(y - 0.5 * Ly, 2);
115: }
116: }
117: DMDAVecRestoreArray(da, v, &va);
118: PetscViewerVTKOpen(comm, filename, FILE_MODE_WRITE, &view);
119: VecView(v, view);
120: PetscViewerDestroy(&view);
121: VecDestroy(&v);
122: DMDestroy(&da);
123: return 0;
124: }
126: /*
127: Write 3D DMDA vector without coordinates in VTK .vts format
129: */
130: PetscErrorCode test_3d_nocoord(const char filename[])
131: {
132: MPI_Comm comm = MPI_COMM_WORLD;
133: const PetscInt M = 10, N = 20, P = 30, dof = 1, sw = 1;
134: const PetscScalar Lx = 1.0, Ly = 1.0, Lz = 1.0;
135: DM da;
136: Vec v;
137: PetscViewer view;
138: DMDALocalInfo info;
139: PetscScalar ***va;
140: PetscInt i, j, k;
142: DMDACreate3d(comm, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_STAR, M, N, P, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, dof, sw, NULL, NULL, NULL, &da);
143: DMSetFromOptions(da);
144: DMSetUp(da);
146: DMDAGetLocalInfo(da, &info);
147: DMCreateGlobalVector(da, &v);
148: DMDAVecGetArray(da, v, &va);
149: for (k = info.zs; k < info.zs + info.zm; k++) {
150: for (j = info.ys; j < info.ys + info.ym; j++) {
151: for (i = info.xs; i < info.xs + info.xm; i++) {
152: PetscScalar x = (Lx * i) / M;
153: PetscScalar y = (Ly * j) / N;
154: PetscScalar z = (Lz * k) / P;
155: va[k][j][i] = PetscPowScalarInt(x - 0.5 * Lx, 2) + PetscPowScalarInt(y - 0.5 * Ly, 2) + PetscPowScalarInt(z - 0.5 * Lz, 2);
156: }
157: }
158: }
159: DMDAVecRestoreArray(da, v, &va);
160: PetscViewerVTKOpen(comm, filename, FILE_MODE_WRITE, &view);
161: VecView(v, view);
162: PetscViewerDestroy(&view);
163: VecDestroy(&v);
164: DMDestroy(&da);
165: return 0;
166: }
168: int main(int argc, char *argv[])
169: {
171: PetscInitialize(&argc, &argv, 0, help);
172: test_3d("3d.vts");
173: test_2d("2d.vts");
174: test_2d_nocoord("2d_nocoord.vts");
175: test_3d_nocoord("3d_nocoord.vts");
176: PetscFinalize();
177: return 0;
178: }
180: /*TEST
182: build:
183: requires: !complex
185: test:
186: nsize: 2
188: TEST*/