Actual source code: ex38.c

  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;

 16:   PetscInitialize(&argc,&argv,0,help);
 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);
 30:   PetscObjectTypeCompareAny((PetscObject)X,&iscuda,VECSEQCUDA,VECMPICUDA,"");
 31:   if (iscuda) { /* trigger a copy of the data on the GPU */
 32:     const PetscScalar *xx;

 34:     VecCUDAGetArrayRead(X,&xx);
 35:     VecCUDARestoreArrayRead(X,&xx);
 36:   }

 38:   VecView(X,viewer);

 40:   idxs[0] = (size - rank - 1)*10 + 5;
 41:   idxs[1] = (size - rank - 1)*10 + 2;
 42:   idxs[2] = (size - rank - 1)*10 + 3;

 44:   ISCreateStride(comm,(rend-rstart)/3+3*(rank>size/2),rstart,1,&is0);
 45:   ISComplement(is0,rstart,rend,&is1);
 46:   ISCreateGeneral(comm,3,idxs,PETSC_USE_POINTER,&is2);

 48:   ISView(is0,viewer);
 49:   ISView(is1,viewer);
 50:   ISView(is2,viewer);

 52:   VecGetSubVector(X,is0,&Y);
 53:   VecGetSubVector(X,is1,&Z);
 54:   VecGetSubVector(X,is2,&W);
 55:   VecView(Y,viewer);
 56:   VecView(Z,viewer);
 57:   VecView(W,viewer);
 58:   VecGetArray(Y,&y);
 59:   y[0] = 1000*(rank+1);
 60:   VecRestoreArray(Y,&y);
 61:   VecGetArray(Z,&z);
 62:   z[0] = -1000*(rank+1);
 63:   VecRestoreArray(Z,&z);
 64:   VecGetArray(W,&w);
 65:   w[0] = -10*(rank+1);
 66:   VecRestoreArray(W,&w);
 67:   VecRestoreSubVector(X,is0,&Y);
 68:   VecRestoreSubVector(X,is1,&Z);
 69:   VecRestoreSubVector(X,is2,&W);
 70:   VecView(X,viewer);

 72:   ISDestroy(&is0);
 73:   ISDestroy(&is1);
 74:   ISDestroy(&is2);
 75:   VecDestroy(&X);
 76:   PetscFinalize();
 77:   return 0;
 78: }

 80: /*TEST

 82:    testset:
 83:       nsize: 3
 84:       output_file: output/ex38_1.out
 85:       filter: grep -v "  type:"
 86:       diff_args: -j
 87:       test:
 88:         suffix: standard
 89:         args: -vec_type standard
 90:       test:
 91:         requires: cuda
 92:         suffix: cuda
 93:         args: -vec_type cuda
 94:       test:
 95:         requires: viennacl
 96:         suffix:  viennacl
 97:         args: -vec_type viennacl
 98:       test:
 99:         requires: kokkos_kernels
100:         suffix: kokkos
101:         args: -vec_type kokkos
102:       test:
103:         requires: hip
104:         suffix: hip
105:         args: -vec_type hip

107: TEST*/