double to int conversion yields strange results

  • Thread starter =?ISO-8859-1?Q?Bj=F8rn_Augestad?=
  • Start date
?

=?ISO-8859-1?Q?Bj=F8rn_Augestad?=

Below is a program which converts a double to an integer in two
different ways, giving me two different values for the int. The basic
expression is 1.0 / (1.0 * 365.0) which should be 365, but one variable
becomes 364 and the other one becomes 365.

Does anyone have any insight to what the problem is?

Thanks in advance.
Bjørn

$ cat d.c
#include <stdio.h>

int main(void)
{
double dd, d = 1.0 / 365.0;
int n, nn;

n = 1.0 / d;
dd = 1.0 / d;
nn = dd;

printf("n==%d nn==%d dd==%f\n", n, nn, dd);
return 0;
}

$ gcc -Wall -O0 -ansi -pedantic -W -Werror -o d d.c
$ ./d
n==364 nn==365 dd==365.000000

$ gcc -v
Reading specs from /usr/lib/gcc/i386-redhat-linux/3.4.2/specs
Configured with: ../configure --prefix=/usr --mandir=/usr/share/man
--infodir=/usr/share/info --enable-shared --enable-threads=posix
--disable-checking --with-system-zlib --enable-__cxa_atexit
--disable-libunwind-exceptions --enable-java-awt=gtk
--host=i386-redhat-linux
Thread model: posix
gcc version 3.4.2 20041017 (Red Hat 3.4.2-6.fc3)
$
 
R

Richard Bos

=?ISO-8859-1?Q?Bj=F8rn_Augestad?= said:
Below is a program which converts a double to an integer in two
different ways, giving me two different values for the int. The basic
expression is 1.0 / (1.0 * 365.0) which should be 365, but one variable
becomes 364 and the other one becomes 365.

How about people reading the FAQ for a change?

<http://www.eskimo.com/~scs/C-faq/q14.1.html>
<http://www.eskimo.com/~scs/C-faq/q14.4.html>
<http://www.eskimo.com/~scs/C-faq/q14.5.html>

Richard
 
M

Michael Mair

Bjørn Augestad said:
Below is a program which converts a double to an integer in two
different ways, giving me two different values for the int. The basic
expression is 1.0 / (1.0 * 365.0) which should be 365, but one variable
becomes 364 and the other one becomes 365.

Does anyone have any insight to what the problem is?

Thanks in advance.
Bjørn

$ cat d.c
#include <stdio.h>

int main(void)
{
double dd, d = 1.0 / 365.0;
int n, nn;

n = 1.0 / d;
dd = 1.0 / d;
nn = dd;

printf("n==%d nn==%d dd==%f\n", n, nn, dd);
return 0;
}

$ gcc -Wall -O0 -ansi -pedantic -W -Werror -o d d.c
$ ./d
n==364 nn==365 dd==365.000000

$ gcc -v
Reading specs from /usr/lib/gcc/i386-redhat-linux/3.4.2/specs
Configured with: ../configure --prefix=/usr --mandir=/usr/share/man
--infodir=/usr/share/info --enable-shared --enable-threads=posix
--disable-checking --with-system-zlib --enable-__cxa_atexit
--disable-libunwind-exceptions --enable-java-awt=gtk
--host=i386-redhat-linux
Thread model: posix
gcc version 3.4.2 20041017 (Red Hat 3.4.2-6.fc3)
$

gcc on x86 architectures may keep certain double variables in
registers where they are effectively long doubles and some may
be stored in memory where they really are doubles (talking about
64 bits vs 80 bits).
Try the option
-ffloat-store
and, if that does not "help", also
-mfpmath=sse -msse2
or whatever is appropriate for you.
This may be inefficient, so you could also use the macro
#define EVAL_AS_DOUBLE(d) (*((volatile double *)&(d)))
on your double objects, whenever you need to be sure that double
really is double.


Cheers
Michael
 
?

=?ISO-8859-1?Q?Bj=F8rn_Augestad?=

