Visual C++ rounds some decimal to double-precision floating-point conversions incorrectly, but it’s not alone; the gcc C compiler and the glibc strtod() function do the same. In this article, I’ll show examples of incorrect conversions in gcc and glibc, and I’ll present a C program that demonstrates the errors.
(To understand my methodology and terminology, please read my article “Incorrectly Rounded Conversions in Visual C++.”)
(Read my December 2013 update to this article: real.c Rounding Is Perfect (GCC Now Converts Correctly).)
(Read my October 2013 update to this article: GCC Conversions Are Incorrect, Architecture or Otherwise.)
(Read my September 2013 update to this article: Correctly Rounded Conversions in GCC and GLIBC.)
My Analysis
I ran my tests on an Intel Core Duo processor running Ubuntu Linux version 10.04, and I used the gcc C compiler version 4.4.3 and the embedded glibc library (eglibc) version 2.11.1. (As far as I can tell, eglibc uses the same conversion code as glibc, so I will refer to eglibc as glibc for the rest of this article.) (Update: I reran my tests on glibc, getting the same results.)
I did automated testing of glibc strtod() string conversions and manual testing of gcc floating-point literal conversions. This is what found in my automated testing of strtod():
- Conversions were done incorrectly only about four ten-thousandths of a percent of the time.
- No conversion was incorrect by more than one unit in the last place (ULP).
- Incorrect conversions included both “hard” and “easy” cases.
- There were incorrect conversions covering all three cases: “halfway”, “greater than halfway”, and “less than halfway.”
My manual testing of gcc was not as comprehensive, of course, but I did manage to find some examples of incorrect conversions.
For the remainder of this article, I will discuss six examples of incorrect conversions — examples uncovered by my automated and manual testing, as well as examples found in Bugzilla bug reports for gcc and glibc.
Incorrect Conversion by Both GCC and GLIBC
Here’s an example that’s converted incorrectly by both gcc and glibc — to the same incorrect value:
Example 1 (Halfway Case)
0.500000000000000166533453693773481063544750213623046875,
which equals 2-1 + 2-53 + 2-54, converts to this binary number:
0.100000000000000000000000000000000000000000000000000011
It has 54 significant bits, with bits 53 and 54 equal to 1. The correctly rounded value — in binary scientific notation — is
1.000000000000000000000000000000000000000000000000001 x 2-1
gcc and glibc strtod() both compute the converted value as
1.0000000000000000000000000000000000000000000000000001 x 2-1
which is one ULP below the correctly rounded value.
(I constructed this example by hand.)
Incorrect Conversions by GLIBC Only
Here are three examples that are converted incorrectly by glibc but converted correctly by gcc:
Example 2 (Halfway Case)
3.518437208883201171875e13 converts to this binary number:
1000000000000000000000000000000000000000000000.00000011
It has 54 significant bits, with bits 53 and 54 equal to 1. The correctly rounded value is
1.000000000000000000000000000000000000000000000000001 x 245
glibc strtod() computes the converted value as
1.0000000000000000000000000000000000000000000000000001 x 245
which is one ULP below the correctly rounded value.
(This example was taken from “Bugzilla Bug 3479: Incorrect rounding in strtod().”)
Example 3 (Greater than Halfway Case)
62.5364939768271845828 converts to this binary number, which has a non-terminating binary fraction portion:
111110.1000100101010111101010110101010011111001010001110000000000001…
Bits 53 and 54 are equal to 1, so the correctly rounded value is
1.11110100010010101011110101011010101001111100101001 x 25
glibc strtod() computes the converted value as
1.1111010001001010101111010101101010100111110010100011 x 25
which is one ULP below the correctly rounded value.
(I discovered this example during random testing.)
Example 4 (Less than Halfway Case)
8.10109172351e-10 converts to this non-terminating binary fraction:
0.00000000000000000000000000000011011110101011100101110101110111100001111110100001100011111111111110…
Significant bit 54 is 0, so it should be rounded down to get the correctly rounded value
1.10111101010111001011101011101111000011111101000011 x 2-31
glibc strtod() computes the converted value as
1.1011110101011100101110101110111100001111110100001101 x 2-31
which is one ULP above the correctly rounded value.
(I discovered this example during random testing.)
Incorrect Conversions by GCC Only
Here are two examples that are converted incorrectly by gcc but converted correctly by glibc:
Example 5 (Halfway Case)
The number
1.50000000000000011102230246251565404236316680908203125
which equals 1 + 2-1 + 2-53, converts to this binary number:
1.10000000000000000000000000000000000000000000000000001
It rounds correctly to 1.1, since bit 53 is 0.
gcc converts it to
1.1000000000000000000000000000000000000000000000000001
which is one ULP above the correctly rounded value.
(I constructed this example by hand.)
Example 6 (Less than Halfway Case)
9007199254740991.4999999999999999999999999999999995 converts to this binary number, which has a non-terminating binary fraction portion:
11111111111111111111111111111111111111111111111111111.011111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111110…
Bit 54 is 0, so it converts correctly to
1.1111111111111111111111111111111111111111111111111111 x 252
gcc converts it to
1.0 x 253
which is one ULP above the correctly rounded value.
(This example was taken from “Bugzilla Bug 21718: real.c rounding not perfect.”)
A C Program that Reports the Behavior
Here’s a C program you can use to confirm incorrect rounding on your own Linux system:
// gcc -o examples examples.c
#include <stdio.h>
#include <stdlib.h>
#define NUM_CONVERT 6
int main (void)
{
unsigned int i;
char decimal[NUM_CONVERT][60] = {
"0.500000000000000166533453693773481063544750213623046875",
"3.518437208883201171875e13",
"62.5364939768271845828",
"8.10109172351e-10",
"1.50000000000000011102230246251565404236316680908203125",
"9007199254740991.4999999999999999999999999999999995"
};
double gcc_conversion[NUM_CONVERT] = {
0.500000000000000166533453693773481063544750213623046875,
3.518437208883201171875e13,
62.5364939768271845828,
8.10109172351e-10,
1.50000000000000011102230246251565404236316680908203125,
9007199254740991.4999999999999999999999999999999995
};
char correct[NUM_CONVERT][30] = {
"0x1.0000000000002p-1",
"0x1.0000000000002p+45",
"0x1.f44abd5aa7ca4p+5",
"0x1.bd5cbaef0fd0cp-31",
"0x1.8p+0",
"0x1.fffffffffffffp+52"
};
for (i=0; i < NUM_CONVERT; i++)
{
printf ("%s\n",decimal[i]);
printf (" Correct = %s\n",correct[i]);
printf (" gcc = %a\n",gcc_conversion[i]);
printf (" strtod = %a\n\n",strtod(decimal[i],NULL));
}
}
This is the output on my system (I shaded the incorrect conversions):
0.500000000000000166533453693773481063544750213623046875 Correct = 0x1.0000000000002p-1 gcc = 0x1.0000000000001p-1 strtod = 0x1.0000000000001p-1 3.518437208883201171875e13 Correct = 0x1.0000000000002p+45 gcc = 0x1.0000000000002p+45 strtod = 0x1.0000000000001p+45 62.5364939768271845828 Correct = 0x1.f44abd5aa7ca4p+5 gcc = 0x1.f44abd5aa7ca4p+5 strtod = 0x1.f44abd5aa7ca3p+5 8.10109172351e-10 Correct = 0x1.bd5cbaef0fd0cp-31 gcc = 0x1.bd5cbaef0fd0cp-31 strtod = 0x1.bd5cbaef0fd0dp-31 1.50000000000000011102230246251565404236316680908203125 Correct = 0x1.8p+0 gcc = 0x1.8000000000001p+0 strtod = 0x1.8p+0 9007199254740991.4999999999999999999999999999999995 Correct = 0x1.fffffffffffffp+52 gcc = 0x1p+53 strtod = 0x1.fffffffffffffp+52
The values for the correct conversions were obtained using David Gay’s strtod() function. I ran a separate program to compute them and then hard coded them in this program. (I could also have constructed the correct hex formats by hand, using my analysis in the examples above.)
If you try this on your system, let me know what happens; I expect the same output. Please include the type of Linux, and version numbers if possible.
Bug Reports
There are existing bug reports for gcc and glibc:
It’s not clear whether these bugs will (or should?) be fixed; if they are, I’d imagine David Gay’s strtod() function would be used to replace the current code.
I would guess this problem affects the other GCC compilers, although I haven’t tested.
GCC/GLIBC Vs. Visual C++
Here are the differences I found between gcc/glibc and Visual C++ in my two analyses:
- glibc strtod() converted incorrectly about two orders of magnitude less frequently than Visual C++ strtod().
- I found no integer values that rounded incorrectly in gcc or glibc, unlike in Visual C++.
- gcc and glibc strtod() can convert the same decimal number differently; I found no such example in Visual C++ and its strtod() function.
- Of the five examples Visual C++ converts incorrectly, only one — example 2 — is converted incorrectly by gcc or glibc (in this case, by both). Of the six examples in this article converted incorrectly by gcc or glibc, two are converted incorrectly by Visual C++: examples 1 and 2.
- When gcc or glibc rounds halfway cases incorrectly, it does so either up or down. I found Visual C++ to consistently round incorrectly down.
- I found an incorrect conversion for a “less than halfway case” in gcc/glibc; I found no such example in Visual C++.
Summary
I’ve demonstrated that a decimal floating-point number is not guaranteed to convert to its closest binary floating-point number. Also, I’ve shown that a decimal floating-point number may be converted to different binary floating-point numbers, depending on which compiler or run-time library is used.
Please check out the discussion of my article on reddit; in particular, see the different results people are getting with different versions of gcc and glibc and when running on different processors.
Hi Rick,
I honestly think your expectations are wrong. Please do educate me if I am the one who’s wrong instead.
Looking at the C99 standard (if you don’t have it, you can download the latest Committee Draft N1124 from [0]), 6.4.4.2 “Floating constants” paragraph 3 says,
“[…] For decimal floating constants […] the result is either the nearest representable value, or the larger or smaller representable value immediately adjacent to the nearest representable value, chosen in an implementation-defined manner. […]”
If you evaluate your test results with these semantics in mind, do you still think the complier / C library is buggy? Ultimately, you did say “[n]o conversion was incorrect by more than one unit in the last place (ULP)”.
Furthermore, you state
“[…] I’ve shown that a decimal floating-point number may be converted to different binary floating-point numbers, depending on which compiler or run-time library is used.”
This too fits in “implementation-defined manner”. (That is, the implementation can chose any allowed way, but it must document its choice).
Cheers,
lacos
[0] http://www.open-std.org/jtc1/sc22/wg14/www/standards.html
Hi lacos,
Yes, I am familiar with that C99 text (in fact, it popped up recently on stackoverflow.com). The bottom line is you’re right — it’s not a bug. I don’t disagree, but that doesn’t mean correct rounding shouldn’t be the goal. After all, there are open “bug” reports in gcc and glibc, which means they think it worthy of consideration. And, there is open source code available that does correct rounding.
BTW, the terms “correct” and “incorrect” might be misleading, but they’re not my terms. A correctly rounded conversion means the closer of two binary floating-point numbers is chosen (I’m assuming the default — “round to nearest even” — but other rounding modes can be applied.)
Thanks for the comment.
clang 1.1.27 does the rounding right (eglibc still fails, that was expected), but I’m worried because it doesn’t mimic gcc behaviour when using gcc_conversion, and that will lead to obscure bugs to someone.
Note: The source program was not modified in any way, gcc in the program refers to clang instead.
esanchma@esanchma-virtualized-linux:~/Escritorio$ clang -o examples-clang examples.c
esanchma@esanchma-virtualized-linux:~/Escritorio$ ./examples-clang
0.500000000000000166533453693773481063544750213623046875
Correct = 0x1.0000000000002p-1
gcc = 0x1.0000000000002p-1
strtod = 0x1.0000000000001p-1
3.518437208883201171875e13
Correct = 0x1.0000000000002p+45
gcc = 0x1.0000000000002p+45
strtod = 0x1.0000000000001p+45
62.5364939768271845828
Correct = 0x1.f44abd5aa7ca4p+5
gcc = 0x1.f44abd5aa7ca4p+5
strtod = 0x1.f44abd5aa7ca3p+5
8.10109172351e-10
Correct = 0x1.bd5cbaef0fd0cp-31
gcc = 0x1.bd5cbaef0fd0cp-31
strtod = 0x1.bd5cbaef0fd0dp-31
1.50000000000000011102230246251565404236316680908203125
Correct = 0x1.8p+0
gcc = 0x1.8p+0
strtod = 0x1.8p+0
9007199254740991.4999999999999999999999999999999995
Correct = 0x1.fffffffffffffp+52
gcc = 0x1.fffffffffffffp+52
strtod = 0x1.fffffffffffffp+52
esanchma@esanchma-virtualized-linux:~/Escritorio$ clang –version
clang version 1.1 (branches/release_27)
Target: i386-pc-linux-gnu
Thread model: posix
eddie,
What’s really interesting is that your results are different than the five other sets of results I have (my results above and the 4 sets of results reported on reddit).
Thanks for the data.
It’s stock ubuntu 10.4 with the default clang in ubuntu multiverse: clang-2.7-0ubuntu1
(1) C99 and IEEE 754-1985 are somewhat lax on conversion requirements, but IEEE 754-2008 requires correct rounding across the entire exponent range and for at least 20 digits in the case of Binary64.
(2) The bug I reported for glibc several years ago actually violates 754-1985.
(3) The code in glibc is perfectly capable of doing the correct conversion. The problem is an ill-advised shortcut that claims n digits (30? I don’t remember now) are sufficient, when in fact the threshold is much larger (it depends on the exponent and may be as large as 752 bits).
Michel.
Michel,
My Example 2 is from your bug report (it’s a halfway case that should round up, not down). Can you elaborate on why it “violates 754-1985”? Thanks.
(BTW, it’s 20 digits that the shortcut uses.)
This reply is a bit late, but I was not following this discussion, and only got
here because of your reply on stds-754 (Re: Conflicts between C standard
and 754-2008).
My bugzilla #3479 violates 754-1985 because it involves fewer than 17
digits and has a middling exponent, within the range where correct rounding
is required.
I’m a bit annoyed that there has been no progress in four years, even though
the fix should be pretty trivial (don’t take the invalid shortcut) — but I’m not
sufficiently plugged in to the glibc maintenance community to provide my own
bug fix.
Michel
Michel,
The example in Bugzilla #3479, 3.518437208883201171875E+013, is 22 digits. Is that required to be correctly rounded?
My mistake! 754-1985 does not require correct rounding for 22 digits.
I misremembered the 22-digit size of the #3479 example, probably
because I’ve also investigated a number of 17-digit (and even 16-digit)
failures on non-Linux systems.
It is still a glibc bug however, because the comments, and the use of
gmp under the covers, clearly suggest that the intention was to perform
correct rounding…
Given that glibc uses gmp to convert at least 20 digits, I suspect that it
does not have bugs within the IEEE 754-1985 correct-rounding range.
(For the record, that range is numbers representable as M*10**N with
integers M and N such that |M| < 1e17 and |N| < 28.)
Of course, by now there is IEEE 754-2008 to contend with! If binary64
is the widest supported FP format, then 20 correctly-rounded digits would
be conforming, provided longer digit strings are first rounded in decimal
to 20 digits, which I suspect is NOT being done. (That's understandable,
because that particular requirement was only formulated in late 2006 if
I remember correctly.)
So perhaps I should torture glibc's strtod() a bit more on post-2008
versions…
Michel.
gcc 4.4.5:
0.500000000000000166533453693773481063544750213623046875
Correct = 0x1.0000000000002p-1
gcc = 0x1.0000000000002p-1
strtod = 0x1.0000000000001p-1
3.518437208883201171875e13
Correct = 0x1.0000000000002p+45
gcc = 0x1.0000000000002p+45
strtod = 0x1.0000000000001p+45
62.5364939768271845828
Correct = 0x1.f44abd5aa7ca4p+5
gcc = 0x1.f44abd5aa7ca4p+5
strtod = 0x1.f44abd5aa7ca3p+5
8.10109172351e-10
Correct = 0x1.bd5cbaef0fd0cp-31
gcc = 0x1.bd5cbaef0fd0cp-31
strtod = 0x1.bd5cbaef0fd0cp-31
1.50000000000000011102230246251565404236316680908203125
Correct = 0x1.8p+0
gcc = 0x1.8p+0
strtod = 0x1.8p+0
9007199254740991.4999999999999999999999999999999995
Correct = 0x1.fffffffffffffp+52
gcc = 0x1.fffffffffffffp+52
strtod = 0x1.fffffffffffffp+52
@Bob,
Thanks for the report. Your results match the results of an AMD64 user on reddit. What system or you on?
Heh, just found this article. Happy to see clang gets it right; I contributed their decimal->binary code 🙂
I submitted the GCC bug 6 years ago; shame there’s little progress.
hmm, quite interestingly, my GCC 4.7.1 x64 on an i5 system running Win7Profx64 (corresponding package using MinGW was downloaded from equation.com) even seems to treat all the doubles as floats:
0.500000000000000166533453693773481063544750213623046875
Correct = 0x1.0000000000002p-1
gcc = 0x1.000000p-1
strtod = 0x1.000000p-1
3.518437208883201171875e13
Correct = 0x1.0000000000002p+45
gcc = 0x1.000000p+45
strtod = 0x1.000000p+45
62.5364939768271845828
Correct = 0x1.f44abd5aa7ca4p+5
gcc = 0x1.f44abdp+5
strtod = 0x1.f44abdp+5
8.10109172351e-10
Correct = 0x1.bd5cbaef0fd0cp-31
gcc = 0x1.bd5cbbp-31
strtod = 0x1.bd5cbbp-31
1.50000000000000011102230246251565404236316680908203125
Correct = 0x1.8p+0
gcc = 0x1.800000p+0
strtod = 0x1.800000p+0
9007199254740991.4999999999999999999999999999999995
Correct = 0x1.fffffffffffffp+52
gcc = 0x1.000000p+53
strtod = 0x2.000000p+52
After having substituted the %a by %f in the printf statement, the results look like this:
0.500000000000000166533453693773481063544750213623046875
Correct = 0x1.0000000000002p-1
gcc = 0.500000
strtod = 0.500000
3.518437208883201171875e13
Correct = 0x1.0000000000002p+45
gcc = 35184372088832.016000
strtod = 35184372088832.016000
62.5364939768271845828
Correct = 0x1.f44abd5aa7ca4p+5
gcc = 62.536494
strtod = 62.536494
8.10109172351e-10
Correct = 0x1.bd5cbaef0fd0cp-31
gcc = 0.000000
strtod = 0.000000
1.50000000000000011102230246251565404236316680908203125
Correct = 0x1.8p+0
gcc = 1.500000
strtod = 1.500000
9007199254740991.4999999999999999999999999999999995
Correct = 0x1.fffffffffffffp+52
gcc = 9007199254740992.000000
strtod = 9007199254740991.000000
—————-
write once, debug everywhere… I could offer results from a non-Intel-driven netbook running WinXPprofx32, but before installing gcc x32 on that device, I’d like to have the double vs. float issue resolved. I’m not a professional C programmer, though…
Regards
Georg
Georg,
My notes say that MinGW does not handle %a correctly. Perhaps this is the problem?
Also, perhaps you can try %.17g instead of %f, to get more digits. From that I can use my converter to see if that matches the correct result (assuming it prints correctly rounded).
Thanks for the output.
A similar test for floats (binary32), with numbers synthesized from the above examples, yields the wrong rounding (0x1p+24) for decimal literal
16777215.4999999999999999999999999999999995,
it must be 0x1.fffffep+23. It fails for GCC, C strtof, and C++ iostream.
Literals I use for testing are:
1) 0.5000000894069671630859375
2) 65536.01171875
3) 14.5382418633098799026
4) 7.803045154438761140010261119e-10
5) 1.500000059604644775390625
6) 16777215.4999999999999999999999999999999995
7) 3.16612277712032069713275043640247858901e+38
8) 5.88828039171069e+37
9) 3.0891034981081734072958e-32
All cases, except case 6, round correctly.
I forgot mention, it was with GCC 5.3.0.
@Waldemar,
How are you testing these? Do the literals have an ‘f’ suffix? If not, you are getting a double rounding error.
OK, I forgot the ‘f’ in the source code literals, and that solved my error issue with GCC.
However, the issues with ‘strtof’, and ‘iostream”s insertion ‘operator>>’ are not solved. And ‘strto(f)’ and strict typed ‘operator>>’ know the size of the destination floating-point value, so, it should perform the correct rounding.
My test code: http://hostcode.sourceforge.net/view/4950
Compile command: g++ -std=c++11 -o main *.cpp
@Waldemar,
My guess is that internally they are converting to double first, then to float. You may want to pursue a bug report. (Please keep us informed.)
It could be, at least for some cases.
Here is a code I’ve used to test rounding of float, double, and long double, with GCC and Clang, adapted from the previous test code:
http://hostcode.sourceforge.net/view/4961
I already made some test, whose results I’ll post in comments to follow.
I’ve had to write my own float to hex string, instead of using the ‘std::hexfloat’ manipulator, because:
1) GCC truncates ‘long double’ and rounds them as plain ‘double’, i.e., it outputs 18446744073709551615.49999999999999995l as 0x1p+64 instead of 0x1.fffffffffffffffep+63.
2) Clang rounds ‘long double’ correctly, but does not output a normalized hexadecimal floating-point string, it prints the above number as 0xf.fffffffffffffffp+60.
Theses are not exhaustive tests, but they evidence some erroneous cases.
g++ -std=c++11 -pedantic -pedantic-errors -Wall -Wextra -Werror -Ofast -UDEBUG -DNDEBUG=1 -lm main.cpp && ./a.out
GNU GCC/G++: 5.3.0
GNU libstdc++: 20151204
GNU new libc: 2.2
Testing ‘float’:
Run-time ‘strtof’ value: literal(16777215.4999999999999999999999999999999995)
– expected = 0x1.fffffep+23
– actual = 0x1.000000p+24
Run-time ‘iostream’ value: literal(16777215.4999999999999999999999999999999995)
– expected = 0x1.fffffep+23
– actual = 0x1.000000p+24
2 failures
Testing ‘double’: 0 failures
Testing ‘long double’: 0 failures
2 failures !!!
RUN FAILED (exit value 2, total time: 386ms)
clang++ -std=c++11 -pedantic -pedantic-errors -Wall -Wextra -Werror -Ofast -UDEBUG -DNDEBUG=1 -lm main.cpp && ./a.out
LLVM Clang: 3.7.1
GNU libstdc++: 20151204
GNU new libc: 2.2
Testing ‘float’:
Run-time ‘strtof’ value: literal(16777215.4999999999999999999999999999999995)
– expected = 0x1.fffffep+23
– actual = 0x1.000000p+24
Run-time ‘iostream’ value: literal(16777215.4999999999999999999999999999999995)
– expected = 0x1.fffffep+23
– actual = 0x1.000000p+24
2 failures
Testing ‘double’: 0 failures
Testing ‘long double’: 0 failures
2 failures !!!
RUN FAILED (exit value 2, total time: 993ms)
Complete output on: http://coliru.stacked-crooked.com/a/26954b871ef5d24a
g++ -std=c++11 -pedantic -pedantic-errors -Wall -Wextra -Werror -Ofast -UDEBUG -DNDEBUG=1 -lm main.cpp && ./a.out
GNU GCC/G++: 5.3.0
GNU libstdc++: 20151204
GNU glibc: 2.15
Testing ‘float’: 6 failures
Testing ‘double’: 8 failures
Testing ‘long double’: 4 failures
18 failures !!!
——————————————
Complete output on: http://coliru.stacked-crooked.com/a/d385d790854c2a96
clang++ -std=c++11 -pedantic -pedantic-errors -Wall -Wextra -Werror -Ofast -UDEBUG -DNDEBUG=1 -lm main.cpp && ./a.out
LLVM Clang: 3.7.0
GNU libstdc++: 20151204
GNU glibc: 2.15
Testing ‘float’: 6 failures
Testing ‘double’: 8 failures
Testing ‘long double’: 4 failures
18 failures !!!
——————————————
Complete output on:
http://coliru.stacked-crooked.com/a/3cca8dd481b79e53
clang++ -stdlib=libc++ -std=c++11 -pedantic -pedantic-errors -Wall -Wextra -Werror -Ofast -UDEBUG -DNDEBUG=1 -lm main.cpp && ./a.out
LLVM Clang: 3.7.0
LLVM libc++: 1101
GNU glibc: 2.15
Testing ‘float’: 4 failures
Testing ‘double’: 6 failures
Testing ‘long double’: 4 failures
14 failures !!!
And this is embarrassing, and confusing…
Complete output on: http://cpp.sh/9mku
g++ -std=c++11 -Wall -Wextra -Wpedantic -O3 main.cpp && ./a.out
GNU GCC/G++: 4.9.2
GNU libstdc++: 20141220
GNU glibc: 2.19
Testing ‘float’: 0 failures
Testing ‘double’: 0 failures
Testing ‘long double’: 0 failures
Exit code: 0 (normal program termination)
So, older version use to yield the correct rounding.
¿Or it’s an architecture issue?