16 ms·

There are only four billion floats, so test them all (2014)

501 points by tape_measure 6y ago

Const-me 6y ago
In addition to many correctness issues, that code is relatively slow. SSE4.1 introduced roundps instruction that's both fast and correct. Unlike AVX2 or FMA3, SSE4.1 support is almost universal by now. In most of my current projects, I check support on startup and refuse to run if SSE4.1 ain't supported. Just don't forget to set the corresponding macros for DirectXMath, Eigen, etc.
ascar 6y ago
Previous submission 2017: https://news.ycombinator.com/item?id=15252426 https://news.ycombinator.com/item?id=15252426 Can't believe this was already 3 years ago. Interesting read! Previous submission 2014: https://news.ycombinator.com/item?id=7135261 https://news.ycombinator.com/item?id=7135261
- tape_measure 6y ago
  Yes, I originally came across this on HN, likely in 2017 based on my surrounding bookmarks. It has stuck with me ever since even though I don't do much programming anymore. At the time, I was a TA in a MATLAB course, so I tried to emphasize that the code runs on a machine and is fundamentally discrete. (.1 + .2 == .3) returning false is an easy way into this concept. 32 bits (and now 64 as @dragontamer calculated) is a huge input space. You could verify 4 int_8 or 2 char[2] functions.
MeteorMarc 6y ago
For python you cannot user the stdlib for this, but it would work for numpy float32.
fallingmeat 6y ago
maybe I dont understand, but is this article only trying to suggest exhaustive testing for single precision, single input functions? obviously this approach is not practical as soon as there are multiple inputs or wider precision.
- pm215 6y ago
  Yes, the article says "Exhaustive testing works brilliantly for functions that take a single float as input. I used this to great effect when rewriting all of the CRT math functions for a game console some years ago. On the other hand, if you have a function that takes multiple floats or a double as input then the search space is too big. In that case a mixture of test cases for suspected problem areas and random testing should work." I think the point is that exhaustive testing is feasible in more situations than you might at first think, and than floating-point arithmetic is a particularly fertile source of edge conditions and corner cases, so it's a good place to go for exhaustive testing where you can.
  - fallingmeat 6y ago
    fertile indeed!!!! haha I should learn to read
- nullc 6y ago
  In many cases it is practical but people wrongly assume it isn't. In cases where its not practical exhausting a subset can still be very powerful. E.g. for a function that takes two floats (x,y) test every value of x input with y={0,-0,1,nan,inf,-1,x,-x,2x,x/2,FLT_EPSILON,FLT_MIN,FLT_MAX,random,random,random...}, and then repeat with x,y swapped. A little knowledge of the domain can help pick better special values, but there is a trade-off in introducing the same blind-spots that permitted bugs in the first place. Still: testing just one argument on all values because that's all you can afford is still a lot better than only testing a couple special values for both arguments. Exhaustive testing, like randomized testing catches errors you weren't even thinking about-- but unlike random testing exhaustion can reliably catch bugs which are extremely rare. "Partial exhaustion" can have some of the same benefits.
  - fallingmeat 6y ago
    mmm, disagree based on my experience. oh also I think people on this thread may be misusing the concept of float epsilon. keep in mind that is not the smallest value of a float.
sega_sai 6y ago
It is all very nice, but if your function has two floating point parameters you have 1.6e19 tests to do...
- andi999 6y ago
  So it only needs 4 Billion times 90s, so is this like never in a million years? (actually it should be only around 12000 year, so you could use a cluster of 48000 Computer and be done in 3 month)
  - Someone 6y ago
    Distribute it over every smartphone on the planet (about 4 billion), and the computations will be done in minutes (of course, that’s dependent on the function being tested). There will be orchestration overhead and transmission latency, but I guess one could get all results in one place in a day, if one had the ability to do the “run it on every smartphone” thing.
- thechao 6y ago
  I have worked at two companies where we brute-force checked binary floating point operators. It’s annoying, error prone, boring, and soooooo satisfying. The Holy Grail^TM is brute force F32-FMA. That one really requires a lot of thought, and can’t be reasonably brute-forced on anyone’s budget.
- ChrisLomont 6y ago
  Then you can prove your algorithm correct using proof libraries like Z3.
- formerly_proven 6y ago
  Totally doable. That's just a few hours using 100k cores.
- TheRealPomax 6y ago
  Worth it if your function is acutely fundamental to the operation of your software because your software runs everyone else's software. You buy some time on a compute cluster with 1024 cores and 4TB of ram, and off you go, after making sure you understand your own code and reduce that number back down to "no it isn't, because two arguments in a given space for fundamental maths operations, taken together, almost certainly don't require combinatorial testing" =)
xchaotic 6y ago
Articles like this make me wonder how everything even works if basic math functions can be so wrong and how a simple concept like floor can have different definitions/ implementations
- marketingPro 6y ago
  Someone correct me if I'm wrong because it's been years since my EE class. Basic math ends at division. Heck if you learn how math works on computers, it feels like a number trick. Those are reliable. The issue comes when you need to track never ending decimal points. Even if every transistor was dedicated to your problem, you'd still need to round. So even if you have 32 bits of decimals, there will be rounding sometimes. The human has to decide what to do. Instead of converting the binary to decimal and calling it done, it's more important to be consistent and predictable. The issues described in the article arent 5/2= 2 billion. The article talks about rounding giving you always even numbers, which you can imagine leads to issues if you are using both float and using logic that involves even and odd numbers. Edit- can someone correct me or add details, I'm very interested but I'm not taking ENG 240 again... :p
  - esperent 6y ago
    Exactly this. The floor function works on binary numbers but has to round them to a decimal place. It's not simple, it just appears to be.
    - bonoboTP 6y ago
      Floor has nothing to do with decimal.
      - erichocean 6y ago
        The period in "3.14" is called a "decimal." The digits after the decimal in "3.14" are enumerated as "decimal places." What each decimal place "means" depends on the numeric base. So, "decimal places" or "rounding to a decimal place" has nothing to do with base 10 (which is what I think you are implying): decimal and "decimal places" exist in all bases, including base 2, and the meaning of each "decimal place" depends on the numeric base (base 10, base 2, etc.).
