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