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

1 files changed, 0 insertions, 251 deletions

diff --git a/libm/double/euclid.c b/libm/double/euclid.c
deleted file mode 100644
index 3a899a6d2..000000000
--- a/libm/double/euclid.c
+++ /dev/null
@@ -1,251 +0,0 @@
-/*							euclid.c
- *
- *	Rational arithmetic routines
- *
- *
- *
- * SYNOPSIS:
- *
- * 
- * typedef struct
- *      {
- *      double n;  numerator
- *      double d;  denominator
- *      }fract;
- *
- * radd( a, b, c )      c = b + a
- * rsub( a, b, c )      c = b - a
- * rmul( a, b, c )      c = b * a
- * rdiv( a, b, c )      c = b / a
- * euclid( &n, &d )     Reduce n/d to lowest terms,
- *                      return greatest common divisor.
- *
- * Arguments of the routines are pointers to the structures.
- * The double precision numbers are assumed, without checking,
- * to be integer valued.  Overflow conditions are reported.
- */
- 
-
-#include <math.h>
-#ifdef ANSIPROT
-extern double fabs ( double );
-extern double floor ( double );
-double euclid( double *, double * );
-#else
-double fabs(), floor(), euclid();
-#endif
-
-extern double MACHEP;
-#define BIG (1.0/MACHEP)
-
-typedef struct
-	{
-	double n; /* numerator */
-	double d; /* denominator */
-	}fract;
-
-/* Add fractions. */
-
-void radd( f1, f2, f3 )
-fract *f1, *f2, *f3;
-{
-double gcd, d1, d2, gcn, n1, n2;
-
-n1 = f1->n;
-d1 = f1->d;
-n2 = f2->n;
-d2 = f2->d;
-if( n1 == 0.0 )
-	{
-	f3->n = n2;
-	f3->d = d2;
-	return;
-	}
-if( n2 == 0.0 )
-	{
-	f3->n = n1;
-	f3->d = d1;
-	return;
-	}
-
-gcd = euclid( &d1, &d2 ); /* common divisors of denominators */
-gcn = euclid( &n1, &n2 ); /* common divisors of numerators */
-/* Note, factoring the numerators
- * makes overflow slightly less likely.
- */
-f3->n = ( n1 * d2 + n2 * d1) * gcn;
-f3->d = d1 * d2 * gcd;
-euclid( &f3->n, &f3->d );
-}
-
-
-/* Subtract fractions. */
-
-void rsub( f1, f2, f3 )
-fract *f1, *f2, *f3;
-{
-double gcd, d1, d2, gcn, n1, n2;
-
-n1 = f1->n;
-d1 = f1->d;
-n2 = f2->n;
-d2 = f2->d;
-if( n1 == 0.0 )
-	{
-	f3->n = n2;
-	f3->d = d2;
-	return;
-	}
-if( n2 == 0.0 )
-	{
-	f3->n = -n1;
-	f3->d = d1;
-	return;
-	}
-
-gcd = euclid( &d1, &d2 );
-gcn = euclid( &n1, &n2 );
-f3->n = (n2 * d1 - n1 * d2) * gcn;
-f3->d = d1 * d2 * gcd;
-euclid( &f3->n, &f3->d );
-}
-
-
-
-
-/* Multiply fractions. */
-
-void rmul( ff1, ff2, ff3 )
-fract *ff1, *ff2, *ff3;
-{
-double d1, d2, n1, n2;
-
-n1 = ff1->n;
-d1 = ff1->d;
-n2 = ff2->n;
-d2 = ff2->d;
-
-if( (n1 == 0.0) || (n2 == 0.0) )
-	{
-	ff3->n = 0.0;
-	ff3->d = 1.0;
-	return;
-	}
-euclid( &n1, &d2 ); /* cross cancel common divisors */
-euclid( &n2, &d1 );
-ff3->n = n1 * n2;
-ff3->d = d1 * d2;
-/* Report overflow. */
-if( (fabs(ff3->n) >= BIG) || (fabs(ff3->d) >= BIG) )
-	{
-	mtherr( "rmul", OVERFLOW );
-	return;
-	}
-/* euclid( &ff3->n, &ff3->d );*/
-}
-
-
-
-/* Divide fractions. */
-
-void rdiv( ff1, ff2, ff3 )
-fract *ff1, *ff2, *ff3;
-{
-double d1, d2, n1, n2;
-
-n1 = ff1->d;	/* Invert ff1, then multiply */
-d1 = ff1->n;
-if( d1 < 0.0 )
-	{ /* keep denominator positive */
-	n1 = -n1;
-	d1 = -d1;
-	}
-n2 = ff2->n;
-d2 = ff2->d;
-if( (n1 == 0.0) || (n2 == 0.0) )
-	{
-	ff3->n = 0.0;
-	ff3->d = 1.0;
-	return;
-	}
-
-euclid( &n1, &d2 ); /* cross cancel any common divisors */
-euclid( &n2, &d1 );
-ff3->n = n1 * n2;
-ff3->d = d1 * d2;
-/* Report overflow. */
-if( (fabs(ff3->n) >= BIG) || (fabs(ff3->d) >= BIG) )
-	{
-	mtherr( "rdiv", OVERFLOW );
-	return;
-	}
-/* euclid( &ff3->n, &ff3->d );*/
-}
-
-
-
-
-
-/* Euclidean algorithm
- *   reduces fraction to lowest terms,
- *   returns greatest common divisor.
- */
-
-
-double euclid( num, den )
-double *num, *den;
-{
-double n, d, q, r;
-
-n = *num; /* Numerator. */
-d = *den; /* Denominator. */
-
-/* Make numbers positive, locally. */
-if( n < 0.0 )
-	n = -n;
-if( d < 0.0 )
-	d = -d;
-
-/* Abort if numbers are too big for integer arithmetic. */
-if( (n >= BIG) || (d >= BIG) )
-	{
-	mtherr( "euclid", OVERFLOW );
-	return(1.0);
-	}
-
-/* Divide by zero, gcd = 1. */
-if(d == 0.0)
-	return( 1.0 );
-
-/* Zero. Return 0/1, gcd = denominator. */
-if(n == 0.0)
-	{
-/*
-	if( *den < 0.0 )
-		*den = -1.0;
-	else
-		*den = 1.0;
-*/
-	*den = 1.0;
-	return( d );
-	}
-
-while( d > 0.5 )
-	{
-/* Find integer part of n divided by d. */
-	q = floor( n/d );
-/* Find remainder after dividing n by d. */
-	r = n - d * q;
-/* The next fraction is d/r. */
-	n = d;
-	d = r;
-	}
-
-if( n < 0.0 )
-	mtherr( "euclid", UNDERFLOW );
-
-*num /= n;
-*den /= n;
-return( n );
-}
-