Actual source code: ex8.c
petsc-3.3-p7 2013-05-11
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>
19: int main(int argc,char **argv)
20: {
22: PetscMPIInt rank;
23: PetscInt i,N,ng,*gindices,rstart,rend,M;
24: PetscScalar one = 1.0;
25: Vec x;
27: PetscInitialize(&argc,&argv,(char *)0,help);
28: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
30: /*
31: Create a parallel vector.
32: - In this case, we specify the size of each processor's local
33: portion, and PETSc computes the global size. Alternatively,
34: PETSc could determine the vector's distribution if we specify
35: just the global size.
36: */
37: VecCreate(PETSC_COMM_WORLD,&x);
38: VecSetSizes(x,rank+1,PETSC_DECIDE);
39: VecSetFromOptions(x);
40: VecGetSize(x,&N);
41: VecSet(x,one);
43: /*
44: Set the local to global ordering for the vector. Each processor
45: generates a list of the global indices for each local index. Note that
46: the local indices are just whatever is convenient for a particular application.
47: In this case we treat the vector as lying on a one dimensional grid and
48: have one ghost point on each end of the blocks owned by each processor.
49: */
51: VecGetSize(x,&M);
52: VecGetOwnershipRange(x,&rstart,&rend);
53: ng = rend - rstart + 2;
54: PetscMalloc(ng*sizeof(PetscInt),&gindices);
55: gindices[0] = rstart - 1;
56: for (i=0; i<ng-1; i++) {
57: gindices[i+1] = gindices[i] + 1;
58: }
59: /* map the first and last point as periodic */
60: if (gindices[0] == -1) gindices[0] = M - 1;
61: if (gindices[ng-1] == M) gindices[ng-1] = 0;
62: {
63: ISLocalToGlobalMapping ltog;
64: ISLocalToGlobalMappingCreate(PETSC_COMM_SELF,ng,gindices,PETSC_COPY_VALUES,<og);
65: VecSetLocalToGlobalMapping(x,ltog);
66: ISLocalToGlobalMappingDestroy(<og);
67: }
68: PetscFree(gindices);
70: /*
71: Set the vector elements.
72: - In this case set the values using the local ordering
73: - Each processor can contribute any vector entries,
74: regardless of which processor "owns" them; any nonlocal
75: contributions will be transferred to the appropriate processor
76: during the assembly process.
77: - In this example, the flag ADD_VALUES indicates that all
78: contributions will be added together.
79: */
80: for (i=0; i<ng; i++) {
81: VecSetValuesLocal(x,1,&i,&one,ADD_VALUES);
82: }
84: /*
85: Assemble vector, using the 2-step process:
86: VecAssemblyBegin(), VecAssemblyEnd()
87: Computations can be done while messages are in transition
88: by placing code between these two statements.
89: */
90: VecAssemblyBegin(x);
91: VecAssemblyEnd(x);
93: /*
94: View the vector; then destroy it.
95: */
96: VecView(x,PETSC_VIEWER_STDOUT_WORLD);
97: VecDestroy(&x);
99: PetscFinalize();
100: return 0;
101: }
102: