summaryrefslogtreecommitdiff
path: root/libm/float/sindgf.c
diff options
context:
space:
mode:
authorEric Andersen <andersen@codepoet.org>2001-05-10 00:40:28 +0000
committerEric Andersen <andersen@codepoet.org>2001-05-10 00:40:28 +0000
commit1077fa4d772832f77a677ce7fb7c2d513b959e3f (patch)
tree579bee13fb0b58d2800206366ec2caecbb15f3fc /libm/float/sindgf.c
parent22358dd7ce7bb49792204b698f01a6f69b9c8e08 (diff)
uClibc now has a math library. muahahahaha!
-Erik
Diffstat (limited to 'libm/float/sindgf.c')
-rw-r--r--libm/float/sindgf.c232
1 files changed, 232 insertions, 0 deletions
diff --git a/libm/float/sindgf.c b/libm/float/sindgf.c
new file mode 100644
index 000000000..a3f5851c8
--- /dev/null
+++ b/libm/float/sindgf.c
@@ -0,0 +1,232 @@
+/* sindgf.c
+ *
+ * Circular sine of angle in degrees
+ *
+ *
+ *
+ * SYNOPSIS:
+ *
+ * float x, y, sindgf();
+ *
+ * y = sindgf( x );
+ *
+ *
+ *
+ * DESCRIPTION:
+ *
+ * Range reduction is into intervals of 45 degrees.
+ *
+ * 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 +-3600 100,000 1.2e-7 3.0e-8
+ *
+ * ERROR MESSAGES:
+ *
+ * message condition value returned
+ * sin total loss x > 2^24 0.0
+ *
+ */
+
+/* cosdgf.c
+ *
+ * Circular cosine of angle in degrees
+ *
+ *
+ *
+ * SYNOPSIS:
+ *
+ * float x, y, cosdgf();
+ *
+ * y = cosdgf( x );
+ *
+ *
+ *
+ * DESCRIPTION:
+ *
+ * Range reduction is into intervals of 45 degrees.
+ *
+ * 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;
+
+/* These are for a 24-bit significand: */
+static float T24M1 = 16777215.;
+
+static float PI180 = 0.0174532925199432957692; /* pi/180 */
+
+float sindgf( float xx )
+{
+float x, y, z;
+long j;
+int sign;
+
+sign = 1;
+x = xx;
+if( xx < 0 )
+ {
+ sign = -1;
+ x = -xx;
+ }
+if( x > T24M1 )
+ {
+ mtherr( "sindgf", TLOSS );
+ return(0.0);
+ }
+j = 0.022222222222222222222 * x; /* integer part of x/45 */
+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;
+ }
+
+x = x - y * 45.0;
+x *= PI180; /* multiply by pi/180 to convert to radians */
+
+z = x * x;
+if( (j==1) || (j==2) )
+ {
+/*
+ y = ((( 2.4462803166E-5 * z
+ - 1.3887580023E-3) * z
+ + 4.1666650433E-2) * z
+ - 4.9999999968E-1) * z
+ + 1.0;
+*/
+
+/* measured relative error in +/- pi/4 is 7.8e-8 */
+ y = (( 2.443315711809948E-005 * z
+ - 1.388731625493765E-003) * z
+ + 4.166664568298827E-002) * 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;
+ }
+
+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 cosdgf( float xx )
+{
+register float x, y, z;
+int j, sign;
+
+/* make argument positive */
+sign = 1;
+x = xx;
+if( x < 0 )
+ x = -x;
+
+if( x > T24M1 )
+ {
+ mtherr( "cosdgf", TLOSS );
+ return(0.0);
+ }
+
+j = 0.02222222222222222222222 * 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;
+
+x = x - y * 45.0; /* x mod 45 degrees */
+x *= PI180; /* multiply by pi/180 to convert to radians */
+
+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 );
+}
+