Actual source code: ex2.c
1: /*$Id: ex2.c,v 1.43 2001/08/07 03:04:42 balay Exp $*/
3: static char help[] = "Tests various 1-dimensional DA routines.\n\n";
5: #include petscda.h
6: #include petscsys.h
10: int main(int argc,char **argv)
11: {
12: int rank,M = 13,ierr,w=1,s=1,wrap=1;
13: DA da;
14: PetscViewer viewer;
15: Vec local,global;
16: PetscScalar value;
17: PetscDraw draw;
18: PetscTruth flg;
20: PetscInitialize(&argc,&argv,(char*)0,help);
22: /* Create viewers */
23: PetscViewerDrawOpen(PETSC_COMM_WORLD,0,"",280,480,600,200,&viewer);
24: PetscViewerDrawGetDraw(viewer,0,&draw);
25: PetscDrawSetDoubleBuffer(draw);
27: /* Readoptions */
28: PetscOptionsGetInt(PETSC_NULL,"-M",&M,PETSC_NULL);
29: PetscOptionsGetInt(PETSC_NULL,"-w",&w,PETSC_NULL);
30: PetscOptionsGetInt(PETSC_NULL,"-s",&s,PETSC_NULL);
32: /* Create distributed array and get vectors */
33: DACreate1d(PETSC_COMM_WORLD,(DAPeriodicType)wrap,M,w,s,PETSC_NULL,&da);
34: DAView(da,viewer);
35: DACreateGlobalVector(da,&global);
36: DACreateLocalVector(da,&local);
38: /* Set global vector; send ghost points to local vectors */
39: value = 1;
40: VecSet(&value,global);
41: DAGlobalToLocalBegin(da,global,INSERT_VALUES,local);
42: DAGlobalToLocalEnd(da,global,INSERT_VALUES,local);
44: /* Scale local vectors according to processor rank; pass to global vector */
45: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
46: value = rank+1;
47: VecScale(&value,local);
48: DALocalToGlobal(da,local,INSERT_VALUES,global);
50: VecView(global,viewer);
51: PetscPrintf(PETSC_COMM_WORLD,"\nGlobal Vector:\n");
52: VecView(global,PETSC_VIEWER_STDOUT_WORLD);
53: PetscPrintf(PETSC_COMM_WORLD,"\n");
55: /* Send ghost points to local vectors */
56: DAGlobalToLocalBegin(da,global,INSERT_VALUES,local);
57: DAGlobalToLocalEnd(da,global,INSERT_VALUES,local);
59: PetscOptionsHasName(PETSC_NULL,"-local_print",&flg);
60: if (flg) {
61: PetscViewer sviewer;
62: PetscSynchronizedPrintf(PETSC_COMM_WORLD,"\nLocal Vector: processor %d\n",rank);
63: PetscViewerGetSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
64: VecView(local,sviewer);
65: PetscViewerRestoreSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
66: PetscSynchronizedFlush(PETSC_COMM_WORLD);
67: }
69: /* Free memory */
70: PetscViewerDestroy(viewer);
71: VecDestroy(global);
72: VecDestroy(local);
73: DADestroy(da);
74: PetscFinalize();
75: return 0;
76: }
77: