One does not simply calculate the absolute value
habr.com
habr.com
Of course, the most famous bit manipulation of floats is the fast inverse square root trick [0]. The bits of a float as an integer literal are proportional to its log2 value (plus a constant), which is handy for computing a reciprocal square root in logspace, which when converted back to a float gets exponentiated, yielding the answer (after a couple more refinement steps). If you actually take the time to manually hammer out the bitwise manipulations, the fast inverse square root is quite simple, but it is nonetheless regarded as black magic, likely due to the totally non-explanatory comments in the Quake source code.
Funnily enough, directly following the fast inverse square root in the Quake source code is a function to compute absolute values by zeroing the float sign bit [1], just as described in the article.
[0] https://en.wikipedia.org/wiki/Fast_inverse_square_root#Alias...
[1] https://github.com/id-Software/Quake-III-Arena/blob/dbe4ddb1...
If a similar article were written about how shifting an integer to the right by one bit is equivalent to dividing it by two, I highly suspect it would attract no attention on this site, aside from a couple dismissive comments grousing about why such a basic fact needs its own article.
I don't know what to take from the fact that Java stdlib didn't use it until Java 18? It does seem like it wasn't obvious to at least someone who knew something.
> If a similar article were written about how shifting an integer to the right by one bit is equivalent to dividing it by two, I highly suspect it would attract no attention on this site,
Nah, it would be full of comments of people saying "well of course I knew that! doesn't everyone who's anyone? I can't believe how many newbies don't though!" :)
Mainstream? I can't think of one. Googling the question yields a stackoverflow question from 2010 that lists the 5 or 6 major architectures (x86, PowerPc, SPARC, Arm,.. ) and about a dozen more as supporting x86. So it's a fairly safe assumption to make.
One area where vanilla IEEE754 isn't used that immediately came to my mind is high performance ML at companies like facebook and Google, the generic range-precision tradeoff that FP offers is inefficient and they ended up designing several new formats with smaller sizes or larger ranges. Googling "Alternatives to FP in ML" or "Custom FP format <company-name>" yields a lot.
Another "niche" area is finance, insurance, taxes and anything nearby, fixed point in those area is, I hear, a must. Because all the usual invariants of high school algebra that FP breaks actually have significance. (imagine the well-known FP flaw "x+epsilon == x, with epsilon non-zero" playing out in these contexts.)
In any case, number portability in a language like C++ is a breeze. You design a myNumber class and overload all relevant operators and that's it. Compiler will make sure to inline the calls if they're a one-line simple floating point operations, but otherwise you can freely change the implementation of the class and the calculations never look different.
In particular, double-entry bookkeeping (standard practice for centuries) uses the accounting equation [1] as a checksum. If the two sides of the equation are ever unequal, then there's been a mistake somewhere that needs to be corrected.
If you add the same value to a large asset account and a small liability account, for example, you can run into problems: The internal rounding performed by floating point arithmetic can mean that the increase to each side of the equation is slightly different, throwing the system out of balance.
Lots of medical devices support floating point format from IEEE-11073 instead. I'd guess there are a significant number of them in existence and they need to be done well.
MSP430-style floats are common in industrial automation and such, since many sensors either use or used to use MSP430s to collect data. And, of course, in that field, legacy carries more inertia than a 50,000-ton press....
Newer DSPs tend to use IEEE-754 floats, but this was not always the case. Older DSPs used every trick in the book to go faster, and "compatibility with standard floating-point formats" was never in their product requirements. So you can guess at what manner of horrors came out of that! Fortunately, their customers are mostly the sort of people who (1) didn't mind changing to standard floats, since it's not like things were stable before; and (2) are very, very comfortable writing this stuff directly in assembly. So they definitely know whether their architecture has an instruction to take the absolute value, needs a short branching program, or just needs a single bit-operation.
Not only that but, at least for Fast Inverse Square Root, it's not advantageous to carry out the bit twiddling since the rsqrtss instruction is faster and more precise, and vectorized instructions are faster still for large strings of floats.
I watched a video[0] a while back about that. Once that coincidence(?) was explained, Q_rsqrt finally made sense.
[0]: https://youtu.be/p8u_k2LIZyo (~20 minutes)
-100…000 = ~100…000 + 1 = 011…111 + 1 = 100…000
where ~ is the bitwise complement. Even worse, in C and C++ specifically, because overflow of signed integers is undefined behavior, so is `abs(INT_MIN)`! To be able to represent the absolute value, you need to either return the corresponding unsigned type (if your language has them), promote to a wider type, or return a tuple `(carry_bit, result)`.Checked returns a Result enum that is an error on overflow, wrapping silently wraps (as in your answer), overflowing wraps but has a boolean indicator whether a wrap has occurred, saturating bounds at the types limits rather than wrapping and unsigned abs always gives the correct result albeit as a different type.
The unqualified abs panics on overflow on debug but silently wraps on release mode, like all Rust arithmetic code does by default.
Regarding signed zero, read this paper from Kahan: https://homes.cs.washington.edu/~ztatlock/599z-17sp/papers/b...
0/(-x) would return -0, for all x > 0
Likewise, underflows on negative numbers would also return -0.
I find this one a little weird, but it seems to be true.
The more normal way to get -0 would be to divide a positive number by negative infinity, or a negative number by positive infinity.
Wouldn't max( x,-x ) cover everything ?Try it and you'll see.
return a >= b ? a : b;
But if x is -0.0, you're back to the same problem that -0.0 >= +0.0 actually turns out to be true and not false, hence the return value will be -0.0.Looking at the actual implementation[0], doubleToRawLongBits is used again here to take care of this edge case. So it would work, but a lot slower than the version given in the article.
0: https://github.com/openjdk/jdk/blob/36e2ddad4d2ef3ce27475af6...
(changing greater equal to just greater). Would the problem then still persist?
The problem is that 0.0 == -0.0. As long as you don't take care of signed zero, you cannot guarantee that abs(x) has no negative sign.
> return Double.longBitsToDouble(
> Double.doubleToRawLongBits(value) & 0x7fffffffffffffffL);
Technically, this will change NaN values into different ones.I don't know if anything does but it could certainly be useful. Imagine if a NaN said exactly where in your code it was created.
I'm not sure how that would be useful. A debugger should have better ways of identifying NaN creation and I don't see anything else using a pointer to code in production. It seems cool though.
It would also lead to problems like: x = 5.0 * y y is NaN with some address as the NaN type. NaN * 5.0 => NaN. Which address (the original or the * 5.0 line) is used?
If you log it, you could trace it back to the line of code at your leisure.
> It would also lead to problems like: x = 5.0 * y y is NaN with some address as the NaN type. NaN * 5.0 => NaN. Which address (the original or the * 5.0 line) is used?
I don't mean to be rude, but isn't the answer to that question extremely obvious?
The behavior that's much more useful and enshrined in the IEEE spec is that you keep the original NaN.
> The behavior that's much more useful and enshrined in the IEEE spec is that you keep the original NaN.
I didn't realize it was in the spec. I just checked, and you're right! That does however seem strange - X was never set by operation NNNN, Y was.
> That does however seem strange - X was never set by operation NNNN, Y was.
If you're using NaN payloads, it would almost never be helpful for 95% of your NaNs to turn into a generic "used a NaN in arithmetic" version.
It makes no sense to bring this up in interview questions. Just about every programmer knows that IEEE 754 contains some quirks. If you want to know the boring details, just look them up in the manual.
For some industries, ensuring your code doesn't fail on corner cases is very important, and questions about quirks show whether a candidate cares about corner cases and their comprehension of them.
Given that the limitations of IEEE754 are the same in in virtually all programming languages, I suspect it is actually really good to use floating point quirks to probe a candidate on how skillful they are at programming.
I roll my eyes at some stack overflow answers for their complete failure to consider even common corner cases like NaN and Infinity.
I have certainly been caught out by my own mistakes in my own code with with integer and floating point errors, and I think I am more careful than most.
The real quirks lie elsewhere, in the hardware implementations and rounding rules (e.g. does FMA have the same precision of multiply-then-add? why does the order of the elements of a sum make a difference? etc.). And once you know these you can treat floating points with familiarity and expect your operations to be predictable to their last bit.
If I were hiring for a developer of FDTD simulations I would expect them to be familiar with things like -0.0, subnormal numbers and NaN types and to be able to learn if necessary how to go around the murky stuff like rounding. Floats shouldn't be mysterious, they are at the end 16/32/64/80/128 bit variables with deterministic operations defined on them.
But as they are corner cases, you are very unlikely to just know them all from memory, so you would need to look them up anyways.
When I need absolute value, I usually do it like that in C++:
return _mm_andnot_pd( _mm_set1_pd( -0.0 ), vec );Of course the answer to "have you covered all weird corner cases" is always 'no', but the article addresses the issue of different NaNs:
> […] However, it turns out that doubleToLongBits is also not entirely trivial, because it canonicalizes NaNs. There are many ways to encode a not-a-number as a double, but only one of them is canonical. These different NaNs are even more similar twins, they can be distinguished neither via Double.compare, nor via any operation. Even string representation is the same. But they look different in computer memory. To avoid surprises, doubleToLongBits converts any NaN to the canonical form, which is encoded in long as 0x7ff8000000000000L. Of course, this procedure adds more conditions, which we do not need here.
> What can we do? It appears that there's another method doubleToRawLongBits. It doesn't do any smart conversion of NaN and just returns exactly the same bit representation: […]
The sign bit is totally independent of the value being encoded, whether that’s a number, q/sNaN, or infinity.
Nitpick: we are not flipping the sign bit, but rather setting it to 0.
chs() ("flip the sign bit") is the other operation that's defined as effectively just an operation on the representation rather than a "do a real floating point computation".
https://github.com/pedrocr/imagepipe/blob/12269f04ac2fa0d165...
The motivation for the article is this change in OpenJDK, where they got a 10% performance improvement by switching from the single-branch implementation to the branchless one:
if (value < 0) return -value; else if (value > 0) return value; else return +0.0
This must be a Java-ism. In most languages I've used you get a division by zero error or NaN.
It exploits the division by zero sign logic of IEEE-754, so that the min/max operations used also handles the edge cases where the ray is parallel to an axis.
IEEE754 section 7.3:
The divideByZero exception shall be signaled if and only if an exact infinite result is defined for an operation on finite operands. The default result of divideByZero shall be an ∞ correctly signed according to the operation: ― For division, when the divisor is zero and the dividend is a finite non-zero number, the sign of the infinity is the exclusive OR of the operands’ signs (see 6.3). ― For logB(0) when logBFormat is a floating-point format, the sign of the infinity is minus (−∞).
You'd get NaN from 0.0/0.0 or division by zero if you did 1/0 as integer division.
Change my mind.
https://people.freebsd.org/~das/kahan86branch.pdf from the father of ieee-754 is illuminating.
And given that we have to implement some sort of approximation to the usual math on a closed subset of the reals (rationals, even), why not extend it with some "numbers" (negative zero, plus/minus infinity, NaNs) that, yes, enable programming convenience?
https://hackernoon.com/negative-zero-bbd5fd790af3
But, bottom line, Kahan (not only the father of IEEE 754, but also, with Golub, father of the SVD!) has thought about it a lot, and though he called the signed zero "a pain in the ass", he also said "There really wasn’t a way around that and you were stuck with it."
1/ (+0) = +∞ but 1/ (–0) = –∞
but not that there be two, say, 3s so that
1/ (3⁺ – 3) = +∞ but 1/ (3⁻ – 3) = –∞
It's not that it's important there be two, just every number already had 2 except for 0.
Curiously, Go gives me a compiler error if I divide constant literals, but assigning 0.0 to a variable and dividing by it does work. On top of that, assigning -0.0 to a float variable leads to it having the value +0.0, which I would call a borderline compiler bug. Both are probably a consequence of Go constant literals not being IEEE-754 floats but some arbitrary precision number types, but it's very unexpected behavior.
Python and C#.
I hope others will take the entire thread into consideration before downvoting.
Or it raises an exception, which once disabled a navy warship for hours:
> Now decommissioned, the USS Yorktown was among the first warhips extensively computerized to reduce crew (by 10% to 374) and costs (by $2.8 million per year).
> On 21 Sept. 1997, the Yorktown was maneuvering off the coast of Cape Charles, VA, when a crewman accidentally ENTERed a blank field into a data base. The blank was treated as a zero and caused a Divide-by-Zero Exception which the data-base program could not handle. It aborted to the operating system, Microsoft Windows NT 4.0, which crashed, bringing down all the ship’s LAN consoles and miniature remote terminals.
> The Yorktown was paralyzed for 2 3/4 hours, unable to control steering, engines or weapons, until the operating system had been re-booted. Fortunately the Yorktown was not in combat nor in crowded shipping lanes.
> If IEEE 754’s default had been in force, the division by zero would have insinuated into the data-base an ∞ and/or NaN, which would have been detected afterwards without a crash.
http://people.eecs.berkeley.edu/~wkahan/Boulder.pdf
https://www.wired.com/1998/07/sunk-by-windows-nt/
EDIT: It was 2 3/4 hours, not 23 hours. I copied/pasted carelessly from the PDF.
1. Interface accepted blank value for a form when a numeric value is required.
2. Application and database permitted a null or coerced a zero when the value is required and cannot be zero. You've got to sanitize your input.
3. Application did not check for or gracefully handle something as common as a division by zero exception. Instead it proceeded to a buffer overrun and crashed.
4. Application crashing apparently caused an OS crash? Yes, okay, it's NT 4.0 and a buffer overrun, but that still should not happen. Protected mode is a thing.
5. There is no redundancy in a system that manages the steering, weapons, and engines of a warship.
6. A warship was unable to restore the single CNC system to working order for 23 hours.
Furthermore... there's no evidence that the software would have behaved more sensibly or more controllably if it had a value of infinity instead of an error. It shouldn't have had a zero in the first place.