Actual source code: ex53.c
1: static const char help[] = "Tests VecShift()\n\n";
3: #include <petscvec.h>
5: static PetscBool PetscIsCloseAtTolScalar(PetscScalar l, PetscScalar r, PetscReal atol, PetscReal rtol)
6: {
7: return (PetscBool)(PetscIsCloseAtTol(PetscRealPart(l), PetscRealPart(r), atol, rtol) && PetscIsCloseAtTol(PetscImaginaryPart(l), PetscImaginaryPart(r), atol, rtol));
8: }
10: static PetscErrorCode CheckVecShift(Vec v, PetscInt n, PetscScalar *array_copy, PetscScalar shift)
11: {
12: const PetscScalar *array;
14: PetscFunctionBegin;
15: for (PetscInt i = 0; i < n; ++i) array_copy[i] += shift;
16: PetscCall(VecShift(v, shift));
17: PetscCall(VecGetArrayRead(v, &array));
18: for (PetscInt i = 0; i < n; ++i) {
19: const PetscScalar actual = array[i], expected = array_copy[i];
21: PetscCheck(PetscIsCloseAtTolScalar(actual, expected, 1e-12, 0.0), PETSC_COMM_SELF, PETSC_ERR_PLIB, "VecShift() returned array[%" PetscInt_FMT "] %g + %gi != expected_array[%" PetscInt_FMT "] %g + %gi", i, (double)PetscRealPart(actual), (double)PetscImaginaryPart(actual), i, (double)PetscRealPart(expected), (double)PetscImaginaryPart(expected));
22: }
23: PetscCall(VecRestoreArrayRead(v, &array));
24: PetscFunctionReturn(PETSC_SUCCESS);
25: }
27: int main(int argc, char **argv)
28: {
29: Vec x;
30: PetscInt n;
31: const PetscScalar *array;
32: PetscScalar *array_copy;
33: PetscReal norm_before, norm_after;
34: PetscBool available;
36: PetscFunctionBeginUser;
37: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
39: PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
40: PetscCall(VecSetSizes(x, PETSC_DECIDE, 10));
41: PetscCall(VecSetFromOptions(x));
43: PetscCall(VecZeroEntries(x));
45: // get a copy of the vectors array, anything we do to the vector via VecShift we will also do
46: // to the copy, and hence they should always match
47: PetscCall(VecGetLocalSize(x, &n));
48: PetscCall(PetscMalloc1(n, &array_copy));
49: PetscCall(VecGetArrayRead(x, &array));
50: PetscCall(PetscArraycpy(array_copy, array, n));
51: PetscCall(VecRestoreArrayRead(x, &array));
53: PetscCall(CheckVecShift(x, n, array_copy, 0.0));
54: PetscCall(CheckVecShift(x, n, array_copy, 1.0));
55: PetscCall(CheckVecShift(x, n, array_copy, -1.0));
56: PetscCall(CheckVecShift(x, n, array_copy, 15.0));
58: PetscCall(VecNorm(x, NORM_2, &norm_before));
59: PetscCall(VecNormAvailable(x, NORM_2, &available, &norm_after));
60: PetscCheck(available, PETSC_COMM_SELF, PETSC_ERR_PLIB, "VecNormAvailable() returned FALSE right after calling VecNorm()");
61: // a shift of zero should not invalidate norms
62: PetscCall(CheckVecShift(x, n, array_copy, 0.0));
63: PetscCall(VecNormAvailable(x, NORM_2, &available, &norm_after));
64: PetscCheck(available, PETSC_COMM_SELF, PETSC_ERR_PLIB, "VecNormAvailable() returned FALSE after calling VecShift() with a shift of 0.0!");
65: // these can be compared with equality as the number should not change *at all*
66: PetscCheck(norm_before == norm_after, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Norms differ before and after calling VecShift() with shift of 0.0: before %g after %g", (double)norm_before, (double)norm_after);
68: PetscCall(PetscFree(array_copy));
69: PetscCall(VecDestroy(&x));
70: PetscCall(PetscFinalize());
71: return 0;
72: }
74: /*TEST
76: testset:
77: output_file: ./output/empty.out
78: nsize: {{1 2}}
79: test:
80: suffix: standard
81: test:
82: requires: defined(PETSC_USE_SHARED_MEMORY)
83: args: -vec_type shared
84: suffix: shared
85: test:
86: requires: viennacl
87: args: -vec_type viennacl
88: suffix: viennacl
89: test:
90: requires: kokkos_kernels
91: args: -vec_type kokkos
92: suffix: kokkos
93: test:
94: requires: cuda
95: args: -vec_type cuda
96: suffix: cuda
97: test:
98: requires: hip
99: args: -vec_type hip
100: suffix: hip
102: TEST*/