throwaway7281 6y ago
I was once in the position to test 3M possible inputs - that took only a few minutes. It's a good strategy to keep in mind.
deleted 6y ago
[deleted]
zaroth 6y ago
> For very small numbers (less than about FLT_EPSILON * 0.25) adding 0.5 gives 0.5 exactly, and this then rounds to zero. Since about 40% of the positive floating-point numbers are smaller than FLT_EPSILON * 0.25 this results in a lot of errors – over 850 million of them! I thought Float Epsilon was the smallest possible number which could be represented by a float? This statement doesn’t seem to make sense...
- Sebb767 6y ago
  I thought that too, but it's actually the difference between 1.0f and the next larger number; so it only represents the smallest difference for the non-dotzero range. See here: http://www.cplusplus.com/reference/cfloat/ http://www.cplusplus.com/reference/cfloat/
- MatmaRex 6y ago
  There are different "epsilon" constants in different languages that mean different things. In C and C++, FLT_EPSILON is "Difference between 1 and the least value greater than 1 that is representable." (http://www.cplusplus.com/reference/cfloat/ http://www.cplusplus.com/reference/cfloat/). You might be thinking of the C# constant.
- neeeeees 6y ago
  FLT_ESPILON is the bound on rounding error for certain operations. Plenty of numbers smaller than it can be represented in most floating point implementations
- jcranmer 6y ago
  When I took numerical analysis, the concept was described as "machine epsilon", which is the smallest number that can be added to 1.0 that is not 1.0. The utility of this metric is that it tells you the best possible relative error, as you can never guarantee a result will be more correct to the true result (in the domain of mathematical real numbers as opposed to machine floating-point numbers) than half of machine epsilon. This concept is often extended into discussions of "units in the last place" (ULPs), which is the magnitude of the value of last bit of the mantissa. So ulp(1.0) = ulp(1.5) = machine epsilon, ulp(2.0) = ulp(3.9) = 2 * machine epsilon, etc. A good floating-point library will document its accuracy in terms of ulps, so that we might say that the maximum error of the exp function is 1 ULP and that of pow might be 6 ULPs.
draw_down 6y ago
If your code works for 1.000000000001 but fails for 1.000000000002, go back.
tialaramex 6y ago
Many years ago now we built something that used a 60-bit truncated hash of URIs. So that's far too small to be comfortable that you won't run into collisions, but it's big enough that for modest sizes (hundreds of millions up to billions) a collision is unlikely. So you can have a slow path for the collision case, so long as that slow path is correct. For unit testing we needed at least one collision so we could see that everything behaves as designed before some real user hits one. We had test data sets with over a billion URIs in them, but none that (as far as we could tell) had a collision. So I just wrote code to mint nonsense URIs for our system based on a small integer seed (like https://foo.example/123456789 https://foo.example/123456789), then span up a system with a bunch of RAM (maybe 8-12GB) to iterate through each seed, storing the hash and seed it was generated from, until it hit a collision after a few gigabytes of RAM, then it spat out the two URIs from the colliding seeds. Bingo, now we had data for unit testing to show that stuff works even in the unlikely case that two distinct URIs have the same hash. Probably a day's work to build that, and then I'd guess several hours every month saved debugging weird collision bugs because our unit tests would now tell us if we'd screwed up and collisions were slipping past.
- seanmcdirmid 6y ago
  Couldn’t you also test the recovery logic with a hasher that yielded collision more often? It’s like take your system, make it worse (use very dumb hashing), and see how it degrades and/or breaks?
  - xxxtentachyon 6y ago
    Or, even simpler, just stub out your hash function and to something that always generates the same value and test with that.
  - rocqua 6y ago
    Seems like that would make it less of an integration test
    - jmholla 6y ago
elvis70 6y ago
"When in doubt, use brute force" -- Ken Thompson
IncRnd 6y ago
> I’ve written in great detail about how to compare floating-point values using an epsilon, but there are times when it is just not appropriate. Sometimes there really is an answer that is correct, and in those cases anything less than perfection is just sloppy. Certain applications will choose fixed point to avoid what appears to be the imprecisions of floating point.
dragontamer 6y ago
Things have gotten even faster since then. * Target: 64-bit space (2^64). * GPUs are more flexible and faster. A Vega64 (a few years old now) has 4096 shaders ("CudaThreads") at 1.6 GHz. It can search the entire 32-bit floating point (or integer) space in less than a second. 1.6 GHz x 4096 shaders == 2^42.5 bitspace per second. * People exhaustively searching a space are relatively patient: and willing to wait months for the results. 60 seconds x 60 minutes x 24 hours x 30 days == 2^21.3 1.6 GHz x 4096 shaders x 60 seconds x 60 minutes x 24 hours x 30 days is roughly equal to 2^63.8. We've done it: a month of GPU effort is all you need to exhaustively search the 64-bit space ------- We've entered a world (~3 years ago), where it is reasonable to brute force the 64-bit space, so long as you're willing to wait about a month for the results. If you're not willing to wait a month, then use 2 GPUs, or 4 GPUs, or a cluster of GPUs like the DGX A100 (or multiple clusters).
- Keyframe 6y ago
  Sounds like a small and relatively cheap farm could do it in hours or minutes even.
  - dragontamer 6y ago
    Hmmm... assuming 4-GPUs per node, and 1200W per node, you'll need 216kW to generate the result in a day. But that also assumes that the operation you're testing is 1-clock tick, which is only true for the most trivial of calculations. Still, 216kWs per day per "test-case clock ticks" shows us that the 2^64 exhaustive bitspace is within a reasonable organization's grasp.
    - dan-robertson 6y ago
      I think your units are wrong. You want 216kW days, not 216kW per day. 216kW days is 18.6Gj.
      - dragontamer 6y ago
        Hmm, "216kW days per test-case clock ticks" does sound more correct from a unit perspective.
known 6y ago
Very difficult due to rounding errors;
vortico 6y ago
Does anyone have code to reduce the search space by a few orders of magnitude to test functions with a float64 argument or two or three float32 arguments in a few seconds? Perhaps by skipping "boring" floats that aren't near endpoints that are unlikely to produce errors that more interesting floats do?
- TheRealPomax 6y ago
  all floats are equally boring.
  - vortico 6y ago
    No, -0.f is more interesting than 4.1045857301e23f. It's more interesting because it's more likely to fail a test where others don't.
    - TheRealPomax 6y ago
      That's a bold assumption: let's make sure it's correct by testing against all floats, instead of only the ones you consider interesting, in the process making any float as boring as any other float again.
- daveFNbuck 6y ago
  That's essentially what fuzzy testing libraries do. You tell the library what types of arguments you want, write your test, and it puts in some known corner cases and a few random values that change each time you run the test (but are saved as new corner cases when the test fails in the better libraries).
- rigtorp 6y ago
  Yes, I've done this using libFuzzer. Supply custom mutation and crossover functions that generates interesting floating point numbers. Then let the coverage guided fuzzer exercise your code. https://rigtorp.se/fuzzing-floating-point-code/ https://rigtorp.se/fuzzing-floating-point-code/
- ramshorns 6y ago
jlebar 6y ago
We did this in XLA. It runs as part of CI, which I always found pretty amazing. https://github.com/tensorflow/tensorflow/blob/bd2608178bca1be95007328fb9e01408a94ba082/tensorflow/compiler/xla/tests/exhaustive_op_test_utils.cc https://github.com/tensorflow/tensorflow/blob/bd2608178bca1b... Found a lot of bugs this way.
NelsonMinar 6y ago
"Several of the ceil functions were implemented by adding 0.5 to the input value and rounding to nearest." What sort of clown would make that mistake and yet think they are capable of writing a floating point library?
- Exuma 6y ago
  Can you give an example of when this wouldn't work? I'm definitely not qualified to write a float library (lol) but I also can't think of any cases where that wouldn't work?
  - brianush1 6y ago
    ceil(5) = round(5.5) = 6 ceil(5) should be 5
    - brucedawson 6y ago
      To expand on this, if you have an integer and you add 0.5 then a round() function has to make a choice, because it is halfway between integers. The IEEE standard says that rounding in general (not just for the round function) should round-to-nearest-even. That means that rounding should go to the closest result, and if there is a tie then it should go to the even number. So, if ceil(x) is implemented as round(x+0.5) then: ceil(4) = round(4.5) = 4 ceil(5) = round(5.5) = 6 The first of those is correct. The second is incorrect.
  - zelphirkalt 6y ago
    This is probably about number precision and loosing digits, when you add something with digits only in the beginning like 0.5 to something small like 0.1000000000000000000000000000013124142579275927597245205 or whatever a barely representable floating point number is. Then the result becomes inprecise and rounding to nearest will result in the wrong value. Such effects and similar should be well known to anyone taking on the responsibility to develop a floating point library.
  - NelsonMinar 6y ago
zelphirkalt 6y ago
Test them all and waste a huge amount of energy and computational resources! I hope people test like that for really critical stuff, but not for your everyday library for your personal startup or e-business. Not even the big players need this.
zzzeek 6y ago
Crap now I have 67 million failures to fix.
ViralBShah 6y ago
We've certainly tested all Float32 inputs to some functions in Julia's math packages. I'm unable to find the relevant issues this very moment, otherwise would certainly post them.
- raphlinus 6y ago
  I suspect https://github.com/JuliaLang/julia/pull/29700 https://github.com/JuliaLang/julia/pull/29700 might be one of them. As the original post states, rounding is (a) amenable to hardware acceleration, (b) easy to get wrong or slow, and (c) possible to test exhaustively.
ChuckMcM 6y ago
I love this. This is another in the "Now that we have computers with gazillions of bits of RAM and many cores and gazillions of cycles, what conventional wisdom is no longer wise?" series :-). I suggested to a client that they use Sqlite3 for their company database rather than a disk based database. And they looked at me shocked, as if I had spoken heresy. I pointed out how many records they could store on a 64GB machine and that is a small data center class machine and it sort of boggled their mind. Add read-only copies, shared from a NAS device gave them perfect scalability in terms of queries per second (their sales folks had their own database that they were updating). It was fun to watch their expression go from horror to wonder to excitement when they finally implemented it.
- michaelmcmillan 6y ago
  Yup! I'm currently running a b2c solution for a client that's grossing ~$3 million/year using Python + sqlite3 behind nginx. I work for another shop that scolds me for even suggesting using something else than k8s. How can a persistent file system possibly scale? Lol. Ironically they're grossing -$3 million/year.
  - tonetheman 6y ago
    I am curious about your case is the traffic mostly read heavy? Ha wish you could share more. Interesting stuff.
