1 | /* Test floating-point arithmetic operations. |
2 | Copyright (C) 1997-2022 Free Software Foundation, Inc. |
3 | This file is part of the GNU C Library. |
4 | |
5 | The GNU C Library is free software; you can redistribute it and/or |
6 | modify it under the terms of the GNU Lesser General Public |
7 | License as published by the Free Software Foundation; either |
8 | version 2.1 of the License, or (at your option) any later version. |
9 | |
10 | The GNU C Library is distributed in the hope that it will be useful, |
11 | but WITHOUT ANY WARRANTY; without even the implied warranty of |
12 | MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU |
13 | Lesser General Public License for more details. |
14 | |
15 | You should have received a copy of the GNU Lesser General Public |
16 | License along with the GNU C Library; if not, see |
17 | <https://www.gnu.org/licenses/>. */ |
18 | #ifndef _GNU_SOURCE |
19 | #define _GNU_SOURCE |
20 | #endif |
21 | #include <math.h> |
22 | #include <stdio.h> |
23 | #include <stdlib.h> |
24 | #include <string.h> |
25 | #include <fenv.h> |
26 | #include <assert.h> |
27 | |
28 | #ifndef ESIZE |
29 | typedef double tocheck_t; |
30 | #define ESIZE 11 |
31 | #define MSIZE 52 |
32 | #define FUNC(x) x |
33 | #endif |
34 | |
35 | #define R_NEAREST 1 |
36 | #define R_ZERO 2 |
37 | #define R_UP 4 |
38 | #define R_DOWN 8 |
39 | #define R_ALL (R_NEAREST|R_ZERO|R_UP|R_DOWN) |
40 | static fenv_t rmodes[4]; |
41 | static const char * const rmnames[4] = |
42 | { "nearest" ,"zero" ,"+Inf" ,"-Inf" }; |
43 | |
44 | typedef union { |
45 | tocheck_t tc; |
46 | unsigned char c[sizeof (tocheck_t)]; |
47 | } union_t; |
48 | |
49 | /* Don't try reading these in a font that doesn't distinguish |
50 | O and zero. */ |
51 | typedef enum { |
52 | P_Z = 0x0, /* 00000...0 */ |
53 | P_000O = 0x1, /* 00011...1 */ |
54 | P_001Z = 0x2, /* 00100...0 */ |
55 | P_00O = 0x3, /* 00111...1 */ |
56 | P_01Z = 0x4, /* 01000...0 */ |
57 | P_010O = 0x5, /* 01011...1 */ |
58 | P_011Z = 0x6, /* 01100...0 */ |
59 | P_0O = 0x7, /* 01111...1 */ |
60 | P_1Z = 0x8, /* 10000...0 */ |
61 | P_100O = 0x9, /* 10011...1 */ |
62 | P_101Z = 0xa, /* 10100...0 */ |
63 | P_10O = 0xb, /* 10111...1 */ |
64 | P_11Z = 0xc, /* 11000...0 */ |
65 | P_110O = 0xd, /* 11011...1 */ |
66 | P_111Z = 0xe, /* 11100...0 */ |
67 | P_O = 0xf, /* 11111...1 */ |
68 | P_Z1 = 0x11, /* 000...001 */ |
69 | P_Z10 = 0x12, /* 000...010 */ |
70 | P_Z11 = 0x13, /* 000...011 */ |
71 | P_0O00 = 0x14, /* 011...100 */ |
72 | P_0O01 = 0x15, /* 011...101 */ |
73 | P_0O0 = 0x16, /* 011...110 */ |
74 | P_1Z1 = 0x19, /* 100...001 */ |
75 | P_1Z10 = 0x1a, /* 100...010 */ |
76 | P_1Z11 = 0x1b, /* 100...011 */ |
77 | P_O00 = 0x1c, /* 111...100 */ |
78 | P_O01 = 0x1d, /* 111...101 */ |
79 | P_O0 = 0x1e, /* 111...110 */ |
80 | P_R = 0x20, /* rrr...rrr */ /* ('r' means random. ) */ |
81 | P_Ro = 0x21, /* rrr...rrr, with odd parity. */ |
82 | P_0R = 0x22, /* 0rr...rrr */ |
83 | P_1R = 0x23, /* 1rr...rrr */ |
84 | P_Rno = 0x24, /* rrr...rrr, but not all ones. */ |
85 | } pattern_t; |
86 | |
87 | static void |
88 | pattern_fill(pattern_t ptn, unsigned char *start, int bitoffset, int count) |
89 | { |
90 | #define bitset(count, value) \ |
91 | start[(count)/8] = (start[(count)/8] & ~(1 << 7-(count)%8) \ |
92 | | (value) << 7-(count)%8) |
93 | int i; |
94 | |
95 | if (ptn >= 0 && ptn <= 0xf) |
96 | { |
97 | /* Patterns between 0 and 0xF have the following format: |
98 | The LSBit is used to fill the last n-3 bits of the pattern; |
99 | The next 3 bits are the first 3 bits of the pattern. */ |
100 | for (i = 0; i < count; i++) |
101 | if (i < 3) |
102 | bitset((bitoffset+i), ptn >> (3-i) & 1); |
103 | else |
104 | bitset((bitoffset+i), ptn >> 0 & 1); |
105 | } |
106 | else if (ptn <= 0x1f) |
107 | { |
108 | /* Patterns between 0x10 and 0x1F have the following format: |
109 | The two LSBits are the last two bits of the pattern; |
110 | The 0x8 bit is the first bit of the pattern; |
111 | The 0x4 bit is used to fill the remainder. */ |
112 | for (i = 0; i < count; i++) |
113 | if (i == 0) |
114 | bitset((bitoffset+i), ptn >> 3 & 1); |
115 | else if (i >= count-2) |
116 | bitset((bitoffset+i), ptn >> (count-1-i) & 1); |
117 | else |
118 | bitset((bitoffset+i), ptn >> 2 & 1); |
119 | } |
120 | else switch (ptn) |
121 | { |
122 | case P_0R: case P_1R: |
123 | assert(count > 0); |
124 | bitset(bitoffset, ptn & 1); |
125 | count--; |
126 | bitoffset++; |
127 | case P_R: |
128 | for (; count > 0; count--, bitoffset++) |
129 | bitset(bitoffset, rand() & 1); |
130 | break; |
131 | case P_Ro: |
132 | { |
133 | int op = 1; |
134 | assert(count > 0); |
135 | for (; count > 1; count--, bitoffset++) |
136 | bitset(bitoffset, op ^= (rand() & 1)); |
137 | bitset(bitoffset, op); |
138 | break; |
139 | } |
140 | case P_Rno: |
141 | { |
142 | int op = 1; |
143 | assert(count > 0); |
144 | for (; count > 1; count--, bitoffset++) |
145 | { |
146 | int r = rand() & 1; |
147 | op &= r; |
148 | bitset(bitoffset, r); |
149 | } |
150 | bitset(bitoffset, rand() & (op ^ 1)); |
151 | break; |
152 | } |
153 | |
154 | default: |
155 | assert(0); |
156 | } |
157 | #undef bitset |
158 | } |
159 | |
160 | static tocheck_t |
161 | pattern(int negative, pattern_t exp, pattern_t mant) |
162 | { |
163 | union_t result; |
164 | #if 0 |
165 | int i; |
166 | #endif |
167 | |
168 | pattern_fill(ptn: negative ? P_O : P_Z, start: result.c, bitoffset: 0, count: 1); |
169 | pattern_fill(ptn: exp, start: result.c, bitoffset: 1, ESIZE); |
170 | pattern_fill(ptn: mant, start: result.c, ESIZE+1, MSIZE); |
171 | #if 0 |
172 | printf("neg=%d exp=%02x mant=%02x: " , negative, exp, mant); |
173 | for (i = 0; i < sizeof (tocheck_t); i++) |
174 | printf("%02x" , result.c[i]); |
175 | printf("\n" ); |
176 | #endif |
177 | return result.tc; |
178 | } |
179 | |
180 | /* Return the closest different tocheck_t to 'x' in the direction of |
181 | 'direction', or 'x' if there is no such value. Assumes 'x' is not |
182 | a NaN. */ |
183 | static tocheck_t |
184 | delta(tocheck_t x, int direction) |
185 | { |
186 | union_t xx; |
187 | int i; |
188 | |
189 | xx.tc = x; |
190 | if (xx.c[0] & 0x80) |
191 | direction = -direction; |
192 | if (direction == +1) |
193 | { |
194 | union_t tx; |
195 | tx.tc = pattern(negative: xx.c[0] >> 7, exp: P_O, mant: P_Z); |
196 | if (memcmp (tx.c, xx.c, sizeof (tocheck_t)) == 0) |
197 | return x; |
198 | } |
199 | for (i = sizeof (tocheck_t)-1; i > 0; i--) |
200 | { |
201 | xx.c[i] += direction; |
202 | if (xx.c[i] != (direction > 0 ? 0 : 0xff)) |
203 | return xx.tc; |
204 | } |
205 | if (direction < 0 && (xx.c[0] & 0x7f) == 0) |
206 | return pattern(negative: ~(xx.c[0] >> 7) & 1, exp: P_Z, mant: P_Z1); |
207 | else |
208 | { |
209 | xx.c[0] += direction; |
210 | return xx.tc; |
211 | } |
212 | } |
213 | |
214 | static int nerrors = 0; |
215 | |
216 | #ifdef FE_ALL_INVALID |
217 | static const int all_exceptions = FE_ALL_INVALID | FE_ALL_EXCEPT; |
218 | #else |
219 | static const int all_exceptions = FE_ALL_EXCEPT; |
220 | #endif |
221 | |
222 | static void |
223 | check_result(int line, const char *rm, tocheck_t expected, tocheck_t actual) |
224 | { |
225 | if (memcmp (&expected, &actual, sizeof (tocheck_t)) != 0) |
226 | { |
227 | unsigned char *ex, *ac; |
228 | size_t i; |
229 | |
230 | printf(format: "%s:%d:round %s:result failed\n" |
231 | " expected result 0x" , __FILE__, line, rm); |
232 | ex = (unsigned char *)&expected; |
233 | ac = (unsigned char *)&actual; |
234 | for (i = 0; i < sizeof (tocheck_t); i++) |
235 | printf(format: "%02x" , ex[i]); |
236 | printf(format: " got 0x" ); |
237 | for (i = 0; i < sizeof (tocheck_t); i++) |
238 | printf(format: "%02x" , ac[i]); |
239 | printf(format: "\n" ); |
240 | nerrors++; |
241 | } |
242 | } |
243 | |
244 | static const struct { |
245 | int except; |
246 | const char *name; |
247 | } excepts[] = { |
248 | #define except_entry(ex) { ex, #ex } , |
249 | #ifdef FE_INEXACT |
250 | except_entry(FE_INEXACT) |
251 | #else |
252 | # define FE_INEXACT 0 |
253 | #endif |
254 | #ifdef FE_DIVBYZERO |
255 | except_entry(FE_DIVBYZERO) |
256 | #else |
257 | # define FE_DIVBYZERO 0 |
258 | #endif |
259 | #ifdef FE_UNDERFLOW |
260 | except_entry(FE_UNDERFLOW) |
261 | #else |
262 | # define FE_UNDERFLOW 0 |
263 | #endif |
264 | #ifdef FE_OVERFLOW |
265 | except_entry(FE_OVERFLOW) |
266 | #else |
267 | # define FE_OVERFLOW 0 |
268 | #endif |
269 | #ifdef FE_INVALID |
270 | except_entry(FE_INVALID) |
271 | #else |
272 | # define FE_INVALID 0 |
273 | #endif |
274 | #ifdef FE_INVALID_SNAN |
275 | except_entry(FE_INVALID_SNAN) |
276 | #else |
277 | # define FE_INVALID_SNAN FE_INVALID |
278 | #endif |
279 | #ifdef FE_INVALID_ISI |
280 | except_entry(FE_INVALID_ISI) |
281 | #else |
282 | # define FE_INVALID_ISI FE_INVALID |
283 | #endif |
284 | #ifdef FE_INVALID_IDI |
285 | except_entry(FE_INVALID_IDI) |
286 | #else |
287 | # define FE_INVALID_IDI FE_INVALID |
288 | #endif |
289 | #ifdef FE_INVALID_ZDZ |
290 | except_entry(FE_INVALID_ZDZ) |
291 | #else |
292 | # define FE_INVALID_ZDZ FE_INVALID |
293 | #endif |
294 | #ifdef FE_INVALID_COMPARE |
295 | except_entry(FE_INVALID_COMPARE) |
296 | #else |
297 | # define FE_INVALID_COMPARE FE_INVALID |
298 | #endif |
299 | #ifdef FE_INVALID_SOFTWARE |
300 | except_entry(FE_INVALID_SOFTWARE) |
301 | #else |
302 | # define FE_INVALID_SOFTWARE FE_INVALID |
303 | #endif |
304 | #ifdef FE_INVALID_SQRT |
305 | except_entry(FE_INVALID_SQRT) |
306 | #else |
307 | # define FE_INVALID_SQRT FE_INVALID |
308 | #endif |
309 | #ifdef FE_INVALID_INTEGER_CONVERSION |
310 | except_entry(FE_INVALID_INTEGER_CONVERSION) |
311 | #else |
312 | # define FE_INVALID_INTEGER_CONVERSION FE_INVALID |
313 | #endif |
314 | }; |
315 | |
316 | static int excepts_missing = 0; |
317 | |
318 | static void |
319 | check_excepts(int line, const char *rm, int expected, int actual) |
320 | { |
321 | if (expected & excepts_missing) |
322 | expected = expected & ~excepts_missing | FE_INVALID_SNAN; |
323 | if ((expected & all_exceptions) != actual) |
324 | { |
325 | size_t i; |
326 | printf(format: "%s:%d:round %s:exceptions failed\n" |
327 | " expected exceptions " , __FILE__, line,rm); |
328 | for (i = 0; i < sizeof (excepts) / sizeof (excepts[0]); i++) |
329 | if (expected & excepts[i].except) |
330 | printf(format: "%s " ,excepts[i].name); |
331 | if ((expected & all_exceptions) == 0) |
332 | printf(format: "- " ); |
333 | printf(format: "got" ); |
334 | for (i = 0; i < sizeof (excepts) / sizeof (excepts[0]); i++) |
335 | if (actual & excepts[i].except) |
336 | printf(format: " %s" ,excepts[i].name); |
337 | if ((actual & all_exceptions) == 0) |
338 | printf(format: "- " ); |
339 | printf(format: ".\n" ); |
340 | nerrors++; |
341 | } |
342 | } |
343 | |
344 | typedef enum { |
345 | B_ADD, B_SUB, B_MUL, B_DIV, B_NEG, B_ABS, B_SQRT |
346 | } op_t; |
347 | typedef struct { |
348 | int line; |
349 | op_t op; |
350 | int a_sgn; |
351 | pattern_t a_exp, a_mant; |
352 | int b_sgn; |
353 | pattern_t b_exp, b_mant; |
354 | int rmode; |
355 | int excepts; |
356 | int x_sgn; |
357 | pattern_t x_exp, x_mant; |
358 | } optest_t; |
359 | static const optest_t optests[] = { |
360 | /* Additions of zero. */ |
361 | {__LINE__,B_ADD, 0,P_Z,P_Z, 0,P_Z,P_Z, R_ALL,0, 0,P_Z,P_Z }, |
362 | {__LINE__,B_ADD, 1,P_Z,P_Z, 0,P_Z,P_Z, R_ALL & ~R_DOWN,0, 0,P_Z,P_Z }, |
363 | {__LINE__,B_ADD, 1,P_Z,P_Z, 0,P_Z,P_Z, R_DOWN,0, 1,P_Z,P_Z }, |
364 | {__LINE__,B_ADD, 1,P_Z,P_Z, 1,P_Z,P_Z, R_ALL,0, 1,P_Z,P_Z }, |
365 | |
366 | /* Additions with NaN. */ |
367 | {__LINE__,B_ADD, 0,P_O,P_101Z, 0,P_Z,P_Z, R_ALL,0, 0,P_O,P_101Z }, |
368 | {__LINE__,B_ADD, 0,P_O,P_01Z, 0,P_Z,P_Z, R_ALL, |
369 | FE_INVALID | FE_INVALID_SNAN, 0,P_O,P_11Z }, |
370 | {__LINE__,B_ADD, 0,P_O,P_Z, 0,P_O,P_0O, R_ALL, |
371 | FE_INVALID | FE_INVALID_SNAN, 0,P_O,P_O }, |
372 | {__LINE__,B_ADD, 0,P_Z,P_Z, 0,P_O,P_11Z, R_ALL,0, 0,P_O,P_11Z }, |
373 | {__LINE__,B_ADD, 0,P_O,P_001Z, 0,P_O,P_001Z, R_ALL, |
374 | FE_INVALID | FE_INVALID_SNAN, 0,P_O,P_101Z }, |
375 | {__LINE__,B_ADD, 0,P_O,P_1Z, 0,P_Z,P_Z, R_ALL,0, 0,P_O,P_1Z }, |
376 | {__LINE__,B_ADD, 0,P_0O,P_Z, 0,P_O,P_10O, R_ALL,0, 0,P_O,P_10O }, |
377 | |
378 | /* Additions with infinity. */ |
379 | {__LINE__,B_ADD, 0,P_O,P_Z, 0,P_Z,P_Z, R_ALL,0, 0,P_O,P_Z }, |
380 | {__LINE__,B_ADD, 0,P_O,P_Z, 1,P_Z,P_Z, R_ALL,0, 0,P_O,P_Z }, |
381 | {__LINE__,B_ADD, 1,P_O,P_Z, 0,P_Z,P_Z, R_ALL,0, 1,P_O,P_Z }, |
382 | {__LINE__,B_ADD, 1,P_O,P_Z, 1,P_Z,P_Z, R_ALL,0, 1,P_O,P_Z }, |
383 | {__LINE__,B_ADD, 0,P_O,P_Z, 0,P_O,P_Z, R_ALL,0, 0,P_O,P_Z }, |
384 | {__LINE__,B_ADD, 1,P_O,P_Z, 1,P_O,P_Z, R_ALL,0, 1,P_O,P_Z }, |
385 | {__LINE__,B_ADD, 0,P_O,P_Z, 1,P_O,P_Z, R_ALL, |
386 | FE_INVALID | FE_INVALID_ISI, 0,P_O,P_1Z }, |
387 | {__LINE__,B_ADD, 1,P_O,P_Z, 0,P_O,P_Z, R_ALL, |
388 | FE_INVALID | FE_INVALID_ISI, 0,P_O,P_1Z }, |
389 | {__LINE__,B_ADD, 0,P_O,P_Z, 0,P_0O,P_Z, R_ALL,0, 0,P_O,P_Z }, |
390 | {__LINE__,B_ADD, 1,P_O,P_Z, 0,P_0O,P_Z, R_ALL,0, 1,P_O,P_Z }, |
391 | {__LINE__,B_ADD, 0,P_O,P_Z, 1,P_0O,P_Z, R_ALL,0, 0,P_O,P_Z }, |
392 | {__LINE__,B_ADD, 1,P_O,P_Z, 1,P_0O,P_Z, R_ALL,0, 1,P_O,P_Z }, |
393 | |
394 | /* Overflow (and zero). */ |
395 | {__LINE__,B_ADD, 0,P_O0,P_Z, 0,P_O0,P_Z, R_NEAREST | R_UP, |
396 | FE_INEXACT | FE_OVERFLOW, 0,P_O,P_Z }, |
397 | {__LINE__,B_ADD, 0,P_O0,P_Z, 0,P_O0,P_Z, R_ZERO | R_DOWN, |
398 | FE_INEXACT | FE_OVERFLOW, 0,P_O0,P_O }, |
399 | {__LINE__,B_ADD, 1,P_O0,P_Z, 1,P_O0,P_Z, R_NEAREST | R_DOWN, |
400 | FE_INEXACT | FE_OVERFLOW, 1,P_O,P_Z }, |
401 | {__LINE__,B_ADD, 1,P_O0,P_Z, 1,P_O0,P_Z, R_ZERO | R_UP, |
402 | FE_INEXACT | FE_OVERFLOW, 1,P_O0,P_O }, |
403 | {__LINE__,B_ADD, 0,P_O0,P_Z, 1,P_O0,P_Z, R_ALL & ~R_DOWN, |
404 | 0, 0,P_Z,P_Z }, |
405 | {__LINE__,B_ADD, 0,P_O0,P_Z, 1,P_O0,P_Z, R_DOWN, |
406 | 0, 1,P_Z,P_Z }, |
407 | |
408 | /* Negation. */ |
409 | {__LINE__,B_NEG, 0,P_Z,P_Z, 0,0,0, R_ALL, 0, 1,P_Z,P_Z }, |
410 | {__LINE__,B_NEG, 1,P_Z,P_Z, 0,0,0, R_ALL, 0, 0,P_Z,P_Z }, |
411 | {__LINE__,B_NEG, 0,P_O,P_Z, 0,0,0, R_ALL, 0, 1,P_O,P_Z }, |
412 | {__LINE__,B_NEG, 1,P_O,P_Z, 0,0,0, R_ALL, 0, 0,P_O,P_Z }, |
413 | {__LINE__,B_NEG, 0,P_O,P_1Z, 0,0,0, R_ALL, 0, 1,P_O,P_1Z }, |
414 | {__LINE__,B_NEG, 1,P_O,P_1Z, 0,0,0, R_ALL, 0, 0,P_O,P_1Z }, |
415 | {__LINE__,B_NEG, 0,P_O,P_01Z, 0,0,0, R_ALL, 0, 1,P_O,P_01Z }, |
416 | {__LINE__,B_NEG, 1,P_O,P_01Z, 0,0,0, R_ALL, 0, 0,P_O,P_01Z }, |
417 | {__LINE__,B_NEG, 0,P_1Z,P_1Z1, 0,0,0, R_ALL, 0, 1,P_1Z,P_1Z1 }, |
418 | {__LINE__,B_NEG, 1,P_1Z,P_1Z1, 0,0,0, R_ALL, 0, 0,P_1Z,P_1Z1 }, |
419 | {__LINE__,B_NEG, 0,P_Z,P_Z1, 0,0,0, R_ALL, 0, 1,P_Z,P_Z1 }, |
420 | {__LINE__,B_NEG, 1,P_Z,P_Z1, 0,0,0, R_ALL, 0, 0,P_Z,P_Z1 }, |
421 | |
422 | /* Absolute value. */ |
423 | {__LINE__,B_ABS, 0,P_Z,P_Z, 0,0,0, R_ALL, 0, 0,P_Z,P_Z }, |
424 | {__LINE__,B_ABS, 1,P_Z,P_Z, 0,0,0, R_ALL, 0, 0,P_Z,P_Z }, |
425 | {__LINE__,B_ABS, 0,P_O,P_Z, 0,0,0, R_ALL, 0, 0,P_O,P_Z }, |
426 | {__LINE__,B_ABS, 1,P_O,P_Z, 0,0,0, R_ALL, 0, 0,P_O,P_Z }, |
427 | {__LINE__,B_ABS, 0,P_O,P_1Z, 0,0,0, R_ALL, 0, 0,P_O,P_1Z }, |
428 | {__LINE__,B_ABS, 1,P_O,P_1Z, 0,0,0, R_ALL, 0, 0,P_O,P_1Z }, |
429 | {__LINE__,B_ABS, 0,P_O,P_01Z, 0,0,0, R_ALL, 0, 0,P_O,P_01Z }, |
430 | {__LINE__,B_ABS, 1,P_O,P_01Z, 0,0,0, R_ALL, 0, 0,P_O,P_01Z }, |
431 | {__LINE__,B_ABS, 0,P_1Z,P_1Z1, 0,0,0, R_ALL, 0, 0,P_1Z,P_1Z1 }, |
432 | {__LINE__,B_ABS, 1,P_1Z,P_1Z1, 0,0,0, R_ALL, 0, 0,P_1Z,P_1Z1 }, |
433 | {__LINE__,B_ABS, 0,P_Z,P_Z1, 0,0,0, R_ALL, 0, 0,P_Z,P_Z1 }, |
434 | {__LINE__,B_ABS, 1,P_Z,P_Z1, 0,0,0, R_ALL, 0, 0,P_Z,P_Z1 }, |
435 | |
436 | /* Square root. */ |
437 | {__LINE__,B_SQRT, 0,P_Z,P_Z, 0,0,0, R_ALL, 0, 0,P_Z,P_Z }, |
438 | {__LINE__,B_SQRT, 1,P_Z,P_Z, 0,0,0, R_ALL, 0, 1,P_Z,P_Z }, |
439 | {__LINE__,B_SQRT, 0,P_O,P_1Z, 0,0,0, R_ALL, 0, 0,P_O,P_1Z }, |
440 | {__LINE__,B_SQRT, 1,P_O,P_1Z, 0,0,0, R_ALL, 0, 1,P_O,P_1Z }, |
441 | {__LINE__,B_SQRT, 0,P_O,P_01Z, 0,0,0, R_ALL, |
442 | FE_INVALID | FE_INVALID_SNAN, 0,P_O,P_11Z }, |
443 | {__LINE__,B_SQRT, 1,P_O,P_01Z, 0,0,0, R_ALL, |
444 | FE_INVALID | FE_INVALID_SNAN, 1,P_O,P_11Z }, |
445 | |
446 | {__LINE__,B_SQRT, 0,P_O,P_Z, 0,0,0, R_ALL, 0, 0,P_O,P_Z }, |
447 | {__LINE__,B_SQRT, 0,P_0O,P_Z, 0,0,0, R_ALL, 0, 0,P_0O,P_Z }, |
448 | |
449 | {__LINE__,B_SQRT, 1,P_O,P_Z, 0,0,0, R_ALL, |
450 | FE_INVALID | FE_INVALID_SQRT, 0,P_O,P_1Z }, |
451 | {__LINE__,B_SQRT, 1,P_1Z,P_1Z1, 0,0,0, R_ALL, |
452 | FE_INVALID | FE_INVALID_SQRT, 0,P_O,P_1Z }, |
453 | {__LINE__,B_SQRT, 1,P_Z,P_Z1, 0,0,0, R_ALL, |
454 | FE_INVALID | FE_INVALID_SQRT, 0,P_O,P_1Z }, |
455 | |
456 | }; |
457 | |
458 | static void |
459 | check_op(void) |
460 | { |
461 | size_t i; |
462 | int j; |
463 | tocheck_t r, a, b, x; |
464 | int raised; |
465 | |
466 | for (i = 0; i < sizeof (optests) / sizeof (optests[0]); i++) |
467 | { |
468 | a = pattern(negative: optests[i].a_sgn, exp: optests[i].a_exp, |
469 | mant: optests[i].a_mant); |
470 | b = pattern(negative: optests[i].b_sgn, exp: optests[i].b_exp, |
471 | mant: optests[i].b_mant); |
472 | x = pattern(negative: optests[i].x_sgn, exp: optests[i].x_exp, |
473 | mant: optests[i].x_mant); |
474 | for (j = 0; j < 4; j++) |
475 | if (optests[i].rmode & 1<<j) |
476 | { |
477 | fesetenv(envp: rmodes+j); |
478 | switch (optests[i].op) |
479 | { |
480 | case B_ADD: r = a + b; break; |
481 | case B_SUB: r = a - b; break; |
482 | case B_MUL: r = a * b; break; |
483 | case B_DIV: r = a / b; break; |
484 | case B_NEG: r = -a; break; |
485 | case B_ABS: r = FUNC(fabs)(x: a); break; |
486 | case B_SQRT: r = FUNC(sqrt)(a); break; |
487 | } |
488 | raised = fetestexcept(excepts: all_exceptions); |
489 | check_result(line: optests[i].line,rm: rmnames[j],expected: x,actual: r); |
490 | check_excepts(line: optests[i].line,rm: rmnames[j], |
491 | expected: optests[i].excepts,actual: raised); |
492 | } |
493 | } |
494 | } |
495 | |
496 | static void |
497 | fail_xr(int line, const char *rm, tocheck_t x, tocheck_t r, tocheck_t xx, |
498 | int xflag) |
499 | { |
500 | size_t i; |
501 | unsigned char *cx, *cr, *cxx; |
502 | |
503 | printf(format: "%s:%d:round %s:fail\n with x=0x" , __FILE__, line,rm); |
504 | cx = (unsigned char *)&x; |
505 | cr = (unsigned char *)&r; |
506 | cxx = (unsigned char *)&xx; |
507 | for (i = 0; i < sizeof (tocheck_t); i++) |
508 | printf(format: "%02x" , cx[i]); |
509 | printf(format: " r=0x" ); |
510 | for (i = 0; i < sizeof (tocheck_t); i++) |
511 | printf(format: "%02x" , cr[i]); |
512 | printf(format: " xx=0x" ); |
513 | for (i = 0; i < sizeof (tocheck_t); i++) |
514 | printf(format: "%02x" , cxx[i]); |
515 | printf(format: " inexact=%d\n" , xflag != 0); |
516 | nerrors++; |
517 | } |
518 | |
519 | static void |
520 | check_sqrt(tocheck_t a) |
521 | { |
522 | int j; |
523 | tocheck_t r0, r1, r2, x0, x1, x2; |
524 | int raised = 0; |
525 | int ok; |
526 | |
527 | for (j = 0; j < 4; j++) |
528 | { |
529 | int excepts; |
530 | |
531 | fesetenv(envp: rmodes+j); |
532 | r1 = FUNC(sqrt)(a); |
533 | excepts = fetestexcept(excepts: all_exceptions); |
534 | fesetenv(FE_DFL_ENV); |
535 | raised |= excepts & ~FE_INEXACT; |
536 | x1 = r1 * r1 - a; |
537 | if (excepts & FE_INEXACT) |
538 | { |
539 | r0 = delta(x: r1,direction: -1); r2 = delta(x: r1,direction: 1); |
540 | switch (1 << j) |
541 | { |
542 | case R_NEAREST: |
543 | x0 = r0 * r0 - a; x2 = r2 * r2 - a; |
544 | ok = fabs(x: x0) >= fabs(x: x1) && fabs(x: x1) <= fabs(x: x2); |
545 | break; |
546 | case R_ZERO: case R_DOWN: |
547 | x2 = r2 * r2 - a; |
548 | ok = x1 <= 0 && x2 >= 0; |
549 | break; |
550 | case R_UP: |
551 | x0 = r0 * r0 - a; |
552 | ok = x1 >= 0 && x0 <= 0; |
553 | break; |
554 | default: |
555 | assert(0); |
556 | } |
557 | } |
558 | else |
559 | ok = x1 == 0; |
560 | if (!ok) |
561 | fail_xr(__LINE__,rm: rmnames[j],x: a,r: r1,xx: x1,xflag: excepts&FE_INEXACT); |
562 | } |
563 | check_excepts(__LINE__,rm: "all" ,expected: 0,actual: raised); |
564 | } |
565 | |
566 | int main(int argc, char **argv) |
567 | { |
568 | int i; |
569 | |
570 | _LIB_VERSION = _IEEE_; |
571 | |
572 | /* Set up environments for rounding modes. */ |
573 | fesetenv(FE_DFL_ENV); |
574 | fesetround(FE_TONEAREST); |
575 | fegetenv(envp: rmodes+0); |
576 | fesetround(FE_TOWARDZERO); |
577 | fegetenv(envp: rmodes+1); |
578 | fesetround(FE_UPWARD); |
579 | fegetenv(envp: rmodes+2); |
580 | fesetround(FE_DOWNWARD); |
581 | fegetenv(envp: rmodes+3); |
582 | |
583 | #if defined(FE_INVALID_SOFTWARE) || defined(FE_INVALID_SQRT) |
584 | /* There's this really stupid feature of the 601... */ |
585 | fesetenv(FE_DFL_ENV); |
586 | feraiseexcept(FE_INVALID_SOFTWARE); |
587 | if (!fetestexcept(FE_INVALID_SOFTWARE)) |
588 | excepts_missing |= FE_INVALID_SOFTWARE; |
589 | fesetenv(FE_DFL_ENV); |
590 | feraiseexcept(FE_INVALID_SQRT); |
591 | if (!fetestexcept(FE_INVALID_SQRT)) |
592 | excepts_missing |= FE_INVALID_SQRT; |
593 | #endif |
594 | |
595 | check_op(); |
596 | for (i = 0; i < 100000; i++) |
597 | check_sqrt(a: pattern(negative: 0, exp: P_Rno, mant: P_R)); |
598 | for (i = 0; i < 100; i++) |
599 | check_sqrt(a: pattern(negative: 0, exp: P_Z, mant: P_R)); |
600 | check_sqrt(a: pattern(negative: 0,exp: P_Z,mant: P_Z1)); |
601 | |
602 | printf(format: "%d errors.\n" , nerrors); |
603 | return nerrors == 0 ? 0 : 1; |
604 | } |
605 | |