Cosine - Maclaurin Series

This is a discussion on Cosine - Maclaurin Series within the C Programming forums, part of the General Programming Boards category; So I decided I wanted to make a C program to converge to the value the function cosine(x) give for ...

  1. #1
    Registered User carrotcake1029's Avatar
    Join Date
    Apr 2008
    Posts
    398

    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("&#37;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;
    }

  2. #2
    Registered User
    Join Date
    Dec 2006
    Location
    Canada
    Posts
    3,183
    without reading the code, double offers much greater precision than float.

  3. #3
    Registered User
    Join Date
    Dec 2006
    Location
    Canada
    Posts
    3,183
    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 &#37; 2) 
              return power_half_exp*power_half_exp*base;
         else
              return power_half_exp*power_half_exp;
    }

  4. #4
    Registered User carrotcake1029's Avatar
    Join Date
    Apr 2008
    Posts
    398
    Quote Originally Posted by cyberfish View Post
    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.

  5. #5
    and the Hat of Guessing tabstop's Avatar
    Join Date
    Nov 2007
    Posts
    14,185
    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.

  6. #6
    Registered User
    Join Date
    Dec 2006
    Location
    Canada
    Posts
    3,183
    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]
    Last edited by cyberfish; 12-05-2008 at 08:07 PM.

  7. #7
    Registered User
    Join Date
    Dec 2006
    Location
    Canada
    Posts
    3,183
    You can also workaround the overflow problem by using a double instead, at the expense of precision.

  8. #8
    Algorithm Dissector iMalc's Avatar
    Join Date
    Dec 2005
    Location
    New Zealand
    Posts
    6,305
    Quote Originally Posted by carrotcake1029 View Post
    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!
    Last edited by iMalc; 12-06-2008 at 12:32 PM.
    My homepage
    Advice: Take only as directed - If symptoms persist, please see your debugger

    Linus Torvalds: "But it clearly is the only right way. The fact that everybody else does it some other way only means that they are wrong"

  9. #9
    and the hat of wrongness Salem's Avatar
    Join Date
    Aug 2001
    Location
    The edge of the known universe
    Posts
    32,558
    > 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;
    If you dance barefoot on the broken glass of undefined behaviour, you've got to expect the occasional cut.
    If at first you don't succeed, try writing your phone number on the exam paper.
    I support http://www.ukip.org/ as the first necessary step to a free Europe.

  10. #10
    Registered User carrotcake1029's Avatar
    Join Date
    Apr 2008
    Posts
    398
    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("&#37;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
    Last edited by carrotcake1029; 12-06-2008 at 07:29 PM.

  11. #11
    Algorithm Dissector iMalc's Avatar
    Join Date
    Dec 2005
    Location
    New Zealand
    Posts
    6,305
    Your factorial calculation (incremental or otherwise) has completely disappeared!
    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.
    My homepage
    Advice: Take only as directed - If symptoms persist, please see your debugger

    Linus Torvalds: "But it clearly is the only right way. The fact that everybody else does it some other way only means that they are wrong"

  12. #12
    King of the Internet Fahrenheit's Avatar
    Join Date
    Oct 2001
    Posts
    128
    Quote Originally Posted by iMalc View Post
    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
    Last edited by Fahrenheit; 12-07-2008 at 09:27 AM.

  13. #13
    Algorithm Dissector iMalc's Avatar
    Join Date
    Dec 2005
    Location
    New Zealand
    Posts
    6,305
    Quote Originally Posted by Fahrenheit View Post
    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.
    My homepage
    Advice: Take only as directed - If symptoms persist, please see your debugger

    Linus Torvalds: "But it clearly is the only right way. The fact that everybody else does it some other way only means that they are wrong"

Popular pages Recent additions subscribe to a feed

Similar Threads

  1. cosine series with out using math.h
    By indrajit_muk in forum C Programming
    Replies: 5
    Last Post: 12-16-2008, 07:17 PM
  2. Cosine fucntion and infinite loop.
    By youareafever in forum C Programming
    Replies: 2
    Last Post: 11-07-2008, 03:45 AM
  3. Implement of a Fast Time Series Evaluation Algorithm
    By BiGreat in forum C Programming
    Replies: 7
    Last Post: 12-04-2007, 01:30 AM
  4. Problem with implementing Liebniz series
    By C_Necessity in forum C Programming
    Replies: 6
    Last Post: 06-15-2005, 12:39 PM
  5. taylor series using only stdio.h
    By faruque in forum C Programming
    Replies: 1
    Last Post: 02-13-2003, 11:50 AM

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21