Actual source code: ex31.c
1: static const char help[] = "Tests VecMaxPointwiseDivide()\n\n";
3: #include <petscvec.h>
5: int main(int argc, char **argv)
6: {
7: Vec x, y;
8: PetscScalar *x_array;
9: PetscInt n;
10: PetscReal max, expected;
12: PetscFunctionBeginUser;
13: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
15: PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
16: PetscCall(VecSetSizes(x, PETSC_DECIDE, 10));
17: PetscCall(VecSetFromOptions(x));
19: PetscCall(VecGetLocalSize(x, &n));
20: PetscCall(VecGetArrayWrite(x, &x_array));
21: for (PetscInt i = 0; i < n; ++i) x_array[i] = (PetscScalar)(i + 1);
22: PetscCall(VecRestoreArrayWrite(x, &x_array));
23: expected = (PetscReal)n;
25: PetscCall(VecDuplicate(x, &y));
27: // check that it works at all
28: PetscCall(VecSet(y, 1.0));
29: PetscCall(VecMaxPointwiseDivide(x, y, &max));
30: PetscCheck(PetscIsCloseAtTol(max, expected, 1e-12, 0.0), PETSC_COMM_SELF, PETSC_ERR_PLIB, "VecMaxPointwiseDivide() returned %g != expected %g for y = 1.0", (double)max, (double)expected);
32: // check that it takes the absolute value
33: PetscCall(VecSet(y, -1.0));
34: PetscCall(VecMaxPointwiseDivide(x, y, &max));
35: PetscCheck(PetscIsCloseAtTol(max, expected, 1e-12, 0.0), PETSC_COMM_SELF, PETSC_ERR_PLIB, "VecMaxPointwiseDivide() returned %g != expected %g for y = -1.0", (double)max, (double)expected);
37: // check that it ignores zero entries in y (treats them as 1.0)
38: PetscCall(VecZeroEntries(y));
39: PetscCall(VecMaxPointwiseDivide(x, y, &max));
40: PetscCheck(PetscIsCloseAtTol(max, expected, 1e-12, 0.0), PETSC_COMM_SELF, PETSC_ERR_PLIB, "VecMaxPointwiseDivide() returned %g != expected %g for y = 0.0", (double)max, (double)expected);
42: PetscCall(VecDestroy(&y));
43: PetscCall(VecDestroy(&x));
44: PetscCall(PetscFinalize());
45: return 0;
46: }
48: /*TEST
50: testset:
51: output_file: ./output/empty.out
52: nsize: {{1 2}}
53: test:
54: suffix: standard
55: test:
56: requires: defined(PETSC_USE_SHARED_MEMORY)
57: args: -vec_type shared
58: suffix: shared
59: test:
60: requires: viennacl
61: args: -vec_type viennacl
62: suffix: viennacl
63: test:
64: requires: kokkos_kernels
65: args: -vec_type kokkos
66: suffix: kokkos
67: test:
68: requires: cuda
69: args: -vec_type cuda
70: suffix: cuda
71: test:
72: requires: hip
73: args: -vec_type hip
74: suffix: hip
76: TEST*/