- dragontamer 6y ago
  https://yourdatafitsinram.net/ https://yourdatafitsinram.net/ The IBM 53-bit quantum compute result from a ton of hard drives + compute power was also intriguing (https://www.ibm.com/blogs/research/2019/10/on-quantum-supremacy/ https://www.ibm.com/blogs/research/2019/10/on-quantum-suprem...). 128PiB of disk space needed for a 54-qubit computation on their system (hard drives). Fun to see where the supercomputers are still needed, but smaller qubits can be emulated on normal computers. 40 qubit calculations probably can be emulated on a 16TB team of hard drives in RAID0.
- kortilla 6y ago
jansan 6y ago
Javascript guy here, can you help me set up that Beowulf cluster so I can start testing my 2^64 cases?
JoeAltmaier 6y ago
My god this issue has been around for too long. In the '80s testing folk learned to shotgun the argument space to test functions, and they've not gone an inch forward since! Or so it seems. A version of the Pentium was released that couldn't divide by 10. And even back then, testing all the mantissas on a simulator would have taken only hours. Instead, they masked, marketed and released a tragically flawed processor. Just test the input space exhaustively, folks! It's a computer; it feels no pain. Work it hard until it fails or works perfectly.
- deleted 6y ago
  [deleted]
- drivebycomment 6y ago
  On most CPUs, there's plenty of 2 64-bit operands, and with SIMD instructions, plenty of CPUs have more than 1K bits total as input to the instruction. Good luck exhaustively testing those.
