11

I have seen various answer here that depicts Strange behavior of pow function in C.
But I Have something different to ask here.

In the below code I have initialized int x = pow(10,2) and int y = pow(10,n) (int n = 2).

In first case it when I print the result it shows 100 and in the other case it comes out to be 99.

I know that pow returns double and it gets truncated on storing in int, but I want to ask why the output comes to be different.

CODE1

#include<stdio.h>
#include<math.h>
int main()
{
     int n = 2;
     int x;
     int y;
     x = pow(10,2);   //Printing Gives Output 100   
     y = pow(10,n);   //Printing Gives Output 99


     printf("%d %d" , x , y);

}

Output : 100 99

Why is the output coming out to be different. ?

My gcc version is 4.9.2

Update :

Code 2

int main()
    {
         int n = 2;
         int x;
         int y;
         x = pow(10,2);   //Printing Gives Output 100   
         y = pow(10,n);   //Printing Gives Output 99
         double k = pow(10,2);
         double l = pow(10,n);


         printf("%d %d\n" , x , y);
         printf("%f %f\n" , k , l);

    }

Output : 100 99
100.000000 100.000000

Update 2 Assembly Instructions FOR CODE1

Generated Assembly Instructions GCC 4.9.2 using gcc -S -masm=intel :

    .LC1:
    .ascii "%d %d\0"
    .text
    .globl  _main
    .def    _main;  .scl    2;  .type   32; .endef
_main:
    push    ebp
    mov ebp, esp
    and esp, -16
    sub esp, 48
    call    ___main
    mov DWORD PTR [esp+44], 2
    mov DWORD PTR [esp+40], 100      //Concerned Line
    fild    DWORD PTR [esp+44]
    fstp    QWORD PTR [esp+8]
    fld QWORD PTR LC0
    fstp    QWORD PTR [esp]
    call    _pow                    //Concerned Line
    fnstcw  WORD PTR [esp+30]
    movzx   eax, WORD PTR [esp+30]
    mov ah, 12
    mov WORD PTR [esp+28], ax
    fldcw   WORD PTR [esp+28]
    fistp   DWORD PTR [esp+36]
    fldcw   WORD PTR [esp+30]
    mov eax, DWORD PTR [esp+36]
    mov DWORD PTR [esp+8], eax
    mov eax, DWORD PTR [esp+40]
    mov DWORD PTR [esp+4], eax
    mov DWORD PTR [esp], OFFSET FLAT:LC1
    call    _printf
    leave
    ret
    .section .rdata,"dr"
    .align 8
LC0:
    .long   0
    .long   1076101120
    .ident  "GCC: (tdm-1) 4.9.2"
    .def    _pow;   .scl    2;  .type   32; .endef
    .def    _printf;    .scl    2;  .type   32; .endef
Suraj Jain
  • 4,463
  • 28
  • 39
