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);
}

3 answers

  • answered 2018-02-13 04:09 user3629249

    after making the corrections, as listed in my comments to the question, the proposed code looks like:

    #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.0; // 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( void )
    {
        double number = atof("6");
        double sineResult = sine(number);
        printf("%1.10f", sineResult);
        return (0);
    }
    

    and the resulting output looks like:

    Loop 0:
    6.000 ^ 3.000 = 216.000
    3.000! = 6.000
    Result = -30.000
    
    Loop 1:
    6.000 ^ 5.000 = 7776.000
    5.000! = 120.000
    Result = 34.800
    
    Loop 2:
    6.000 ^ 7.000 = 279936.000
    7.000! = 5040.000
    Result = -20.743
    
    Loop 3:
    6.000 ^ 9.000 = 10077696.000
    9.000! = 362880.000
    Result = 7.029
    
    Loop 4:
    6.000 ^ 11.000 = 362797056.000
    11.000! = 39916800.000
    Result = -2.060
    
    Loop 5:
    6.000 ^ 13.000 = 13060694016.000
    13.000! = 6227020800.000
    Result = 0.037
    
    Loop 6:
    6.000 ^ 15.000 = 470184984576.000
    15.000! = 1307674368000.000
    Result = -0.322
    
    Loop 7:
    6.000 ^ 17.000 = 16926659444736.000
    17.000! = 355687428096000.000
    Result = -0.275
    
    Loop 8:
    6.000 ^ 19.000 = 609359740010496.000
    19.000! = 121645100408832000.000
    Result = -0.280
    
    Loop 9:
    6.000 ^ 21.000 = 21936950640377856.000
    21.000! = 51090942171709440000.000
    Result = -0.279
    
    -0.2793866930
    

  • answered 2018-02-13 17:32 VladP

    Taylor expansion has an error that depends on the argument scope as well as the order of the Taylor expansion. I believe that you have overreached the boundaries of the argument. See here for more examples: www.dotancohen.com/eng/taylor-sine.php

  • answered 2018-02-13 18:49 Marco

    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!