Cosine - Maclaurin Series

So I decided I wanted to make a C program to converge to the value the function cosine(x) give for any given x. Having used the Maclaurin Series in the past for Calculus, I knew this should be an easy implementation.

Using Maclaurin series in this case, you do a summation from 0 to infinity. Well obviously you can't make it go to infinity, so I want to iterate as many times for accuracy as possible.

Well it seems that if I set my iterations past 10, my numbers go real wacky. I need to set it higher because if I run 0 through cosine(0), I naturally get 1. But if I run 2pi through, I get 0.996515. Can anyone see if I am overflowing anything or doing something wrong?

Also, I know I am no perfect programmer, but if you will critique, please also help me with my original problem, because this is no final draft.

Code:

#include <stdio.h> float power (float base, float exponent); unsigned long long factorial (int input); float cosine (float input); float pi = 3.141592; int main (void) { printf("%f\n", cosine(2*pi)); return 0; } float power (float base, float exponent) { int i; float product = 1.0; for (i = 0; i < exponent; i++) product *= base; return product; } unsigned long long factorial (int input) { int i; unsigned long long product = 1; for (i = 0; i < input; i++) product *= input - i; return product; } float cosine (float input) { unsigned int i; float sum = 0.0; for (i = 0; i < 10; i++) sum += power(-1, i) * power(input, 2*i) / factorial(2*i); return sum; }

without reading the code, double offers much greater precision than float.

after reading your code (btw when I said without reading... I meant without reading in detail :)), why have the exponent of your power function as a float if it's going to be an integer all the time?

Also, the power function can be more efficiently implemented as a recursive function (someone on the forum pointed it out to me a few months ago)

Code:

double power(double base, int exp) { if (exp == 0) return 1; if (exp == 1) return base; double power_half_exp = power(base, exp / 2); if (exp % 2) return power_half_exp*power_half_exp*base; else return power_half_exp*power_half_exp; }

Quote:

Originally Posted by cyberfish

why have the exponent of your power function as a float if it's going to be an integer all the time?

So whats pi^2? An int?

And I am not looking for efficiency.

I tried changing to double, and nothing changed.

20! may overflow, depending on long long.

It would probably be better to get from one term to the next (make it negative, multiply by x*x, divide by i+1 and i+2), as that will be more accurate.

Ah okay nvm, remembered the algorithm wrong.

I would change everything to double nevertheless.

[edit]
Oh and, your power function only works for integer exponents (that's why I assumed it's only going to take int exponents). Try power(1.5, 1.5).

The standard library has a function "pow()" for that in math.h
[/edit]

You can also workaround the overflow problem by using a double instead, at the expense of precision.

Quote:

Originally Posted by carrotcake1029

So I decided I wanted to make a C program to converge to the value the function cosine(x) give for any given x. Having used the Maclaurin Series in the past for Calculus, I knew this should be an easy implementation.

Using Maclaurin series in this case, you do a summation from 0 to infinity. Well obviously you can't make it go to infinity, so I want to iterate as many times for accuracy as possible.

Well it seems that if I set my iterations past 10, my numbers go real wacky. I need to set it higher because if I run 0 through cosine(0), I naturally get 1. But if I run 2pi through, I get 0.996515. Can anyone see if I am overflowing anything or doing something wrong?

Also, I know I am no perfect programmer, but if you will critique, please also help me with my original problem, because this is no final draft.

The approximation of trig functions via successive polynomials only works over the range 0 to pi/2 iirc. You need to map the value you input into that range. Start by taking the fmod of the value with 2pi, then fold it acording to the function you're approximating from there, and negate the result if necessary. I.e. for cos after the fmod, you check which quadrant it is in and subtract it from pi or whatever in order to map the value to the range o to pi/2. Even the built-in sin and cos do this folding.

You should also use M_PI rather than providing a less accurate version of this constant yourself. Just include the right header file.

Don't forget that once you've calculated n!, you can get (n+1)! by simply multiplying n! by (n+1). I.e. You don't need to throw away the previous factorial you calculated and start over!

> power(input, 2*i) / factorial(2*i);
If you think about it (look at the expression), you should be able to incrementally calculate the next term without having to calculate one increasingly large number, to be divided by another increasingly large number.

> power(-1, i)
If you just want to toggle between -1 and +1, just do
toggle = -toggle;

Thank you for all your replies. Salem, thanks for bringing me to that realization. I couldn't see it until I worked it out on paper. Here is the final code. Thanks.

Post any potential problems you see. Thanks.

Code:

#include <stdio.h> #include <math.h> double cosine (double input); int main (void) { printf("%f\n", cosine(M_PI)); return 0; } double cosine (double input) { int i, j, toggle = 1; double angle, product, sum = 0.0; angle = input; while (angle >= (2*M_PI)) angle -= (2*M_PI); for (i = 0; i < 100; i++) { product = 1.0; for (j = (2*i); j > 0 ; j--) product *= angle / j; sum += product * toggle; toggle = -toggle; } return sum; }

Edit: How can I avoid it producing -0.00000 for 3pi/2? lol

Your factorial calculation (incremental or otherwise) has completely disappeared! :eek:
You also need to do more folding of the input angle to the range 0 .. pi/2. The first step is to map it to 0 .. 2pi, but after that you have to fold it furthur, but instead of just modding it, you have to take into account where the points of symmetry of the curve are, reflecting either the input or the output somehow in the x-axis or y-axis as needed. Draw yourself a diagram and stare at it to figure this out.

Quote:

Originally Posted by iMalc

The approximation of trig functions via successive polynomials only works over the range 0 to pi/2 iirc.

taylor series for sin and cos are entire

Quote:

Originally Posted by Fahrenheit

taylor series for sin and cos are entire

I stand corrected, my statement was too general.
In either case though, you're always going to get a more accurate result when folding the range towards the neighbourhood of zero, given a finite number of terms.