Actual source code: ex8.c
petsc-3.11.4 2019-09-28
2: static char help[] = "Demonstrates using a local ordering to set values into a parallel vector.\n\n";
4: /*T
5: Concepts: vectors^assembling vectors with local ordering;
6: Processors: n
7: T*/
9: /*
10: Include "petscvec.h" so that we can use vectors. Note that this file
11: automatically includes:
12: petscsys.h - base PETSc routines petscis.h - index sets
13: petscviewer.h - viewers
14: */
15: #include <petscvec.h>
17: int main(int argc,char **argv)
18: {
20: PetscMPIInt rank;
21: PetscInt i,ng,*gindices,rstart,rend,M;
22: PetscScalar one = 1.0;
23: Vec x;
25: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
26: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
28: /*
29: Create a parallel vector.
30: - In this case, we specify the size of each processor's local
31: portion, and PETSc computes the global size. Alternatively,
32: PETSc could determine the vector's distribution if we specify
33: just the global size.
34: */
35: VecCreate(PETSC_COMM_WORLD,&x);
36: VecSetSizes(x,rank+1,PETSC_DECIDE);
37: VecSetFromOptions(x);
38: VecSet(x,one);
40: /*
41: Set the local to global ordering for the vector. Each processor
42: generates a list of the global indices for each local index. Note that
43: the local indices are just whatever is convenient for a particular application.
44: In this case we treat the vector as lying on a one dimensional grid and
45: have one ghost point on each end of the blocks owned by each processor.
46: */
48: VecGetSize(x,&M);
49: VecGetOwnershipRange(x,&rstart,&rend);
50: ng = rend - rstart + 2;
51: PetscMalloc1(ng,&gindices);
52: gindices[0] = rstart - 1;
53: for (i=0; i<ng-1; i++) gindices[i+1] = gindices[i] + 1;
54: /* map the first and last point as periodic */
55: if (gindices[0] == -1) gindices[0] = M - 1;
56: if (gindices[ng-1] == M) gindices[ng-1] = 0;
57: {
58: ISLocalToGlobalMapping ltog;
59: ISLocalToGlobalMappingCreate(PETSC_COMM_SELF,1,ng,gindices,PETSC_COPY_VALUES,<og);
60: VecSetLocalToGlobalMapping(x,ltog);
61: ISLocalToGlobalMappingDestroy(<og);
62: }
63: PetscFree(gindices);
65: /*
66: Set the vector elements.
67: - In this case set the values using the local ordering
68: - Each processor can contribute any vector entries,
69: regardless of which processor "owns" them; any nonlocal
70: contributions will be transferred to the appropriate processor
71: during the assembly process.
72: - In this example, the flag ADD_VALUES indicates that all
73: contributions will be added together.
74: */
75: for (i=0; i<ng; i++) {
76: VecSetValuesLocal(x,1,&i,&one,ADD_VALUES);
77: }
79: /*
80: Assemble vector, using the 2-step process:
81: VecAssemblyBegin(), VecAssemblyEnd()
82: Computations can be done while messages are in transition
83: by placing code between these two statements.
84: */
85: VecAssemblyBegin(x);
86: VecAssemblyEnd(x);
88: /*
89: View the vector; then destroy it.
90: */
91: VecView(x,PETSC_VIEWER_STDOUT_WORLD);
92: VecDestroy(&x);
94: PetscFinalize();
95: return ierr;
96: }
98: /*TEST
100: test:
101: nsize: 4
103: TEST*/