1 static char help[] = "Scatters from a parallel vector to a sequential vector. \n\ 2 uses block index sets\n\n"; 3 4 #include <petscvec.h> 5 6 int main(int argc, char **argv) 7 { 8 PetscInt bs = 1, n = 5, i, low; 9 PetscInt ix0[3] = {5, 7, 9}, iy0[3] = {1, 2, 4}, ix1[3] = {2, 3, 4}, iy1[3] = {0, 1, 3}; 10 PetscMPIInt size, rank; 11 PetscScalar *array; 12 Vec x, y; 13 IS isx, isy; 14 VecScatter ctx; 15 16 PetscFunctionBeginUser; 17 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 18 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 19 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 20 21 PetscCheck(size >= 2, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "Must run more than one processor"); 22 23 PetscCall(PetscOptionsGetInt(NULL, NULL, "-bs", &bs, NULL)); 24 n = bs * n; 25 26 /* Create vector x over shared memory */ 27 PetscCall(VecCreate(PETSC_COMM_WORLD, &x)); 28 PetscCall(VecSetSizes(x, n, PETSC_DECIDE)); 29 PetscCall(VecSetFromOptions(x)); 30 31 PetscCall(VecGetOwnershipRange(x, &low, NULL)); 32 PetscCall(VecGetArray(x, &array)); 33 for (i = 0; i < n; i++) array[i] = (PetscScalar)(i + low); 34 PetscCall(VecRestoreArray(x, &array)); 35 36 /* Create a sequential vector y */ 37 PetscCall(VecCreateSeq(PETSC_COMM_SELF, n, &y)); 38 PetscCall(VecSet(y, 0.0)); 39 40 /* Create two index sets */ 41 if (rank == 0) { 42 PetscCall(ISCreateBlock(PETSC_COMM_SELF, bs, 3, ix0, PETSC_COPY_VALUES, &isx)); 43 PetscCall(ISCreateBlock(PETSC_COMM_SELF, bs, 3, iy0, PETSC_COPY_VALUES, &isy)); 44 } else { 45 PetscCall(ISCreateBlock(PETSC_COMM_SELF, bs, 3, ix1, PETSC_COPY_VALUES, &isx)); 46 PetscCall(ISCreateBlock(PETSC_COMM_SELF, bs, 3, iy1, PETSC_COPY_VALUES, &isy)); 47 } 48 49 if (rank == 10) { 50 PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n[%d] isx:\n", rank)); 51 PetscCall(ISView(isx, PETSC_VIEWER_STDOUT_SELF)); 52 } 53 54 PetscCall(VecScatterCreate(x, isx, y, isy, &ctx)); 55 PetscCall(VecScatterSetFromOptions(ctx)); 56 57 /* Test forward vecscatter */ 58 PetscCall(VecScatterBegin(ctx, x, y, ADD_VALUES, SCATTER_FORWARD)); 59 PetscCall(VecScatterEnd(ctx, x, y, ADD_VALUES, SCATTER_FORWARD)); 60 if (rank == 0) { 61 PetscCall(PetscPrintf(PETSC_COMM_SELF, "[%d] y:\n", rank)); 62 PetscCall(VecView(y, PETSC_VIEWER_STDOUT_SELF)); 63 } 64 65 /* Test reverse vecscatter */ 66 PetscCall(VecScale(y, -1.0)); 67 if (rank) PetscCall(VecScale(y, 1.0 / (size - 1))); 68 69 PetscCall(VecScatterBegin(ctx, y, x, ADD_VALUES, SCATTER_REVERSE)); 70 PetscCall(VecScatterEnd(ctx, y, x, ADD_VALUES, SCATTER_REVERSE)); 71 PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD)); 72 73 /* Free spaces */ 74 PetscCall(VecScatterDestroy(&ctx)); 75 PetscCall(ISDestroy(&isx)); 76 PetscCall(ISDestroy(&isy)); 77 PetscCall(VecDestroy(&x)); 78 PetscCall(VecDestroy(&y)); 79 PetscCall(PetscFinalize()); 80 return 0; 81 } 82 83 /*TEST 84 85 test: 86 nsize: 3 87 88 TEST*/ 89