1 static const char help[] = "Tests VecPointwiseMaxAbs()\n\n";
2
3 #include <petscvec.h>
4
TestPointwiseMaxAbs(Vec result,Vec x,Vec y,Vec ref)5 static PetscErrorCode TestPointwiseMaxAbs(Vec result, Vec x, Vec y, Vec ref)
6 {
7 PetscInt n;
8 const PetscScalar *array, *ref_array;
9
10 PetscFunctionBegin;
11 PetscCall(VecPointwiseMaxAbs(result, x, y));
12 PetscCall(VecGetLocalSize(result, &n));
13 PetscCall(VecGetArrayRead(result, &array));
14 PetscCall(VecGetArrayRead(ref, &ref_array));
15 for (PetscInt i = 0; i < n; ++i) {
16 const PetscReal expected = PetscAbsScalar(ref_array[i]);
17 const PetscReal actual_real = PetscRealPart(array[i]), actual_imag = PetscImaginaryPart(array[i]);
18
19 PetscCheck(actual_imag == (PetscReal)0.0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "VecPointwiseMaxAbs() did not properly take absolute value, imaginary part %g != 0.0", (double)actual_imag);
20 PetscCheck(actual_real >= (PetscReal)0.0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "VecPointwiseMaxAbs() did not properly take absolute value, real part %g < 0.0", (double)actual_real);
21 PetscCheck(actual_real == expected, PETSC_COMM_SELF, PETSC_ERR_PLIB, "VecShift() returned array[%" PetscInt_FMT "] %g + %gi != expected_array[%" PetscInt_FMT "] %g + 0.0i", i, (double)actual_real, (double)actual_imag, i, (double)expected);
22 }
23 PetscCall(VecRestoreArrayRead(ref, &ref_array));
24 PetscCall(VecRestoreArrayRead(result, &array));
25 PetscFunctionReturn(PETSC_SUCCESS);
26 }
27
main(int argc,char ** argv)28 int main(int argc, char **argv)
29 {
30 Vec x, y, z;
31
32 PetscFunctionBeginUser;
33 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
34
35 PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
36 PetscCall(VecSetSizes(x, PETSC_DECIDE, 10));
37 PetscCall(VecSetFromOptions(x));
38 PetscCall(VecDuplicate(x, &y));
39 PetscCall(VecDuplicate(x, &z));
40
41 PetscCall(VecSet(x, 0.0));
42 PetscCall(VecSet(y, 10.0));
43
44 // Basic correctness tests, z should always match abs(y) exactly
45 PetscCall(TestPointwiseMaxAbs(z, x, y, y));
46 PetscCall(VecSet(x, 1.0));
47 PetscCall(TestPointwiseMaxAbs(z, x, y, y));
48 PetscCall(VecSet(x, -1.0));
49 PetscCall(TestPointwiseMaxAbs(z, x, y, y));
50 PetscCall(VecSet(y, -10.0));
51 PetscCall(TestPointwiseMaxAbs(z, x, y, y));
52
53 // Test that it works if x and y are the same vector
54 PetscCall(VecSet(x, 0.0));
55 PetscCall(TestPointwiseMaxAbs(z, x, x, x));
56 PetscCall(VecSet(x, 1.0));
57 PetscCall(TestPointwiseMaxAbs(z, x, x, x));
58 PetscCall(VecSet(x, -1.0));
59 PetscCall(TestPointwiseMaxAbs(z, x, x, x));
60
61 // Test that it works if z is one of x or y
62 PetscCall(VecSet(z, 0.0));
63 PetscCall(VecSet(x, 0.0));
64 PetscCall(TestPointwiseMaxAbs(z, x, z, x));
65 PetscCall(VecSet(x, 1.0));
66 PetscCall(TestPointwiseMaxAbs(z, z, x, x));
67 PetscCall(VecSet(x, -10.0));
68 PetscCall(TestPointwiseMaxAbs(z, x, z, x));
69
70 // Test that it works if all vectors are the same
71 PetscCall(VecSet(z, 0.0));
72 PetscCall(TestPointwiseMaxAbs(z, z, z, z));
73 PetscCall(VecSet(z, 1.0));
74 PetscCall(TestPointwiseMaxAbs(z, z, z, z));
75 PetscCall(VecSet(z, -1.0));
76 PetscCall(TestPointwiseMaxAbs(z, z, z, z));
77
78 PetscCall(VecDestroy(&x));
79 PetscCall(VecDestroy(&y));
80 PetscCall(VecDestroy(&z));
81 PetscCall(PetscFinalize());
82 return 0;
83 }
84
85 /*TEST
86
87 testset:
88 output_file: output/empty.out
89 nsize: {{1 2}}
90 test:
91 suffix: standard
92 test:
93 requires: defined(PETSC_USE_SHARED_MEMORY)
94 args: -vec_type shared
95 suffix: shared
96 test:
97 requires: viennacl
98 args: -vec_type viennacl
99 suffix: viennacl
100 test:
101 requires: kokkos_kernels
102 args: -vec_type kokkos
103 suffix: kokkos
104 test:
105 requires: cuda
106 args: -vec_type cuda
107 suffix: cuda
108 test:
109 requires: hip
110 args: -vec_type hip
111 suffix: hip
112
113 TEST*/
114