Actual source code: ex39.c
1: static char help[] = "Tests PetscIsCloseAtTol() routine.\n";
3: #include <petscsys.h>
5: PETSC_INTERN PetscReal zero;
6: PetscReal zero = 0;
7: PETSC_INTERN PetscReal zero2;
8: PetscReal zero2 = 0;
10: #define CALL(call) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "%s -> %s\n", #call, (call) ? "True" : "False"))
12: int main(int argc, char **argv)
13: {
14: PetscReal eps = PETSC_MACHINE_EPSILON;
15: PetscReal neg_zero = PetscRealConstant(-0.0);
16: PetscReal pos_zero = PetscRealConstant(+0.0);
17: PetscReal neg_one = PetscRealConstant(-1.0);
18: PetscReal pos_one = PetscRealConstant(+1.0);
19: PetscReal neg_inf = neg_one / zero; /* -inf */
20: PetscReal pos_inf = pos_one / zero; /* +inf */
21: PetscReal x_nan = zero2 / zero; /* NaN */ /* some compilers may optimize out zero/zero and set x_nan = 1! */
23: PetscFunctionBeginUser;
24: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
26: CALL(PetscIsCloseAtTol(pos_zero, neg_zero, 0, 0));
27: CALL(PetscIsCloseAtTol(pos_one, pos_one, 0, 0));
28: CALL(PetscIsCloseAtTol(pos_one, neg_one, 0, 0));
29: CALL(PetscIsCloseAtTol(pos_one, neg_one, 0, 2));
30: CALL(PetscIsCloseAtTol(pos_one, neg_one, 2, 0));
32: CALL(PetscIsCloseAtTol(pos_one + eps, pos_one, 0, 0));
33: CALL(PetscIsCloseAtTol(pos_one - eps, pos_one, 0, 0));
34: CALL(PetscIsCloseAtTol(pos_one + eps, pos_one, 0, 0));
35: CALL(PetscIsCloseAtTol(pos_one - eps, pos_one, 0, 0));
37: CALL(PetscIsCloseAtTol(pos_one + eps, pos_one, 0, eps));
38: CALL(PetscIsCloseAtTol(pos_one - eps, pos_one, 0, eps));
39: CALL(PetscIsCloseAtTol(pos_one + eps, pos_one, eps, 0));
40: CALL(PetscIsCloseAtTol(pos_one - eps, pos_one, eps, 0));
42: CALL(PetscIsCloseAtTol(pos_one + 2 * eps, pos_one, eps, 0));
43: CALL(PetscIsCloseAtTol(pos_one - 2 * eps, pos_one, eps, 0));
44: CALL(PetscIsCloseAtTol(pos_one + 2 * eps, pos_one, 0, eps));
45: CALL(PetscIsCloseAtTol(pos_one - 2 * eps, pos_one, 0, eps));
47: CALL(PetscIsCloseAtTol(neg_inf, neg_zero, 2, 2));
48: CALL(PetscIsCloseAtTol(neg_inf, pos_zero, 2, 2));
49: CALL(PetscIsCloseAtTol(neg_inf, neg_one, 2, 2));
50: CALL(PetscIsCloseAtTol(neg_inf, pos_one, 2, 2));
51: CALL(PetscIsCloseAtTol(neg_inf, neg_inf, 2, 2));
52: CALL(PetscIsCloseAtTol(neg_inf, pos_inf, 2, 2));
53: CALL(PetscIsCloseAtTol(neg_inf, x_nan, 2, 2));
55: CALL(PetscIsCloseAtTol(pos_inf, neg_zero, 2, 2));
56: CALL(PetscIsCloseAtTol(pos_inf, pos_zero, 2, 2));
57: CALL(PetscIsCloseAtTol(pos_inf, neg_one, 2, 2));
58: CALL(PetscIsCloseAtTol(pos_inf, pos_one, 2, 2));
59: CALL(PetscIsCloseAtTol(pos_inf, neg_inf, 2, 2));
60: CALL(PetscIsCloseAtTol(pos_inf, pos_inf, 2, 2));
61: CALL(PetscIsCloseAtTol(pos_inf, x_nan, 2, 2));
63: CALL(PetscIsCloseAtTol(x_nan, neg_zero, 2, 2));
64: CALL(PetscIsCloseAtTol(x_nan, pos_zero, 2, 2));
65: CALL(PetscIsCloseAtTol(x_nan, neg_one, 2, 2));
66: CALL(PetscIsCloseAtTol(x_nan, pos_one, 2, 2));
67: CALL(PetscIsCloseAtTol(x_nan, neg_inf, 2, 2));
68: CALL(PetscIsCloseAtTol(x_nan, pos_inf, 2, 2));
69: CALL(PetscIsCloseAtTol(x_nan, x_nan, 2, 2));
71: PetscCall(PetscFinalize());
72: return 0;
73: }
75: /*TEST
77: test:
78: args: -fp_trap 0
79: output_file: output/ex39.out
81: TEST*/