[comp.os.minix] Trig-function for libfp.a

hall@cod.NOSC.MIL (Robert R. Hall) (08/03/89)

Now that Peter S. Housel has given Minix floating point capabilites,
here is a couple of library module to do trigonometery problems.

Nice work Peter, that was a superior piece of work. Well done.

Robert R. Hall
hall@nosc.mil

echo x - arctan.c
sed '/^X/s///' > arctan.c << '/'
X/*
X		ARCTAN.C
X*/
X
X#define FALSE 0
X#define TRUE !FALSE
X#define PI 3.141592654
X#define HALF_PI 1.570796327
X
Xstatic double 
Xseries(x)
X   double          x;
X
X{
X   static double   coef[8] = {
X			 0.9999993, -0.3332986, 0.1994654, -0.1390853,
X		  9.642004e-2, -5.59098e-2, 2.186123e-2, -4.05406e-3};
X
X   double          x2, ans;
X   int             i;
X
X   x2 = x * x;
X   ans = coef[7];
X   for (i = 6; i >= 0; i--)
X      ans = ans * x2 + coef[i];
X   ans = x * ans;
X   return (ans);
X}
X
Xdouble 
Xatan(x)
X   double          x;
X
X{
X   double          ans;
X
X   if (x > 1.0)
X      ans = -(series(1.0 / x) - HALF_PI);
X   else if (x < -1.0)
X      ans = -(series(1.0 / x) + HALF_PI);
X   else
X      ans = series(x);
X   return (ans);
X}
X
X/*
X			ATAN2
X	   A four quadrant arctan function
X	   returns arctan of y/x in range from pi to -pi
X*/
X
Xdouble 
Xatan2(y, x)
X   double          y, x;
X
X{
X   double          ans;
X   int             x_neg, y_neg;
X
X   if ((x == 0.0) && (y == 0.0))
X      ans = 0.0;
X   else
X   {
X      if (x < 0.0)
X      {
X	 x_neg = TRUE;
X	 x = -x;
X      } else
X	 x_neg = FALSE;
X      if (y < 0.0)
X      {
X	 y_neg = TRUE;
X	 y = -y;
X      } else
X	 y_neg = FALSE;
X      if (y < x)
X	 ans = series(y / x);
X      else
X	 ans = HALF_PI - series(x / y);
X      if (x_neg)
X	 ans = PI - ans;
X      if (y_neg)
X	 ans = -ans;
X   }
X   return (ans);
X}
/
echo x - sincos.c
sed '/^X/s///' > sincos.c << '/'
X/*
X		SINCOS.C
X*/
X
X#define PI 3.141592654
X#define HALF_PI 1.570796327
X#define TWO_PI 6.283185308
X
Xstatic double 
Xseries(x)
X   double          x;
X
X{
X   static double   coef[5] = {
X	     1.0, -0.1666665, 8.333026e-3, -1.980741e-4, 2.601887e-6};
X   double          x2, ans;
X   int             i;
X
X   x2 = x * x;
X   ans = coef[4];
X   for (i = 3; i >= 0; i--)
X      ans = ans * x2 + coef[i];
X   return (ans * x);
X}
X
Xdouble 
Xsin(x)
X   double          x;
X
X{
X
X   if ((x > 1.0e5) || (x < -1.0e5))	/* needed to guarantee that
X					   while loop will terminate
X					   for extremely large angles */
X      x = 0;
X   while (x > PI)
X      x = x - TWO_PI;
X   while (x < -PI)
X      x = x + TWO_PI;
X   if (x > HALF_PI)
X      x = PI - x;
X   if (x < -HALF_PI)
X      x = -(PI + x);
X   return (series(x));
X}
X
Xdouble 
Xcos(x)
X   double          x;
X
X{
X   return (sin(x + HALF_PI));
X}
/
echo x - test_trig.c
sed '/^X/s///' > test_trig.c << '/'
X#include <math.h>
X
Xmain()
X{
X   double          x, y, z, theta, pi, rad2deg, deg2rad;
X
X   pi = 4.0 * atan(1.0);
X   rad2deg = 180.0 / pi;
X   deg2rad = pi / 180.0;
X
X   printf("  Theta     cosine     sine     arctan\n\n");
X   for (theta = -360.0; theta <= 360.0; theta += 15.0)
X   {
X      x = cos(deg2rad * theta);
X      y = sin(deg2rad * theta);
X      z = rad2deg * atan2(y, x);
X      printf(" %7.2f  %8.5f  %8.5f  %7.2f\n", theta, x, y, z);
X   }
X}
/