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