diff --git a/src/libc/acos.src b/src/libc/acosf.src similarity index 100% rename from src/libc/acos.src rename to src/libc/acosf.src diff --git a/src/libc/acosh.c b/src/libc/acoshf.c similarity index 64% rename from src/libc/acosh.c rename to src/libc/acoshf.c index 34ef3a415..b42f45646 100644 --- a/src/libc/acosh.c +++ b/src/libc/acoshf.c @@ -6,8 +6,3 @@ float acoshf(float x) } double acosh(double) __attribute__((alias("acoshf"))); - -long double acoshl(long double x) -{ - return logl(x + sqrtl(x * x - 1)); -} diff --git a/src/libc/acoshl.c b/src/libc/acoshl.c new file mode 100644 index 000000000..d30cdad42 --- /dev/null +++ b/src/libc/acoshl.c @@ -0,0 +1,6 @@ +#include + +long double acoshl(long double x) +{ + return logl(x + sqrtl(x * x - 1)); +} diff --git a/src/libc/asin.c b/src/libc/asin.c deleted file mode 100644 index 03f647fdb..000000000 --- a/src/libc/asin.c +++ /dev/null @@ -1,37 +0,0 @@ -/* Copyright (c) 2000-2008 Zilog, Inc. */ - -/** - * asin(arg) return the arcsin of arg. - * arctan is called after appropriate range reduction. - */ - -#include -#include - -#define pio2 1.57079632679490 - -float _asinf_c(float arg) { - float sign, temp; - - sign = 1.; - if(arg < 0) { - arg = -arg; - sign = -1.; - } - - if(arg > 1.) { - errno = EDOM; - return(0.); - } - - temp = sqrt(1. - arg*arg); - if(arg > 0.7) { - temp = pio2 - atan(temp/arg); - } else { - temp = atan(arg/temp); - } - - return(sign*temp); -} - -double _asin_c(double) __attribute__((alias("_asinf_c"))); diff --git a/src/libc/asinf.c b/src/libc/asinf.c new file mode 100644 index 000000000..8a24fe5a8 --- /dev/null +++ b/src/libc/asinf.c @@ -0,0 +1,37 @@ +/* Copyright (c) 2000-2008 Zilog, Inc. */ + +/** + * asin(arg) return the arcsin of arg. + * arctan is called after appropriate range reduction. + */ + +#include +#include + +#define pio2 1.57079632679490 + +float _asinf_c(float arg) { + float sign, temp; + + sign = 1.; + if(arg < 0) { + arg = -arg; + sign = -1.; + } + + if(arg > 1.) { + errno = EDOM; + return(0.); + } + + temp = sqrtf(1. - arg*arg); + if(arg > 0.7) { + temp = pio2 - atanf(temp/arg); + } else { + temp = atanf(arg/temp); + } + + return(sign*temp); +} + +double _asin_c(double) __attribute__((alias("_asinf_c"))); diff --git a/src/libc/asin.src b/src/libc/asinf.src similarity index 100% rename from src/libc/asin.src rename to src/libc/asinf.src diff --git a/src/libc/asinh.c b/src/libc/asinhf.c similarity index 64% rename from src/libc/asinh.c rename to src/libc/asinhf.c index deb268d2b..d3c41bde7 100644 --- a/src/libc/asinh.c +++ b/src/libc/asinhf.c @@ -6,8 +6,3 @@ float asinhf(float x) } double asinh(double) __attribute__((alias("asinhf"))); - -long double asinhl(long double x) -{ - return logl(x + sqrtl(x * x + 1)); -} diff --git a/src/libc/asinhl.c b/src/libc/asinhl.c new file mode 100644 index 000000000..603ad0889 --- /dev/null +++ b/src/libc/asinhl.c @@ -0,0 +1,6 @@ +#include + +long double asinhl(long double x) +{ + return logl(x + sqrtl(x * x + 1)); +} diff --git a/src/libc/atan.c b/src/libc/atan.c deleted file mode 100644 index 31eb395ab..000000000 --- a/src/libc/atan.c +++ /dev/null @@ -1,84 +0,0 @@ -/* Copyright (c) 2000-2008 Zilog, Inc. */ - -/** - * floating-point arctangent - * - * atan returns the value of the arctangent of its - * argument in the range [-pi/2,pi/2]. - * - * atan2 returns the arctangent of arg1/arg2 - * in the range [-pi,pi]. - * - * there are no error returns. - * - * coefficients are #5077 from Hart & Cheney. (19.56D) - */ -#include - -#define sq2p1 2.41421356237309e0 -#define sq2m1 0.414213562373095e0 -#define pio2 1.57079632679489e0 -#define pio4 0.785398163397448e0 -#define p4 0.161536412982230e2 -#define p3 0.268425481955040e3 -#define p2 0.115302935154049e4 -#define p1 0.178040631643320e4 -#define p0 0.896785974036639e3 -#define q4 0.589569705084446e2 -#define q3 0.536265374031215e3 -#define q2 0.166678381488163e4 -#define q1 0.207933497444541e4 -#define q0 0.896785974036639e3 - -/** - * atan makes its argument positive and - * calls the inner routine satan. - */ - -float _atanf_c(float arg) { - float satan(float); - - if(arg>0) { - return(satan(arg)); - } else { - return(-satan(-arg)); - } -} - -double _atan_c(double) __attribute__((alias("_atanf_c"))); - -/** - * atan2 discovers what quadrant the angle - * is in and calls atan. - */ - - -/** - * xatan evaluates a series valid in the - * range [-0.414...,+0.414...]. - */ - -static float xatan(float arg) { - float argsq; - float value; - - argsq = arg*arg; - value = ((((p4*argsq + p3)*argsq + p2)*argsq + p1)*argsq + p0); - value = value/(((((argsq + q4)*argsq + q3)*argsq + q2)*argsq + q1)*argsq + q0); - return(value*arg); -} - -/** - * satan reduces its argument (known to be positive) - * to the range [0,0.414...] and calls xatan. - */ - -float satan(float arg) { - if(arg < sq2m1) { - return(xatan(arg)); - } else if(arg > sq2p1) { - return(pio2 - xatan(1.0/arg)); - } else { - return(pio4 + xatan((arg-1.0)/(arg+1.0))); - } -} diff --git a/src/libc/atan2.c b/src/libc/atan2.c deleted file mode 100644 index 046e4eb72..000000000 --- a/src/libc/atan2.c +++ /dev/null @@ -1,28 +0,0 @@ -/* Copyright (c) 2000-2008 Zilog, Inc. */ -#include - -#define pio2 1.57079632679489e0 - -float _atan2f_c(float arg1,float arg2) { - float satan(float); - - if((arg1+arg2)==arg1) { - if(arg1 >= 0.) { - return(pio2); - } else { - return(-pio2); - } - } else if(arg2 < 0.) { - if(arg1 >= 0.) { - return(pio2+pio2 - satan(-arg1/arg2)); - } else { - return(-pio2-pio2 + satan(arg1/arg2)); - } - } else if(arg1 > 0) { - return(satan(arg1/arg2)); - } else { - return(-satan(-arg1/arg2)); - } -} - -double _atan2_c(double, double) __attribute__((alias("_atan2f_c"))); diff --git a/src/libc/atan2f.c b/src/libc/atan2f.c new file mode 100644 index 000000000..48957946c --- /dev/null +++ b/src/libc/atan2f.c @@ -0,0 +1,28 @@ +/* Copyright (c) 2000-2008 Zilog, Inc. */ +#include + +#define pio2 1.57079632679489e0 + +float _atan2f_c(float arg1,float arg2) { + float satan(float); + + if((arg1+arg2)==arg1) { + if(arg1 >= 0.) { + return(pio2); + } else { + return(-pio2); + } + } else if(arg2 < 0.) { + if(arg1 >= 0.) { + return(pio2+pio2 - satan(-arg1/arg2)); + } else { + return(-pio2-pio2 + satan(arg1/arg2)); + } + } else if(arg1 > 0) { + return(satan(arg1/arg2)); + } else { + return(-satan(-arg1/arg2)); + } +} + +double _atan2_c(double, double) __attribute__((alias("_atan2f_c"))); diff --git a/src/libc/atan2.src b/src/libc/atan2f.src similarity index 100% rename from src/libc/atan2.src rename to src/libc/atan2f.src diff --git a/src/libc/atanf.c b/src/libc/atanf.c new file mode 100644 index 000000000..152d21b03 --- /dev/null +++ b/src/libc/atanf.c @@ -0,0 +1,84 @@ +/* Copyright (c) 2000-2008 Zilog, Inc. */ + +/** + * floating-point arctangent + * + * atan returns the value of the arctangent of its + * argument in the range [-pi/2,pi/2]. + * + * atan2 returns the arctangent of arg1/arg2 + * in the range [-pi,pi]. + * + * there are no error returns. + * + * coefficients are #5077 from Hart & Cheney. (19.56D) + */ +#include + +#define sq2p1 2.41421356237309e0 +#define sq2m1 0.414213562373095e0 +#define pio2 1.57079632679489e0 +#define pio4 0.785398163397448e0 +#define p4 0.161536412982230e2 +#define p3 0.268425481955040e3 +#define p2 0.115302935154049e4 +#define p1 0.178040631643320e4 +#define p0 0.896785974036639e3 +#define q4 0.589569705084446e2 +#define q3 0.536265374031215e3 +#define q2 0.166678381488163e4 +#define q1 0.207933497444541e4 +#define q0 0.896785974036639e3 + +/** + * atan makes its argument positive and + * calls the inner routine satan. + */ + +float _atanf_c(float arg) { + float satan(float); + + if(arg>0) { + return(satan(arg)); + } else { + return(-satan(-arg)); + } +} + +double _atan_c(double) __attribute__((alias("_atanf_c"))); + +/** + * atan2 discovers what quadrant the angle + * is in and calls atan. + */ + + +/** + * xatan evaluates a series valid in the + * range [-0.414...,+0.414...]. + */ + +static float xatan(float arg) { + float argsq; + float value; + + argsq = arg*arg; + value = ((((p4*argsq + p3)*argsq + p2)*argsq + p1)*argsq + p0); + value = value/(((((argsq + q4)*argsq + q3)*argsq + q2)*argsq + q1)*argsq + q0); + return(value*arg); +} + +/** + * satan reduces its argument (known to be positive) + * to the range [0,0.414...] and calls xatan. + */ + +float satan(float arg) { + if(arg < sq2m1) { + return(xatan(arg)); + } else if(arg > sq2p1) { + return(pio2 - xatan(1.0/arg)); + } else { + return(pio4 + xatan((arg-1.0)/(arg+1.0))); + } +} diff --git a/src/libc/atan.src b/src/libc/atanf.src similarity index 100% rename from src/libc/atan.src rename to src/libc/atanf.src diff --git a/src/libc/atanh.c b/src/libc/atanhf.c similarity index 63% rename from src/libc/atanh.c rename to src/libc/atanhf.c index 8ae497138..c62c075dc 100644 --- a/src/libc/atanh.c +++ b/src/libc/atanhf.c @@ -6,8 +6,3 @@ float atanhf(float x) } double atanh(double) __attribute__((alias("atanhf"))); - -long double atanhl(long double x) -{ - return .5l * logl((1 + x) / (1 - x)); -} diff --git a/src/libc/atanhl.c b/src/libc/atanhl.c new file mode 100644 index 000000000..e3c0a756c --- /dev/null +++ b/src/libc/atanhl.c @@ -0,0 +1,6 @@ +#include + +long double atanhl(long double x) +{ + return .5L * logl((1 + x) / (1 - x)); +} diff --git a/src/libc/bsearch.c b/src/libc/bsearch.c index e43eb543a..77cbb8c7d 100644 --- a/src/libc/bsearch.c +++ b/src/libc/bsearch.c @@ -1,11 +1,12 @@ /************************************************************************/ -/* */ -/* Copyright (C)1987-2008 by */ -/* Zilog, Inc. */ -/* */ -/* San Jose, California */ -/* */ +/* */ +/* Copyright (C)1987-2008 by */ +/* Zilog, Inc. */ +/* */ +/* San Jose, California */ +/* */ /************************************************************************/ + #include #include @@ -14,62 +15,63 @@ * bsearch - binary search * * Inputs: -* key - key to search for -* base - base of array to be sorted -* num - number of elements -* width - size of each element -* comp - pointer to function for comparison +* key - key to search for +* base - base of array to be sorted +* num - number of elements +* width - size of each element +* comp - pointer to function for comparison * * Returns: -* nothing +* nothing * *************************************************/ -void *bsearch(void *keyp, void *ptr, size_t num, size_t width, - int (*comp)(const void *, const void *)) -{ - char *key = keyp; - char *base = ptr; - unsigned int mid; - unsigned int low; - unsigned int high; - int d; - unsigned int pmid; - char *addr; +void *bsearch( + void *keyp, void *ptr, size_t num, size_t width, + int (*comp)(const void *, const void *) +) { + char *key = keyp; + char *base = ptr; + unsigned int mid; + unsigned int low; + unsigned int high; + int d; + unsigned int pmid; + char *addr; - /* make high the nearest power of two for */ - /* efficiency and to ensure we always */ - /* terminate properly. */ - /* Note that high is always higher */ - /* than it should be so that we will */ - /* not fail to find the last entry in the */ - /* table. */ + /* make high the nearest power of two for */ + /* efficiency and to ensure we always */ + /* terminate properly. */ + /* Note that high is always higher */ + /* than it should be so that we will */ + /* not fail to find the last entry in the */ + /* table. */ - high = 0x0001; - while (high <= num) - high <<= 1; - low = 0; - mid = 0; + high = 0x0001; + while (high <= num) + high <<= 1; + low = 0; + mid = 0; - /* begin the search */ + /* begin the search */ - for(;;) { - pmid = mid; - mid = ((high - low) >> 1) + low; + for(;;) { + pmid = mid; + mid = ((high - low) >> 1) + low; - if (pmid == mid) /* terminate because we're */ - return(NULL); /* oscilating. */ + if (pmid == mid) /* terminate because we're */ + return(NULL); /* oscilating. */ - if (mid >= num) { /* we're above the array. */ - high = mid; /* pretend element is larger*/ - continue; /* than the key. */ - } + if (mid >= num) { /* we're above the array. */ + high = mid; /* pretend element is larger*/ + continue; /* than the key. */ + } - d = (*comp)(key,addr = base + mid * width); - if (d == 0) /* we found it */ - return(addr); - if (d < 0) /* key is less than mid, */ - high = mid; /* set high to mid. */ - else /* key is greater than mid, */ - low = mid; /* set low to mid. */ - } + d = (*comp)(key,addr = base + mid * width); + if (d == 0) /* we found it */ + return(addr); + if (d < 0) /* key is less than mid, */ + high = mid; /* set high to mid. */ + else /* key is greater than mid, */ + low = mid; /* set low to mid. */ + } } diff --git a/src/libc/cbrt.c b/src/libc/cbrtf.c similarity index 100% rename from src/libc/cbrt.c rename to src/libc/cbrtf.c diff --git a/src/libc/cos.src b/src/libc/cosf.src similarity index 100% rename from src/libc/cos.src rename to src/libc/cosf.src diff --git a/src/libc/cosh.c b/src/libc/cosh.c deleted file mode 100644 index fe39354b0..000000000 --- a/src/libc/cosh.c +++ /dev/null @@ -1,51 +0,0 @@ -/************************************************************************/ -/* */ -/* Copyright (C)1987-2008 by */ -/* Zilog, Inc. */ -/* */ -/* San Jose, California */ -/* */ -/************************************************************************/ -/* - sinh(arg) returns the hyperbolic sine of its floating- - point argument. - - The exponential function is called for arguments - greater in magnitude than 0.5. - - A series is used for arguments smaller in magnitude than 0.5. - The coefficients are #2029 from Hart & Cheney. (20.36D) - - cosh(arg) is computed from the exponential function for - all arguments. -*/ - -#include - -#define p0 -0.630767364049772e+6 -#define p1 -0.899127202203951e+5 -#define p2 -0.289421135598956e+4 -#define p3 -0.263056321339750e+2 -#define q0 -0.630767364049772e+6 -#define q1 0.152151737879002e+5 -#define q2 -0.173678953558234e+3 - -float _coshf_c(float arg) { - float val; - - if(arg < 0) { - arg = -arg; - } - - val = expf(arg); - - if(arg > 21.) { - return val/2; - } - - val += expf(-arg); - val /= 2; - return val; -} - -double _cosh_c(double) __attribute__((alias("_coshf_c"))); diff --git a/src/libc/coshf.c b/src/libc/coshf.c new file mode 100644 index 000000000..d6516299e --- /dev/null +++ b/src/libc/coshf.c @@ -0,0 +1,51 @@ +/************************************************************************/ +/* */ +/* Copyright (C)1987-2008 by */ +/* Zilog, Inc. */ +/* */ +/* San Jose, California */ +/* */ +/************************************************************************/ +/* + sinh(arg) returns the hyperbolic sine of its floating- + point argument. + + The exponential function is called for arguments + greater in magnitude than 0.5. + + A series is used for arguments smaller in magnitude than 0.5. + The coefficients are #2029 from Hart & Cheney. (20.36D) + + cosh(arg) is computed from the exponential function for + all arguments. +*/ + +#include + +#define p0 -0.630767364049772e+6 +#define p1 -0.899127202203951e+5 +#define p2 -0.289421135598956e+4 +#define p3 -0.263056321339750e+2 +#define q0 -0.630767364049772e+6 +#define q1 0.152151737879002e+5 +#define q2 -0.173678953558234e+3 + +float _coshf_c(float arg) { + float val; + + if(arg < 0) { + arg = -arg; + } + + val = expf(arg); + + if(arg > 21.) { + return val/2; + } + + val += expf(-arg); + val /= 2; + return val; +} + +double _cosh_c(double) __attribute__((alias("_coshf_c"))); diff --git a/src/libc/cosh.src b/src/libc/coshf.src similarity index 100% rename from src/libc/cosh.src rename to src/libc/coshf.src diff --git a/src/libc/difftime.c b/src/libc/difftime.c index c98ec889a..f25b312e8 100644 --- a/src/libc/difftime.c +++ b/src/libc/difftime.c @@ -2,6 +2,6 @@ double difftime(time_t end, time_t beginning) { - /* assuming typedef unsigned long time_t */ + /* assuming typedef unsigned long time_t */ return (double)((signed long)(end - beginning)); } diff --git a/src/libc/erfc.c b/src/libc/erfc.c deleted file mode 100644 index a6ca51406..000000000 --- a/src/libc/erfc.c +++ /dev/null @@ -1,17 +0,0 @@ -#include - -float erfcf(float x) -{ - static const float p = 0.47047f, a1 = 0.3480242f, a2 = -0.0958798f, a3 = 0.7478556f; - const float t = 1.0f / (1.0f + p * fabsf(x)); - return copysignf(t * (a1 + t * (a2 + t * a3)) * expf(-x * x), x); -} - -double erfc(double) __attribute__((alias("erfcf"))); - -long double erfcl(long double x) -{ - static const long double p = 0.3275911l, a1 = 0.254829592l, a2 = -0.284496736l, a3 = 1.421413741l, a4 = -1.453152027l, a5 = 1.061405429l; - const long double t = 1.0l / (1.0l + p * fabsl(x)); - return copysignl(t * (a1 + t * (a2 + t * (a3 + t * (a4 + t * a5)))) * expl(-x * x), x); -} diff --git a/src/libc/erfcf.c b/src/libc/erfcf.c new file mode 100644 index 000000000..42a38bc15 --- /dev/null +++ b/src/libc/erfcf.c @@ -0,0 +1,22 @@ +#include + +/** + * Algorithm from: + * https://en.wikipedia.org/wiki/Error_function#Approximation_with_elementary_functions + */ +float erfcf(float x) +{ + static const float + p = 0.47047f, + a1 = 0.3480242f, + a2 = -0.0958798f, + a3 = 0.7478556f; + const float t = 1.0f / (1.0f + p * fabsf(x)); + float ret = t * (a1 + t * (a2 + t * a3)) * expf(-x * x); + if (signbit(x)) { + ret = 2.0f - ret; + } + return ret; +} + +double erfc(double) __attribute__((alias("erfcf"))); diff --git a/src/libc/erfcl.c b/src/libc/erfcl.c new file mode 100644 index 000000000..d2ce8eac5 --- /dev/null +++ b/src/libc/erfcl.c @@ -0,0 +1,22 @@ +#include + +/** + * Algorithm from: + * https://en.wikipedia.org/wiki/Error_function#Approximation_with_elementary_functions + */ +long double erfcl(long double x) +{ + static const long double + p = 0.3275911L, + a1 = 0.254829592L, + a2 = -0.284496736L, + a3 = 1.421413741L, + a4 = -1.453152027L, + a5 = 1.061405429L; + const long double t = 1.0L / (1.0L + p * fabsl(x)); + long double ret = t * (a1 + t * (a2 + t * (a3 + t * (a4 + t * a5)))) * expl(-x * x); + if (signbit(x)) { + ret = 2.0L - ret; + } + return ret; +} diff --git a/src/libc/erf.c b/src/libc/erff.c similarity index 65% rename from src/libc/erf.c rename to src/libc/erff.c index 74bb7ae33..421688f2b 100644 --- a/src/libc/erf.c +++ b/src/libc/erff.c @@ -6,8 +6,3 @@ float erff(float x) } double erf(double) __attribute__((alias("erff"))); - -long double erfl(long double x) -{ - return 1 - erfcl(x); -} diff --git a/src/libc/erfl.c b/src/libc/erfl.c new file mode 100644 index 000000000..a7027954f --- /dev/null +++ b/src/libc/erfl.c @@ -0,0 +1,6 @@ +#include + +long double erfl(long double x) +{ + return 1 - erfcl(x); +} diff --git a/src/libc/exp.c b/src/libc/exp.c deleted file mode 100644 index 7616cd492..000000000 --- a/src/libc/exp.c +++ /dev/null @@ -1,50 +0,0 @@ -/************************************************************************/ -/* */ -/* Copyright (C) 2000-2008 Zilog, Inc. */ -/* */ -/************************************************************************/ -/* - exp returns the exponential function of its - floating-point argument. - - The coefficients are #1069 from Hart and Cheney. (22.35D) -*/ - -#include -#include - -#define p0 0.20803843466947e7 -#define p1 0.30286971697440e5 -#define p2 0.60614853300611e2 -#define q0 0.60027203602388e7 -#define q1 0.32772515180829e6 -#define q2 0.17492876890931e4 -#define log2e 1.44269504088896 -#define sqrt2 1.41421356237310 -#define maxf 10000 - -float _expf_c(float arg) { - float fraction; - float temp1, temp2, xsq; - int ent; - - if ( arg == 0.0 ){ - return 1.0; - } - if ( arg < -maxf ){ - return 0.0; - } - if ( arg > maxf ){ - errno = ERANGE; - return HUGE_VALF; - } - arg *= log2e; - ent = floor( arg ); - fraction = arg - ent - 0.5; - xsq = fraction * fraction; - temp1 = ((p2 * xsq + p1) * xsq + p0) * fraction; - temp2 = ((1.0 * xsq + q2) * xsq + q1) * xsq + q0; - return ldexp( sqrt2 * (temp2+temp1) / (temp2-temp1), ent ); -} - -double _exp_c(double) __attribute__((alias("_expf_c"))); diff --git a/src/libc/exp2.c b/src/libc/exp2f.c similarity index 100% rename from src/libc/exp2.c rename to src/libc/exp2f.c diff --git a/src/libc/expf.c b/src/libc/expf.c new file mode 100644 index 000000000..d4dc3da00 --- /dev/null +++ b/src/libc/expf.c @@ -0,0 +1,50 @@ +/************************************************************************/ +/* */ +/* Copyright (C) 2000-2008 Zilog, Inc. */ +/* */ +/************************************************************************/ +/* + exp returns the exponential function of its + floating-point argument. + + The coefficients are #1069 from Hart and Cheney. (22.35D) +*/ + +#include +#include + +#define p0 0.20803843466947e7 +#define p1 0.30286971697440e5 +#define p2 0.60614853300611e2 +#define q0 0.60027203602388e7 +#define q1 0.32772515180829e6 +#define q2 0.17492876890931e4 +#define log2e 1.44269504088896 +#define sqrt2 1.41421356237310 +#define maxf 10000 + +float _expf_c(float arg) { + float fraction; + float temp1, temp2, xsq; + int ent; + + if ( arg == 0.0 ){ + return 1.0; + } + if ( arg < -maxf ){ + return 0.0; + } + if ( arg > maxf ){ + errno = ERANGE; + return HUGE_VALF; + } + arg *= log2e; + ent = floor( arg ); + fraction = arg - ent - 0.5; + xsq = fraction * fraction; + temp1 = ((p2 * xsq + p1) * xsq + p0) * fraction; + temp2 = ((1.0 * xsq + q2) * xsq + q1) * xsq + q0; + return ldexp( sqrt2 * (temp2+temp1) / (temp2-temp1), ent ); +} + +double _exp_c(double) __attribute__((alias("_expf_c"))); diff --git a/src/libc/exp.src b/src/libc/expf.src similarity index 100% rename from src/libc/exp.src rename to src/libc/expf.src diff --git a/src/libc/expm1.c b/src/libc/expm1f.c similarity index 100% rename from src/libc/expm1.c rename to src/libc/expm1f.c diff --git a/src/libc/floorf.c b/src/libc/floorf.c index 582058c6e..a5e5ea5fd 100644 --- a/src/libc/floorf.c +++ b/src/libc/floorf.c @@ -1,25 +1,26 @@ /************************************************************************/ -/* */ -/* Copyright (C) 2000-2008 Zilog, Inc. */ -/* */ -/* San Jose, California */ -/* */ +/* */ +/* Copyright (C) 2000-2008 Zilog, Inc. */ +/* */ +/* San Jose, California */ +/* */ /************************************************************************/ #include float _floorf_c(float d) { - float fraction; + float fraction; - if (d<0.0) { - d = -d; - fraction = modff(d, &d); - if (fraction != 0.0) - d += 1; - d = -d; - } else { - fraction = modff(d, &d); - } - return(d); + if (d<0.0) { + d = -d; + fraction = modff(d, &d); + if (fraction != 0.0) { + d += 1; + } + d = -d; + } else { + fraction = modff(d, &d); + } + return(d); } double _floor_c(double) __attribute__((alias("_floorf_c"))); diff --git a/src/libc/frexpf.c b/src/libc/frexpf.c index 4a549a74a..4fe24b048 100644 --- a/src/libc/frexpf.c +++ b/src/libc/frexpf.c @@ -5,9 +5,9 @@ #if 1 /************************************************************************/ /* */ -/* Copyright (C) 2000-2008 Zilog, Inc. */ +/* Copyright (C) 2000-2008 Zilog, Inc. */ /* */ -/* San Jose, California */ +/* San Jose, California */ /* */ /************************************************************************/ diff --git a/src/libc/lgamma.c b/src/libc/lgamma.c deleted file mode 100644 index 3fb425833..000000000 --- a/src/libc/lgamma.c +++ /dev/null @@ -1,37 +0,0 @@ -/* - * gamma.c - public domain implementation of function tgamma(3m) - * reference - Haruhiko Okumura: C-gengo niyoru saishin algorithm jiten - * (New Algorithm handbook in C language) (Gijyutsu hyouronsha, Tokyo, 1991) [in Japanese] - * http://oku.edu.mie-u.ac.jp/~okumura/algo/ -*/ - -#include -#include - -#define N 8 - -#define B0 1 /* Bernoulli numbers */ -#define B1 (-1.0 / 2.0) -#define B2 ( 1.0 / 6.0) -#define B4 (-1.0 / 30.0) -#define B6 ( 1.0 / 42.0) -#define B8 (-1.0 / 30.0) -#define B10 ( 5.0 / 66.0) -#define B12 (-691.0 / 2730.0) -#define B14 ( 7.0 / 6.0) -#define B16 (-3617.0 / 510.0) - -float lgammaf(float x) { /* the natural logarithm of the Gamma function. */ - float v, w; - - v = 1; - while (x < N) { v *= x; x++; } - w = 1 / (x * x); - return ((((((((B16 / (16 * 15)) * w + (B14 / (14 * 13))) * w - + (B12 / (12 * 11))) * w + (B10 / (10 * 9))) * w - + (B8 / ( 8 * 7))) * w + (B6 / ( 6 * 5))) * w - + (B4 / ( 4 * 3))) * w + (B2 / ( 2 * 1))) / x - + 0.5 * M_LOG_2M_PI - log(v) - x + (x - 0.5) * log(x); -} - -double lgamma(double) __attribute__((alias("lgammaf"))); diff --git a/src/libc/lgammaf.c b/src/libc/lgammaf.c new file mode 100644 index 000000000..07fb48c37 --- /dev/null +++ b/src/libc/lgammaf.c @@ -0,0 +1,49 @@ +/* + * gamma.c - public domain implementation of function tgamma(3m) + * reference - Haruhiko Okumura: C-gengo niyoru saishin algorithm jiten + * (New Algorithm handbook in C language) (Gijyutsu hyouronsha, Tokyo, 1991) [in Japanese] + * http://oku.edu.mie-u.ac.jp/~okumura/algo/ +*/ + +#include +#include + +#define N 8 + +#define B0 1.0 /* Bernoulli numbers */ +#define B1 (-1.0 / 2.0) +#define B2 ( 1.0 / 6.0) +#define B4 (-1.0 / 30.0) +#define B6 ( 1.0 / 42.0) +#define B8 (-1.0 / 30.0) +#define B10 ( 5.0 / 66.0) +#define B12 (-691.0 / 2730.0) +#define B14 ( 7.0 / 6.0) +#define B16 (-3617.0 / 510.0) + +float lgammaf(float x) { /* the natural logarithm of the Gamma function. */ + float v, w; + v = 1.0; + + /** + * This loop will take forever to terminate if `x < -100.0f`, so we have a + * maximum iteration count to ensure that the loop will terminate in a + * reasonable amount of time. `v` should overflow when `x < -33.0f` + */ + const int maximum_iter = 36 + N; + for (int iter = 0; iter < maximum_iter; iter++) { + if (x < (float)N) { + break; + } + v *= x; + x++; + } + w = 1.0 / (x * x); + return ((((((((B16 / (16.0 * 15.0)) * w + (B14 / (14.0 * 13.0))) * w + + (B12 / (12.0 * 11.0))) * w + (B10 / (10.0 * 9.0))) * w + + (B8 / ( 8.0 * 7.0))) * w + (B6 / ( 6.0 * 5.0))) * w + + (B4 / ( 4.0 * 3.0))) * w + (B2 / ( 2.0 * 1.0))) / x + + 0.5 * (float)M_LOG_2M_PI - logf(v) - x + (x - 0.5) * logf(x); +} + +double lgamma(double) __attribute__((alias("lgammaf"))); diff --git a/src/libc/llrint.c b/src/libc/llrintf.c similarity index 100% rename from src/libc/llrint.c rename to src/libc/llrintf.c diff --git a/src/libc/llround.c b/src/libc/llroundf.c similarity index 100% rename from src/libc/llround.c rename to src/libc/llroundf.c diff --git a/src/libc/log.c b/src/libc/log.c deleted file mode 100644 index f6b735215..000000000 --- a/src/libc/log.c +++ /dev/null @@ -1,53 +0,0 @@ -/************************************************************************/ -/* */ -/* Copyright (C)1999-2008 by Zilog, Inc. */ -/* */ -/************************************************************************/ -/* - log returns the natural logarithm of its floating - point argument. - - The coefficients are #2705 from Hart & Cheney. (19.38D) - - It calls frexp. -*/ - -#include -#include - -#define log2 0.693147180559945e0 -#define ln10 2.30258509299404 -#define sqrto2 0.707106781186548e0 -#define p0 -0.240139179559211e2 -#define p1 0.309572928215377e2 -#define p2 -0.963769093368687e1 -#define p3 0.421087371217980e0 -#define q0 -0.120069589779605e2 -#define q1 0.194809660700890e2 -#define q2 -0.891110902798312e1 - -float _logf_c(float arg) -{ - double x,z, zsq, temp; - int exp; - - if (arg <= 0.0) { - errno = EDOM; - return -HUGE_VALF; - } - x = frexpf(arg, & exp); - if ( x < sqrto2 ){ - x *= 2; - exp--; - } - - z = (x-1)/(x+1); - zsq = z*z; - - temp = ((p3*zsq + p2)*zsq + p1)*zsq + p0; - temp = temp/(((1.0*zsq + q2)*zsq + q1)*zsq + q0); - temp = temp*z + exp*log2; - return temp; -} - -double _log_c(double) __attribute__((alias("_logf_c"))); diff --git a/src/libc/log10.c b/src/libc/log10f.c similarity index 73% rename from src/libc/log10.c rename to src/libc/log10f.c index 2eef9e58b..36c815861 100644 --- a/src/libc/log10.c +++ b/src/libc/log10f.c @@ -2,7 +2,7 @@ float _log10f_c(float x) { - return logf(x) / M_LN10; + return logf(x) * (float)(M_LOG10E); } double _log10_c(double) __attribute__((alias("_log10f_c"))); diff --git a/src/libc/log10.src b/src/libc/log10f.src similarity index 100% rename from src/libc/log10.src rename to src/libc/log10f.src diff --git a/src/libc/log1p.c b/src/libc/log1pf.c similarity index 93% rename from src/libc/log1p.c rename to src/libc/log1pf.c index 99f596665..f4123ead8 100644 --- a/src/libc/log1p.c +++ b/src/libc/log1pf.c @@ -1,7 +1,7 @@ #include float log1pf(float x) { - if (fabs(x) <= 0.125) { + if (fabsf(x) <= 0.125) { // pade(series(ln(1+x),x=0,6,polynom),x,5,3) // (-57*x**2-90*x)/(x**3-21*x**2-102*x-90) // relative error less than 1e-7 diff --git a/src/libc/log2.c b/src/libc/log2f.c similarity index 70% rename from src/libc/log2.c rename to src/libc/log2f.c index 6ffe86aa3..0e3d17103 100644 --- a/src/libc/log2.c +++ b/src/libc/log2f.c @@ -2,7 +2,7 @@ float log2f(float x) { - return logf(x) * (float)(1 / M_LN2); + return logf(x) * (float)(M_LOG2E); } double log2(double) __attribute__((alias("log2f"))); diff --git a/src/libc/logf.c b/src/libc/logf.c new file mode 100644 index 000000000..6bff1d1a2 --- /dev/null +++ b/src/libc/logf.c @@ -0,0 +1,53 @@ +/************************************************************************/ +/* */ +/* Copyright (C) 1999-2008 by Zilog, Inc. */ +/* */ +/************************************************************************/ +/* + log returns the natural logarithm of its floating + point argument. + + The coefficients are #2705 from Hart & Cheney. (19.38D) + + It calls frexp. +*/ + +#include +#include + +#define ln2 0.693147180559945e0 +#define ln10 2.30258509299404 +#define sqrto2 0.707106781186548e0 +#define p0 -0.240139179559211e2 +#define p1 0.309572928215377e2 +#define p2 -0.963769093368687e1 +#define p3 0.421087371217980e0 +#define q0 -0.120069589779605e2 +#define q1 0.194809660700890e2 +#define q2 -0.891110902798312e1 + +float _logf_c(float arg) +{ + float x, z, zsq, temp; + int expon; + + if (arg <= 0.0) { + errno = EDOM; + return -HUGE_VALF; + } + x = frexpf(arg, & expon); + if ( x < sqrto2 ){ + x *= 2; + expon--; + } + + z = (x-1)/(x+1); + zsq = z*z; + + temp = ((p3*zsq + p2)*zsq + p1)*zsq + p0; + temp = temp/(((1.0*zsq + q2)*zsq + q1)*zsq + q0); + temp = temp*z + expon*ln2; + return temp; +} + +double _log_c(double) __attribute__((alias("_logf_c"))); diff --git a/src/libc/log.src b/src/libc/logf.src similarity index 100% rename from src/libc/log.src rename to src/libc/logf.src diff --git a/src/libc/lrint.c b/src/libc/lrintf.c similarity index 100% rename from src/libc/lrint.c rename to src/libc/lrintf.c diff --git a/src/libc/lround.c b/src/libc/lroundf.c similarity index 100% rename from src/libc/lround.c rename to src/libc/lroundf.c diff --git a/src/libc/nextafterf.c b/src/libc/nextafterf.c new file mode 100644 index 000000000..65f4e5ee8 --- /dev/null +++ b/src/libc/nextafterf.c @@ -0,0 +1,58 @@ +#include +#include +#include + +typedef union F32_pun { + float flt; + uint32_t bin; +} F32_pun; + +#define Float32_pos_subnorm_min UINT32_C(0x00000001) +#define Float32_neg_subnorm_min UINT32_C(0x80000001) + +float nextafterf(float x, float y) { + F32_pun arg_x, arg_y; + arg_x.flt = x; + arg_y.flt = y; + + if (isnan(y)) { + return y; + } + if (isnan(x)) { + return x; + } + if (arg_x.bin == arg_y.bin) { + return y; + } + + if (iszero(x)) { + if (iszero(y)) { + // special case where `+0.0 --> -0.0` and `-0.0 --> +0.0` + return y; + } + feraiseexcept(FE_INEXACT | FE_UNDERFLOW); + F32_pun ret; + ret.bin = signbit(y) ? Float32_neg_subnorm_min : Float32_pos_subnorm_min; + return ret.flt; + } + if (isless(x, y) != signbit(x)) { + // Towards positive/negative infinity + arg_x.bin++; + } else { + // Towards negative/positive zero + arg_x.bin--; + } + if (isnormal(arg_x.flt)) { + return arg_x.flt; + } + if (isinf(arg_x.flt)) { + // overflow to infinity + feraiseexcept(FE_INEXACT | FE_OVERFLOW); + } else { + // result is subnormal or zero + feraiseexcept(FE_INEXACT | FE_UNDERFLOW); + } + return arg_x.flt; +} + +double nextafter(double, double) __attribute__((alias("nextafterf"))); diff --git a/src/libc/nextafter.c b/src/libc/nextafterl.c similarity index 52% rename from src/libc/nextafter.c rename to src/libc/nextafterl.c index 056755aa0..b478c3c3a 100644 --- a/src/libc/nextafter.c +++ b/src/libc/nextafterl.c @@ -2,63 +2,6 @@ #include #include -typedef union F32_pun { - float flt; - uint32_t bin; -} F32_pun; - -#define Float32_pos_subnorm_min UINT32_C(0x00000001) -#define Float32_neg_subnorm_min UINT32_C(0x80000001) - -float nextafterf(float x, float y) { - F32_pun arg_x, arg_y; - arg_x.flt = x; - arg_y.flt = y; - - if (isnan(y)) { - return y; - } - if (isnan(x)) { - return x; - } - if (arg_x.bin == arg_y.bin) { - return y; - } - - if (iszero(x)) { - if (iszero(y)) { - // special case where `+0.0 --> -0.0` and `-0.0 --> +0.0` - return y; - } - feraiseexcept(FE_INEXACT | FE_UNDERFLOW); - F32_pun ret; - ret.bin = signbit(y) ? Float32_neg_subnorm_min : Float32_pos_subnorm_min; - return ret.flt; - } - if (isless(x, y) != signbit(x)) { - // Towards positive/negative infinity - arg_x.bin++; - } else { - // Towards negative/positive zero - arg_x.bin--; - } - if (isnormal(arg_x.flt)) { - return arg_x.flt; - } - if (isinf(arg_x.flt)) { - // overflow to infinity - feraiseexcept(FE_INEXACT | FE_OVERFLOW); - } else { - // result is subnormal or zero - feraiseexcept(FE_INEXACT | FE_UNDERFLOW); - } - return arg_x.flt; -} - -double nextafter(double, double) __attribute__((alias("nextafterf"))); - - - typedef union F64_pun { long double flt; uint64_t bin; diff --git a/src/libc/nextupdown.c b/src/libc/nextupdownf.c similarity index 52% rename from src/libc/nextupdown.c rename to src/libc/nextupdownf.c index f171c1407..89d1358c5 100644 --- a/src/libc/nextupdown.c +++ b/src/libc/nextupdownf.c @@ -37,6 +37,8 @@ float nextupf(float x) { return arg_x.flt; } +double nextup(double) __attribute__((alias("nextupf"))); + float nextdownf(float x) { F32_pun arg_x; arg_x.flt = x; @@ -62,67 +64,4 @@ float nextdownf(float x) { return arg_x.flt; } -double nextup(double) __attribute__((alias("nextupf"))); double nextdown(double) __attribute__((alias("nextdownf"))); - - - -typedef union F64_pun { - long double flt; - uint64_t bin; -} F64_pun; - -#define Float64_pos_subnorm_min UINT64_C(0x0000000000000001) -#define Float64_neg_subnorm_min UINT64_C(0x8000000000000001) -#define Float64_pos_inf UINT64_C(0x7FF0000000000000) -#define Float64_neg_inf UINT64_C(0xFFF0000000000000) - -long double nextupl(long double x) { - F64_pun arg_x; - arg_x.flt = x; - - if (isnan(x) || arg_x.bin == Float64_pos_inf) { - // return unmodified - return x; - } - - if (iszero(x)) { - F64_pun ret; - ret.bin = Float64_pos_subnorm_min; - return ret.flt; - } - - if (signbit(x)) { - // Towards negative zero - arg_x.bin--; - } else { - // Towards positive infinity - arg_x.bin++; - } - return arg_x.flt; -} - -long double nextdownl(long double x) { - F64_pun arg_x; - arg_x.flt = x; - - if (isnan(x) || arg_x.bin == Float64_neg_inf) { - // return unmodified - return x; - } - - if (iszero(x)) { - F64_pun ret; - ret.bin = Float64_neg_subnorm_min; - return ret.flt; - } - - if (signbit(x)) { - // Towards negative infinity - arg_x.bin++; - } else { - // Towards positive zero - arg_x.bin--; - } - return arg_x.flt; -} diff --git a/src/libc/nextupdownl.c b/src/libc/nextupdownl.c new file mode 100644 index 000000000..f2b98cce3 --- /dev/null +++ b/src/libc/nextupdownl.c @@ -0,0 +1,63 @@ +#include +#include +#include + +typedef union F64_pun { + long double flt; + uint64_t bin; +} F64_pun; + +#define Float64_pos_subnorm_min UINT64_C(0x0000000000000001) +#define Float64_neg_subnorm_min UINT64_C(0x8000000000000001) +#define Float64_pos_inf UINT64_C(0x7FF0000000000000) +#define Float64_neg_inf UINT64_C(0xFFF0000000000000) + +long double nextupl(long double x) { + F64_pun arg_x; + arg_x.flt = x; + + if (isnan(x) || arg_x.bin == Float64_pos_inf) { + // return unmodified + return x; + } + + if (iszero(x)) { + F64_pun ret; + ret.bin = Float64_pos_subnorm_min; + return ret.flt; + } + + if (signbit(x)) { + // Towards negative zero + arg_x.bin--; + } else { + // Towards positive infinity + arg_x.bin++; + } + return arg_x.flt; +} + +long double nextdownl(long double x) { + F64_pun arg_x; + arg_x.flt = x; + + if (isnan(x) || arg_x.bin == Float64_neg_inf) { + // return unmodified + return x; + } + + if (iszero(x)) { + F64_pun ret; + ret.bin = Float64_neg_subnorm_min; + return ret.flt; + } + + if (signbit(x)) { + // Towards negative infinity + arg_x.bin++; + } else { + // Towards positive zero + arg_x.bin--; + } + return arg_x.flt; +} diff --git a/src/libc/pow.c b/src/libc/pow.c deleted file mode 100644 index ebaf6bd3b..000000000 --- a/src/libc/pow.c +++ /dev/null @@ -1,34 +0,0 @@ -/************************************************************************/ -/* */ -/* Copyright (C)1987-2008 by */ -/* Zilog, Inc. */ -/* */ -/* San Jose, California */ -/* */ -/************************************************************************/ - -#include -#include - -float _powf_c(float arg1, float arg2) { - float result; - long temp; - - if ( arg1 > 0.0 ){ - return expf( arg2 * log( arg1 ) ); - } - if ( arg1 < 0.0 ){ - temp = arg2; - if ( temp == arg2 ){ - result = expf( arg2 * log( -arg1 ) ); - return temp & 1 ? -result : result; - } - errno = EDOM; - } - if ( arg2 <= 0.0 ){ - errno = EDOM; - } - return 0.0; -} - -double _pow_c(double, double) __attribute__((alias("_powf_c"))); diff --git a/src/libc/powf.c b/src/libc/powf.c new file mode 100644 index 000000000..f3392fd08 --- /dev/null +++ b/src/libc/powf.c @@ -0,0 +1,34 @@ +/************************************************************************/ +/* */ +/* Copyright (C)1987-2008 by */ +/* Zilog, Inc. */ +/* */ +/* San Jose, California */ +/* */ +/************************************************************************/ + +#include +#include + +float _powf_c(float arg1, float arg2) { + float result; + long temp; + + if ( arg1 > 0.0 ){ + return expf( arg2 * logf( arg1 ) ); + } + if ( arg1 < 0.0 ){ + temp = arg2; + if ( temp == arg2 ){ + result = expf( arg2 * logf( -arg1 ) ); + return temp & 1 ? -result : result; + } + errno = EDOM; + } + if ( arg2 <= 0.0 ){ + errno = EDOM; + } + return 0.0; +} + +double _pow_c(double, double) __attribute__((alias("_powf_c"))); diff --git a/src/libc/pow.src b/src/libc/powf.src similarity index 100% rename from src/libc/pow.src rename to src/libc/powf.src diff --git a/src/libc/qsort.c b/src/libc/qsort.c index 348fc6d5e..335a064ec 100644 --- a/src/libc/qsort.c +++ b/src/libc/qsort.c @@ -1,11 +1,12 @@ /************************************************************************/ -/* */ -/* Copyright (C)1987-2008 by */ -/* Zilog, Inc. */ -/* */ -/* San Jose, California */ -/* */ +/* */ +/* Copyright (C)1987-2008 by */ +/* Zilog, Inc. */ +/* */ +/* San Jose, California */ +/* */ /************************************************************************/ + #include #include @@ -14,89 +15,90 @@ * qsort - Quick sort * * Inputs: -* base - base of array to be sorted -* nel - number of elements -* size - size of each element -* compar - pointer to function for comparison +* base - base of array to be sorted +* nel - number of elements +* size - size of each element +* compar - pointer to function for comparison * * Returns: -* nothing +* nothing * *************************************************/ static void swapmem(char *a,char *b,size_t size) { - register char t; - register size_t i; + register char t; + register size_t i; - for (i=0;il = i; - sp->r = r; - ++sp; - } - r = j; /* continue sorting left partition */ - } - else { - if (base < j) { /* stack request for left partition */ - sp->l = base; - sp->r = j; - ++sp; - } - base = i; /* continue sorting right partition */ - } - } while (base < r); + sp = stack; + r = base + (nel-1)*size; + for (;;) { + do { + x = base + (r-base)/size/2 * size; + i = base; + j = r; + do { + while ((*compar)(i,x) < 0) + i += size; + while ((*compar)(x,j) < 0) + j -= size; + if (i < j) { + swapmem(i, j, size); + if (i == x) + x = j; + else if (j == x) + x = i; + } + if (i <= j) { + i += size; + j -= size; + } + } while (i <= j); + if (j-base < r-i) { + if (i < r) { /* stack request for right partition */ + sp->l = i; + sp->r = r; + ++sp; + } + r = j; /* continue sorting left partition */ + } + else { + if (base < j) { /* stack request for left partition */ + sp->l = base; + sp->r = j; + ++sp; + } + base = i; /* continue sorting right partition */ + } + } while (base < r); - if (sp <= stack) - break; - --sp; - base = sp->l; - r = sp->r; - } + if (sp <= stack) + break; + --sp; + base = sp->l; + r = sp->r; + } } diff --git a/src/libc/remainder.c b/src/libc/remainderf.c similarity index 100% rename from src/libc/remainder.c rename to src/libc/remainderf.c diff --git a/src/libc/remquo.c b/src/libc/remquof.c similarity index 100% rename from src/libc/remquo.c rename to src/libc/remquof.c diff --git a/src/libc/rint.c b/src/libc/rintf.c similarity index 100% rename from src/libc/rint.c rename to src/libc/rintf.c diff --git a/src/libc/roundeven.c b/src/libc/roundeven.c index 957a401a9..b2c8e7318 100644 --- a/src/libc/roundeven.c +++ b/src/libc/roundeven.c @@ -5,7 +5,7 @@ float roundevenf(float x) float i, f = modff(x, &i); if (!isgreaterequal(f, .5f)) return i; if (f == .5f) { - float half = ldexp(i, -1); + float half = ldexpf(i, -1); if (truncf(half) == half) return i; } return signbit(i) ? i - 1 : i + 1; diff --git a/src/libc/signbitf.src b/src/libc/signbitf.src new file mode 100644 index 000000000..3be7e06f9 --- /dev/null +++ b/src/libc/signbitf.src @@ -0,0 +1,14 @@ + assume adl=1 + + section .text + + public __signbitf + +; bool _signbitf(float) +__signbitf: + ld hl, 6 + add hl, sp + ld a, (hl) + rla + sbc a, a + ret diff --git a/src/libc/signbit.src b/src/libc/signbitl.src similarity index 51% rename from src/libc/signbit.src rename to src/libc/signbitl.src index e48ae18fe..120a32844 100644 --- a/src/libc/signbit.src +++ b/src/libc/signbitl.src @@ -2,16 +2,7 @@ section .text - public __signbitf, __signbitl - -; bool _signbitf(float) -__signbitf: - ld hl, 6 - add hl, sp - ld a, (hl) - rla - sbc a, a - ret + public __signbitl ; bool _signbitl(long double) __signbitl: diff --git a/src/libc/sin.c b/src/libc/sin.c deleted file mode 100644 index 2395c99e1..000000000 --- a/src/libc/sin.c +++ /dev/null @@ -1,65 +0,0 @@ -/************************************************************************/ -/* */ -/* Copyright (C) 1999-2008 by Zilog, Inc. */ -/* */ -/************************************************************************/ -/* - C program for floating point sin/cos. - Calls modf. - There are no error exits. - Coefficients are #3370 from Hart & Cheney (18.80D). -*/ -#include - - -#define twoopi 0.636619772367581 -#define p0 0.135788409787738e8 -#define p1 -0.494290810090284e7 -#define p2 0.440103053537527e6 -#define p3 -0.138472724998245e5 -#define p4 0.145968840666577e3 -#define q0 0.864455865292253e7 -#define q1 0.408179225234330e6 -#define q2 0.946309610153821e4 -#define q3 0.132653490878614e3 - -float sinus(float arg, int quad) -{ - float e, f; - int k; - float ysq; - float x,y; - float temp1, temp2; - - x = arg; - if(x<0) { - x = -x; - quad = quad + 2; - } - x = x*twoopi; /*underflow?*/ - if(x>32764){ - y = modff(x,&e); - e = e + quad; - modff(0.25*e,&f); - quad = e - 4*f; - }else{ - k = x; - y = x - k; - quad = (quad + k) & 03; - } - if (quad & 01) - y = 1-y; - if(quad > 1) - y = -y; - - ysq = y*y; - temp1 = ((((p4*ysq+p3)*ysq+p2)*ysq+p1)*ysq+p0)*y; - temp2 = ((((ysq+q3)*ysq+q2)*ysq+q1)*ysq+q0); - return(temp1/temp2); -} - -float _sinf_c(float arg) { - return sinus(arg, 0); -} - -double _sin_c(double) __attribute__((alias("_sinf_c"))); diff --git a/src/libc/sinf.c b/src/libc/sinf.c new file mode 100644 index 000000000..b8544f5c6 --- /dev/null +++ b/src/libc/sinf.c @@ -0,0 +1,64 @@ +/************************************************************************/ +/* */ +/* Copyright (C) 1999-2008 by Zilog, Inc. */ +/* */ +/************************************************************************/ +/* + C program for floating point sin/cos. + Calls modf. + There are no error exits. + Coefficients are #3370 from Hart & Cheney (18.80D). +*/ +#include + +#define two_over_pi 0.636619772367581 +#define p0 0.135788409787738e8 +#define p1 -0.494290810090284e7 +#define p2 0.440103053537527e6 +#define p3 -0.138472724998245e5 +#define p4 0.145968840666577e3 +#define q0 0.864455865292253e7 +#define q1 0.408179225234330e6 +#define q2 0.946309610153821e4 +#define q3 0.132653490878614e3 + +float sinus(float arg, int quad) +{ + float e, f; + int k; + float ysq; + float x,y; + float temp1, temp2; + + x = arg; + if(x<0) { + x = -x; + quad = quad + 2; + } + x = x*two_over_pi; /*underflow?*/ + if(x>32764){ + y = modff(x,&e); + e = e + quad; + modff(0.25*e,&f); + quad = e - 4*f; + }else{ + k = x; + y = x - k; + quad = (quad + k) & 03; + } + if (quad & 01) + y = 1-y; + if(quad > 1) + y = -y; + + ysq = y*y; + temp1 = ((((p4*ysq+p3)*ysq+p2)*ysq+p1)*ysq+p0)*y; + temp2 = ((((ysq+q3)*ysq+q2)*ysq+q1)*ysq+q0); + return(temp1/temp2); +} + +float _sinf_c(float arg) { + return sinus(arg, 0); +} + +double _sin_c(double) __attribute__((alias("_sinf_c"))); diff --git a/src/libc/sin.src b/src/libc/sinf.src similarity index 100% rename from src/libc/sin.src rename to src/libc/sinf.src diff --git a/src/libc/sinh.c b/src/libc/sinh.c deleted file mode 100644 index 079c4c657..000000000 --- a/src/libc/sinh.c +++ /dev/null @@ -1,61 +0,0 @@ -/************************************************************************/ -/* */ -/* Copyright (C)1987-2008 by */ -/* Zilog, Inc. */ -/* */ -/* San Jose, California */ -/* */ -/************************************************************************/ -/* - sinh(arg) returns the hyperbolic sine of its floating- - point argument. - - The exponential function is called for arguments - greater in magnitude than 0.5. - - A series is used for arguments smaller in magnitude than 0.5. - The coefficients are #2029 from Hart & Cheney. (20.36D) - - cosh(arg) is computed from the exponential function for - all arguments. -*/ - -#include - -#define p0 -0.630767364049772e+6 -#define p1 -0.899127202203951e+5 -#define p2 -0.289421135598956e+4 -#define p3 -0.263056321339750e+2 -#define q0 -0.630767364049772e+6 -#define q1 0.152151737879002e+5 -#define q2 -0.173678953558234e+3 - -float _sinhf_c(float arg) { - float temp, argsq; - register int sign; - - sign = 1; - if(arg < 0) { - arg = -arg; - sign = -1; - } - - if(arg > 21.) { - temp = expf(arg)/2; - if (sign>0) - return(temp); - else - return(-temp); - } - - if(arg > 0.5) { - return(sign*(expf(arg) - expf(-arg))/2); - } - - argsq = arg*arg; - temp = (((p3*argsq+p2)*argsq+p1)*argsq+p0)*arg; - temp /= (((argsq+q2)*argsq+q1)*argsq+q0); - return(sign*temp); -} - -double _sinh_c(double, double *) __attribute__((alias("_sinhf_c"))); diff --git a/src/libc/sinhf.c b/src/libc/sinhf.c new file mode 100644 index 000000000..e94fada98 --- /dev/null +++ b/src/libc/sinhf.c @@ -0,0 +1,61 @@ +/************************************************************************/ +/* */ +/* Copyright (C)1987-2008 by */ +/* Zilog, Inc. */ +/* */ +/* San Jose, California */ +/* */ +/************************************************************************/ +/* + sinh(arg) returns the hyperbolic sine of its floating- + point argument. + + The exponential function is called for arguments + greater in magnitude than 0.5. + + A series is used for arguments smaller in magnitude than 0.5. + The coefficients are #2029 from Hart & Cheney. (20.36D) + + cosh(arg) is computed from the exponential function for + all arguments. +*/ + +#include + +#define p0 -0.630767364049772e+6 +#define p1 -0.899127202203951e+5 +#define p2 -0.289421135598956e+4 +#define p3 -0.263056321339750e+2 +#define q0 -0.630767364049772e+6 +#define q1 0.152151737879002e+5 +#define q2 -0.173678953558234e+3 + +float _sinhf_c(float arg) { + float temp, argsq; + register int sign; + + sign = 1; + if(arg < 0) { + arg = -arg; + sign = -1; + } + + if(arg > 21.) { + temp = expf(arg)/2; + if (sign>0) + return(temp); + else + return(-temp); + } + + if(arg > 0.5) { + return(sign*(expf(arg) - expf(-arg))/2); + } + + argsq = arg*arg; + temp = (((p3*argsq+p2)*argsq+p1)*argsq+p0)*arg; + temp /= (((argsq+q2)*argsq+q1)*argsq+q0); + return(sign*temp); +} + +double _sinh_c(double, double *) __attribute__((alias("_sinhf_c"))); diff --git a/src/libc/sinh.src b/src/libc/sinhf.src similarity index 100% rename from src/libc/sinh.src rename to src/libc/sinhf.src diff --git a/src/libc/strtof.c b/src/libc/strtof.c index 7ba41b1a4..912cfdd24 100644 --- a/src/libc/strtof.c +++ b/src/libc/strtof.c @@ -1,10 +1,10 @@ /************************************************************************/ -/* */ -/* Copyright (C)1987-2008 by */ -/* Zilog, Inc. */ -/* */ -/* San Jose, California */ -/* */ +/* */ +/* Copyright (C)1987-2008 by */ +/* Zilog, Inc. */ +/* */ +/* San Jose, California */ +/* */ /************************************************************************/ #include @@ -24,18 +24,17 @@ typedef union F32_pun { /************************************************* * -* strtod - string to double conversion +* strtof - string to float conversion * * Inputs: -* str - pointer to the character string -* endptr - pointer to pointer to char to -* put address of first char past -* the end of the string -- or NULL +* str - pointer to the character string +* endptr - pointer to pointer to char to +* put address of first char past +* the end of the string -- or NULL * Returns: -* the value of the number +* the value of the number * *************************************************/ - /** * @remarks `*str >= '0' && *str <= '9'` is smaller than calls to `isdigit(*str)` * @todo Add support for INF INFINITY NAN NAN(...) diff --git a/src/libc/strtold.c b/src/libc/strtold.c index 219cb1c61..ac4e86b7b 100644 --- a/src/libc/strtold.c +++ b/src/libc/strtold.c @@ -1,12 +1,12 @@ #if 1 /************************************************************************/ -/* */ -/* Copyright (C)1987-2008 by */ -/* Zilog, Inc. */ -/* */ -/* San Jose, California */ -/* */ +/* */ +/* Copyright (C)1987-2008 by */ +/* Zilog, Inc. */ +/* */ +/* San Jose, California */ +/* */ /************************************************************************/ #include @@ -29,12 +29,12 @@ typedef union F64_pun { * strtold - string to long double conversion * * Inputs: -* str - pointer to the character string -* endptr - pointer to pointer to char to -* put address of first char past -* the end of the string -- or NULL +* str - pointer to the character string +* endptr - pointer to pointer to char to +* put address of first char past +* the end of the string -- or NULL * Returns: -* the value of the number +* the value of the number * *************************************************/ diff --git a/src/libc/tan.c b/src/libc/tan.c deleted file mode 100644 index 768d3d96e..000000000 --- a/src/libc/tan.c +++ /dev/null @@ -1,79 +0,0 @@ -/************************************************************************/ -/* */ -/* Copyright (C)1987-2008 by */ -/* Zilog, Inc. */ -/* */ -/* San Jose, California */ -/* */ -/************************************************************************/ -/* - floating point tangent - - A series is used after range reduction. - Coefficients are #4285 from Hart & Cheney. (19.74D) -*/ - -#include -#include - -#define invpi 1.27323954473516 -#define p0 -0.130682026475483e+5 -#define p1 0.105597090171495e+4 -#define p2 -0.155068565348327e+2 -#define p3 0.342255438724100e-1 -#define p4 0.338663864267717e-4 -#define q0 -0.166389523894712e+5 -#define q1 0.476575136291648e+4 -#define q2 -0.155503316403171e+3 - -float _tanf_c(float arg) -{ - float sign, temp, e, x, xsq; - int flag, i; - - flag = 0; - sign = 1.; - if(arg < 0.){ - arg = -arg; - sign = -1.; - } - arg = arg*invpi; /*overflow?*/ - x = modff(arg, &e); - i = e; - switch(i%4) { - case 1: - x = 1. - x; - flag = 1; - break; - - case 2: - sign = - sign; - flag = 1; - break; - - case 3: - x = 1. - x; - sign = - sign; - break; - - case 0: - break; - } - - xsq = x*x; - temp = ((((p4*xsq+p3)*xsq+p2)*xsq+p1)*xsq+p0)*x; - temp = temp/(((1.0*xsq+q2)*xsq+q1)*xsq+q0); - - if(flag == 1) { - if(temp == 0.) { - errno = ERANGE; - if (sign>0) - return(HUGE_VALF); - return(-HUGE_VALF); - } - temp = 1./temp; - } - return(sign*temp); -} - -double _tan_c(double) __attribute__((alias("_tanf_c"))); diff --git a/src/libc/tanf.c b/src/libc/tanf.c new file mode 100644 index 000000000..eff712c77 --- /dev/null +++ b/src/libc/tanf.c @@ -0,0 +1,80 @@ +/************************************************************************/ +/* */ +/* Copyright (C)1987-2008 by */ +/* Zilog, Inc. */ +/* */ +/* San Jose, California */ +/* */ +/************************************************************************/ +/* + floating point tangent + + A series is used after range reduction. + Coefficients are #4285 from Hart & Cheney. (19.74D) +*/ + +#include +#include + +#define four_mul_invpi 1.27323954473516 + +#define p0 -0.130682026475483e+5 +#define p1 0.105597090171495e+4 +#define p2 -0.155068565348327e+2 +#define p3 0.342255438724100e-1 +#define p4 0.338663864267717e-4 +#define q0 -0.166389523894712e+5 +#define q1 0.476575136291648e+4 +#define q2 -0.155503316403171e+3 + +float _tanf_c(float arg) +{ + float sign, temp, e, x, xsq; + int flag, i; + + flag = 0; + sign = 1.; + if(arg < 0.){ + arg = -arg; + sign = -1.; + } + arg = arg*four_mul_invpi; /*overflow?*/ + x = modff(arg, &e); + i = e; + switch(i%4) { + case 1: + x = 1. - x; + flag = 1; + break; + + case 2: + sign = - sign; + flag = 1; + break; + + case 3: + x = 1. - x; + sign = - sign; + break; + + case 0: + break; + } + + xsq = x*x; + temp = ((((p4*xsq+p3)*xsq+p2)*xsq+p1)*xsq+p0)*x; + temp = temp/(((1.0*xsq+q2)*xsq+q1)*xsq+q0); + + if(flag == 1) { + if(temp == 0.) { + errno = ERANGE; + if (sign>0) + return(HUGE_VALF); + return(-HUGE_VALF); + } + temp = 1./temp; + } + return(sign*temp); +} + +double _tan_c(double) __attribute__((alias("_tanf_c"))); diff --git a/src/libc/tan.src b/src/libc/tanf.src similarity index 100% rename from src/libc/tan.src rename to src/libc/tanf.src diff --git a/src/libc/tanh.c b/src/libc/tanh.c deleted file mode 100644 index fd867cb83..000000000 --- a/src/libc/tanh.c +++ /dev/null @@ -1,35 +0,0 @@ -/************************************************************************/ -/* */ -/* Copyright (C)1987-2008 by */ -/* Zilog, Inc. */ -/* */ -/* San Jose, California */ -/* */ -/************************************************************************/ -/* - tanh(arg) computes the hyperbolic tangent of its floating - point argument. - - sinh and cosh are called except for large arguments, which - would cause overflow improperly. -*/ - -#include - -float _tanhf_c(float arg) -{ - float sign; - - sign = 1.; - if(arg < 0.){ - arg = -arg; - sign = -1.; - } - - if(arg > 21.) - return(sign); - - return(sign*sinhf(arg)/coshf(arg)); -} - -double _tanh_c(double) __attribute__((alias("_tanhf_c"))); diff --git a/src/libc/tanhf.c b/src/libc/tanhf.c new file mode 100644 index 000000000..bef020ccd --- /dev/null +++ b/src/libc/tanhf.c @@ -0,0 +1,35 @@ +/************************************************************************/ +/* */ +/* Copyright (C)1987-2008 by */ +/* Zilog, Inc. */ +/* */ +/* San Jose, California */ +/* */ +/************************************************************************/ +/* + tanh(arg) computes the hyperbolic tangent of its floating + point argument. + + sinh and cosh are called except for large arguments, which + would cause overflow improperly. +*/ + +#include + +float _tanhf_c(float arg) +{ + float sign; + + sign = 1.; + if(arg < 0.){ + arg = -arg; + sign = -1.; + } + + if(arg > 21.) + return(sign); + + return(sign*sinhf(arg)/coshf(arg)); +} + +double _tanh_c(double) __attribute__((alias("_tanhf_c"))); diff --git a/src/libc/tanh.src b/src/libc/tanhf.src similarity index 100% rename from src/libc/tanh.src rename to src/libc/tanhf.src diff --git a/src/libc/tgamma.c b/src/libc/tgammaf.c similarity index 76% rename from src/libc/tgamma.c rename to src/libc/tgammaf.c index 589efa91f..b05a006dc 100644 --- a/src/libc/tgamma.c +++ b/src/libc/tgammaf.c @@ -11,19 +11,19 @@ float tgammaf(float x) { /* Gamma function */ if (x == 0.0) { /* Pole Error */ errno = ERANGE; - return 1/x < 0 ? -HUGE_VAL : HUGE_VAL; + return signbit(x) ? -HUGE_VALF : HUGE_VALF; } if (x < 0) { int sign; - static float zero = 0.0; + static float zero = 0.0; float i, f; f = modff(-x, &i); if (f == 0.0) { /* Domain Error */ errno = EDOM; - return zero/zero; + return zero/zero; /* probably better to return NAN here */ } sign = (fmodf(i, 2.0) != 0.0) ? 1 : -1; - return sign * M_PI / (sinf(M_PI * f) * expf(lgammaf(1 - x))); + return (sign * M_PI) / (sinf(M_PI * f) * expf(lgammaf(1 - x))); } return expf(lgammaf(x)); } diff --git a/src/linker.mk b/src/linker.mk index 2e0b0b7dc..d99b6108d 100644 --- a/src/linker.mk +++ b/src/linker.mk @@ -63,8 +63,8 @@ linker_script: $(STATIC_FILES) $(LINKED_FILES) $(SHARED_FILES) $(Q)$(call APPEND_FILES,source ,ce,$(sort $(CE_FILES))) $(Q)$(call APPEND,if HAS_LIBC) $(Q)$(call APPEND_FILES, source ,libc,$(sort $(LIBC_FILES))) + $(Q)$(call APPEND_FILES, source ,softfloat,$(sort $(SOFTFLOAT_FILES))) $(Q)$(call APPEND,end if) $(Q)$(call APPEND,if HAS_LIBCXX) $(Q)$(call APPEND_FILES, source ,libcxx,$(sort $(LIBCXX_FILES))) - $(Q)$(call APPEND_FILES, source ,softfloat,$(sort $(SOFTFLOAT_FILES))) $(Q)$(call APPEND,end if)