xref: /freebsd/contrib/llvm-project/libc/src/__support/FPUtil/generic/add_sub.h (revision bb722a7d0f1642bff6487f943ad0427799a6e5bf)
1 //===-- Add and subtract IEEE 754 floating-point numbers --------*- C++ -*-===//
2 //
3 // Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
4 // See https://llvm.org/LICENSE.txt for license information.
5 // SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
6 //
7 //===----------------------------------------------------------------------===//
8 
9 #ifndef LLVM_LIBC_SRC___SUPPORT_FPUTIL_GENERIC_ADD_SUB_H
10 #define LLVM_LIBC_SRC___SUPPORT_FPUTIL_GENERIC_ADD_SUB_H
11 
12 #include "hdr/fenv_macros.h"
13 #include "src/__support/CPP/algorithm.h"
14 #include "src/__support/CPP/bit.h"
15 #include "src/__support/CPP/type_traits.h"
16 #include "src/__support/FPUtil/BasicOperations.h"
17 #include "src/__support/FPUtil/FEnvImpl.h"
18 #include "src/__support/FPUtil/FPBits.h"
19 #include "src/__support/FPUtil/cast.h"
20 #include "src/__support/FPUtil/dyadic_float.h"
21 #include "src/__support/FPUtil/rounding_mode.h"
22 #include "src/__support/macros/attributes.h"
23 #include "src/__support/macros/config.h"
24 #include "src/__support/macros/optimization.h"
25 
26 namespace LIBC_NAMESPACE_DECL {
27 namespace fputil::generic {
28 
29 template <bool IsSub, typename OutType, typename InType>
30 LIBC_INLINE cpp::enable_if_t<cpp::is_floating_point_v<OutType> &&
31                                  cpp::is_floating_point_v<InType> &&
32                                  sizeof(OutType) <= sizeof(InType),
33                              OutType>
add_or_sub(InType x,InType y)34 add_or_sub(InType x, InType y) {
35   using OutFPBits = FPBits<OutType>;
36   using OutStorageType = typename OutFPBits::StorageType;
37   using InFPBits = FPBits<InType>;
38   using InStorageType = typename InFPBits::StorageType;
39 
40   constexpr int GUARD_BITS_LEN = 3;
41   constexpr int RESULT_FRACTION_LEN = InFPBits::FRACTION_LEN + GUARD_BITS_LEN;
42   constexpr int RESULT_MANTISSA_LEN = RESULT_FRACTION_LEN + 1;
43 
44   using DyadicFloat =
45       DyadicFloat<cpp::bit_ceil(static_cast<size_t>(RESULT_MANTISSA_LEN))>;
46 
47   InFPBits x_bits(x);
48   InFPBits y_bits(y);
49 
50   bool is_effectively_add = (x_bits.sign() == y_bits.sign()) != IsSub;
51 
52   if (LIBC_UNLIKELY(x_bits.is_inf_or_nan() || y_bits.is_inf_or_nan() ||
53                     x_bits.is_zero() || y_bits.is_zero())) {
54     if (x_bits.is_nan() || y_bits.is_nan()) {
55       if (x_bits.is_signaling_nan() || y_bits.is_signaling_nan())
56         raise_except_if_required(FE_INVALID);
57 
58       if (x_bits.is_quiet_nan()) {
59         InStorageType x_payload = x_bits.get_mantissa();
60         x_payload >>= InFPBits::FRACTION_LEN - OutFPBits::FRACTION_LEN;
61         return OutFPBits::quiet_nan(x_bits.sign(),
62                                     static_cast<OutStorageType>(x_payload))
63             .get_val();
64       }
65 
66       if (y_bits.is_quiet_nan()) {
67         InStorageType y_payload = y_bits.get_mantissa();
68         y_payload >>= InFPBits::FRACTION_LEN - OutFPBits::FRACTION_LEN;
69         return OutFPBits::quiet_nan(y_bits.sign(),
70                                     static_cast<OutStorageType>(y_payload))
71             .get_val();
72       }
73 
74       return OutFPBits::quiet_nan().get_val();
75     }
76 
77     if (x_bits.is_inf()) {
78       if (y_bits.is_inf()) {
79         if (!is_effectively_add) {
80           raise_except_if_required(FE_INVALID);
81           return OutFPBits::quiet_nan().get_val();
82         }
83 
84         return OutFPBits::inf(x_bits.sign()).get_val();
85       }
86 
87       return OutFPBits::inf(x_bits.sign()).get_val();
88     }
89 
90     if (y_bits.is_inf())
91       return OutFPBits::inf(y_bits.sign()).get_val();
92 
93     if (x_bits.is_zero()) {
94       if (y_bits.is_zero()) {
95         switch (quick_get_round()) {
96         case FE_DOWNWARD:
97           return OutFPBits::zero(Sign::NEG).get_val();
98         default:
99           return OutFPBits::zero(Sign::POS).get_val();
100         }
101       }
102 
103       // volatile prevents Clang from converting tmp to OutType and then
104       // immediately back to InType before negating it, resulting in double
105       // rounding.
106       volatile InType tmp = y;
107       if constexpr (IsSub)
108         tmp = -tmp;
109       return cast<OutType>(tmp);
110     }
111 
112     if (y_bits.is_zero())
113       return cast<OutType>(x);
114   }
115 
116   InType x_abs = x_bits.abs().get_val();
117   InType y_abs = y_bits.abs().get_val();
118 
119   if (x_abs == y_abs && !is_effectively_add) {
120     switch (quick_get_round()) {
121     case FE_DOWNWARD:
122       return OutFPBits::zero(Sign::NEG).get_val();
123     default:
124       return OutFPBits::zero(Sign::POS).get_val();
125     }
126   }
127 
128   Sign result_sign = Sign::POS;
129 
130   if (x_abs > y_abs) {
131     result_sign = x_bits.sign();
132   } else if (x_abs < y_abs) {
133     if (is_effectively_add)
134       result_sign = y_bits.sign();
135     else if (y_bits.is_pos())
136       result_sign = Sign::NEG;
137   } else if (is_effectively_add) {
138     result_sign = x_bits.sign();
139   }
140 
141   InFPBits max_bits(cpp::max(x_abs, y_abs));
142   InFPBits min_bits(cpp::min(x_abs, y_abs));
143 
144   InStorageType result_mant;
145 
146   if (max_bits.is_subnormal()) {
147     // min_bits must be subnormal too.
148 
149     if (is_effectively_add)
150       result_mant = max_bits.get_mantissa() + min_bits.get_mantissa();
151     else
152       result_mant = max_bits.get_mantissa() - min_bits.get_mantissa();
153 
154     result_mant <<= GUARD_BITS_LEN;
155   } else {
156     InStorageType max_mant = max_bits.get_explicit_mantissa() << GUARD_BITS_LEN;
157     InStorageType min_mant = min_bits.get_explicit_mantissa() << GUARD_BITS_LEN;
158 
159     int alignment = (max_bits.get_biased_exponent() - max_bits.is_normal()) -
160                     (min_bits.get_biased_exponent() - min_bits.is_normal());
161 
162     InStorageType aligned_min_mant =
163         min_mant >> cpp::min(alignment, RESULT_MANTISSA_LEN);
164     bool aligned_min_mant_sticky;
165 
166     if (alignment <= GUARD_BITS_LEN)
167       aligned_min_mant_sticky = false;
168     else if (alignment > InFPBits::FRACTION_LEN + GUARD_BITS_LEN)
169       aligned_min_mant_sticky = true;
170     else
171       aligned_min_mant_sticky =
172           (static_cast<InStorageType>(
173               min_mant << (InFPBits::STORAGE_LEN - alignment))) != 0;
174 
175     InStorageType min_mant_sticky(static_cast<int>(aligned_min_mant_sticky));
176 
177     if (is_effectively_add)
178       result_mant = max_mant + (aligned_min_mant | min_mant_sticky);
179     else
180       result_mant = max_mant - (aligned_min_mant | min_mant_sticky);
181   }
182 
183   int result_exp = max_bits.get_explicit_exponent() - RESULT_FRACTION_LEN;
184   DyadicFloat result(result_sign, result_exp, result_mant);
185   return result.template as<OutType, /*ShouldSignalExceptions=*/true>();
186 }
187 
188 template <typename OutType, typename InType>
189 LIBC_INLINE cpp::enable_if_t<cpp::is_floating_point_v<OutType> &&
190                                  cpp::is_floating_point_v<InType> &&
191                                  sizeof(OutType) <= sizeof(InType),
192                              OutType>
add(InType x,InType y)193 add(InType x, InType y) {
194   return add_or_sub</*IsSub=*/false, OutType>(x, y);
195 }
196 
197 template <typename OutType, typename InType>
198 LIBC_INLINE cpp::enable_if_t<cpp::is_floating_point_v<OutType> &&
199                                  cpp::is_floating_point_v<InType> &&
200                                  sizeof(OutType) <= sizeof(InType),
201                              OutType>
sub(InType x,InType y)202 sub(InType x, InType y) {
203   return add_or_sub</*IsSub=*/true, OutType>(x, y);
204 }
205 
206 } // namespace fputil::generic
207 } // namespace LIBC_NAMESPACE_DECL
208 
209 #endif // LLVM_LIBC_SRC___SUPPORT_FPUTIL_GENERIC_ADD_SUB_H
210