Calculating the sine function for ourselves:

  sin(x) =   x / 1!
           - x*x*x / 3!
           + x*x*x*x*x / 5!
           - x*x*x*x*x*x*x / 7!
           + x*x*x*x*x*x*x*x*x / 9!
         and so on, the more terms we have, the more
         accurate the result is

(yes, there is no point in writing the division by 1! as 1!
 is 1. but it makes the repeating pattern clearer)

So we need to add up values over a range controlling the
powers and factorials. As written above the range would be
1 to 9, in steps of 2.

Apart from the sine function that does the adding over a range,
we need three functions:

  double power(const double x, const int p);   for x to the power of p
  double factorial(const int n);               (results too big for ints)
  double plus_or_minus(const int n)            so that plus_or_minus(3) is -1,
                                               plus_or_minus(5) = 1,
                                               plus_or_minus(7) = -1
                                               and so on.

the function that brings everything together would be

  double sine(const double x, const int first_p, const int last_p);

       x is the angle we want the sine of, first_p and last_p are used
       to run through the controlling values 1, 3, 5, 7, 9, ...

and to keep users happy:

  double sine(const double x)
  { return sine(x, 1, 39); }

(it is OK to have many functions with the same name so long as you can
 tell them apart by their parameters only. Either having a different number
 of parameters, or having parameters with clearly different types)

The following is what we got in class on Thursday 18th September, 
including two different versions of the power function and two
different methods for controlling the repetition in the sine function.


#include <library.h>


// IN CLASS I had a particular value for last_p built into the power
// function itself, so there was no third parameter but there was a 
// "result so far" parameter, I'll put both versions here

double power1(const double x, const int p)
{ if (p == 0)
    return 1.0;
  return x * power1(x, p - 1); }


double power2(const double sofar, const double apower, const int p)
{ if (p == 0)
    return sofar;
  if (p % 2 == 0)
    return power2(sofar, apower * apower, p / 2);
  else
    return power2(sofar * apower, apower * apower, p / 2); }

double power2(const double x, const int p)
{ return power2(1.0, x, p); }


double factorial(const int n)
{ if (n == 0)
    return 1.0;
  return n * factorial(n - 1); }


int plus_or_minus(const int n)
{ return 2 - n % 4; }


double sine1(const double x, const int p, const int last_p)
{ if (p > last_p)
    return 0.0;
  return plus_or_minus(p) * power1(x, p) / factorial(p) + sine1(x, p + 2, last_p); }

double sine1(const double x)
{ return sine1(x, 1, 65); }


double sine2(const double sofar, const double x, const int p)
{ const int max_power = 65;
  if (p > max_power)
    return sofar;
  return sine2(sofar + plus_or_minus(p) * power2(x, p) / factorial(p), x, p + 2); }

double sine2(const double x)
{ return sine2(0.0, x, 1); }


void test(const string s, const double x)
{ const double our_answer = sine2(x),
               correct_answer = sin(x);
  cout << s << " = " << x << ": " 
       << our_answer << " should be " << correct_answer
       << " error " << fabs(our_answer - correct_answer)
       << "\n"; }

const double pi = acos(-1.0);    // this is the most accurate calculation on most hardware

int main()
{ test("zero", 0.0);
  test("30 degrees", pi / 6.0);
  test("45 degrees", pi / 4.0);
  test("90 degrees", pi / 2.0); }