wtn 6y ago
In the early 1980s, truncated floating point calculations caused the Vancouver Stock Exchange index value to drift down for about 23 months, until it was shown as less than half the correct value: https://en.wikipedia.org/wiki/Vancouver_Stock_Exchange#Rounding_errors_on_its_Index_price https://en.wikipedia.org/wiki/Vancouver_Stock_Exchange#Round...
peter303 6y ago
I recall hackers building standard password hashes of all eight character ascii combinations. Then create a reverse lookup table. But thats 400 trillion entries for stand 68 character ascii password subset. 400 trillion hashes is not a big computation in the blockchain era.
- SirYandi 6y ago
  More info: https://en.m.wikipedia.org/wiki/Rainbow_table https://en.m.wikipedia.org/wiki/Rainbow_table
hinkley 6y ago
I decided to build a scatter plot to demonstrate what I was 95% sure was a data loss bug (preproduction, thankfully) in the initializer for a CSPRNG. I figured plotting a few bytes if output (256x256) would do nicely. I plotted all possible two byte inputs just to make sure it was working, and that ran instantly, so I tried 3 bytes and that was still pretty fast, so I did 4 while I was at lunch or in a meeting. Sure enough, big black splotches from dropping some of the seed data on the floor. Please to be fixing. Old farts know a lot of things have been tried before and how those worked out. But we also think doing a complex operation 4 billion times is going to be irretrievably slow. Always test your assumptions.
uyt 6y ago
I was reading about an implementation of Barrett reduction recently: https://en.wikipedia.org/wiki/Barrett_reduction https://en.wikipedia.org/wiki/Barrett_reduction Specifically one that uses 80-bit long double for calculating `(__int128_t)a * b % M` from this github: https://github.com/kth-competitive-programming/kactl/blob/master/content/number-theory/ModMulLL.h https://github.com/kth-competitive-programming/kactl/blob/ma... But how do you test the correctness of such code? You can't just stress test all possible product of 64 bit numbers! In the end you still gotta go back to good old fashion proofs: https://github.com/kth-competitive-programming/kactl/blob/master/doc/modmul-proof.pdf https://github.com/kth-competitive-programming/kactl/blob/ma...
tagh 6y ago
Makes sense when you're refactoring an existing function with one numeric argument, and which evaluates rapidly.
CGamesPlay 6y ago
1. Given the thoroughness of the write-up, I'm going to assume that the author did not find that any of the discrepancies in the implementations were actually bugs in the reference implementation. He goes as far as to say that he reimplemented the reference implantation so this seems like a safe assumption. Still, it would be nice to confirm. Given that other implementations seem so shoddy, it wouldn't be surprising to find a bug in the reference. 2. "Conventional wisdom Nazis" Maybe the author was just bored of writing at this point, but some justification for why bitwise comparing floats was a correct decision would have been nice. Is it because there's no math involving unrepresentable intermediates? Because there's no integration happening, so no way for errors to accumulate?
- brucedawson 6y ago
  "Sometimes there really is an answer that is correct, and in those cases anything less than perfection is just sloppy." Anytime there is a correct answer, and it can be determined exactly (+, -, *, /, sqrt) you should expect nothing less than perfection. ceil, round, and floor fall into this category fairly self-evidently, I think.
alkonaut 6y ago
How do you do this for a function with more than one or two inputs? If the domain is 4B to the power five, you are out of luck. Very few of us need to implement "floor" or any other function that takes 1 float input.
jart 6y ago
"Seymour Cray didn't care that 81.0/3.0 did not give exactly 27.0 on the CDC 6000 class machines; and he was universally respected for making the fastest machines around." ──Linus Torvalds
lentil_soup 6y ago
"Round to nearest even is the default IEEE rounding mode. This means that 5.5 rounds to 6, and 6.5 also rounds to 6" Does anyone know the reason behind this? why the nearest even and not a hard cut like "always round up with >= 0.5"?