CuriosGuy
  • 181
  • 1
  • 9
  • 1
    I have no problem https://ideone.com/sbaCit – Ed Heal Feb 10 '17 at 16:52
  • 2
    I get 100 both times. What is your build command? – Trevor Hickey Feb 10 '17 at 16:52
  • @EugeneSh. Please read it again , it has different answer n = 2 pow(10,n) and pow(10,2) , in first it is 99 and in second with 100. – CuriosGuy Feb 10 '17 at 16:53
  • Cannot reproduce (VS2013 on x64). – barak manos Feb 10 '17 at 16:53
  • Compiling it with Apple LLVM version 8.0.0 (clang-800.0.42.1), I obtain 100 and 100. With and without `#include`. – Ortomala Lokni Feb 10 '17 at 16:54
  • Can you print the resulting `double` values and see any differences? Just as an experiment. – Eugene Sh. Feb 10 '17 at 17:06
  • BTW, you can look for the source code of `pow` corresponding to your gcc version, duplicate it and run with the debugger to see what's going on. – Eugene Sh. Feb 10 '17 at 17:11
  • @CuriosGuy Just googled version of your compiler and it looks like it can be a bug (if you use optimisation) See my answer for details – Yuriy Ivaskevych Feb 10 '17 at 17:16
  • It's more likely glibc (or your C standard library) than gcc. gcc only generates the call to `pow` and links to it; the library itself is the likely cause of rounding errors. – trent Feb 10 '17 at 17:27
  • The exact same problem had already been found here: http://stackoverflow.com/questions/19126809/c-i-got-different-results-with-pow10-2-and-pow10-j-j-2?noredirect=1&lq=1, but with another compiler. Since the compiler is so important in this case, maybe this question shouldn't be considered a duplicate. – Fabio says Reinstate Monica Feb 10 '17 at 17:44
  • 1
    @FabioTurati Actually looks like exact duplicate, but I would rather close that one as a duplicate of this one, as it is higher quality... – Eugene Sh. Feb 10 '17 at 17:49
  • Simple solution: don't mix floating point functions with integers if you need exact results. – too honest for this site Feb 10 '17 at 17:50
  • @Olaf Did you read the question, i am asking why it is happening what happens behind the scene atleast read question first.. – CuriosGuy Feb 10 '17 at 19:47
  • 3
    This question has an answer here: http://stackoverflow.com/questions/32353301/will-intpown-m-be-wrong-for-some-positive-integers-n-m. The C standard says the precision of math functions like pow() is implementation defined. The glibc library fixed bugs due to rounding mode recently, pow(1.01,1.1) would crash (!) in one case. Discussed here: https://sourceware.org/bugzilla/show_bug.cgi?id=3976 Fixed here: https://github.com/lattera/glibc/commit/b7cd39e8f8c5cf2844f20eb03f545d19c4c25987. – amdn Feb 10 '17 at 21:59
  • typo... previous comment should read "accuracy" not "precision" – amdn Feb 10 '17 at 22:41
  • @amdn Sir , What does this mean std::pow(int, int) – CuriosGuy Feb 11 '17 at 03:58
  • 1
    The answer to your question about what does `__gnu_cxx::__promote_2::__value>::__type, __gnu_cxx::__promote::__value>::__type>::__type std::pow(int, int)` mean is roughly that the implementation calls a template function and converts (promotes) your arguments from their original type to the type required by the `std::pow(int, int)`. It isn't immediately obvious why the conversions are deemed necessary. – Jonathan Leffler Feb 11 '17 at 04:06
  • @JonathanLeffler That means there is pow function that takes int argument , how can it be ? What is meant by template function ? – CuriosGuy Feb 11 '17 at 04:08
  • 1
    Amongst other things, it means you're compiling your code with a C++ compiler, which makes the C tag on the question odd. How are you abusing the system so that you're compiling it as C++? You must either be using `g++` as the compiler, or possibly using `gcc` but specifying a file name with an extension like `.cpp` or `.cxx` — or perhaps `.C` (upper-case letter C) as the extension, and `gcc` is interpreting that as a C++ source file and running appropriate compiler phases for C++ instead of C. – Jonathan Leffler Feb 11 '17 at 04:12
  • I have given link, i have compiled using gcc 6.3 in godbolt.org – CuriosGuy Feb 11 '17 at 04:13
  • 1
    The IDE says it is C++ source, so your code is being compiled as C++. – Jonathan Leffler Feb 11 '17 at 04:15
  • @JonathanLeffler I noticed now , godbolt is doing that , wait – CuriosGuy Feb 11 '17 at 04:15
  • @JonathanLeffler I dont have gcc 6.3 , but my friends codeblock has it , maybe that is why his pow function is giving correct answer , why is pow giving wrong answer in pow(10,n) and pow(10,2) it is giving 100. How can i generate assembly online for C, ? – CuriosGuy Feb 11 '17 at 04:18
  • You should read the [documentation](https://gcc.gnu.org/onlinedocs/), but the answer is `gcc -S`. – Jonathan Leffler Feb 11 '17 at 04:19
  • 1
    @JonathanLeffler I used this answer http://stackoverflow.com/a/25781924/7546499 , and now godbolt is treating code as C https://godbolt.org/g/LCslhx – CuriosGuy Feb 11 '17 at 04:25
  • Re. "update 3", it means you used a C++ compiler, which is irrelevant and should be removed since this is a C-only question. You should also remove the C++ parts of "update 2".' – M.M Feb 11 '17 at 04:51
  • @M.M Done , Sir , can you tell me how to generate the complete floating point number of pow(10,n) just so that i could see where it fell short – CuriosGuy Feb 11 '17 at 04:54
  • Same Problem Here http://stackoverflow.com/questions/18155883/strange-behaviour-of-the-pow-function – Suraj Jain Feb 19 '17 at 07:47
  • Compilers can constant fold away calculations sometimes and this can lead to differences compared to when the calculation is done at run-time, usually using `-fno-builtin` will disable. Also see [log(10.0) can compile but log(0.0) cannot?](http://stackoverflow.com/q/24294578/1708801) and [Inconsistent strcmp() return value when passing strings as pointers or as literals](http://stackoverflow.com/q/27751221/1708801) – Shafik Yaghmour Feb 21 '17 at 05:07

5 Answers5

7

I know that pow returns double and it gets truncated on storing in int, but I want to ask why the output comes to be different.

You must first, if you haven't already, divest yourself of the idea that floating-point numbers are in any way sensible or predictable. double only approximates real numbers and almost anything you do with a double is likely to be an approximation to the actual result.

That said, as you have realized, pow(10, n) resulted in a value like 99.99999999999997, which is an approximation accurate to 15 significant figures. And then you told it to truncate to the largest integer less than that, so it threw away most of those.

(Aside: there is rarely a good reason to convert a double to an int. Usually you should either format it for display with something like sprintf("%.0f", x), which does rounding correctly, or use the floor function, which can handle floating-point numbers that may be out of the range of an int. If neither of those suit your purpose, like in currency or date calculations, possibly you should not be using floating point numbers at all.)

There are two weird things going on here. First, why is pow(10, n) inaccurate? 10, 2, and 100 are all precisely representable as double. The best answer I can offer is that the C standard library you are using has a bug. (The compiler and the standard library, which I assume are gcc and glibc, are developed on different release schedules and by different teams. If pow is returning inaccurate results, that is probably a bug in glibc, not gcc.)

In the comments on your question, amdn found a glibc bug to do with FP rounding that might be related and another Q&A that goes into more detail about why this happens and how it's not a violation of the C standard. chux's answer also addresses this. (C doesn't require implementation of IEEE 754, but even if it did, pow isn't required to use correct rounding.) I will still call this a glibc bug, because it's an undesirable property.

(It's also conceivable, though unlikely, that your processor's FPU is wrong.)

Second, why is pow(10, n) different from pow(10, 2)? This one is far easier. gcc optimizes away function calls for which the result can be calculated at compile time, so pow(10, 2) is almost certainly being optimized to 100.0. If you look at the generated assembly code, you will find only one call to pow.

The GCC manual, section 6.59 describes which standard library functions may be treated in this way (follow the link for the full list):

The remaining functions are provided for optimization purposes.

With the exception of built-ins that have library equivalents such as the standard C library functions discussed below, or that expand to library calls, GCC built-in functions are always expanded inline and thus do not have corresponding entry points and their address cannot be obtained. Attempting to use them in an expression other than a function call results in a compile-time error.

[...]

The ISO C90 functions abort, abs, acos, asin, atan2, atan, calloc, ceil, cosh, cos, exit, exp, fabs, floor, fmod, fprintf, fputs, frexp, fscanf, isalnum, isalpha, iscntrl, isdigit, isgraph, islower, isprint, ispunct, isspace, isupper, isxdigit, tolower, toupper, labs, ldexp, log10, log, malloc, memchr, memcmp, memcpy, memset, modf, pow, printf, putchar, puts, scanf, sinh, sin, snprintf, sprintf, sqrt, sscanf, strcat, strchr, strcmp, strcpy, strcspn, strlen, strncat, strncmp, strncpy, strpbrk, strrchr, strspn, strstr, tanh, tan, vfprintf, vprintf and vsprintf are all recognized as built-in functions unless -fno-builtin is specified (or -fno-builtin-function is specified for an individual function).

So it would seem you can disable this behavior with -fno-builtin-pow.

Community
  • 1
  • 1
trent
  • 25,033
  • 7
  • 51
  • 90
  • As you pointed in the comments, the math library is precompiled, so I doubt `pow(10,2)` can be computed in the compile-time, as the compiler is just not aware of the implementation of `pow`. – Eugene Sh. Feb 10 '17 at 17:29
  • @EugeneSh. Try it yourself. Compiling with `gcc -S`, only one call to `pow` is in the assembly output. My x86 isn't quite good enough to figure out what happens to the other one, but I assume it gets optimized to the equivalent of `100.0`. – trent Feb 10 '17 at 17:32
  • @trentcl It is completely legit theoretically. But I wonder how the compiler knows what it should be optimized to not knowing what `pow` is doing in practice... Unless the `pow` function is hardcoded in the compiler itself. – Eugene Sh. Feb 10 '17 at 17:34
  • @trentcl Aha! Now you have nailed it. I would give you +5 if I could :) – Eugene Sh. Feb 10 '17 at 17:50
  • Wait a second - does this mean that when you call `pow(10, 2)` the compiler understands that it can be optimized away and replaced with the result (so far so good), but then to calculate the result it uses its own implementation of `pow`, which is built into it, and might be different from the one in the standard library which is called at run time for `pow(10, n)`, as in this case? It sounds crazily dangerous! I hope the fact that the 2 libraries aren't the same in that version of gcc is recognized as a bug! – Fabio says Reinstate Monica Feb 10 '17 at 18:25
  • 1
    @FabioTurati Apparently, the danger of such an optimization is exactly demonstrated by the question :) – Eugene Sh. Feb 10 '17 at 19:12
  • *[...] because squaring 10 should not invoke floating-point rounding errors.* Why not? `pow(x, y)` is typically implemented as `exp(log(x) * y)` where `exp` is e^x and `log` is the natural logarithm. – nwellnhof Feb 10 '17 at 20:25
  • +1 @trentcl excellent answer. I've added a comment to the OP regarding a bug fix to glibc, feel free to include in your answer. – amdn Feb 10 '17 at 22:10
  • @nwellnhof A fair point. I meant "should" in a more general sense. The C standard technically allows 2.0 + 2.0 to equal 4.000001, but surely it shouldn't! However, even in IEEE 754, it's not required for `pow` to have correct rounding as for `+ - * / sqrt`, so you are right. (I'm still pretty sure this is a bug though.) – trent Feb 10 '17 at 22:27
  • @trentcl Can You See The Updated Question what does this mean std::pow(int, int) – CuriosGuy Feb 11 '17 at 04:00
  • @trentcl: The fact that the exponent 2.0 "is precisely representable as a double" is totally irrelevant. 0.5 is also precisely representable as a double, but I think we would all be surprised if `pow(2.0, 0.5)` returned the precisely correct value. :) While some standard libraries do check to see if the exponent is a small integer and attempt to use that information to minimize inaccuracy (even making the result precise), `pow` is under no obligation to do so. A good implementation will probably do better than `pow(a, b) => exp(a, log(b))`, but that would be legal, too.... – rici Feb 11 '17 at 04:19
  • ... You are 100% correct that the compiler is allowed to compute standard library functions at compile time if it can deduce the arguments, and the standard explicitly allows compile-time computations to be different from run-time computations, even for ordinary floating-point arithmetic. The difference should be minor, maybe only an ULP, but it can happen. The real problem with this code is that the floor operation (and implicit double->integer conversion) introduces numerical instability, and should be avoided... – rici Feb 11 '17 at 04:23
  • ... floating-point arithmetic is underspecified in the C standard, but IEEE-749 does a much better job and the requirements it provides are reasonable and sufficient for programmers to analyze code and deduce consequences. If you start with the hypothesis that floating point arithmetic is inherently unpredictable, you will be much less likely to learn how to predict it. :) – rici Feb 11 '17 at 04:25
  • Floating point numbers are predictable. If they weren't then they would be completely useless. – M.M Feb 11 '17 at 04:49
  • @trentcl Does your answer means the same http://stackoverflow.com/a/19127020/7546499 , it is talking about constant folding ? – CuriosGuy Feb 11 '17 at 05:03
  • @CuriosGuy: "constant folding" is another term for compile-time computation, so that is the same thing. The standard library is part of the standard, and users are not allowed to define functions with the same name as standard library functions if they have #include'd the standard headers. So the compiler is allowed to compute standard math functions at compile time. – rici Feb 11 '17 at 05:08
  • @rici So the compiler has another implementation of pow, the answer we are commenting on is talking about built in pow function. – CuriosGuy Feb 11 '17 at 05:10
  • What does this line mean " GCC built-in functions are always expanded inline and thus do not have corresponding entry points and their address cannot be obtained. Attempting to use them in an expression other than a function call results in a compile-time error." – CuriosGuy Feb 11 '17 at 05:14
  • @CuriosGuy You should ask that as a separate question! – trent Feb 11 '17 at 14:37
  • @rici From a C perspective, floating-point math *is* inherently unpredictable, since it may do different things on different processors (theoretically, even without recompilation). For a given stackup, assuming IEEE 754 and only using `+ - * / sqrt`, sure, that's predictable. – trent Feb 11 '17 at 14:43
  • @trentcl You have this line in your answer , and i did not understand that "GCC built-in functions are always expanded inline and thus do not have corresponding entry points and their address cannot be obtained. Attempting to use them in an expression other than a function call results in a compile-time error." – CuriosGuy Feb 11 '17 at 15:10
  • @curiosguy: i thought I answered that but my comment either vanished or never was registered. Anyway, it means what it says: you cannot take the address of an inline function because it doesn't exist as an entity. This has nothing to do with your question and tentcl's quote of that line is also irrelevant. – rici Feb 11 '17 at 17:22
