1 #ifndef __ULPH_H__
2 #define __ULPH_H__
4 #include <assert.h>
6 struct ulp {
7         double ulp_max;
8         double ulp_min;
9         double ulp_avg;
10         size_t ulp_skipped;
12         long double ulp_maxl;
13         long double ulp_minl;
14         long double ulp_avgl;
15         size_t ulp_skippedl;
16 };
18 struct ulp_complex {
19         struct ulp ulp_real;
20         struct ulp ulp_imag;
21         struct ulp ulp_total;
22 };
24 /*******************************************************************************
25  *                              Real arithmetic
26  ******************************************************************************/
27 #define ULP_INIT(u)                             \
28         assert(u);                              \
29         (u)->ulp_max      = -DBL_MAX;           \
30         (u)->ulp_min      =  DBL_MAX;           \
31         (u)->ulp_avg      =  0.0;               \
32         (u)->ulp_skipped  =  0;                 \
33         (u)->ulp_maxl     = -LDBL_MAX;          \
34         (u)->ulp_minl     =  LDBL_MAX;          \
35         (u)->ulp_avgl     =  0.0;               \
36         (u)->ulp_skippedl =  0;
38 #define ULP_UPDATE(u, n)                        \
39 assert(u);                                      \
40 do {                                            \
41         if ((n) > (u)->ulp_max)                 \
42                   (u)->ulp_max = n;             \
43         if ((n) < (u)->ulp_min)                 \
44                   (u)->ulp_min = n;             \
45         (u)->ulp_avg += n;                      \
46 } while(0)
48 #define ULP_UPDATEL(u, nl)                      \
49 assert(u);                                      \
50 do {                                            \
51         if ((nl) > (u)->ulp_maxl)               \
52                    (u)->ulp_maxl = nl;          \
53         if ((nl) < (u)->ulp_minl)               \
54                    (u)->ulp_minl = nl;          \
55         (u)->ulp_avgl += nl;                    \
56 } while(0)
59 /*******************************************************************************
60  *                              Complex arithmetic
61  ******************************************************************************/
62 #define ULP_COMPLEX_INIT(uc)                                            \
63         ULP_INIT(&uc->ulp_real);                                        \
64         ULP_INIT(&uc->ulp_imag);                                        \
65         ULP_INIT(&uc->ulp_total);
67 #define ULP_COMPLEX_UPDATE(uc, n)                                       \
68 assert(uc);                                                             \
69 ULP_UPDATE(&uc->ulp_real, creal(n));                                    \
70 ULP_UPDATE(&uc->ulp_imag, cimag(n));                                    \
71 ULP_UPDATE(&uc->ulp_total, fabs(creal(n)) + fabs(cimag(n)));
73 #define ULP_COMPLEX_UPDATEL(uc, n)                                      \
74 assert(uc);                                                             \
75 ULP_UPDATEL(&uc->ulp_real, creall(n));                                  \
76 ULP_UPDATEL(&uc->ulp_imag, cimagl(n));                                  \
77 ULP_UPDATEL(&uc->ulp_total, fabsl(creall(n)) + fabsl(cimagl(n)));
79 double calculp_double(double computed, double exact);
80 long double calculp_long_double(long double computedl, long double exactl);
82 void getfunctionulp_real(const struct fentry *f, struct ulp *u);
83 void getfunctionulp_complex(const struct fentry *f, struct ulp_complex *u);
85 void printulps_double(struct ulp u);
86 void printulps_long_double(struct ulp u);
87 void printulps_double_complex(struct ulp_complex u);
88 void printulps_long_double_complex(struct ulp_complex u);
90 #define NITERATIONS (1 * 1000)
92 #endif  /* ! __ULP_H__ */