501 points by tape_measure 6y ago

In addition to many correctness issues, that code is relatively slow. SSE4.1 introduced roundps instruction that's both fast and correct. Unlike AVX2 or FMA3, SSE4.1 support is almost universal by now. In most of my current projects, I check support on startup and refuse to run if SSE4.1 ain't supported. Just don't forget to set the corresponding macros for DirectXMath, Eigen, etc.

ascar 6y ago

Previous submission 2017: https://news.ycombinator.com/item?id=15252426 https://news.ycombinator.com/item?id=15252426 Can't believe this was already 3 years ago. Interesting read! Previous submission 2014: https://news.ycombinator.com/item?id=7135261 https://news.ycombinator.com/item?id=7135261

tape_measure 6y ago
Yes, I originally came across this on HN, likely in 2017 based on my surrounding bookmarks. It has stuck with me ever since even though I don't do much programming anymore. At the time, I was a TA in a MATLAB course, so I tried to emphasize that the code runs on a machine and is fundamentally discrete. (.1 + .2 == .3) returning false is an easy way into this concept. 32 bits (and now 64 as @dragontamer calculated) is a huge input space. You could verify 4 int_8 or 2 char[2] functions.

MeteorMarc 6y ago

For python you cannot user the stdlib for this, but it would work for numpy float32.

fallingmeat 6y ago

maybe I dont understand, but is this article only trying to suggest exhaustive testing for single precision, single input functions? obviously this approach is not practical as soon as there are multiple inputs or wider precision.

pm215 6y ago
Yes, the article says "Exhaustive testing works brilliantly for functions that take a single float as input. I used this to great effect when rewriting all of the CRT math functions for a game console some years ago. On the other hand, if you have a function that takes multiple floats or a double as input then the search space is too big. In that case a mixture of test cases for suspected problem areas and random testing should work." I think the point is that exhaustive testing is feasible in more situations than you might at first think, and than floating-point arithmetic is a particularly fertile source of edge conditions and corner cases, so it's a good place to go for exhaustive testing where you can.
- fallingmeat 6y ago
  fertile indeed!!!! haha I should learn to read
nullc 6y ago
In many cases it is practical but people wrongly assume it isn't. In cases where its not practical exhausting a subset can still be very powerful. E.g. for a function that takes two floats (x,y) test every value of x input with y={0,-0,1,nan,inf,-1,x,-x,2x,x/2,FLT_EPSILON,FLT_MIN,FLT_MAX,random,random,random...}, and then repeat with x,y swapped. A little knowledge of the domain can help pick better special values, but there is a trade-off in introducing the same blind-spots that permitted bugs in the first place. Still: testing just one argument on all values because that's all you can afford is still a lot better than only testing a couple special values for both arguments. Exhaustive testing, like randomized testing catches errors you weren't even thinking about-- but unlike random testing exhaustion can reliably catch bugs which are extremely rare. "Partial exhaustion" can have some of the same benefits.
- fallingmeat 6y ago
  mmm, disagree based on my experience. oh also I think people on this thread may be misusing the concept of float epsilon. keep in mind that is not the smallest value of a float.

sega_sai 6y ago

It is all very nice, but if your function has two floating point parameters you have 1.6e19 tests to do...

andi999 6y ago
So it only needs 4 Billion times 90s, so is this like never in a million years? (actually it should be only around 12000 year, so you could use a cluster of 48000 Computer and be done in 3 month)
- Someone 6y ago
  Distribute it over every smartphone on the planet (about 4 billion), and the computations will be done in minutes (of course, that’s dependent on the function being tested). There will be orchestration overhead and transmission latency, but I guess one could get all results in one place in a day, if one had the ability to do the “run it on every smartphone” thing.
thechao 6y ago
I have worked at two companies where we brute-force checked binary floating point operators. It’s annoying, error prone, boring, and soooooo satisfying. The Holy Grail^TM is brute force F32-FMA. That one really requires a lot of thought, and can’t be reasonably brute-forced on anyone’s budget.
ChrisLomont 6y ago
Then you can prove your algorithm correct using proof libraries like Z3.
formerly_proven 6y ago
Totally doable. That's just a few hours using 100k cores.
TheRealPomax 6y ago
Worth it if your function is acutely fundamental to the operation of your software because your software runs everyone else's software. You buy some time on a compute cluster with 1024 cores and 4TB of ram, and off you go, after making sure you understand your own code and reduce that number back down to "no it isn't, because two arguments in a given space for fundamental maths operations, taken together, almost certainly don't require combinatorial testing" =)

xchaotic 6y ago

Articles like this make me wonder how everything even works if basic math functions can be so wrong and how a simple concept like floor can have different definitions/ implementations

marketingPro 6y ago
Someone correct me if I'm wrong because it's been years since my EE class. Basic math ends at division. Heck if you learn how math works on computers, it feels like a number trick. Those are reliable. The issue comes when you need to track never ending decimal points. Even if every transistor was dedicated to your problem, you'd still need to round. So even if you have 32 bits of decimals, there will be rounding sometimes. The human has to decide what to do. Instead of converting the binary to decimal and calling it done, it's more important to be consistent and predictable. The issues described in the article arent 5/2= 2 billion. The article talks about rounding giving you always even numbers, which you can imagine leads to issues if you are using both float and using logic that involves even and odd numbers. Edit- can someone correct me or add details, I'm very interested but I'm not taking ENG 240 again... :p
- esperent 6y ago
  Exactly this. The floor function works on binary numbers but has to round them to a decimal place. It's not simple, it just appears to be.
  - bonoboTP 6y ago
    Floor has nothing to do with decimal.
    - erichocean 6y ago
      The period in "3.14" is called a "decimal." The digits after the decimal in "3.14" are enumerated as "decimal places." What each decimal place "means" depends on the numeric base. So, "decimal places" or "rounding to a decimal place" has nothing to do with base 10 (which is what I think you are implying): decimal and "decimal places" exist in all bases, including base 2, and the meaning of each "decimal place" depends on the numeric base (base 10, base 2, etc.).

