Actual source code: ex2.c
1: /*$Id: ex2.c,v 1.15 2001/08/07 03:04:45 balay Exp $*/
3: static char help[] = "Tests DAGlobalToNaturalAllCreate() using contour plotting for 2d DAs.nn";
5: #include petscda.h
6: #include petscsys.h
8: int main(int argc,char **argv)
9: {
10: int i,j,rank,M = 10,N = 8,m = PETSC_DECIDE,n = PETSC_DECIDE,ierr;
11: PetscTruth flg;
12: DA da;
13: PetscViewer viewer;
14: Vec localall,global;
15: PetscScalar value,*vlocal;
16: DAPeriodicType ptype = DA_NONPERIODIC;
17: DAStencilType stype = DA_STENCIL_BOX;
18: VecScatter tolocalall,fromlocalall;
20: PetscInitialize(&argc,&argv,(char*)0,help);
21: PetscViewerDrawOpen(PETSC_COMM_WORLD,0,"",300,0,300,300,&viewer);
23: /* Read options */
24: PetscOptionsGetInt(PETSC_NULL,"-M",&M,PETSC_NULL);
25: PetscOptionsGetInt(PETSC_NULL,"-N",&N,PETSC_NULL);
26: PetscOptionsGetInt(PETSC_NULL,"-m",&m,PETSC_NULL);
27: PetscOptionsGetInt(PETSC_NULL,"-n",&n,PETSC_NULL);
28: PetscOptionsHasName(PETSC_NULL,"-star_stencil",&flg);
29: if (flg) stype = DA_STENCIL_STAR;
31: /* Create distributed array and get vectors */
32: DACreate2d(PETSC_COMM_WORLD,ptype,stype,
33: M,N,m,n,1,1,PETSC_NULL,PETSC_NULL,&da);
34: DACreateGlobalVector(da,&global);
35: VecCreateSeq(PETSC_COMM_SELF,M*N,&localall);
37: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
38: value = 5.0*rank;
39: VecSet(&value,global);
41: VecView(global,viewer);
43: /*
44: Create Scatter from global DA parallel vector to local vector that
45: contains all entries
46: */
47: DAGlobalToNaturalAllCreate(da,&tolocalall);
48: DANaturalAllToGlobalCreate(da,&fromlocalall);
50: VecScatterBegin(global,localall,INSERT_VALUES,SCATTER_FORWARD,tolocalall);
51: VecScatterEnd(global,localall,INSERT_VALUES,SCATTER_FORWARD,tolocalall);
53: VecGetArray(localall,&vlocal);
54: for (j=0; j<N; j++) {
55: for (i=0; i<M; i++) {
56: *vlocal++ += i + j*M;
57: }
58: }
59: VecRestoreArray(localall,&vlocal);
61: /* scatter back to global vector */
62: VecScatterBegin(localall,global,INSERT_VALUES,SCATTER_FORWARD,fromlocalall);
63: VecScatterEnd(localall,global,INSERT_VALUES,SCATTER_FORWARD,fromlocalall);
65: VecView(global,viewer);
67: /* Free memory */
68: VecScatterDestroy(tolocalall);
69: VecScatterDestroy(fromlocalall);
70: PetscViewerDestroy(viewer);
71: VecDestroy(localall);
72: VecDestroy(global);
73: DADestroy(da);
74: PetscFinalize();
75: return 0;
76: }
77: