Approximating sin(x)

Taylor Series

Probably the one most are familiar with, This is a really understandable way of expressing sin(x) in terms of a polynomial
\[ \sin(x) = \frac{x}{1!} - \frac{x^3}{3!} + \frac{x^5}{5!} \cdots \]In the above graph, you can visually compare the accuracy of the approximation with the sine function for increasing number of terms.

It shows that even for 4 or 3 terms, you can get a very viable approximation for \(x \in [-\pi/2,\pi/2]\)

In the way I understand it, The polynomial aims to match the nth derivative of the function with its own derivative, at a certain point
From its error graph, we can also see that its error sharply increases at the ends
This nature of the graph remains even when the interval is cut to \([-\pi/8,\pi/8]\)

Taking the general form of the polynomial to be \(P(x) = C_0x^0 + C_1x^1 + C_2x^2 \cdots \)

Bhaskara's Approximation

This is a formula for calculating sin(x) discovered by 7th centuary Indian mathematician Baskara,
$$\sin(\theta) \approx \frac{16x(\pi-x)}{5\pi^2-4x(\pi-x)}$$ And it's so accurate in the range [0,pi] that its graph is coincident with sin(x).

The plot above is of the error of this function w/ sin(x), which as you can infer is pretty low

And it should be faster than calculating using Taylor's series

A derivation for the function
Wikipedia Article

Minmax Polynomial

A minmax polynomail gives a approximation of a function while minimizing the maximum error as much as possible.

Exactly how its obtained via the Remez Algorithm, is something I haven't really understood yet so I can't say anything more

An implementation in C is given here:

static double minimaxsin(double x)
{
    static const
    double a0 =  1.0,
           a1 = -1.666666666640169148537065260055e-1,
           a2 =  8.333333316490113523036717102793e-3,
           a3 = -1.984126600659171392655484413285e-4,
           a4 =  2.755690114917374804474016589137e-6,
           a5 = -2.502845227292692953118686710787e-8,
           a6 =  1.538730635926417598443354215485e-10;
    double x2 = x * x;
    /*
        a0x + a1x^3 + a2x^5 + a3x^7...
    */
    return x * (a0 + x2 * (a1 + x2 * (a2 + x2
             * (a3 + x2 * (a4 + x2 * (a5 + x2 * a6))))));
}	

/*
 *      3. sin(x) is approximated by a polynomial of degree 13 on
 *         [0,pi/4]
 *                               3            13
 *              sin(x) ~ x + S1*x + ... + S6*x
 *         where
 *
 *      |sin(x)         2     4     6     8     10     12  |     -58
 *      |----- - (1+S1*x +S2*x +S3*x +S4*x +S5*x  +S6*x   )| <= 2
 *      |  x                                               |
 *
 *      4. sin(x+y) = sin(x) + sin'(x')*y
 *                  ~ sin(x) + (1-x*x/2)*y
 *         For better accuracy, let
 *                   3      2      2      2      2
 *              r = x *(S2+x *(S3+x *(S4+x *(S5+x *S6))))
 *         then                   3    2
 *              sin(x) = x + (S1*x + (x *(r-y/2)+y))
*/

static const double
S1  = -1.66666666666666324348e-01, /* 0xBFC55555, 0x55555549 */
S2  =  8.33333333332248946124e-03, /* 0x3F811111, 0x1110F8A6 */
S3  = -1.98412698298579493134e-04, /* 0xBF2A01A0, 0x19C161D5 */
S4  =  2.75573137070700676789e-06, /* 0x3EC71DE3, 0x57B1FE7D */
S5  = -2.50507602534068634195e-08, /* 0xBE5AE5E6, 0x8A2B9CEB */
S6  =  1.58969099521155010221e-10; /* 0x3DE5D93A, 0x5ACFD57C */

double __sin(double x, double dx, int is_dx_not_zero)
{
	double_t x2,x3,x4,r;

	x2 = x*x;
	x4 = x2*x2;
	r = S2 + x2*(S3 + x2*S4) + x2*x4*(S5 + x2*S6);
	x3 = x2*x;
	if (is_dx_not_zero == 0)
        // This unwraps to a taylor expansion of 6 terms
		return x + x3*(S1 + x2*r);
	else
        
		return x - ((x2*(0.5*dx - x3*r) - dx) - x3*S1);
}	

Conclusion

The results from my very naive implementation are as follows, This is without any effort made to optimize the code and most is copied from here and there so the comparison isn't entirely fair,
I also added the SIMD AVX sin function from libmvec for fun too :)

Standard Library:
        took 3266 ms
Taylor Series (7T):
        took 1212 ms
Bhaskara Formula:
        took 1127 ms
Minmax Polynomial (7T):
        took 1210 ms
Libmvec (4xD):
        took 565 ms