mailing list of musl libc
 help / color / mirror / Atom feed
* [musl] Re: Potentially incorrect musl scalbn results on AArch64?
       [not found] <CAEeofcgZbgsP6D-idQXFv6qqQw-d0kVeBHiCi2xqEeSoD-JcZw@mail.gmail.com>
@ 2021-02-26 11:43 ` Szabolcs Nagy
  2021-02-26 12:45   ` Alexander Richardson
  0 siblings, 1 reply; 2+ messages in thread
From: Szabolcs Nagy @ 2021-02-26 11:43 UTC (permalink / raw)
  To: Alexander Richardson; +Cc: musl

* Alexander Richardson <Alexander.Richardson@cl.cam.ac.uk> [2021-02-25 10:27:11 +0000]:
> Hello,
> 
> I've recently been tracking down testsuite failures on FreeBSD aarch64
> and as part of this updated the FreeBSD scalbn* implementations to use
> the musl versions. However, two of the scalbn tests are failing on
> non-x86 architectures (https://godbolt.org/z/rax7f6)
> For example, scalbn(1, -1023) returns
> "1.1125369292536006915451e-308"/0x0.8p-1022 on x86, but if I run the
> tests on aarch64 I get 0 instead.

i added musl list on cc

i cannot reproduce your issue (i.e. the c code works for me on
all targets as is)

one issue can be that if freebsd incorrectly sets the fpu on
aarch64 into flush-subnormals-to-zero mode.

or a clang compiler bug (which we have seen before wrt floating
point optimizations, although not wrong results, only wrong fenv)

> I'm not particularly familiar with floating-point calculations, but it
> appears to me that this could be caused by x86's extended precision
> during calculations?
> If I cast the result to (long double) on aarch64 prior to the
> multiplication, I get the expected result on AArch64 (but that's
> obviously slow and won't work on architectures where long double ==
> double).
> I've attached the current workaround, but I'm sure there is a better
> solution to this. Or possibly the test is incorrect and 0 is a
> perfectly valid result?
> 
> Kind regards,
> Alex
> 
> 
> diff --git a/lib/msun/src/s_scalbn.c b/lib/msun/src/s_scalbn.c
> index 219cd8f0c989..0d344840862f 100644
> --- a/lib/msun/src/s_scalbn.c
> +++ b/lib/msun/src/s_scalbn.c
> @@ -29,6 +29,19 @@ double scalbn(double x, int n)
>         }
>         u.i = (uint64_t)(0x3ff+n)<<52;
>         x = y * u.f;
> +#if !defined(__amd64__) && !defined(__i386__)
> +       /*
> +        * x86 performs the multiplication with higher precision, but on
> +        * non-x86 architectures we might get 0 instead of a tiny value. To work
> +        * around this problem perform the multiplication with float128 (slow).
> +        * TODO: This doesn't work on e.g. MIPS where long double == double.
> +        */
> +       if (x == 0.) {
> +               x = (long double)y * u.f;
> +               /* fprintf(stderr, "\ttrying again: %a/%a\n", x,
> (double)((long double)y * u.f)); */
> +               return x;
> +       }
> +#endif
>         return x;
>  }

^ permalink raw reply	[flat|nested] 2+ messages in thread

* [musl] Re: Potentially incorrect musl scalbn results on AArch64?
  2021-02-26 11:43 ` [musl] Re: Potentially incorrect musl scalbn results on AArch64? Szabolcs Nagy
@ 2021-02-26 12:45   ` Alexander Richardson
  0 siblings, 0 replies; 2+ messages in thread
From: Alexander Richardson @ 2021-02-26 12:45 UTC (permalink / raw)
  To: Szabolcs Nagy; +Cc: musl

