Actual source code: ex47.c
petsc-3.13.6 2020-09-29
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;
21: PetscErrorCode ierr;
23: 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);
24: DMSetFromOptions(da);
25: DMSetUp(da);
27: DMDASetUniformCoordinates(da,0.0,Lx,0.0,Ly,0.0,Lz);
28: DMDAGetLocalInfo(da,&info);
29: DMCreateGlobalVector(da,&v);
30: DMDAVecGetArray(da,v,&va);
31: for (k=info.zs; k<info.zs+info.zm; k++) {
32: for (j=info.ys; j<info.ys+info.ym; j++) {
33: for (i=info.xs; i<info.xs+info.xm; i++) {
34: PetscScalar x = (Lx*i)/M;
35: PetscScalar y = (Ly*j)/N;
36: PetscScalar z = (Lz*k)/P;
37: va[k][j][i] = PetscPowScalarInt(x-0.5*Lx,2)+PetscPowScalarInt(y-0.5*Ly,2)+PetscPowScalarInt(z-0.5*Lz,2);
38: }
39: }
40: }
41: DMDAVecRestoreArray(da,v,&va);
42: PetscViewerVTKOpen(comm,filename,FILE_MODE_WRITE,&view);
43: VecView(v,view);
44: PetscViewerDestroy(&view);
45: VecDestroy(&v);
46: DMDestroy(&da);
47: return 0;
48: }
51: /*
52: Write 2D DMDA vector with coordinates in VTK .vts format
54: */
55: PetscErrorCode test_2d(const char filename[])
56: {
57: MPI_Comm comm = MPI_COMM_WORLD;
58: const PetscInt M=10,N=20,dof=1,sw=1;
59: const PetscScalar Lx=1.0,Ly=1.0,Lz=1.0;
60: DM da;
61: Vec v;
62: PetscViewer view;
63: DMDALocalInfo info;
64: PetscScalar **va;
65: PetscInt i,j;
66: PetscErrorCode ierr;
68: DMDACreate2d(comm,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_STAR, M,N,PETSC_DECIDE,PETSC_DECIDE,dof,sw,NULL,NULL,&da);
69: DMSetFromOptions(da);
70: DMSetUp(da);
71: DMDASetUniformCoordinates(da,0.0,Lx,0.0,Ly,0.0,Lz);
72: DMDAGetLocalInfo(da,&info);
73: DMCreateGlobalVector(da,&v);
74: DMDAVecGetArray(da,v,&va);
75: for (j=info.ys; j<info.ys+info.ym; j++) {
76: for (i=info.xs; i<info.xs+info.xm; i++) {
77: PetscScalar x = (Lx*i)/M;
78: PetscScalar y = (Ly*j)/N;
79: va[j][i] = PetscPowScalarInt(x-0.5*Lx,2)+PetscPowScalarInt(y-0.5*Ly,2);
80: }
81: }
82: DMDAVecRestoreArray(da,v,&va);
83: PetscViewerVTKOpen(comm,filename,FILE_MODE_WRITE,&view);
84: VecView(v,view);
85: PetscViewerDestroy(&view);
86: VecDestroy(&v);
87: DMDestroy(&da);
88: return 0;
89: }
92: /*
93: Write 2D DMDA vector without coordinates in VTK .vts format
95: */
96: PetscErrorCode test_2d_nocoord(const char filename[])
97: {
98: MPI_Comm comm = MPI_COMM_WORLD;
99: const PetscInt M=10,N=20,dof=1,sw=1;
100: const PetscScalar Lx=1.0,Ly=1.0;
101: DM da;
102: Vec v;
103: PetscViewer view;
104: DMDALocalInfo info;
105: PetscScalar **va;
106: PetscInt i,j;
107: PetscErrorCode ierr;
109: DMDACreate2d(comm,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_STAR,M,N,PETSC_DECIDE,PETSC_DECIDE,dof,sw,NULL,NULL,&da);
110: DMSetFromOptions(da);
111: DMSetUp(da);
112: DMDAGetLocalInfo(da,&info);
113: DMCreateGlobalVector(da,&v);
114: DMDAVecGetArray(da,v,&va);
115: for (j=info.ys; j<info.ys+info.ym; j++) {
116: for (i=info.xs; i<info.xs+info.xm; i++) {
117: PetscScalar x = (Lx*i)/M;
118: PetscScalar y = (Ly*j)/N;
119: va[j][i] = PetscPowScalarInt(x-0.5*Lx,2)+PetscPowScalarInt(y-0.5*Ly,2);
120: }
121: }
122: DMDAVecRestoreArray(da,v,&va);
123: PetscViewerVTKOpen(comm,filename,FILE_MODE_WRITE,&view);
124: VecView(v,view);
125: PetscViewerDestroy(&view);
126: VecDestroy(&v);
127: DMDestroy(&da);
128: return 0;
129: }
132: /*
133: Write 3D DMDA vector without coordinates in VTK .vts format
135: */
136: PetscErrorCode test_3d_nocoord(const char filename[])
137: {
138: MPI_Comm comm = MPI_COMM_WORLD;
139: const PetscInt M=10,N=20,P=30,dof=1,sw=1;
140: const PetscScalar Lx=1.0,Ly=1.0,Lz=1.0;
141: DM da;
142: Vec v;
143: PetscViewer view;
144: DMDALocalInfo info;
145: PetscScalar ***va;
146: PetscInt i,j,k;
147: PetscErrorCode ierr;
149: 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);
150: DMSetFromOptions(da);
151: DMSetUp(da);
153: DMDAGetLocalInfo(da,&info);
154: DMCreateGlobalVector(da,&v);
155: DMDAVecGetArray(da,v,&va);
156: for (k=info.zs; k<info.zs+info.zm; k++) {
157: for (j=info.ys; j<info.ys+info.ym; j++) {
158: for (i=info.xs; i<info.xs+info.xm; i++) {
159: PetscScalar x = (Lx*i)/M;
160: PetscScalar y = (Ly*j)/N;
161: PetscScalar z = (Lz*k)/P;
162: va[k][j][i] = PetscPowScalarInt(x-0.5*Lx,2)+PetscPowScalarInt(y-0.5*Ly,2)+PetscPowScalarInt(z-0.5*Lz,2);
163: }
164: }
165: }
166: DMDAVecRestoreArray(da,v,&va);
167: PetscViewerVTKOpen(comm,filename,FILE_MODE_WRITE,&view);
168: VecView(v,view);
169: PetscViewerDestroy(&view);
170: VecDestroy(&v);
171: DMDestroy(&da);
172: return 0;
173: }
175: int main(int argc, char *argv[])
176: {
179: PetscInitialize(&argc,&argv,0,help);if (ierr) return ierr;
180: test_3d("3d.vts");
181: test_2d("2d.vts");
182: test_2d_nocoord("2d_nocoord.vts");
183: test_3d_nocoord("3d_nocoord.vts");
184: PetscFinalize();
185: return ierr;
186: }
189: /*TEST
191: build:
192: requires: !complex
194: test:
195: nsize: 2
197: TEST*/