Actual source code: ex28.c

  1: static char help[] = "Tests repeated VecDotBegin()/VecDotEnd().\n\n";

  3: #include <petscvec.h>
  4: #define CheckError(a, b, tol) \
  5:   do { \
  6:     if (!PetscIsCloseAtTol(a, b, 0, tol)) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Real error at line %d, tol %g: %s %g %s %g diff %g\n", __LINE__, (double)tol, #a, (double)(a), #b, (double)(b), (double)((a) - (b)))); \
  7:   } while (0)

  9: #define CheckErrorScalar(a, b, tol) \
 10:   do { \
 11:     if (!PetscIsCloseAtTol(PetscRealPart(a), PetscRealPart(b), 0, tol)) { \
 12:       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Real error at line %d, tol %g: %s %g %s %g diff %g\n", __LINE__, (double)tol, #a, (double)PetscRealPart(a), #b, (double)PetscRealPart(b), (double)PetscRealPart((a) - (b)))); \
 13:     } \
 14:     if (!PetscIsCloseAtTol(PetscImaginaryPart(a), PetscImaginaryPart(b), 0, PETSC_SMALL)) { \
 15:       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Imag error at line %d, tol %g: %s %g %s %g diff %g\n", __LINE__, (double)tol, #a, (double)PetscImaginaryPart(a), #b, (double)PetscImaginaryPart(b), (double)PetscImaginaryPart((a) - (b)))); \
 16:     } \
 17:   } while (0)

 19: int main(int argc, char **argv)
 20: {
 21:   PetscInt    n = 25, i, row0 = 0;
 22:   PetscScalar two = 2.0, result1, result2, results[40], value, ten = 10.0;
 23:   PetscScalar result1a, result2a;
 24:   PetscReal   result3, result4, result[2], result3a, result4a, resulta[2];
 25:   Vec         x, y, vecs[40];
 26:   PetscReal   tol = PETSC_SMALL;

 28:   PetscFunctionBeginUser;
 29:   PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));

 31:   /* create vectors */
 32:   PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
 33:   PetscCall(VecSetSizes(x, n, PETSC_DECIDE));
 34:   PetscCall(VecSetFromOptions(x));
 35:   PetscCall(VecDuplicate(x, &y));
 36:   PetscCall(VecSetRandom(x, NULL));
 37:   PetscCall(VecViewFromOptions(x, NULL, "-x_view"));
 38:   PetscCall(VecSet(y, two));

 40:   /*
 41:         Test mixing dot products and norms that require sums
 42:   */
 43:   result1 = result2 = 0.0;
 44:   result3 = result4 = 0.0;
 45:   PetscCall(VecDotBegin(x, y, &result1));
 46:   PetscCall(VecDotBegin(y, x, &result2));
 47:   PetscCall(VecNormBegin(y, NORM_2, &result3));
 48:   PetscCall(VecNormBegin(x, NORM_1, &result4));
 49:   PetscCall(PetscCommSplitReductionBegin(PetscObjectComm((PetscObject)x)));
 50:   PetscCall(VecDotEnd(x, y, &result1));
 51:   PetscCall(VecDotEnd(y, x, &result2));
 52:   PetscCall(VecNormEnd(y, NORM_2, &result3));
 53:   PetscCall(VecNormEnd(x, NORM_1, &result4));

 55:   PetscCall(VecDot(x, y, &result1a));
 56:   PetscCall(VecDot(y, x, &result2a));
 57:   PetscCall(VecNorm(y, NORM_2, &result3a));
 58:   PetscCall(VecNorm(x, NORM_1, &result4a));

 60:   CheckErrorScalar(result1, result1a, tol);
 61:   CheckErrorScalar(result2, result2a, tol);
 62:   CheckError(result3, result3a, tol);
 63:   CheckError(result4, result4a, tol);

 65:   /*
 66:         Test norms that only require abs
 67:   */
 68:   result1 = result2 = 0.0;
 69:   result3 = result4 = 0.0;
 70:   PetscCall(VecNormBegin(y, NORM_MAX, &result3));
 71:   PetscCall(VecNormBegin(x, NORM_MAX, &result4));
 72:   PetscCall(VecNormEnd(y, NORM_MAX, &result3));
 73:   PetscCall(VecNormEnd(x, NORM_MAX, &result4));

 75:   PetscCall(VecNorm(x, NORM_MAX, &result4a));
 76:   PetscCall(VecNorm(y, NORM_MAX, &result3a));
 77:   CheckError(result3, result3a, tol);
 78:   CheckError(result4, result4a, tol);

 80:   /*
 81:         Tests dot,  max, 1, norm
 82:   */
 83:   result1 = result2 = 0.0;
 84:   result3 = result4 = 0.0;
 85:   PetscCall(VecSetValues(x, 1, &row0, &ten, INSERT_VALUES));
 86:   PetscCall(VecAssemblyBegin(x));
 87:   PetscCall(VecAssemblyEnd(x));

 89:   PetscCall(VecDotBegin(x, y, &result1));
 90:   PetscCall(VecDotBegin(y, x, &result2));
 91:   PetscCall(VecNormBegin(x, NORM_MAX, &result3));
 92:   PetscCall(VecNormBegin(x, NORM_1, &result4));
 93:   PetscCall(VecDotEnd(x, y, &result1));
 94:   PetscCall(VecDotEnd(y, x, &result2));
 95:   PetscCall(VecNormEnd(x, NORM_MAX, &result3));
 96:   PetscCall(VecNormEnd(x, NORM_1, &result4));

 98:   PetscCall(VecDot(x, y, &result1a));
 99:   PetscCall(VecDot(y, x, &result2a));
