Inverse square root. This function is exact to rounding for all rounding modes and thus generally more accurate than directly computing 1 / sqrt(\a arg) in half-precision, in addition to also being faster. \param arg function argument \return reciprocal of square root of \a arg \exception FE_INVALID for signaling NaN and negative arguments \exception FE_INEXACT according to rounding
| 3802 | /// reciprocal of square root of \a arg \exception FE_INVALID for signaling NaN |
| 3803 | /// and negative arguments \exception FE_INEXACT according to rounding |
| 3804 | inline half rsqrt(half arg) { |
| 3805 | #ifdef HALF_ARITHMETIC_TYPE |
| 3806 | return half( |
| 3807 | detail::binary, |
| 3808 | detail::float2half<half::round_style>( |
| 3809 | detail::internal_t(1) / |
| 3810 | std::sqrt(detail::half2float<detail::internal_t>(arg.data_)))); |
| 3811 | #else |
| 3812 | unsigned int abs = arg.data_ & 0x7FFF, bias = 0x4000; |
| 3813 | if (!abs || arg.data_ >= 0x7C00) |
| 3814 | return half(detail::binary, (abs > 0x7C00) ? detail::signal(arg.data_) |
| 3815 | : (arg.data_ > 0x8000) ? detail::invalid() |
| 3816 | : !abs ? detail::pole(arg.data_ & 0x8000) |
| 3817 | : 0); |
| 3818 | for (; abs < 0x400; abs <<= 1, bias -= 0x400) |
| 3819 | ; |
| 3820 | unsigned int frac = (abs += bias) & 0x7FF; |
| 3821 | if (frac == 0x400) return half(detail::binary, 0x7A00 - (abs >> 1)); |
| 3822 | if ((half::round_style == std::round_to_nearest && |
| 3823 | (frac == 0x3FE || frac == 0x76C)) || |
| 3824 | (half::round_style != std::round_to_nearest && |
| 3825 | (frac == 0x15A || frac == 0x3FC || frac == 0x401 || frac == 0x402 || |
| 3826 | frac == 0x67B))) |
| 3827 | return pow(arg, half(detail::binary, 0xB800)); |
| 3828 | detail::uint32 f = 0x17376 - abs, mx = (abs & 0x3FF) | 0x400, |
| 3829 | my = ((f >> 1) & 0x3FF) | 0x400, mz = my * my; |
| 3830 | int expy = (f >> 11) - 31, expx = 32 - (abs >> 10), i = mz >> 21; |
| 3831 | for (mz = 0x60000000 - (((mz >> i) * mx) >> (expx - 2 * expy - i)); |
| 3832 | mz < 0x40000000; mz <<= 1, --expy) |
| 3833 | ; |
| 3834 | i = (my *= mz >> 10) >> 31; |
| 3835 | expy += i; |
| 3836 | my = (my >> (20 + i)) + 1; |
| 3837 | i = (mz = my * my) >> 21; |
| 3838 | for (mz = 0x60000000 - (((mz >> i) * mx) >> (expx - 2 * expy - i)); |
| 3839 | mz < 0x40000000; mz <<= 1, --expy) |
| 3840 | ; |
| 3841 | i = (my *= (mz >> 10) + 1) >> 31; |
| 3842 | return half(detail::binary, |
| 3843 | detail::fixed2half<half::round_style, 30, false, false, true>( |
| 3844 | my >> i, expy + i + 14)); |
| 3845 | #endif |
| 3846 | } |
| 3847 | |
| 3848 | /// Cubic root. |
| 3849 | /// This function is exact to rounding for all rounding modes. |