Michael said:
Bjørn Augestad said:
Below is a program which converts a double to an integer in two
different ways, giving me two different values for the int. The basic
expression is 1.0 / (1.0 * 365.0) which should be 365, but one
variable becomes 364 and the other one becomes 365.

Does anyone have any insight to what the problem is?
[snip]


gcc on x86 architectures may keep certain double variables in
registers where they are effectively long doubles and some may
be stored in memory where they really are doubles (talking about
64 bits vs 80 bits).
Try the option
-ffloat-store
and, if that does not "help", also
-mfpmath=sse -msse2
or whatever is appropriate for you.
This may be inefficient, so you could also use the macro
#define EVAL_AS_DOUBLE(d) (*((volatile double *)&(d)))
on your double objects, whenever you need to be sure that double
really is double.

Thanks, Michael.
-msse2 did the trick,-ffloat-store and -mfpmath didn't work.

What an odd bug that was. I had been using -march=prescott and
everything was working fine. Then I removed that option for portability
reasons and then the code started to break in strange ways.

I'll rewrite the code to avoid the problems completely.
Bjørn
 
R

Richard Bos

=?ISO-8859-1?Q?Bj=F8rn_Augestad?= said:
Thanks, but why do you think this applies to my question?

_Never_ assume that floating point computations are exact. You are
obviously computing 364.999999999 in one case, and 365.000000000001 in
the other. If you'd read the FAQ, you'd have realised this.

Richard
 
M

Michael Mair

Richard said:
_Never_ assume that floating point computations are exact. You are
obviously computing 364.999999999 in one case, and 365.000000000001 in
the other. If you'd read the FAQ, you'd have realised this.

However, this is a gcc bug due to the incredibly *** x86 architecture.
Bjorn got the semantics right, so even though pointing people to the
FAQ is hardly ever wrong, his problem was of a different nature and
could not be solved within the program because gcc left the realms of
standard C.


Cheers
Michael
 
M

Martin Ambuhl

Bjørn Augestad said:
Below is a program which converts a double to an integer in two
different ways, giving me two different values for the int. The basic
expression is 1.0 / (1.0 * 365.0) which should be 365, but one variable
becomes 364 and the other one becomes 365.

Does anyone have any insight to what the problem is?

That you, like many clueless people, have not learned how floating point
works, have not bothered to check the FAQ before posting, have not
checked the archives for countless threads on the same subject -- some
in the last week -- started by other clueless wonders, or even followed
the newsgroup for a week before posting. Will September *never* end?
 
C

CBFalconer

Bjørn Augestad said:
Below is a program which converts a double to an integer in two
different ways, giving me two different values for the int. The
basic expression is 1.0 / (1.0 * 365.0) which should be 365, but
one variable becomes 364 and the other one becomes 365.

Does anyone have any insight to what the problem is?

$ cat d.c
#include <stdio.h>

int main(void)
{
double dd, d = 1.0 / 365.0;
int n, nn;

n = 1.0 / d;
dd = 1.0 / d;
nn = dd;

printf("n==%d nn==%d dd==%f\n", n, nn, dd);
return 0;
}

After all the huffing and puffing about the FAQ, nobody has told
him what the fundamental cause is. The x86 does floating point
calculations with an 80 bit operand, which allows for extra
precision. Those operands are in the FP processor internal stack,
and stay there until needed. When the computation 1.0 / d is
directly converted to an integer, the extra bits are present. When
it is first stored in a 64 bit double, they have to be removed, and
thus the final conversion to an integer is different.

There are various options to control this, but fundamentally the
FAQ is correct - floating point values are imprecise. Optimization
level will probably also affect the result.
 
L

Lawrence Kirby

On Fri, 11 Feb 2005 13:58:12 +0100, Michael Mair wrote:

....
However, this is a gcc bug due to the incredibly *** x86 architecture.
Bjorn got the semantics right, so even though pointing people to the
FAQ is hardly ever wrong, his problem was of a different nature and
could not be solved within the program because gcc left the realms of
standard C.

