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:   PetscFunctionBeginUser;
 17:   PetscCall(PetscInitialize(&argc, &argv, 0, help));
 18:   comm   = PETSC_COMM_WORLD;
 19:   viewer = PETSC_VIEWER_STDOUT_WORLD;
 20:   PetscCallMPI(MPI_Comm_size(comm, &size));
 21:   PetscCallMPI(MPI_Comm_rank(comm, &rank));

 23:   PetscCall(VecCreate(comm, &X));
 24:   PetscCall(VecSetSizes(X, 10, PETSC_DETERMINE));
 25:   PetscCall(VecSetFromOptions(X));
 26:   PetscCall(VecGetOwnershipRange(X, &rstart, &rend));

 28:   PetscCall(VecGetArray(X, &x));
 29:   for (i = 0; i < rend - rstart; i++) x[i] = rstart + i;
 30:   PetscCall(VecRestoreArray(X, &x));
 31:   PetscCall(PetscObjectTypeCompareAny((PetscObject)X, &iscuda, VECSEQCUDA, VECMPICUDA, ""));
 32:   if (iscuda) { /* trigger a copy of the data on the GPU */
 33:     const PetscScalar *xx;

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

 39:   PetscCall(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:   PetscCall(ISCreateStride(comm, (rend - rstart) / 3 + 3 * (rank > size / 2), rstart, 1, &is0));
 46:   PetscCall(ISComplement(is0, rstart, rend, &is1));
 47:   PetscCall(ISCreateGeneral(comm, 3, idxs, PETSC_USE_POINTER, &is2));

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

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

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

 81: /*TEST

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

108: TEST*/