Actual source code: ex39.c
petsc-3.8.4 2018-03-24
2: static char help[] = "Tests mirror boundary conditions in 1-d.\n\n";
4: #include <petscdm.h>
5: #include <petscdmda.h>
7: int main(int argc,char **argv)
8: {
10: PetscInt M = 6,stencil_width = 1, dof = 1,m,xstart,i,j;
11: DM da;
12: Vec global,local;
13: PetscScalar **vglobal;
15: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
16: PetscOptionsGetInt(NULL,0,"-stencil_width",&stencil_width,0);
17: PetscOptionsGetInt(NULL,0,"-dof",&dof,0);
19: DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_MIRROR,M,dof,stencil_width,NULL,&da);
20: DMSetFromOptions(da);
21: DMSetUp(da);
22: DMDAGetCorners(da,&xstart,0,0,&m,0,0);
24: DMCreateGlobalVector(da,&global);
25: DMDAVecGetArrayDOF(da,global,&vglobal);
26: for (i=xstart; i<xstart+m; i++) {
27: for (j=0; j<dof; j++) {
28: vglobal[i][j] = 100*(i+1) + j;
29: }
30: }
31: DMDAVecRestoreArrayDOF(da,global,&vglobal);
33: DMCreateLocalVector(da,&local);
34: DMGlobalToLocalBegin(da,global,INSERT_VALUES,local);
35: DMGlobalToLocalEnd(da,global,INSERT_VALUES,local);
37: PetscSequentialPhaseBegin(PETSC_COMM_WORLD,1);
38: VecView(local,PETSC_VIEWER_STDOUT_SELF);
39: PetscSequentialPhaseEnd(PETSC_COMM_WORLD,1);
40: VecView(global,PETSC_VIEWER_STDOUT_WORLD);
42: DMDestroy(&da);
43: VecDestroy(&local);
44: VecDestroy(&global);
46: PetscFinalize();
47: return ierr;
48: }