More discussions with Ranko re: Atan follow.
First I asked Ranko how he figured out the original function:
"The approximation I suggested x/(1+0.28*x^2) is very close to the best rational approximation of the arctan function. But what is the best rational approximation? A rational function x/(1 + ax^2) is the best approximation if the constant a is selected so that the difference arctan(x) - x/(1+ax^2) is equi-oscillating on [-1,1]
If you draw the graph of arctan(x) - x/(1+0.28 * x^2) you will notice that it has alternatingly a minimum of -0.00415 at -1, then a positive maximum of 0.00488 between -1 and 0, then a negative minimum of -0.00488 between 0 and 1 and finally, a positive maximum of 0.00415 at the point 1.
For a rational function to be the best approximation to the arctan function, these four values must be equal.
When you carry out the necessary computation you find out that
for a = 0.280872 these four values are
-0.0046800276, +0.0046800931,-0.0046800931,+0.004680276
These numbers differ in the 8th decimal place!
This is, of course, only important to mathematicians.
The value a = 0.28 is accurate enough for all computer related applications. If |x| <= 1 we can use the approximation x/(1 + 0.28 * x^2), For |x|>=1 we can use the approximation pi/2 - x/(x^2 + 0.28). In both cases the error is smaller than 0.005!
It is quite possible, however, that in some applications an
accuracy of 6 decimal digits or more
is needed. In that case we would use the function (x + a*x^3)/(1+b*x^2
+ c*x^4)
or something similar."
I then wondered about the speed of the rational approxinmation compared to the table lookup method, and sent Ranko the following code:
I was somewhat suprised at the results of the speed tests Ranko carried out. The rational approximation method actually beats the table lookup method in terms of speed!
The times for 1 million calls on Ranko's machine are:
rational approximation | 88 ticks |
table lookup | 112 ticks |
MathLib's atan | 156 ticks |
The maximum error of each method in degrees:
rational approximation | 0.0004 degrees |
table lookup | 2.5 degrees |
MathLib's atan | Accuracy not known, but considered accurate |
We know that the table lookup methods accuracy can be improved with interpolation, but quite frankly why bother?
r(x) = (x + 0.0.43157974*x^3)/(1 + 0.76443945*x^2 + 0.05831938*x^4)
For -1 <= x <= 1 we have
|arctan(x) - r(x)| < 0.00000642
This means clearly that r(x) gives at least 4 correct decimal digits of arctan(x) for any x in the interval [-1,1]. For example, arctan(1) = pi/4 = 0.785398163 and r(1) = 0.785391746.
Notice that the first three digits of coefficients in both approximations coincide. But the difference is profound. With the new r(x), the difference d(x) = arctan(x) - r(x) is equi-oscillating at 8 points:
d(-1) |
- 0.00000641741 |
d(-0.911087) |
+ 0.00000641589 |
d(-0.679890) |
- 0.00000641194 |
d(-0.360013) |
+ 0.00000641687 |
d(0.360013) |
- 0.00000641687 |
d(0.679890) |
+ 0.00000641194 |
d(0.911087) |
- 0.00000641589 |
d(1) |
- 0.00000641741 |
In mathematics, we say that r(x) is the best rational approximation
to arctan(x) of degree (3,4) on [-1,1]. Theoretically, the digits
on the right hand side should have the same absolute value. Here,
they are equal up to the 9th digit."