Skip site navigation (1)Skip section navigation (2)
Date:      Sat, 14 Dec 1996 22:25:06 +1100
From:      Bruce Evans <bde@zeta.org.au>
To:        current@FreeBSD.org, dmaddox@scsn.net
Subject:   Re: Something broken in libm
Message-ID:  <199612141125.WAA02738@godzilla.zeta.org.au>

next in thread | raw e-mail | index | archive | help
>    There is definitely something broken in libm when HAVE_FPU is
>defined...  More specifically, the problem seems to be in the
>__ieee754_pow() function defined in src/lib/msun/src/e_pow.c.  I
>have not been able to nail down the exact problem, but if anyone
>else is interested --

Try this fix.  __ieee754_pow() only uses two functions that change
when HAVE_FPU is defined, __ieee754_sqrt() and scalbn().
__iee754_sqrt() is to simple to be be broken, but scalbn() left
garbage on the FPU stack.

---
diff -c2 s_scalbn.S~ s_scalbn.S
*** s_scalbn.S~	Fri Aug 19 21:14:30 1994
--- s_scalbn.S	Sat Dec 14 21:36:34 1996
***************
*** 42,44 ****
--- 42,45 ----
  	fldl	4(%esp)
  	fscale
+ 	fstp	%st(1)
  	ret
---

It's surprising that this problem hasn't been noticed before.
scalbn() is also used for log10() and remainders mod pi/2.

Some packages ignore all floating point errors (see fpsetmask(3)).
FPU stack overflow should not be ignored, but unfortunately there
is no way to trap it without also trapping for invalid operands.
When the stack overflows, pushing further operands onto it pushes
NaNs instead of the operands, so nothing much can be computed.
E.g., the following program termintes after 1 iteration instead
of 1e6 iterations:

main()
{
    double d;

    /* Fill up FPU stack. */
    asm("fldz");
    asm("fldz");
    asm("fldz");
    asm("fldz");
    asm("fldz");
    asm("fldz");
    asm("fldz");
    asm("fldz");

    for (d = 0; d < 1e6; )
	++d;
}

Bruce



Want to link to this message? Use this URL: <https://mail-archive.FreeBSD.org/cgi/mid.cgi?199612141125.WAA02738>