Actual source code: ex38.c
petsc-3.12.2 2019-11-22
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,W;
9: PetscMPIInt rank,size;
10: PetscInt i,rstart,rend,idxs[3];
11: PetscScalar *x;
12: PetscViewer viewer;
13: IS is0,is1,is2;
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: idxs[0] = (size - rank - 1)*10 + 5;
32: idxs[1] = (size - rank - 1)*10 + 2;
33: idxs[2] = (size - rank - 1)*10 + 3;
35: ISCreateStride(comm,(rend-rstart)/3+3*(rank>size/2),rstart,1,&is0);
36: ISComplement(is0,rstart,rend,&is1);
37: ISCreateGeneral(comm,3,idxs,PETSC_USE_POINTER,&is2);
39: ISView(is0,viewer);
40: ISView(is1,viewer);
41: ISView(is2,viewer);
43: VecGetSubVector(X,is0,&Y);
44: VecGetSubVector(X,is1,&Z);
45: VecGetSubVector(X,is2,&W);
46: VecView(Y,viewer);
47: VecView(Z,viewer);
48: VecView(W,viewer);
49: VecRestoreSubVector(X,is0,&Y);
50: VecRestoreSubVector(X,is1,&Z);
51: VecRestoreSubVector(X,is2,&W);
53: ISDestroy(&is0);
54: ISDestroy(&is1);
55: ISDestroy(&is2);
56: VecDestroy(&X);
57: PetscFinalize();
58: return ierr;
59: }
62: /*TEST
64: testset:
65: nsize: 3
66: output_file: output/ex38_1.out
67: filter: grep -v " type:"
68: test:
69: suffix: standard
70: args: -vec_type standard
71: test:
72: requires: cuda
73: suffix: cuda
74: args: -vec_type cuda
75: test:
76: requires: viennacl
77: suffix: viennacl
78: args: -vec_type viennacl
80: TEST*/