100:   PetscCall(VecNorm(x, NORM_MAX, &result3a));
101:   PetscCall(VecNorm(x, NORM_1, &result4a));

103:   CheckErrorScalar(result1, result1a, tol);
104:   CheckErrorScalar(result2, result2a, tol);
105:   CheckError(result3, result3a, tol);
106:   CheckError(result4, result4a, tol);

108:   /*
109:        tests 1_and_2 norm
110:   */
111:   PetscCall(VecNormBegin(x, NORM_MAX, &result3));
112:   PetscCall(VecNormBegin(x, NORM_1_AND_2, result));
113:   PetscCall(VecNormBegin(y, NORM_MAX, &result4));
114:   PetscCall(VecNormEnd(x, NORM_MAX, &result3));
115:   PetscCall(VecNormEnd(x, NORM_1_AND_2, result));
116:   PetscCall(VecNormEnd(y, NORM_MAX, &result4));

118:   PetscCall(VecNorm(x, NORM_MAX, &result3a));
119:   PetscCall(VecNorm(x, NORM_1_AND_2, resulta));
120:   PetscCall(VecNorm(y, NORM_MAX, &result4a));

122:   CheckError(result3, result3a, tol);
123:   CheckError(result4, result4a, tol);
124:   CheckError(result[0], resulta[0], tol);
125:   CheckError(result[1], resulta[1], tol);

127:   PetscCall(VecDestroy(&x));
128:   PetscCall(VecDestroy(&y));

130:   /*
131:        Tests computing a large number of operations that require
132:     allocating a larger data structure internally
133:   */
134:   for (i = 0; i < 40; i++) {
135:     PetscCall(VecCreate(PETSC_COMM_WORLD, vecs + i));
136:     PetscCall(VecSetSizes(vecs[i], PETSC_DECIDE, n));
137:     PetscCall(VecSetFromOptions(vecs[i]));
138:     value = (PetscReal)i;
139:     PetscCall(VecSet(vecs[i], value));
140:   }
141:   for (i = 0; i < 39; i++) PetscCall(VecDotBegin(vecs[i], vecs[i + 1], results + i));
142:   for (i = 0; i < 39; i++) {
143:     PetscScalar expected = 25.0 * i * (i + 1);
144:     PetscCall(VecDotEnd(vecs[i], vecs[i + 1], results + i));
145:     CheckErrorScalar(results[i], expected, tol);
146:   }
147:   for (i = 0; i < 40; i++) PetscCall(VecDestroy(&vecs[i]));

149:   PetscCall(PetscFinalize());
150:   return 0;
151: }

153: /*TEST

155:    test:
156:       nsize: 3

158:    testset:
159:      nsize: 3
160:      output_file: output/ex28_1.out

162:      test:
163:         suffix: 2
164:         args: -splitreduction_async -options_left no

166:      test:
167:         suffix: 2_cuda
168:         args: -splitreduction_async -vec_type cuda
169:         requires: cuda

171:      test:
172:         suffix: cuda
173:         args: -vec_type cuda
174:         requires: cuda

176:      test:
177:         suffix: 2_hip
178:         args: -splitreduction_async -vec_type hip
179:         requires: hip

181:      test:
182:         suffix: hip
183:         args: -vec_type hip
184:         requires: hip

186:      test:
187:         suffix: 2_kokkos
188:         args: -splitreduction_async -vec_type kokkos -options_left no
189:         requires: kokkos_kernels

191:      test:
192:         suffix: kokkos
193:         args: -vec_type kokkos
194:         requires: kokkos_kernels
195: TEST*/