1 ///////////////////////////////////////////////////////////////
2 //  Copyright 2012-20 John Maddock.
3 //  Copyright 2019-20 Christopher Kormanyos.
4 //  Copyright 2019-20 Madhur Chauhan.
5 //  Distributed under the Boost Software License, Version 1.0.
6 //  (See accompanying file LICENSE_1_0.txt or copy at https://www.boost.org/LICENSE_1_0.txt
7 //
8 // Comparison operators for cpp_int_backend:
9 //
10 #ifndef BOOST_MP_CPP_INT_MUL_HPP
11 #define BOOST_MP_CPP_INT_MUL_HPP
12 
13 #include <limits>
14 #include <boost/multiprecision/detail/endian.hpp>
15 #include <boost/multiprecision/detail/assert.hpp>
16 #include <boost/multiprecision/integer.hpp>
17 
18 namespace boost { namespace multiprecision { namespace backends {
19 
20 #ifdef _MSC_VER
21 #pragma warning(push)
22 #pragma warning(disable : 4127) // conditional expression is constant
23 #endif
24 //
25 // Multiplication by a single limb:
26 //
27 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, unsigned MinBits2, unsigned MaxBits2, cpp_integer_type SignType2, cpp_int_check_type Checked2, class Allocator2>
28 inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && !is_trivial_cpp_int<cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const cpp_int_backend<MinBits2,MaxBits2,SignType2,Checked2,Allocator2> & a,const limb_type & val)29 eval_multiply(
30     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&       result,
31     const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>& a,
32     const limb_type&                                                            val) noexcept((is_non_throwing_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value))
33 {
34    if (!val)
35    {
36       result = static_cast<limb_type>(0);
37       return;
38    }
39    if ((void*)&a != (void*)&result)
40       result.resize(a.size(), a.size());
41    double_limb_type                                                                                  carry = 0;
42    typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_pointer       p     = result.limbs();
43    typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_pointer       pe    = result.limbs() + result.size();
44    typename cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>::const_limb_pointer pa    = a.limbs();
45    while (p != pe)
46    {
47       carry += static_cast<double_limb_type>(*pa) * static_cast<double_limb_type>(val);
48 #ifdef __MSVC_RUNTIME_CHECKS
49       *p = static_cast<limb_type>(carry & ~static_cast<limb_type>(0));
50 #else
51       *p = static_cast<limb_type>(carry);
52 #endif
53       carry >>= cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits;
54       ++p, ++pa;
55    }
56    if (carry)
57    {
58       unsigned i = result.size();
59       result.resize(i + 1, i + 1);
60       if (result.size() > i)
61          result.limbs()[i] = static_cast<limb_type>(carry);
62    }
63    result.sign(a.sign());
64    if (is_fixed_precision<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value)
65       result.normalize();
66 }
67 
68 //
69 // resize_for_carry forces a resize of the underlying buffer only if a previous request
70 // for "required" elements could possibly have failed, *and* we have checking enabled.
71 // This will cause an overflow error inside resize():
72 //
73 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1>
resize_for_carry(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> &,unsigned)74 inline BOOST_MP_CXX14_CONSTEXPR void resize_for_carry(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& /*result*/, unsigned /*required*/) {}
75 
76 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, class Allocator1>
resize_for_carry(cpp_int_backend<MinBits1,MaxBits1,SignType1,checked,Allocator1> & result,unsigned required)77 inline BOOST_MP_CXX14_CONSTEXPR void resize_for_carry(cpp_int_backend<MinBits1, MaxBits1, SignType1, checked, Allocator1>& result, unsigned required)
78 {
79    if (result.size() < required)
80       result.resize(required, required);
81 }
82 //
83 // Minimum number of limbs required for Karatsuba to be worthwhile:
84 //
85 #ifdef BOOST_MP_KARATSUBA_CUTOFF
86 const size_t karatsuba_cutoff = BOOST_MP_KARATSUBA_CUTOFF;
87 #else
88 const size_t karatsuba_cutoff = 40;
89 #endif
90 //
91 // Core (recursive) Karatsuba multiplication, all the storage required is allocated upfront and
92 // passed down the stack in this routine.  Note that all the cpp_int_backend's must be the same type
93 // and full variable precision.  Karatsuba really doesn't play nice with fixed-size integers.  If necessary
94 // fixed precision integers will get aliased as variable-precision types before this is called.
95 //
96 template <unsigned MinBits, unsigned MaxBits, cpp_int_check_type Checked, class Allocator>
multiply_karatsuba(cpp_int_backend<MinBits,MaxBits,signed_magnitude,Checked,Allocator> & result,const cpp_int_backend<MinBits,MaxBits,signed_magnitude,Checked,Allocator> & a,const cpp_int_backend<MinBits,MaxBits,signed_magnitude,Checked,Allocator> & b,typename cpp_int_backend<MinBits,MaxBits,signed_magnitude,Checked,Allocator>::scoped_shared_storage & storage)97 inline void multiply_karatsuba(
98     cpp_int_backend<MinBits, MaxBits, signed_magnitude, Checked, Allocator>&       result,
99     const cpp_int_backend<MinBits, MaxBits, signed_magnitude, Checked, Allocator>& a,
100     const cpp_int_backend<MinBits, MaxBits, signed_magnitude, Checked, Allocator>& b,
101     typename cpp_int_backend<MinBits, MaxBits, signed_magnitude, Checked, Allocator>::scoped_shared_storage& storage)
102 {
103    using cpp_int_type = cpp_int_backend<MinBits, MaxBits, signed_magnitude, Checked, Allocator>;
104 
105    unsigned as = a.size();
106    unsigned bs = b.size();
107    //
108    // Termination condition: if either argument is smaller than karatsuba_cutoff
109    // then schoolboy multiplication will be faster:
110    //
111    if ((as < karatsuba_cutoff) || (bs < karatsuba_cutoff))
112    {
113       eval_multiply(result, a, b);
114       return;
115    }
116    //
117    // Partitioning size: split the larger of a and b into 2 halves
118    //
119    unsigned n  = (as > bs ? as : bs) / 2 + 1;
120    //
121    // Partition a and b into high and low parts.
122    // ie write a, b as a = a_h * 2^n + a_l, b = b_h * 2^n + b_l
123    //
124    // We could copy the high and low parts into new variables, but we'll
125    // use aliasing to reference the internal limbs of a and b.  There is one wart here:
126    // if a and b are mismatched in size, then n may be larger than the smaller
127    // of a and b.  In that situation the high part is zero, and we have no limbs
128    // to alias, so instead alias a local variable.
129    // This raises 2 questions:
130    // * Is this the best way to partition a and b?
131    // * Since we have one high part zero, the arithmetic simplifies considerably,
132    //   so should we have a special routine for this?
133    //
134    unsigned          sz = (std::min)(as, n);
135    const cpp_int_type a_l(a.limbs(), 0, sz);
136 
137    sz = (std::min)(bs, n);
138    const cpp_int_type b_l(b.limbs(), 0, sz);
139 
140    limb_type          zero = 0;
141    const cpp_int_type a_h(as > n ? a.limbs() + n : &zero, 0, as > n ? as - n : 1);
142    const cpp_int_type b_h(bs > n ? b.limbs() + n : &zero, 0, bs > n ? bs - n : 1);
143    //
144    // The basis for the Karatsuba algorithm is as follows:
145    //
146    // let                x = a_h * b_ h
147    //                    y = a_l * b_l
148    //                    z = (a_h + a_l)*(b_h + b_l) - x - y
149    // and therefore  a * b = x * (2 ^ (2 * n))+ z * (2 ^ n) + y
150    //
151    // Begin by allocating our temporaries, these alias the memory already allocated in the shared storage:
152    //
153    cpp_int_type t1(storage, 2 * n + 2);
154    cpp_int_type t2(storage, n + 1);
155    cpp_int_type t3(storage, n + 1);
156    //
157    // Now we want:
158    //
159    // result = | a_h*b_h  | a_l*b_l |
160    // (bits)              <-- 2*n -->
161    //
162    // We create aliases for the low and high parts of result, and multiply directly into them:
163    //
164    cpp_int_type result_low(result.limbs(), 0, 2 * n);
165    cpp_int_type result_high(result.limbs(), 2 * n, result.size() - 2 * n);
166    //
167    // low part of result is a_l * b_l:
168    //
169    multiply_karatsuba(result_low, a_l, b_l, storage);
170    //
171    // We haven't zeroed out memory in result, so set to zero any unused limbs,
172    // if a_l and b_l have mostly random bits then nothing happens here, but if
173    // one is zero or nearly so, then a memset might be faster... it's not clear
174    // that it's worth the extra logic though (and is darn hard to measure
175    // what the "average" case is).
176    //
177    for (unsigned i = result_low.size(); i < 2 * n; ++i)
178       result.limbs()[i] = 0;
179    //
180    // Set the high part of result to a_h * b_h:
181    //
182    multiply_karatsuba(result_high, a_h, b_h, storage);
183    for (unsigned i = result_high.size() + 2 * n; i < result.size(); ++i)
184       result.limbs()[i] = 0;
185    //
186    // Now calculate (a_h+a_l)*(b_h+b_l):
187    //
188    add_unsigned(t2, a_l, a_h);
189    add_unsigned(t3, b_l, b_h);
190    multiply_karatsuba(t1, t2, t3, storage); // t1 = (a_h+a_l)*(b_h+b_l)
191    //
192    // There is now a slight deviation from Karatsuba, we want to subtract
193    // a_l*b_l + a_h*b_h from t1, but rather than use an addition and a subtraction
194    // plus one temporary, we'll use 2 subtractions.  On the minus side, a subtraction
195    // is on average slightly slower than an addition, but we save a temporary (ie memory)
196    // and also hammer the same piece of memory over and over rather than 2 disparate
197    // memory regions.  Overall it seems to be a slight win.
198    //
199    subtract_unsigned(t1, t1, result_high);
200    subtract_unsigned(t1, t1, result_low);
201    //
202    // The final step is to left shift t1 by n bits and add to the result.
203    // Rather than do an actual left shift, we can simply alias the result
204    // and add to the alias:
205    //
206    cpp_int_type result_alias(result.limbs(), n, result.size() - n);
207    add_unsigned(result_alias, result_alias, t1);
208    //
209    // Free up storage for use by sister branches to this one:
210    //
211    storage.deallocate(t1.capacity() + t2.capacity() + t3.capacity());
212 
213    result.normalize();
214 }
215 
karatsuba_storage_size(unsigned s)216 inline unsigned karatsuba_storage_size(unsigned s)
217 {
218    //
219    // This estimates how much memory we will need based on
220    // s-limb multiplication.  In an ideal world the number of limbs
221    // would halve with each recursion, and our storage requirements
222    // would be 4s in the limit, and rather less in practice since
223    // we bail out long before we reach one limb.  In the real world
224    // we don't quite halve s in each recursion, so this is an heuristic
225    // which over-estimates how much we need.  We could compute an exact
226    // value, but it would be rather time consuming.
227    //
228    return 5 * s;
229 }
230 //
231 // There are 2 entry point routines for Karatsuba multiplication:
232 // one for variable precision types, and one for fixed precision types.
233 // These are responsible for allocating all the storage required for the recursive
234 // routines above, and are always at the outermost level.
235 //
236 // Normal variable precision case comes first:
237 //
238 template <unsigned MinBits, unsigned MaxBits, cpp_integer_type SignType, cpp_int_check_type Checked, class Allocator>
239 inline typename std::enable_if<!is_fixed_precision<cpp_int_backend<MinBits, MaxBits, SignType, Checked, Allocator> >::value>::type
setup_karatsuba(cpp_int_backend<MinBits,MaxBits,SignType,Checked,Allocator> & result,const cpp_int_backend<MinBits,MaxBits,SignType,Checked,Allocator> & a,const cpp_int_backend<MinBits,MaxBits,SignType,Checked,Allocator> & b)240 setup_karatsuba(
241    cpp_int_backend<MinBits, MaxBits, SignType, Checked, Allocator>& result,
242    const cpp_int_backend<MinBits, MaxBits, SignType, Checked, Allocator>& a,
243    const cpp_int_backend<MinBits, MaxBits, SignType, Checked, Allocator>& b)
244 {
245    unsigned as = a.size();
246    unsigned bs = b.size();
247    unsigned s = as > bs ? as : bs;
248    unsigned storage_size = karatsuba_storage_size(s);
249    if (storage_size < 300)
250    {
251       //
252       // Special case: if we don't need too much memory, we can use stack based storage
253       // and save a call to the allocator, this allows us to use Karatsuba multiply
254       // at lower limb counts than would otherwise be possible:
255       //
256       limb_type limbs[300];
257       typename cpp_int_backend<MinBits, MaxBits, SignType, Checked, Allocator>::scoped_shared_storage storage(limbs, storage_size);
258       multiply_karatsuba(result, a, b, storage);
259    }
260    else
261    {
262       typename cpp_int_backend<MinBits, MaxBits, SignType, Checked, Allocator>::scoped_shared_storage storage(result.allocator(), storage_size);
263       multiply_karatsuba(result, a, b, storage);
264    }
265 }
266 
267 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, unsigned MinBits2, unsigned MaxBits2, cpp_integer_type SignType2, cpp_int_check_type Checked2, class Allocator2, unsigned MinBits3, unsigned MaxBits3, cpp_integer_type SignType3, cpp_int_check_type Checked3, class Allocator3>
268 inline typename std::enable_if<is_fixed_precision<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value || is_fixed_precision<cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2> >::value || is_fixed_precision<cpp_int_backend<MinBits3, MaxBits3, SignType3, Checked3, Allocator3> >::value>::type
setup_karatsuba(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const cpp_int_backend<MinBits2,MaxBits2,SignType2,Checked2,Allocator2> & a,const cpp_int_backend<MinBits3,MaxBits3,SignType3,Checked3,Allocator3> & b)269 setup_karatsuba(
270     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&       result,
271     const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>& a,
272     const cpp_int_backend<MinBits3, MaxBits3, SignType3, Checked3, Allocator3>& b)
273 {
274    //
275    // Now comes the fixed precision case.
276    // In fact Karatsuba doesn't really work with fixed precision since the logic
277    // requires that we calculate all the bits of the result (especially in the
278    // temporaries used internally).  So... we'll convert all the arguments
279    // to variable precision types by aliasing them, this also
280    // reduce the number of template instantations:
281    //
282    using variable_precision_type = cpp_int_backend<0, 0, signed_magnitude, unchecked, std::allocator<limb_type> >;
283    variable_precision_type a_t(a.limbs(), 0, a.size()), b_t(b.limbs(), 0, b.size());
284    unsigned as = a.size();
285    unsigned bs = b.size();
286    unsigned s = as > bs ? as : bs;
287    unsigned sz = as + bs;
288    unsigned storage_size = karatsuba_storage_size(s);
289 
290    if (!is_fixed_precision<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value || (sz * sizeof(limb_type) * CHAR_BIT <= MaxBits1))
291    {
292       // Result is large enough for all the bits of the result, so we can use aliasing:
293       result.resize(sz, sz);
294       variable_precision_type t(result.limbs(), 0, result.size());
295       typename variable_precision_type::scoped_shared_storage storage(t.allocator(), storage_size);
296       multiply_karatsuba(t, a_t, b_t, storage);
297       result.resize(t.size(), t.size());
298    }
299    else
300    {
301       //
302       // Not enough bit in result for the answer, so we must use a temporary
303       // and then truncate (ie modular arithmetic):
304       //
305       typename variable_precision_type::scoped_shared_storage storage(variable_precision_type::allocator_type(), sz + storage_size);
306       variable_precision_type t(storage, sz);
307       multiply_karatsuba(t, a_t, b_t, storage);
308       //
309       // If there is truncation, and result is a checked type then this will throw:
310       //
311       result = t;
312    }
313 }
314 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, unsigned MinBits2, unsigned MaxBits2, cpp_integer_type SignType2, cpp_int_check_type Checked2, class Allocator2, unsigned MinBits3, unsigned MaxBits3, cpp_integer_type SignType3, cpp_int_check_type Checked3, class Allocator3>
315 inline typename std::enable_if<!is_fixed_precision<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && !is_fixed_precision<cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2> >::value && !is_fixed_precision<cpp_int_backend<MinBits3, MaxBits3, SignType3, Checked3, Allocator3> >::value>::type
setup_karatsuba(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const cpp_int_backend<MinBits2,MaxBits2,SignType2,Checked2,Allocator2> & a,const cpp_int_backend<MinBits3,MaxBits3,SignType3,Checked3,Allocator3> & b)316 setup_karatsuba(
317    cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result,
318    const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>& a,
319    const cpp_int_backend<MinBits3, MaxBits3, SignType3, Checked3, Allocator3>& b)
320 {
321    //
322    // Variable precision, mixed arguments, just alias and forward:
323    //
324    using variable_precision_type = cpp_int_backend<0, 0, signed_magnitude, unchecked, std::allocator<limb_type> >;
325    variable_precision_type a_t(a.limbs(), 0, a.size()), b_t(b.limbs(), 0, b.size());
326    unsigned as = a.size();
327    unsigned bs = b.size();
328    unsigned s = as > bs ? as : bs;
329    unsigned sz = as + bs;
330    unsigned storage_size = karatsuba_storage_size(s);
331 
332    result.resize(sz, sz);
333    variable_precision_type t(result.limbs(), 0, result.size());
334    typename variable_precision_type::scoped_shared_storage storage(t.allocator(), storage_size);
335    multiply_karatsuba(t, a_t, b_t, storage);
336 }
337 
338 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, unsigned MinBits2, unsigned MaxBits2, cpp_integer_type SignType2, cpp_int_check_type Checked2, class Allocator2, unsigned MinBits3, unsigned MaxBits3, cpp_integer_type SignType3, cpp_int_check_type Checked3, class Allocator3>
339 inline BOOST_MP_CXX14_CONSTEXPR void
eval_multiply_comba(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const cpp_int_backend<MinBits2,MaxBits2,SignType2,Checked2,Allocator2> & a,const cpp_int_backend<MinBits3,MaxBits3,SignType3,Checked3,Allocator3> & b)340 eval_multiply_comba(
341     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&       result,
342     const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>& a,
343     const cpp_int_backend<MinBits3, MaxBits3, SignType3, Checked3, Allocator3>& b) noexcept((is_non_throwing_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value))
344 {
345    //
346    // see PR #182
347    // Comba Multiplier - based on Paul Comba's
348    // Exponentiation cryptosystems on the IBM PC, 1990
349    //
350    int as                                                                                         = a.size(),
351        bs                                                                                         = b.size(),
352        rs                                                                                         = result.size();
353    typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_pointer pr = result.limbs();
354 
355    double_limb_type carry    = 0,
356                     temp     = 0;
357    limb_type      overflow   = 0;
358    const unsigned limb_bits  = sizeof(limb_type) * CHAR_BIT;
359    const bool     must_throw = rs < as + bs - 1;
360    for (int r = 0, lim = (std::min)(rs, as + bs - 1); r < lim; ++r, overflow = 0)
361    {
362       int i = r >= as ? as - 1 : r,
363           j = r - i,
364           k = i < bs - j ? i + 1 : bs - j; // min(i+1, bs-j);
365 
366       typename cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>::const_limb_pointer pa = a.limbs() + i;
367       typename cpp_int_backend<MinBits3, MaxBits3, SignType3, Checked3, Allocator3>::const_limb_pointer pb = b.limbs() + j;
368 
369       temp = carry;
370       carry += static_cast<double_limb_type>(*(pa)) * (*(pb));
371       overflow += carry < temp;
372       for (--k; k; k--)
373       {
374          temp = carry;
375          carry += static_cast<double_limb_type>(*(--pa)) * (*(++pb));
376          overflow += carry < temp;
377       }
378       *(pr++) = static_cast<limb_type>(carry);
379       carry   = (static_cast<double_limb_type>(overflow) << limb_bits) | (carry >> limb_bits);
380    }
381    if (carry || must_throw)
382    {
383       resize_for_carry(result, as + bs);
384       if ((int)result.size() >= as + bs)
385          *pr = static_cast<limb_type>(carry);
386    }
387 }
388 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, unsigned MinBits2, unsigned MaxBits2, cpp_integer_type SignType2, cpp_int_check_type Checked2, class Allocator2, unsigned MinBits3, unsigned MaxBits3, cpp_integer_type SignType3, cpp_int_check_type Checked3, class Allocator3>
389 inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && !is_trivial_cpp_int<cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2> >::value && !is_trivial_cpp_int<cpp_int_backend<MinBits3, MaxBits3, SignType3, Checked3, Allocator3> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const cpp_int_backend<MinBits2,MaxBits2,SignType2,Checked2,Allocator2> & a,const cpp_int_backend<MinBits3,MaxBits3,SignType3,Checked3,Allocator3> & b)390 eval_multiply(
391     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&       result,
392     const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>& a,
393     const cpp_int_backend<MinBits3, MaxBits3, SignType3, Checked3, Allocator3>& b)
394    noexcept((is_non_throwing_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value
395       && (karatsuba_cutoff * sizeof(limb_type) * CHAR_BIT > MaxBits1)
396       && (karatsuba_cutoff * sizeof(limb_type)* CHAR_BIT > MaxBits2)
397       && (karatsuba_cutoff * sizeof(limb_type)* CHAR_BIT > MaxBits3)))
398 {
399    // Uses simple (O(n^2)) multiplication when the limbs are less
400    // otherwise switches to karatsuba algorithm based on experimental value (~40 limbs)
401    //
402    // Trivial cases first:
403    //
404    unsigned                                                                                          as = a.size();
405    unsigned                                                                                          bs = b.size();
406    typename cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>::const_limb_pointer pa = a.limbs();
407    typename cpp_int_backend<MinBits3, MaxBits3, SignType3, Checked3, Allocator3>::const_limb_pointer pb = b.limbs();
408    if (as == 1)
409    {
410       bool s = b.sign() != a.sign();
411       if (bs == 1)
412       {
413          result = static_cast<double_limb_type>(*pa) * static_cast<double_limb_type>(*pb);
414       }
415       else
416       {
417          limb_type l = *pa;
418          eval_multiply(result, b, l);
419       }
420       result.sign(s);
421       return;
422    }
423    if (bs == 1)
424    {
425       bool      s = b.sign() != a.sign();
426       limb_type l = *pb;
427       eval_multiply(result, a, l);
428       result.sign(s);
429       return;
430    }
431 
432    if ((void*)&result == (void*)&a)
433    {
434       cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> t(a);
435       eval_multiply(result, t, b);
436       return;
437    }
438    if ((void*)&result == (void*)&b)
439    {
440       cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> t(b);
441       eval_multiply(result, a, t);
442       return;
443    }
444 
445    constexpr const double_limb_type limb_max = ~static_cast<limb_type>(0u);
446    constexpr const double_limb_type double_limb_max = ~static_cast<double_limb_type>(0u);
447 
448    result.resize(as + bs, as + bs - 1);
449 #ifndef BOOST_MP_NO_CONSTEXPR_DETECTION
450    if (!BOOST_MP_IS_CONST_EVALUATED(as) && (as >= karatsuba_cutoff && bs >= karatsuba_cutoff))
451 #else
452    if (as >= karatsuba_cutoff && bs >= karatsuba_cutoff)
453 #endif
454    {
455       setup_karatsuba(result, a, b);
456       //
457       // Set the sign of the result:
458       //
459       result.sign(a.sign() != b.sign());
460       return;
461    }
462    typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_pointer pr = result.limbs();
463    static_assert(double_limb_max - 2 * limb_max >= limb_max * limb_max, "failed limb size sanity check");
464 
465 #ifndef BOOST_MP_NO_CONSTEXPR_DETECTION
466    if (BOOST_MP_IS_CONST_EVALUATED(as))
467    {
468       for (unsigned i = 0; i < result.size(); ++i)
469          pr[i] = 0;
470    }
471    else
472 #endif
473    std::memset(pr, 0, result.size() * sizeof(limb_type));
474 
475 #if defined(BOOST_MP_COMBA)
476        //
477        // Comba Multiplier might not be efficient because of less efficient assembly
478        // by the compiler as of 09/01/2020 (DD/MM/YY). See PR #182
479        // Till then this will lay dormant :(
480        //
481        eval_multiply_comba(result, a, b);
482 #else
483 
484    double_limb_type carry = 0;
485    for (unsigned i = 0; i < as; ++i)
486    {
487       BOOST_MP_ASSERT(result.size() > i);
488       unsigned inner_limit = !is_fixed_precision<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value ? bs : (std::min)(result.size() - i, bs);
489       unsigned j           = 0;
490       for (; j < inner_limit; ++j)
491       {
492          BOOST_MP_ASSERT(i + j < result.size());
493 #if (!defined(__GLIBCXX__) && !defined(__GLIBCPP__)) || !BOOST_WORKAROUND(BOOST_GCC_VERSION, <= 50100)
494          BOOST_MP_ASSERT(!std::numeric_limits<double_limb_type>::is_specialized || ((std::numeric_limits<double_limb_type>::max)() - carry >
495                                                                                  static_cast<double_limb_type>(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::max_limb_value) * static_cast<double_limb_type>(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::max_limb_value)));
496 #endif
497          carry += static_cast<double_limb_type>(pa[i]) * static_cast<double_limb_type>(pb[j]);
498          BOOST_MP_ASSERT(!std::numeric_limits<double_limb_type>::is_specialized || ((std::numeric_limits<double_limb_type>::max)() - carry >= pr[i + j]));
499          carry += pr[i + j];
500 #ifdef __MSVC_RUNTIME_CHECKS
501          pr[i + j] = static_cast<limb_type>(carry & ~static_cast<limb_type>(0));
502 #else
503          pr[i + j] = static_cast<limb_type>(carry);
504 #endif
505          carry >>= cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::limb_bits;
506          BOOST_MP_ASSERT(carry <= (cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::max_limb_value));
507       }
508       if (carry)
509       {
510          resize_for_carry(result, i + j + 1); // May throw if checking is enabled
511          if (i + j < result.size())
512 #ifdef __MSVC_RUNTIME_CHECKS
513             pr[i + j] = static_cast<limb_type>(carry & ~static_cast<limb_type>(0));
514 #else
515             pr[i + j] = static_cast<limb_type>(carry);
516 #endif
517       }
518       carry = 0;
519    }
520 #endif // ifdef(BOOST_MP_COMBA) ends
521 
522    result.normalize();
523    //
524    // Set the sign of the result:
525    //
526    result.sign(a.sign() != b.sign());
527 }
528 
529 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, unsigned MinBits2, unsigned MaxBits2, cpp_integer_type SignType2, cpp_int_check_type Checked2, class Allocator2>
530 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && !is_trivial_cpp_int<cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const cpp_int_backend<MinBits2,MaxBits2,SignType2,Checked2,Allocator2> & a)531 eval_multiply(
532     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&       result,
533     const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>& a)
534    noexcept((noexcept(eval_multiply(std::declval<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>&>()))))
535 {
536    eval_multiply(result, result, a);
537 }
538 
539 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1>
540 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const limb_type & val)541 eval_multiply(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result, const limb_type& val)
542    noexcept((noexcept(eval_multiply(std::declval<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const limb_type&>()))))
543 {
544    eval_multiply(result, result, val);
545 }
546 
547 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, unsigned MinBits2, unsigned MaxBits2, cpp_integer_type SignType2, cpp_int_check_type Checked2, class Allocator2>
548 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && !is_trivial_cpp_int<cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const cpp_int_backend<MinBits2,MaxBits2,SignType2,Checked2,Allocator2> & a,const double_limb_type & val)549 eval_multiply(
550     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&       result,
551     const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>& a,
552     const double_limb_type&                                                     val)
553    noexcept(
554       (noexcept(eval_multiply(std::declval<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>&>(), std::declval<const limb_type&>())))
555       && (noexcept(eval_multiply(std::declval<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>&>(), std::declval<const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>())))
556    )
557 {
558    if (val <= (std::numeric_limits<limb_type>::max)())
559    {
560       eval_multiply(result, a, static_cast<limb_type>(val));
561    }
562    else
563    {
564 #if BOOST_MP_ENDIAN_LITTLE_BYTE && !defined(BOOST_MP_TEST_NO_LE)
565       cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> t(val);
566 #else
567       cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> t;
568       t = val;
569 #endif
570       eval_multiply(result, a, t);
571    }
572 }
573 
574 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1>
575 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const double_limb_type & val)576 eval_multiply(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result, const double_limb_type& val)
577    noexcept((noexcept(eval_multiply(std::declval<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const double_limb_type&>()))))
578 {
579    eval_multiply(result, result, val);
580 }
581 
582 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, unsigned MinBits2, unsigned MaxBits2, cpp_integer_type SignType2, cpp_int_check_type Checked2, class Allocator2>
583 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && !is_trivial_cpp_int<cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const cpp_int_backend<MinBits2,MaxBits2,SignType2,Checked2,Allocator2> & a,const signed_limb_type & val)584 eval_multiply(
585     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&       result,
586     const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>& a,
587     const signed_limb_type&                                                     val)
588    noexcept((noexcept(eval_multiply(std::declval<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>&>(), std::declval<const limb_type&>()))))
589 {
590    if (val > 0)
591       eval_multiply(result, a, static_cast<limb_type>(val));
592    else
593    {
594       eval_multiply(result, a, static_cast<limb_type>(boost::multiprecision::detail::unsigned_abs(val)));
595       result.negate();
596    }
597 }
598 
599 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1>
600 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const signed_limb_type & val)601 eval_multiply(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result, const signed_limb_type& val)
602    noexcept((noexcept(eval_multiply(std::declval<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const limb_type&>()))))
603 {
604    eval_multiply(result, result, val);
605 }
606 
607 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, unsigned MinBits2, unsigned MaxBits2, cpp_integer_type SignType2, cpp_int_check_type Checked2, class Allocator2>
608 inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && !is_trivial_cpp_int<cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const cpp_int_backend<MinBits2,MaxBits2,SignType2,Checked2,Allocator2> & a,const signed_double_limb_type & val)609 eval_multiply(
610     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&       result,
611     const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>& a,
612     const signed_double_limb_type&                                              val)
613    noexcept(
614    (noexcept(eval_multiply(std::declval<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>&>(), std::declval<const limb_type&>())))
615       && (noexcept(eval_multiply(std::declval<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>&>(), std::declval<const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>())))
616    )
617 {
618    if (val > 0)
619    {
620       if (val <= (std::numeric_limits<limb_type>::max)())
621       {
622          eval_multiply(result, a, static_cast<limb_type>(val));
623          return;
624       }
625    }
626    else if (val >= -static_cast<signed_double_limb_type>((std::numeric_limits<limb_type>::max)()))
627    {
628       eval_multiply(result, a, static_cast<limb_type>(boost::multiprecision::detail::unsigned_abs(val)));
629       result.negate();
630       return;
631    }
632 #if BOOST_MP_ENDIAN_LITTLE_BYTE && !defined(BOOST_MP_TEST_NO_LE)
633    cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> t(val);
634 #else
635    cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> t;
636    t = val;
637 #endif
638    eval_multiply(result, a, t);
639 }
640 
641 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1>
642 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<!is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const signed_double_limb_type & val)643 eval_multiply(cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result, const signed_double_limb_type& val)
644    noexcept(
645    (noexcept(eval_multiply(std::declval<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const limb_type&>())))
646    && (noexcept(eval_multiply(std::declval<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>(), std::declval<const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&>())))
647    )
648 {
649    eval_multiply(result, result, val);
650 }
651 
652 //
653 // Now over again for trivial cpp_int's:
654 //
655 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1>
656 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<
657     is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && (is_signed_number<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value || is_signed_number<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value)>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & o)658 eval_multiply(
659     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&       result,
660     const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& o) noexcept((is_non_throwing_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value))
661 {
662    *result.limbs() = detail::checked_multiply(*result.limbs(), *o.limbs(), typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::checked_type());
663    result.sign(result.sign() != o.sign());
664    result.normalize();
665 }
666 
667 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1>
668 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<
669     is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && is_unsigned_number<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & o)670 eval_multiply(
671     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&       result,
672     const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& o) noexcept((is_non_throwing_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value))
673 {
674    *result.limbs() = detail::checked_multiply(*result.limbs(), *o.limbs(), typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::checked_type());
675    result.normalize();
676 }
677 
678 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1>
679 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<
680     is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && (is_signed_number<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value || is_signed_number<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value)>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & a,const cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & b)681 eval_multiply(
682     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&       result,
683     const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a,
684     const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& b) noexcept((is_non_throwing_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value))
685 {
686    *result.limbs() = detail::checked_multiply(*a.limbs(), *b.limbs(), typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::checked_type());
687    result.sign(a.sign() != b.sign());
688    result.normalize();
689 }
690 
691 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1>
692 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<
693     is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && is_unsigned_number<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,const cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & a,const cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & b)694 eval_multiply(
695     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&       result,
696     const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& a,
697     const cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& b) noexcept((is_non_throwing_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value))
698 {
699    *result.limbs() = detail::checked_multiply(*a.limbs(), *b.limbs(), typename cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>::checked_type());
700    result.normalize();
701 }
702 
703 //
704 // Special routines for multiplying two integers to obtain a multiprecision result:
705 //
706 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1>
707 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<
708     !is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,signed_double_limb_type a,signed_double_limb_type b)709 eval_multiply(
710     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result,
711     signed_double_limb_type a, signed_double_limb_type b)
712 {
713    constexpr const signed_double_limb_type mask = ~static_cast<limb_type>(0);
714    constexpr const unsigned limb_bits = sizeof(limb_type) * CHAR_BIT;
715 
716    bool s = false;
717    if (a < 0)
718    {
719       a = -a;
720       s = true;
721    }
722    if (b < 0)
723    {
724       b = -b;
725       s = !s;
726    }
727    double_limb_type w = a & mask;
728    double_limb_type x = a >> limb_bits;
729    double_limb_type y = b & mask;
730    double_limb_type z = b >> limb_bits;
731 
732    result.resize(4, 4);
733    limb_type* pr = result.limbs();
734 
735    double_limb_type carry = w * y;
736 #ifdef __MSVC_RUNTIME_CHECKS
737    pr[0] = static_cast<limb_type>(carry & ~static_cast<limb_type>(0));
738    carry >>= limb_bits;
739    carry += w * z + x * y;
740    pr[1] = static_cast<limb_type>(carry & ~static_cast<limb_type>(0));
741    carry >>= limb_bits;
742    carry += x * z;
743    pr[2] = static_cast<limb_type>(carry & ~static_cast<limb_type>(0));
744    pr[3] = static_cast<limb_type>(carry >> limb_bits);
745 #else
746    pr[0] = static_cast<limb_type>(carry);
747    carry >>= limb_bits;
748    carry += w * z + x * y;
749    pr[1] = static_cast<limb_type>(carry);
750    carry >>= limb_bits;
751    carry += x * z;
752    pr[2] = static_cast<limb_type>(carry);
753    pr[3] = static_cast<limb_type>(carry >> limb_bits);
754 #endif
755    result.sign(s);
756    result.normalize();
757 }
758 
759 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1>
760 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<
761     !is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,double_limb_type a,double_limb_type b)762 eval_multiply(
763     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result,
764     double_limb_type a, double_limb_type b)
765 {
766    constexpr const signed_double_limb_type mask = ~static_cast<limb_type>(0);
767    constexpr const unsigned limb_bits = sizeof(limb_type) * CHAR_BIT;
768 
769    double_limb_type w = a & mask;
770    double_limb_type x = a >> limb_bits;
771    double_limb_type y = b & mask;
772    double_limb_type z = b >> limb_bits;
773 
774    result.resize(4, 4);
775    limb_type* pr = result.limbs();
776 
777    double_limb_type carry = w * y;
778 #ifdef __MSVC_RUNTIME_CHECKS
779    pr[0] = static_cast<limb_type>(carry & ~static_cast<limb_type>(0));
780    carry >>= limb_bits;
781    carry += w * z;
782    pr[1] = static_cast<limb_type>(carry & ~static_cast<limb_type>(0));
783    carry >>= limb_bits;
784    pr[2] = static_cast<limb_type>(carry & ~static_cast<limb_type>(0));
785    carry = x * y + pr[1];
786    pr[1] = static_cast<limb_type>(carry & ~static_cast<limb_type>(0));
787    carry >>= limb_bits;
788    carry += pr[2] + x * z;
789    pr[2] = static_cast<limb_type>(carry & ~static_cast<limb_type>(0));
790    pr[3] = static_cast<limb_type>(carry >> limb_bits);
791 #else
792    pr[0] = static_cast<limb_type>(carry);
793    carry >>= limb_bits;
794    carry += w * z;
795    pr[1] = static_cast<limb_type>(carry);
796    carry >>= limb_bits;
797    pr[2] = static_cast<limb_type>(carry);
798    carry = x * y + pr[1];
799    pr[1] = static_cast<limb_type>(carry);
800    carry >>= limb_bits;
801    carry += pr[2] + x * z;
802    pr[2] = static_cast<limb_type>(carry);
803    pr[3] = static_cast<limb_type>(carry >> limb_bits);
804 #endif
805    result.sign(false);
806    result.normalize();
807 }
808 
809 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1,
810           unsigned MinBits2, unsigned MaxBits2, cpp_integer_type SignType2, cpp_int_check_type Checked2, class Allocator2>
811 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<
812     !is_trivial_cpp_int<cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::value && is_trivial_cpp_int<cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2> >::value && is_trivial_cpp_int<cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2> >::value>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,cpp_int_backend<MinBits2,MaxBits2,SignType2,Checked2,Allocator2> const & a,cpp_int_backend<MinBits2,MaxBits2,SignType2,Checked2,Allocator2> const & b)813 eval_multiply(
814     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>&       result,
815     cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2> const& a,
816     cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2> const& b)
817 {
818    using canonical_type = typename boost::multiprecision::detail::canonical<typename cpp_int_backend<MinBits2, MaxBits2, SignType2, Checked2, Allocator2>::local_limb_type, cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1> >::type;
819    eval_multiply(result, static_cast<canonical_type>(*a.limbs()), static_cast<canonical_type>(*b.limbs()));
820    result.sign(a.sign() != b.sign());
821 }
822 
823 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, class SI>
824 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<boost::multiprecision::detail::is_signed<SI>::value && boost::multiprecision::detail::is_integral<SI>::value && (sizeof(SI) <= sizeof(signed_double_limb_type) / 2)>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,SI a,SI b)825 eval_multiply(
826     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result,
827     SI a, SI b)
828 {
829    result = static_cast<signed_double_limb_type>(a) * static_cast<signed_double_limb_type>(b);
830 }
831 
832 template <unsigned MinBits1, unsigned MaxBits1, cpp_integer_type SignType1, cpp_int_check_type Checked1, class Allocator1, class UI>
833 BOOST_MP_FORCEINLINE BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<boost::multiprecision::detail::is_unsigned<UI>::value && (sizeof(UI) <= sizeof(signed_double_limb_type) / 2)>::type
eval_multiply(cpp_int_backend<MinBits1,MaxBits1,SignType1,Checked1,Allocator1> & result,UI a,UI b)834 eval_multiply(
835     cpp_int_backend<MinBits1, MaxBits1, SignType1, Checked1, Allocator1>& result,
836     UI a, UI b)
837 {
838    result = static_cast<double_limb_type>(a) * static_cast<double_limb_type>(b);
839 }
840 
841 #ifdef _MSC_VER
842 #pragma warning(pop)
843 #endif
844 
845 }}} // namespace boost::multiprecision::backends
846 
847 #endif
848