Actual source code: ex44.c

petsc-main 2021-04-13
Report Typos and Errors

  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*/