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

Jonathan Wakely <redi at gcc dot gnu.org> changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
     Ever confirmed|0                           |1
           Keywords|                            |wrong-code
   Last reconfirmed|                            |2025-08-20
             Status|UNCONFIRMED                 |NEW

--- Comment #1 from Jonathan Wakely <redi at gcc dot gnu.org> ---
The standard says "The effect of instantiating the primary template of complex
for any type that is not a cv-unqualified floating-point type (6.9.2) is
unspecified."

But we can do better anyway.

The purpose of the divisions by __s is to avoid overflow in the x*x+y+y
expression, but maybe we should assume no overflow happens, and only do that if
needed:

  template<typename _Tp>
    inline _Tp
    __complex_abs(const complex<_Tp>& __z)
    {
      _Tp __x = __z.real();
      _Tp __y = __z.imag();

      _Tp __xx, __yy, __sum;
      if (!__builtin_mul_overflow(__x, __x, &__xx))
        if (!__builtin_mul_overflow(__y, __y, &__yy))
          if (!__builtin_add_overflow(__xx, __yy, &__sum))
            return std::sqrt(__sum);

      // as before ...

This would only use the truncating integer divisions for cases where overflow
would happen.

We could also use double precision FP arithmetic for the cases where overflow
would happen:

  // 26.2.7/3 abs(__z):  Returns the magnitude of __z.
  template<typename _Tp>
    inline _Tp
    __complex_abs(const complex<_Tp>& __z)
    {
      _Tp __x = __z.real();
      _Tp __y = __z.imag();

      _Tp __xx, __yy, __sum;
      if (!__builtin_mul_overflow(__x, __x, &__xx))
        if (!__builtin_mul_overflow(__y, __y, &__yy))
          if (!__builtin_add_overflow(__xx, __yy, &__sum))
            return std::sqrt(__sum);

      const double __s = std::max(abs(__x), abs(__y));
      if (__s == _Tp())  // well ...
        return __s;
      double __x_s = __x / __s;
      double __y_s = __y / __s;
      return __s * sqrt(__x_s * __x_s + __y_s * __y_s);
    }

This would be consistent with the general rule in the standard that integer
arguments are promoted to double for calls to <cmath> functions such as
std::pow.

Reply via email to