diff options
author | Eric Andersen <andersen@codepoet.org> | 2001-11-22 14:04:29 +0000 |
---|---|---|
committer | Eric Andersen <andersen@codepoet.org> | 2001-11-22 14:04:29 +0000 |
commit | 7ce331c01ce6eb7b3f5c715a38a24359da9c6ee2 (patch) | |
tree | 3a7e8476e868ae15f4da1b7ce26b2db6f434468c /libm/float/sinf.c | |
parent | c117dd5fb183afb1a4790a6f6110d88704be6bf8 (diff) |
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
Diffstat (limited to 'libm/float/sinf.c')
-rw-r--r-- | libm/float/sinf.c | 283 |
1 files changed, 0 insertions, 283 deletions
diff --git a/libm/float/sinf.c b/libm/float/sinf.c deleted file mode 100644 index 2f1bb45b8..000000000 --- a/libm/float/sinf.c +++ /dev/null @@ -1,283 +0,0 @@ -/* sinf.c - * - * Circular sine - * - * - * - * SYNOPSIS: - * - * float x, y, sinf(); - * - * y = sinf( x ); - * - * - * - * DESCRIPTION: - * - * Range reduction is into intervals of pi/4. The reduction - * error is nearly eliminated by contriving an extended precision - * modular arithmetic. - * - * Two polynomial approximating functions are employed. - * Between 0 and pi/4 the sine is approximated by - * x + x**3 P(x**2). - * Between pi/4 and pi/2 the cosine is represented as - * 1 - x**2 Q(x**2). - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * IEEE -4096,+4096 100,000 1.2e-7 3.0e-8 - * IEEE -8192,+8192 100,000 3.0e-7 3.0e-8 - * - * ERROR MESSAGES: - * - * message condition value returned - * sin total loss x > 2^24 0.0 - * - * Partial loss of accuracy begins to occur at x = 2^13 - * = 8192. Results may be meaningless for x >= 2^24 - * The routine as implemented flags a TLOSS error - * for x >= 2^24 and returns 0.0. - */ - -/* cosf.c - * - * Circular cosine - * - * - * - * SYNOPSIS: - * - * float x, y, cosf(); - * - * y = cosf( x ); - * - * - * - * DESCRIPTION: - * - * Range reduction is into intervals of pi/4. The reduction - * error is nearly eliminated by contriving an extended precision - * modular arithmetic. - * - * Two polynomial approximating functions are employed. - * Between 0 and pi/4 the cosine is approximated by - * 1 - x**2 Q(x**2). - * Between pi/4 and pi/2 the sine is represented as - * x + x**3 P(x**2). - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * IEEE -8192,+8192 100,000 3.0e-7 3.0e-8 - */ - -/* -Cephes Math Library Release 2.2: June, 1992 -Copyright 1985, 1987, 1988, 1992 by Stephen L. Moshier -Direct inquiries to 30 Frost Street, Cambridge, MA 02140 -*/ - - -/* Single precision circular sine - * test interval: [-pi/4, +pi/4] - * trials: 10000 - * peak relative error: 6.8e-8 - * rms relative error: 2.6e-8 - */ -#include <math.h> - - -static float FOPI = 1.27323954473516; - -extern float PIO4F; -/* Note, these constants are for a 32-bit significand: */ -/* -static float DP1 = 0.7853851318359375; -static float DP2 = 1.30315311253070831298828125e-5; -static float DP3 = 3.03855025325309630e-11; -static float lossth = 65536.; -*/ - -/* These are for a 24-bit significand: */ -static float DP1 = 0.78515625; -static float DP2 = 2.4187564849853515625e-4; -static float DP3 = 3.77489497744594108e-8; -static float lossth = 8192.; -static float T24M1 = 16777215.; - -static float sincof[] = { --1.9515295891E-4, - 8.3321608736E-3, --1.6666654611E-1 -}; -static float coscof[] = { - 2.443315711809948E-005, --1.388731625493765E-003, - 4.166664568298827E-002 -}; - -float sinf( float xx ) -{ -float *p; -float x, y, z; -register unsigned long j; -register int sign; - -sign = 1; -x = xx; -if( xx < 0 ) - { - sign = -1; - x = -xx; - } -if( x > T24M1 ) - { - mtherr( "sinf", TLOSS ); - return(0.0); - } -j = FOPI * x; /* integer part of x/(PI/4) */ -y = j; -/* map zeros to origin */ -if( j & 1 ) - { - j += 1; - y += 1.0; - } -j &= 7; /* octant modulo 360 degrees */ -/* reflect in x axis */ -if( j > 3) - { - sign = -sign; - j -= 4; - } - -if( x > lossth ) - { - mtherr( "sinf", PLOSS ); - x = x - y * PIO4F; - } -else - { -/* Extended precision modular arithmetic */ - x = ((x - y * DP1) - y * DP2) - y * DP3; - } -/*einits();*/ -z = x * x; -if( (j==1) || (j==2) ) - { -/* measured relative error in +/- pi/4 is 7.8e-8 */ -/* - y = (( 2.443315711809948E-005 * z - - 1.388731625493765E-003) * z - + 4.166664568298827E-002) * z * z; -*/ - p = coscof; - y = *p++; - y = y * z + *p++; - y = y * z + *p++; - y *= z * z; - y -= 0.5 * z; - y += 1.0; - } -else - { -/* Theoretical relative error = 3.8e-9 in [-pi/4, +pi/4] */ -/* - y = ((-1.9515295891E-4 * z - + 8.3321608736E-3) * z - - 1.6666654611E-1) * z * x; - y += x; -*/ - p = sincof; - y = *p++; - y = y * z + *p++; - y = y * z + *p++; - y *= z * x; - y += x; - } -/*einitd();*/ -if(sign < 0) - y = -y; -return( y); -} - - -/* Single precision circular cosine - * test interval: [-pi/4, +pi/4] - * trials: 10000 - * peak relative error: 8.3e-8 - * rms relative error: 2.2e-8 - */ - -float cosf( float xx ) -{ -float x, y, z; -int j, sign; - -/* make argument positive */ -sign = 1; -x = xx; -if( x < 0 ) - x = -x; - -if( x > T24M1 ) - { - mtherr( "cosf", TLOSS ); - return(0.0); - } - -j = FOPI * x; /* integer part of x/PIO4 */ -y = j; -/* integer and fractional part modulo one octant */ -if( j & 1 ) /* map zeros to origin */ - { - j += 1; - y += 1.0; - } -j &= 7; -if( j > 3) - { - j -=4; - sign = -sign; - } - -if( j > 1 ) - sign = -sign; - -if( x > lossth ) - { - mtherr( "cosf", PLOSS ); - x = x - y * PIO4F; - } -else -/* Extended precision modular arithmetic */ - x = ((x - y * DP1) - y * DP2) - y * DP3; - -z = x * x; - -if( (j==1) || (j==2) ) - { - y = (((-1.9515295891E-4 * z - + 8.3321608736E-3) * z - - 1.6666654611E-1) * z * x) - + x; - } -else - { - y = (( 2.443315711809948E-005 * z - - 1.388731625493765E-003) * z - + 4.166664568298827E-002) * z * z; - y -= 0.5 * z; - y += 1.0; - } -if(sign < 0) - y = -y; -return( y ); -} - |