throwaway7281 6y ago

I was once in the position to test 3M possible inputs - that took only a few minutes. It's a good strategy to keep in mind.

deleted 6y ago

[deleted]

zaroth 6y ago

> For very small numbers (less than about FLT_EPSILON * 0.25) adding 0.5 gives 0.5 exactly, and this then rounds to zero. Since about 40% of the positive floating-point numbers are smaller than FLT_EPSILON * 0.25 this results in a lot of errors – over 850 million of them! I thought Float Epsilon was the smallest possible number which could be represented by a float? This statement doesn’t seem to make sense...

Sebb767 6y ago
I thought that too, but it's actually the difference between 1.0f and the next larger number; so it only represents the smallest difference for the non-dotzero range. See here: http://www.cplusplus.com/reference/cfloat/ http://www.cplusplus.com/reference/cfloat/
MatmaRex 6y ago
There are different "epsilon" constants in different languages that mean different things. In C and C++, FLT_EPSILON is "Difference between 1 and the least value greater than 1 that is representable." (http://www.cplusplus.com/reference/cfloat/ http://www.cplusplus.com/reference/cfloat/). You might be thinking of the C# constant.
neeeeees 6y ago
FLT_ESPILON is the bound on rounding error for certain operations. Plenty of numbers smaller than it can be represented in most floating point implementations
jcranmer 6y ago
When I took numerical analysis, the concept was described as "machine epsilon", which is the smallest number that can be added to 1.0 that is not 1.0. The utility of this metric is that it tells you the best possible relative error, as you can never guarantee a result will be more correct to the true result (in the domain of mathematical real numbers as opposed to machine floating-point numbers) than half of machine epsilon. This concept is often extended into discussions of "units in the last place" (ULPs), which is the magnitude of the value of last bit of the mantissa. So ulp(1.0) = ulp(1.5) = machine epsilon, ulp(2.0) = ulp(3.9) = 2 * machine epsilon, etc. A good floating-point library will document its accuracy in terms of ulps, so that we might say that the maximum error of the exp function is 1 ULP and that of pow might be 6 ULPs.

draw_down 6y ago

If your code works for 1.000000000001 but fails for 1.000000000002, go back.

tialaramex 6y ago

Many years ago now we built something that used a 60-bit truncated hash of URIs. So that's far too small to be comfortable that you won't run into collisions, but it's big enough that for modest sizes (hundreds of millions up to billions) a collision is unlikely. So you can have a slow path for the collision case, so long as that slow path is correct. For unit testing we needed at least one collision so we could see that everything behaves as designed before some real user hits one. We had test data sets with over a billion URIs in them, but none that (as far as we could tell) had a collision. So I just wrote code to mint nonsense URIs for our system based on a small integer seed (like https://foo.example/123456789 https://foo.example/123456789), then span up a system with a bunch of RAM (maybe 8-12GB) to iterate through each seed, storing the hash and seed it was generated from, until it hit a collision after a few gigabytes of RAM, then it spat out the two URIs from the colliding seeds. Bingo, now we had data for unit testing to show that stuff works even in the unlikely case that two distinct URIs have the same hash. Probably a day's work to build that, and then I'd guess several hours every month saved debugging weird collision bugs because our unit tests would now tell us if we'd screwed up and collisions were slipping past.

seanmcdirmid 6y ago
Couldn’t you also test the recovery logic with a hasher that yielded collision more often? It’s like take your system, make it worse (use very dumb hashing), and see how it degrades and/or breaks?
- xxxtentachyon 6y ago
  Or, even simpler, just stub out your hash function and to something that always generates the same value and test with that.
- rocqua 6y ago
  Seems like that would make it less of an integration test
  - jmholla 6y ago

elvis70 6y ago

"When in doubt, use brute force" -- Ken Thompson

IncRnd 6y ago

> I’ve written in great detail about how to compare floating-point values using an epsilon, but there are times when it is just not appropriate. Sometimes there really is an answer that is correct, and in those cases anything less than perfection is just sloppy. Certain applications will choose fixed point to avoid what appears to be the imprecisions of floating point.

dragontamer 6y ago

Things have gotten even faster since then. * Target: 64-bit space (2^64). * GPUs are more flexible and faster. A Vega64 (a few years old now) has 4096 shaders ("CudaThreads") at 1.6 GHz. It can search the entire 32-bit floating point (or integer) space in less than a second. 1.6 GHz x 4096 shaders == 2^42.5 bitspace per second. * People exhaustively searching a space are relatively patient: and willing to wait months for the results. 60 seconds x 60 minutes x 24 hours x 30 days == 2^21.3 1.6 GHz x 4096 shaders x 60 seconds x 60 minutes x 24 hours x 30 days is roughly equal to 2^63.8. We've done it: a month of GPU effort is all you need to exhaustively search the 64-bit space ------- We've entered a world (~3 years ago), where it is reasonable to brute force the 64-bit space, so long as you're willing to wait about a month for the results. If you're not willing to wait a month, then use 2 GPUs, or 4 GPUs, or a cluster of GPUs like the DGX A100 (or multiple clusters).

Keyframe 6y ago
Sounds like a small and relatively cheap farm could do it in hours or minutes even.
- dragontamer 6y ago
  Hmmm... assuming 4-GPUs per node, and 1200W per node, you'll need 216kW to generate the result in a day. But that also assumes that the operation you're testing is 1-clock tick, which is only true for the most trivial of calculations. Still, 216kWs per day per "test-case clock ticks" shows us that the 2^64 exhaustive bitspace is within a reasonable organization's grasp.
  - dan-robertson 6y ago
    I think your units are wrong. You want 216kW days, not 216kW per day. 216kW days is 18.6Gj.
    - dragontamer 6y ago
      Hmm, "216kW days per test-case clock ticks" does sound more correct from a unit perspective.

known 6y ago

Very difficult due to rounding errors;

vortico 6y ago

Does anyone have code to reduce the search space by a few orders of magnitude to test functions with a float64 argument or two or three float32 arguments in a few seconds? Perhaps by skipping "boring" floats that aren't near endpoints that are unlikely to produce errors that more interesting floats do?

TheRealPomax 6y ago
all floats are equally boring.
- vortico 6y ago
  No, -0.f is more interesting than 4.1045857301e23f. It's more interesting because it's more likely to fail a test where others don't.
  - TheRealPomax 6y ago
    That's a bold assumption: let's make sure it's correct by testing against all floats, instead of only the ones you consider interesting, in the process making any float as boring as any other float again.
daveFNbuck 6y ago
That's essentially what fuzzy testing libraries do. You tell the library what types of arguments you want, write your test, and it puts in some known corner cases and a few random values that change each time you run the test (but are saved as new corner cases when the test fails in the better libraries).
rigtorp 6y ago
Yes, I've done this using libFuzzer. Supply custom mutation and crossover functions that generates interesting floating point numbers. Then let the coverage guided fuzzer exercise your code. https://rigtorp.se/fuzzing-floating-point-code/ https://rigtorp.se/fuzzing-floating-point-code/
ramshorns 6y ago

jlebar 6y ago

We did this in XLA. It runs as part of CI, which I always found pretty amazing. https://github.com/tensorflow/tensorflow/blob/bd2608178bca1be95007328fb9e01408a94ba082/tensorflow/compiler/xla/tests/exhaustive_op_test_utils.cc https://github.com/tensorflow/tensorflow/blob/bd2608178bca1b... Found a lot of bugs this way.

NelsonMinar 6y ago

"Several of the ceil functions were implemented by adding 0.5 to the input value and rounding to nearest." What sort of clown would make that mistake and yet think they are capable of writing a floating point library?

Exuma 6y ago
Can you give an example of when this wouldn't work? I'm definitely not qualified to write a float library (lol) but I also can't think of any cases where that wouldn't work?
- brianush1 6y ago
  ceil(5) = round(5.5) = 6 ceil(5) should be 5
  - brucedawson 6y ago
    To expand on this, if you have an integer and you add 0.5 then a round() function has to make a choice, because it is halfway between integers. The IEEE standard says that rounding in general (not just for the round function) should round-to-nearest-even. That means that rounding should go to the closest result, and if there is a tie then it should go to the even number. So, if ceil(x) is implemented as round(x+0.5) then: ceil(4) = round(4.5) = 4 ceil(5) = round(5.5) = 6 The first of those is correct. The second is incorrect.
- zelphirkalt 6y ago
  This is probably about number precision and loosing digits, when you add something with digits only in the beginning like 0.5 to something small like 0.1000000000000000000000000000013124142579275927597245205 or whatever a barely representable floating point number is. Then the result becomes inprecise and rounding to nearest will result in the wrong value. Such effects and similar should be well known to anyone taking on the responsibility to develop a floating point library.
- NelsonMinar 6y ago

zelphirkalt 6y ago

Test them all and waste a huge amount of energy and computational resources! I hope people test like that for really critical stuff, but not for your everyday library for your personal startup or e-business. Not even the big players need this.

zzzeek 6y ago

Crap now I have 67 million failures to fix.

ViralBShah 6y ago

We've certainly tested all Float32 inputs to some functions in Julia's math packages. I'm unable to find the relevant issues this very moment, otherwise would certainly post them.

raphlinus 6y ago
I suspect https://github.com/JuliaLang/julia/pull/29700 https://github.com/JuliaLang/julia/pull/29700 might be one of them. As the original post states, rounding is (a) amenable to hardware acceleration, (b) easy to get wrong or slow, and (c) possible to test exhaustively.

ChuckMcM 6y ago

I love this. This is another in the "Now that we have computers with gazillions of bits of RAM and many cores and gazillions of cycles, what conventional wisdom is no longer wise?" series :-). I suggested to a client that they use Sqlite3 for their company database rather than a disk based database. And they looked at me shocked, as if I had spoken heresy. I pointed out how many records they could store on a 64GB machine and that is a small data center class machine and it sort of boggled their mind. Add read-only copies, shared from a NAS device gave them perfect scalability in terms of queries per second (their sales folks had their own database that they were updating). It was fun to watch their expression go from horror to wonder to excitement when they finally implemented it.

