2
$\begingroup$

In a nutshell, the Moler and Morrison algorithm is a fast method for calculating euclidean distances in a numerically stable way by using reflections instead of the pythagorean theorem.

In order to reflect a point $(x,y), 0\le y\le x,$ towards the $x$-axis in a norm-preserving way, Moler and Morrison use the line through $(0,0)$ and $(x,\frac{y}{2})$ as an approximate angle-bisector for the reflection of $(x,y)$; the convergence rate is cubic, as can be seen from the taylor series of $$\arctan\left(\frac{\tan(x)}{2}\right) = \frac{x}{2} + \frac{x^3}{8} + \frac{x^5}{32} + \frac{11 x^7}{1920} + \frac{25 x^9}{96768} - \frac{3443 x^{11}}{9676800} + O(x^{12})$$

Now, out of curiosity, I checked, what would happen, if the angle-bisector were approximated by the line through $(0,0)$ and $(x+\frac{y}{2}\frac{y}{x},\frac{y}{2})$, i.e. the line through the origin and the midpoint on the normal to $(x,y)$; well, somewhat surprising to me, the taylor series of $$\arctan\left(\frac{\frac{\sin(x)}{2}}{\cos(x)+\frac{\sin(x)}{2}\tan(x)}\right) = \frac{x}{2} - \frac{x^3}{8} - \frac{x^5}{32} - \frac{11 x^7}{1920} - \frac{25 x^9}{96768} + \frac{3443 x^{11}}{9676800} + O(x^{13})$$ and, investigating the obvious idea of taking as the $x$-coordinate of the point, through which the approximate angle-bisector should pass, the arithmetic mean of those two points, one gets $$\arctan\left(\frac{\frac{\sin(x)}{2}}{2\frac{\cos(x)}{2}+\frac{\sin(x)}{4}\tan(x)}\right) = \frac{x}{2} - \frac{x^5}{32} - \frac{5 x^7}{384} - \frac{x^9}{288} - \frac{289 x^{11}}{387072} + O(x^{13})$$ so the error of approximating the angle-bisector has reduced from $O(x^3)$ to $O(x^5)$.

Question:
would an algorithm that uses the better approximation of the angle-bisector really be an improvement over the Moler Morrison algorithm, when taking into account the additional operations that are necessary for obtaining the better approximation for angle-bisector?

I know, that an answer will depend on the specific circumstances; therefore I am not looking for a yes/no decision, but rather a list of answers that address the specifics of different situations, e.g. datatypes or coding schemes like IEEE 754 or BCD for real implementations or, for different theoretical models of computation.

$\endgroup$

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.