xref: /petsc/src/vec/vec/tutorials/ex9.c (revision e08b1d6d0faae6eca507e20c9d3498f81719d047)
1 
2 static char help[] = "Demonstrates use of VecCreateGhost().\n\n";
3 
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 */
11 
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>
19 
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;
27 
28   PetscFunctionBeginUser;
29   PetscCall(PetscInitialize(&argc,&argv,(char*)0,help));
30   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank));
31   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size));
32   PetscCheck(size == 2,PETSC_COMM_WORLD,PETSC_ERR_WRONG_MPI_SIZE,"Must run example with two processors");
33 
34   /*
35      Construct a two dimensional graph connecting nlocal degrees of
36      freedom per processor. From this we will generate the global
37      indices of needed ghost values
38 
39      For simplicity we generate the entire graph on each processor:
40      in real application the graph would stored in parallel, but this
41      example is only to demonstrate the management of ghost padding
42      with VecCreateGhost().
43 
44      In this example we consider the vector as representing
45      degrees of freedom in a one dimensional grid with periodic
46      boundary conditions.
47 
48         ----Processor  1---------  ----Processor 2 --------
49          0    1   2   3   4    5    6    7   8   9   10   11
50                                |----|
51          |-------------------------------------------------|
52 
53   */
54 
55   if (rank == 0) {
56     ifrom[0] = 11; ifrom[1] = 6;
57   } else {
58     ifrom[0] = 0;  ifrom[1] = 5;
59   }
60 
61   /*
62      Create the vector with two slots for ghost points. Note that both
63      the local vector (lx) and the global vector (gx) share the same
64      array for storing vector values.
65   */
66   PetscCall(PetscOptionsHasName(NULL,NULL,"-allocate",&flg));
67   PetscCall(PetscOptionsHasName(NULL,NULL,"-vecmpisetghost",&flg2));
68   PetscCall(PetscOptionsHasName(NULL,NULL,"-minvalues",&flg3));
69   if (flg) {
70     PetscCall(PetscMalloc1(nlocal+nghost,&tarray));
71     PetscCall(VecCreateGhostWithArray(PETSC_COMM_WORLD,nlocal,PETSC_DECIDE,nghost,ifrom,tarray,&gxs));
72   } else if (flg2) {
73     PetscCall(VecCreate(PETSC_COMM_WORLD,&gxs));
74     PetscCall(VecSetType(gxs,VECMPI));
75     PetscCall(VecSetSizes(gxs,nlocal,PETSC_DECIDE));
76     PetscCall(VecMPISetGhost(gxs,nghost,ifrom));
77   } else {
78     PetscCall(VecCreateGhost(PETSC_COMM_WORLD,nlocal,PETSC_DECIDE,nghost,ifrom,&gxs));
79   }
80 
81   /*
82       Test VecDuplicate()
83   */
84   PetscCall(VecDuplicate(gxs,&gx));
85   PetscCall(VecDestroy(&gxs));
86 
87   /*
88      Access the local representation
89   */
90   PetscCall(VecGhostGetLocalForm(gx,&lx));
91 
92   /*
93      Set the values from 0 to 12 into the "global" vector
94   */
95   PetscCall(VecGetOwnershipRange(gx,&rstart,&rend));
96   for (i=rstart; i<rend; i++) {
97     value = (PetscScalar) i;
98     PetscCall(VecSetValues(gx,1,&i,&value,INSERT_VALUES));
99   }
100   PetscCall(VecAssemblyBegin(gx));
101   PetscCall(VecAssemblyEnd(gx));
102 
103   PetscCall(VecGhostUpdateBegin(gx,INSERT_VALUES,SCATTER_FORWARD));
104   PetscCall(VecGhostUpdateEnd(gx,INSERT_VALUES,SCATTER_FORWARD));
105 
106   /*
107      Print out each vector, including the ghost padding region.
108   */
109   PetscCall(VecGetArray(lx,&array));
110   for (i=0; i<nlocal+nghost; i++) {
111     PetscCall(PetscSynchronizedPrintf(PETSC_COMM_WORLD,"%" PetscInt_FMT " %g\n",i,(double)PetscRealPart(array[i])));
112   }
113   PetscCall(VecRestoreArray(lx,&array));
114   PetscCall(PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT));
115   PetscCall(VecGhostRestoreLocalForm(gx,&lx));
116 
117   /* Another test that sets ghost values and then accumulates onto the owning processors using MIN_VALUES */
118   if (flg3) {
119     if (rank == 0)PetscCall(PetscSynchronizedPrintf(PETSC_COMM_WORLD,"\nTesting VecGhostUpdate with MIN_VALUES\n"));
120     PetscCall(VecGhostGetLocalForm(gx,&lx));
121     PetscCall(VecGetArray(lx,&array));
122     for (i=0; i<nghost; i++) array[nlocal+i] = rank ? (PetscScalar)4 : (PetscScalar)8;
123     PetscCall(VecRestoreArray(lx,&array));
124     PetscCall(VecGhostRestoreLocalForm(gx,&lx));
125 
126     PetscCall(VecGhostUpdateBegin(gx,MIN_VALUES,SCATTER_REVERSE));
127     PetscCall(VecGhostUpdateEnd(gx,MIN_VALUES,SCATTER_REVERSE));
128 
129     PetscCall(VecGhostGetLocalForm(gx,&lx));
130     PetscCall(VecGetArray(lx,&array));
131 
132     for (i=0; i<nlocal+nghost; i++) {
133       PetscCall(PetscSynchronizedPrintf(PETSC_COMM_WORLD,"%" PetscInt_FMT " %g\n",i,(double)PetscRealPart(array[i])));
134     }
135     PetscCall(VecRestoreArray(lx,&array));
136     PetscCall(PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT));
137     PetscCall(VecGhostRestoreLocalForm(gx,&lx));
138   }
139 
140   PetscCall(VecDestroy(&gx));
141 
142   if (flg) PetscCall(PetscFree(tarray));
143   PetscCall(PetscFinalize());
144   return 0;
145 }
146 
147 /*TEST
148 
149      test:
150        nsize: 2
151 
152      test:
153        suffix: 2
154        nsize: 2
155        args: -allocate
156        output_file: output/ex9_1.out
157 
158      test:
159        suffix: 3
160        nsize: 2
161        args: -vecmpisetghost
162        output_file: output/ex9_1.out
163 
164      test:
165        suffix: 4
166        nsize: 2
167        args: -minvalues
168        output_file: output/ex9_2.out
169        requires: !complex
170 
171 TEST*/
172