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