Actual source code: ex28.c
2: static char help[] = "Tests repeated VecDotBegin()/VecDotEnd().\n\n";
4: #include <petscvec.h>
5: #define CheckError(a, b, tol) \
6: do { \
7: 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)))); \
8: } while (0)
10: #define CheckErrorScalar(a, b, tol) \
11: do { \
12: if (!PetscIsCloseAtTol(PetscRealPart(a), PetscRealPart(b), 0, tol)) { \
13: 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)))); \
14: } \
15: if (!PetscIsCloseAtTol(PetscImaginaryPart(a), PetscImaginaryPart(b), 0, PETSC_SMALL)) { \
16: 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)))); \
17: } \
18: } while (0)
20: int main(int argc, char **argv)
21: {
22: PetscInt n = 25, i, row0 = 0;
23: PetscScalar two = 2.0, result1, result2, results[40], value, ten = 10.0;
24: PetscScalar result1a, result2a;
25: PetscReal result3, result4, result[2], result3a, result4a, resulta[2];
26: Vec x, y, vecs[40];
27: PetscReal tol = PETSC_SMALL;
29: PetscFunctionBeginUser;
30: PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
32: /* create vectors */
33: PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
34: PetscCall(VecSetSizes(x, n, PETSC_DECIDE));
35: PetscCall(VecSetFromOptions(x));
36: PetscCall(VecDuplicate(x, &y));
37: PetscCall(VecSetRandom(x, NULL));
38: PetscCall(VecViewFromOptions(x, NULL, "-x_view"));
39: PetscCall(VecSet(y, two));
41: /*
42: Test mixing dot products and norms that require sums
43: */
44: result1 = result2 = 0.0;
45: result3 = result4 = 0.0;
46: PetscCall(VecDotBegin(x, y, &result1));
47: PetscCall(VecDotBegin(y, x, &result2));
48: PetscCall(VecNormBegin(y, NORM_2, &result3));
49: PetscCall(VecNormBegin(x, NORM_1, &result4));
50: PetscCall(PetscCommSplitReductionBegin(PetscObjectComm((PetscObject)x)));
51: PetscCall(VecDotEnd(x, y, &result1));
52: PetscCall(VecDotEnd(y, x, &result2));
53: PetscCall(VecNormEnd(y, NORM_2, &result3));
54: PetscCall(VecNormEnd(x, NORM_1, &result4));
56: PetscCall(VecDot(x, y, &result1a));
57: PetscCall(VecDot(y, x, &result2a));
58: PetscCall(VecNorm(y, NORM_2, &result3a));
59: PetscCall(VecNorm(x, NORM_1, &result4a));
61: CheckErrorScalar(result1, result1a, tol);
62: CheckErrorScalar(result2, result2a, tol);
63: CheckError(result3, result3a, tol);
64: CheckError(result4, result4a, tol);
66: /*
67: Test norms that only require abs
68: */
69: result1 = result2 = 0.0;
70: result3 = result4 = 0.0;
71: PetscCall(VecNormBegin(y, NORM_MAX, &result3));
72: PetscCall(VecNormBegin(x, NORM_MAX, &result4));
73: PetscCall(VecNormEnd(y, NORM_MAX, &result3));
74: PetscCall(VecNormEnd(x, NORM_MAX, &result4));
76: PetscCall(VecNorm(x, NORM_MAX, &result4a));
77: PetscCall(VecNorm(y, NORM_MAX, &result3a));
78: CheckError(result3, result3a, tol);
79: CheckError(result4, result4a, tol);
81: /*
82: Tests dot, max, 1, norm
83: */
84: result1 = result2 = 0.0;
85: result3 = result4 = 0.0;
86: PetscCall(VecSetValues(x, 1, &row0, &ten, INSERT_VALUES));
87: PetscCall(VecAssemblyBegin(x));
88: PetscCall(VecAssemblyEnd(x));
90: PetscCall(VecDotBegin(x, y, &result1));
91: PetscCall(VecDotBegin(y, x, &result2));
92: PetscCall(VecNormBegin(x, NORM_MAX, &result3));
93: PetscCall(VecNormBegin(x, NORM_1, &result4));
94: PetscCall(VecDotEnd(x, y, &result1));
95: PetscCall(VecDotEnd(y, x, &result2));
96: PetscCall(VecNormEnd(x, NORM_MAX, &result3));
97: PetscCall(VecNormEnd(x, NORM_1, &result4));
99: PetscCall(VecDot(x, y, &result1a));
100: PetscCall(VecDot(y, x, &result2a));
101: PetscCall(VecNorm(x, NORM_MAX, &result3a));
102: PetscCall(VecNorm(x, NORM_1, &result4a));
104: CheckErrorScalar(result1, result1a, tol);
105: CheckErrorScalar(result2, result2a, tol);
106: CheckError(result3, result3a, tol);
107: CheckError(result4, result4a, tol);
109: /*
110: tests 1_and_2 norm
111: */
112: PetscCall(VecNormBegin(x, NORM_MAX, &result3));
113: PetscCall(VecNormBegin(x, NORM_1_AND_2, result));
114: PetscCall(VecNormBegin(y, NORM_MAX, &result4));
115: PetscCall(VecNormEnd(x, NORM_MAX, &result3));
116: PetscCall(VecNormEnd(x, NORM_1_AND_2, result));
117: PetscCall(VecNormEnd(y, NORM_MAX, &result4));
119: PetscCall(VecNorm(x, NORM_MAX, &result3a));
120: PetscCall(VecNorm(x, NORM_1_AND_2, resulta));
121: PetscCall(VecNorm(y, NORM_MAX, &result4a));
123: CheckError(result3, result3a, tol);
124: CheckError(result4, result4a, tol);
125: CheckError(result[0], resulta[0], tol);
126: CheckError(result[1], resulta[1], tol);
128: PetscCall(VecDestroy(&x));
129: PetscCall(VecDestroy(&y));
131: /*
132: Tests computing a large number of operations that require
133: allocating a larger data structure internally
134: */
135: for (i = 0; i < 40; i++) {
136: PetscCall(VecCreate(PETSC_COMM_WORLD, vecs + i));
137: PetscCall(VecSetSizes(vecs[i], PETSC_DECIDE, n));
138: PetscCall(VecSetFromOptions(vecs[i]));
139: value = (PetscReal)i;
140: PetscCall(VecSet(vecs[i], value));
141: }
142: for (i = 0; i < 39; i++) PetscCall(VecDotBegin(vecs[i], vecs[i + 1], results + i));
143: for (i = 0; i < 39; i++) {
144: PetscScalar expected = 25.0 * i * (i + 1);
145: PetscCall(VecDotEnd(vecs[i], vecs[i + 1], results + i));
146: CheckErrorScalar(results[i], expected, tol);
147: }
148: for (i = 0; i < 40; i++) PetscCall(VecDestroy(&vecs[i]));
150: PetscCall(PetscFinalize());
151: return 0;
152: }
154: /*TEST
156: test:
157: nsize: 3
159: testset:
160: nsize: 3
161: output_file: output/ex28_1.out
163: test:
164: suffix: 2
165: args: -splitreduction_async
167: test:
168: suffix: 2_cuda
169: args: -splitreduction_async -vec_type cuda
170: requires: cuda
172: test:
173: suffix: cuda
174: args: -vec_type cuda
175: requires: cuda
177: test:
178: suffix: 2_hip
179: args: -splitreduction_async -vec_type hip
180: requires: hip
182: test:
183: suffix: hip
184: args: -vec_type hip
185: requires: hip
187: test:
188: suffix: 2_kokkos
189: args: -splitreduction_async -vec_type kokkos
190: requires: kokkos_kernels
192: test:
193: suffix: kokkos
194: args: -vec_type kokkos
195: requires: kokkos_kernels
196: TEST*/