Search code examples
cfactorialtrigonometrytaylor-series

Creating a custom sine function


I tried to create a custom sine function using c and the Taylor Series for calculating sin with 10 terms in the series, but I'm getting the wrong results when I try to find the sine(x) where x > 6.

It works well for -5 < x < 5, but anything out of that range isn't producing the correct results.

I expect sin(10) to return something close to -0.5440, but get 1418.0269775391

I've put everything in a single file so it's easier.

#include <stdio.h>
#include <stdlib.h>

double factorial(double n);
double power(double n, double pow);
double sine(double n);

// This is supposed to all go in a .c file and reference the .h stuff above
// This is the actual implementation of the functions declared above
double factorial(double n) {
    // 0! = 1 so just return it
    if(n == 0) {
        return 1;
    }
    // Recursively call factorial with n-1 until n == 0
    return n * (factorial(n - 1)); 
}


double power(double n, double power) {
    double result = n;
    // Loop as many times as the power and just multiply itself power amount of times
    for(int i = 1; i < power; i++) {
        result = n * result;
    }
    return result;
}

double sine(double n) {
    double result = n;
    double coefficent = 3; // Increment this by 2 each loop
    for(int i = 0; i < 10; i++) { // Change 10 to go out to more/less terms
        double pow = power(n, coefficent);
        double frac = factorial(coefficent);
        printf("Loop %d:\n%2.3f ^ %2.3f = %2.3f\n", i, n, coefficent, pow);
        printf("%2.3f! = %2.3f\n", coefficent, frac);

        // Switch between adding/subtracting
        if(i % 2 == 0) { // If the index of the loop is divided by 2, the index is even, so subtract
            result = result - (pow/frac); // x - ((x^3)/(3!)) - ((x^5)/(5!))...
        } else {
            result = result + (pow/frac); // x - ((x^3)/(3!)) + ((x^5)/(5!))...
        }
        coefficent = coefficent + 2;
        printf("Result = %2.3f\n\n", result);
    }
    return result;
}

// main starting point. This is suppossed to #include "functions.c" which contain the above functions in it
int main(int argc, char** argv) {

    double number = atof(argv[1]); // argv[1] = "6"
    double sineResult = sine(number);
    printf("%1.10f", sineResult);
    return (0);
}

Solution

  • As I already said in Python: Calculate sine/cosine with a precision of up to 1 million digits

    The real Taylor expansion centered in x0 is:

    where Rn is the Lagrange Remainder

    Note that Rn grows fast as soon as x moves away from the center x0.

    Since you are implementing the Maclaurin series (Taylor series centered in 0) and not the general Taylor series, your function will give really wrong results when trying to calculate sin(x) for big values of x.

    So before the for loop in your sine() function you must reduce the domain to at least [-pi, pi]... better if you reduce it to [0, pi] and take advantage of sine's parity.

    To fix your code you'll need fmod() from math.h, so you can do:

    #include <math.h>
    
    // Your code
    
    double sine (double n) {
        // Define PI
        const double my_pi = 3.14159265358979323846;
        // Sine's period is 2*PI
        n = fmod(n, 2 * my_pi);
        // Any negative angle can be brought back
        // to it's equivalent positive angle
        if (n < 0) {
            n = 2 * my_pi - n;
        }
        // Sine is an odd function...
        // let's take advantage of it.
        char sign = 1;
        if (n > my_pi) {
            n -= my_pi;
            sign = -1;
        }
        // Now n is in range [0, PI].
    
        // The rest of your function is fine
    
        return sign * result;
    }
    

    Now if you really hate math.h module, you can implement your own fmod() like this,

    double fmod(double a, double b)
    {
        double frac = a / b;
        int floor = frac > 0 ? (int)frac : (int)(frac - 0.9999999999999999);
        return (a - b * floor);
    }
    

    Try it online!