1/* 2 * ==================================================== 3 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. 4 * 5 * Developed at SunPro, a Sun Microsystems, Inc. business. 6 * Permission to use, copy, modify, and distribute this 7 * software is freely granted, provided that this notice 8 * is preserved. 9 * ==================================================== 10 */ 11 12/* 13 * from: @(#)fdlibm.h 5.1 93/09/24 14 * $FreeBSD$ 15 */ 16 17#ifndef _MATH_PRIVATE_H_ 18#define _MATH_PRIVATE_H_ 19 20#include <sys/types.h> 21#include <machine/endian.h> 22 23#include <stdint.h> 24typedef uint32_t u_int32_t; 25typedef uint64_t u_int64_t; 26 27/* 28 * The original fdlibm code used statements like: 29 * n0 = ((*(int*)&one)>>29)^1; * index of high word * 30 * ix0 = *(n0+(int*)&x); * high word of x * 31 * ix1 = *((1-n0)+(int*)&x); * low word of x * 32 * to dig two 32 bit words out of the 64 bit IEEE floating point 33 * value. That is non-ANSI, and, moreover, the gcc instruction 34 * scheduler gets it wrong. We instead use the following macros. 35 * Unlike the original code, we determine the endianness at compile 36 * time, not at run time; I don't see much benefit to selecting 37 * endianness at run time. 38 */ 39 40/* 41 * A union which permits us to convert between a double and two 32 bit 42 * ints. 43 */ 44 45#ifdef __arm__ 46#if defined(__VFP_FP__) || defined(__ARM_EABI__) 47#define IEEE_WORD_ORDER BYTE_ORDER 48#else 49#define IEEE_WORD_ORDER BIG_ENDIAN 50#endif 51#else /* __arm__ */ 52#define IEEE_WORD_ORDER BYTE_ORDER 53#endif 54 55#if IEEE_WORD_ORDER == BIG_ENDIAN 56 57typedef union 58{ 59 double value; 60 struct 61 { 62 u_int32_t msw; 63 u_int32_t lsw; 64 } parts; 65 struct 66 { 67 u_int64_t w; 68 } xparts; 69} ieee_double_shape_type; 70 71#endif 72 73#if IEEE_WORD_ORDER == LITTLE_ENDIAN 74 75typedef union 76{ 77 double value; 78 struct 79 { 80 u_int32_t lsw; 81 u_int32_t msw; 82 } parts; 83 struct 84 { 85 u_int64_t w; 86 } xparts; 87} ieee_double_shape_type; 88 89#endif 90 91/* Get two 32 bit ints from a double. */ 92 93#define EXTRACT_WORDS(ix0,ix1,d) \ 94do { \ 95 ieee_double_shape_type ew_u; \ 96 ew_u.value = (d); \ 97 (ix0) = ew_u.parts.msw; \ 98 (ix1) = ew_u.parts.lsw; \ 99} while (0) 100 101/* Get a 64-bit int from a double. */ 102#define EXTRACT_WORD64(ix,d) \ 103do { \ 104 ieee_double_shape_type ew_u; \ 105 ew_u.value = (d); \ 106 (ix) = ew_u.xparts.w; \ 107} while (0) 108 109/* Get the more significant 32 bit int from a double. */ 110 111#define GET_HIGH_WORD(i,d) \ 112do { \ 113 ieee_double_shape_type gh_u; \ 114 gh_u.value = (d); \ 115 (i) = gh_u.parts.msw; \ 116} while (0) 117 118/* Get the less significant 32 bit int from a double. */ 119 120#define GET_LOW_WORD(i,d) \ 121do { \ 122 ieee_double_shape_type gl_u; \ 123 gl_u.value = (d); \ 124 (i) = gl_u.parts.lsw; \ 125} while (0) 126 127/* Set a double from two 32 bit ints. */ 128 129#define INSERT_WORDS(d,ix0,ix1) \ 130do { \ 131 ieee_double_shape_type iw_u; \ 132 iw_u.parts.msw = (ix0); \ 133 iw_u.parts.lsw = (ix1); \ 134 (d) = iw_u.value; \ 135} while (0) 136 137/* Set a double from a 64-bit int. */ 138#define INSERT_WORD64(d,ix) \ 139do { \ 140 ieee_double_shape_type iw_u; \ 141 iw_u.xparts.w = (ix); \ 142 (d) = iw_u.value; \ 143} while (0) 144 145/* Set the more significant 32 bits of a double from an int. */ 146 147#define SET_HIGH_WORD(d,v) \ 148do { \ 149 ieee_double_shape_type sh_u; \ 150 sh_u.value = (d); \ 151 sh_u.parts.msw = (v); \ 152 (d) = sh_u.value; \ 153} while (0) 154 155/* Set the less significant 32 bits of a double from an int. */ 156 157#define SET_LOW_WORD(d,v) \ 158do { \ 159 ieee_double_shape_type sl_u; \ 160 sl_u.value = (d); \ 161 sl_u.parts.lsw = (v); \ 162 (d) = sl_u.value; \ 163} while (0) 164 165/* 166 * A union which permits us to convert between a float and a 32 bit 167 * int. 168 */ 169 170typedef union 171{ 172 float value; 173 /* FIXME: Assumes 32 bit int. */ 174 unsigned int word; 175} ieee_float_shape_type; 176 177/* Get a 32 bit int from a float. */ 178 179#define GET_FLOAT_WORD(i,d) \ 180do { \ 181 ieee_float_shape_type gf_u; \ 182 gf_u.value = (d); \ 183 (i) = gf_u.word; \ 184} while (0) 185 186/* Set a float from a 32 bit int. */ 187 188#define SET_FLOAT_WORD(d,i) \ 189do { \ 190 ieee_float_shape_type sf_u; \ 191 sf_u.word = (i); \ 192 (d) = sf_u.value; \ 193} while (0) 194 195/* 196 * Get expsign and mantissa as 16 bit and 64 bit ints from an 80 bit long 197 * double. 198 */ 199 200#define EXTRACT_LDBL80_WORDS(ix0,ix1,d) \ 201do { \ 202 union IEEEl2bits ew_u; \ 203 ew_u.e = (d); \ 204 (ix0) = ew_u.xbits.expsign; \ 205 (ix1) = ew_u.xbits.man; \ 206} while (0) 207 208/* 209 * Get expsign and mantissa as one 16 bit and two 64 bit ints from a 128 bit 210 * long double. 211 */ 212 213#define EXTRACT_LDBL128_WORDS(ix0,ix1,ix2,d) \ 214do { \ 215 union IEEEl2bits ew_u; \ 216 ew_u.e = (d); \ 217 (ix0) = ew_u.xbits.expsign; \ 218 (ix1) = ew_u.xbits.manh; \ 219 (ix2) = ew_u.xbits.manl; \ 220} while (0) 221 222/* Get expsign as a 16 bit int from a long double. */ 223 224#define GET_LDBL_EXPSIGN(i,d) \ 225do { \ 226 union IEEEl2bits ge_u; \ 227 ge_u.e = (d); \ 228 (i) = ge_u.xbits.expsign; \ 229} while (0) 230 231/* 232 * Set an 80 bit long double from a 16 bit int expsign and a 64 bit int 233 * mantissa. 234 */ 235 236#define INSERT_LDBL80_WORDS(d,ix0,ix1) \ 237do { \ 238 union IEEEl2bits iw_u; \ 239 iw_u.xbits.expsign = (ix0); \ 240 iw_u.xbits.man = (ix1); \ 241 (d) = iw_u.e; \ 242} while (0) 243 244/* 245 * Set a 128 bit long double from a 16 bit int expsign and two 64 bit ints 246 * comprising the mantissa. 247 */ 248 249#define INSERT_LDBL128_WORDS(d,ix0,ix1,ix2) \ 250do { \ 251 union IEEEl2bits iw_u; \ 252 iw_u.xbits.expsign = (ix0); \ 253 iw_u.xbits.manh = (ix1); \ 254 iw_u.xbits.manl = (ix2); \ 255 (d) = iw_u.e; \ 256} while (0) 257 258/* Set expsign of a long double from a 16 bit int. */ 259 260#define SET_LDBL_EXPSIGN(d,v) \ 261do { \ 262 union IEEEl2bits se_u; \ 263 se_u.e = (d); \ 264 se_u.xbits.expsign = (v); \ 265 (d) = se_u.e; \ 266} while (0) 267 268#ifdef __i386__ 269/* Long double constants are broken on i386. */ 270#define LD80C(m, ex, v) { \ 271 .xbits.man = __CONCAT(m, ULL), \ 272 .xbits.expsign = (0x3fff + (ex)) | ((v) < 0 ? 0x8000 : 0), \ 273} 274#else 275/* The above works on non-i386 too, but we use this to check v. */ 276#define LD80C(m, ex, v) { .e = (v), } 277#endif 278 279#ifdef FLT_EVAL_METHOD 280/* 281 * Attempt to get strict C99 semantics for assignment with non-C99 compilers. 282 */ 283#if FLT_EVAL_METHOD == 0 || __GNUC__ == 0 284#define STRICT_ASSIGN(type, lval, rval) ((lval) = (rval)) 285#else 286#define STRICT_ASSIGN(type, lval, rval) do { \ 287 volatile type __lval; \ 288 \ 289 if (sizeof(type) >= sizeof(long double)) \ 290 (lval) = (rval); \ 291 else { \ 292 __lval = (rval); \ 293 (lval) = __lval; \ 294 } \ 295} while (0) 296#endif 297#endif /* FLT_EVAL_METHOD */ 298 299/* Support switching the mode to FP_PE if necessary. */ 300#if defined(__i386__) && !defined(NO_FPSETPREC) 301#define ENTERI() \ 302 long double __retval; \ 303 fp_prec_t __oprec; \ 304 \ 305 if ((__oprec = fpgetprec()) != FP_PE) \ 306 fpsetprec(FP_PE) 307#define RETURNI(x) do { \ 308 __retval = (x); \ 309 if (__oprec != FP_PE) \ 310 fpsetprec(__oprec); \ 311 RETURNF(__retval); \ 312} while (0) 313#else 314#define ENTERI(x) 315#define RETURNI(x) RETURNF(x) 316#endif 317 318/* Default return statement if hack*_t() is not used. */ 319#define RETURNF(v) return (v) 320 321/* 322 * 2sum gives the same result as 2sumF without requiring |a| >= |b| or 323 * a == 0, but is slower. 324 */ 325#define _2sum(a, b) do { \ 326 __typeof(a) __s, __w; \ 327 \ 328 __w = (a) + (b); \ 329 __s = __w - (a); \ 330 (b) = ((a) - (__w - __s)) + ((b) - __s); \ 331 (a) = __w; \ 332} while (0) 333 334/* 335 * 2sumF algorithm. 336 * 337 * "Normalize" the terms in the infinite-precision expression a + b for 338 * the sum of 2 floating point values so that b is as small as possible 339 * relative to 'a'. (The resulting 'a' is the value of the expression in 340 * the same precision as 'a' and the resulting b is the rounding error.) 341 * |a| must be >= |b| or 0, b's type must be no larger than 'a's type, and 342 * exponent overflow or underflow must not occur. This uses a Theorem of 343 * Dekker (1971). See Knuth (1981) 4.2.2 Theorem C. The name "TwoSum" 344 * is apparently due to Skewchuk (1997). 345 * 346 * For this to always work, assignment of a + b to 'a' must not retain any 347 * extra precision in a + b. This is required by C standards but broken 348 * in many compilers. The brokenness cannot be worked around using 349 * STRICT_ASSIGN() like we do elsewhere, since the efficiency of this 350 * algorithm would be destroyed by non-null strict assignments. (The 351 * compilers are correct to be broken -- the efficiency of all floating 352 * point code calculations would be destroyed similarly if they forced the 353 * conversions.) 354 * 355 * Fortunately, a case that works well can usually be arranged by building 356 * any extra precision into the type of 'a' -- 'a' should have type float_t, 357 * double_t or long double. b's type should be no larger than 'a's type. 358 * Callers should use these types with scopes as large as possible, to 359 * reduce their own extra-precision and efficiciency problems. In 360 * particular, they shouldn't convert back and forth just to call here. 361 */ 362#ifdef DEBUG 363#define _2sumF(a, b) do { \ 364 __typeof(a) __w; \ 365 volatile __typeof(a) __ia, __ib, __r, __vw; \ 366 \ 367 __ia = (a); \ 368 __ib = (b); \ 369 assert(__ia == 0 || fabsl(__ia) >= fabsl(__ib)); \ 370 \ 371 __w = (a) + (b); \ 372 (b) = ((a) - __w) + (b); \ 373 (a) = __w; \ 374 \ 375 /* The next 2 assertions are weak if (a) is already long double. */ \ 376 assert((long double)__ia + __ib == (long double)(a) + (b)); \ 377 __vw = __ia + __ib; \ 378 __r = __ia - __vw; \ 379 __r += __ib; \ 380 assert(__vw == (a) && __r == (b)); \ 381} while (0) 382#else /* !DEBUG */ 383#define _2sumF(a, b) do { \ 384 __typeof(a) __w; \ 385 \ 386 __w = (a) + (b); \ 387 (b) = ((a) - __w) + (b); \ 388 (a) = __w; \ 389} while (0) 390#endif /* DEBUG */ 391 392/* 393 * Set x += c, where x is represented in extra precision as a + b. 394 * x must be sufficiently normalized and sufficiently larger than c, 395 * and the result is then sufficiently normalized. 396 * 397 * The details of ordering are that |a| must be >= |c| (so that (a, c) 398 * can be normalized without extra work to swap 'a' with c). The details of 399 * the normalization are that b must be small relative to the normalized 'a'. 400 * Normalization of (a, c) makes the normalized c tiny relative to the 401 * normalized a, so b remains small relative to 'a' in the result. However, 402 * b need not ever be tiny relative to 'a'. For example, b might be about 403 * 2**20 times smaller than 'a' to give about 20 extra bits of precision. 404 * That is usually enough, and adding c (which by normalization is about 405 * 2**53 times smaller than a) cannot change b significantly. However, 406 * cancellation of 'a' with c in normalization of (a, c) may reduce 'a' 407 * significantly relative to b. The caller must ensure that significant 408 * cancellation doesn't occur, either by having c of the same sign as 'a', 409 * or by having |c| a few percent smaller than |a|. Pre-normalization of 410 * (a, b) may help. 411 * 412 * This is is a variant of an algorithm of Kahan (see Knuth (1981) 4.2.2 413 * exercise 19). We gain considerable efficiency by requiring the terms to 414 * be sufficiently normalized and sufficiently increasing. 415 */ 416#define _3sumF(a, b, c) do { \ 417 __typeof(a) __tmp; \ 418 \ 419 __tmp = (c); \ 420 _2sumF(__tmp, (a)); \ 421 (b) += (a); \ 422 (a) = __tmp; \ 423} while (0) 424 425/* 426 * Common routine to process the arguments to nan(), nanf(), and nanl(). 427 */ 428void _scan_nan(uint32_t *__words, int __num_words, const char *__s); 429 430#ifdef _COMPLEX_H 431 432/* 433 * C99 specifies that complex numbers have the same representation as 434 * an array of two elements, where the first element is the real part 435 * and the second element is the imaginary part. 436 */ 437typedef union { 438 float complex f; 439 float a[2]; 440} float_complex; 441typedef union { 442 double complex f; 443 double a[2]; 444} double_complex; 445typedef union { 446 long double complex f; 447 long double a[2]; 448} long_double_complex; 449#define REALPART(z) ((z).a[0]) 450#define IMAGPART(z) ((z).a[1]) 451 452/* 453 * Inline functions that can be used to construct complex values. 454 * 455 * The C99 standard intends x+I*y to be used for this, but x+I*y is 456 * currently unusable in general since gcc introduces many overflow, 457 * underflow, sign and efficiency bugs by rewriting I*y as 458 * (0.0+I)*(y+0.0*I) and laboriously computing the full complex product. 459 * In particular, I*Inf is corrupted to NaN+I*Inf, and I*-0 is corrupted 460 * to -0.0+I*0.0. 461 */ 462static __inline float complex 463cpackf(float x, float y) 464{ 465 float_complex z; 466 467 REALPART(z) = x; 468 IMAGPART(z) = y; 469 return (z.f); 470} 471 472static __inline double complex 473cpack(double x, double y) 474{ 475 double_complex z; 476 477 REALPART(z) = x; 478 IMAGPART(z) = y; 479 return (z.f); 480} 481 482static __inline long double complex 483cpackl(long double x, long double y) 484{ 485 long_double_complex z; 486 487 REALPART(z) = x; 488 IMAGPART(z) = y; 489 return (z.f); 490} 491#endif /* _COMPLEX_H */ 492 493#ifdef __GNUCLIKE_ASM 494 495/* Asm versions of some functions. */ 496 497#ifdef __amd64__ 498static __inline int 499irint(double x) 500{ 501 int n; 502 503 asm("cvtsd2si %1,%0" : "=r" (n) : "x" (x)); 504 return (n); 505} 506#define HAVE_EFFICIENT_IRINT 507#endif 508 509#ifdef __i386__ 510static __inline int 511irint(double x) 512{ 513 int n; 514 515 asm("fistl %0" : "=m" (n) : "t" (x)); 516 return (n); 517} 518#define HAVE_EFFICIENT_IRINT 519#endif 520 521#if defined(__amd64__) || defined(__i386__) 522static __inline int 523irintl(long double x) 524{ 525 int n; 526 527 asm("fistl %0" : "=m" (n) : "t" (x)); 528 return (n); 529} 530#define HAVE_EFFICIENT_IRINTL 531#endif 532 533#endif /* __GNUCLIKE_ASM */ 534 535#ifdef DEBUG 536#if defined(__amd64__) || defined(__i386__) 537#define breakpoint() asm("int $3") 538#else 539#include <signal.h> 540 541#define breakpoint() raise(SIGTRAP) 542#endif 543#endif 544 545/* Write a pari script to test things externally. */ 546#ifdef DOPRINT 547#include <stdio.h> 548 549#ifndef DOPRINT_SWIZZLE 550#define DOPRINT_SWIZZLE 0 551#endif 552 553#ifdef DOPRINT_LD80 554 555#define DOPRINT_START(xp) do { \ 556 uint64_t __lx; \ 557 uint16_t __hx; \ 558 \ 559 /* Hack to give more-problematic args. */ \ 560 EXTRACT_LDBL80_WORDS(__hx, __lx, *xp); \ 561 __lx ^= DOPRINT_SWIZZLE; \ 562 INSERT_LDBL80_WORDS(*xp, __hx, __lx); \ 563 printf("x = %.21Lg; ", (long double)*xp); \ 564} while (0) 565#define DOPRINT_END1(v) \ 566 printf("y = %.21Lg; z = 0; show(x, y, z);\n", (long double)(v)) 567#define DOPRINT_END2(hi, lo) \ 568 printf("y = %.21Lg; z = %.21Lg; show(x, y, z);\n", \ 569 (long double)(hi), (long double)(lo)) 570 571#elif defined(DOPRINT_D64) 572 573#define DOPRINT_START(xp) do { \ 574 uint32_t __hx, __lx; \ 575 \ 576 EXTRACT_WORDS(__hx, __lx, *xp); \ 577 __lx ^= DOPRINT_SWIZZLE; \ 578 INSERT_WORDS(*xp, __hx, __lx); \ 579 printf("x = %.21Lg; ", (long double)*xp); \ 580} while (0) 581#define DOPRINT_END1(v) \ 582 printf("y = %.21Lg; z = 0; show(x, y, z);\n", (long double)(v)) 583#define DOPRINT_END2(hi, lo) \ 584 printf("y = %.21Lg; z = %.21Lg; show(x, y, z);\n", \ 585 (long double)(hi), (long double)(lo)) 586 587#elif defined(DOPRINT_F32) 588 589#define DOPRINT_START(xp) do { \ 590 uint32_t __hx; \ 591 \ 592 GET_FLOAT_WORD(__hx, *xp); \ 593 __hx ^= DOPRINT_SWIZZLE; \ 594 SET_FLOAT_WORD(*xp, __hx); \ 595 printf("x = %.21Lg; ", (long double)*xp); \ 596} while (0) 597#define DOPRINT_END1(v) \ 598 printf("y = %.21Lg; z = 0; show(x, y, z);\n", (long double)(v)) 599#define DOPRINT_END2(hi, lo) \ 600 printf("y = %.21Lg; z = %.21Lg; show(x, y, z);\n", \ 601 (long double)(hi), (long double)(lo)) 602 603#else /* !DOPRINT_LD80 && !DOPRINT_D64 (LD128 only) */ 604 605#ifndef DOPRINT_SWIZZLE_HIGH 606#define DOPRINT_SWIZZLE_HIGH 0 607#endif 608 609#define DOPRINT_START(xp) do { \ 610 uint64_t __lx, __llx; \ 611 uint16_t __hx; \ 612 \ 613 EXTRACT_LDBL128_WORDS(__hx, __lx, __llx, *xp); \ 614 __llx ^= DOPRINT_SWIZZLE; \ 615 __lx ^= DOPRINT_SWIZZLE_HIGH; \ 616 INSERT_LDBL128_WORDS(*xp, __hx, __lx, __llx); \ 617 printf("x = %.36Lg; ", (long double)*xp); \ 618} while (0) 619#define DOPRINT_END1(v) \ 620 printf("y = %.36Lg; z = 0; show(x, y, z);\n", (long double)(v)) 621#define DOPRINT_END2(hi, lo) \ 622 printf("y = %.36Lg; z = %.36Lg; show(x, y, z);\n", \ 623 (long double)(hi), (long double)(lo)) 624 625#endif /* DOPRINT_LD80 */ 626 627#else /* !DOPRINT */ 628#define DOPRINT_START(xp) 629#define DOPRINT_END1(v) 630#define DOPRINT_END2(hi, lo) 631#endif /* DOPRINT */ 632 633#define RETURNP(x) do { \ 634 DOPRINT_END1(x); \ 635 RETURNF(x); \ 636} while (0) 637#define RETURNPI(x) do { \ 638 DOPRINT_END1(x); \ 639 RETURNI(x); \ 640} while (0) 641#define RETURN2P(x, y) do { \ 642 DOPRINT_END2((x), (y)); \ 643 RETURNF((x) + (y)); \ 644} while (0) 645#define RETURN2PI(x, y) do { \ 646 DOPRINT_END2((x), (y)); \ 647 RETURNI((x) + (y)); \ 648} while (0) 649#ifdef STRUCT_RETURN 650#define RETURNSP(rp) do { \ 651 if (!(rp)->lo_set) \ 652 RETURNP((rp)->hi); \ 653 RETURN2P((rp)->hi, (rp)->lo); \ 654} while (0) 655#define RETURNSPI(rp) do { \ 656 if (!(rp)->lo_set) \ 657 RETURNPI((rp)->hi); \ 658 RETURN2PI((rp)->hi, (rp)->lo); \ 659} while (0) 660#endif 661#define SUM2P(x, y) ({ \ 662 const __typeof (x) __x = (x); \ 663 const __typeof (y) __y = (y); \ 664 \ 665 DOPRINT_END2(__x, __y); \ 666 __x + __y; \ 667}) 668 669/* 670 * ieee style elementary functions 671 * 672 * We rename functions here to improve other sources' diffability 673 * against fdlibm. 674 */ 675#define __ieee754_sqrt sqrt 676#define __ieee754_acos acos 677#define __ieee754_acosh acosh 678#define __ieee754_log log 679#define __ieee754_log2 log2 680#define __ieee754_atanh atanh 681#define __ieee754_asin asin 682#define __ieee754_atan2 atan2 683#define __ieee754_exp exp 684#define __ieee754_cosh cosh 685#define __ieee754_fmod fmod 686#define __ieee754_pow pow 687#define __ieee754_lgamma lgamma 688#define __ieee754_gamma gamma 689#define __ieee754_lgamma_r lgamma_r 690#define __ieee754_gamma_r gamma_r 691#define __ieee754_log10 log10 692#define __ieee754_sinh sinh 693#define __ieee754_hypot hypot 694#define __ieee754_j0 j0 695#define __ieee754_j1 j1 696#define __ieee754_y0 y0 697#define __ieee754_y1 y1 698#define __ieee754_jn jn 699#define __ieee754_yn yn 700#define __ieee754_remainder remainder 701#define __ieee754_scalb scalb 702#define __ieee754_sqrtf sqrtf 703#define __ieee754_acosf acosf 704#define __ieee754_acoshf acoshf 705#define __ieee754_logf logf 706#define __ieee754_atanhf atanhf 707#define __ieee754_asinf asinf 708#define __ieee754_atan2f atan2f 709#define __ieee754_expf expf 710#define __ieee754_coshf coshf 711#define __ieee754_fmodf fmodf 712#define __ieee754_powf powf 713#define __ieee754_lgammaf lgammaf 714#define __ieee754_gammaf gammaf 715#define __ieee754_lgammaf_r lgammaf_r 716#define __ieee754_gammaf_r gammaf_r 717#define __ieee754_log10f log10f 718#define __ieee754_log2f log2f 719#define __ieee754_sinhf sinhf 720#define __ieee754_hypotf hypotf 721#define __ieee754_j0f j0f 722#define __ieee754_j1f j1f 723#define __ieee754_y0f y0f 724#define __ieee754_y1f y1f 725#define __ieee754_jnf jnf 726#define __ieee754_ynf ynf 727#define __ieee754_remainderf remainderf 728#define __ieee754_scalbf scalbf 729 730/* fdlibm kernel function */ 731int __kernel_rem_pio2(double*,double*,int,int,int); 732 733/* double precision kernel functions */ 734#ifndef INLINE_REM_PIO2 735int __ieee754_rem_pio2(double,double*); 736#endif 737double __kernel_sin(double,double,int); 738double __kernel_cos(double,double); 739double __kernel_tan(double,double,int); 740double __ldexp_exp(double,int); 741#ifdef _COMPLEX_H 742double complex __ldexp_cexp(double complex,int); 743#endif 744 745/* float precision kernel functions */ 746#ifndef INLINE_REM_PIO2F 747int __ieee754_rem_pio2f(float,double*); 748#endif 749#ifndef INLINE_KERNEL_SINDF 750float __kernel_sindf(double); 751#endif 752#ifndef INLINE_KERNEL_COSDF 753float __kernel_cosdf(double); 754#endif 755#ifndef INLINE_KERNEL_TANDF 756float __kernel_tandf(double,int); 757#endif 758float __ldexp_expf(float,int); 759#ifdef _COMPLEX_H 760float complex __ldexp_cexpf(float complex,int); 761#endif 762 763/* long double precision kernel functions */ 764long double __kernel_sinl(long double, long double, int); 765long double __kernel_cosl(long double, long double); 766long double __kernel_tanl(long double, long double, int); 767 768#endif /* !_MATH_PRIVATE_H_ */ 769