https://gcc.gnu.org/bugzilla/show_bug.cgi?id=89991

--- Comment #20 from Jonathan Wakely <redi at gcc dot gnu.org> ---
(In reply to Steve Kargl from comment #16)
> If Andrew is correct and a builtin is called,

Wasn't that me, not Andrew?

> you might find
> my results if you use -fno-builtins (check spelling).

No, same results. Calling __builtin_csqrt doesn't necessarily mean GCC
evaluates it. Without optimisation it still generates a call to csqrt from
libc.

> Looking at ./libstdc++-v3/include/std/complex, one finds.
> 
>   // 26.2.8/13 sqrt(__z): Returns the complex square root of __z.
>   //                     The branch cut is on the negative axis.
>   template<typename _Tp>
>     complex<_Tp>
>     __complex_sqrt(const complex<_Tp>& __z)
>     {
>       _Tp __x = __z.real();
>       _Tp __y = __z.imag();
> 
>       if (__x == _Tp())
>         {
>           _Tp __t = sqrt(abs(__y) / 2);
>           return complex<_Tp>(__t, __y < _Tp() ? -__t : __t);
>         }
>       else
>         {
>           _Tp __t = sqrt(2 * (std::abs(__z) + abs(__x)));
>           _Tp __u = __t / 2;
>           return __x > _Tp()
>             ? complex<_Tp>(__u, __y / __t)
>             : complex<_Tp>(abs(__y) / __t, __y < _Tp() ? -__u : __u);
>         }
>     }
> 
> Doesn't this gets the wrong answer for __y = -0 (as -0 < 0 is false)?

Yes, but that code shouldn't be used for modern targets ...

(In reply to kargl from comment #19)
> I get the expected.  So, if you're on a system that has
> _GLIBCXX_USE_C99_COMPLEX, you won't see the bug.

Wow, why isn't libstdc++ using the C99 <math.h> functions on FreeBSD?

I'll have to look into that.

> It is likely that everywhere that a construct of the
> form __y < _Tp() ? -__u : __u appear, it needs to use
> copysign.

That won't always work, because the generic functions should really only get
used when _Tp is not one of float, double or long double. And in that case
there might be no copysign for the type.

For float, double and long double we should be using the libc routines. So the
bug is that FreeBSD isn't using them.

The _original_ bug report is for std::pow though, and on Ubuntu, which does use
glibc. Comment 9 needs more analysis.

Reply via email to