Thursday, April 03, 2014

Decimal Floating Point Arithmetic

If you're doing things like financial calculations, you have to be careful about using conventional binary floating point because it can't represent decimal fractions exactly.

One approach is to use "scaled" numbers, e.g. represent your dollar amount in cents or hundredths of cents so you are always working in integers. And it requires big integers, 32 bits is only about 9 decimal digits and the 52 bits of double floats is about 15. You really need 64 bit integers which are about 19 digits. (10 bits ~ 3 decimal digits) But that still doesn't give you the ability to deal with general purpose floating point.

So Suneido has always had a decimal floating point numeric type. (Internally, for performance, it also uses plain integers when possible.) Another advantage of a decimal type is that it is simple and quick to convert to and from string form.

Back when I first wrote Suneido (~ 15 years ago) there were no 64 bit integers in C++ compilers ("long" was 32 bits) and 32 bits wasn't sufficient precision. So I had to use multiple values to hold the coefficient. Since I had to use multiple integers anyway, to simplify overflow (by using 32 bit ints for intermediate results) I used four 16 bit ints, each one holding four decimal digits for an overall precision of 16 decimal digits. (To simplify "shifting" the exponent is in terms of the 16 bit ints, i.e. it jumps 4 decimals at a time. This "granularity" causes problems with precision. Depending on the exponent, in the worst case you get as few as 10 decimal digits of precision.)

Of course, having to use multiple integers and trying to get decent performance complicated the code, especially division. I won't claim it's the greatest code, but nevertheless it's worked reasonably well for a long time.

When I implemented jSuneido, I used Java's BigDecimal. Because of the different implementation there were a few minor differences, but they were mostly edge cases that didn't matter in practical usage. (Unfortunately I had made the external dump format for numbers mirror cSuneido's internal representation so it's a little awkward converting to and from BigDecimals.)

Recently, we've started to run into issues with using API's that deal with 64 bit integers, because we don't have enough precision to store them. In jSuneido it would be easy to bump up the BigDecimal precision to 20 digits. In theory I could do the same with cSuneido, but unfortunately, the code is fairly specific to the current precision. e.g. loops are unrolled. The thought of making this change is not pleasant :-(

The other problem is that some of the code assumes that you can convert to and from 64 bit integers losslessly. But 20 decimal digits won't always fit in a 64 bit integer.

Now that we have 64 bit integer types, the obvious answer seems to be to use a 64 bit integer for the coefficient. This will be faster and simpler than using multiple small integers, and probably faster than BigDecimal since it handles arbitrary precision. And if I used the same approach in both cSuneido and jSuneido this would ensure consistent results.

Since I'm in the middle of playing with Go, I figured I'd try writing a Go version first. It should be relatively easy to port to C++ and Java if I decide to.

It took me a couple of days to write it. One of the challenges is detecting overflow when calculating with 64 bit integers, since you don't have a larger type to use for intermediate calculations. Hacker's Delight provided a few tips for this. Another useful reference was General Decimal Arithmetic.

It's about 500 lines for add, subtract, multiply, divide, and conversion to and from strings. (That's about half the size of the cSuneido C++ code.) Since I'm new to Go, it may not be the most idiomatic code. And I have only done basic testing and refactoring. "float10" isn't the greatest name. Maybe "decimal" or even "dec"? (in keeping with Go's predilection for short names) I'm open to suggestions...

I chose to pass and return by value rather than by pointer. I'm not sure if this is the best choice for a 10 byte struct. Would it be faster to pass by pointer? Returning by pointer forces heap allocation for intermediate results which isn't ideal. Pass and return by value is a good fit for immutable values which are my preference.

Go makes it really easy to benchmark so I checked the speed of division (the slowest operation). Micro-benchmarks are always dubious, but it gave me a rough idea. It showed about 400 ns per divide (on my iMac). I don't have comparable benchmarks for cSuneido or jSuneido, but that seems pretty good. I'm pretty sure it's better than cSuneido. (Of course, it's nowhere near as fast as native binary floating point done in hardware. The same benchmark with float64 gives about 7 ns per divide, although this is so small that it's even less likely to be accurate.)

As far as evaluating Go, so far I like it. Of course, it's well suited to low level code like this. Sublime Text + GoSublime works well. (The only issues have been with learning Sublime since I haven't used it much.) I might have broken out the debugger a couple of times if I'd been working in Java or C++, but I get the impression the debugger story for Go isn't that great. I managed easily enough with old school prints :-) I plan to give Eclipse + GoClipse a try at some point since I'm already familiar with Eclipse.

The code is embedded below but it's probably easier to read (or download) on GitHub.


Anonymous said...

Based on just a quick browse of the source code, it looks pretty good. Two things popped out at me. One is that elements with common type can be elided, so

func Add(x, y Float10) Float10

rather than

func Add(x Float10, y Float10) Float10

Second, I don't like the Nlz function name (coming from someone who doesn't do decimal math). Maybe LeadZero, or LZero or something? Nlz sounds like some form of algorithms. LeadZero combined with the signature is pretty clear, though maybe too long for the simplicity. LZero at least hilights the zero rather than the N

Andrew McKinlay said...

Thanks for the review.

Yeah, I still forget about eliding common types, and reading it still takes me a minute to figure out what it means. My first response is to think there's a type missing. I'll get used to it.

"Nlz" (Number of Leading Zeroes) is a somewhat standard name (e.g. in Hacker's Delight) I agree it's a bit cryptic if you're not familiar with it.

Andrew McKinlay said...

I wonder why gofmt doesn't clean up the common types? Does that mean there isn't a consensus on whether that's "good" or not?

Could you clean this up with a fmt rewrite rule?