Re: Strange FP behaviour

["Tim Prince" <tprince at computer dot org>]

Evidently, my_sqr() is returning effectively long double, and main() is
rounding the second instance to double, but not the first.  As I'm sure
you're aware, it is impractical to avoid this effect in general on an
extended precision architecture like x87, although it is surprising for
it to occur in such a simple case.

When I tried this case, running with newlib, where printf() doesn't
support long double, the displayed result has been cast to double, and a
known bug in the display of the difference is triggered, so there are
questions of library implementation as well.

Tim Prince
----- Original Message -----
From: <michael.hiegemann@chkra.mail.abb.com>
To: <gcc-bugs@gcc.gnu.org>
Sent: Tuesday, July 11, 2000 2:30 AM
Subject: Strange FP behaviour


>
>
> Hi folks,
>
> My system is SuSE 6.4 (i386) - Kernel 2.2.14, P90 on ASUS P/I-P55TP4XE
board,
> libc.so.6
>
> The following code (fp_bug.c) seems to be compiled incorrect:
>
> #include <stdio.h>
> #include <math.h>
>
>
> double my_sqr (double y)
> {
>      return y * y;
> }
>
>
> int main ()
> {
>      double y, result, result2, delta;
>
>
>      y       = 6.03856090847925908e-01;
>
>      result  = (0.5 - my_sqr (y)) + 0.5;
>
>      result2 = my_sqr (y);
>      result2 = (0.5 - result2) + 0.5;
>
>      delta   = fabs (result - result2);
>
>      printf ("result  = %24.16e\n", result);
>      printf ("result2 = %24.16e\n", result2);
>      printf ("delta   = %24.16e\n", delta);
>
>      return 1;
> }
>
>
> The program was compiled via
>
> gcc -v fp_bug.c
>
> and the compiler output is
>
>
> Reading specs from /usr/lib/gcc-lib/i486-suse-linux/2.95.2/specs
> gcc version 2.95.2 19991024 (release)
>  /usr/lib/gcc-lib/i486-suse-linux/2.95.2/cpp -lang-c -v -D__GNUC__=2
> -D__GNUC_MINOR__=95 -D__ELF__ -Dunix -D__i386__ -Dlinux -D__ELF__ -D__
unix__
> -D__i386__ -D__linux__ -D__unix -D__linux -Asystem(posix) -Wall -pedan
tic
> -Acpu(i386) -Amachine(i386) -Di386 -D__i386 -D__i386__ -Di486 -D__i486
> -D__i486__ fp_bug.c /tmp/cc2Auxw0.i
> GNU CPP version 2.95.2 19991024 (release) (i386 Linux/ELF)
> #include "..." search starts here:
> #include <...> search starts here:
>  /usr/local/include
>  /usr/lib/gcc-lib/i486-suse-linux/2.95.2/include
>  /usr/include
> End of search list.
> The following default directories have been omitted from the search
path:
>  /usr/include/g++
>
/usr/lib/gcc-lib/i486-suse-linux/2.95.2/../../../../i486-suse-linux/incl
ude
> End of omitted list.
>  /usr/lib/gcc-lib/i486-suse-linux/2.95.2/cc1
/tmp/cc2Auxw0.i -quiet -dumpbase
> fp_bug.c -Wall -pedantic -version -o /tmp/ccgR5is8.s
> GNU C version 2.95.2 19991024 (release) (i486-suse-linux) compiled by
GNU C
> version 2.95.2 19991024 (release).
>  /usr/i486-suse-linux/bin/as -V -Qy -o /tmp/ccxSoHuw.o /tmp/ccgR5is8.s
> GNU assembler version 2.9.5 (i486-suse-linux) using BFD version
2.9.5.0.24
>  /usr/lib/gcc-lib/i486-suse-linux/2.95.2/collect2 -m
elf_i386 -dynamic-linker
> /lib/ld-linux.so.2 /usr/lib/crt1.o /usr/lib/crti.o
> /usr/lib/gcc-lib/i486-suse-linux/2.95.2/crtbegin.o
> -L/usr/lib/gcc-lib/i486-suse-linux/2.95.2 -L/usr/i486-suse-linux/lib
> /tmp/ccxSoHuw.o -lgcc -lc -lgcc
/usr/lib/gcc-lib/i486-suse-linux/2.95.2/crtend.o
>  /usr/lib/crtn.o
>
>
> The program a.out generates the following output:
>
> result  =  6.3535782154586140e-01
> result2 =  6.3535782154586151e-01
> delta   =  1.1102230246251565e-16
>
>
> A friend tested it with his PII/400 system, ASUS P2B-B board, OS SuSE
linux 6.4,
>  Kernel 2.2.14.
> He got the same program output, the compiler info is virtually the
same and the
> assembler output is identical.
>
> We expect the delta to be exactly zero. Could you have a look at this?
>
>
> Greetings from Switzerland,
>
> Michael Hiegemann and Christoph Schilling
>
>
>
>
> For sake of completeness, the assembler output (gcc -S fp_bug.c) is
>
>
>      .file     "fp_bug.c"
>      .version  "01.01"
> gcc2_compiled.:
> .text
>      .align 16
> .globl my_sqr
>      .type      my_sqr,@function
> my_sqr:
>      pushl %ebp
>      movl %esp,%ebp
>      fldl 8(%ebp)
>      fmull 8(%ebp)
>      jmp .L2
>      .p2align 4,,7
> .L2:
>      movl %ebp,%esp
>      popl %ebp
>      ret
> .Lfe1:
>      .size      my_sqr,.Lfe1-my_sqr
> .section  .rodata
> .LC2:
>      .string   "result  = %24.16e\n"
> .LC3:
>      .string   "result2 = %24.16e\n"
> .LC4:
>      .string   "delta   = %24.16e\n"
>      .align 8
> .LC0:
>      .long 0x235d4f7,0x3fe352ca
>      .align 8
> .LC1:
>      .long 0x0,0x3fe00000
> .text
>      .align 16
> .globl main
>      .type      main,@function
> main:
>      pushl %ebp
>      movl %esp,%ebp
>      subl $40,%esp
>      fldl .LC0
>      fstpl -8(%ebp)
>      addl $-8,%esp
>      fldl -8(%ebp)
>      subl $8,%esp
>      fstpl (%esp)
>      call my_sqr
>      addl $16,%esp
>      fldl .LC1
>      fsubp %st,%st(1)
>      fldl .LC1
>      faddp %st,%st(1)
>      fstpl -16(%ebp)
>      addl $-8,%esp
>      fldl -8(%ebp)
>      subl $8,%esp
>      fstpl (%esp)
>      call my_sqr
>      addl $16,%esp
>      fstpl -24(%ebp)
>      fldl .LC1
>      fsubl -24(%ebp)
>      fldl .LC1
>      faddp %st,%st(1)
>      fstpl -24(%ebp)
>      fldl -16(%ebp)
>      fsubl -24(%ebp)
>      fabs
>      fstpl -32(%ebp)
>      addl $-4,%esp
>      fldl -16(%ebp)
>      subl $8,%esp
>      fstpl (%esp)
>      pushl $.LC2
>      call printf
>      addl $16,%esp
>      addl $-4,%esp
>      fldl -24(%ebp)
>      subl $8,%esp
>      fstpl (%esp)
>      pushl $.LC3
>      call printf
>      addl $16,%esp
>      addl $-4,%esp
>      fldl -32(%ebp)
>      subl $8,%esp
>      fstpl (%esp)
>      pushl $.LC4
>      call printf
>      addl $16,%esp
>      movl $1,%eax
>      jmp .L3
>      .p2align 4,,7
> .L3:
>      movl %ebp,%esp
>      popl %ebp
>      ret
> .Lfe2:
>      .size      main,.Lfe2-main
>      .ident    "GCC: (GNU) 2.95.2 19991024 (release)"
>
>