Actual source code: ex44.c


  2: #include <petscvec.h>

  4: static char help[] = "Tests vecScatter Sequential to Sequential for (CUDA) vectors\n\
  5:   -m # : the size of the vectors\n                                      \
  6:   -n # : the numer of indices (with n<=m)\n                             \
  7:   -toFirst # : the starting index of the output vector for strided scatters\n \
  8:   -toStep # : the step size into the output vector for strided scatters\n \
  9:   -fromFirst # : the starting index of the input vector for strided scatters\n\
 10:   -fromStep # : the step size into the input vector for strided scatters\n\n";

 12: int main(int argc, char *argv[])
 13: {
 14:   Vec         X, Y;
 15:   PetscInt    m, n, i, n1, n2;
 16:   PetscInt    toFirst, toStep, fromFirst, fromStep;
 17:   PetscInt   *idx, *idy;
 18:   PetscBool   flg;
 19:   IS          toISStrided, fromISStrided, toISGeneral, fromISGeneral;
 20:   VecScatter  vscatSStoSS, vscatSStoSG, vscatSGtoSS, vscatSGtoSG;
 21:   ScatterMode mode;
 22:   InsertMode  addv;

 25:   PetscInitialize(&argc, &argv, 0, help);
 26:   PetscOptionsGetInt(NULL, NULL, "-m", &m, &flg);
 27:   if (!flg) m = 100;

 29:   PetscOptionsGetInt(NULL, NULL, "-n", &n, &flg);
 30:   if (!flg) n = 30;

 32:   PetscOptionsGetInt(NULL, NULL, "-toFirst", &toFirst, &flg);
 33:   if (!flg) toFirst = 3;

 35:   PetscOptionsGetInt(NULL, NULL, "-toStep", &toStep, &flg);
 36:   if (!flg) toStep = 3;

 38:   PetscOptionsGetInt(NULL, NULL, "-fromFirst", &fromFirst, &flg);
 39:   if (!flg) fromFirst = 2;

 41:   PetscOptionsGetInt(NULL, NULL, "-fromStep", &fromStep, &flg);
 42:   if (!flg) fromStep = 2;

 44:   if (n > m) {
 45:     PetscPrintf(PETSC_COMM_WORLD, "The vector sizes are %" PetscInt_FMT ". The number of elements being scattered is %" PetscInt_FMT "\n", m, n);
 46:     PetscPrintf(PETSC_COMM_WORLD, "Adjust the parameters such that m>=n\n");
 47:   } else if (toFirst + (n - 1) * toStep >= m) {
 48:     PetscPrintf(PETSC_COMM_WORLD, "The vector sizes are %" PetscInt_FMT ". The number of elements being scattered is %" PetscInt_FMT "\n", m, n);
 49:     PetscPrintf(PETSC_COMM_WORLD, "For the Strided Scatter, toFirst=%" PetscInt_FMT " and toStep=%" PetscInt_FMT ".\n", toFirst, toStep);
 50:     PetscPrintf(PETSC_COMM_WORLD, "This produces an index (toFirst+(n-1)*toStep)>=m\n");
 51:     PetscPrintf(PETSC_COMM_WORLD, "Adjust the parameterrs accordingly with -m, -n, -toFirst, or -toStep\n");
 52:   } else if (fromFirst + (n - 1) * fromStep >= m) {
 53:     PetscPrintf(PETSC_COMM_WORLD, "The vector sizes are %" PetscInt_FMT ". The number of elements being scattered is %" PetscInt_FMT "\n", m, n);
 54:     PetscPrintf(PETSC_COMM_WORLD, "For the Strided Scatter, fromFirst=%" PetscInt_FMT " and fromStep=%" PetscInt_FMT ".\n", fromFirst, toStep);
 55:     PetscPrintf(PETSC_COMM_WORLD, "This produces an index (fromFirst+(n-1)*fromStep)>=m\n");
 56:     PetscPrintf(PETSC_COMM_WORLD, "Adjust the parameterrs accordingly with -m, -n, -fromFirst, or -fromStep\n");
 57:   } else {
 58:     PetscPrintf(PETSC_COMM_WORLD, "m=%" PetscInt_FMT "\tn=%" PetscInt_FMT "\tfromFirst=%" PetscInt_FMT "\tfromStep=%" PetscInt_FMT "\ttoFirst=%" PetscInt_FMT "\ttoStep=%" PetscInt_FMT "\n", m, n, fromFirst, fromStep, toFirst, toStep);
 59:     PetscPrintf(PETSC_COMM_WORLD, "fromFirst+(n-1)*fromStep=%" PetscInt_FMT "\ttoFirst+(n-1)*toStep=%" PetscInt_FMT "\n", fromFirst + (n - 1) * fromStep, toFirst + (n - 1) * toStep);

 61:     /* Build the vectors */
 62:     VecCreate(PETSC_COMM_WORLD, &Y);
 63:     VecSetSizes(Y, m, PETSC_DECIDE);
 64:     VecCreate(PETSC_COMM_WORLD, &X);
 65:     VecSetSizes(X, m, PETSC_DECIDE);

 67:     VecSetFromOptions(Y);
 68:     VecSetFromOptions(X);
 69:     VecSet(X, 2.0);
 70:     VecSet(Y, 1.0);

 72:     /* Build the strided index sets */
 73:     ISCreate(PETSC_COMM_WORLD, &toISStrided);
 74:     ISCreate(PETSC_COMM_WORLD, &fromISStrided);
 75:     ISSetType(toISStrided, ISSTRIDE);
 76:     ISSetType(fromISStrided, ISSTRIDE);
 77:     ISStrideSetStride(fromISStrided, n, fromFirst, fromStep);
 78:     ISStrideSetStride(toISStrided, n, toFirst, toStep);

 80:     /* Build the general index sets */
 81:     PetscMalloc1(n, &idx);
 82:     PetscMalloc1(n, &idy);
 83:     for (i = 0; i < n; i++) {
 84:       idx[i] = i % m;
 85:       idy[i] = (i + m) % m;
 86:     }
 87:     n1 = n;
 88:     n2 = n;
 89:     PetscSortRemoveDupsInt(&n1, idx);
 90:     PetscSortRemoveDupsInt(&n2, idy);

 92:     ISCreateGeneral(PETSC_COMM_WORLD, n1, idx, PETSC_COPY_VALUES, &toISGeneral);
 93:     ISCreateGeneral(PETSC_COMM_WORLD, n2, idy, PETSC_COPY_VALUES, &fromISGeneral);

 95:     /* set the mode and the insert/add parameter */
 96:     mode = SCATTER_FORWARD;
 97:     addv = ADD_VALUES;

 99:     /* VecScatter : Seq Strided to Seq Strided */
100:     VecScatterCreate(X, fromISStrided, Y, toISStrided, &vscatSStoSS);
101:     VecScatterBegin(vscatSStoSS, X, Y, addv, mode);
102:     VecScatterEnd(vscatSStoSS, X, Y, addv, mode);
103:     VecScatterDestroy(&vscatSStoSS);

105:     /* VecScatter : Seq General to Seq Strided */
106:     VecScatterCreate(Y, fromISGeneral, X, toISStrided, &vscatSGtoSS);
107:     VecScatterBegin(vscatSGtoSS, Y, X, addv, mode);
108:     VecScatterEnd(vscatSGtoSS, Y, X, addv, mode);
109:     VecScatterDestroy(&vscatSGtoSS);

111:     /* VecScatter : Seq General to Seq General */
112:     VecScatterCreate(X, fromISGeneral, Y, toISGeneral, &vscatSGtoSG);
113:     VecScatterBegin(vscatSGtoSG, X, Y, addv, mode);
114:     VecScatterEnd(vscatSGtoSG, X, Y, addv, mode);
115:     VecScatterDestroy(&vscatSGtoSG);

117:     /* VecScatter : Seq Strided to Seq General */
118:     VecScatterCreate(Y, fromISStrided, X, toISGeneral, &vscatSStoSG);
119:     VecScatterBegin(vscatSStoSG, Y, X, addv, mode);
120:     VecScatterEnd(vscatSStoSG, Y, X, addv, mode);
121:     VecScatterDestroy(&vscatSStoSG);

123:     /* view the results */
124:     VecView(Y, PETSC_VIEWER_STDOUT_WORLD);

126:     /* Cleanup */
127:     VecDestroy(&X);
128:     VecDestroy(&Y);
129:     ISDestroy(&toISStrided);
130:     ISDestroy(&fromISStrided);
131:     ISDestroy(&toISGeneral);
132:     ISDestroy(&fromISGeneral);
133:     PetscFree(idx);
134:     PetscFree(idy);
135:   }
136:   PetscFinalize();
137:   return 0;
138: }

140: /*TEST

142:    test:
143:       suffix: cuda
144:       args: -vec_type cuda
145:       requires: cuda

147: TEST*/