2

Why is the output coming out to be different. ? (in the updated appended code)

We do not know the values are that different.

When comparing the textual out of int/double, be sure to print the double with sufficient precision to see if it is 100.000000 or just near 100.000000 or in hex to remove all doubt.

printf("%d %d\n" , x , y);
// printf("%f %f\n" , k , l);
// Is it the FP number just less than 100?
printf("%.17e %.17e\n" , k , l);  // maybe 9.99999999999999858e+01
printf("%a %a\n" , k , l);        // maybe 0x1.8ffffffffffff0000p+6

Why is the output coming out to be different. ? (in the original code)

C does not specify the accuracy of most <math.h> functions. The following are all compliant results.

// Higher quality functions return 100.0
pow(10,2) --> 100.0   
// Lower quality and/or faster one may return nearby results
pow(10,2) --> 100.0000000000000142...
pow(10,2) --> 99.9999999999999857...

Assigning a floating point (FP) number to an int simple drops the fraction regardless of how close the fraction is to 1.0

When converting FP to an integer, better to control the conversion and round to cope with minor computational differences.

// long int lround(double x);
long i = lround(pow(10.0,2.0));
Community
  • 1
  • 1
chux - Reinstate Monica
  • 143,097
  • 13
  • 135
  • 256
0

You're not the first to find this. Here's a discussion form 2013: pow() cast to integer, unexpected result

