```// Uses first 10 terms of Taylor Series double my_cos(double x) { double sum = 0; for (int n = 0; n < 10; n++) { sum += pow(-1, n) * pow(x, 2 * n) / factorial(2 * n); } return sum; }```
double cos2(double x, int n) {
double numerator = 1.0;
double denominator = 1.0;
double pi = 3.14159265358979323846;
int i;
for (i = 1; i <= n; ++i) {
numerator *= pow((x - n * pi / 2.0), 2);
denominator *= pow((n * pi / 2.0), 2);
}
return numerator / denominator;
}(Edit for format)
I’ve never worked directly on floating point libraries like this, but having dealt with similar stuff, I suspect if this were in production it would end up peppered with IFDEF to support various implementation differences across hardware types, even if it is compliant with eg IEEE standards.
Taylor series are essentially “Give me a polynomial p(x). At x=a, I want the value to be b, the first derivative to be c, the second derivative to be d, etc. etc.”. (Usually, a=0, also known as a MacLaurin series.) This gives you a _great_ approximation in the immediate vicinity of a, but usually a pretty horrible one the further away you get; in a sense, you spend all your approximation energy on that single point.
Minimax approximations are “Give me a polynomial p(x). Over the area x=[a,b], I want it to never stray far from f(x); give me the one that minimizes the worst-case distance”. (The distance function is usually either the absolute error |p(x)-f(x)|, or the relative error |p(x)-f(x)|/|f(x)|, but you can have an arbitrary weight function in some implementations.) So you consider everything in that area to be important, not just one point.
Taylor series are simple to calculate (which can help if you need to calculate the coefficients on the-fly, or by hand), good for analysis, and easy to understand, so most engineers will have heard of them in calculus. The Remez algorithm, on the other hand, essensially requires a computer, and is significantly more complicated to understand and implement. (Remez is not the only way to generate minimax polynomials; Chebychev polynomials come close, and I believe Sollya uses some other algorithm.) Both will converge to the actual function giving infinite terms, but minimax polynomials will be _much_ closer in most real-life situations. Also, if you are really pressed for it, you can have rational minimax polynomials (p(x)/q(x)) that approximate tricky functions _really_ well even at low degree, but of course at the cost of a division.
There's no real tricks or anything to this code, it's just a polynomial evaluation of a Taylor series (which you probably learned, and subsequently forgot, from calculus), with some steps done using techniques for getting higher precision out of a floating-point value, and (in some other functions in libm, not here) occasional use of the bit representation of floating-point values. It's standard techniques, but in a field that is incredibly niche and where there is very little value in delving into if you don't have a need for it.