Actual source code: ex41.c
2: static char help[] = "Tests mirror boundary conditions in 3-d.\n\n";
4: #include <petscdm.h>
5: #include <petscdmda.h>
7: int main(int argc,char **argv)
8: {
9: PetscInt M = 2, N = 3, P = 4,stencil_width = 1, dof = 1,m,n,p,xstart,ystart,zstart,i,j,k,c;
10: DM da;
11: Vec global,local;
12: PetscScalar ****vglobal;
13: PetscViewer sview;
14: PetscScalar sum;
16: PetscInitialize(&argc,&argv,(char*)0,help);
17: PetscOptionsGetInt(NULL,0,"-stencil_width",&stencil_width,0);
18: PetscOptionsGetInt(NULL,0,"-dof",&dof,0);
20: DMDACreate3d(PETSC_COMM_WORLD,DM_BOUNDARY_MIRROR,DM_BOUNDARY_MIRROR,DM_BOUNDARY_MIRROR,DMDA_STENCIL_STAR,M,N,P,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,dof,stencil_width,NULL,NULL,NULL,&da);
21: DMSetFromOptions(da);
22: DMSetUp(da);
23: DMDAGetCorners(da,&xstart,&ystart,&zstart,&m,&n,&p);
25: DMCreateGlobalVector(da,&global);
26: DMDAVecGetArrayDOF(da,global,&vglobal);
27: for (k=zstart; k<zstart+p; k++) {
28: for (j=ystart; j<ystart+n; j++) {
29: for (i=xstart; i<xstart+m; i++) {
30: for (c=0; c<dof; c++) {
31: vglobal[k][j][i][c] = 1000*k + 100*j + 10*i + c;
32: }
33: }
34: }
35: }
36: DMDAVecRestoreArrayDOF(da,global,&vglobal);
38: DMCreateLocalVector(da,&local);
39: DMGlobalToLocalBegin(da,global,ADD_VALUES,local);
40: DMGlobalToLocalEnd(da,global,ADD_VALUES,local);
42: VecSum(local,&sum);
43: PetscSynchronizedPrintf(PETSC_COMM_WORLD,"sum %g\n",(double)PetscRealPart(sum));
44: PetscSynchronizedFlush(PETSC_COMM_WORLD,stdout);
45: PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sview);
46: VecView(local,sview);
47: PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sview);
48: PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD);
50: DMDestroy(&da);
51: VecDestroy(&local);
52: VecDestroy(&global);
54: PetscFinalize();
55: return 0;
56: }
58: /*TEST
60: test:
62: test:
63: suffix: 2
64: nsize: 3
66: TEST*/