I'm speculating that the assembly code produced by the tcc guys is causing the second value to be rounded down after calculating a result that is REALLY close to 100. Like mikijov said in that historic post, looks like the bug has been fixed.

Community
  • 1
  • 1
Kellen
  • 1
  • 1
  • So, it's something with your compiler or your standard library. I could tr to recreate the problem on my university's cse machines... would that help? – Kellen Feb 10 '17 at 17:06
  • He says it's GCC 4.9, which is close to 3 years old now. – AntonH Feb 10 '17 at 17:09
  • I could not reproduce the error, I did not modify your code. I am using gcc version 4.9.4 you gotta update. the programmers who make gcc made a bit of a mistake. – Kellen Feb 10 '17 at 17:12
  • Please see the question again. – CuriosGuy Feb 11 '17 at 05:07
0

As others have mentioned, Code 2 returns 99 due to floating point truncation. The reason why Code 1 returns a different and correct answer is because of a libc optimization.

When the power is a small positive integer, it is more efficient to perform the operation as repeated multiplication. The simpler path removes roundoff. Since this is inlined you don't see function calls being made.

doron
  • 27,972
  • 12
  • 65
  • 103
-1

You've fooled it into thinking that the inputs are real and so it gives an approximate answer, which happens to be slightly under 100, e.g. 99.999999 that is then truncated to 99.

Malcolm McLean
  • 6,258
  • 1
  • 17
  • 18
  • 1
    Please read question again , I know about floating inaccuracy thing it is different question . – CuriosGuy Feb 10 '17 at 17:03
  • If you give pow an integer constant, it obviously special-cases it and inlines or somehow calculates with exact arithmetic. If you pass an integer variable, it will be promoted to a double, and must go down the approximation route. – Malcolm McLean Feb 10 '17 at 19:10
  • still your answer is short of many things.. update it to contain all you know about topic including this comment it tells a lot.. and was this question bad... 3 downvots – CuriosGuy Feb 10 '17 at 19:20