4

Consider this code

x = 0.109354682484;
IntegerPart[x/(x/2)]
(* 1 *)

Precision[x]
(* MachinePrecision *)

Why does it give 1 ?

Version number: 9.0 on Mac 10.9.2

screenshot

enter image description here

Update:

If we use an undefined variable, IntegerPart[x0/(x0/2)] gives 2. Since Mathematica never gives warnings about this x0, I'm assuming for any x0 it is true.

If we calculate the same integer part using fortran, we get 2 instead of 1.

program main
implicit none

real(8):: x=0.109354682484
real(8):: y=1.4

write(*,*) int(x/(x/2))
write(*,*) nint(x/(x/2))
write(*,*) int(y)

end program main

compiled with ifort -O0 main.f90

output of above fortran code is

       2
       2
       1

according to here, int is a fortran intrinsic function that calculate the integer part.

Is this a bug?

xslittlegrass
  • 27,549
  • 9
  • 97
  • 186

1 Answers1

7

The explanation is interesting here. I tried the same in C++, and worked a bit extra to make sure the compiler won't optimize away the divisions (looking at the assembly output, it may optimize it away if you're not careful). Indeed, I get 2 with C++.

And here's why:

C++ does the equivalent of 

IntegerPart@Divide[x, Divide[x, 2]]

(* ==> 2 *)

while in Mathematica you're computing the equivalent of

x*(1/((1/2)*x))

which gives a different result.

If I do the calculation as (x* (1/ (0.5*x))) in C++, then the result is less than 2 and the integer part is 1.

Relevant reading:

Community
  • 1
Szabolcs
  • 234,956
  • 30
  • 623
  • 1,263
  • Thanks! I never knew that a/b is interpreted as a*(1/b). That makes me worried a little bit about all the codes I have written. Actually I caught on this because of a strange behavior of my package, and it took me a whole morning to track down to problem to IntegerPart I used in that package. – xslittlegrass Apr 01 '14 at 18:54
  • 1
    The observed behavior is wrong under IEEE-754 as well as generic binary floating-point implementations. The result of x/2 or 0.5 * 2 has the same mantissa as x, and the exponent is decreased by one. Dividing two numbers having the same mantissa then results in a mantissa of 1 (stored implicitly), and the exponent is the difference in operand exponents. The only way for a binary implementation to not get exactly 2 for the result of x/(x/2) or x/(0.5*x) is if the scaled value becomes denormal. – Ben Voigt Apr 01 '14 at 23:13
  • @BenVoigt But as I said it does not compute x/(0.5*x). It computes x * (1/(0.5*x)). I don't know what's correct by the standard, but the result is precisely the same as what I get from C++, which just uses a few simple floating point instructions (based on the asm output). Mathematica also uses hardware floating point arithmetic by default. – Szabolcs Apr 02 '14 at 04:55
  • @BenVoigt The reason why it does the arithmetic in this order is explained in the post I linked to. – Szabolcs Apr 02 '14 at 04:57
  • Oh right, I thought about only one reciprocal and there are two. – Ben Voigt Apr 02 '14 at 06:39