22// Distributed under the MIT software license, see the accompanying
33// file COPYING or http://www.opensource.org/licenses/mit-license.php.
44
5+ #include < arith_uint256.h>
56#include < util/feefrac.h>
67#include < test/fuzz/FuzzedDataProvider.h>
78#include < test/fuzz/fuzz.h>
1314
1415namespace {
1516
16- /* * Compute a * b, represented in 4x32 bits, highest limb first. */
17- std::array<uint32_t , 4 > Mul128 (uint64_t a, uint64_t b)
18- {
19- std::array<uint32_t , 4 > ret{0 , 0 , 0 , 0 };
20-
21- /* * Perform ret += v << (32 * pos), at 128-bit precision. */
22- auto add_fn = [&](uint64_t v, int pos) {
23- uint64_t accum{0 };
24- for (int i = 0 ; i + pos < 4 ; ++i) {
25- // Add current value at limb pos in ret.
26- accum += ret[3 - pos - i];
27- // Add low or high half of v.
28- if (i == 0 ) accum += v & 0xffffffff ;
29- if (i == 1 ) accum += v >> 32 ;
30- // Store lower half of result in limb pos in ret.
31- ret[3 - pos - i] = accum & 0xffffffff ;
32- // Leave carry in accum.
33- accum >>= 32 ;
34- }
35- // Make sure no overflow.
36- assert (accum == 0 );
37- };
38-
39- // Multiply the 4 individual limbs (schoolbook multiply, with base 2^32).
40- add_fn ((a & 0xffffffff ) * (b & 0xffffffff ), 0 );
41- add_fn ((a >> 32 ) * (b & 0xffffffff ), 1 );
42- add_fn ((a & 0xffffffff ) * (b >> 32 ), 1 );
43- add_fn ((a >> 32 ) * (b >> 32 ), 2 );
44- return ret;
45- }
17+ /* * The maximum absolute value of an int64_t, as an arith_uint256 (2^63). */
18+ const auto MAX_ABS_INT64 = arith_uint256{1 } << 63 ;
4619
47- /* comparison helper for std::array */
48- std::strong_ordering compare_arrays (const std::array<uint32_t , 4 >& a, const std::array<uint32_t , 4 >& b) {
49- for (size_t i = 0 ; i < a.size (); ++i) {
50- if (a[i] != b[i]) return a[i] <=> b[i];
20+ /* * Construct an arith_uint256 whose value equals abs(x). */
21+ arith_uint256 Abs256 (int64_t x)
22+ {
23+ if (x >= 0 ) {
24+ // For positive numbers, pass through the value.
25+ return arith_uint256{static_cast <uint64_t >(x)};
26+ } else if (x > std::numeric_limits<int64_t >::min ()) {
27+ // For negative numbers, negate first.
28+ return arith_uint256{static_cast <uint64_t >(-x)};
29+ } else {
30+ // Special case for x == -2^63 (for which -x results in integer overflow).
31+ return MAX_ABS_INT64;
5132 }
52- return std::strong_ordering::equal;
5333}
5434
5535std::strong_ordering MulCompare (int64_t a1, int64_t a2, int64_t b1, int64_t b2)
@@ -59,23 +39,14 @@ std::strong_ordering MulCompare(int64_t a1, int64_t a2, int64_t b1, int64_t b2)
5939 int sign_b = (b1 == 0 ? 0 : b1 < 0 ? -1 : 1 ) * (b2 == 0 ? 0 : b2 < 0 ? -1 : 1 );
6040 if (sign_a != sign_b) return sign_a <=> sign_b;
6141
62- // Compute absolute values.
63- uint64_t abs_a1 = static_cast <uint64_t >(a1), abs_a2 = static_cast <uint64_t >(a2);
64- uint64_t abs_b1 = static_cast <uint64_t >(b1), abs_b2 = static_cast <uint64_t >(b2);
65- // Use (~x + 1) instead of the equivalent (-x) to silence the linter; mod 2^64 behavior is
66- // intentional here.
67- if (a1 < 0 ) abs_a1 = ~abs_a1 + 1 ;
68- if (a2 < 0 ) abs_a2 = ~abs_a2 + 1 ;
69- if (b1 < 0 ) abs_b1 = ~abs_b1 + 1 ;
70- if (b2 < 0 ) abs_b2 = ~abs_b2 + 1 ;
42+ // Compute absolute values of products.
43+ auto mul_abs_a = Abs256 (a1) * Abs256 (a2), mul_abs_b = Abs256 (b1) * Abs256 (b2);
7144
7245 // Compute products of absolute values.
73- auto mul_abs_a = Mul128 (abs_a1, abs_a2);
74- auto mul_abs_b = Mul128 (abs_b1, abs_b2);
7546 if (sign_a < 0 ) {
76- return compare_arrays ( mul_abs_b, mul_abs_a) ;
47+ return mul_abs_b <=> mul_abs_a;
7748 } else {
78- return compare_arrays ( mul_abs_a, mul_abs_b) ;
49+ return mul_abs_a <=> mul_abs_b;
7950 }
8051}
8152
0 commit comments