michaelmcmillan 6y ago
Yup! I'm currently running a b2c solution for a client that's grossing ~$3 million/year using Python + sqlite3 behind nginx. I work for another shop that scolds me for even suggesting using something else than k8s. How can a persistent file system possibly scale? Lol. Ironically they're grossing -$3 million/year.
- tonetheman 6y ago
  I am curious about your case is the traffic mostly read heavy? Ha wish you could share more. Interesting stuff.
dragontamer 6y ago
https://yourdatafitsinram.net/ https://yourdatafitsinram.net/ The IBM 53-bit quantum compute result from a ton of hard drives + compute power was also intriguing (https://www.ibm.com/blogs/research/2019/10/on-quantum-supremacy/ https://www.ibm.com/blogs/research/2019/10/on-quantum-suprem...). 128PiB of disk space needed for a 54-qubit computation on their system (hard drives). Fun to see where the supercomputers are still needed, but smaller qubits can be emulated on normal computers. 40 qubit calculations probably can be emulated on a 16TB team of hard drives in RAID0.
kortilla 6y ago

jansan 6y ago

Javascript guy here, can you help me set up that Beowulf cluster so I can start testing my 2^64 cases?

JoeAltmaier 6y ago

My god this issue has been around for too long. In the '80s testing folk learned to shotgun the argument space to test functions, and they've not gone an inch forward since! Or so it seems. A version of the Pentium was released that couldn't divide by 10. And even back then, testing all the mantissas on a simulator would have taken only hours. Instead, they masked, marketed and released a tragically flawed processor. Just test the input space exhaustively, folks! It's a computer; it feels no pain. Work it hard until it fails or works perfectly.

