Taylor series in C with functions by dividing first then multiplying

Question

I made a program that calculates the respective Taylor series, but the professor wants us to change it up by dividing first then multiplying. This is my code and I'm using the divide() for the first part. My power() and fact() are just references to my old code.

#include <stdio.h>
#include <math.h>

int main()
{
double  x, mySin(), myCos(), myExp();
char    more;

do      {
    printf ("\n\t\t\t\t\tInput x: ");
    scanf  ("%lf", &x);
    printf ("\n\t\t\t\t\t\t\tLibraryResult  \t MyResult");
    printf ("\n\t\t\tSin(%6.2f)\t      %9.6f\t\t%9.6f", x, sin(x), mySin(x));
    printf ("\n\t\t\tCos(%6.2f)\t      %9.6f\t\t%9.6f", x, cos(x), myCos(x));
    printf ("\n\t\t\tExp(%6.2f)\t      %9.6f\t\t%9.6f", x, exp(x), myExp(x));
    printf ("\n\n\t\t\t\t\tDo more (Y/N) ?");
    scanf  ("%s", &more);
} while (more == 'y' || more == 'Y');
}

double mySin(double x)
{
double sum = 0., divide();
int    i, sign = 1;

for (i = 0; i < 30; i++, sign = - sign)
    sum = sum + sign * divide(x) ; //power(x, 2 * i + 1) / fact(2 * i + 1);

return sum;
}

double myCos(double x)
{
double sum = 0., divide();
int    i, sign = 1;

for (i = 0; i < 30; i++, sign = - sign)
    sum = sum + sign * divide(x);//power(x, 2 * i) / fact(2 * i);

return sum;

}

double myExp(double x)
{
double sum = 0., divide();
int i;

for (i= 0; i < 30; i++)
    sum = sum + divide(x); //power(x, i) / fact(i);


return sum;
}


double divide(int n, double x)
{
int i;
double div = 1.;

for (i = 1; i < n; i++)
    div =  x / i;

return div;

}

/*double fact(int n)
{
int i;
double prod = 1.;

for (i = 1; i <= n; i++)
    prod = prod * i;

return prod;
}

double power (double x, int n)
{
int i;
double prod = 1.;

for (i = 0; i < n; i++)
    prod = prod * x;

return prod;
}
*/

I've been messing around with the code and when in the divide(), using i < n yields no result when I enter x. However, using i < 30 does, but the calculations are wrong. I'm trying to figure this out on my own so just a pointer would be helpful! Thanks!

Answer 1

I figured I would drop you an answer, since your issue is systemic.

First, pull the functions out into forward declarations. No more:

double  x, mySin(), myCos(), myExp();

Instead put this above int main():

double divide(int n, double x);
double mySin(double x);
double myCos(double x);
double myExp(double x);

This also provides you a great opportunity to document what these do. This will actually solve your problem! (Yes, even though the documentation is just a comment):

// Returns x^n / n!
double divide(int n, double x);

Looking at that comment, suddenly mySin makes no sense:

sum = sum + sign * divide(x); // this is sum += sign * x^? / ?!

So, correcting mySin would look something like this: (note that pulling logic out of the loop statement is generally a good idea)

// The series for sin(x) starts with x^1/1 and continues by increasing by 2
for (i = 1; i < 30; i += 2) {
    sum += sign * divide(i,x); // += saves you from writing sum + ...
    sign = -sign;
}

This should produce the correct output, but it doesn't. Unfortunately you also have a problem in your divide(). Keeping its documentation in mind, lets take a look at it:

double divide(int n, double x)
{
int i;
double div = 1.;

for (i = 1; i < n; i++)
    div =  x / i;

return div;

}

This is very close to correct, its just that the equals sign isn't quite right. We want to multiply by x/i each iteration.

for (i = 1; i < n; i++)
    div *= x / i;

Is much closer, but lets count how many times that will loop. Its n-1 times (if we started at 0 it would be n times, but we'd have a divide by 0 error). Let's fix that:

for (i = 1; i <= n; i++)
    div *= x / i;

Now this is technically doing division at the same time as multiplication, which will produce the most accurate results (by keeping div near 'normal' numbers). If you divide first you will start using really small numbers and if you multiply first you will start using really big numbers. Either way you will lose accuracy (since double's are most accurate on moderately sized numbers). If you do split that loop apart it turns into:

divide loop:

for (i = 1; i <= n; i++)
    div /= i;

multiply loop:

for (i = 1; i <= n; i++)
    div *= x;

For most inputs, it won't make much difference at all. However if you have 'extreme' values you will end up killing your doubles. For example, when multiplying first, if you use x=1000, you'll out-run the space of the double before you actually start dividing. The result will be nan, since the double just can't handle it. Similarly if you pump up the iterations (instead of 30, do 30000000) it will obviously take much longer (about 100000 times longer of course), but if you also divide first you'll end up zeroing out your double before you get to multiplying. As such you will get a result of 0 (Note: I did not test this, this is theory). It is also entirely possible that the compiler sees the optimization opportunity and takes it for you.