From 7ce331c01ce6eb7b3f5c715a38a24359da9c6ee2 Mon Sep 17 00:00:00 2001 From: Eric Andersen Date: Thu, 22 Nov 2001 14:04:29 +0000 Subject: Totally rework the math library, this time based on the MacOs X math library (which is itself based on the math lib from FreeBSD). -Erik --- libm/Makefile | 62 +- libm/README | 50 +- libm/ceilfloor.c | 179 ++ libm/double/Makefile | 114 - libm/double/README.txt | 5845 ---------------------------------------------- libm/double/acos.c | 58 - libm/double/acosh.c | 167 -- libm/double/airy.c | 965 -------- libm/double/arcdot.c | 110 - libm/double/asin.c | 324 --- libm/double/asinh.c | 165 -- libm/double/atan.c | 393 ---- libm/double/atanh.c | 156 -- libm/double/bdtr.c | 263 --- libm/double/bernum.c | 74 - libm/double/beta.c | 201 -- libm/double/btdtr.c | 64 - libm/double/cbrt.c | 142 -- libm/double/chbevl.c | 82 - libm/double/chdtr.c | 200 -- libm/double/cheby.c | 149 -- libm/double/clog.c | 1043 --------- libm/double/cmplx.c | 461 ---- libm/double/coil.c | 63 - libm/double/const.c | 252 -- libm/double/cosh.c | 83 - libm/double/cpmul.c | 104 - libm/double/dawsn.c | 392 ---- libm/double/dcalc.c | 1512 ------------ libm/double/dcalc.h | 77 - libm/double/dtestvec.c | 543 ----- libm/double/ei.c | 1062 --------- libm/double/eigens.c | 181 -- libm/double/ellie.c | 148 -- libm/double/ellik.c | 148 -- libm/double/ellpe.c | 195 -- libm/double/ellpj.c | 171 -- libm/double/ellpk.c | 234 -- libm/double/eltst.c | 37 - libm/double/euclid.c | 251 -- libm/double/exp.c | 203 -- libm/double/exp10.c | 223 -- libm/double/exp2.c | 183 -- libm/double/expn.c | 208 -- libm/double/fabs.c | 56 - libm/double/fac.c | 263 --- libm/double/fdtr.c | 237 -- libm/double/fftr.c | 237 -- libm/double/floor.c | 531 ----- libm/double/fltest.c | 272 --- libm/double/fltest2.c | 18 - libm/double/fltest3.c | 259 -- libm/double/fresnl.c | 515 ---- libm/double/gamma.c | 685 ------ libm/double/gdtr.c | 130 -- libm/double/gels.c | 232 -- libm/double/hyp2f1.c | 460 ---- libm/double/hyperg.c | 386 --- libm/double/i0.c | 397 ---- libm/double/i1.c | 402 ---- libm/double/igam.c | 210 -- libm/double/igami.c | 187 -- libm/double/incbet.c | 409 ---- libm/double/incbi.c | 313 --- libm/double/isnan.c | 237 -- libm/double/iv.c | 116 - libm/double/j0.c | 543 ----- libm/double/j1.c | 515 ---- libm/double/jn.c | 133 -- libm/double/jv.c | 884 ------- libm/double/k0.c | 333 --- libm/double/k1.c | 335 --- libm/double/kn.c | 255 -- libm/double/kolmogorov.c | 243 -- libm/double/levnsn.c | 82 - libm/double/log.c | 341 --- libm/double/log10.c | 250 -- libm/double/log2.c | 348 --- libm/double/lrand.c | 86 - libm/double/lsqrt.c | 85 - libm/double/ltstd.c | 469 ---- libm/double/minv.c | 61 - libm/double/mod2pi.c | 122 - libm/double/monot.c | 308 --- libm/double/mtherr.c | 102 - libm/double/mtransp.c | 61 - libm/double/mtst.c | 464 ---- libm/double/nbdtr.c | 222 -- libm/double/ndtr.c | 481 ---- libm/double/ndtri.c | 417 ---- libm/double/noncephes.c | 127 - libm/double/paranoia.c | 2156 ----------------- libm/double/pdtr.c | 184 -- libm/double/planck.c | 223 -- libm/double/polevl.c | 97 - libm/double/polmisc.c | 309 --- libm/double/polrt.c | 227 -- libm/double/polylog.c | 467 ---- libm/double/polyn.c | 471 ---- libm/double/polyr.c | 533 ----- libm/double/pow.c | 756 ------ libm/double/powi.c | 186 -- libm/double/psi.c | 201 -- libm/double/revers.c | 156 -- libm/double/rgamma.c | 209 -- libm/double/round.c | 79 - libm/double/setprec.c | 10 - libm/double/shichi.c | 599 ----- libm/double/sici.c | 675 ------ libm/double/simpsn.c | 81 - libm/double/simq.c | 180 -- libm/double/sin.c | 387 --- libm/double/sincos.c | 364 --- libm/double/sindg.c | 308 --- libm/double/sinh.c | 148 -- libm/double/spence.c | 205 -- libm/double/sqrt.c | 178 -- libm/double/stdtr.c | 225 -- libm/double/struve.c | 312 --- libm/double/tan.c | 304 --- libm/double/tandg.c | 267 --- libm/double/tanh.c | 141 -- libm/double/time-it.c | 38 - libm/double/unity.c | 138 -- libm/double/yn.c | 114 - libm/double/zeta.c | 189 -- libm/double/zetac.c | 599 ----- libm/e_acos.c | 111 + libm/e_acosh.c | 69 + libm/e_asin.c | 120 + libm/e_atan2.c | 130 ++ libm/e_atanh.c | 74 + libm/e_cosh.c | 93 + libm/e_exp.c | 167 ++ libm/e_fmod.c | 140 ++ libm/e_gamma.c | 34 + libm/e_gamma_r.c | 33 + libm/e_hypot.c | 128 + libm/e_j0.c | 487 ++++ libm/e_j1.c | 486 ++++ libm/e_jn.c | 281 +++ libm/e_lgamma.c | 34 + libm/e_lgamma_r.c | 316 +++ libm/e_log.c | 146 ++ libm/e_log10.c | 98 + libm/e_pow.c | 308 +++ libm/e_rem_pio2.c | 183 ++ libm/e_remainder.c | 80 + libm/e_scalb.c | 55 + libm/e_sinh.c | 86 + libm/e_sqrt.c | 453 ++++ libm/float/Makefile | 59 - libm/float/README.txt | 4721 ------------------------------------- libm/float/acoshf.c | 97 - libm/float/airyf.c | 377 --- libm/float/asinf.c | 186 -- libm/float/asinhf.c | 88 - libm/float/atanf.c | 190 -- libm/float/atanhf.c | 92 - libm/float/bdtrf.c | 247 -- libm/float/betaf.c | 122 - libm/float/cbrtf.c | 119 - libm/float/chbevlf.c | 86 - libm/float/chdtrf.c | 210 -- libm/float/clogf.c | 669 ------ libm/float/cmplxf.c | 407 ---- libm/float/constf.c | 20 - libm/float/coshf.c | 67 - libm/float/dawsnf.c | 168 -- libm/float/ellief.c | 115 - libm/float/ellikf.c | 113 - libm/float/ellpef.c | 105 - libm/float/ellpjf.c | 161 -- libm/float/ellpkf.c | 128 - libm/float/exp10f.c | 115 - libm/float/exp2f.c | 116 - libm/float/expf.c | 122 - libm/float/expnf.c | 207 -- libm/float/facf.c | 106 - libm/float/fdtrf.c | 214 -- libm/float/floorf.c | 526 ----- libm/float/fresnlf.c | 173 -- libm/float/gammaf.c | 423 ---- libm/float/gdtrf.c | 144 -- libm/float/hyp2f1f.c | 442 ---- libm/float/hypergf.c | 384 --- libm/float/i0f.c | 160 -- libm/float/i1f.c | 177 -- libm/float/igamf.c | 223 -- libm/float/igamif.c | 112 - libm/float/incbetf.c | 424 ---- libm/float/incbif.c | 197 -- libm/float/ivf.c | 114 - libm/float/j0f.c | 228 -- libm/float/j0tst.c | 43 - libm/float/j1f.c | 211 -- libm/float/jnf.c | 124 - libm/float/jvf.c | 848 ------- libm/float/k0f.c | 175 -- libm/float/k1f.c | 174 -- libm/float/knf.c | 252 -- libm/float/log10f.c | 129 - libm/float/log2f.c | 129 - libm/float/logf.c | 128 - libm/float/mtherr.c | 99 - libm/float/nantst.c | 54 - libm/float/nbdtrf.c | 141 -- libm/float/ndtrf.c | 281 --- libm/float/ndtrif.c | 186 -- libm/float/pdtrf.c | 188 -- libm/float/polevlf.c | 99 - libm/float/polynf.c | 520 ----- libm/float/powf.c | 338 --- libm/float/powif.c | 156 -- libm/float/powtst.c | 41 - libm/float/psif.c | 153 -- libm/float/rgammaf.c | 130 -- libm/float/setprec.c | 10 - libm/float/shichif.c | 212 -- libm/float/sicif.c | 279 --- libm/float/sindgf.c | 232 -- libm/float/sinf.c | 283 --- libm/float/sinhf.c | 87 - libm/float/spencef.c | 135 -- libm/float/sqrtf.c | 140 -- libm/float/stdtrf.c | 154 -- libm/float/struvef.c | 315 --- libm/float/tandgf.c | 206 -- libm/float/tanf.c | 192 -- libm/float/tanhf.c | 88 - libm/float/ynf.c | 120 - libm/float/zetacf.c | 266 --- libm/float/zetaf.c | 175 -- libm/fp_private.h | 112 + libm/fpmacros.c | 239 ++ libm/frexpldexp.c | 73 + libm/k_cos.c | 96 + libm/k_rem_pio2.c | 320 +++ libm/k_sin.c | 79 + libm/k_standard.c | 782 +++++++ libm/k_tan.c | 131 ++ libm/ldouble/Makefile | 122 - libm/ldouble/README.txt | 3502 --------------------------- libm/ldouble/acoshl.c | 167 -- libm/ldouble/arcdotl.c | 108 - libm/ldouble/asinhl.c | 156 -- libm/ldouble/asinl.c | 249 -- libm/ldouble/atanhl.c | 163 -- libm/ldouble/atanl.c | 376 --- libm/ldouble/bdtrl.c | 260 --- libm/ldouble/btdtrl.c | 68 - libm/ldouble/cbrtl.c | 143 -- libm/ldouble/chdtrl.c | 200 -- libm/ldouble/clogl.c | 720 ------ libm/ldouble/cmplxl.c | 461 ---- libm/ldouble/coshl.c | 89 - libm/ldouble/econst.c | 96 - libm/ldouble/ehead.h | 45 - libm/ldouble/elliel.c | 146 -- libm/ldouble/ellikl.c | 148 -- libm/ldouble/ellpel.c | 173 -- libm/ldouble/ellpjl.c | 164 -- libm/ldouble/ellpkl.c | 203 -- libm/ldouble/exp10l.c | 192 -- libm/ldouble/exp2l.c | 166 -- libm/ldouble/expl.c | 183 -- libm/ldouble/fdtrl.c | 237 -- libm/ldouble/floorl.c | 432 ---- libm/ldouble/flrtstl.c | 104 - libm/ldouble/fltestl.c | 265 --- libm/ldouble/gammal.c | 764 ------ libm/ldouble/gdtrl.c | 130 -- libm/ldouble/gelsl.c | 240 -- libm/ldouble/ieee.c | 4182 --------------------------------- libm/ldouble/igamil.c | 193 -- libm/ldouble/igaml.c | 220 -- libm/ldouble/incbetl.c | 406 ---- libm/ldouble/incbil.c | 305 --- libm/ldouble/isnanl.c | 186 -- libm/ldouble/j0l.c | 541 ----- libm/ldouble/j1l.c | 551 ----- libm/ldouble/jnl.c | 130 -- libm/ldouble/lcalc.c | 1484 ------------ libm/ldouble/lcalc.h | 79 - libm/ldouble/ldrand.c | 175 -- libm/ldouble/log10l.c | 319 --- libm/ldouble/log2l.c | 302 --- libm/ldouble/logl.c | 292 --- libm/ldouble/lparanoi.c | 2348 ------------------- libm/ldouble/monotl.c | 307 --- libm/ldouble/mtherr.c | 102 - libm/ldouble/mtstl.c | 521 ----- libm/ldouble/nantst.c | 61 - libm/ldouble/nbdtrl.c | 197 -- libm/ldouble/ndtril.c | 416 ---- libm/ldouble/ndtrl.c | 473 ---- libm/ldouble/pdtrl.c | 184 -- libm/ldouble/polevll.c | 182 -- libm/ldouble/powil.c | 164 -- libm/ldouble/powl.c | 739 ------ libm/ldouble/sinhl.c | 150 -- libm/ldouble/sinl.c | 342 --- libm/ldouble/sqrtl.c | 172 -- libm/ldouble/stdtrl.c | 225 -- libm/ldouble/tanhl.c | 129 - libm/ldouble/tanl.c | 279 --- libm/ldouble/testvect.c | 497 ---- libm/ldouble/unityl.c | 128 - libm/ldouble/wronkl.c | 67 - libm/ldouble/ynl.c | 113 - libm/logb.c | 104 + libm/math_private.h | 231 ++ libm/rndint.c | 627 +++++ libm/s_asinh.c | 65 + libm/s_atan.c | 139 ++ libm/s_cbrt.c | 93 + libm/s_ceil.c | 82 + libm/s_copysign.c | 40 + libm/s_cos.c | 82 + libm/s_erf.c | 314 +++ libm/s_expm1.c | 228 ++ libm/s_fabs.c | 35 + libm/s_finite.c | 35 + libm/s_floor.c | 83 + libm/s_frexp.c | 61 + libm/s_ilogb.c | 51 + libm/s_ldexp.c | 34 + libm/s_lib_version.c | 39 + libm/s_log1p.c | 173 ++ libm/s_logb.c | 44 + libm/s_matherr.c | 30 + libm/s_modf.c | 85 + libm/s_nextafter.c | 79 + libm/s_rint.c | 88 + libm/s_scalbn.c | 66 + libm/s_signgam.c | 3 + libm/s_significand.c | 34 + libm/s_sin.c | 82 + libm/s_tan.c | 76 + libm/s_tanh.c | 86 + libm/scalb.c | 87 + libm/sign.c | 58 + libm/w_acos.c | 43 + libm/w_acosh.c | 42 + libm/w_asin.c | 44 + libm/w_atan2.c | 42 + libm/w_atanh.c | 47 + libm/w_cabs.c | 20 + libm/w_cosh.c | 42 + libm/w_drem.c | 15 + libm/w_exp.c | 53 + libm/w_fmod.c | 43 + libm/w_gamma.c | 49 + libm/w_gamma_r.c | 46 + libm/w_hypot.c | 43 + libm/w_j0.c | 69 + libm/w_j1.c | 70 + libm/w_jn.c | 92 + libm/w_lgamma.c | 49 + libm/w_lgamma_r.c | 46 + libm/w_log.c | 43 + libm/w_log10.c | 46 + libm/w_pow.c | 61 + libm/w_remainder.c | 42 + libm/w_scalb.c | 62 + libm/w_sinh.c | 42 + libm/w_sqrt.c | 42 + 367 files changed, 10699 insertions(+), 91124 deletions(-) create mode 100644 libm/ceilfloor.c delete mode 100644 libm/double/Makefile delete mode 100644 libm/double/README.txt delete mode 100644 libm/double/acos.c delete mode 100644 libm/double/acosh.c delete mode 100644 libm/double/airy.c delete mode 100644 libm/double/arcdot.c delete mode 100644 libm/double/asin.c delete mode 100644 libm/double/asinh.c delete mode 100644 libm/double/atan.c delete mode 100644 libm/double/atanh.c delete mode 100644 libm/double/bdtr.c delete mode 100644 libm/double/bernum.c delete mode 100644 libm/double/beta.c delete mode 100644 libm/double/btdtr.c delete mode 100644 libm/double/cbrt.c delete mode 100644 libm/double/chbevl.c delete mode 100644 libm/double/chdtr.c delete mode 100644 libm/double/cheby.c delete mode 100644 libm/double/clog.c delete mode 100644 libm/double/cmplx.c delete mode 100644 libm/double/coil.c delete mode 100644 libm/double/const.c delete mode 100644 libm/double/cosh.c delete mode 100644 libm/double/cpmul.c delete mode 100644 libm/double/dawsn.c delete mode 100644 libm/double/dcalc.c delete mode 100644 libm/double/dcalc.h delete mode 100644 libm/double/dtestvec.c delete mode 100644 libm/double/ei.c delete mode 100644 libm/double/eigens.c delete mode 100644 libm/double/ellie.c delete mode 100644 libm/double/ellik.c delete mode 100644 libm/double/ellpe.c delete mode 100644 libm/double/ellpj.c delete mode 100644 libm/double/ellpk.c delete mode 100644 libm/double/eltst.c delete mode 100644 libm/double/euclid.c delete mode 100644 libm/double/exp.c delete mode 100644 libm/double/exp10.c delete mode 100644 libm/double/exp2.c delete mode 100644 libm/double/expn.c delete mode 100644 libm/double/fabs.c delete mode 100644 libm/double/fac.c delete mode 100644 libm/double/fdtr.c delete mode 100644 libm/double/fftr.c delete mode 100644 libm/double/floor.c delete mode 100644 libm/double/fltest.c delete mode 100644 libm/double/fltest2.c delete mode 100644 libm/double/fltest3.c delete mode 100644 libm/double/fresnl.c delete mode 100644 libm/double/gamma.c delete mode 100644 libm/double/gdtr.c delete mode 100644 libm/double/gels.c delete mode 100644 libm/double/hyp2f1.c delete mode 100644 libm/double/hyperg.c delete mode 100644 libm/double/i0.c delete mode 100644 libm/double/i1.c delete mode 100644 libm/double/igam.c delete mode 100644 libm/double/igami.c delete mode 100644 libm/double/incbet.c delete mode 100644 libm/double/incbi.c delete mode 100644 libm/double/isnan.c delete mode 100644 libm/double/iv.c delete mode 100644 libm/double/j0.c delete mode 100644 libm/double/j1.c delete mode 100644 libm/double/jn.c delete mode 100644 libm/double/jv.c delete mode 100644 libm/double/k0.c delete mode 100644 libm/double/k1.c delete mode 100644 libm/double/kn.c delete mode 100644 libm/double/kolmogorov.c delete mode 100644 libm/double/levnsn.c delete mode 100644 libm/double/log.c delete mode 100644 libm/double/log10.c delete mode 100644 libm/double/log2.c delete mode 100644 libm/double/lrand.c delete mode 100644 libm/double/lsqrt.c delete mode 100644 libm/double/ltstd.c delete mode 100644 libm/double/minv.c delete mode 100644 libm/double/mod2pi.c delete mode 100644 libm/double/monot.c delete mode 100644 libm/double/mtherr.c delete mode 100644 libm/double/mtransp.c delete mode 100644 libm/double/mtst.c delete mode 100644 libm/double/nbdtr.c delete mode 100644 libm/double/ndtr.c delete mode 100644 libm/double/ndtri.c delete mode 100644 libm/double/noncephes.c delete mode 100644 libm/double/paranoia.c delete mode 100644 libm/double/pdtr.c delete mode 100644 libm/double/planck.c delete mode 100644 libm/double/polevl.c delete mode 100644 libm/double/polmisc.c delete mode 100644 libm/double/polrt.c delete mode 100644 libm/double/polylog.c delete mode 100644 libm/double/polyn.c delete mode 100644 libm/double/polyr.c delete mode 100644 libm/double/pow.c delete mode 100644 libm/double/powi.c delete mode 100644 libm/double/psi.c delete mode 100644 libm/double/revers.c delete mode 100644 libm/double/rgamma.c delete mode 100644 libm/double/round.c delete mode 100644 libm/double/setprec.c delete mode 100644 libm/double/shichi.c delete mode 100644 libm/double/sici.c delete mode 100644 libm/double/simpsn.c delete mode 100644 libm/double/simq.c delete mode 100644 libm/double/sin.c delete mode 100644 libm/double/sincos.c delete mode 100644 libm/double/sindg.c delete mode 100644 libm/double/sinh.c delete mode 100644 libm/double/spence.c delete mode 100644 libm/double/sqrt.c delete mode 100644 libm/double/stdtr.c delete mode 100644 libm/double/struve.c delete mode 100644 libm/double/tan.c delete mode 100644 libm/double/tandg.c delete mode 100644 libm/double/tanh.c delete mode 100644 libm/double/time-it.c delete mode 100644 libm/double/unity.c delete mode 100644 libm/double/yn.c delete mode 100644 libm/double/zeta.c delete mode 100644 libm/double/zetac.c create mode 100644 libm/e_acos.c create mode 100644 libm/e_acosh.c create mode 100644 libm/e_asin.c create mode 100644 libm/e_atan2.c create mode 100644 libm/e_atanh.c create mode 100644 libm/e_cosh.c create mode 100644 libm/e_exp.c create mode 100644 libm/e_fmod.c create mode 100644 libm/e_gamma.c create mode 100644 libm/e_gamma_r.c create mode 100644 libm/e_hypot.c create mode 100644 libm/e_j0.c create mode 100644 libm/e_j1.c create mode 100644 libm/e_jn.c create mode 100644 libm/e_lgamma.c create mode 100644 libm/e_lgamma_r.c create mode 100644 libm/e_log.c create mode 100644 libm/e_log10.c create mode 100644 libm/e_pow.c create mode 100644 libm/e_rem_pio2.c create mode 100644 libm/e_remainder.c create mode 100644 libm/e_scalb.c create mode 100644 libm/e_sinh.c create mode 100644 libm/e_sqrt.c delete mode 100644 libm/float/Makefile delete mode 100644 libm/float/README.txt delete mode 100644 libm/float/acoshf.c delete mode 100644 libm/float/airyf.c delete mode 100644 libm/float/asinf.c delete mode 100644 libm/float/asinhf.c delete mode 100644 libm/float/atanf.c delete mode 100644 libm/float/atanhf.c delete mode 100644 libm/float/bdtrf.c delete mode 100644 libm/float/betaf.c delete mode 100644 libm/float/cbrtf.c delete mode 100644 libm/float/chbevlf.c delete mode 100644 libm/float/chdtrf.c delete mode 100644 libm/float/clogf.c delete mode 100644 libm/float/cmplxf.c delete mode 100644 libm/float/constf.c delete mode 100644 libm/float/coshf.c delete mode 100644 libm/float/dawsnf.c delete mode 100644 libm/float/ellief.c delete mode 100644 libm/float/ellikf.c delete mode 100644 libm/float/ellpef.c delete mode 100644 libm/float/ellpjf.c delete mode 100644 libm/float/ellpkf.c delete mode 100644 libm/float/exp10f.c delete mode 100644 libm/float/exp2f.c delete mode 100644 libm/float/expf.c delete mode 100644 libm/float/expnf.c delete mode 100644 libm/float/facf.c delete mode 100644 libm/float/fdtrf.c delete mode 100644 libm/float/floorf.c delete mode 100644 libm/float/fresnlf.c delete mode 100644 libm/float/gammaf.c delete mode 100644 libm/float/gdtrf.c delete mode 100644 libm/float/hyp2f1f.c delete mode 100644 libm/float/hypergf.c delete mode 100644 libm/float/i0f.c delete mode 100644 libm/float/i1f.c delete mode 100644 libm/float/igamf.c delete mode 100644 libm/float/igamif.c delete mode 100644 libm/float/incbetf.c delete mode 100644 libm/float/incbif.c delete mode 100644 libm/float/ivf.c delete mode 100644 libm/float/j0f.c delete mode 100644 libm/float/j0tst.c delete mode 100644 libm/float/j1f.c delete mode 100644 libm/float/jnf.c delete mode 100644 libm/float/jvf.c delete mode 100644 libm/float/k0f.c delete mode 100644 libm/float/k1f.c delete mode 100644 libm/float/knf.c delete mode 100644 libm/float/log10f.c delete mode 100644 libm/float/log2f.c delete mode 100644 libm/float/logf.c delete mode 100644 libm/float/mtherr.c delete mode 100644 libm/float/nantst.c delete mode 100644 libm/float/nbdtrf.c delete mode 100644 libm/float/ndtrf.c delete mode 100644 libm/float/ndtrif.c delete mode 100644 libm/float/pdtrf.c delete mode 100644 libm/float/polevlf.c delete mode 100644 libm/float/polynf.c delete mode 100644 libm/float/powf.c delete mode 100644 libm/float/powif.c delete mode 100644 libm/float/powtst.c delete mode 100644 libm/float/psif.c delete mode 100644 libm/float/rgammaf.c delete mode 100644 libm/float/setprec.c delete mode 100644 libm/float/shichif.c delete mode 100644 libm/float/sicif.c delete mode 100644 libm/float/sindgf.c delete mode 100644 libm/float/sinf.c delete mode 100644 libm/float/sinhf.c delete mode 100644 libm/float/spencef.c delete mode 100644 libm/float/sqrtf.c delete mode 100644 libm/float/stdtrf.c delete mode 100644 libm/float/struvef.c delete mode 100644 libm/float/tandgf.c delete mode 100644 libm/float/tanf.c delete mode 100644 libm/float/tanhf.c delete mode 100644 libm/float/ynf.c delete mode 100644 libm/float/zetacf.c delete mode 100644 libm/float/zetaf.c create mode 100644 libm/fp_private.h create mode 100644 libm/fpmacros.c create mode 100644 libm/frexpldexp.c create mode 100644 libm/k_cos.c create mode 100644 libm/k_rem_pio2.c create mode 100644 libm/k_sin.c create mode 100644 libm/k_standard.c create mode 100644 libm/k_tan.c delete mode 100644 libm/ldouble/Makefile delete mode 100644 libm/ldouble/README.txt delete mode 100644 libm/ldouble/acoshl.c delete mode 100644 libm/ldouble/arcdotl.c delete mode 100644 libm/ldouble/asinhl.c delete mode 100644 libm/ldouble/asinl.c delete mode 100644 libm/ldouble/atanhl.c delete mode 100644 libm/ldouble/atanl.c delete mode 100644 libm/ldouble/bdtrl.c delete mode 100644 libm/ldouble/btdtrl.c delete mode 100644 libm/ldouble/cbrtl.c delete mode 100644 libm/ldouble/chdtrl.c delete mode 100644 libm/ldouble/clogl.c delete mode 100644 libm/ldouble/cmplxl.c delete mode 100644 libm/ldouble/coshl.c delete mode 100644 libm/ldouble/econst.c delete mode 100644 libm/ldouble/ehead.h delete mode 100644 libm/ldouble/elliel.c delete mode 100644 libm/ldouble/ellikl.c delete mode 100644 libm/ldouble/ellpel.c delete mode 100644 libm/ldouble/ellpjl.c delete mode 100644 libm/ldouble/ellpkl.c delete mode 100644 libm/ldouble/exp10l.c delete mode 100644 libm/ldouble/exp2l.c delete mode 100644 libm/ldouble/expl.c delete mode 100644 libm/ldouble/fdtrl.c delete mode 100644 libm/ldouble/floorl.c delete mode 100644 libm/ldouble/flrtstl.c delete mode 100644 libm/ldouble/fltestl.c delete mode 100644 libm/ldouble/gammal.c delete mode 100644 libm/ldouble/gdtrl.c delete mode 100644 libm/ldouble/gelsl.c delete mode 100644 libm/ldouble/ieee.c delete mode 100644 libm/ldouble/igamil.c delete mode 100644 libm/ldouble/igaml.c delete mode 100644 libm/ldouble/incbetl.c delete mode 100644 libm/ldouble/incbil.c delete mode 100644 libm/ldouble/isnanl.c delete mode 100644 libm/ldouble/j0l.c delete mode 100644 libm/ldouble/j1l.c delete mode 100644 libm/ldouble/jnl.c delete mode 100644 libm/ldouble/lcalc.c delete mode 100644 libm/ldouble/lcalc.h delete mode 100644 libm/ldouble/ldrand.c delete mode 100644 libm/ldouble/log10l.c delete mode 100644 libm/ldouble/log2l.c delete mode 100644 libm/ldouble/logl.c delete mode 100644 libm/ldouble/lparanoi.c delete mode 100644 libm/ldouble/monotl.c delete mode 100644 libm/ldouble/mtherr.c delete mode 100644 libm/ldouble/mtstl.c delete mode 100644 libm/ldouble/nantst.c delete mode 100644 libm/ldouble/nbdtrl.c delete mode 100644 libm/ldouble/ndtril.c delete mode 100644 libm/ldouble/ndtrl.c delete mode 100644 libm/ldouble/pdtrl.c delete mode 100644 libm/ldouble/polevll.c delete mode 100644 libm/ldouble/powil.c delete mode 100644 libm/ldouble/powl.c delete mode 100644 libm/ldouble/sinhl.c delete mode 100644 libm/ldouble/sinl.c delete mode 100644 libm/ldouble/sqrtl.c delete mode 100644 libm/ldouble/stdtrl.c delete mode 100644 libm/ldouble/tanhl.c delete mode 100644 libm/ldouble/tanl.c delete mode 100644 libm/ldouble/testvect.c delete mode 100644 libm/ldouble/unityl.c delete mode 100644 libm/ldouble/wronkl.c delete mode 100644 libm/ldouble/ynl.c create mode 100644 libm/logb.c create mode 100644 libm/math_private.h create mode 100644 libm/rndint.c create mode 100644 libm/s_asinh.c create mode 100644 libm/s_atan.c create mode 100644 libm/s_cbrt.c create mode 100644 libm/s_ceil.c create mode 100644 libm/s_copysign.c create mode 100644 libm/s_cos.c create mode 100644 libm/s_erf.c create mode 100644 libm/s_expm1.c create mode 100644 libm/s_fabs.c create mode 100644 libm/s_finite.c create mode 100644 libm/s_floor.c create mode 100644 libm/s_frexp.c create mode 100644 libm/s_ilogb.c create mode 100644 libm/s_ldexp.c create mode 100644 libm/s_lib_version.c create mode 100644 libm/s_log1p.c create mode 100644 libm/s_logb.c create mode 100644 libm/s_matherr.c create mode 100644 libm/s_modf.c create mode 100644 libm/s_nextafter.c create mode 100644 libm/s_rint.c create mode 100644 libm/s_scalbn.c create mode 100644 libm/s_signgam.c create mode 100644 libm/s_significand.c create mode 100644 libm/s_sin.c create mode 100644 libm/s_tan.c create mode 100644 libm/s_tanh.c create mode 100644 libm/scalb.c create mode 100644 libm/sign.c create mode 100644 libm/w_acos.c create mode 100644 libm/w_acosh.c create mode 100644 libm/w_asin.c create mode 100644 libm/w_atan2.c create mode 100644 libm/w_atanh.c create mode 100644 libm/w_cabs.c create mode 100644 libm/w_cosh.c create mode 100644 libm/w_drem.c create mode 100644 libm/w_exp.c create mode 100644 libm/w_fmod.c create mode 100644 libm/w_gamma.c create mode 100644 libm/w_gamma_r.c create mode 100644 libm/w_hypot.c create mode 100644 libm/w_j0.c create mode 100644 libm/w_j1.c create mode 100644 libm/w_jn.c create mode 100644 libm/w_lgamma.c create mode 100644 libm/w_lgamma_r.c create mode 100644 libm/w_log.c create mode 100644 libm/w_log10.c create mode 100644 libm/w_pow.c create mode 100644 libm/w_remainder.c create mode 100644 libm/w_scalb.c create mode 100644 libm/w_sinh.c create mode 100644 libm/w_sqrt.c (limited to 'libm') diff --git a/libm/Makefile b/libm/Makefile index 5813ee9e3..b5ac92f80 100644 --- a/libm/Makefile +++ b/libm/Makefile @@ -25,31 +25,43 @@ include $(TOPDIR)Rules.mak LIBM=libm.a LIBM_SHARED=libm.so LIBM_SHARED_FULLNAME=libm-$(MAJOR_VERSION).$(MINOR_VERSION).so +TARGET_CC= $(TOPDIR)extra/gcc-uClibc/$(TARGET_ARCH)-uclibc-gcc +TARGET_CFLAGS+=-D_IEEE_LIBM -D_ISOC99_SOURCE -D_SVID_SOURCE -DIRS= -ifeq ($(strip $(HAS_LIBM_FLOAT)),true) - DIRS+=float +ifeq ($(strip $(DO_C89_ONLY)),true) +CSRC = FIXME +else +CSRC = e_acos.c e_acosh.c e_asin.c e_atan2.c e_atanh.c e_cosh.c\ + e_exp.c e_fmod.c e_gamma.c e_gamma_r.c e_hypot.c e_j0.c\ + e_j1.c e_jn.c e_lgamma.c e_lgamma_r.c e_log.c e_log10.c\ + e_pow.c e_remainder.c e_rem_pio2.c e_scalb.c e_sinh.c\ + e_sqrt.c k_cos.c k_rem_pio2.c k_sin.c k_standard.c k_tan.c\ + s_asinh.c s_atan.c s_cbrt.c s_ceil.c s_copysign.c s_cos.c\ + s_erf.c s_expm1.c s_fabs.c s_finite.c s_floor.c s_frexp.c\ + s_ilogb.c s_ldexp.c s_lib_version.c s_log1p.c s_logb.c\ + s_matherr.c s_modf.c s_nextafter.c s_rint.c s_scalbn.c\ + s_signgam.c s_significand.c s_sin.c s_tan.c s_tanh.c\ + w_acos.c w_acosh.c w_asin.c w_atan2.c w_atanh.c w_cabs.c\ + w_cosh.c w_drem.c w_exp.c w_fmod.c w_gamma.c w_gamma_r.c\ + w_hypot.c w_j0.c w_j1.c w_jn.c w_lgamma.c w_lgamma_r.c\ + w_log.c w_log10.c w_pow.c w_remainder.c w_scalb.c w_sinh.c\ + w_sqrt.c ceilfloor.c fpmacros.c frexpldexp.c logb.c rndint.c\ + scalb.c sign.c endif -ifeq ($(strip $(HAS_LIBM_DOUBLE)),true) - DIRS+=double -endif -ifeq ($(strip $(HAS_LIBM_LONG_DOUBLE)),true) - DIRS+=ldouble -endif -ALL_SUBDIRS = float double ldouble +COBJS=$(patsubst %.c,%.o, $(CSRC)) +OBJS=$(COBJS) + -all: $(LIBM) -$(LIBM): subdirs +all: $(OBJS) $(LIBM) + +$(LIBM): ar-target @if [ -f $(LIBM) ] ; then \ install -d $(TOPDIR)lib; \ rm -f $(TOPDIR)lib/$(LIBM); \ install -m 644 $(LIBM) $(TOPDIR)lib; \ fi; -tags: - ctags -R - shared: all if [ -f $(LIBM) ] ; then \ $(TARGET_CC) $(TARGET_LDFLAGS) -nostdlib -shared -o $(LIBM_SHARED_FULLNAME) \ @@ -61,18 +73,18 @@ shared: all (cd $(TOPDIR)lib; ln -sf $(LIBM_SHARED_FULLNAME) $(LIBM_SHARED).$(MAJOR_VERSION)); \ fi; -subdirs: $(patsubst %, _dir_%, $(DIRS)) -subdirs_clean: $(patsubst %, _dirclean_%, $(ALL_SUBDIRS)) - -$(patsubst %, _dir_%, $(DIRS)) : dummy - $(MAKE) -C $(patsubst _dir_%, %, $@) +ar-target: $(OBJS) + $(AR) $(ARFLAGS) $(LIBM) $(OBJS) -$(patsubst %, _dirclean_%, $(ALL_SUBDIRS)) : dummy - $(MAKE) -C $(patsubst _dirclean_%, %, $@) clean +$(COBJS): %.o : %.c + $(TARGET_CC) $(TARGET_CFLAGS) -c $< -o $@ + $(STRIPTOOL) -x -R .note -R .comment $*.o -clean: subdirs_clean - rm -f *.[oa] *~ core $(LIBM_SHARED)* $(LIBM_SHARED_FULLNAME)* +$(OBJ): Makefile -.PHONY: dummy +tags: + ctags -R +clean: + rm -f *.[oa] *~ core $(LIBM_SHARED)* $(LIBM_SHARED_FULLNAME)* diff --git a/libm/README b/libm/README index 023e46846..c275d1b9a 100644 --- a/libm/README +++ b/libm/README @@ -1,42 +1,16 @@ -The actual routines included in this math library are derived almost -exclusively from the Cephes Mathematical Library, which "is copyrighted by the -author [and] may be used freely but ... comes with no support or guarantee" +The routines included in this math library are derived from the +math library for Apple's MacOS X/Darwin math library, which was +itself swiped from FreeBSD. The original copyright information +is as follows: -It has been ported to fit into uClibc and generally behave -by Erik Andersen , - 5 May, 2001 + Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. --------------------------------------------------- + Developed at SunPro, a Sun Microsystems, Inc. business. + Permission to use, copy, modify, and distribute this + software is freely granted, provided that this notice + is preserved. - Some software in this archive may be from the book _Methods and -Programs for Mathematical Functions_ (Prentice-Hall, 1989) or -from the Cephes Mathematical Library, a commercial product. In -either event, it is copyrighted by the author. What you see here -may be used freely but it comes with no support or guarantee. +It has been ported to work with uClibc and generally behave +by Erik Andersen + 22 May, 2001 - The two known misprints in the book are repaired here in the -source listings for the gamma function and the incomplete beta -integral. - - - Stephen L. Moshier - moshier@world.std.com - --------------------------------------------------- - -19 November 1992 - -ZIP archive constructed and index compiled. - -To reconstruct the original directory structure, use the -d switch: - - C:\CEPHES>pkunzip -d cephes - -This archive includes all the programs in the /netlib/cephes directory -on research.att.com as of 17 Nov 92. The file "index" will tell you in -what directory and file each function can be found. If there is -something else mentioned in cephes.doc that you need, you can check -research.att.com to see whether it has been added. Failing that, you -can contact Stephen Moshier. - - Jim Van Zandt diff --git a/libm/ceilfloor.c b/libm/ceilfloor.c new file mode 100644 index 000000000..9607435c3 --- /dev/null +++ b/libm/ceilfloor.c @@ -0,0 +1,179 @@ +#if defined(__ppc__) +/******************************************************************************* +* * +* File ceilfloor.c, * +* Function ceil(x) and floor(x), * +* Implementation of ceil and floor for the PowerPC. * +* * +* Copyright © 1991 Apple Computer, Inc. All rights reserved. * +* * +* Written by Ali Sazegari, started on November 1991, * +* * +* based on math.h, library code for Macintoshes with a 68881/68882 * +* by Jim Thomas. * +* * +* W A R N I N G: This routine expects a 64 bit double model. * +* * +* December 03 1992: first rs6000 port. * +* July 14 1993: comment changes and addition of #pragma fenv_access. * +* May 06 1997: port of the ibm/taligent ceil and floor routines. * +* April 11 2001: first port to os x using gcc. * +* June 13 2001: replaced __setflm with in-line assembly * +* * +*******************************************************************************/ + +#if !defined(__ppc__) +#define asm(x) +#endif + +static const double twoTo52 = 4503599627370496.0; +static const unsigned long signMask = 0x80000000ul; + +typedef union + { + struct { +#if defined(__BIG_ENDIAN__) + unsigned long int hi; + unsigned long int lo; +#else + unsigned long int lo; + unsigned long int hi; +#endif + } words; + double dbl; + } DblInHex; + +/******************************************************************************* +* Functions needed for the computation. * +*******************************************************************************/ + +/******************************************************************************* +* Ceil(x) returns the smallest integer not less than x. * +*******************************************************************************/ + +double ceil ( double x ) + { + DblInHex xInHex,OldEnvironment; + register double y; + register unsigned long int xhi; + register int target; + + xInHex.dbl = x; + xhi = xInHex.words.hi & 0x7fffffffUL; // xhi is the high half of |x| + target = ( xInHex.words.hi < signMask ); + + if ( xhi < 0x43300000ul ) +/******************************************************************************* +* Is |x| < 2.0^52? * +*******************************************************************************/ + { + if ( xhi < 0x3ff00000ul ) +/******************************************************************************* +* Is |x| < 1.0? * +*******************************************************************************/ + { + if ( ( xhi | xInHex.words.lo ) == 0ul ) // zero x is exact case + return ( x ); + else + { // inexact case + asm ("mffs %0" : "=f" (OldEnvironment.dbl)); + OldEnvironment.words.lo |= 0x02000000ul; + asm ("mtfsf 255,%0" : /*NULLOUT*/ : /*IN*/ "f" ( OldEnvironment.dbl )); + if ( target ) + return ( 1.0 ); + else + return ( -0.0 ); + } + } +/******************************************************************************* +* Is 1.0 < |x| < 2.0^52? * +*******************************************************************************/ + if ( target ) + { + y = ( x + twoTo52 ) - twoTo52; // round at binary pt. + if ( y < x ) + return ( y + 1.0 ); + else + return ( y ); + } + + else + { + y = ( x - twoTo52 ) + twoTo52; // round at binary pt. + if ( y < x ) + return ( y + 1.0 ); + else + return ( y ); + } + } +/******************************************************************************* +* |x| >= 2.0^52 or x is a NaN. * +*******************************************************************************/ + return ( x ); + } + +/******************************************************************************* +* Floor(x) returns the largest integer not greater than x. * +*******************************************************************************/ + +double floor ( double x ) + { + DblInHex xInHex,OldEnvironment; + register double y; + register unsigned long int xhi; + register long int target; + + xInHex.dbl = x; + xhi = xInHex.words.hi & 0x7fffffffUL; // xhi is the high half of |x| + target = ( xInHex.words.hi < signMask ); + + if ( xhi < 0x43300000ul ) +/******************************************************************************* +* Is |x| < 2.0^52? * +*******************************************************************************/ + { + if ( xhi < 0x3ff00000ul ) +/******************************************************************************* +* Is |x| < 1.0? * +*******************************************************************************/ + { + if ( ( xhi | xInHex.words.lo ) == 0ul ) // zero x is exact case + return ( x ); + else + { // inexact case + asm ("mffs %0" : "=f" (OldEnvironment.dbl)); + OldEnvironment.words.lo |= 0x02000000ul; + asm ("mtfsf 255,%0" : /*NULLOUT*/ : /*IN*/ "f" ( OldEnvironment.dbl )); + if ( target ) + return ( 0.0 ); + else + return ( -1.0 ); + } + } +/******************************************************************************* +* Is 1.0 < |x| < 2.0^52? * +*******************************************************************************/ + if ( target ) + { + y = ( x + twoTo52 ) - twoTo52; // round at binary pt. + if ( y > x ) + return ( y - 1.0 ); + else + return ( y ); + } + + else + { + y = ( x - twoTo52 ) + twoTo52; // round at binary pt. + if ( y > x ) + return ( y - 1.0 ); + else + return ( y ); + } + } +/******************************************************************************* +* |x| >= 2.0^52 or x is a NaN. * +*******************************************************************************/ + return ( x ); + } +#endif /* __ppc__ */ diff --git a/libm/double/Makefile b/libm/double/Makefile deleted file mode 100644 index a53b44d2e..000000000 --- a/libm/double/Makefile +++ /dev/null @@ -1,114 +0,0 @@ -# Makefile for uClibc's math library -# Copyright (C) 2001 by Lineo, inc. -# -# This math library is derived primarily from the Cephes Math Library, -# copyright by Stephen L. Moshier -# -# This program is free software; you can redistribute it and/or modify it under -# the terms of the GNU Library General Public License as published by the Free -# Software Foundation; either version 2 of the License, or (at your option) any -# later version. -# -# This program is distributed in the hope that it will be useful, but WITHOUT -# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS -# FOR A PARTICULAR PURPOSE. See the GNU Library General Public License for more -# details. -# -# You should have received a copy of the GNU Library General Public License -# along with this program; if not, write to the Free Software Foundation, Inc., -# 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA -# - -TOPDIR=../../ -include $(TOPDIR)Rules.mak - -LIBM=../libm.a -TARGET_CC= $(TOPDIR)extra/gcc-uClibc/$(TARGET_ARCH)-uclibc-gcc - -CSRC=acosh.c airy.c asin.c asinh.c atan.c atanh.c bdtr.c beta.c \ - btdtr.c cbrt.c chbevl.c chdtr.c clog.c cmplx.c const.c \ - cosh.c dawsn.c ei.c ellie.c ellik.c ellpe.c ellpj.c ellpk.c \ - exp.c exp10.c exp2.c expn.c fabs.c fac.c fdtr.c \ - fresnl.c gamma.c gdtr.c hyp2f1.c hyperg.c i0.c i1.c igami.c incbet.c \ - incbi.c igam.c isnan.c iv.c j0.c j1.c jn.c jv.c k0.c k1.c kn.c kolmogorov.c \ - log.c log2.c log10.c lrand.c nbdtr.c ndtr.c ndtri.c pdtr.c planck.c \ - polevl.c polmisc.c polylog.c polyn.c pow.c powi.c psi.c rgamma.c round.c \ - shichi.c sici.c sin.c sindg.c sinh.c spence.c stdtr.c struve.c \ - tan.c tandg.c tanh.c unity.c yn.c zeta.c zetac.c \ - sqrt.c floor.c setprec.c mtherr.c noncephes.c - -COBJS=$(patsubst %.c,%.o, $(CSRC)) - - -OBJS=$(COBJS) - -all: $(OBJS) $(LIBM) - -$(LIBM): ar-target - -ar-target: $(OBJS) - $(AR) $(ARFLAGS) $(LIBM) $(OBJS) - -$(COBJS): %.o : %.c - $(TARGET_CC) $(TARGET_CFLAGS) -c $< -o $@ - $(STRIPTOOL) -x -R .note -R .comment $*.o - -$(OBJ): Makefile - -clean: - rm -f *.[oa] *~ core - - - -#----------------------------------------- - -#all: libmd.a mtst dtestvec monot dcalc paranoia - -time-it: time-it.o - $(TARGET_CC) -o time-it time-it.o - -time-it.o: time-it.c - $(TARGET_CC) -O2 -c time-it.c - -dcalc: dcalc.o libmd.a - $(TARGET_CC) -o dcalc dcalc.o libmd.a - -mtst: mtst.o libmd.a - $(TARGET_CC) -v -o mtst mtst.o libmd.a - -mtst.o: mtst.c - $(TARGET_CC) -O2 -Wall -c mtst.c - -dtestvec: dtestvec.o libmd.a - $(TARGET_CC) -o dtestvec dtestvec.o libmd.a - -dtestvec.o: dtestvec.c - $(TARGET_CC) -g -c dtestvec.c - -monot: monot.o libmd.a - $(TARGET_CC) -o monot monot.o libmd.a - -monot.o: monot.c - $(TARGET_CC) -g -c monot.c - -paranoia: paranoia.o setprec.o libmd.a - $(TARGET_CC) -o paranoia paranoia.o setprec.o libmd.a - -paranoia.o: paranoia.c - $(TARGET_CC) $(TARGET_CFLAGS) -Wno-implicit -c paranoia.c - -libmd.a: $(OBJS) $(INCS) - $(AR) rv libmd.a $(OBJS) - -#clean: -# rm -f *.o -# rm -f mtst -# rm -f paranoia -# rm -f dcalc -# rm -f dtestvec -# rm -f monot -# rm -f libmd.a -# rm -f time-it -# rm -f dtestvec - - diff --git a/libm/double/README.txt b/libm/double/README.txt deleted file mode 100644 index f2cb6c3dc..000000000 --- a/libm/double/README.txt +++ /dev/null @@ -1,5845 +0,0 @@ -/* acosh.c - * - * Inverse hyperbolic cosine - * - * - * - * SYNOPSIS: - * - * double x, y, acosh(); - * - * y = acosh( x ); - * - * - * - * DESCRIPTION: - * - * Returns inverse hyperbolic cosine of argument. - * - * If 1 <= x < 1.5, a rational approximation - * - * sqrt(z) * P(z)/Q(z) - * - * where z = x-1, is used. Otherwise, - * - * acosh(x) = log( x + sqrt( (x-1)(x+1) ). - * - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC 1,3 30000 4.2e-17 1.1e-17 - * IEEE 1,3 30000 4.6e-16 8.7e-17 - * - * - * ERROR MESSAGES: - * - * message condition value returned - * acosh domain |x| < 1 NAN - * - */ - -/* airy.c - * - * Airy function - * - * - * - * SYNOPSIS: - * - * double x, ai, aip, bi, bip; - * int airy(); - * - * airy( x, _&ai, _&aip, _&bi, _&bip ); - * - * - * - * DESCRIPTION: - * - * Solution of the differential equation - * - * y"(x) = xy. - * - * The function returns the two independent solutions Ai, Bi - * and their first derivatives Ai'(x), Bi'(x). - * - * Evaluation is by power series summation for small x, - * by rational minimax approximations for large x. - * - * - * - * ACCURACY: - * Error criterion is absolute when function <= 1, relative - * when function > 1, except * denotes relative error criterion. - * For large negative x, the absolute error increases as x^1.5. - * For large positive x, the relative error increases as x^1.5. - * - * Arithmetic domain function # trials peak rms - * IEEE -10, 0 Ai 10000 1.6e-15 2.7e-16 - * IEEE 0, 10 Ai 10000 2.3e-14* 1.8e-15* - * IEEE -10, 0 Ai' 10000 4.6e-15 7.6e-16 - * IEEE 0, 10 Ai' 10000 1.8e-14* 1.5e-15* - * IEEE -10, 10 Bi 30000 4.2e-15 5.3e-16 - * IEEE -10, 10 Bi' 30000 4.9e-15 7.3e-16 - * DEC -10, 0 Ai 5000 1.7e-16 2.8e-17 - * DEC 0, 10 Ai 5000 2.1e-15* 1.7e-16* - * DEC -10, 0 Ai' 5000 4.7e-16 7.8e-17 - * DEC 0, 10 Ai' 12000 1.8e-15* 1.5e-16* - * DEC -10, 10 Bi 10000 5.5e-16 6.8e-17 - * DEC -10, 10 Bi' 7000 5.3e-16 8.7e-17 - * - */ - -/* asin.c - * - * Inverse circular sine - * - * - * - * SYNOPSIS: - * - * double x, y, asin(); - * - * y = asin( x ); - * - * - * - * DESCRIPTION: - * - * Returns radian angle between -pi/2 and +pi/2 whose sine is x. - * - * A rational function of the form x + x**3 P(x**2)/Q(x**2) - * is used for |x| in the interval [0, 0.5]. If |x| > 0.5 it is - * transformed by the identity - * - * asin(x) = pi/2 - 2 asin( sqrt( (1-x)/2 ) ). - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC -1, 1 40000 2.6e-17 7.1e-18 - * IEEE -1, 1 10^6 1.9e-16 5.4e-17 - * - * - * ERROR MESSAGES: - * - * message condition value returned - * asin domain |x| > 1 NAN - * - */ - /* acos() - * - * Inverse circular cosine - * - * - * - * SYNOPSIS: - * - * double x, y, acos(); - * - * y = acos( x ); - * - * - * - * DESCRIPTION: - * - * Returns radian angle between 0 and pi whose cosine - * is x. - * - * Analytically, acos(x) = pi/2 - asin(x). However if |x| is - * near 1, there is cancellation error in subtracting asin(x) - * from pi/2. Hence if x < -0.5, - * - * acos(x) = pi - 2.0 * asin( sqrt((1+x)/2) ); - * - * or if x > +0.5, - * - * acos(x) = 2.0 * asin( sqrt((1-x)/2) ). - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC -1, 1 50000 3.3e-17 8.2e-18 - * IEEE -1, 1 10^6 2.2e-16 6.5e-17 - * - * - * ERROR MESSAGES: - * - * message condition value returned - * asin domain |x| > 1 NAN - */ - -/* asinh.c - * - * Inverse hyperbolic sine - * - * - * - * SYNOPSIS: - * - * double x, y, asinh(); - * - * y = asinh( x ); - * - * - * - * DESCRIPTION: - * - * Returns inverse hyperbolic sine of argument. - * - * If |x| < 0.5, the function is approximated by a rational - * form x + x**3 P(x)/Q(x). Otherwise, - * - * asinh(x) = log( x + sqrt(1 + x*x) ). - * - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC -3,3 75000 4.6e-17 1.1e-17 - * IEEE -1,1 30000 3.7e-16 7.8e-17 - * IEEE 1,3 30000 2.5e-16 6.7e-17 - * - */ - -/* atan.c - * - * Inverse circular tangent - * (arctangent) - * - * - * - * SYNOPSIS: - * - * double x, y, atan(); - * - * y = atan( x ); - * - * - * - * DESCRIPTION: - * - * Returns radian angle between -pi/2 and +pi/2 whose tangent - * is x. - * - * Range reduction is from three intervals into the interval - * from zero to 0.66. The approximant uses a rational - * function of degree 4/5 of the form x + x**3 P(x)/Q(x). - * - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC -10, 10 50000 2.4e-17 8.3e-18 - * IEEE -10, 10 10^6 1.8e-16 5.0e-17 - * - */ - /* atan2() - * - * Quadrant correct inverse circular tangent - * - * - * - * SYNOPSIS: - * - * double x, y, z, atan2(); - * - * z = atan2( y, x ); - * - * - * - * DESCRIPTION: - * - * Returns radian angle whose tangent is y/x. - * Define compile time symbol ANSIC = 1 for ANSI standard, - * range -PI < z <= +PI, args (y,x); else ANSIC = 0 for range - * 0 to 2PI, args (x,y). - * - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * IEEE -10, 10 10^6 2.5e-16 6.9e-17 - * See atan.c. - * - */ - -/* atanh.c - * - * Inverse hyperbolic tangent - * - * - * - * SYNOPSIS: - * - * double x, y, atanh(); - * - * y = atanh( x ); - * - * - * - * DESCRIPTION: - * - * Returns inverse hyperbolic tangent of argument in the range - * MINLOG to MAXLOG. - * - * If |x| < 0.5, the rational form x + x**3 P(x)/Q(x) is - * employed. Otherwise, - * atanh(x) = 0.5 * log( (1+x)/(1-x) ). - * - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC -1,1 50000 2.4e-17 6.4e-18 - * IEEE -1,1 30000 1.9e-16 5.2e-17 - * - */ - -/* bdtr.c - * - * Binomial distribution - * - * - * - * SYNOPSIS: - * - * int k, n; - * double p, y, bdtr(); - * - * y = bdtr( k, n, p ); - * - * DESCRIPTION: - * - * Returns the sum of the terms 0 through k of the Binomial - * probability density: - * - * k - * -- ( n ) j n-j - * > ( ) p (1-p) - * -- ( j ) - * j=0 - * - * The terms are not summed directly; instead the incomplete - * beta integral is employed, according to the formula - * - * y = bdtr( k, n, p ) = incbet( n-k, k+1, 1-p ). - * - * The arguments must be positive, with p ranging from 0 to 1. - * - * ACCURACY: - * - * Tested at random points (a,b,p), with p between 0 and 1. - * - * a,b Relative error: - * arithmetic domain # trials peak rms - * For p between 0.001 and 1: - * IEEE 0,100 100000 4.3e-15 2.6e-16 - * See also incbet.c. - * - * ERROR MESSAGES: - * - * message condition value returned - * bdtr domain k < 0 0.0 - * n < k - * x < 0, x > 1 - */ - /* bdtrc() - * - * Complemented binomial distribution - * - * - * - * SYNOPSIS: - * - * int k, n; - * double p, y, bdtrc(); - * - * y = bdtrc( k, n, p ); - * - * DESCRIPTION: - * - * Returns the sum of the terms k+1 through n of the Binomial - * probability density: - * - * n - * -- ( n ) j n-j - * > ( ) p (1-p) - * -- ( j ) - * j=k+1 - * - * The terms are not summed directly; instead the incomplete - * beta integral is employed, according to the formula - * - * y = bdtrc( k, n, p ) = incbet( k+1, n-k, p ). - * - * The arguments must be positive, with p ranging from 0 to 1. - * - * ACCURACY: - * - * Tested at random points (a,b,p). - * - * a,b Relative error: - * arithmetic domain # trials peak rms - * For p between 0.001 and 1: - * IEEE 0,100 100000 6.7e-15 8.2e-16 - * For p between 0 and .001: - * IEEE 0,100 100000 1.5e-13 2.7e-15 - * - * ERROR MESSAGES: - * - * message condition value returned - * bdtrc domain x<0, x>1, n 1 - */ - -/* beta.c - * - * Beta function - * - * - * - * SYNOPSIS: - * - * double a, b, y, beta(); - * - * y = beta( a, b ); - * - * - * - * DESCRIPTION: - * - * - - - * | (a) | (b) - * beta( a, b ) = -----------. - * - - * | (a+b) - * - * For large arguments the logarithm of the function is - * evaluated using lgam(), then exponentiated. - * - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC 0,30 1700 7.7e-15 1.5e-15 - * IEEE 0,30 30000 8.1e-14 1.1e-14 - * - * ERROR MESSAGES: - * - * message condition value returned - * beta overflow log(beta) > MAXLOG 0.0 - * a or b <0 integer 0.0 - * - */ - -/* btdtr.c - * - * Beta distribution - * - * - * - * SYNOPSIS: - * - * double a, b, x, y, btdtr(); - * - * y = btdtr( a, b, x ); - * - * - * - * DESCRIPTION: - * - * Returns the area from zero to x under the beta density - * function: - * - * - * x - * - - - * | (a+b) | | a-1 b-1 - * P(x) = ---------- | t (1-t) dt - * - - | | - * | (a) | (b) - - * 0 - * - * - * This function is identical to the incomplete beta - * integral function incbet(a, b, x). - * - * The complemented function is - * - * 1 - P(1-x) = incbet( b, a, x ); - * - * - * ACCURACY: - * - * See incbet.c. - * - */ - -/* cbrt.c - * - * Cube root - * - * - * - * SYNOPSIS: - * - * double x, y, cbrt(); - * - * y = cbrt( x ); - * - * - * - * DESCRIPTION: - * - * Returns the cube root of the argument, which may be negative. - * - * Range reduction involves determining the power of 2 of - * the argument. A polynomial of degree 2 applied to the - * mantissa, and multiplication by the cube root of 1, 2, or 4 - * approximates the root to within about 0.1%. Then Newton's - * iteration is used three times to converge to an accurate - * result. - * - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC -10,10 200000 1.8e-17 6.2e-18 - * IEEE 0,1e308 30000 1.5e-16 5.0e-17 - * - */ - -/* chbevl.c - * - * Evaluate Chebyshev series - * - * - * - * SYNOPSIS: - * - * int N; - * double x, y, coef[N], chebevl(); - * - * y = chbevl( x, coef, N ); - * - * - * - * DESCRIPTION: - * - * Evaluates the series - * - * N-1 - * - ' - * y = > coef[i] T (x/2) - * - i - * i=0 - * - * of Chebyshev polynomials Ti at argument x/2. - * - * Coefficients are stored in reverse order, i.e. the zero - * order term is last in the array. Note N is the number of - * coefficients, not the order. - * - * If coefficients are for the interval a to b, x must - * have been transformed to x -> 2(2x - b - a)/(b-a) before - * entering the routine. This maps x from (a, b) to (-1, 1), - * over which the Chebyshev polynomials are defined. - * - * If the coefficients are for the inverted interval, in - * which (a, b) is mapped to (1/b, 1/a), the transformation - * required is x -> 2(2ab/x - b - a)/(b-a). If b is infinity, - * this becomes x -> 4a/x - 1. - * - * - * - * SPEED: - * - * Taking advantage of the recurrence properties of the - * Chebyshev polynomials, the routine requires one more - * addition per loop than evaluating a nested polynomial of - * the same degree. - * - */ - -/* chdtr.c - * - * Chi-square distribution - * - * - * - * SYNOPSIS: - * - * double df, x, y, chdtr(); - * - * y = chdtr( df, x ); - * - * - * - * DESCRIPTION: - * - * Returns the area under the left hand tail (from 0 to x) - * of the Chi square probability density function with - * v degrees of freedom. - * - * - * inf. - * - - * 1 | | v/2-1 -t/2 - * P( x | v ) = ----------- | t e dt - * v/2 - | | - * 2 | (v/2) - - * x - * - * where x is the Chi-square variable. - * - * The incomplete gamma integral is used, according to the - * formula - * - * y = chdtr( v, x ) = igam( v/2.0, x/2.0 ). - * - * - * The arguments must both be positive. - * - * - * - * ACCURACY: - * - * See igam(). - * - * ERROR MESSAGES: - * - * message condition value returned - * chdtr domain x < 0 or v < 1 0.0 - */ - /* chdtrc() - * - * Complemented Chi-square distribution - * - * - * - * SYNOPSIS: - * - * double v, x, y, chdtrc(); - * - * y = chdtrc( v, x ); - * - * - * - * DESCRIPTION: - * - * Returns the area under the right hand tail (from x to - * infinity) of the Chi square probability density function - * with v degrees of freedom: - * - * - * inf. - * - - * 1 | | v/2-1 -t/2 - * P( x | v ) = ----------- | t e dt - * v/2 - | | - * 2 | (v/2) - - * x - * - * where x is the Chi-square variable. - * - * The incomplete gamma integral is used, according to the - * formula - * - * y = chdtr( v, x ) = igamc( v/2.0, x/2.0 ). - * - * - * The arguments must both be positive. - * - * - * - * ACCURACY: - * - * See igamc(). - * - * ERROR MESSAGES: - * - * message condition value returned - * chdtrc domain x < 0 or v < 1 0.0 - */ - /* chdtri() - * - * Inverse of complemented Chi-square distribution - * - * - * - * SYNOPSIS: - * - * double df, x, y, chdtri(); - * - * x = chdtri( df, y ); - * - * - * - * - * DESCRIPTION: - * - * Finds the Chi-square argument x such that the integral - * from x to infinity of the Chi-square density is equal - * to the given cumulative probability y. - * - * This is accomplished using the inverse gamma integral - * function and the relation - * - * x/2 = igami( df/2, y ); - * - * - * - * - * ACCURACY: - * - * See igami.c. - * - * ERROR MESSAGES: - * - * message condition value returned - * chdtri domain y < 0 or y > 1 0.0 - * v < 1 - * - */ - -/* clog.c - * - * Complex natural logarithm - * - * - * - * SYNOPSIS: - * - * void clog(); - * cmplx z, w; - * - * clog( &z, &w ); - * - * - * - * DESCRIPTION: - * - * Returns complex logarithm to the base e (2.718...) of - * the complex argument x. - * - * If z = x + iy, r = sqrt( x**2 + y**2 ), - * then - * w = log(r) + i arctan(y/x). - * - * The arctangent ranges from -PI to +PI. - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC -10,+10 7000 8.5e-17 1.9e-17 - * IEEE -10,+10 30000 5.0e-15 1.1e-16 - * - * Larger relative error can be observed for z near 1 +i0. - * In IEEE arithmetic the peak absolute error is 5.2e-16, rms - * absolute error 1.0e-16. - */ - -/* cexp() - * - * Complex exponential function - * - * - * - * SYNOPSIS: - * - * void cexp(); - * cmplx z, w; - * - * cexp( &z, &w ); - * - * - * - * DESCRIPTION: - * - * Returns the exponential of the complex argument z - * into the complex result w. - * - * If - * z = x + iy, - * r = exp(x), - * - * then - * - * w = r cos y + i r sin y. - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC -10,+10 8700 3.7e-17 1.1e-17 - * IEEE -10,+10 30000 3.0e-16 8.7e-17 - * - */ - /* csin() - * - * Complex circular sine - * - * - * - * SYNOPSIS: - * - * void csin(); - * cmplx z, w; - * - * csin( &z, &w ); - * - * - * - * DESCRIPTION: - * - * If - * z = x + iy, - * - * then - * - * w = sin x cosh y + i cos x sinh y. - * - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC -10,+10 8400 5.3e-17 1.3e-17 - * IEEE -10,+10 30000 3.8e-16 1.0e-16 - * Also tested by csin(casin(z)) = z. - * - */ - /* ccos() - * - * Complex circular cosine - * - * - * - * SYNOPSIS: - * - * void ccos(); - * cmplx z, w; - * - * ccos( &z, &w ); - * - * - * - * DESCRIPTION: - * - * If - * z = x + iy, - * - * then - * - * w = cos x cosh y - i sin x sinh y. - * - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC -10,+10 8400 4.5e-17 1.3e-17 - * IEEE -10,+10 30000 3.8e-16 1.0e-16 - */ - /* ctan() - * - * Complex circular tangent - * - * - * - * SYNOPSIS: - * - * void ctan(); - * cmplx z, w; - * - * ctan( &z, &w ); - * - * - * - * DESCRIPTION: - * - * If - * z = x + iy, - * - * then - * - * sin 2x + i sinh 2y - * w = --------------------. - * cos 2x + cosh 2y - * - * On the real axis the denominator is zero at odd multiples - * of PI/2. The denominator is evaluated by its Taylor - * series near these points. - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC -10,+10 5200 7.1e-17 1.6e-17 - * IEEE -10,+10 30000 7.2e-16 1.2e-16 - * Also tested by ctan * ccot = 1 and catan(ctan(z)) = z. - */ - /* ccot() - * - * Complex circular cotangent - * - * - * - * SYNOPSIS: - * - * void ccot(); - * cmplx z, w; - * - * ccot( &z, &w ); - * - * - * - * DESCRIPTION: - * - * If - * z = x + iy, - * - * then - * - * sin 2x - i sinh 2y - * w = --------------------. - * cosh 2y - cos 2x - * - * On the real axis, the denominator has zeros at even - * multiples of PI/2. Near these points it is evaluated - * by a Taylor series. - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC -10,+10 3000 6.5e-17 1.6e-17 - * IEEE -10,+10 30000 9.2e-16 1.2e-16 - * Also tested by ctan * ccot = 1 + i0. - */ - /* casin() - * - * Complex circular arc sine - * - * - * - * SYNOPSIS: - * - * void casin(); - * cmplx z, w; - * - * casin( &z, &w ); - * - * - * - * DESCRIPTION: - * - * Inverse complex sine: - * - * 2 - * w = -i clog( iz + csqrt( 1 - z ) ). - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * DEC -10,+10