Search code examples
c++doubleprecisiontaylor-series

C++ - Estimating cos(x) using Taylor Series Approximation


To get some more practice in C++, I decided to do some basic math functions without the aid of the math library. I've made a power and factorial function and they seem to work well. However, I'm having lots of problems regarding my Taylor Series cosine function.

Wikipedia Cosine Taylor Series

It outputs a good approximation at cos(1), cos(2), and begins losing precision at cos(3) and cos(4). Beyond that, its answer becomes completely wrong. The following are results from ./a.out

Input an angle in radians, output will be its cosine
1
Output is: 0.540302

Input an angle in radians, output will be its cosine
2
Output is: -0.415873

Input an angle in radians, output will be its cosine
3
Output is: -0.974777

Input an angle in radians, output will be its cosine
4
Output is: -0.396825                       <-------------Should be approx. -0.654

Input an angle in radians, output will be its cosine
5
Output is: 2.5284                          <-------------Should be approx.  0.284

Here is the complete source code:

#include <iostream>
#include <iomanip>

using std::cout;
using std::cin;
using std::endl;

int factorial(int factorial_input) {

    int original_input = factorial_input;
    int loop_length = factorial_input - 1;

    if(factorial_input == 1 || factorial_input == 0) {

        return 1;
    }

    for(int i=1; i != loop_length; i++) {

        factorial_input = factorial_input - 1;

        original_input = original_input * factorial_input;

    }

    return original_input;
}

double power(double base_input, double exponent_input) {

    double power_output = base_input;

    if(exponent_input == 0) {

        return 1;
    }

    if(base_input == 0) {

        return 0;
    }

    for(int i=0; i < exponent_input -1; i++){

        power_output = power_output * base_input;

    }
    return power_output;

}

double cos(double user_input) {

    double sequence[5] = { 0 };  //The container for each generated elemement.
    double cos_value = 0;        //The final output.
    double variable_x = 0;       //The user input x, being raised to the power 2n

    int alternating_one = 0;     //The (-1) that is being raised to the nth power,so switches back and forth from -1 to 1 
    int factorial_denom = 0;     //Factorial denominator (2n)!
    int loop_lim = sizeof(sequence)/sizeof(double);            //The upper limit of the series (where to stop), depends on size of sequence. Bigger is more precision.

    for(int n=0; n < loop_lim; n++) {

        alternating_one = power(-1, n);
        variable_x = power(user_input, (n*2));
        factorial_denom = factorial((n*2));

        sequence[n] =  alternating_one * variable_x / factorial_denom;
        cout << "Element[" << n << "] is: " << sequence[n] << endl;     //Prints out the value of each element for debugging.
    }

    //This loop sums together all the elements of the sequence.
    for(int i=0; i < loop_lim; i++) {                   

        cos_value = cos_value + sequence[i];

    }

    return cos_value;
}

int main() {

    double user_input = 0;
    double cos_output;

    cout << "Input an angle in radians, output will be its cosine" << endl;
    cin >> user_input;
    cos_output = cos(user_input);

    cout << "Output is: " << cos_output << endl;

}

At five iterations, my function should maintain accuracy until after around x > 4.2 according to this graph on Desmos:

Desmos Graph

Also, when I set the series up to use 20 iterations or more (it generates smaller and smaller numbers which should make the answer more precise), the elements start acting very unpredictable. This is the ./a.out with the sequence debugger on so that we may see what each element contains. The input is 1.

Input an angle in radians, output will be its cosine
1
Element[0] is: 1
Element[1] is: -0.5
Element[2] is: 0.0416667
Element[3] is: -0.00138889
Element[4] is: 2.48016e-05
Element[5] is: -2.75573e-07
Element[6] is: 2.08768e-09
Element[7] is: -7.81894e-10
Element[8] is: 4.98955e-10
Element[9] is: 1.11305e-09
Element[10] is: -4.75707e-10
Element[11] is: 1.91309e-09
Element[12] is: -1.28875e-09
Element[13] is: 5.39409e-10
Element[14] is: -7.26886e-10
Element[15] is: -7.09579e-10
Element[16] is: -4.65661e-10
Element[17] is: -inf
Element[18] is: inf
Element[19] is: -inf
Output is: -nan

Can anyone point out what things I'm doing wrong and what I should be doing better? I'm new to C++ so I still have a lot of misconceptions. Thank you so much for taking the time to read this!


Solution

  • You have the following problems:


    In the graph you are showing in the picture k is included in the sum, while you are excluding it in your code. Therefore k=5 in the Desmos graph is equal to double sequence[6] = { 0 } in your code.

    This fixes the output for user_input = 4.

    For user_input = 5 you can then compare to the graph to see that it gives a similar result as well (which is already far off of the true value)


    Then you will have bugs for larger number of terms, because the factorial function outputs int, but the factorial grows so quickly that it will go out-of-range of the values int can hold quickly and also quickly out-of-range of any integer type. You should return double and let original_input be double as well, if you want to support a somewhat (though not much) larger input range.


    In power you take the exponent as double, but work with it as if it was an integer. In particular you use it for the limit of loop iterations. That will only work correctly as long as the values are small enough to be exactly representable by double. As soon as the values become larger, the number of loop iterations will become inexact.

    Use int as second parameter to power instead.


    If one were to implement cos with this approach, one would normally use cos symmetry first, to reduce the range to something smaller, e.g. [0,pi/2] first, by using e.g. that cos(x + 2pi) = cos(x) and cos(x+pi) = - cos(x) and cos(-x) = cos(x), etc.