Wednesday, January 02, 2008

Slightly Less Double (Im)precision

Anonymous writes

You may want to try evaluating the polynomial using a different form. For example, given the polynomial: A*x^3 + B*x^2 + C*x + D one is often tempted to just enter it "as is". However it can also be expressed this way: ((A*x+B)*x+C)*x+D Generally speaking polynomials often work better using floating point this alternate way. Some programming languages/systems know this and convert the expression before evaluating it.

...and he is absolutely right. It makes the code of heck of a lot more efficient, too, because it eliminates the calls to Math.pow. That being said, it does not completely fix the problem. The polynomial line is a lot smoother, and the fitting algorithm yields a lower degree polynomial for the same mean error, but I still think the results are too fuzzy compared to the higher precision floating point. Here's a graph to show the difference:

Compared to the previous result:

Further improvement could probably be obtained by taking a close look at the QR Decomposition algorithm used to do the least-squares fitting.

In my opinion, the problem here is not so much the double-prevision floating points are bad. They are not. For many applications, especially with carefully crafted algorithms, they are great. They are certainly much higher performance than their higher-precision object-oriented kin. I'll warn you: Don't do high precision matrix operations with your MacBook in your lap - it gets really hot. The double-precision version finishes in a blink. It also takes a several orders of magnitude longer than using doubles. The problem is that, as an abstraction for real numbers, doubles are extremely leaky. Of course, this could be extended to any fixed-prevision, floating-point representation of numbers, depending the application.

Basically, I think in most applications doubles represent a premature optimization. Higher-precision numbers should be used by default, and then the precision reduced in order to improve performance, rather than doubles being used by default and then higher-precision numbers being considered if the programmer realizes that he has a problem due to lack of precision.

Of course, the problem is I don't know how precise is enough, because it depends entirely on the application. I'm tempted so say that, whenever possible, exact representations should be used. I've done a little research into it, and I'll do some more. There's tons of papers on the subject, but everything I've encountered so far seems to require intimate knowledge of floating points, the algorithm using the, and possibly the data being fed into the algorithm. That could help with library functions, such as matrix decomposition and solving, which could automatically scale up the precision of their internal operations in order to meet the expected resulting precision of the calling code, but that would still be leaky for "user-implemented" algorithms. What I really want is something that will work in the general case with reasonable performance characteristics, which can then be tuned for specific applications by reducing precision.

Sphere: Related Content

1 comment:

pstickne said...

Indeed. It's an amazing difference what simply carefully ordering floating point operations can make... take a completely upfront approach to solving a system of linear equation vs one that also uses pivots in the underlying matrix.