Re: ccmath::copysign/signbit

From: Gero Peterhoff (g.peterhoff_at_[hidden])
Date: 2023-04-22 09:03:15

Hi Matt,
I don't think you understood what I was getting at. I mean:

> Gero,
> 
>> 1) copysign uses abs, which is wrong. copysign/signbit may *only* change/check the sign (https://en.cppreference.com/w/cpp/numeric/math/copysign https://en.cppreference.com/w/cpp/numeric/math/signbit)!
> 
> I think you are conflating signbit and copysign here. Copysign is to return the magnitude of the first argument (which by definition is the absolute value) with the sign of the second argument. We use signbit in the implementation of this because copysign can be used to manipulate the signs of NaNs.
> 

Your copysign uses abs (https://github.com/boostorg/math/blob/develop/include/boost/math/ccmath/copysign.hpp). This is wrong, because it cannot be guaranteed (abs) that *only* the sign is manipulated. e.g:

using T = boost::float64_t;
constexpr T
   mag = std::numeric_limits<T>::signaling_NaN(),
   sgn = -42,
   a = std::copysign(mag, sgn),
   b = boost::math::ccmath::copysign(mag, sgn);

std::cout << "std\t" << std::bit_cast<int64_t>(a) << "\n";
std::cout << "cc\t"  << std::bit_cast<int64_t>(b) << "\n";

std	-3377699720527872
cc	-2251799813685248

The correct procedure is (as with signbit):

template <typename T>
constexpr Type copysign_impl(const T mag, const T sgn) noexcept
{
      if constexpr (std::is_same_v<T, float>)
      {
         auto bits = BOOST_MATH_BIT_CAST(IEEEf2bits, mag);
        bits.sign = signbit_impl(sgn);
         return BOOST_MATH_BIT_CAST(T, bits);
      }
      else if constexpr (std::is_same_v<T, double>)
      {
        auto bits = BOOST_MATH_BIT_CAST(IEEEd2bits, mag);
        bits.sign = signbit_impl(sgn);
         return BOOST_MATH_BIT_CAST(T, bits);
      }
      #ifndef BOOST_MATH_UNSUPPORTED_LONG_DOUBLE
      else if constexpr (std::is_same_v<T, long double>)
      {
        auto bits = BOOST_MATH_BIT_CAST(IEEEl2bits, mag);
        bits.sign = signbit_impl(sgn);
         return BOOST_MATH_BIT_CAST(T, bits);
      }
      #endif
      else
      {
          // error handling
      }
}


>> 2) In signbit.hpp you implement struct IEEEfN2bits. But these are incomplete:
>> - long double can be 16/12/10 bytes depending on platform/compiler option. You do not take this into account.
>> - Correct types are already available in boost/math/cstdfloat/cstdfloat_types.hpp.
>>
> 
> You will find implementations of IEEEl2bits here https://github.com/boostorg/math/blob/develop/include/boost/math/ccmath/signbit.hpp#L94 <https://github.com/boostorg/math/blob/develop/include/boost/math/ccmath/signbit.hpp#L94>. They account for 64, 80, and 128 bit long doubles based on LDBL_MANT_DIG and LDBL_MAX_EXP. If for instance you use compiler flags on GCC to alter the size of a long double those macros are changed accordingly.
>

This is incomplete. In case long double == boost::float80_t (LDBL_MANT_DIG==64 && LDBL_MAX_EXP==16384) the storagesize (sizeof(long double)) may be different depending on platform/system/compiler option/etc, see e.g. https://gcc.gnu.org/onlinedocs/gcc/x86-Options.html section "-m128bit-long-double":
- x86_64 sizeof(long double) == 16
- x86_32 sizeof(long double) == 12
- if necessary actually sizeof(long double) == 10
Your 80 bit long double IEEEl2bits has *always* a size of 16 bytes, so BOOST_MATH_BIT_CAST fails if the storagesize does not fit. So for this case 3 different variants are needed:
sizeof(long double) == 16 -> 48 padding bits
sizeof(long double) == 12 -> 16 padding bits
sizeof(long double) == 10 -> no padding


I hope I was able to clear this up.

thx
Gero