Actual source code: ex9.c
2: static char help[] = "Demonstrates use of VecCreateGhost().\n\n";
4: /*
5: Description: Ghost padding is one way to handle local calculations that
6: involve values from other processors. VecCreateGhost() provides
7: a way to create vectors with extra room at the end of the vector
8: array to contain the needed ghost values from other processors,
9: vector computations are otherwise unaffected.
10: */
12: /*
13: Include "petscvec.h" so that we can use vectors. Note that this file
14: automatically includes:
15: petscsys.h - base PETSc routines petscis.h - index sets
16: petscviewer.h - viewers
17: */
18: #include <petscvec.h>
20: int main(int argc,char **argv)
21: {
22: PetscMPIInt rank,size;
23: PetscInt nlocal = 6,nghost = 2,ifrom[2],i,rstart,rend;
24: PetscBool flg,flg2,flg3;
25: PetscScalar value,*array,*tarray=0;
26: Vec lx,gx,gxs;
28: PetscInitialize(&argc,&argv,(char*)0,help);
29: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
30: MPI_Comm_size(PETSC_COMM_WORLD,&size);
33: /*
34: Construct a two dimensional graph connecting nlocal degrees of
35: freedom per processor. From this we will generate the global
36: indices of needed ghost values
38: For simplicity we generate the entire graph on each processor:
39: in real application the graph would stored in parallel, but this
40: example is only to demonstrate the management of ghost padding
41: with VecCreateGhost().
43: In this example we consider the vector as representing
44: degrees of freedom in a one dimensional grid with periodic
45: boundary conditions.
47: ----Processor 1--------- ----Processor 2 --------
48: 0 1 2 3 4 5 6 7 8 9 10 11
49: |----|
50: |-------------------------------------------------|
52: */
54: if (rank == 0) {
55: ifrom[0] = 11; ifrom[1] = 6;
56: } else {
57: ifrom[0] = 0; ifrom[1] = 5;
58: }
60: /*
61: Create the vector with two slots for ghost points. Note that both
62: the local vector (lx) and the global vector (gx) share the same
63: array for storing vector values.
64: */
65: PetscOptionsHasName(NULL,NULL,"-allocate",&flg);
66: PetscOptionsHasName(NULL,NULL,"-vecmpisetghost",&flg2);
67: PetscOptionsHasName(NULL,NULL,"-minvalues",&flg3);
68: if (flg) {
69: PetscMalloc1(nlocal+nghost,&tarray);
70: VecCreateGhostWithArray(PETSC_COMM_WORLD,nlocal,PETSC_DECIDE,nghost,ifrom,tarray,&gxs);
71: } else if (flg2) {
72: VecCreate(PETSC_COMM_WORLD,&gxs);
73: VecSetType(gxs,VECMPI);
74: VecSetSizes(gxs,nlocal,PETSC_DECIDE);
75: VecMPISetGhost(gxs,nghost,ifrom);
76: } else {
77: VecCreateGhost(PETSC_COMM_WORLD,nlocal,PETSC_DECIDE,nghost,ifrom,&gxs);
78: }
80: /*
81: Test VecDuplicate()
82: */
83: VecDuplicate(gxs,&gx);
84: VecDestroy(&gxs);
86: /*
87: Access the local representation
88: */
89: VecGhostGetLocalForm(gx,&lx);
91: /*
92: Set the values from 0 to 12 into the "global" vector
93: */
94: VecGetOwnershipRange(gx,&rstart,&rend);
95: for (i=rstart; i<rend; i++) {
96: value = (PetscScalar) i;
97: VecSetValues(gx,1,&i,&value,INSERT_VALUES);
98: }
99: VecAssemblyBegin(gx);
100: VecAssemblyEnd(gx);
102: VecGhostUpdateBegin(gx,INSERT_VALUES,SCATTER_FORWARD);
103: VecGhostUpdateEnd(gx,INSERT_VALUES,SCATTER_FORWARD);
105: /*
106: Print out each vector, including the ghost padding region.
107: */
108: VecGetArray(lx,&array);
109: for (i=0; i<nlocal+nghost; i++) {
110: PetscSynchronizedPrintf(PETSC_COMM_WORLD,"%" PetscInt_FMT " %g\n",i,(double)PetscRealPart(array[i]));
111: }
112: VecRestoreArray(lx,&array);
113: PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT);
114: VecGhostRestoreLocalForm(gx,&lx);
116: /* Another test that sets ghost values and then accumulates onto the owning processors using MIN_VALUES */
117: if (flg3) {
118: if (rank == 0)PetscSynchronizedPrintf(PETSC_COMM_WORLD,"\nTesting VecGhostUpdate with MIN_VALUES\n");
119: VecGhostGetLocalForm(gx,&lx);
120: VecGetArray(lx,&array);
121: for (i=0; i<nghost; i++) array[nlocal+i] = rank ? (PetscScalar)4 : (PetscScalar)8;
122: VecRestoreArray(lx,&array);
123: VecGhostRestoreLocalForm(gx,&lx);
125: VecGhostUpdateBegin(gx,MIN_VALUES,SCATTER_REVERSE);
126: VecGhostUpdateEnd(gx,MIN_VALUES,SCATTER_REVERSE);
128: VecGhostGetLocalForm(gx,&lx);
129: VecGetArray(lx,&array);
131: for (i=0; i<nlocal+nghost; i++) {
132: PetscSynchronizedPrintf(PETSC_COMM_WORLD,"%" PetscInt_FMT " %g\n",i,(double)PetscRealPart(array[i]));
133: }
134: VecRestoreArray(lx,&array);
135: PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT);
136: VecGhostRestoreLocalForm(gx,&lx);
137: }
139: VecDestroy(&gx);
141: if (flg) PetscFree(tarray);
142: PetscFinalize();
143: return 0;
144: }
146: /*TEST
148: test:
149: nsize: 2
151: test:
152: suffix: 2
153: nsize: 2
154: args: -allocate
155: output_file: output/ex9_1.out
157: test:
158: suffix: 3
159: nsize: 2
160: args: -vecmpisetghost
161: output_file: output/ex9_1.out
163: test:
164: suffix: 4
165: nsize: 2
166: args: -minvalues
167: output_file: output/ex9_2.out
168: requires: !complex
170: TEST*/