Actual source code: ex12.c
petsc-3.11.3 2019-06-26
1: static char help[] = "Test DMStag 2d star stencil\n\n";
2: #include <petscdm.h>
3: #include <petscdmstag.h>
5: int main(int argc,char **argv)
6: {
7: PetscErrorCode ierr;
8: DM dm;
9: Vec vec,vecLocal1,vecLocal2;
10: PetscScalar *a,***a1,***a2,expected,sum;
11: PetscInt startx,starty,nx,ny,i,j,d,is,js,dof0,dof1,dof2,dofTotal,stencilWidth,ngx,ngy;
12: DMBoundaryType boundaryTypex,boundaryTypey;
13: PetscMPIInt rank;
15: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
16: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
17: dof0 = 1;
18: dof1 = 1;
19: dof2 = 1;
20: stencilWidth = 2;
21: DMStagCreate2d(PETSC_COMM_WORLD,DM_BOUNDARY_PERIODIC,DM_BOUNDARY_PERIODIC,4,4,PETSC_DECIDE,PETSC_DECIDE,dof0,dof1,dof2,DMSTAG_STENCIL_STAR,stencilWidth,NULL,NULL,&dm);
22: DMSetFromOptions(dm);
23: DMSetUp(dm);
24: DMStagGetDOF(dm,&dof0,&dof1,&dof2,NULL);
25: dofTotal = dof0 + 2*dof1 + dof2;
26: DMStagGetStencilWidth(dm,&stencilWidth);
28: DMCreateLocalVector(dm,&vecLocal1);
29: VecDuplicate(vecLocal1,&vecLocal2);
31: DMCreateGlobalVector(dm,&vec);
32: VecSet(vec,1.0);
33: VecSet(vecLocal1,0.0);
34: DMGlobalToLocalBegin(dm,vec,INSERT_VALUES,vecLocal1);
35: DMGlobalToLocalEnd(dm,vec,INSERT_VALUES,vecLocal1);
37: DMStagGetCorners(dm,&startx,&starty,NULL,&nx,&ny,NULL,NULL,NULL,NULL);
38: DMStagVecGetArrayDOFRead(dm,vecLocal1,&a1);
39: DMStagVecGetArrayDOF(dm,vecLocal2,&a2);
40: for (j=starty; j<starty + ny; ++j) {
41: for (i=startx; i<startx + nx; ++i) {
42: for (d=0; d<dofTotal; ++d) {
43: if (a1[j][i][d] != 1.0) {
44: PetscPrintf(PETSC_COMM_SELF,"[%d] Unexpected value %g (expecting %g)\n",rank,a1[j][i][d],1.0);
45: }
46: a2[j][i][d] = 0.0;
47: for (js = -stencilWidth; js <= stencilWidth; ++js) {
48: a2[j][i][d] += a1[j+js][i][d];
49: }
50: for (is = -stencilWidth; is <= stencilWidth; ++is) {
51: a2[j][i][d] += a1[j][i+is][d];
52: }
53: a2[j][i][d] -= a1[j][i][d];
54: }
55: }
56: }
57: DMStagVecRestoreArrayDOFRead(dm,vecLocal1,&a1);
58: DMStagVecRestoreArrayDOF(dm,vecLocal2,&a2);
60: DMLocalToGlobalBegin(dm,vecLocal2,INSERT_VALUES,vec);
61: DMLocalToGlobalEnd(dm,vecLocal2,INSERT_VALUES,vec);
63: /* For the all-periodic case, some additional checks */
64: DMStagGetBoundaryTypes(dm,&boundaryTypex,&boundaryTypey,NULL);
65: if (boundaryTypex == DM_BOUNDARY_PERIODIC && boundaryTypey == DM_BOUNDARY_PERIODIC) {
67: DMStagGetGhostCorners(dm,NULL,NULL,NULL,&ngx,&ngy,NULL);
68: expected = (ngx*ngy - 4*stencilWidth*stencilWidth)*dofTotal;
69: VecSum(vecLocal1,&sum);
70: if (sum != expected) {
71: PetscPrintf(PETSC_COMM_SELF,"[%d] Unexpected sum of local entries %g (expected %g)\n",rank,sum,expected);
72: }
74: VecGetArray(vec,&a);
75: expected = 1 + 4*stencilWidth;
76: for (i=0; i<ny*nx*dofTotal; ++i) {
77: if (a[i] != expected) {
78: PetscPrintf(PETSC_COMM_SELF,"[%d] Unexpected value %g (expecting %g)\n",rank,a[i],expected);
79: }
80: }
81: VecRestoreArray(vec,&a);
82: }
84: VecDestroy(&vec);
85: VecDestroy(&vecLocal1);
86: VecDestroy(&vecLocal2);
87: DMDestroy(&dm);
88: PetscFinalize();
89: return ierr;
90: }
92: /*TEST
94: test:
95: suffix: 1
96: nsize: 4
97: args: -stag_ranks_x 2 -stag_ranks_y 2 -stag_stencil_width 1
99: test:
100: suffix: 2
101: nsize: 6
102: args: -stag_ranks_x 3 -stag_ranks_y 2 -stag_dof_0 2 -stag_grid_x 6
104: test:
105: suffix: 3
106: nsize: 4
107: args: -stag_dof_0 3 -stag_dof_1 2 -stag_dof_2 4 2 -stag_stencil_width 3 -stag_grid_x 6 -stag_grid_y 6
109: test:
110: suffix: 4
111: nsize: 4
112: args: -stag_stencil_width 1 -stag_grid_x 2 -stag_grid_y 2 -stag_boundary_type_x ghosted
114: test:
115: suffix: 5
116: nsize: 4
117: args: -stag_stencil_width 1 -stag_grid_x 2 -stag_grid_y 2 -stag_boundary_type_y ghosted
119: test:
120: suffix: 6
121: nsize: 4
122: args: -stag_stencil_width 1 -stag_grid_x 3 -stag_grid_y 2 -stag_boundary_type_x ghosted -stag_boundary_type_y ghosted
124: test:
125: suffix: 7
126: nsize: 4
127: args: -stag_stencil_width 1 -stag_grid_x 2 -stag_grid_y 2 -stag_boundary_type_y ghosted
129: test:
130: suffix: 8
131: nsize: 6
132: args: -stag_stencil_width 1 -stag_grid_y 2 -stag_grid_x 19 -stag_boundary_type_y ghosted -stag_ranks_x 6
133: TEST*/