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