How so? While gcc does have conformance issues in the area of floating
point I don't see any problem with the output of this program. In
particular C doesn't require n and nn to be set to the same value.
In fact the difference between them, where an intermediate value is
stored in a double object in one case and not in the other, is a classic
example of where the values are allowed to differ.

Lawrence
 
J

jacob navia

This happens because gcc sets the default precision to full
precision in the x86, i.e. 80 bits precision.

Using the lcc-win32 compiler with your program I obtain the
same result as gcc.

If I add however a call to a function to set the precision
of all operations to 64 bits only, I get the correct
result:
#include <stdio.h>
#include <fenv.h>
int main(void)
{
double dd, d = 1.0 / 365.0;
int n, nn;

DoublePrecision(); // Set the machine to 64 bits
// precision *only*
n = (double)(1.0 / d);
dd = 1.0 / d;
nn = dd;

printf("n==%d nn==%d dd==%f\n", n, nn, dd);
return 0;
}
Result:

n==365 nn==365 dd==365.000000

If you do not use the lcc-win32 compiler, gcc will
swallow this assembly:
.text
.globl _DoublePrecision
_DoublePrecision:
push %eax
fnstcw (%esp)
movl (%esp),%eax
btr $8,(%esp)
orw $0x200,(%esp)
fldcw (%esp)
popl %ecx
movb %ah,%al
andl $3,%eax
ret

You compile it with gcc and you set the machine precision
at the beginning of the prograM.

Maybe gcc offers a functional equivalent, I do not know.
 
M

Michael Mair

Lawrence said:
On Fri, 11 Feb 2005 13:58:12 +0100, Michael Mair wrote:

...



How so? While gcc does have conformance issues in the area of floating
point I don't see any problem with the output of this program. In
particular C doesn't require n and nn to be set to the same value.
In fact the difference between them, where an intermediate value is
stored in a double object in one case and not in the other, is a classic
example of where the values are allowed to differ.

Sorry, I shot off the other mail thinking of the following situation:

Q: If the code had been
double dd, d = 1.0 / 365.0;
int n, nn;

n = (double) (1.0 / d); /*difference*/
dd = 1.0 / d;
nn = dd;
I would expect n==nn. Is there a loop hole in the standard
s.th. this could be not the case?

TIA
Michael
 
M

Mark McIntyre

Even though it says, a==b is wrong, some _simple_ like

if(i == 1.0) {

seems to work. Is it due to the nature that 1.0 can be encoded "good enough"
in FP regs?

Possibly, tho there's no guarantee. It might work today, and fail
tomorrow., or work on one platform but not on another. Or always work
perfectly on every platform you're able to test on *except* the one your
new customer is going to use but hasn't told you about yet....
 
A

Andrey Tarasevich

Michael said:
However, this is a gcc bug due to the incredibly *** x86 architecture.
...

This has absolutely nothing to do with x86 architecture. These issues
are caused by certain properties of the floating-point representation
(natural to virtually any architecture) and by some inherent properties
of binary notation system itself. For example, floating-point binary
notation of a decimal value '0.1' has fractional part of infinite
length, which means that this value cannot be represented precisely by
any architecture using traditional floating-point representation.

In OP's case, 1 is divided by 365. The expected resultant value is
1/365, which is

0.(000000001011001110001100111110011011)

in binary notation. It is a periodic fraction of infinite length. It
cannot be ever represented precisely in floating-point format of finite
length on any architecture. Ever. Absolutely impossible. And it has
nothing to do with x86 and/or GCC.

There are alternative ways to represent fractional numbers, which can at
least partially solve these issues, by they have their own drawbacks.
 
C

Chris Torek

This happens because gcc sets the default precision to full
precision in the x86, i.e. 80 bits precision.

Actually, it does not: it leaves the precision bits alone. Someone
or something else is responsible for setting them initially.
If I add however a call to a function to set the precision
of all operations to 64 bits only, I get the correct result ...

"Correct" according to some other standard than Standard C,
anyway. (Standard C is pretty lax about FP results.)

Unfortunately, setting the precision control word on the x86 is
not sufficient in general. It will solve this particular problem,
but not other problems on the x86.
If you do not use the lcc-win32 compiler, gcc will
swallow this assembly:

(I assume you mean gas -- gcc itself is a C compiler, not an
assembler. :) )

