1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
|
/*
* fp-test-log2.c - test QEMU's softfloat log2
*
* Copyright (C) 2020, Linaro, Ltd.
*
* License: GNU GPL, version 2 or later.
* See the COPYING file in the top-level directory.
*/
#ifndef HW_POISON_H
#error Must define HW_POISON_H to work around TARGET_* poisoning
#endif
#include "qemu/osdep.h"
#include "qemu/cutils.h"
#include <math.h>
#include "fpu/softfloat.h"
typedef union {
double d;
float64 i;
} ufloat64;
static int errors;
static void compare(ufloat64 test, ufloat64 real, ufloat64 soft, bool exact)
{
int msb;
uint64_t ulp = UINT64_MAX;
if (real.i == soft.i) {
return;
}
msb = 63 - __builtin_clzll(real.i ^ soft.i);
if (msb < 52) {
if (real.i > soft.i) {
ulp = real.i - soft.i;
} else {
ulp = soft.i - real.i;
}
}
/* glibc allows 3 ulp error in its libm-test-ulps; allow 4 here */
if (!exact && ulp <= 4) {
return;
}
printf("test: %016" PRIx64 " %+.13a\n"
" sf: %016" PRIx64 " %+.13a\n"
"libm: %016" PRIx64 " %+.13a\n",
test.i, test.d, soft.i, soft.d, real.i, real.d);
if (msb == 63) {
printf("Error in sign!\n\n");
} else if (msb >= 52) {
printf("Error in exponent: %d\n\n",
(int)(soft.i >> 52) - (int)(real.i >> 52));
} else {
printf("Error in fraction: %" PRIu64 " ulp\n\n", ulp);
}
if (++errors == 20) {
exit(1);
}
}
int main(int ac, char **av)
{
ufloat64 test, real, soft;
float_status qsf = {0};
int i;
set_float_rounding_mode(float_round_nearest_even, &qsf);
test.d = 0.0;
real.d = -__builtin_inf();
soft.i = float64_log2(test.i, &qsf);
compare(test, real, soft, true);
test.d = 1.0;
real.d = 0.0;
soft.i = float64_log2(test.i, &qsf);
compare(test, real, soft, true);
test.d = 2.0;
real.d = 1.0;
soft.i = float64_log2(test.i, &qsf);
compare(test, real, soft, true);
test.d = 4.0;
real.d = 2.0;
soft.i = float64_log2(test.i, &qsf);
compare(test, real, soft, true);
test.d = 0x1p64;
real.d = 64.0;
soft.i = float64_log2(test.i, &qsf);
compare(test, real, soft, true);
test.d = __builtin_inf();
real.d = __builtin_inf();
soft.i = float64_log2(test.i, &qsf);
compare(test, real, soft, true);
for (i = 0; i < 10000; ++i) {
test.d = drand48() + 1.0; /* [1.0, 2.0) */
real.d = log2(test.d);
soft.i = float64_log2(test.i, &qsf);
compare(test, real, soft, false);
test.d = drand48() * 100; /* [0.0, 100) */
real.d = log2(test.d);
soft.i = float64_log2(test.i, &qsf);
compare(test, real, soft, false);
}
return 0;
}
|