deleted 6y ago
[deleted]
drivebycomment 6y ago
On most CPUs, there's plenty of 2 64-bit operands, and with SIMD instructions, plenty of CPUs have more than 1K bits total as input to the instruction. Good luck exhaustively testing those.

wtn 6y ago

In the early 1980s, truncated floating point calculations caused the Vancouver Stock Exchange index value to drift down for about 23 months, until it was shown as less than half the correct value: https://en.wikipedia.org/wiki/Vancouver_Stock_Exchange#Rounding_errors_on_its_Index_price https://en.wikipedia.org/wiki/Vancouver_Stock_Exchange#Round...

peter303 6y ago

I recall hackers building standard password hashes of all eight character ascii combinations. Then create a reverse lookup table. But thats 400 trillion entries for stand 68 character ascii password subset. 400 trillion hashes is not a big computation in the blockchain era.

SirYandi 6y ago
More info: https://en.m.wikipedia.org/wiki/Rainbow_table https://en.m.wikipedia.org/wiki/Rainbow_table

hinkley 6y ago

I decided to build a scatter plot to demonstrate what I was 95% sure was a data loss bug (preproduction, thankfully) in the initializer for a CSPRNG. I figured plotting a few bytes if output (256x256) would do nicely. I plotted all possible two byte inputs just to make sure it was working, and that ran instantly, so I tried 3 bytes and that was still pretty fast, so I did 4 while I was at lunch or in a meeting. Sure enough, big black splotches from dropping some of the seed data on the floor. Please to be fixing. Old farts know a lot of things have been tried before and how those worked out. But we also think doing a complex operation 4 billion times is going to be irretrievably slow. Always test your assumptions.

uyt 6y ago

I was reading about an implementation of Barrett reduction recently: https://en.wikipedia.org/wiki/Barrett_reduction https://en.wikipedia.org/wiki/Barrett_reduction Specifically one that uses 80-bit long double for calculating `(__int128_t)a * b % M` from this github: https://github.com/kth-competitive-programming/kactl/blob/master/content/number-theory/ModMulLL.h https://github.com/kth-competitive-programming/kactl/blob/ma... But how do you test the correctness of such code? You can't just stress test all possible product of 64 bit numbers! In the end you still gotta go back to good old fashion proofs: https://github.com/kth-competitive-programming/kactl/blob/master/doc/modmul-proof.pdf https://github.com/kth-competitive-programming/kactl/blob/ma...

tagh 6y ago

Makes sense when you're refactoring an existing function with one numeric argument, and which evaluates rapidly.

CGamesPlay 6y ago

1. Given the thoroughness of the write-up, I'm going to assume that the author did not find that any of the discrepancies in the implementations were actually bugs in the reference implementation. He goes as far as to say that he reimplemented the reference implantation so this seems like a safe assumption. Still, it would be nice to confirm. Given that other implementations seem so shoddy, it wouldn't be surprising to find a bug in the reference. 2. "Conventional wisdom Nazis" Maybe the author was just bored of writing at this point, but some justification for why bitwise comparing floats was a correct decision would have been nice. Is it because there's no math involving unrepresentable intermediates? Because there's no integration happening, so no way for errors to accumulate?

brucedawson 6y ago
"Sometimes there really is an answer that is correct, and in those cases anything less than perfection is just sloppy." Anytime there is a correct answer, and it can be determined exactly (+, -, *, /, sqrt) you should expect nothing less than perfection. ceil, round, and floor fall into this category fairly self-evidently, I think.

alkonaut 6y ago

How do you do this for a function with more than one or two inputs? If the domain is 4B to the power five, you are out of luck. Very few of us need to implement "floor" or any other function that takes 1 float input.

jart 6y ago

"Seymour Cray didn't care that 81.0/3.0 did not give exactly 27.0 on the CDC 6000 class machines; and he was universally respected for making the fastest machines around." ──Linus Torvalds

lentil_soup 6y ago

"Round to nearest even is the default IEEE rounding mode. This means that 5.5 rounds to 6, and 6.5 also rounds to 6" Does anyone know the reason behind this? why the nearest even and not a hard cut like "always round up with >= 0.5"?