[snippage]
Maybe gcc offers a functional equivalent, I do not know.

GCC has a way to do the trick with inline assembly, but this is
off-topic, and does not solve the problem anyway, because the x86
FPU is obnoxious. (Exercise for the x86 student: what 64-bit result
do you expect for 1.0e200 * 1.0e200? What 80-bit result do you
expect for 1.0e200L * 1.0e200L? What 80-bit result will you get
in the FPU stack for the first expression, no matter how you fiddle
with the control word? Does it matter if overflows do not occur
when they should?)
 
M

Michael Mair

Andrey said:
This has absolutely nothing to do with x86 architecture. These issues
are caused by certain properties of the floating-point representation
(natural to virtually any architecture) and by some inherent properties
of binary notation system itself. For example, floating-point binary
notation of a decimal value '0.1' has fractional part of infinite
length, which means that this value cannot be represented precisely by
any architecture using traditional floating-point representation.

In OP's case, 1 is divided by 365. The expected resultant value is
1/365, which is

0.(000000001011001110001100111110011011)

in binary notation. It is a periodic fraction of infinite length. It
cannot be ever represented precisely in floating-point format of finite
length on any architecture. Ever. Absolutely impossible. And it has
nothing to do with x86 and/or GCC.

There are alternative ways to represent fractional numbers, which can at
least partially solve these issues, by they have their own drawbacks.

Please have a look at my reply to Lawrence Kirby's reply.
I was thinking about the wrong situation and thought the OP
had used a cast -- in which case I would have expected the same
result but would not have gotten it with gcc.

-Michael
 
J

jacob navia

Chris said:
GCC has a way to do the trick with inline assembly, but this is
off-topic, and does not solve the problem anyway, because the x86
FPU is obnoxious. (Exercise for the x86 student: what 64-bit result
do you expect for 1.0e200 * 1.0e200? What 80-bit result do you
expect for 1.0e200L * 1.0e200L? What 80-bit result will you get
in the FPU stack for the first expression, no matter how you fiddle
with the control word? Does it matter if overflows do not occur
when they should?)

Overflow will be detected only when WRITING the numbers...

This precision setting will be used when writing/reading numbers
into the FPU, not when doing operations, as it should...

That is why the precision setting doesn't provoke an
overflow when multiplying 1e200*1e200: the number fits
in 80 bits with a dynamic range until around 1.18e4932,
only when you attempt to write the result into RAM
does the CPU notice... Ohh shh... I was in overflow!

jacob
 
M

Michael Mair

Christian said:
Michael Mair wrote:



How does this differ from "n = 1.0 / d"?

1.0/d may have been computed with higher precision and the result
may be there in higher precision, too. By telling the compiler
explicitly that I want the result converted to double (before
it is converted to int), I would expect the excess precision
to disappear.
In fact, if I try to calculate ***_EPSILON with gcc
and test
fltepsilon = 1.0;
while(1.0F < 1.0+fltepsilon)
fltepsilon /= 2.0;
fltepsilon *= 2.0;
I find fltepsilon==LDBL_EPSILON but if I insert a cast
while(1.0F < (float)(1.0+fltepsilon))
I arrive at the right result. The same does not hold for double
and DBL_EPSILON unless I force gcc to do so.

In my understanding, inserting the cast should lead to the same
result as storing into an intermediate variable of the type we cast
to and using this variable for the next operation.

Maybe I understand something wrong but with the cast I would
expect the original example to work as expected by the OP on a
conforming implementation.


Cheers
Michael
 

Ask a Question

Want to reply to this thread or ask your own question?

You'll need to choose a username for the site, which only take a couple of moments. After that, you can post your question and our members will help you out.

Ask a Question

Members online

Forum statistics

Threads
473,744
Messages
2,569,484
Members
44,904
Latest member
HealthyVisionsCBDPrice

Latest Threads

Top