1///////////////////////////////////////////////////////////////
2// Copyright 2012 John Maddock. Distributed under the Boost
3// Software License, Version 1.0. (See accompanying file
4// LICENSE_1_0.txt or copy at https://www.boost.org/LICENSE_1_0.txt
5
6#ifndef BOOST_MP_DETAIL_INTEGER_OPS_HPP
7#define BOOST_MP_DETAIL_INTEGER_OPS_HPP
8
9#include <boost/multiprecision/number.hpp>
10#include <boost/multiprecision/detail/no_exceptions_support.hpp>
11
12namespace boost { namespace multiprecision {
13
14namespace default_ops {
15
16template <class Backend>
17inline BOOST_MP_CXX14_CONSTEXPR void eval_qr(const Backend& x, const Backend& y, Backend& q, Backend& r)
18{
19 eval_divide(q, x, y);
20 eval_modulus(r, x, y);
21}
22
23template <class Backend, class Integer>
24inline BOOST_MP_CXX14_CONSTEXPR Integer eval_integer_modulus(const Backend& x, Integer val)
25{
26 BOOST_MP_USING_ABS
27 using default_ops::eval_convert_to;
28 using default_ops::eval_modulus;
29 using int_type = typename boost::multiprecision::detail::canonical<Integer, Backend>::type;
30 Backend t;
31 eval_modulus(t, x, static_cast<int_type>(val));
32 Integer result(0);
33 eval_convert_to(&result, t);
34 return abs(result);
35}
36
37template <class B>
38inline BOOST_MP_CXX14_CONSTEXPR void eval_gcd(B& result, const B& a, const B& b)
39{
40 using default_ops::eval_get_sign;
41 using default_ops::eval_is_zero;
42 using default_ops::eval_lsb;
43
44 std::ptrdiff_t shift(0);
45
46 B u(a), v(b);
47
48 int s = eval_get_sign(u);
49
50 /* GCD(0,x) := x */
51 if (s < 0)
52 {
53 u.negate();
54 }
55 else if (s == 0)
56 {
57 result = v;
58 return;
59 }
60 s = eval_get_sign(v);
61 if (s < 0)
62 {
63 v.negate();
64 }
65 else if (s == 0)
66 {
67 result = u;
68 return;
69 }
70
71 /* Let shift := lg K, where K is the greatest power of 2
72 dividing both u and v. */
73
74 std::size_t us = eval_lsb(u);
75 std::size_t vs = eval_lsb(v);
76 shift = static_cast<std::ptrdiff_t>((std::min)(a: us, b: vs));
77 eval_right_shift(u, us);
78 eval_right_shift(v, vs);
79
80 do
81 {
82 /* Now u and v are both odd, so diff(u, v) is even.
83 Let u = min(u, v), v = diff(u, v)/2. */
84 s = u.compare(v);
85 if (s > 0)
86 u.swap(v);
87 if (s == 0)
88 break;
89 eval_subtract(v, u);
90 vs = eval_lsb(v);
91 eval_right_shift(v, vs);
92 } while (true);
93
94 result = u;
95 eval_left_shift(result, shift);
96}
97
98template <class B>
99inline BOOST_MP_CXX14_CONSTEXPR void eval_lcm(B& result, const B& a, const B& b)
100{
101 using ui_type = typename std::tuple_element<0, typename B::unsigned_types>::type;
102 B t;
103 eval_gcd(t, a, b);
104
105 if (eval_is_zero(t))
106 {
107 result = static_cast<ui_type>(0);
108 }
109 else
110 {
111 eval_divide(result, a, t);
112 eval_multiply(result, b);
113 }
114 if (eval_get_sign(result) < 0)
115 result.negate();
116}
117
118} // namespace default_ops
119
120template <class Backend, expression_template_option ExpressionTemplates>
121inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<number_category<Backend>::value == number_kind_integer>::type
122divide_qr(const number<Backend, ExpressionTemplates>& x, const number<Backend, ExpressionTemplates>& y,
123 number<Backend, ExpressionTemplates>& q, number<Backend, ExpressionTemplates>& r)
124{
125 using default_ops::eval_qr;
126 eval_qr(x.backend(), y.backend(), q.backend(), r.backend());
127}
128
129template <class Backend, expression_template_option ExpressionTemplates, class tag, class A1, class A2, class A3, class A4>
130inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<number_category<Backend>::value == number_kind_integer>::type
131divide_qr(const number<Backend, ExpressionTemplates>& x, const multiprecision::detail::expression<tag, A1, A2, A3, A4>& y,
132 number<Backend, ExpressionTemplates>& q, number<Backend, ExpressionTemplates>& r)
133{
134 divide_qr(x, number<Backend, ExpressionTemplates>(y), q, r);
135}
136
137template <class tag, class A1, class A2, class A3, class A4, class Backend, expression_template_option ExpressionTemplates>
138inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<number_category<Backend>::value == number_kind_integer>::type
139divide_qr(const multiprecision::detail::expression<tag, A1, A2, A3, A4>& x, const number<Backend, ExpressionTemplates>& y,
140 number<Backend, ExpressionTemplates>& q, number<Backend, ExpressionTemplates>& r)
141{
142 divide_qr(number<Backend, ExpressionTemplates>(x), y, q, r);
143}
144
145template <class tag, class A1, class A2, class A3, class A4, class tagb, class A1b, class A2b, class A3b, class A4b, class Backend, expression_template_option ExpressionTemplates>
146inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<number_category<Backend>::value == number_kind_integer>::type
147divide_qr(const multiprecision::detail::expression<tag, A1, A2, A3, A4>& x, const multiprecision::detail::expression<tagb, A1b, A2b, A3b, A4b>& y,
148 number<Backend, ExpressionTemplates>& q, number<Backend, ExpressionTemplates>& r)
149{
150 divide_qr(number<Backend, ExpressionTemplates>(x), number<Backend, ExpressionTemplates>(y), q, r);
151}
152
153template <class Backend, expression_template_option ExpressionTemplates, class Integer>
154inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<boost::multiprecision::detail::is_integral<Integer>::value && (number_category<Backend>::value == number_kind_integer), Integer>::type
155integer_modulus(const number<Backend, ExpressionTemplates>& x, Integer val)
156{
157 using default_ops::eval_integer_modulus;
158 return eval_integer_modulus(x.backend(), val);
159}
160
161template <class tag, class A1, class A2, class A3, class A4, class Integer>
162inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<boost::multiprecision::detail::is_integral<Integer>::value && (number_category<typename multiprecision::detail::expression<tag, A1, A2, A3, A4>::result_type>::value == number_kind_integer), Integer>::type
163integer_modulus(const multiprecision::detail::expression<tag, A1, A2, A3, A4>& x, Integer val)
164{
165 using result_type = typename multiprecision::detail::expression<tag, A1, A2, A3, A4>::result_type;
166 return integer_modulus(result_type(x), val);
167}
168
169template <class Backend, expression_template_option ExpressionTemplates>
170inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<number_category<Backend>::value == number_kind_integer, std::size_t>::type
171lsb(const number<Backend, ExpressionTemplates>& x)
172{
173 using default_ops::eval_lsb;
174 return eval_lsb(x.backend());
175}
176
177template <class tag, class A1, class A2, class A3, class A4>
178inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<number_category<typename multiprecision::detail::expression<tag, A1, A2, A3, A4>::result_type>::value == number_kind_integer, std::size_t>::type
179lsb(const multiprecision::detail::expression<tag, A1, A2, A3, A4>& x)
180{
181 using number_type = typename multiprecision::detail::expression<tag, A1, A2, A3, A4>::result_type;
182 number_type n(x);
183 using default_ops::eval_lsb;
184 return eval_lsb(n.backend());
185}
186
187template <class Backend, expression_template_option ExpressionTemplates>
188inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<number_category<Backend>::value == number_kind_integer, std::size_t>::type
189msb(const number<Backend, ExpressionTemplates>& x)
190{
191 using default_ops::eval_msb;
192 return eval_msb(x.backend());
193}
194
195template <class tag, class A1, class A2, class A3, class A4>
196inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<number_category<typename multiprecision::detail::expression<tag, A1, A2, A3, A4>::result_type>::value == number_kind_integer, std::size_t>::type
197msb(const multiprecision::detail::expression<tag, A1, A2, A3, A4>& x)
198{
199 using number_type = typename multiprecision::detail::expression<tag, A1, A2, A3, A4>::result_type;
200 number_type n(x);
201 using default_ops::eval_msb;
202 return eval_msb(n.backend());
203}
204
205template <class Backend, expression_template_option ExpressionTemplates>
206inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<number_category<Backend>::value == number_kind_integer, bool>::type
207bit_test(const number<Backend, ExpressionTemplates>& x, std::size_t index)
208{
209 using default_ops::eval_bit_test;
210 return eval_bit_test(x.backend(), index);
211}
212
213template <class tag, class A1, class A2, class A3, class A4>
214inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<number_category<typename multiprecision::detail::expression<tag, A1, A2, A3, A4>::result_type>::value == number_kind_integer, bool>::type
215bit_test(const multiprecision::detail::expression<tag, A1, A2, A3, A4>& x, std::size_t index)
216{
217 using number_type = typename multiprecision::detail::expression<tag, A1, A2, A3, A4>::result_type;
218 number_type n(x);
219 using default_ops::eval_bit_test;
220 return eval_bit_test(n.backend(), index);
221}
222
223template <class Backend, expression_template_option ExpressionTemplates>
224inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<number_category<Backend>::value == number_kind_integer, number<Backend, ExpressionTemplates>&>::type
225bit_set(number<Backend, ExpressionTemplates>& x, std::size_t index)
226{
227 using default_ops::eval_bit_set;
228 eval_bit_set(x.backend(), index);
229 return x;
230}
231
232template <class Backend, expression_template_option ExpressionTemplates>
233inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<number_category<Backend>::value == number_kind_integer, number<Backend, ExpressionTemplates>&>::type
234bit_unset(number<Backend, ExpressionTemplates>& x, std::size_t index)
235{
236 using default_ops::eval_bit_unset;
237 eval_bit_unset(x.backend(), index);
238 return x;
239}
240
241template <class Backend, expression_template_option ExpressionTemplates>
242inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<number_category<Backend>::value == number_kind_integer, number<Backend, ExpressionTemplates>&>::type
243bit_flip(number<Backend, ExpressionTemplates>& x, std::size_t index)
244{
245 using default_ops::eval_bit_flip;
246 eval_bit_flip(x.backend(), index);
247 return x;
248}
249
250namespace default_ops {
251
252//
253// Within powm, we need a type with twice as many digits as the argument type, define
254// a traits class to obtain that type:
255//
256template <class Backend>
257struct double_precision_type
258{
259 using type = Backend;
260};
261
262//
263// If the exponent is a signed integer type, then we need to
264// check the value is positive:
265//
266template <class Backend>
267inline BOOST_MP_CXX14_CONSTEXPR void check_sign_of_backend(const Backend& v, const std::integral_constant<bool, true>)
268{
269 if (eval_get_sign(v) < 0)
270 {
271 BOOST_MP_THROW_EXCEPTION(std::runtime_error("powm requires a positive exponent."));
272 }
273}
274template <class Backend>
275inline BOOST_MP_CXX14_CONSTEXPR void check_sign_of_backend(const Backend&, const std::integral_constant<bool, false>) {}
276//
277// Calculate (a^p)%c:
278//
279template <class Backend>
280BOOST_MP_CXX14_CONSTEXPR void eval_powm(Backend& result, const Backend& a, const Backend& p, const Backend& c)
281{
282 using default_ops::eval_bit_test;
283 using default_ops::eval_get_sign;
284 using default_ops::eval_modulus;
285 using default_ops::eval_multiply;
286 using default_ops::eval_right_shift;
287
288 using double_type = typename double_precision_type<Backend>::type ;
289 using ui_type = typename boost::multiprecision::detail::canonical<unsigned char, double_type>::type;
290
291 check_sign_of_backend(p, std::integral_constant<bool, std::numeric_limits<number<Backend> >::is_signed>());
292
293 double_type x, y(a), b(p), t;
294 x = ui_type(1u);
295
296 while (eval_get_sign(b) > 0)
297 {
298 if (eval_bit_test(b, 0))
299 {
300 eval_multiply(t, x, y);
301 eval_modulus(x, t, c);
302 }
303 eval_multiply(t, y, y);
304 eval_modulus(y, t, c);
305 eval_right_shift(b, ui_type(1));
306 }
307 Backend x2(x);
308 eval_modulus(result, x2, c);
309}
310
311template <class Backend, class Integer>
312BOOST_MP_CXX14_CONSTEXPR void eval_powm(Backend& result, const Backend& a, const Backend& p, Integer c)
313{
314 using double_type = typename double_precision_type<Backend>::type ;
315 using ui_type = typename boost::multiprecision::detail::canonical<unsigned char, double_type>::type;
316 using i1_type = typename boost::multiprecision::detail::canonical<Integer, double_type>::type ;
317 using i2_type = typename boost::multiprecision::detail::canonical<Integer, Backend>::type ;
318
319 using default_ops::eval_bit_test;
320 using default_ops::eval_get_sign;
321 using default_ops::eval_modulus;
322 using default_ops::eval_multiply;
323 using default_ops::eval_right_shift;
324
325 check_sign_of_backend(p, std::integral_constant<bool, std::numeric_limits<number<Backend> >::is_signed>());
326
327 if (eval_get_sign(p) < 0)
328 {
329 BOOST_MP_THROW_EXCEPTION(std::runtime_error("powm requires a positive exponent."));
330 }
331
332 double_type x, y(a), b(p), t;
333 x = ui_type(1u);
334
335 while (eval_get_sign(b) > 0)
336 {
337 if (eval_bit_test(b, 0))
338 {
339 eval_multiply(t, x, y);
340 eval_modulus(x, t, static_cast<i1_type>(c));
341 }
342 eval_multiply(t, y, y);
343 eval_modulus(y, t, static_cast<i1_type>(c));
344 eval_right_shift(b, ui_type(1));
345 }
346 Backend x2(x);
347 eval_modulus(result, x2, static_cast<i2_type>(c));
348}
349
350template <class Backend, class Integer>
351BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<boost::multiprecision::detail::is_unsigned<Integer>::value >::type eval_powm(Backend& result, const Backend& a, Integer b, const Backend& c)
352{
353 using double_type = typename double_precision_type<Backend>::type ;
354 using ui_type = typename boost::multiprecision::detail::canonical<unsigned char, double_type>::type;
355
356 using default_ops::eval_bit_test;
357 using default_ops::eval_get_sign;
358 using default_ops::eval_modulus;
359 using default_ops::eval_multiply;
360 using default_ops::eval_right_shift;
361
362 double_type x, y(a), t;
363 x = ui_type(1u);
364
365 while (b > 0)
366 {
367 if (b & 1)
368 {
369 eval_multiply(t, x, y);
370 eval_modulus(x, t, c);
371 }
372 eval_multiply(t, y, y);
373 eval_modulus(y, t, c);
374 b >>= 1;
375 }
376 Backend x2(x);
377 eval_modulus(result, x2, c);
378}
379
380template <class Backend, class Integer>
381BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<boost::multiprecision::detail::is_signed<Integer>::value && boost::multiprecision::detail::is_integral<Integer>::value>::type eval_powm(Backend& result, const Backend& a, Integer b, const Backend& c)
382{
383 if (b < 0)
384 {
385 BOOST_MP_THROW_EXCEPTION(std::runtime_error("powm requires a positive exponent."));
386 }
387 eval_powm(result, a, static_cast<typename boost::multiprecision::detail::make_unsigned<Integer>::type>(b), c);
388}
389
390template <class Backend, class Integer1, class Integer2>
391BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<boost::multiprecision::detail::is_unsigned<Integer1>::value >::type eval_powm(Backend& result, const Backend& a, Integer1 b, Integer2 c)
392{
393 using double_type = typename double_precision_type<Backend>::type ;
394 using ui_type = typename boost::multiprecision::detail::canonical<unsigned char, double_type>::type;
395 using i1_type = typename boost::multiprecision::detail::canonical<Integer1, double_type>::type ;
396 using i2_type = typename boost::multiprecision::detail::canonical<Integer2, Backend>::type ;
397
398 using default_ops::eval_bit_test;
399 using default_ops::eval_get_sign;
400 using default_ops::eval_modulus;
401 using default_ops::eval_multiply;
402 using default_ops::eval_right_shift;
403
404 double_type x, y(a), t;
405 x = ui_type(1u);
406
407 while (b > 0)
408 {
409 if (b & 1)
410 {
411 eval_multiply(t, x, y);
412 eval_modulus(x, t, static_cast<i1_type>(c));
413 }
414 eval_multiply(t, y, y);
415 eval_modulus(y, t, static_cast<i1_type>(c));
416 b >>= 1;
417 }
418 Backend x2(x);
419 eval_modulus(result, x2, static_cast<i2_type>(c));
420}
421
422template <class Backend, class Integer1, class Integer2>
423BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<boost::multiprecision::detail::is_signed<Integer1>::value && boost::multiprecision::detail::is_integral<Integer1>::value>::type eval_powm(Backend& result, const Backend& a, Integer1 b, Integer2 c)
424{
425 if (b < 0)
426 {
427 BOOST_MP_THROW_EXCEPTION(std::runtime_error("powm requires a positive exponent."));
428 }
429 eval_powm(result, a, static_cast<typename boost::multiprecision::detail::make_unsigned<Integer1>::type>(b), c);
430}
431
432struct powm_func
433{
434 template <class T, class U, class V>
435 BOOST_MP_CXX14_CONSTEXPR void operator()(T& result, const T& b, const U& p, const V& m) const
436 {
437 eval_powm(result, b, p, m);
438 }
439 template <class R, class T, class U, class V>
440 BOOST_MP_CXX14_CONSTEXPR void operator()(R& result, const T& b, const U& p, const V& m) const
441 {
442 T temp;
443 eval_powm(temp, b, p, m);
444 result = std::move(temp);
445 }
446};
447
448} // namespace default_ops
449
450template <class T, class U, class V>
451inline BOOST_MP_CXX14_CONSTEXPR typename std::enable_if<
452 (number_category<T>::value == number_kind_integer) &&
453 (is_number<T>::value || is_number_expression<T>::value) &&
454 (is_number<U>::value || is_number_expression<U>::value || boost::multiprecision::detail::is_integral<U>::value) &&
455 (is_number<V>::value || is_number_expression<V>::value || boost::multiprecision::detail::is_integral<V>::value),
456 typename std::conditional<
457 is_no_et_number<T>::value,
458 T,
459 typename std::conditional<
460 is_no_et_number<U>::value,
461 U,
462 typename std::conditional<
463 is_no_et_number<V>::value,
464 V,
465 detail::expression<detail::function, default_ops::powm_func, T, U, V> >::type>::type>::type>::type
466powm(const T& b, const U& p, const V& mod)
467{
468 return detail::expression<detail::function, default_ops::powm_func, T, U, V>(
469 default_ops::powm_func(), b, p, mod);
470}
471
472}} // namespace boost::multiprecision
473
474#endif
475

source code of boost/libs/multiprecision/include/boost/multiprecision/detail/integer_ops.hpp