Actual source code: ex44.c

      2: static char help[] = "Test VecConcatenate both in serial and parallel.\n";

      4: #include <petscvec.h>

      6: int main(int argc,char **args)
      7: {
      8:   Vec                *x, x_test, y, y_test;
      9:   IS                 *x_is;
     10:   VecScatter         y_to_x, x_to_y;
     11:   PetscInt           i, j, nx, shift, x_size, y_size, *idx;
     12:   PetscScalar        *vals;
     13:   PetscBool          flg, x_equal, y_equal;
     14:   PetscErrorCode     ierr;

     16:   PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
     17:   PetscOptionsGetInt(NULL,NULL,"-nx",&nx,&flg);
     18:   if (!flg) nx = 3;

     20:   y_size = 0;
     21:   shift = 0;
     22:   PetscMalloc1(nx, &x);
     23:   for (i=0; i<nx; i++) {
     24:     x_size = 2*(i + 1);
     25:     y_size += x_size;
     26:     VecCreate(PETSC_COMM_WORLD, &x[i]);
     27:     VecSetSizes(x[i], PETSC_DECIDE, x_size);
     28:     VecSetFromOptions(x[i]);
     29:     VecSetUp(x[i]);
     30:     PetscMalloc2(x_size, &idx, x_size, &vals);
     31:     for (j=0; j<x_size; j++) {
     32:       idx[j] = j;
     33:       vals[j] = (PetscScalar)(shift + j + 1);
     34:     }
     35:     shift += x_size;
     36:     VecSetValues(x[i], x_size, (const PetscInt*)idx, (const PetscScalar*)vals, INSERT_VALUES);
     37:     VecAssemblyBegin(x[i]);
     38:     VecAssemblyEnd(x[i]);
     39:     PetscFree2(idx, vals);
     40:     PetscPrintf(PETSC_COMM_WORLD, "Original X[%i] vector\n", i);
     41:     VecView(x[i], PETSC_VIEWER_STDOUT_WORLD);
     42:   }
     43:   VecCreate(PETSC_COMM_WORLD, &y);
     44:   VecSetSizes(y, PETSC_DECIDE, y_size);
     45:   VecSetFromOptions(y);
     46:   VecSetUp(y);
     47:   PetscMalloc2(y_size, &idx, y_size, &vals);
     48:   for (j=0; j<y_size; j++) {
     49:     idx[j] = j;
     50:     vals[j] = (PetscScalar)(j + 1);
     51:   }
     52:   VecSetValues(y, y_size, (const PetscInt*)idx, (const PetscScalar*)vals, INSERT_VALUES);
     53:   VecAssemblyBegin(y);
     54:   VecAssemblyEnd(y);
     55:   PetscFree2(idx, vals);
     56:   PetscPrintf(PETSC_COMM_WORLD, "Expected Y vector\n");
     57:   VecView(y, PETSC_VIEWER_STDOUT_WORLD);
     58:   PetscPrintf(PETSC_COMM_WORLD, "------------------------------------------------------------\n");

     60:   /* ---------- base VecConcatenate() test ----------- */
     61:   VecConcatenate(nx, (const Vec*)x, &y_test, &x_is);
     62:   PetscPrintf(PETSC_COMM_WORLD, "Testing VecConcatenate() for Y = [X[1], X[2], ...]\n");
     63:   VecView(y_test, PETSC_VIEWER_STDOUT_WORLD);
     64:   y_equal = PETSC_FALSE;
     65:   VecEqual(y_test, y, &y_equal);
     66:   if (!y_equal) {
     67:     PetscPrintf(PETSC_COMM_WORLD, "  FAIL\n");
     68:   } else {
     69:     PetscPrintf(PETSC_COMM_WORLD, "  PASS\n");
     70:   }
     71:   PetscPrintf(PETSC_COMM_WORLD, "------------------------------------------------------------\n");

     73:   /* ---------- test VecConcatenate() without IS (checks for dangling memory from IS) ----------- */
     74:   VecDestroy(&y_test);
     75:   VecConcatenate(nx, (const Vec*)x, &y_test, NULL);
     76:   PetscPrintf(PETSC_COMM_WORLD, "Testing VecConcatenate() for Y = [X[1], X[2], ...] w/o IS\n");
     77:   VecView(y_test, PETSC_VIEWER_STDOUT_WORLD);
     78:   y_equal = PETSC_FALSE;
     79:   VecEqual(y_test, y, &y_equal);
     80:   if (!y_equal) {
     81:     PetscPrintf(PETSC_COMM_WORLD, "  FAIL\n");
     82:   } else {
     83:     PetscPrintf(PETSC_COMM_WORLD, "  PASS\n");
     84:   }
     85:   PetscPrintf(PETSC_COMM_WORLD, "------------------------------------------------------------\n");

     87:   /* ---------- using index sets on expected Y instead of concatenated Y ----------- */
     88:   for (i=0; i<nx; i++) {
     89:     VecGetSubVector(y, x_is[i], &x_test);
     90:     PetscPrintf(PETSC_COMM_WORLD, "Testing index set for X[%i] component\n", i);
     91:     VecView(x_test, PETSC_VIEWER_STDOUT_WORLD);
     92:     x_equal = PETSC_FALSE;
     93:     VecEqual(x_test, x[i], &x_equal);
     94:     if (!x_equal) {
     95:       PetscPrintf(PETSC_COMM_WORLD, "  FAIL\n");
     96:     } else {
     97:       PetscPrintf(PETSC_COMM_WORLD, "  PASS\n");
     98:     }
     99:     VecRestoreSubVector(y, x_is[i], &x_test);
    100:   }
    101:   PetscPrintf(PETSC_COMM_WORLD, "------------------------------------------------------------\n");

    103:   /* ---------- using VecScatter to communicate data from Y to X[i] for all i ----------- */
    104:   for (i=0; i<nx; i++) {
    105:     VecDuplicate(x[i], &x_test);
    106:     VecZeroEntries(x_test);
    107:     VecScatterCreate(y, x_is[i], x[i], NULL, &y_to_x);
    108:     VecScatterBegin(y_to_x, y, x_test, INSERT_VALUES, SCATTER_FORWARD);
    109:     VecScatterEnd(y_to_x, y, x_test, INSERT_VALUES, SCATTER_FORWARD);
    110:     VecScatterDestroy(&y_to_x);
    111:     PetscPrintf(PETSC_COMM_WORLD, "Testing VecScatter for Y -> X[%i]\n", i);
    112:     VecView(x_test, PETSC_VIEWER_STDOUT_WORLD);
    113:     x_equal = PETSC_FALSE;
    114:     VecEqual(x_test, x[i], &x_equal);
    115:     if (!x_equal) {
    116:       PetscPrintf(PETSC_COMM_WORLD, "  FAIL\n");
    117:     } else {
    118:       PetscPrintf(PETSC_COMM_WORLD, "  PASS\n");
    119:     }
    120:     VecDestroy(&x_test);
    121:   }
    122:   PetscPrintf(PETSC_COMM_WORLD, "------------------------------------------------------------\n");

    124:   /* ---------- using VecScatter to communicate data from X[i] to Y for all i ----------- */
    125:   VecZeroEntries(y_test);
    126:   for (i=0; i<nx; i++) {
    127:     VecScatterCreate(x[i], NULL, y, x_is[i], &x_to_y);
    128:     VecScatterBegin(x_to_y, x[i], y_test, INSERT_VALUES, SCATTER_FORWARD);
    129:     VecScatterEnd(x_to_y, x[i], y_test, INSERT_VALUES, SCATTER_FORWARD);
    130:     VecScatterDestroy(&x_to_y);
    131:   }
    132:   PetscPrintf(PETSC_COMM_WORLD, "Testing VecScatter for X[:] -> Y\n");
    133:   VecView(y_test, PETSC_VIEWER_STDOUT_WORLD);
    134:   y_equal = PETSC_FALSE;
    135:   VecEqual(y_test, y, &y_equal);
    136:   if (!y_equal) {
    137:     PetscPrintf(PETSC_COMM_WORLD, "  FAIL\n");
    138:   } else {
    139:     PetscPrintf(PETSC_COMM_WORLD, "  PASS\n");
    140:   }
    141:   VecDestroy(&y_test);

    143:   VecDestroy(&y);
    144:   for (i=0; i<nx; i++) {
    145:     VecDestroy(&x[i]);
    146:     ISDestroy(&x_is[i]);
    147:   }
    148:   PetscFree(x);
    149:   PetscFree(x_is);
    150:   PetscFinalize();
    151:   return ierr;
    152: }

    154: /*TEST

    156:     test:
    157:         suffix: serial

    159:     test:
    160:         suffix: parallel
    161:         nsize: 2

    163:     test:
    164:         suffix: cuda
    165:         nsize: 2
    166:         args: -vec_type cuda
    167:         requires: cuda

    169:     test:
    170:         suffix: uneven
    171:         nsize: 5

    173: TEST*/