1 /* Native implementation of soft float functions */ 2 #include <math.h> 3 4 #if (defined(CONFIG_BSD) && !defined(__APPLE__) && !defined(__GLIBC__)) \ 5 || defined(CONFIG_SOLARIS) 6 #include <ieeefp.h> 7 #define fabsf(f) ((float)fabs(f)) 8 #else 9 #include <fenv.h> 10 #endif 11 12 #if defined(__OpenBSD__) || defined(__NetBSD__) 13 #include <sys/param.h> 14 #endif 15 16 /* 17 * Define some C99-7.12.3 classification macros and 18 * some C99-.12.4 for Solaris systems OS less than 10, 19 * or Solaris 10 systems running GCC 3.x or less. 20 * Solaris 10 with GCC4 does not need these macros as they 21 * are defined in <iso/math_c99.h> with a compiler directive 22 */ 23 #if defined(CONFIG_SOLARIS) && \ 24 ((CONFIG_SOLARIS_VERSION <= 9 ) || \ 25 ((CONFIG_SOLARIS_VERSION == 10) && (__GNUC__ < 4))) \ 26 || (defined(__OpenBSD__) && (OpenBSD < 200811)) 27 /* 28 * C99 7.12.3 classification macros 29 * and 30 * C99 7.12.14 comparison macros 31 * 32 * ... do not work on Solaris 10 using GNU CC 3.4.x. 33 * Try to workaround the missing / broken C99 math macros. 34 */ 35 #if defined(__OpenBSD__) 36 #define unordered(x, y) (isnan(x) || isnan(y)) 37 #endif 38 39 #ifdef __NetBSD__ 40 #ifndef isgreater 41 #define isgreater(x, y) __builtin_isgreater(x, y) 42 #endif 43 #ifndef isgreaterequal 44 #define isgreaterequal(x, y) __builtin_isgreaterequal(x, y) 45 #endif 46 #ifndef isless 47 #define isless(x, y) __builtin_isless(x, y) 48 #endif 49 #ifndef islessequal 50 #define islessequal(x, y) __builtin_islessequal(x, y) 51 #endif 52 #ifndef isunordered 53 #define isunordered(x, y) __builtin_isunordered(x, y) 54 #endif 55 #endif 56 57 58 #define isnormal(x) (fpclass(x) >= FP_NZERO) 59 #define isgreater(x, y) ((!unordered(x, y)) && ((x) > (y))) 60 #define isgreaterequal(x, y) ((!unordered(x, y)) && ((x) >= (y))) 61 #define isless(x, y) ((!unordered(x, y)) && ((x) < (y))) 62 #define islessequal(x, y) ((!unordered(x, y)) && ((x) <= (y))) 63 #define isunordered(x,y) unordered(x, y) 64 #endif 65 66 #if defined(__sun__) && !defined(CONFIG_NEEDS_LIBSUNMATH) 67 68 #ifndef isnan 69 # define isnan(x) \ 70 (sizeof (x) == sizeof (long double) ? isnan_ld (x) \ 71 : sizeof (x) == sizeof (double) ? isnan_d (x) \ 72 : isnan_f (x)) 73 static inline int isnan_f (float x) { return x != x; } 74 static inline int isnan_d (double x) { return x != x; } 75 static inline int isnan_ld (long double x) { return x != x; } 76 #endif 77 78 #ifndef isinf 79 # define isinf(x) \ 80 (sizeof (x) == sizeof (long double) ? isinf_ld (x) \ 81 : sizeof (x) == sizeof (double) ? isinf_d (x) \ 82 : isinf_f (x)) 83 static inline int isinf_f (float x) { return isnan (x - x); } 84 static inline int isinf_d (double x) { return isnan (x - x); } 85 static inline int isinf_ld (long double x) { return isnan (x - x); } 86 #endif 87 #endif 88 89 typedef float float32; 90 typedef double float64; 91 #ifdef FLOATX80 92 typedef long double floatx80; 93 #endif 94 95 typedef union { 96 float32 f; 97 uint32_t i; 98 } float32u; 99 typedef union { 100 float64 f; 101 uint64_t i; 102 } float64u; 103 #ifdef FLOATX80 104 typedef union { 105 floatx80 f; 106 struct { 107 uint64_t low; 108 uint16_t high; 109 } i; 110 } floatx80u; 111 #endif 112 113 /*---------------------------------------------------------------------------- 114 | Software IEC/IEEE floating-point rounding mode. 115 *----------------------------------------------------------------------------*/ 116 #if (defined(CONFIG_BSD) && !defined(__APPLE__) && !defined(__GLIBC__)) \ 117 || defined(CONFIG_SOLARIS) 118 #if defined(__OpenBSD__) 119 #define FE_RM FP_RM 120 #define FE_RP FP_RP 121 #define FE_RZ FP_RZ 122 #endif 123 enum { 124 float_round_nearest_even = FP_RN, 125 float_round_down = FP_RM, 126 float_round_up = FP_RP, 127 float_round_to_zero = FP_RZ 128 }; 129 #else 130 enum { 131 float_round_nearest_even = FE_TONEAREST, 132 float_round_down = FE_DOWNWARD, 133 float_round_up = FE_UPWARD, 134 float_round_to_zero = FE_TOWARDZERO 135 }; 136 #endif 137 138 typedef struct float_status { 139 int float_rounding_mode; 140 #ifdef FLOATX80 141 int floatx80_rounding_precision; 142 #endif 143 } float_status; 144 145 void set_float_rounding_mode(int val STATUS_PARAM); 146 #ifdef FLOATX80 147 void set_floatx80_rounding_precision(int val STATUS_PARAM); 148 #endif 149 150 /*---------------------------------------------------------------------------- 151 | Software IEC/IEEE integer-to-floating-point conversion routines. 152 *----------------------------------------------------------------------------*/ 153 float32 int32_to_float32( int STATUS_PARAM); 154 float32 uint32_to_float32( unsigned int STATUS_PARAM); 155 float64 int32_to_float64( int STATUS_PARAM); 156 float64 uint32_to_float64( unsigned int STATUS_PARAM); 157 #ifdef FLOATX80 158 floatx80 int32_to_floatx80( int STATUS_PARAM); 159 #endif 160 #ifdef FLOAT128 161 float128 int32_to_float128( int STATUS_PARAM); 162 #endif 163 float32 int64_to_float32( int64_t STATUS_PARAM); 164 float32 uint64_to_float32( uint64_t STATUS_PARAM); 165 float64 int64_to_float64( int64_t STATUS_PARAM); 166 float64 uint64_to_float64( uint64_t v STATUS_PARAM); 167 #ifdef FLOATX80 168 floatx80 int64_to_floatx80( int64_t STATUS_PARAM); 169 #endif 170 #ifdef FLOAT128 171 float128 int64_to_float128( int64_t STATUS_PARAM); 172 #endif 173 174 /*---------------------------------------------------------------------------- 175 | Software IEC/IEEE single-precision conversion constants. 176 *----------------------------------------------------------------------------*/ 177 #define float32_zero (0.0) 178 #define float32_one (1.0) 179 #define float32_ln2 (0.6931471) 180 #define float32_pi (3.1415926) 181 #define float32_half (0.5) 182 183 /*---------------------------------------------------------------------------- 184 | Software IEC/IEEE single-precision conversion routines. 185 *----------------------------------------------------------------------------*/ 186 int float32_to_int32( float32 STATUS_PARAM); 187 int float32_to_int32_round_to_zero( float32 STATUS_PARAM); 188 unsigned int float32_to_uint32( float32 a STATUS_PARAM); 189 unsigned int float32_to_uint32_round_to_zero( float32 a STATUS_PARAM); 190 int64_t float32_to_int64( float32 STATUS_PARAM); 191 int64_t float32_to_int64_round_to_zero( float32 STATUS_PARAM); 192 float64 float32_to_float64( float32 STATUS_PARAM); 193 #ifdef FLOATX80 194 floatx80 float32_to_floatx80( float32 STATUS_PARAM); 195 #endif 196 #ifdef FLOAT128 197 float128 float32_to_float128( float32 STATUS_PARAM); 198 #endif 199 200 /*---------------------------------------------------------------------------- 201 | Software IEC/IEEE single-precision operations. 202 *----------------------------------------------------------------------------*/ 203 float32 float32_round_to_int( float32 STATUS_PARAM); 204 INLINE float32 float32_add( float32 a, float32 b STATUS_PARAM) 205 { 206 return a + b; 207 } 208 INLINE float32 float32_sub( float32 a, float32 b STATUS_PARAM) 209 { 210 return a - b; 211 } 212 INLINE float32 float32_mul( float32 a, float32 b STATUS_PARAM) 213 { 214 return a * b; 215 } 216 INLINE float32 float32_div( float32 a, float32 b STATUS_PARAM) 217 { 218 return a / b; 219 } 220 float32 float32_rem( float32, float32 STATUS_PARAM); 221 float32 float32_sqrt( float32 STATUS_PARAM); 222 INLINE int float32_eq_quiet( float32 a, float32 b STATUS_PARAM) 223 { 224 return a == b; 225 } 226 INLINE int float32_le( float32 a, float32 b STATUS_PARAM) 227 { 228 return a <= b; 229 } 230 INLINE int float32_lt( float32 a, float32 b STATUS_PARAM) 231 { 232 return a < b; 233 } 234 INLINE int float32_eq( float32 a, float32 b STATUS_PARAM) 235 { 236 return a <= b && a >= b; 237 } 238 INLINE int float32_le_quiet( float32 a, float32 b STATUS_PARAM) 239 { 240 return islessequal(a, b); 241 } 242 INLINE int float32_lt_quiet( float32 a, float32 b STATUS_PARAM) 243 { 244 return isless(a, b); 245 } 246 INLINE int float32_unordered( float32 a, float32 b STATUS_PARAM) 247 { 248 return isunordered(a, b); 249 } 250 INLINE int float32_unordered_quiet( float32 a, float32 b STATUS_PARAM) 251 { 252 return isunordered(a, b); 253 } 254 int float32_compare( float32, float32 STATUS_PARAM ); 255 int float32_compare_quiet( float32, float32 STATUS_PARAM ); 256 int float32_is_signaling_nan( float32 ); 257 int float32_is_quiet_nan( float32 ); 258 int float32_is_any_nan( float32 ); 259 260 INLINE float32 float32_abs(float32 a) 261 { 262 return fabsf(a); 263 } 264 265 INLINE float32 float32_chs(float32 a) 266 { 267 return -a; 268 } 269 270 INLINE float32 float32_is_infinity(float32 a) 271 { 272 return fpclassify(a) == FP_INFINITE; 273 } 274 275 INLINE float32 float32_is_neg(float32 a) 276 { 277 float32u u; 278 u.f = a; 279 return u.i >> 31; 280 } 281 282 INLINE float32 float32_is_zero(float32 a) 283 { 284 return fpclassify(a) == FP_ZERO; 285 } 286 287 INLINE float32 float32_scalbn(float32 a, int n STATUS_PARAM) 288 { 289 return scalbnf(a, n); 290 } 291 292 /*---------------------------------------------------------------------------- 293 | Software IEC/IEEE double-precision conversion constants. 294 *----------------------------------------------------------------------------*/ 295 #define float64_zero (0.0) 296 #define float64_one (1.0) 297 #define float64_ln2 (0.693147180559945) 298 #define float64_pi (3.141592653589793) 299 #define float64_half (0.5) 300 301 /*---------------------------------------------------------------------------- 302 | Software IEC/IEEE double-precision conversion routines. 303 *----------------------------------------------------------------------------*/ 304 int float64_to_int32( float64 STATUS_PARAM ); 305 int float64_to_int32_round_to_zero( float64 STATUS_PARAM ); 306 unsigned int float64_to_uint32( float64 STATUS_PARAM ); 307 unsigned int float64_to_uint32_round_to_zero( float64 STATUS_PARAM ); 308 int64_t float64_to_int64( float64 STATUS_PARAM ); 309 int64_t float64_to_int64_round_to_zero( float64 STATUS_PARAM ); 310 uint64_t float64_to_uint64( float64 STATUS_PARAM ); 311 uint64_t float64_to_uint64_round_to_zero( float64 STATUS_PARAM ); 312 float32 float64_to_float32( float64 STATUS_PARAM ); 313 #ifdef FLOATX80 314 floatx80 float64_to_floatx80( float64 STATUS_PARAM ); 315 #endif 316 #ifdef FLOAT128 317 float128 float64_to_float128( float64 STATUS_PARAM ); 318 #endif 319 320 /*---------------------------------------------------------------------------- 321 | Software IEC/IEEE double-precision operations. 322 *----------------------------------------------------------------------------*/ 323 float64 float64_round_to_int( float64 STATUS_PARAM ); 324 float64 float64_trunc_to_int( float64 STATUS_PARAM ); 325 INLINE float64 float64_add( float64 a, float64 b STATUS_PARAM) 326 { 327 return a + b; 328 } 329 INLINE float64 float64_sub( float64 a, float64 b STATUS_PARAM) 330 { 331 return a - b; 332 } 333 INLINE float64 float64_mul( float64 a, float64 b STATUS_PARAM) 334 { 335 return a * b; 336 } 337 INLINE float64 float64_div( float64 a, float64 b STATUS_PARAM) 338 { 339 return a / b; 340 } 341 float64 float64_rem( float64, float64 STATUS_PARAM ); 342 float64 float64_sqrt( float64 STATUS_PARAM ); 343 INLINE int float64_eq_quiet( float64 a, float64 b STATUS_PARAM) 344 { 345 return a == b; 346 } 347 INLINE int float64_le( float64 a, float64 b STATUS_PARAM) 348 { 349 return a <= b; 350 } 351 INLINE int float64_lt( float64 a, float64 b STATUS_PARAM) 352 { 353 return a < b; 354 } 355 INLINE int float64_eq( float64 a, float64 b STATUS_PARAM) 356 { 357 return a <= b && a >= b; 358 } 359 INLINE int float64_le_quiet( float64 a, float64 b STATUS_PARAM) 360 { 361 return islessequal(a, b); 362 } 363 INLINE int float64_lt_quiet( float64 a, float64 b STATUS_PARAM) 364 { 365 return isless(a, b); 366 367 } 368 INLINE int float64_unordered( float64 a, float64 b STATUS_PARAM) 369 { 370 return isunordered(a, b); 371 } 372 INLINE int float64_unordered_quiet( float64 a, float64 b STATUS_PARAM) 373 { 374 return isunordered(a, b); 375 } 376 int float64_compare( float64, float64 STATUS_PARAM ); 377 int float64_compare_quiet( float64, float64 STATUS_PARAM ); 378 int float64_is_signaling_nan( float64 ); 379 int float64_is_any_nan( float64 ); 380 int float64_is_quiet_nan( float64 ); 381 382 INLINE float64 float64_abs(float64 a) 383 { 384 return fabs(a); 385 } 386 387 INLINE float64 float64_chs(float64 a) 388 { 389 return -a; 390 } 391 392 INLINE float64 float64_is_infinity(float64 a) 393 { 394 return fpclassify(a) == FP_INFINITE; 395 } 396 397 INLINE float64 float64_is_neg(float64 a) 398 { 399 float64u u; 400 u.f = a; 401 return u.i >> 63; 402 } 403 404 INLINE float64 float64_is_zero(float64 a) 405 { 406 return fpclassify(a) == FP_ZERO; 407 } 408 409 INLINE float64 float64_scalbn(float64 a, int n STATUS_PARAM) 410 { 411 return scalbn(a, n); 412 } 413 414 #ifdef FLOATX80 415 416 /*---------------------------------------------------------------------------- 417 | Software IEC/IEEE extended double-precision conversion constants. 418 *----------------------------------------------------------------------------*/ 419 #define floatx80_zero (0.0L) 420 #define floatx80_one (1.0L) 421 #define floatx80_ln2 (0.69314718055994530943L) 422 #define floatx80_pi (3.14159265358979323851L) 423 #define floatx80_half (0.5L) 424 425 /*---------------------------------------------------------------------------- 426 | Software IEC/IEEE extended double-precision conversion routines. 427 *----------------------------------------------------------------------------*/ 428 int floatx80_to_int32( floatx80 STATUS_PARAM ); 429 int floatx80_to_int32_round_to_zero( floatx80 STATUS_PARAM ); 430 int64_t floatx80_to_int64( floatx80 STATUS_PARAM); 431 int64_t floatx80_to_int64_round_to_zero( floatx80 STATUS_PARAM); 432 float32 floatx80_to_float32( floatx80 STATUS_PARAM ); 433 float64 floatx80_to_float64( floatx80 STATUS_PARAM ); 434 #ifdef FLOAT128 435 float128 floatx80_to_float128( floatx80 STATUS_PARAM ); 436 #endif 437 438 /*---------------------------------------------------------------------------- 439 | Software IEC/IEEE extended double-precision operations. 440 *----------------------------------------------------------------------------*/ 441 floatx80 floatx80_round_to_int( floatx80 STATUS_PARAM ); 442 INLINE floatx80 floatx80_add( floatx80 a, floatx80 b STATUS_PARAM) 443 { 444 return a + b; 445 } 446 INLINE floatx80 floatx80_sub( floatx80 a, floatx80 b STATUS_PARAM) 447 { 448 return a - b; 449 } 450 INLINE floatx80 floatx80_mul( floatx80 a, floatx80 b STATUS_PARAM) 451 { 452 return a * b; 453 } 454 INLINE floatx80 floatx80_div( floatx80 a, floatx80 b STATUS_PARAM) 455 { 456 return a / b; 457 } 458 floatx80 floatx80_rem( floatx80, floatx80 STATUS_PARAM ); 459 floatx80 floatx80_sqrt( floatx80 STATUS_PARAM ); 460 INLINE int floatx80_eq_quiet( floatx80 a, floatx80 b STATUS_PARAM) 461 { 462 return a == b; 463 } 464 INLINE int floatx80_le( floatx80 a, floatx80 b STATUS_PARAM) 465 { 466 return a <= b; 467 } 468 INLINE int floatx80_lt( floatx80 a, floatx80 b STATUS_PARAM) 469 { 470 return a < b; 471 } 472 INLINE int floatx80_eq( floatx80 a, floatx80 b STATUS_PARAM) 473 { 474 return a <= b && a >= b; 475 } 476 INLINE int floatx80_le_quiet( floatx80 a, floatx80 b STATUS_PARAM) 477 { 478 return islessequal(a, b); 479 } 480 INLINE int floatx80_lt_quiet( floatx80 a, floatx80 b STATUS_PARAM) 481 { 482 return isless(a, b); 483 484 } 485 INLINE int floatx80_unordered( floatx80 a, floatx80 b STATUS_PARAM) 486 { 487 return isunordered(a, b); 488 } 489 INLINE int floatx80_unordered_quiet( floatx80 a, floatx80 b STATUS_PARAM) 490 { 491 return isunordered(a, b); 492 } 493 int floatx80_compare( floatx80, floatx80 STATUS_PARAM ); 494 int floatx80_compare_quiet( floatx80, floatx80 STATUS_PARAM ); 495 int floatx80_is_signaling_nan( floatx80 ); 496 int floatx80_is_quiet_nan( floatx80 ); 497 int floatx80_is_any_nan( floatx80 ); 498 499 INLINE floatx80 floatx80_abs(floatx80 a) 500 { 501 return fabsl(a); 502 } 503 504 INLINE floatx80 floatx80_chs(floatx80 a) 505 { 506 return -a; 507 } 508 509 INLINE floatx80 floatx80_is_infinity(floatx80 a) 510 { 511 return fpclassify(a) == FP_INFINITE; 512 } 513 514 INLINE floatx80 floatx80_is_neg(floatx80 a) 515 { 516 floatx80u u; 517 u.f = a; 518 return u.i.high >> 15; 519 } 520 521 INLINE floatx80 floatx80_is_zero(floatx80 a) 522 { 523 return fpclassify(a) == FP_ZERO; 524 } 525 526 INLINE floatx80 floatx80_scalbn(floatx80 a, int n STATUS_PARAM) 527 { 528 return scalbnl(a, n); 529 } 530 531 #endif 532
