Actual source code: ex43.c
1: static char help[] = "Tests VecMDot(),VecDot(),VecMTDot(), and VecTDot()\n";
3: #include <petscvec.h>
4: #include <petscdevice.h>
6: static PetscErrorCode TestRepeatedMDot(void)
7: {
8: const PetscStreamType streams[] = {PETSC_STREAM_DEFAULT};
9: const PetscInt counts[] = {2, 8, 9};
10: PetscDeviceContext saved, dctx;
11: PetscInt n = 8 * 1024 * 1024, N;
12: PetscMPIInt rank;
13: PetscBool empty = PETSC_FALSE;
15: PetscFunctionBeginUser;
16: PetscCall(PetscOptionsGetInt(NULL, NULL, "-mdot_size", &n, NULL));
17: PetscCall(PetscOptionsGetBool(NULL, NULL, "-mdot_empty_rank", &empty, NULL));
18: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
19: if (empty && rank) n = 0;
20: PetscCall(PetscDeviceContextGetCurrentContext(&saved));
21: for (size_t s = 0; s < PETSC_STATIC_ARRAY_LENGTH(streams); ++s) {
22: Vec x, y, ys[9];
23: PetscScalar *a;
24: PetscScalar dots[9];
25: const PetscScalar *ar;
27: PetscCall(PetscDeviceContextDuplicate(saved, &dctx));
28: PetscCall(PetscDeviceContextSetStreamType(dctx, streams[s]));
29: PetscCall(PetscDeviceContextSetUp(dctx));
30: PetscCall(PetscDeviceContextSetCurrentContext(dctx));
31: PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
32: PetscCall(VecSetSizes(x, n, PETSC_DECIDE));
33: PetscCall(VecSetFromOptions(x));
34: PetscCall(VecGetSize(x, &N));
35: PetscCall(VecDuplicate(x, &y));
36: PetscCall(VecSetPinnedMemoryMin(y, 0));
37: PetscCall(VecSet(x, 1));
38: PetscCall(VecGetArrayWrite(y, &a));
39: for (PetscInt i = 0; i < n; ++i) a[i] = 1;
40: PetscCall(VecRestoreArrayWrite(y, &a));
41: PetscCall(VecGetArrayReadAndMemType(y, &ar, NULL));
42: PetscCall(VecRestoreArrayReadAndMemType(y, &ar));
43: for (PetscInt j = 0; j < 9; ++j) ys[j] = y;
44: for (size_t c = 0; c < PETSC_STATIC_ARRAY_LENGTH(counts); ++c) {
45: PetscInt nv = counts[c];
47: // Warm allocations and kernels before making the host copy newer.
48: PetscCall(VecMDot(x, nv, ys, dots));
49: for (PetscInt pass = 0; pass < 3; ++pass) {
50: #if PetscDefined(USE_COMPLEX)
51: PetscScalar value = PetscCMPLX(pass + 2, 1);
52: #else
53: PetscScalar value = pass + 2;
54: #endif
55: PetscScalar expected = N * PetscConj(value);
57: PetscCall(VecGetArrayWrite(y, &a));
58: for (PetscInt i = 0; i < n; ++i) a[i] = value;
59: PetscCall(VecRestoreArrayWrite(y, &a));
60: // The first call needs an upload; the second uses the current device copy.
61: for (PetscInt repeat = 0; repeat < 2; ++repeat) {
62: PetscCall(VecMDot(x, nv, ys, dots));
63: for (PetscInt j = 0; j < nv; ++j)
64: PetscCheck(PetscAbsScalar(dots[j] - expected) <= 100 * PETSC_MACHINE_EPSILON * PetscMax(PetscAbsScalar(expected), 1), PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Incorrect repeated-input dot product %" PetscInt_FMT " (nv %" PetscInt_FMT ", stream type %d, repeat %" PetscInt_FMT "): error %g", j, nv, (int)streams[s], repeat, (double)PetscAbsScalar(dots[j] - expected));
65: }
66: }
67: }
68: PetscCall(VecDestroy(&y));
69: PetscCall(VecDestroy(&x));
70: PetscCall(PetscDeviceContextSetCurrentContext(saved));
71: PetscCall(PetscDeviceContextDestroy(&dctx));
72: }
73: PetscFunctionReturn(PETSC_SUCCESS);
74: }
76: int main(int argc, char **argv)
77: {
78: Vec *V, t;
79: PetscInt i, j, reps, n = 15, k = 6;
80: PetscRandom rctx;
81: PetscScalar *val_dot, *val_mdot, *tval_dot, *tval_mdot;
82: PetscBool test_repeated = PETSC_FALSE;
84: PetscFunctionBeginUser;
85: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
86: PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL));
87: PetscCall(PetscOptionsGetInt(NULL, NULL, "-k", &k, NULL));
88: PetscCall(PetscOptionsGetBool(NULL, NULL, "-test_repeated_mdot", &test_repeated, NULL));
89: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test with %" PetscInt_FMT " random vectors of length %" PetscInt_FMT "\n", k, n));
90: PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rctx));
91: PetscCall(PetscRandomSetFromOptions(rctx));
92: #if PetscDefined(USE_COMPLEX)
93: PetscCall(PetscRandomSetInterval(rctx, -1. + 4. * PETSC_i, 1. + 5. * PETSC_i));
94: #else
95: PetscCall(PetscRandomSetInterval(rctx, -1., 1.));
96: #endif
97: PetscCall(VecCreate(PETSC_COMM_WORLD, &t));
98: PetscCall(VecSetSizes(t, n, PETSC_DECIDE));
99: PetscCall(VecSetFromOptions(t));
100: PetscCall(VecDuplicateVecs(t, k, &V));
101: PetscCall(VecSetRandom(t, rctx));
102: PetscCall(VecViewFromOptions(t, NULL, "-t_view"));
103: PetscCall(PetscMalloc1(k, &val_dot));
104: PetscCall(PetscMalloc1(k, &val_mdot));
105: PetscCall(PetscMalloc1(k, &tval_dot));
106: PetscCall(PetscMalloc1(k, &tval_mdot));
107: for (i = 0; i < k; i++) PetscCall(VecSetRandom(V[i], rctx));
108: for (reps = 0; reps < 20; reps++) {
109: for (i = 1; i < k; i++) {
110: PetscCall(VecMDot(t, i, V, val_mdot));
111: PetscCall(VecMTDot(t, i, V, tval_mdot));
112: for (j = 0; j < i; j++) {
113: PetscCall(VecDot(t, V[j], &val_dot[j]));
114: PetscCall(VecTDot(t, V[j], &tval_dot[j]));
115: }
116: /* Check result */
117: for (j = 0; j < i; j++) {
118: if (PetscAbsScalar(val_mdot[j] - val_dot[j]) / PetscAbsScalar(val_dot[j]) > 1e-5) {
119: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "[TEST FAILED] i=%" PetscInt_FMT ", j=%" PetscInt_FMT ", val_mdot[j]=%g, val_dot[j]=%g\n", i, j, (double)PetscAbsScalar(val_mdot[j]), (double)PetscAbsScalar(val_dot[j])));
120: break;
121: }
122: if (PetscAbsScalar(tval_mdot[j] - tval_dot[j]) / PetscAbsScalar(tval_dot[j]) > 1e-5) {
123: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "[TEST FAILED] i=%" PetscInt_FMT ", j=%" PetscInt_FMT ", tval_mdot[j]=%g, tval_dot[j]=%g\n", i, j, (double)PetscAbsScalar(tval_mdot[j]), (double)PetscAbsScalar(tval_dot[j])));
124: break;
125: }
126: }
127: }
128: }
129: if (test_repeated) PetscCall(TestRepeatedMDot());
130: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test completed successfully!\n"));
131: PetscCall(PetscFree(val_dot));
132: PetscCall(PetscFree(val_mdot));
133: PetscCall(PetscFree(tval_dot));
134: PetscCall(PetscFree(tval_mdot));
135: PetscCall(VecDestroyVecs(k, &V));
136: PetscCall(VecDestroy(&t));
137: PetscCall(PetscRandomDestroy(&rctx));
138: PetscCall(PetscFinalize());
139: return 0;
140: }
142: /*TEST
144: test:
146: testset:
147: output_file: output/ex43_1.out
149: test:
150: suffix: cuda
151: args: -vec_type cuda -random_type curand
152: requires: cuda
154: test:
155: suffix: kokkos
156: args: -vec_type kokkos
157: requires: kokkos_kernels
159: test:
160: suffix: hip
161: args: -vec_type hip
162: requires: hip
164: test:
165: suffix: repeated_cuda
166: nsize: {{1 2}}
167: args: -vec_type cuda -test_repeated_mdot
168: requires: cuda
170: test:
171: suffix: repeated_cuda_empty
172: nsize: 2
173: args: -vec_type cuda -test_repeated_mdot -mdot_empty_rank
174: requires: cuda
176: test:
177: suffix: repeated_hip
178: nsize: {{1 2}}
179: args: -vec_type hip -test_repeated_mdot
180: requires: hip
182: test:
183: suffix: repeated_hip_empty
184: nsize: 2
185: args: -vec_type hip -test_repeated_mdot -mdot_empty_rank
186: requires: hip
187: TEST*/