systfloat.c revision 1.7 1 /* $NetBSD: systfloat.c,v 1.7 2004/04/15 19:01:57 matt Exp $ */
2
3 /* This is a derivative work. */
4
5 /*-
6 * Copyright (c) 2001 The NetBSD Foundation, Inc.
7 * All rights reserved.
8 *
9 * This code is derived from software contributed to The NetBSD Foundation
10 * by Ross Harvey.
11 *
12 * Redistribution and use in source and binary forms, with or without
13 * modification, are permitted provided that the following conditions
14 * are met:
15 * 1. Redistributions of source code must retain the above copyright
16 * notice, this list of conditions and the following disclaimer.
17 * 2. Redistributions in binary form must reproduce the above copyright
18 * notice, this list of conditions and the following disclaimer in the
19 * documentation and/or other materials provided with the distribution.
20 * 3. All advertising materials mentioning features or use of this software
21 * must display the following acknowledgement:
22 * This product includes software developed by the NetBSD
23 * Foundation, Inc. and its contributors.
24 * 4. Neither the name of The NetBSD Foundation nor the names of its
25 * contributors may be used to endorse or promote products derived
26 * from this software without specific prior written permission.
27 *
28 * THIS SOFTWARE IS PROVIDED BY THE NETBSD FOUNDATION, INC. AND CONTRIBUTORS
29 * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
30 * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
31 * PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE FOUNDATION OR CONTRIBUTORS
32 * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
33 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
34 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
35 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
36 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
37 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
38 * POSSIBILITY OF SUCH DAMAGE.
39 */
40
41 /*
42 ===============================================================================
43
44 This C source file is part of TestFloat, Release 2a, a package of programs
45 for testing the correctness of floating-point arithmetic complying to the
46 IEC/IEEE Standard for Floating-Point.
47
48 Written by John R. Hauser. More information is available through the Web
49 page `http://HTTP.CS.Berkeley.EDU/~jhauser/arithmetic/TestFloat.html'.
50
51 THIS SOFTWARE IS DISTRIBUTED AS IS, FOR FREE. Although reasonable effort
52 has been made to avoid it, THIS SOFTWARE MAY CONTAIN FAULTS THAT WILL AT
53 TIMES RESULT IN INCORRECT BEHAVIOR. USE OF THIS SOFTWARE IS RESTRICTED TO
54 PERSONS AND ORGANIZATIONS WHO CAN AND WILL TAKE FULL RESPONSIBILITY FOR ANY
55 AND ALL LOSSES, COSTS, OR OTHER PROBLEMS ARISING FROM ITS USE.
56
57 Derivative works are acceptable, even for commercial purposes, so long as
58 (1) they include prominent notice that the work is derivative, and (2) they
59 include prominent notice akin to these four paragraphs for those parts of
60 this code that are retained.
61
62 ===============================================================================
63 */
64
65 #include <sys/cdefs.h>
66 #ifndef __lint
67 __RCSID("$NetBSD: systfloat.c,v 1.7 2004/04/15 19:01:57 matt Exp $");
68 #endif
69
70 #include <math.h>
71 #include <ieeefp.h>
72 #include "milieu.h"
73 #include "softfloat.h"
74 #include "systfloat.h"
75 #include "systflags.h"
76 #include "systmodes.h"
77
78 typedef union {
79 float32 f32;
80 float f;
81 } union32;
82 typedef union {
83 float64 f64;
84 double d;
85 } union64;
86 #if defined( FLOATX80 ) && defined( LONG_DOUBLE_IS_FLOATX80 )
87 typedef union {
88 floatx80 fx80;
89 long double ld;
90 } unionx80;
91 #endif
92 #if defined( FLOAT128 ) && defined( LONG_DOUBLE_IS_FLOAT128 )
93 typedef union {
94 float128 f128;
95 long double ld;
96 } union128;
97 #endif
98
99 fp_except
100 syst_float_flags_clear(void)
101 {
102 return fpsetsticky(0)
103 & (FP_X_IMP | FP_X_UFL | FP_X_OFL | FP_X_DZ | FP_X_INV);
104 }
105
106 void
107 syst_float_set_rounding_mode(fp_rnd direction)
108 {
109 fpsetround(direction);
110 fpsetmask(0);
111 }
112
113 float32 syst_int32_to_float32( int32 a )
114 {
115 const union32 uz = { .f = a };
116
117 return uz.f32;
118
119 }
120
121 float64 syst_int32_to_float64( int32 a )
122 {
123 const union64 uz = { .d = a };
124
125 return uz.f64;
126
127 }
128
129 #if defined( FLOATX80 ) && defined( LONG_DOUBLE_IS_FLOATX80 )
130
131 floatx80 syst_int32_to_floatx80( int32 a )
132 {
133 const unionx80 uz = { .ld = a };
134
135 return uz.fx80;
136
137 }
138
139 #endif
140
141 #if defined( FLOAT128 ) && defined( LONG_DOUBLE_IS_FLOAT128 )
142
143 float128 syst_int32_to_float128( int32 a )
144 {
145 const union128 uz = { .ld = a };
146
147 return uz.f128;
148
149 }
150
151 #endif
152
153 #ifdef BITS64
154
155 float32 syst_int64_to_float32( int64 a )
156 {
157 const union32 uz = { .f = a };
158
159 return uz.f32;
160 }
161
162 float64 syst_int64_to_float64( int64 a )
163 {
164 const union64 uz = { .d = a };
165
166 return uz.f64;
167 }
168
169 #if defined( FLOATX80 ) && defined( LONG_DOUBLE_IS_FLOATX80 )
170
171 floatx80 syst_int64_to_floatx80( int64 a )
172 {
173 const unionx80 uz = { .ld = a };
174
175 return uz.fx80;
176 }
177
178 #endif
179
180 #if defined( FLOAT128 ) && defined( LONG_DOUBLE_IS_FLOAT128 )
181
182 float128 syst_int64_to_float128( int64 a )
183 {
184 const union128 uz = { .ld = a };
185
186 return uz.f128;
187 }
188
189 #endif
190
191 #endif
192
193 int32 syst_float32_to_int32_round_to_zero( float32 a )
194 {
195 const union32 uz = { .f32 = a };
196
197 return uz.f;
198
199 }
200
201 #ifdef BITS64
202
203 int64 syst_float32_to_int64_round_to_zero( float32 a )
204 {
205 const union32 uz = { .f32 = a };
206
207 return uz.f;
208 }
209
210 #endif
211
212 float64 syst_float32_to_float64( float32 a )
213 {
214 const union32 ua = { .f32 = a };
215 union64 uz;
216
217 uz.d = ua.f;
218 return uz.f64;
219
220 }
221
222 #if defined( FLOATX80 ) && defined( LONG_DOUBLE_IS_FLOATX80 )
223
224 floatx80 syst_float32_to_floatx80( float32 a )
225 {
226 const union32 ua = { .f32 = a };
227 unionx80 uz;
228
229 uz.ld = ua.f;
230 return uz.fx80;
231 }
232
233 #endif
234
235 #if defined( FLOAT128 ) && defined( LONG_DOUBLE_IS_FLOAT128 )
236
237 float128 syst_float32_to_float128( float32 a )
238 {
239 const union32 ua = { .f32 = a };
240 union128 ub;
241
242 ub.ld = ua.f;
243 return ub.f128;
244 }
245
246 #endif
247
248 float32 syst_float32_add( float32 a, float32 b )
249 {
250 const union32 ua = { .f32 = a }, ub = { .f32 = b };
251 union32 uz;
252
253 uz.f = ua.f + ub.f;
254 return uz.f32;
255 }
256
257 float32 syst_float32_sub( float32 a, float32 b )
258 {
259 const union32 ua = { .f32 = a }, ub = { .f32 = b };
260 union32 uz;
261
262 uz.f = ua.f - ub.f;
263 return uz.f32;
264 }
265
266 float32 syst_float32_mul( float32 a, float32 b )
267 {
268 const union32 ua = { .f32 = a }, ub = { .f32 = b };
269 union32 uz;
270
271 uz.f = ua.f * ub.f;
272 return uz.f32;
273 }
274
275 float32 syst_float32_div( float32 a, float32 b )
276 {
277 const union32 ua = { .f32 = a }, ub = { .f32 = b };
278 union32 uz;
279
280 uz.f = ua.f / ub.f;
281 return uz.f32;
282 }
283
284 flag syst_float32_eq( float32 a, float32 b )
285 {
286 const union32 ua = { .f32 = a }, ub = { .f32 = b };
287
288 return ua.f == ub.f;
289 }
290
291 flag syst_float32_le( float32 a, float32 b )
292 {
293 const union32 ua = { .f32 = a }, ub = { .f32 = b };
294
295 return ua.f <= ub.f;
296 }
297
298 flag syst_float32_lt( float32 a, float32 b )
299 {
300 const union32 ua = { .f32 = a }, ub = { .f32 = b };
301
302 return ua.f < ub.f;
303 }
304
305 int32 syst_float64_to_int32_round_to_zero( float64 a )
306 {
307 const union64 uz = { .f64 = a };
308
309 return uz.d;
310 }
311
312 #ifdef BITS64
313
314 int64 syst_float64_to_int64_round_to_zero( float64 a )
315 {
316 const union64 uz = { .f64 = a };
317
318 return uz.d;
319 }
320
321 #endif
322
323 float32 syst_float64_to_float32( float64 a )
324 {
325 const union64 ua = { .f64 = a };
326 union32 uz;
327
328 uz.f = ua.d;
329 return uz.f32;
330 }
331
332 #if defined( FLOATX80 ) && defined( LONG_DOUBLE_IS_FLOATX80 )
333
334 floatx80 syst_float64_to_floatx80( float64 a )
335 {
336 const union64 ua = { .f64 = a };
337 unionx80 u;
338
339 u.ld = ua.d;
340 return u.fx80;
341 }
342
343 #endif
344
345 #if defined( FLOAT128 ) && defined( LONG_DOUBLE_IS_FLOAT128 )
346
347 float128 syst_float64_to_float128( float64 a )
348 {
349 const union64 ua = { .f64 = a };
350 union128 uz;
351
352 uz.ld = ua.d;
353 return uz.f128;
354 }
355
356 #endif
357
358 float64 syst_float64_add( float64 a, float64 b )
359 {
360 const union64 ua = { .f64 = a }, ub = { .f64 = b };
361 union64 uz;
362
363 uz.d = ua.d + ub.d;
364 return uz.f64;
365 }
366
367 float64 syst_float64_sub( float64 a, float64 b )
368 {
369 const union64 ua = { .f64 = a }, ub = { .f64 = b };
370 union64 uz;
371
372 uz.d = ua.d - ub.d;
373 return uz.f64;
374 }
375
376 float64 syst_float64_mul( float64 a, float64 b )
377 {
378 const union64 ua = { .f64 = a }, ub = { .f64 = b };
379 union64 uz;
380
381 uz.d = ua.d * ub.d;
382 return uz.f64;
383 }
384
385 float64 syst_float64_div( float64 a, float64 b )
386 {
387 const union64 ua = { .f64 = a }, ub = { .f64 = b };
388 union64 uz;
389
390 uz.d = ua.d / ub.d;
391 return uz.f64;
392 }
393
394 float64 syst_float64_sqrt( float64 a )
395 {
396 const union64 ua = { .f64 = a };
397 union64 uz;
398
399 uz.d = sqrt(ua.d);
400 return uz.f64;
401 }
402
403 flag syst_float64_eq( float64 a, float64 b )
404 {
405 const union64 ua = { .f64 = a }, ub = { .f64 = b };
406
407 return ua.d == ub.d;
408 }
409
410 flag syst_float64_le( float64 a, float64 b )
411 {
412 const union64 ua = { .f64 = a }, ub = { .f64 = b };
413
414 return ua.d <= ub.d;
415 }
416
417 flag syst_float64_lt( float64 a, float64 b )
418 {
419 const union64 ua = { .f64 = a }, ub = { .f64 = b };
420
421 return ua.d < ub.d;
422 }
423
424 #if defined( FLOATX80 ) && defined( LONG_DOUBLE_IS_FLOATX80 )
425
426 int32 syst_floatx80_to_int32_round_to_zero( floatx80 a )
427 {
428 const unionx80 uz = { .fx80 = a };
429
430 return uz.ld;
431 }
432
433 #ifdef BITS64
434
435 int64 syst_floatx80_to_int64_round_to_zero( floatx80 a )
436 {
437 const unionx80 uz = { .fx80 = a };
438
439 return uz.ld;
440 }
441
442 #endif
443
444 float32 syst_floatx80_to_float32( floatx80 a )
445 {
446 const unionx80 ua = { .fx80 = a };
447 union32 uz;
448
449 uz.f = ua.ld;
450 return uz.f32;
451 }
452
453 float64 syst_floatx80_to_float64( floatx80 a )
454 {
455 const unionx80 ua = { .fx80 = a };
456 union64 uz;
457
458 uz.d = ua.ld;
459 return uz.f64;
460 }
461
462 floatx80 syst_floatx80_add( floatx80 a, floatx80 b )
463 {
464 const unionx80 ua = { .fx80 = a }, ub = { .fx80 = b };
465 unionx80 uz;
466
467 uz.ld = ua.ld + ub.ld;
468 return uz.fx80;
469 }
470
471 floatx80 syst_floatx80_sub( floatx80 a, floatx80 b )
472 {
473 const unionx80 ua = { .fx80 = a }, ub = { .fx80 = b };
474 unionx80 uz;
475
476 uz.ld = ua.ld - ub.ld;
477 return uz.fx80;
478 }
479
480 floatx80 syst_floatx80_mul( floatx80 a, floatx80 b )
481 {
482 const unionx80 ua = { .fx80 = a }, ub = { .fx80 = b };
483 unionx80 uz;
484
485 uz.ld = ua.ld * ub.ld;
486 return uz.fx80;
487 }
488
489 floatx80 syst_floatx80_div( floatx80 a, floatx80 b )
490 {
491 const unionx80 ua = { .fx80 = a }, ub = { .fx80 = b };
492 unionx80 uz;
493
494 uz.ld = ua.ld / ub.ld;
495 return uz.fx80;
496 }
497
498 flag syst_floatx80_eq( floatx80 a, floatx80 b )
499 {
500 const unionx80 ua = { .fx80 = a }, ub = { .fx80 = b };
501
502 return ua.ld == ub.ld;
503 }
504
505 flag syst_floatx80_le( floatx80 a, floatx80 b )
506 {
507 const unionx80 ua = { .fx80 = a }, ub = { .fx80 = b };
508
509 return ua.ld <= ub.ld;
510 }
511
512 flag syst_floatx80_lt( floatx80 a, floatx80 b )
513 {
514 const unionx80 ua = { .fx80 = a }, ub = { .fx80 = b };
515
516 return ua.ld < ub.ld;
517 }
518
519 #endif
520
521 #if defined( FLOAT128 ) && defined( LONG_DOUBLE_IS_FLOAT128 )
522
523 int32 syst_float128_to_int32_round_to_zero( float128 a )
524 {
525 const union128 ua = { .f128 = a };
526
527 return ua.ld;
528 }
529
530 #ifdef BITS64
531
532 int64 syst_float128_to_int64_round_to_zero( float128 a )
533 {
534 const union128 ua = { .f128 = a };
535
536 return ua.ld;
537 }
538
539 #endif
540
541 float32 syst_float128_to_float32( float128 a )
542 {
543 const union128 ua = { .f128 = a };
544 union32 uz;
545
546 uz.f = ua.ld;
547 return uz.f32;
548
549 }
550
551 float64 syst_float128_to_float64( float128 a )
552 {
553 const union128 ua = { .f128 = a };
554 union64 uz;
555
556 uz.d = ua.ld;
557 return uz.f64;
558 }
559
560 float128 syst_float128_add( float128 a, float128 b )
561 {
562 const union128 ua = { .f128 = a }, ub = { .f128 = b };
563 union128 uz;
564
565 uz.ld = ua.ld + ub.ld;
566 return uz.f128;
567
568 }
569
570 float128 syst_float128_sub( float128 a, float128 b )
571 {
572 const union128 ua = { .f128 = a }, ub = { .f128 = b };
573 union128 uz;
574
575 uz.ld = ua.ld - ub.ld;
576 return uz.f128;
577 }
578
579 float128 syst_float128_mul( float128 a, float128 b )
580 {
581 const union128 ua = { .f128 = a }, ub = { .f128 = b };
582 union128 uz;
583
584 uz.ld = ua.ld * ub.ld;
585 return uz.f128;
586 }
587
588 float128 syst_float128_div( float128 a, float128 b )
589 {
590 const union128 ua = { .f128 = a }, ub = { .f128 = b };
591 union128 uz;
592
593 uz.ld = ua.ld / ub.ld;
594 return uz.f128;
595 }
596
597 flag syst_float128_eq( float128 a, float128 b )
598 {
599 const union128 ua = { .f128 = a }, ub = { .f128 = b };
600
601 return ua.ld == ub.ld;
602 }
603
604 flag syst_float128_le( float128 a, float128 b )
605 {
606 const union128 ua = { .f128 = a }, ub = { .f128 = b };
607
608 return ua.ld <= ub.ld;
609 }
610
611 flag syst_float128_lt( float128 a, float128 b )
612 {
613 const union128 ua = { .f128 = a }, ub = { .f128 = b };
614
615 return ua.ld < ub.ld;
616 }
617
618 #endif
619