"For JPL's highest accuracy calculations, which are for interplanetary navigation, we use 3.141592653589793" (15 digits).
How Many Decimals of Pi Do We Really Need? : https://www.jpl.nasa.gov/edu/news/2016/3/16/how-many-decimal...
The precise thing they are good at is dealing with number in a wide range of magnitudes. Where as fixed point numbers can not be used if the magnitudes vary wildly.
You can only use fixed point arithmetic if you know that every intermediate calculation you make will take place in a specific range of precision. E.g. your base precision might be millimeters, so a 32 bit fixed point number is exact up to one millimeter, but can at maximum contain a distance of 2^32-1 millimeters, so around 4.3 billion millimeters. But again you have to keep in mind that this is the maximum value for every intermediate result. E.g. when calculating the distance between two point in 3D space you need a power of 3, so every value you calculate the power of needs to have a value of less than the third root of 4.3 billion.
This makes fixed point arithmetic very hard to correctly implement and requires a very deep analysis of the system, to make sure that that the arithmetic is correct.
I guess the GP comment was discussing that, with this new measurement of pi, we now have enough precision (in pi) to reference a point this small on an object this far away. Once you account for all the other uncertainties in referencing that point, as you mentioned, all that precision in one dimension of the measurement is completely meaningless.
It still feels weird that you'd use an arithmetic with guaranteed imprecision in a field like this, but I can definitely see that, as long as you constrain the scales, it's more than enough.
e.g, with two decimal digits: (2.83 * 0.10) = 0.283, which is stored as 0.28.
You can not in general avoid cancellation, even claiming that is ridiculous. WTF are you even saying.
Or you can just use something like herbie that thinks about it for you: https://herbie.uwplse.org/
Sometimes there are ways to mitigate this, sometimes there aren't. Sometimes you need to precondition, sometimes you need to rearrange, sometimes you need a different algorithm, sometimes you need to normalize, sometimes you need to use a different arithmetic and so on.
For solving linear systems alone, there are definitely thousands of papers dealing with the problems arising from this. For every single algorithm you write and for all data which comes into that algorithm, you need a careful analysis if you want to exclude the potential of significant numerical errors.
Your comment makes it seem like this is a small problem, where you can just look at an algorithm for a time and fix it, this is literally a hundred year research project in numerics.
> For solving linear systems alone, there are definitely thousands of papers dealing with the problems arising from this. For every single algorithm you write and for all data which comes into that algorithm, you need a careful analysis if you want to exclude the potential of significant numerical errors.
It sounds like we agree that cancellation is avoidable with some analysis, and there are hundreds of techniques you can use to deal with it, but mostly it's the ~5 you listed there. And as you suggest, I don't believe this is nearly as significant a problem in the general case as you think it is. A careful error analysis is possible if you care (and if ever you cared, it would be on a spacecraft), and far easier in floating point than in many other number systems, including fixed point number systems.
Numeric systems that truly fix cancellation are incredibly big and heavy, and cannot usually be used for real-time calculations in a generic form. Fixed point certainly doesn't fix cancellation - it introduces precision loss issues on every operation you do that causes a number to go down in magnitude. It is actually harder to design systems in fixed point that avoid massive precision losses than it is in floating point, and the error analysis is much more substantial.
My original comment was about manned space flight in particular. If your application is relatively generic I think it is completely okay, if you are aware of it and mitigate the most pressing issues.
>Numeric systems that truly fix cancellation are incredibly big and heavy, and cannot usually be used for real-time calculations in a generic form.
You can use interval arithmetic, which guarantees that you at least know when cancellation has occurred. Interval arithmetic is fast enough for real time, although it has its own significant drawbacks.
> It is actually harder to design systems in fixed point that avoid massive precision losses than it is in floating point, and the error analysis is much more substantial.
Absolutely. My point was, that a manned space craft, might just be the point to do it.
Everything we have been talking about relates to space flight. In fact, with humans on board, you can afford to be a lot less precise, because they can work around most numerical issues by hand. The Apollo guidance computers, for example, were prone to occasional instances of gimbal lock and numerical instability, and the astronauts just fixed it.
> You can use interval arithmetic, which guarantees that you at least know when cancellation has occurred. Interval arithmetic is fast enough for real time, although it has its own significant drawbacks.
Interval arithmetic does not prevent cancellation. It's just two floating point calculations, both of which are actually less precise than the one you would do otherwise (you don't use default rounding for interval arithmetic, you round the bottom down and the top up). You do know when things have been canceled, but you know that in a floating point calculation anyway if you have done the error analysis.
My overall point here is that NASA isn't missing anything by using floating point instead of using other weird or exotic arithmetic systems. Double-precision floating point combined with a rudimentary error analysis and some algebra is good enough for pretty much everything, and you may not be able to do better at all with fixed point. Designing fixed point algorithms also depends on a very careful analysis of interval ranges and precisions, and often gets you nothing over just using "double", where the error analysis is easier anyway.
If you need to do better than double, there's also double-double arithmetic for your hard parts, which is a similar speed to interval arithmetic and doubles the precision you get beyond double.
There is no general way to mitigate that, you can use numerically superior algorithms or screen your inputs, but these only help in specific cases. There is no general way to avoid this, every algorithm needs to be treated specifically.
You don't. "-" is exact for fixed point unless the operation falls outside the range of valid values.
In floating point, almost every operator (other than subtraction) has precision of the full width of the mantissa minus 0.5 ULPs. All operators are not guaranteed to always be exact, but they are far more precise on average than equivalent operators in fixed point.
Cancellation isn't an issue of exactness, it's an issue of precision.
E.g. (a-b)*c, which is the common example for cancellation, if a and b are very close, can have an unbounded error compared to the result in the real numbers, in floating point. Since all operations besides "/" are exact in fixed point, no error can be introduced by this operation in fixed point (if all operations are representable).
Claiming that fixed and floating point are suffering the same way is just wrong.
For example, suppose your fixed point format is "the integers" and your floating point format has 6 significant digits: if you have real-valued a = 100000.5 and b = 100001.9, both number systems will round a to 100001 and b to 100002. In both cases, (b - a) will be 1 while (b - a) should be 1.4 if done in the reals. That rounding problem exists in fixed point just as much as in floating point. In both systems, the operation that causes the cancellation is itself an exact calculation, but the issue is that it's not precise. Fixed point will just give you 1 in the register while floating point will add a bunch of spurious trailing zeros. Floating point can represent 1.4, though, while fixed point can't. If a and b were represented exactly (a = 100001 and b = 100002 in the reals), there would be no problem in either number system.
The only times that you get better cancellation behavior are when you have more precision to the initial results, which when comparing double precision float to 64-bit fixed point comes when your operands in fixed point have their MSB at the 53rd position or above. That only happens when your dynamic range is so deeply limited that you can't do much math.
When you are thinking about cancellation numerically, exact is a red herring. Precise is what you want to think about.
It depends on how many bits you have. For example compare f128 with Q8.8. Which one do you think would give better astronomical calculation results?
Double-double would be about 31 digits, and quad precision would get you 34.
Single-precision gets you a bit more than 7 digits.
"The 53-bit significand precision gives from 15 to 17 significant decimal digits precision (2−53 ≈ 1.11 × 10−16). If a decimal string with at most 15 significant digits is converted to the IEEE 754 double-precision format, giving a normal number, and then converted back to a decimal string with the same number of digits" (from Wikipedia)
Simulating physical systems to extremely high precision (e.g. more than double precision) in general seems pointless in most situations because of those effects.
By the way, the first Ariane 5 launch blew up because of floating point error, specifically an overflow when converting a 64-bit float to an int. So be careful with floats!
The general rule of thumb in numerical analysis is you need roughly twice the working precision as the output precision. Double-precision floating point has ~16 decimal digits of precision, which means the output should generally be good for ~8 decimal digits; with single-precision, you have ~7 decimal digits of working precision, or about ~3-4 decimal digits of output precision.
In other words, a 32-bit floating-point number doesn't leave with enough useful precision for many cases, whereas a 64-bit floating-point number is good for most use cases.
> That goes back to the Intel 8087 chip, the floating-point coprocessor for the IBM PC. A double-precision real in the 8087 provided ~15 digits of accuracy, because that's the way Berkeley floating-point expert William Kahan designed its number representation. This representation was standardized and became the IEEE 754 floating point standard that almost everyone uses now.
It predates 8087! VAX had 64-bit floats with similar precision to IEEE 754 double precision. There's probably even older uses of 64-ish-bit floating-point types, but my knowledge of computers in the 60's and 70's is pretty poor. I fully expect you'd see similar results on those computers, though: you need enough decimal digits for working precision, and word-sized floating point numbers are just too small to have enough.
The 8087 itself doesn't use double precision types, it uses 80-bit types internally, which have 64 bits of mantissa (or ~19 decimal digits), although the reason for the 80-bit type is primarily to get higher precision for intermediate results on implementing transcendental functions.
I don’t think Kahan had a direct part in the design of the 8087. https://en.wikipedia.org/wiki/Intel_8087#Design_and_developm... agrees, saying “Palmer credited William Kahan's writings on floating point as a significant influence on their design.”
https://ieeemilestones.ethw.org/w/images/7/7f/Wk_an_intervie...
Also see https://math.berkeley.edu/news/congratulations-professor-wil...
And Kahan's Turing award: "During a long and productive relationship with Intel he specified the design for its floating-point arithmetic on several chips starting with the 8087" https://amturing.acm.org/award_winners/kahan_1023746.cfm
IOW, the likely reason for not storing more mantissa bits than that is someone involved in designing the 8087 determined that even NASA doesn't need more precision than that.
“atan(1) * 4”
casts to double?
- I wonder if this cast is always correct in C [ie.: math.h], no matter the datatype and/or the number base?
Floating point arithmetic is deterministic. As long as it is implemented as specified atan(1) has to give the floating point number which is the closest approximation to the real number pi/4 (in the current rounding mode), the multiplication by 4 means that precision can be lost and potentially your result is no longer the closest possible approximation to pi.
What? This is not true at all. The standards specifies adherence to IEEE 754 arithmetic.
You can read the standard here: https://www.open-std.org/jtc1/sc22/wg14/www/docs/n1570.pdf
Page 507 for adherence to number formats. Page 517 for atan adhering to IEEE 754 specification for the functions defined therein, which guarantees best possible results for individual operations.
Any C implementation where atan gives a result which is inconsistent with IEEE 754 specification does not adhere to the standard.
> also, no OS provided libm produces correctly rounded results for all inputs.
Every IEEE 754 conforming library does adhere to the best possible rounding guarantee. If you have any evidence to the contrary that would be a disaster and should be reported to the vendor of that library ASAP.
Can you provide some function and some input which violates the IEEE 754 guarantee together with the specific library and version? Or are you just making stuff up?
Even beyond transcendental functions, 754 isn't deterministic in practice because implementations have choices that aren't always equivalent. Using FMA vs separate multiplication and addition leads to different results in real programs, even though both methods are individually deterministic.
But then it doesn't conform to the standard. It is pretty unambiguous on that point.
From Section 9.2:
"A conforming operation shall return results correctly rounded for the applicable rounding direction for all operands in its domain."
I do not see how two conforming implementations can differ in results.
>Using FMA vs separate multiplication and addition leads to different results in real programs, even though both methods are individually deterministic.
Obviously. I never claimed that the arithmetic was invariant under transformations which change floating point operations, but are equivalent for real numbers. That would be ridiculous.
Is there actually an example of two programs performing identical operations under the same environment that give different results where both implementations conform to the standard?
>Even beyond transcendental functions, 754 isn't deterministic in practice because implementations have choices that aren't always equivalent.
Could you give an example? Where are implementations allowed to differ? And are these cases relevant, in the sense that identical operations lead to differing results? Or do they just relate to error handling and signaling.
> 9. Recommended operations
> Clause 5 completely specifies the operations required for all supported arithmetic formats. This clause specifies additional operations, recommended for all supported arithmetic formats.
Hyperbolic tan is in the list of recommended functions, and yet: https://github.com/numpy/numpy/issues/9187
Who cares? The C standard for math.h requires these functions to be present as specified. They are specified to round correctly, the C standard specifies them to be present as specified, therefore the C standard specifies them as present and correctly rounded. I literally quoted the relevant sections, there are no conforming C specification which give different results.
>Hyperbolic tan is in the list of recommended functions, and yet: https://github.com/numpy/numpy/issues/9187
Any evidence whatsoever that this is caused by two differing implementations of tanh, which BOTH conform to the IEEE 754 standard?
Everyone is free to write their own tanh, it is totally irrelevant what numpy gives, unless there are calls to two standard confirming tanh function which for the same datatype produce different results.
Forgive me, but I cannot see that in the document sections you point out. The closest I can see is F.10-3, on page 517, but my reading of that is that it only applies to the Special cases (i.e values in Section 9.2.1), not the full domain.
In fact, my reading of F.10-10 (page 518) suggests that a conforming implementation does not even have to honor the rounding mode.
I'm not aware of any libm implementations that will guarantee correct rounding across all inputs for all types. I'm aware of a few libm's that will guarantee that for floats (e.g. rlibm: https://people.cs.rutgers.edu/~sn349/rlibm/ ), but these are not common.
Genuinely a bit shocked by this.
That's not the current version of C. The best document right now is https://www.open-std.org/jtc1/sc22/wg14/www/docs/n3220.pdf (which is C23 with one or two editorial fixes).
Note that Annex §F.3¶20 says
> However, correct rounding, which ISO/IEC 60559 specifies for its operations, is not required for the C functions in the table.
which proceeds to list most of the functions in IEEE 754-2008/2019 section 9.2.
A C compiler which claims conformance to IEEE 754 need not correctly round those functions. Most C libraries do not correctly round all of these functions: https://members.loria.fr/PZimmermann/papers/accuracy.pdf
(in summary, llvm-libc correctly rounds all functions, as it's explicitly borrowing from one of the correctly-rounded efforts; of the other implementations, Intel's library usually gets the closest, but not always).
(Insert long rant about icc enabling daz/ftz by default.)
It is a bit shift to the right, so where do the new bits come from? Why would the two new bits be the correct ones?
This is about the approximation to pi not the approximation to float(atan(1))*4, it is exact (but irrelevant) for the later, for the former you loose two bits, so you have a 25% chance of correctly rounding towards pi.