Actual source code: ex11.c
petsc-3.4.2 2013-07-02
2: static char help[] = "Tests various 2-dimensional DMDA routines.\n\n";
4: #include <petscdmda.h>
5: #include <petscdraw.h>
9: int main(int argc,char **argv)
10: {
11: PetscInt M = 5,N = 4,dof=1,s=1,bx=0,by=0,i,n,j,k,m,cnt,wrap;
13: DM da;
14: PetscViewer viewer;
15: Vec local,locala,global,coors;
16: PetscScalar *xy,*alocal;
17: PetscDraw draw;
18: char fname[16];
20: PetscInitialize(&argc,&argv,(char*)0,help);
22: /* Create viewers */
23: PetscViewerDrawOpen(PETSC_COMM_WORLD,0,"",PETSC_DECIDE,PETSC_DECIDE,600,200,&viewer);
24: PetscViewerDrawGetDraw(viewer,0,&draw);
25: PetscDrawSetDoubleBuffer(draw);
27: /* Read options */
28: PetscOptionsGetInt(NULL,"-M",&M,NULL);
29: PetscOptionsGetInt(NULL,"-N",&N,NULL);
30: PetscOptionsGetInt(NULL,"-dof",&dof,NULL);
31: PetscOptionsGetInt(NULL,"-s",&s,NULL);
32: PetscOptionsGetInt(NULL,"-periodic_x",&wrap,NULL);
33: PetscOptionsGetInt(NULL,"-periodic_y",&wrap,NULL);
35: /* Create distributed array and get vectors */
36: DMDACreate2d(PETSC_COMM_WORLD,(DMDABoundaryType)bx,(DMDABoundaryType)by,DMDA_STENCIL_BOX,M,N,PETSC_DECIDE,
37: PETSC_DECIDE,dof,s,NULL,NULL,&da);
38: DMDASetUniformCoordinates(da,0.0,1.0,0.0,1.0,0.0,0.0);
39: for (i=0; i<dof; i++) {
40: sprintf(fname,"Field %d",(int)i);
41: DMDASetFieldName(da,i,fname);
42: }
44: DMView(da,viewer);
45: DMCreateGlobalVector(da,&global);
46: DMCreateLocalVector(da,&local);
47: DMCreateLocalVector(da,&locala);
48: DMGetCoordinates(da,&coors);
49: VecGetArray(coors,&xy);
51: VecView(coors,PETSC_VIEWER_STDOUT_SELF);
53: /* Set values into local vectors */
54: VecGetArray(local,&alocal);
55: DMDAGetGhostCorners(da,0,0,0,&m,&n,0);
56: n = n/dof;
57: for (k=0; k<dof; k++) {
58: cnt = 0;
59: for (j=0; j<n; j++) {
60: for (i=0; i<m; i++) {
61: alocal[k+dof*cnt] = PetscSinScalar(2.0*PETSC_PI*(k+1)*xy[2*cnt]);
62: cnt++;
63: }
64: }
65: }
66: VecRestoreArray(local,&alocal);
67: VecRestoreArray(coors,&xy);
69: DMLocalToGlobalBegin(da,local,INSERT_VALUES,global);
70: DMLocalToGlobalEnd(da,local,INSERT_VALUES,global);
72: VecView(global,viewer);
74: /* Send ghost points to local vectors */
75: DMGlobalToLocalBegin(da,global,INSERT_VALUES,locala);
76: DMGlobalToLocalEnd(da,global,INSERT_VALUES,locala);
78: /* Free memory */
79: PetscViewerDestroy(&viewer);
80: VecDestroy(&global);
81: VecDestroy(&local);
82: VecDestroy(&locala);
83: DMDestroy(&da);
84: PetscFinalize();
85: return 0;
86: }