A bit more BeOS compatible math.h functions.

git-svn-id: file:///srv/svn/repos/haiku/trunk/current@6962 a95241bf-73f2-0310-859d-f6bbb57e9c96
This commit is contained in:
Axel Dörfler
2004-03-13 18:48:55 +00:00
parent 1e37a9ac80
commit f7f6185398
4 changed files with 163 additions and 153 deletions
+50 -13
View File
@@ -1,4 +1,5 @@
/* Copyright (c) 1991, 1993 /*-
* Copyright (c) 1991, 1993
* The Regents of the University of California. All rights reserved. * The Regents of the University of California. All rights reserved.
* *
* Redistribution and use in source and binary forms, with or without * Redistribution and use in source and binary forms, with or without
@@ -28,20 +29,19 @@
* LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
* OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
* SUCH DAMAGE. * SUCH DAMAGE.
*
* $FreeBSD: src/lib/libc/i386/gen/isinf.c,v 1.6 1999/08/27 23:59:21 peter Exp $
*/ */
#include <sys/types.h>
#include <math.h> #include <math.h>
/* ToDo: on BeOS/x86, these are defined as __isnan(), and __isinf() - the
* isnan()/isinf() are only macros - we will have to reflect this here.
*/
int isnan(double); // double
int isinf(double);
int int
isnan(double d) __isnan(double d)
{ {
register struct IEEEdp { register struct IEEEdp {
unsigned manl : 32; unsigned manl : 32;
@@ -55,14 +55,51 @@ isnan(double d)
int int
isinf(double d) __isinf(double d)
{ {
register struct IEEEdp { register struct IEEEdp {
u_int manl : 32; unsigned manl : 32;
u_int manh : 20; unsigned manh : 20;
u_int exp : 11; unsigned exp : 11;
u_int sign : 1; unsigned sign : 1;
} *p = (struct IEEEdp *)&d; } *p = (struct IEEEdp *)&d;
return p->exp == 2047 && !p->manh && !p->manl; return p->exp == 2047 && !p->manh && !p->manl;
} }
// #pragma mark -
// float
int
__isnanf(float f)
{
// ToDo: fix implementation!
return __isnan((double)f);
}
int
__isinff(float f)
{
return __isinf((double)f);
}
// #pragma mark -
// long double
int
__isnanl(long double d)
{
return __isnan((double)d);
}
int
__isinfl(long double d)
{
return __isinf((double)d);
}
+45 -15
View File
@@ -33,22 +33,15 @@
* $FreeBSD: src/lib/libc/i386/gen/isinf.c,v 1.6 1999/08/27 23:59:21 peter Exp $ * $FreeBSD: src/lib/libc/i386/gen/isinf.c,v 1.6 1999/08/27 23:59:21 peter Exp $
*/ */
#if defined(LIBC_RCS) && !defined(lint)
static const char rcsid[] = "$FreeBSD: src/lib/libc/i386/gen/isinf.c,v 1.6 1999/08/27 23:59:21 peter Exp $";
#endif /* LIBC_RCS and not lint */
#include <ktypes.h>
#include <math.h> #include <math.h>
/* ToDo: on BeOS/x86, these are defined as __isnan(), and __isinf() - the
* isnan()/isinf() are only macros - we will have to reflect this here.
*/
int isnan(double); // double
int isinf(double);
int int
isnan(double d) __isnan(double d)
{ {
register struct IEEEdp { register struct IEEEdp {
unsigned manl : 32; unsigned manl : 32;
@@ -62,14 +55,51 @@ isnan(double d)
int int
isinf(double d) __isinf(double d)
{ {
register struct IEEEdp { register struct IEEEdp {
u_int manl : 32; unsigned manl : 32;
u_int manh : 20; unsigned manh : 20;
u_int exp : 11; unsigned exp : 11;
u_int sign : 1; unsigned sign : 1;
} *p = (struct IEEEdp *)&d; } *p = (struct IEEEdp *)&d;
return p->exp == 2047 && !p->manh && !p->manl; return p->exp == 2047 && !p->manh && !p->manl;
} }
// #pragma mark -
// float
int
__isnanf(float f)
{
// ToDo: fix implementation!
return __isnan((double)f);
}
int
__isinff(float f)
{
return __isinf((double)f);
}
// #pragma mark -
// long double
int
__isnanl(long double d)
{
return __isnan((double)d);
}
int
__isinfl(long double d)
{
return __isinf((double)d);
}
+66 -56
View File
@@ -83,7 +83,9 @@
* from decimal to binary accurately enough to produce the hexadecimal values * from decimal to binary accurately enough to produce the hexadecimal values
* shown. * shown.
*/ */
#include <math.h>
#include <mathimpl.h> #include <mathimpl.h>
vc(r2p1hi, 2.4142135623730950345E0 ,8279,411a,ef32,99fc, 2, .9A827999FCEF32) vc(r2p1hi, 2.4142135623730950345E0 ,8279,411a,ef32,99fc, 2, .9A827999FCEF32)
@@ -95,72 +97,74 @@ ic(r2p1lo, 1.2537167179050217666E-16 , -53, 1.21165F626CDD5)
ic(sqrt2, 1.4142135623730951455E0 , 0, 1.6A09E667F3BCD) ic(sqrt2, 1.4142135623730951455E0 , 0, 1.6A09E667F3BCD)
#ifdef vccast #ifdef vccast
#define r2p1hi vccast(r2p1hi) # define r2p1hi vccast(r2p1hi)
#define r2p1lo vccast(r2p1lo) # define r2p1lo vccast(r2p1lo)
#define sqrt2 vccast(sqrt2) # define sqrt2 vccast(sqrt2)
#endif #endif
double double
hypot(double x, double y) hypot(double x, double y)
{ {
static double const zero=0; static double const zero = 0;
static double const one=1; static double const one = 1;
static double const small=1.0E-18; /* fl(1+small)==1 */ static double const small = 1.0E-18; /* fl(1+small)==1 */
static int const ibig=30; /* fl(1+2**(2*ibig))==1 */ static int const ibig = 30; /* fl(1+2**(2*ibig))==1 */
double t;
double r;
int exp;
if(finite(x)) { if (finite(x)) {
if(finite(y)) { if (finite(y)) {
x=copysign(x,one); double t;
y=copysign(y,one); double r;
if(y > x) { int e;
t=x; x=y; y=t;
}
if(x == zero) return(zero);
if(y == zero) return(x);
exp= logb(x);
if(exp-(int)logb(y) > ibig ) {
/* raise inexact flag and return |x| */
(void volatile)(one+small);
return(x);
}
/* start computing sqrt(x^2 + y^2) */ x = copysign(x, one);
r=x-y; y = copysign(y, one);
if(r>y) { /* x/y > 2 */ if (y > x) {
r=x/y; t = x; x = y; y = t;
r=r+sqrt(one+r*r); }
} else { /* 1 <= x/y <= 2 */
r/=y; t=r*(r+2.0);
r+=t/(sqrt2+sqrt(2.0+t));
r+=r2p1lo; r+=r2p1hi;
}
r=y/r; if (x == zero)
return(x+r); return zero;
if (y == zero)
return x;
} e = logb(x);
if (e - (int)logb(y) > ibig ) {
/* raise inexact flag and return |x| */
(void volatile)(one + small);
return x;
}
else if(y==y) /* y is +-INF */ /* start computing sqrt(x^2 + y^2) */
return(copysign(y,one)); r = x - y;
else if (r > y) { /* x/y > 2 */
return(y); /* y is NaN and x is finite */ r = x / y;
r = r + sqrt(one + r*r);
} else { /* 1 <= x/y <= 2 */
r /= y;
t = r * (r + 2.0);
r += t / (sqrt2 + sqrt(2.0 + t));
r += r2p1lo;
r += r2p1hi;
}
} else if(x==x) { r = y / r;
return (copysign(x,one)); /* x is +-INF */ return x + r;
} else if(finite(y)) { } else if (y == y) { /* y is +-INF */
return(x); /* x is NaN, y is finite */ return copysign(y, one);
#if !defined(vax)&&!defined(tahoe) } else
} else if(y!=y) { return y; /* y is NaN and x is finite */
return(y); /* x and y is NaN */ } else if (x == x) {
#endif /* !defined(vax)&&!defined(tahoe) */ return copysign(x, one); /* x is +-INF */
} else { } else if (finite(y)) {
return(copysign(y,one)); /* y is INF */ return x; /* x is NaN, y is finite */
} } else if (y != y) {
return y; /* x and y is NaN */
} else
return copysign(y, one); /* y is INF */
} }
/* CABS(Z) /* CABS(Z)
* RETURN THE ABSOLUTE VALUE OF THE COMPLEX NUMBER Z = X + iY * RETURN THE ABSOLUTE VALUE OF THE COMPLEX NUMBER Z = X + iY
* DOUBLE PRECISION (VAX D format 56 bits, IEEE DOUBLE 53 BITS) * DOUBLE PRECISION (VAX D format 56 bits, IEEE DOUBLE 53 BITS)
@@ -174,17 +178,23 @@ hypot(double x, double y)
* cabs(z) = hypot(x,y) . * cabs(z) = hypot(x,y) .
*/ */
// ToDo: complex functions are currently missing - BeOS doesn't have
// a header for them, but libroot.so seems to define all functions
#if 0
double double
cabs(struct complex z) cabs(struct complex z)
{ {
return hypot(z.x,z.y); return hypot(z.x, z.y);
} }
double double
z_abs(struct complex const *z) z_abs(struct complex const *z)
{ {
return hypot(z->x,z->y); return hypot(z->x, z->y);
} }
#endif
/* A faster but less accurate version of cabs(x,y) */ /* A faster but less accurate version of cabs(x,y) */
#if 0 #if 0
+2 -69
View File
@@ -59,19 +59,10 @@
* fmod(x,0), fmod(INF,y) are invalid operations and NaN is returned. * fmod(x,0), fmod(INF,y) are invalid operations and NaN is returned.
* *
*/ */
#if !defined(vax) && !defined(tahoe)
extern int isnan(),finite();
#endif /* !defined(vax) && !defined(tahoe) */
extern double frexp(),ldexp(),fabs();
#ifdef TEST_FMOD
static double
_fmod(x,y)
#else /* TEST_FMOD */
double double
fmod(x,y) fmod(double x, double y)
#endif /* TEST_FMOD */
double x,y;
{ {
int ir,iy; int ir,iy;
double r,w; double r,w;
@@ -93,61 +84,3 @@ double x,y;
} }
return x >= (double)0 ? r : -r; return x >= (double)0 ? r : -r;
} }
#ifdef TEST_FMOD
extern long random();
extern double fmod();
#define NTEST 10000
#define NCASES 3
static int nfail = 0;
static void
doit(x,y)
double x,y;
{
double ro = fmod(x,y),rn = _fmod(x,y);
if (ro != rn) {
(void)printf(" x = 0x%08.8x %08.8x (%24.16e)\n",x,x);
(void)printf(" y = 0x%08.8x %08.8x (%24.16e)\n",y,y);
(void)printf(" fmod = 0x%08.8x %08.8x (%24.16e)\n",ro,ro);
(void)printf("_fmod = 0x%08.8x %08.8x (%24.16e)\n",rn,rn);
(void)printf("\n");
}
}
main()
{
register int i,cases;
double x,y;
srandom(12345);
for (i = 0; i < NTEST; i++) {
x = (double)random();
y = (double)random();
for (cases = 0; cases < NCASES; cases++) {
switch (cases) {
case 0:
break;
case 1:
y = (double)1/y; break;
case 2:
x = (double)1/x; break;
default:
abort(); break;
}
doit(x,y);
doit(x,-y);
doit(-x,y);
doit(-x,-y);
}
}
if (nfail)
(void)printf("Number of failures: %d (out of a total of %d)\n",
nfail,NTEST*NCASES*4);
else
(void)printf("No discrepancies were found\n");
exit(0);
}
#endif /* TEST_FMOD */