Actual source code: ex38.c
petsc-3.9.0 2018-04-07
1: static const char help[] = "Test VecGetSubVector()\n\n";
3: #include <petscvec.h>
5: int main(int argc, char *argv[])
6: {
7: MPI_Comm comm;
8: Vec X,Y,Z;
9: PetscMPIInt rank,size;
10: PetscInt i,rstart,rend;
11: PetscScalar *x;
12: PetscViewer viewer;
13: IS is0,is1;
16: PetscInitialize(&argc,&argv,0,help);if (ierr) return ierr;
17: comm = PETSC_COMM_WORLD;
18: viewer = PETSC_VIEWER_STDOUT_WORLD;
19: MPI_Comm_size(comm,&size);
20: MPI_Comm_rank(comm,&rank);
22: VecCreate(comm,&X);
23: VecSetSizes(X,10,PETSC_DETERMINE);
24: VecSetFromOptions(X);
25: VecGetOwnershipRange(X,&rstart,&rend);
27: VecGetArray(X,&x);
28: for (i=0; i<rend-rstart; i++) x[i] = rstart+i;
29: VecRestoreArray(X,&x);
31: ISCreateStride(comm,(rend-rstart)/3+3*(rank>size/2),rstart,1,&is0);
32: ISComplement(is0,rstart,rend,&is1);
34: ISView(is0,viewer);
35: ISView(is1,viewer);
37: VecGetSubVector(X,is0,&Y);
38: VecGetSubVector(X,is1,&Z);
39: VecView(Y,viewer);
40: VecView(Z,viewer);
41: VecRestoreSubVector(X,is0,&Y);
42: VecRestoreSubVector(X,is1,&Z);
44: ISDestroy(&is0);
45: ISDestroy(&is1);
46: VecDestroy(&X);
47: PetscFinalize();
48: return ierr;
49: }
52: /*TEST
54: test:
55: nsize: 3
57: TEST*/