|
67 | 67 | #include <universal/number/cfloat/cfloat.hpp> |
68 | 68 | #include <universal/number/bfloat16/bfloat16.hpp> |
69 | 69 | #include <universal/number/elreal/elreal.hpp> |
70 | | -#include <universal/verification/dyadic_exact.hpp> |
| 70 | +#include <universal/verification/elreal_oracle.hpp> |
71 | 71 | #include <universal/verification/test_suite.hpp> |
72 | 72 |
|
73 | 73 | namespace { |
74 | 74 |
|
75 | 75 | using namespace sw::universal; |
| 76 | +// exact_real / exact_block / exact_blocks / exact_value + ZBCL_EXACT_WINDOW: the |
| 77 | +// shared exact-dyadic oracle this file pioneered (#1022), now consolidated (#1035). |
| 78 | +using namespace sw::universal::elreal_oracle; |
76 | 79 |
|
77 | | -// Exact value of a binary floating-point value v as a dyadic rational, with NO |
78 | | -// precision loss at any significand width. |
79 | | -// |
80 | | -// Two regimes, chosen at compile time: |
81 | | -// |
82 | | -// * p = digits <= 53: double represents v exactly, so from_double is exact and |
83 | | -// cheap. Covers float(24), double(53), half(11), bfloat16(8), cfloat<24,5>(19), |
84 | | -// cfloat<32,8>(24) -- every <= 53-bit host. |
85 | | -// |
86 | | -// * p > 53 with the Universal FP API (cfloat quad and up): read the encoding |
87 | | -// DIRECTLY -- sign, scale, and the fbits stored fraction bits -- and assemble |
88 | | -// value = (-1)^sign * (2^fbits + F) * 2^(scale - fbits), fbits = p - 1. |
89 | | -// Reading the encoding directly keeps this oracle self-contained: it depends |
90 | | -// on no cfloat math function, only on the bit layout. (Historically it also |
91 | | -// side-stepped two wide-precision cfloat bugs that an earlier frexp/floor |
92 | | -// extraction tripped over -- cfloat floor mis-handling large integers, #1026, |
93 | | -// and cfloat frexp's non-std [1,2) fraction, #1027; both now fixed. That |
94 | | -// earlier extraction had passed #1023's quad sweeps only because they fed |
95 | | -// double-derived <= 53-bit q128 values, silently mis-extracting genuine |
96 | | -// 113-bit ones.) The bit-based path is verified consistent across widths and |
97 | | -// against an independent 2x-wider cfloat product. |
98 | | -template <typename T> |
99 | | -dyadic exact_real(T v) { |
100 | | - if (v == T(0)) return dyadic(); |
101 | | - if constexpr (std::numeric_limits<T>::digits <= 53) { |
102 | | - return dyadic::from_double(static_cast<double>(v)); |
103 | | - } else if constexpr (has_universal_fp_api_v<T>) { |
104 | | - constexpr int fbits = std::numeric_limits<T>::digits - 1; |
105 | | - assert(v.isnormal() && "exact_real bit path expects a normal value"); |
106 | | - dyadic::bigint F(0); |
107 | | - for (int i = 0; i < fbits; ++i) { |
108 | | - if (v.test(static_cast<unsigned>(i))) { |
109 | | - dyadic::bigint bit(1); bit <<= i; F = F + bit; |
110 | | - } |
111 | | - } |
112 | | - dyadic::bigint M(1); M <<= fbits; M = M + F; // hidden bit + fraction |
113 | | - return dyadic(v.sign() ? -M : M, v.scale() - fbits); |
114 | | - } else { |
115 | | - static_assert(has_universal_fp_api_v<T>, |
116 | | - "exact_real: wide (>53-bit) native hosts are not supported yet; add a " |
117 | | - "std::frexp-based extraction (std frexp is correct for native types) " |
118 | | - "when such a host is introduced."); |
119 | | - return dyadic(); |
120 | | - } |
121 | | -} |
122 | | - |
123 | | -// exact value of a single block as a dyadic rational (value(b) = v * 2^exp). |
124 | | -// Shares no code with the block/EFT/threeAdd/add algorithms under test. |
125 | | -template <typename FpType> |
126 | | -dyadic exact_block(const block<FpType>& b) { |
127 | | - if (b.is_zero_block()) return dyadic(); |
128 | | - dyadic d = exact_real(b.v); |
129 | | - d.scale += b.exp; // multiply by 2^exp exactly (value = v * 2^exp) |
130 | | - return d; |
131 | | -} |
132 | | - |
133 | | -template <typename FpType> |
134 | | -dyadic exact_blocks(const std::vector<block<FpType>>& bs) { |
135 | | - dyadic acc; // 0 |
136 | | - for (const auto& b : bs) acc = acc + exact_block(b); |
137 | | - return acc; |
138 | | -} |
139 | | - |
140 | | -// Number of ZBCL blocks forced when computing an exact value. Finite sums of the |
141 | | -// test inputs settle well within this; kept identical to addition.cpp's window |
142 | | -// so the two files agree on the "exact value" of the same stream. |
143 | | -constexpr std::size_t ZBCL_EXACT_WINDOW = 32; |
144 | | - |
145 | | -template <typename FpType> |
146 | | -dyadic exact_zbcl(const ZBCL<FpType>& z) { |
147 | | - return exact_blocks(z.take(ZBCL_EXACT_WINDOW)); |
148 | | -} |
| 80 | +// The exact-dyadic helpers (exact_real / exact_block / exact_blocks / exact_value) |
| 81 | +// and ZBCL_EXACT_WINDOW that this oracle pioneered (#1022) now live in the shared |
| 82 | +// header <universal/verification/elreal_oracle.hpp> (#1035) and are brought into |
| 83 | +// scope by the using-directive above. The EFT/threeAdd/add reference checks below |
| 84 | +// build on them but share no code with the algorithms under test. |
149 | 85 |
|
150 | 86 | // value of a block in double, and its ulp -- used only for the bfloat16 |
151 | 87 | // truncation bound (double has ~45 bits of headroom over bfloat16's 8). |
@@ -225,7 +161,7 @@ int check_add_exact(double a, double b, const std::string& tag) { |
225 | 161 | auto za = from_native<FpType>(a); |
226 | 162 | auto zb = from_native<FpType>(b); |
227 | 163 | auto z = add(za, zb); |
228 | | - if ((exact_zbcl(za) + exact_zbcl(zb)) != exact_zbcl(z)) { |
| 164 | + if ((exact_value(za) + exact_value(zb)) != exact_value(z)) { |
229 | 165 | std::cout << tag << " add value WRONG: a=" << a << " b=" << b << "\n"; |
230 | 166 | return 1; |
231 | 167 | } |
|
0 commit comments