Actual source code: ex38.c
petsc-3.14.4 2021-02-03
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,*y,*w,*z;
12: PetscViewer viewer;
13: IS is0,is1,is2;
14: PetscBool iscuda;
17: PetscInitialize(&argc,&argv,0,help);if (ierr) return ierr;
18: comm = PETSC_COMM_WORLD;
19: viewer = PETSC_VIEWER_STDOUT_WORLD;
20: MPI_Comm_size(comm,&size);
21: MPI_Comm_rank(comm,&rank);
23: VecCreate(comm,&X);
24: VecSetSizes(X,10,PETSC_DETERMINE);
25: VecSetFromOptions(X);
26: VecGetOwnershipRange(X,&rstart,&rend);
28: VecGetArray(X,&x);
29: for (i=0; i<rend-rstart; i++) x[i] = rstart+i;
30: VecRestoreArray(X,&x);
31: PetscObjectTypeCompareAny((PetscObject)X,&iscuda,VECSEQCUDA,VECMPICUDA,"");
32: if (iscuda) { /* trigger a copy of the data on the GPU */
33: const PetscScalar *xx;
35: VecCUDAGetArrayRead(X,&xx);
36: VecCUDARestoreArrayRead(X,&xx);
37: }
39: VecView(X,viewer);
41: idxs[0] = (size - rank - 1)*10 + 5;
42: idxs[1] = (size - rank - 1)*10 + 2;
43: idxs[2] = (size - rank - 1)*10 + 3;
45: ISCreateStride(comm,(rend-rstart)/3+3*(rank>size/2),rstart,1,&is0);
46: ISComplement(is0,rstart,rend,&is1);
47: ISCreateGeneral(comm,3,idxs,PETSC_USE_POINTER,&is2);
49: ISView(is0,viewer);
50: ISView(is1,viewer);
51: ISView(is2,viewer);
53: VecGetSubVector(X,is0,&Y);
54: VecGetSubVector(X,is1,&Z);
55: VecGetSubVector(X,is2,&W);
56: VecView(Y,viewer);
57: VecView(Z,viewer);
58: VecView(W,viewer);
59: VecGetArray(Y,&y);
60: y[0] = 1000*(rank+1);
61: VecRestoreArray(Y,&y);
62: VecGetArray(Z,&z);
63: z[0] = -1000*(rank+1);
64: VecRestoreArray(Z,&z);
65: VecGetArray(W,&w);
66: w[0] = -10*(rank+1);
67: VecRestoreArray(W,&w);
68: VecRestoreSubVector(X,is0,&Y);
69: VecRestoreSubVector(X,is1,&Z);
70: VecRestoreSubVector(X,is2,&W);
71: VecView(X,viewer);
73: ISDestroy(&is0);
74: ISDestroy(&is1);
75: ISDestroy(&is2);
76: VecDestroy(&X);
77: PetscFinalize();
78: return ierr;
79: }
82: /*TEST
84: testset:
85: nsize: 3
86: output_file: output/ex38_1.out
87: filter: grep -v " type:"
88: diff_args: -j
89: test:
90: suffix: standard
91: args: -vec_type standard
92: test:
93: requires: cuda
94: suffix: cuda
95: args: -vec_type cuda
96: test:
97: requires: viennacl
98: suffix: viennacl
99: args: -vec_type viennacl
100: test:
101: requires: kokkos_kernels
102: suffix: kokkos
103: args: -vec_type kokkos
105: TEST*/