Actual source code: ex11.c


  2: static char help[] = "Scatters between parallel vectors. \n\
  3: uses block index sets\n\n";

  5: #include <petscvec.h>

  7: int main(int argc, char **argv)
  8: {
  9:   PetscInt       bs = 1, n = 5, N, i, low;
 10:   PetscInt       ix0[3] = {5, 7, 9}, iy0[3] = {1, 2, 4}, ix1[3] = {2, 3, 1}, iy1[3] = {0, 3, 9};
 11:   PetscMPIInt    size, rank;
 12:   PetscScalar   *array;
 13:   Vec            x, x1, y;
 14:   IS             isx, isy;
 15:   VecScatter     ctx;
 16:   VecScatterType type;
 17:   PetscBool      flg;

 20:   PetscInitialize(&argc, &argv, (char *)0, help);
 21:   MPI_Comm_size(PETSC_COMM_WORLD, &size);
 22:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);


 26:   PetscOptionsGetInt(NULL, NULL, "-bs", &bs, NULL);
 27:   n = bs * n;

 29:   /* Create vector x over shared memory */
 30:   VecCreate(PETSC_COMM_WORLD, &x);
 31:   VecSetSizes(x, n, PETSC_DECIDE);
 32:   VecSetFromOptions(x);

 34:   VecGetOwnershipRange(x, &low, NULL);
 35:   VecGetArray(x, &array);
 36:   for (i = 0; i < n; i++) array[i] = (PetscScalar)(i + low);
 37:   VecRestoreArray(x, &array);
 38:   /* VecView(x,PETSC_VIEWER_STDOUT_WORLD); */

 40:   /* Test some vector functions */
 41:   VecAssemblyBegin(x);
 42:   VecAssemblyEnd(x);

 44:   VecGetSize(x, &N);
 45:   VecGetLocalSize(x, &n);

 47:   VecDuplicate(x, &x1);
 48:   VecCopy(x, x1);
 49:   VecEqual(x, x1, &flg);

 52:   VecScale(x1, 2.0);
 53:   VecSet(x1, 10.0);
 54:   /* VecView(x1,PETSC_VIEWER_STDOUT_WORLD); */

 56:   /* Create vector y over shared memory */
 57:   VecCreate(PETSC_COMM_WORLD, &y);
 58:   VecSetSizes(y, n, PETSC_DECIDE);
 59:   VecSetFromOptions(y);
 60:   VecGetArray(y, &array);
 61:   for (i = 0; i < n; i++) array[i] = -(PetscScalar)(i + 100 * rank);
 62:   VecRestoreArray(y, &array);
 63:   VecAssemblyBegin(y);
 64:   VecAssemblyEnd(y);
 65:   /* VecView(y,PETSC_VIEWER_STDOUT_WORLD); */

 67:   /* Create two index sets */
 68:   if (rank == 0) {
 69:     ISCreateBlock(PETSC_COMM_SELF, bs, 3, ix0, PETSC_COPY_VALUES, &isx);
 70:     ISCreateBlock(PETSC_COMM_SELF, bs, 3, iy0, PETSC_COPY_VALUES, &isy);
 71:   } else {
 72:     ISCreateBlock(PETSC_COMM_SELF, bs, 3, ix1, PETSC_COPY_VALUES, &isx);
 73:     ISCreateBlock(PETSC_COMM_SELF, bs, 3, iy1, PETSC_COPY_VALUES, &isy);
 74:   }

 76:   if (rank == 10) {
 77:     PetscPrintf(PETSC_COMM_SELF, "\n[%d] isx:\n", rank);
 78:     ISView(isx, PETSC_VIEWER_STDOUT_SELF);
 79:     PetscPrintf(PETSC_COMM_SELF, "\n[%d] isy:\n", rank);
 80:     ISView(isy, PETSC_VIEWER_STDOUT_SELF);
 81:   }

 83:   /* Create Vector scatter */
 84:   VecScatterCreate(x, isx, y, isy, &ctx);
 85:   VecScatterSetFromOptions(ctx);
 86:   VecScatterGetType(ctx, &type);
 87:   PetscPrintf(PETSC_COMM_WORLD, "scatter type %s\n", type);

 89:   /* Test forward vecscatter */
 90:   VecSet(y, 0.0);
 91:   VecScatterBegin(ctx, x, y, ADD_VALUES, SCATTER_FORWARD);
 92:   VecScatterEnd(ctx, x, y, ADD_VALUES, SCATTER_FORWARD);
 93:   PetscPrintf(PETSC_COMM_WORLD, "\nSCATTER_FORWARD y:\n");
 94:   VecView(y, PETSC_VIEWER_STDOUT_WORLD);

 96:   /* Test reverse vecscatter */
 97:   VecSet(x, 0.0);
 98:   VecSet(y, 0.0);
 99:   VecGetOwnershipRange(y, &low, NULL);
100:   VecGetArray(y, &array);
101:   for (i = 0; i < n; i++) array[i] = (PetscScalar)(i + low);
102:   VecRestoreArray(y, &array);
103:   VecScatterBegin(ctx, y, x, ADD_VALUES, SCATTER_REVERSE);
104:   VecScatterEnd(ctx, y, x, ADD_VALUES, SCATTER_REVERSE);
105:   PetscPrintf(PETSC_COMM_WORLD, "\nSCATTER_REVERSE x:\n");
106:   VecView(x, PETSC_VIEWER_STDOUT_WORLD);

108:   /* Free objects */
109:   VecScatterDestroy(&ctx);
110:   ISDestroy(&isx);
111:   ISDestroy(&isy);
112:   VecDestroy(&x);
113:   VecDestroy(&x1);
114:   VecDestroy(&y);
115:   PetscFinalize();
116:   return 0;
117: }

119: /*TEST

121:    test:
122:       nsize: 2

124: TEST*/