Python rounds float values by converting them to string and then back
github.com
github.com
From the comment in protobuf source (which does the same thing as Python), mentioned in the Twitter thread:
(...) An arguably better strategy would be to use the algorithm described in "How to Print Floating-Point Numbers Accurately" by Steele & White, e.g. as implemented by David M. Gay's dtoa(). It turns out, however, that the following implementation is about as fast as DMG's code. Furthermore, DMG's code locks mutexes, which means it will not scale well on multi-core machines. DMG's code is slightly more accurate (in that it will never use more digits than necessary), but this is probably irrelevant for most users.
Rob Pike and Ken Thompson also have an implementation of dtoa() in third_party/fmt/fltfmt.cc. Their implementation is similar to this one in that it makes guesses and then uses strtod() to check them. (...)
https://github.com/protocolbuffers/protobuf/blob/ed4321d1cb3...
The C/C++ standards do not require formatting to round correctly or even be portable. I recently had an issue where a developer used this method to round floats for display, and there were differences on PC and on Mac. It literally rounded something like 18.25 to 18.2 on one platform and 18.3 on the other. This led to all sorts of other bugs as some parts of the program used text to transmit data, which ended up in weird states.
The culprit was this terrible method. If you want anything approaching consistency or predictability, do not use formatting to round floating point numbers. Pick a numerically stable method, which will be much faster of done correctly.
Coincidentally, C/C++ do not require any of their formatting and parsing routines to round-trip floating point values correctly (except the newly added hex formatted floats which are a direct binary representation, and some newly added function allowing an obscure trick I do not recall at the moment... )
The linked-to method uses PyOS_snprintf(). Its documentation at https://docs.python.org/3/c-api/conversion.html says:
"""PyOS_snprintf() and PyOS_vsnprintf() wrap the Standard C library functions snprintf() and vsnprintf(). Their purpose is to guarantee consistent behavior in corner cases, which the Standard C functions do not."""
The python wrapper also does not specify in this case, so you should not use them for rounding, or you will have the same problem. No where on that page does it specify proper rounding will be cross platform.
Simply do it with floats. There are perfectly good, numerically stable, fast rounding methods, that avoid all this nonsense.
It’s obvious that PyOS_snprintf is not a standard library function
snprintf() wrappers. If the platform has vsnprintf, we use it, else we
emulate it in a half-hearted way. Even if the platform has it, we wrap
it because platforms differ in what vsnprintf does in case the buffer
is too small:
The "the best ways to go about it" comment links to the protobuf code, which also uses snprintf.
The (what I think is the) relevant C99 spec at http://www.open-std.org/jtc1/sc22/WG14/www/docs/n1256.pdf says:
> For e, E, f, F, g, and G conversions, if the number of significant decimal digits is at most DECIMAL_DIG, then the result should be correctly rounded. If the number of significant decimal digits is more than DECIMAL_DIG but the source value is exactly representable with DECIMAL_DIG digits, then the result should be an exact representation with trailing zeros. Otherwise, the source value is bounded by two adjacent decimal strings L<U, both having DECIMAL_DIG significant digits; the value of the resultant decimal string D should satisfy L ≤ D ≤ U, with the extra stipulation that the error should have a correct sign for the current rounding direction.
So either 1) "The C/C++ standards do not require formatting to round correctly or even be portable.", in which case Python and protobuf are doing it wrong and somehow this issue was never detected, or 2) The C/C++ standards do require correct rounding, but the case described by ChrisLomont didn't quite meet the spec requirements to get precision and rounding modes to match across platforms. Or 3), I don't know what I'm talking about.
Here's [1] where you can query the current floating-point environment in C: "Specifics about this type depend on the library implementation".
Here's [2] where you can set some rounding modes in C++: "Additional rounding modes may be supported by an implementation.". Note this does not have by default bankers rounding which is used to make many scientific calculations more stable (lowers errors and drift in accumulated calculations). Many platforms do this by default, but it's not in the standard.
You can chase down this rabbit hole. I (and several others) did during the issue on the last project, and got to where it was well-known in numerics circles that this is not a well-defined process in C/C++. If it were, printing and parsing should round-trip, and it does not before a recent C++ addition, and now it only is guaranteed in a special case.
Other good ways could trade-off edge case comprehensiveness for performance or whatever. That doesn’t make this way less good.
The issue here is you don't know what a library that formats a float does, and is the function is not specified clearly (as in C/C++), you have zero way of knowing what you will get.
Thus I said to do it yourself, using proper numerics.
For what it's worth, it looks like different standard libraries make different choices on whether float->string conversion cares about the current rounding mode.
Consider that a float's internal representation is in base 2, and you're trying to round in base 10. Even if you didn't use a string, I'd assume you'd have to create an array of ints that contain base 10 digits in order to do the rounding, unless there are some weird math tricks that can be employed that can avoid you having to process all base 2 digits. And an array of ints isn't all that computationally different from a string.
(Realistically, calling wordexp should just abort the program. Now I actually want to make a hacked up musl that aborts in all the various "libc functions no one should ever use" and see how far I get into a Ubuntu boot..)
It's somewhat messier than I remember, because it uses csh as the first choice and falls back to sh.
Perhaps from the perspective of an end user running things from a shell. Generally speaking though, shelling out from within a program is not ideal.
/* XXX this is _not_ designed to be fast */laughing so hard :')
And imagine the debug errors:
>perl error X
"But I'm just calling libc ?!?"
$ otool -L /usr/bin/perl
/usr/bin/perl:
/System/Library/Frameworks/CoreFoundation.framework/Versions/A/CoreFoundation (compatibility version 150.0.0, current version 1663.0.0)
/usr/lib/libSystem.B.dylib (compatibility version 1.0.0, current version 1281.0.0)I always implemented round to a specific digit based on the built-in roundss/roundsd functions which are native x86-64 assembler instructions (i.e. https://www.felixcloutier.com/x86/roundsd).
I do not understand why this would not be preferable to the string method.
float round( float x, int digits, int base) { float factor = pow( base, digits ); return roundss( x * factor ) / factor; }
I guess this has the effect of not working for numbers near the edge of it's range.
One could check this and fall back to the string method. Or alternatively use higher precision doubles internally:
float round( float x, int digits, int base ) { double factor = pow( base, digits ); return (float)( roundsd( x * factor ) / factor ); }
But then what do you do if you have a double rounded and want to maintain all precision? I think there is likely some way to do that by somehow unpacking the double into a manual mantissa and exponent each of which are doubles and doing this manually - or maybe using some type of float128 library (https://www.boost.org/doc/libs/1_63_0/libs/multiprecision/do...)...
But changing this implementation now could cause slight differences and if someone was rounding then hashing this type of changes could be horrible if not behind some type of opt-in.
I don’t know of truly fast algorithms for string to float, although I improved upon our CRT’s performance by 40%.
Ryu is more than 100x slower than something like
rval = floor(100*val+0.5)/100.0
(which is not quite right due to numerical issues, but close, and illustrates the idea).
Formatting, to get a rounded float, is terribly slow.
I also consider simultaneously deterministic, fast, portable rounding or binary formats to be important for decentralized delegation of computation, say a TrueBit style platform:
Death by 1000 papercuts....
High level languages are also annoying in that they don't provide great support for working with them as binary denominated rationals, e.g. there's no round_base2 in python, and there's no hex float representation in printf.
https://retrocomputing.stackexchange.com/questions/7810/why-...
Humans generally are afraid of new words (especially weird sounding ones) and often will assume that the subject is complex and might intimidate them.
But unknown words can have extremely simple meanings, or be synonyms of already known words.
By asserting: "It's actually pretty simple" You give them confidence that there's not reason to be afraid of the topic or of the words.
Whenever I've helped older people with technology they've never used before (a new tablet or similar), if I started off with any suggestion that it's less than simple, they'll almost certainly frame the problem scope in their mind as difficult, and give up, because they're already exhibiting some animosity toward learning a new thing.
If instead you phrase it as "this is really easy, let me show you how...", you short-circuit this process by framing their expectations differently, and that little bit of extra confidence ("this is easy") can help them through the learning process.
I've found you can't simply show them, either. It's almost better to say "It's easy" and then go through the process, because it's absolutely necessary to establish expectations first. They're already afraid of it (it's new), so doing something to get their guard down can go a long way toward helping them explore on their own. I tried this experiment with my mother, and some weeks later she'd have a problem and discover the solution herself specifically because she was convinced it was easy to do. This can (and will) backfire if you're not careful about how you do it, but I've had far more success using this tool than other techniques individually (e.g. writing down instructions).
This doesn't broadly apply to areas outside education and support (or even to all areas in education), but for simple things that people may express an irrational fear over, it works and it works well. A good teacher will use this technique successfully with their students, so if you're teaching someone, use it!
If someone has some anxiety about not understanding something, telling them it's actually pretty simple can just reinforce the framing they already have going in that maybe they're too dumb to get it.
I've found it's usually better to acknowledge that it's a little difficult or otherwise totally normal not to already know / have grasped the thing in question.
I think the intent is in the right place when saying "it's actually pretty simple" -- you want to provide optimism. The approach I like is along the lines of "this part is a little tricky, so let's break it down."
> It's actually pretty simple. We'll be looking at something called the "D-Wave P-500", which is a version of the P500 chip for quantum computers.
> It's basically a single bit computer, but with more than 500 qubits. Which means that our "real number" will have more numbers than the number of qubits that are available. That's really important.
> Quantum computers are theoretically able to do more things than just solve equations. For example, the way that a quantum computer uses energy from an electron to solve a classical math problem, or the way that it can break a complex calculation into smaller bits of information that each can solve on its own, is very different from how computers currently work.
> But I am not suggesting that a quantum computer can be used to solve more abstract problems. Because that would be crazy.
> But to give an example of what it could do, imagine doing a number crunching function that was 10× faster than a classical chip, and that had some really useful, and practical things that would be interesting to try.
Because Talk to Transformer is trained on real-world data, this supports the hypothesis that the phrase "It's actually pretty simple" is often followed by an unintelligible and highly technical explanation.
Computers can only natively store integers,
so they need some way of representing decimal numbers.
> Computers can only natively store integers, so they need some way of representing decimal numbers. This representation comes with some degree of inaccuracy. . . Why does this happen? It's actually pretty simple. When you have a base 10 system (like ours), it can only express fractions that use a prime factor of the base. . .
Perhaps a more formally correct way to put it is that integers (and natural numbers) are the only numbers that a computer can manage in a way that behaves reasonably similarly to the corresponding mathematic set. Specifically, as long as you stick to their numeric range, computer ints behave like a group, just like real integers do. But IEEE floats break the definition of a field in every which way, so they're really not a great match for the rationals.
That said, you could represent the rationals as a pair of integers, and that would be better-behaved, and some programming languages do do that. But I'm not aware of an ISA that supports it directly.
Perhaps it would be better phrased as "computers can only store whole numbers, so they need some way to represent other information, including integers, rational numbers and approximations to complex numbers."
In some cases, rounding is performed for the primary purpose of displaying a number as a string, in which case it can't be any less complicated than the string conversion function itself.
import math
x = 0.49999999999999994
print(x-0.5)
print(math.floor(x+0.5))
-5.551115123125783e-17
1
Yet another reminder that floating point calculations are approximations, and not exact.
julia> using Benchmarktools
julia> @btime floor(0.49999999999999997+0.5)
1.0
julia> @btime floor(0.49999999999999997+BigFloat(0.5))
0.0
julia> @code_native floor(0.49999999999999997+BigFloat(0.5))
.text
; ┌ @ floatfuncs.jl:152 within `floor'
subq $24, %rsp
movq %rsi, 16(%rsp)
movq (%rsi), %rax
; │┌ @ floatfuncs.jl:152 within `#floor#543'
movq %rax, (%rsp)
movabsq $jl_system_image_data, %rax
movq %rax, 8(%rsp)
movabsq $japi1_round_16005, %rax
movabsq $jl_system_image_data, %rdi
movq %rsp, %rsi
movl $2, %edx
callq *%rax
; │└
addq $24, %rsp
retq
nopw %cs:(%rax,%rax)
; └
julia> @code_native floor(0.49999999999999997+0.5)
.text
; ┌ @ floatfuncs.jl:152 within `floor'
; │┌ @ floatfuncs.jl:152 within `#floor#543'
; ││┌ @ floatfuncs.jl:152 within `round'
vroundsd $9, %xmm0, %xmm0, %xmm0
; │└└
retq
nopw (%rax,%rax)
; └floating-point calculations are absolutely exact. What's not is conversions between decimals and floating-point.
Are you sure about that? What about 1.00000000000001^2 (using eg 64 bit double)?
Or just explicitly check for 0.5 - 1ulp as a special corner case.
Is there a phrase for the ratio between the frequency of an apparent archetype of a bug/feature and the real-world occurrences of said bug/feature? If not then perhaps the "Fudderson-Hypeman ratio" in honor of its namesakes.
For example, I'm sure every C programmer on here has their favored way to quickly demo what bugs may come from C's null-delimited strings. But even though C programmers are quick to cite that deficiency, I'd bet there's a greater occurrence of C string bugs in the wild. Thus we get a relatively low Fudderson-Hypeman ratio.
On the other hand: "0.1 + 0.2 != 0.3"? I'm just thinking back through the mailing list and issue tracker for a realtime DSP environment that uses single-precision floats exclusively as the numeric data type. My first approximation is that there are significantly more didactic quotes of that example than reports of problems due to the class of bugs that archetype represents.
Does anyone have some real-world data to trump my rank speculation? (Keep in mind that simply replying with more didactic examples will raise the Fudderson-Hypeman ratio.)
I remember even writing a program that tested every possible floating point number (must have only been 32 bit). I think I used ctypes and interpreted every binary combination of 32 bits as a float, turned it into a string, then back and checked equality. A lot of them were NaN.
I seem to recall that ~0.5% of the IEEE 32 bit float space is NaN.
That means there's 25 bits that we can change while still being either a NaN or an inf. But two of those values are infs, so we need to remove them. Divide that by the entire range and we have (2^25 - 2) / 2^32 = 16777215/2147483648, or about 0.78124995%.
Example of implementing it the sane way: https://github.com/numpy/numpy/blob/75ea05fc0af60c685e6c071d...
Every step of this function is complex and expensive, especially printing a float as a decimal is very complex. And round is routinely used in a tight loop.
>>> round(56294995342131.5, 2)
56294995342131.5
>>> round(56294995342131.5, 3)
56294995342131.5
>>> np.round(56294995342131.5, 2)
56294995342131.5
>>> np.round(56294995342131.5, 3)
56294995342131.51You kind of got me thinking now. The decimal representation of a number is really a string representation (in the sense of a certain sequence of characters). Hence rounding to a certain decimal is essentially a string operation. You can of course do it by (say) dividing by 10^whatever or something else in some numerical fashion, but the more I think about it, the more natural it is to just think of the whole thing as a string.
Or you could flip it around and consider that the string manipulation can also be described numerically so whether you consider the operation as a string operation or a numerical operation is sort irrelevant. It's just a point of view.
It is incorrect to speak of "the" decimal representation of a number, as many numbers have non-unique decimal representations, the most famous example being 1.000...=0.999...
That's true and I'll concede that point, but it's not really relevant to what I said. That just means some numbers have different string representations that represent the same object. That doesn't really contradict anything in my post except the use of the definite article.
It would previously scale up, round (ceil/floor really) then scale down. That turned out to induce severe precision issues: https://bugs.python.org/issue1869
The only silly part of ieee754 2008 is the fact that they specified two representations (DPD, championed by IBM, and BID, championed by Intel) with no way to tell them apart.
CPython rounds float values by converting them to string and then back
Jython instead uses BigDecimal::doubleValue: https://github.com/jythontools/jython/blob/b9ff520f4f6523120...
But as another comment noted, BigDecimal::doubleValue can pull a similar trick: https://news.ycombinator.com/item?id=20818586
%timeit round(12335423552.33, -6)
%timeit int(12335423552.33 / 1_000_000.0) * 1_000_000.0
500 ns ± 3.88 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)
219 ns ± 1.1 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)