On Fri, 26 Feb 2021 at 11:43, Szabolcs Nagy <nsz@port70.net> wrote:
>
> * Alexander Richardson <Alexander.Richardson@cl.cam.ac.uk> [2021-02-25 10:27:11 +0000]:
> > Hello,
> >
> > I've recently been tracking down testsuite failures on FreeBSD aarch64
> > and as part of this updated the FreeBSD scalbn* implementations to use
> > the musl versions. However, two of the scalbn tests are failing on
> > non-x86 architectures (https://godbolt.org/z/rax7f6)
> > For example, scalbn(1, -1023) returns
> > "1.1125369292536006915451e-308"/0x0.8p-1022 on x86, but if I run the
> > tests on aarch64 I get 0 instead.
>
> i added musl list on cc
>
> i cannot reproduce your issue (i.e. the c code works for me on
> all targets as is)
>
> one issue can be that if freebsd incorrectly sets the fpu on
> aarch64 into flush-subnormals-to-zero mode.
>
> or a clang compiler bug (which we have seen before wrt floating
> point optimizations, although not wrong results, only wrong fenv)
>

Thanks very much for that suggestion! Turns out that as of
https://cgit.freebsd.org/src/commit/?id=65618fdda0f272a823e6701966421bdca0efa301
FreeBSD sets the flush-subnormals-to-zero flag on startup so this is a
FreeBSD issue and I can confirm that the code works as expected when I
clear the flag.
I've submitted a possible fix to FreeBSD in https://reviews.freebsd.org/D28938.

Alex

> > I'm not particularly familiar with floating-point calculations, but it
> > appears to me that this could be caused by x86's extended precision
> > during calculations?
> > If I cast the result to (long double) on aarch64 prior to the
> > multiplication, I get the expected result on AArch64 (but that's
> > obviously slow and won't work on architectures where long double ==
> > double).
> > I've attached the current workaround, but I'm sure there is a better
> > solution to this. Or possibly the test is incorrect and 0 is a
> > perfectly valid result?
> >
> > Kind regards,
> > Alex
> >
> >
> > diff --git a/lib/msun/src/s_scalbn.c b/lib/msun/src/s_scalbn.c
> > index 219cd8f0c989..0d344840862f 100644
> > --- a/lib/msun/src/s_scalbn.c
> > +++ b/lib/msun/src/s_scalbn.c
> > @@ -29,6 +29,19 @@ double scalbn(double x, int n)
> >         }
> >         u.i = (uint64_t)(0x3ff+n)<<52;
> >         x = y * u.f;
> > +#if !defined(__amd64__) && !defined(__i386__)
> > +       /*
> > +        * x86 performs the multiplication with higher precision, but on
> > +        * non-x86 architectures we might get 0 instead of a tiny value. To work
> > +        * around this problem perform the multiplication with float128 (slow).
> > +        * TODO: This doesn't work on e.g. MIPS where long double == double.
> > +        */
> > +       if (x == 0.) {
> > +               x = (long double)y * u.f;
> > +               /* fprintf(stderr, "\ttrying again: %a/%a\n", x,
> > (double)((long double)y * u.f)); */
> > +               return x;
> > +       }
> > +#endif
> >         return x;
> >  }

^ permalink raw reply	[flat|nested] 2+ messages in thread

end of thread, other threads:[~2021-02-26 12:49 UTC | newest]

Thread overview: 2+ messages (download: mbox.gz / follow: Atom feed)
-- links below jump to the message on this page --
     [not found] <CAEeofcgZbgsP6D-idQXFv6qqQw-d0kVeBHiCi2xqEeSoD-JcZw@mail.gmail.com>
2021-02-26 11:43 ` [musl] Re: Potentially incorrect musl scalbn results on AArch64? Szabolcs Nagy
2021-02-26 12:45   ` Alexander Richardson

mailing list of musl libc

This inbox may be cloned and mirrored by anyone:

	git clone --mirror http://inbox.vuxu.org/musl

	# If you have public-inbox 1.1+ installed, you may
	# initialize and index your mirror using the following commands:
	public-inbox-init -V1 musl musl/ http://inbox.vuxu.org/musl \
		musl@inbox.vuxu.org
	public-inbox-index musl

Example config snippet for mirrors.
Newsgroup available over NNTP:
	nntp://inbox.vuxu.org/vuxu.archive.musl


code repositories for the project(s) associated with this inbox:

	https://git.vuxu.org/mirror/musl/

AGPL code for this site: git clone https://public-inbox.org/public-inbox.git