Direct-BT v3.3.0-1-gc2d430c
Direct-BT - Direct Bluetooth Programming.
big_int.hpp
Go to the documentation of this file.
1/**
2 * big_int_t (this file)
3 * (c) 2024 Gothel Software e.K.
4 *
5 * Includes code from: Botan BigInt (bigint.h, ..)
6 * (C) 1999-2011,2012,2014,2019 Jack Lloyd
7 * 2007 FlexSecure
8 *
9 * Includes code from: Botan Division Algorithms (divide.h)
10 * (C) 1999-2007,2012,2018,2021 Jack Lloyd
11 *
12 * jaulib including this code is released under the MIT License (see COPYING)
13 * Botan itself is released under the Simplified BSD License (see COPYING)
14 */
15
16#ifndef JAU_BIG_INT_HPP_
17#define JAU_BIG_INT_HPP_
18
19#include <cstdint>
20#include <cassert>
21
23#include <jau/byte_util.hpp>
24#include <jau/string_util.hpp>
25
26namespace jau::mp {
27
28 /** \addtogroup Integer
29 *
30 * @{
31 */
32
33 /**
34 * Arbitrary precision integer type
35 *
36 * @anchor bigint_storage_format
37 * ### Local storage format
38 * Internally the big integer is stored in an array of mp_word_t ordered little-endian alike,
39 * with the least significant word at the array-bottom and most significant word at the array-top.
40 *
41 * The mp_word_t itself is stored in jau::endian::native!
42 */
43 class BigInt {
44 public:
45 /**
46 * Sign symbol definitions for positive and negative numbers
47 */
48 enum sign_t { negative = 0, positive = 1 };
49
50 public:
51 BigInt() noexcept = default;
52
53 BigInt(const BigInt& o) noexcept = default;
54
55 ~BigInt() noexcept = default;
56
57 /**
58 * Create a 0-value big_int
59 */
60 static BigInt zero() { return BigInt(); }
61
62 /**
63 * Create a 1-value big_int
64 */
65 static BigInt one() { return BigInt::from_word(1); }
66
67 /**
68 * Create big_int from an unsigned 64 bit integer
69 * @param n initial value of this big_int
70 */
71 static BigInt from_u64(uint64_t n) {
72 BigInt bn;
73 if( 64 == mp_word_bits ) {
74 bn.set_word_at(0, n);
75 } else {
76 bn.set_word_at(1, static_cast<mp_word_t>(n >> 32));
77 bn.set_word_at(0, static_cast<mp_word_t>(n));
78 }
79 return bn;
80 }
81
82 /**
83 * Create big_int from a mp_word_t (limb)
84 * @param n initial value of this big_int
85 */
87 BigInt bn;
88 bn.set_word_at(0, n);
89 return bn;
90 }
91
92 /**
93 * Create big_int from a signed 32 bit integer
94 * @param n initial value of this big_int
95 */
96 static BigInt from_s32(int32_t n) {
97 if(n >= 0) {
98 return BigInt::from_u64(static_cast<uint64_t>(n));
99 } else {
100 return -BigInt::from_u64(static_cast<uint64_t>(-n));
101 }
102 }
103
104 /**
105 * Create big_int of specified size, all zeros
106 * @param n size of the internal register in words
107 */
108 static BigInt with_capacity(size_t n) {
109 BigInt bn;
110 bn.grow_to(n);
111 return bn;
112 }
113
114 /**
115 * Create a power of two
116 * @param n the power of two to create
117 * @return big_int_t representing 2^n
118 */
119 static BigInt power_of_2(size_t n) {
120 BigInt b;
121 b.set_bit(n);
122 return b;
123 }
124
125 /**
126 * Create big_int_t from an unsigned 64 bit integer
127 * @param n initial value of this big_int_t
128 */
129 BigInt(uint64_t n) {
130 if( 64 == mp_word_bits ) {
131 m_data.set_word_at(0, n);
132 } else {
133 m_data.set_word_at(1, static_cast<mp_word_t>(n >> 32));
134 m_data.set_word_at(0, static_cast<mp_word_t>(n));
135 }
136 }
137
138 /**
139 * Construct a big_int_t from a string encoded as hexadecimal or decimal.
140 *
141 * Both number bases may lead a `-`, denoting a negative number.
142 *
143 * Hexadecimal is detected by a leading `0x`.
144 */
145 BigInt(const std::string& str) {
146 size_t markers = 0;
147 bool is_negative = false;
148
149 if(str.length() > 0 && str[0] == '-')
150 {
151 markers += 1;
152 is_negative = true;
153 }
154
155 if(str.length() > markers + 2 && str[markers ] == '0' &&
156 str[markers + 1] == 'x')
157 {
158 markers += 2;
159 *this = hex_decode(cast_char_ptr_to_uint8(str.data()) + markers, str.length() - markers, lb_endian_t::big);
160 } else {
161 *this = dec_decode(cast_char_ptr_to_uint8(str.data()) + markers, str.length() - markers);
162 }
163
165 else set_sign(positive);
166 }
167
168 /**
169 * Create a big_int_t from an integer in a byte array with given byte_len,
170 * considering the given byte_order.
171 *
172 * The value is stored in the local storage format, see \ref bigint_storage_format
173 *
174 * @param buf the byte array holding the value
175 * @param byte_len size of buf in bytes
176 * @param littleEndian
177 */
178 BigInt(const uint8_t buf[], size_t byte_len, const lb_endian_t byte_order) {
179 binary_decode(buf, byte_len, byte_order);
180 }
181
182 BigInt(std::vector<mp_word_t>&& other_reg) noexcept {
183 this->swap_reg(other_reg);
184 }
185
186 BigInt(BigInt&& other) noexcept {
187 this->swap(other);
188 }
189
190 BigInt& operator=(const BigInt& r) = default;
191
192 BigInt& operator=(BigInt&& other) noexcept {
193 if(this != &other) {
194 this->swap(other);
195 }
196 return *this;
197 }
198
199 /**
200 * Swap this value with another
201 * @param other big_int to swap values with
202 */
203 void swap(BigInt& other) noexcept {
204 m_data.swap(other.m_data);
205 std::swap(m_signedness, other.m_signedness);
206 }
207
208 private:
209 void swap_reg(std::vector<mp_word_t>& reg) noexcept {
210 m_data.swap(reg);
211 // sign left unchanged
212 }
213
214 public:
215 /** Unary negation operator, returns new negative instance of this. */
216 BigInt operator-() const noexcept { return BigInt(*this).flip_sign(); }
217
218 /**
219 * @param n the offset to get a byte from
220 * @result byte at offset n
221 */
222 uint8_t byte_at(size_t n) const noexcept {
223 return get_byte_var_be(sizeof(mp_word_t) - (n % sizeof(mp_word_t)) - 1,
224 word_at(n / sizeof(mp_word_t)));
225 }
226
227 /**
228 * Return the mp_word_t at a specified position of the internal register
229 * @param n position in the register
230 * @return value at position n
231 */
232 mp_word_t word_at(size_t n) const noexcept {
233 return m_data.get_word_at(n);
234 }
235
236 /**
237 * Return a const pointer to the register
238 * @result a pointer to the start of the internal register
239 */
240 const mp_word_t* data() const { return m_data.const_data(); }
241
242 /**
243 * Tests if the sign of the integer is negative
244 * @result true, iff the integer has a negative sign
245 */
246 bool is_negative() const noexcept { return sign() == negative; }
247
248 /**
249 * Tests if the sign of the integer is positive
250 * @result true, iff the integer has a positive sign
251 */
252 bool is_positive() const noexcept { return sign() == positive; }
253
254 /**
255 * Return the sign of the integer
256 * @result the sign of the integer
257 */
258 sign_t sign() const noexcept { return m_signedness; }
259
260 /**
261 * @result the opposite sign of the represented integer value
262 */
263 sign_t reverse_sign() const noexcept {
264 if(sign() == positive) {
265 return negative;
266 }
267 return positive;
268 }
269
270 /**
271 * Flip the sign of this big_int
272 */
273 BigInt& flip_sign() noexcept {
274 return set_sign(reverse_sign());
275 }
276
277 /**
278 * Set sign of the integer
279 * @param sign new Sign to set
280 */
282 if(sign == negative && is_zero()) {
283 sign = positive;
284 }
285 m_signedness = sign;
286 return *this;
287 }
288
289 /** Returns absolute (positive) value of this instance */
290 BigInt abs() const noexcept {
291 return BigInt(*this).set_sign(positive);
292 }
293
294 /**
295 * Give size of internal register
296 * @result size of internal register in words
297 */
298 size_t size() const noexcept { return m_data.size(); }
299
300 /**
301 * Return how many words we need to hold this value
302 * @result significant words of the represented integer value
303 */
304 size_t sig_words() const noexcept { return m_data.sig_words(); }
305
306 /** Returns byte length of this integer */
307 size_t bytes() const { return jau::round_up(bits(), 8U) / 8U; }
308
309 /** Returns bit length of this integer */
310 size_t bits() const noexcept {
311 const size_t words = sig_words();
312 if(words == 0) {
313 return 0;
314 }
315 const size_t full_words = (words - 1) * mp_word_bits;
316 const size_t top_bits = mp_word_bits - top_bits_free();
317 return full_words + top_bits;
318 }
319
320 /**
321 * Zeroize the big_int. The size of the underlying register is not
322 * modified.
323 */
324 void clear() { m_data.set_to_zero(); m_signedness = positive; }
325
326 /**
327 * Compares this instance against other considering both sign() value
328 * Returns
329 * - -1 if this < other
330 * - 0 if this == other
331 * - 1 if this > other
332 */
333 int compare(const BigInt& b) const noexcept {
334 return cmp(b, true);
335 }
336
337 /**
338 * Test if the integer has an even value
339 * @result true if the integer is even, false otherwise
340 */
341 bool is_even() const noexcept { return get_bit(0) == 0; }
342
343 /**
344 * Test if the integer has an odd value
345 * @result true if the integer is odd, false otherwise
346 */
347 bool is_odd() const noexcept { return get_bit(0) == 1; }
348
349 /**
350 * Test if the integer is not zero
351 * @result true if the integer is non-zero, false otherwise
352 */
353 bool is_nonzero() const noexcept { return !is_zero(); }
354
355 /**
356 * Test if the integer is zero
357 * @result true if the integer is zero, false otherwise
358 */
359 bool is_zero() const noexcept {
360 return sig_words() == 0;
361 }
362
363 /**
364 * Return bit value at specified position
365 * @param n the bit offset to test
366 * @result true, if the bit at position n is set, false otherwise
367 */
368 bool get_bit(size_t n) const noexcept {
369 return (word_at(n / mp_word_bits) >> (n % mp_word_bits)) & 1;
370 }
371
372 /**
373 * Set bit at specified position
374 * @param n bit position to set
375 */
376 void set_bit(size_t n) noexcept {
377 conditionally_set_bit(n, true);
378 }
379
380 /**
381 * Conditionally set bit at specified position. Note if set_it is
382 * false, nothing happens, and if the bit is already set, it
383 * remains set.
384 *
385 * @param n bit position to set
386 * @param set_it if the bit should be set
387 */
388 void conditionally_set_bit(size_t n, bool set_it) noexcept {
389 const size_t which = n / mp_word_bits;
390 const mp_word_t mask = static_cast<mp_word_t>(set_it) << (n % mp_word_bits);
391 m_data.set_word_at(which, word_at(which) | mask);
392 }
393
394 /** ! operator, returns `true` iff this is zero, otherwise false. */
395 bool operator !() const noexcept { return is_zero(); }
396
397 bool operator==(const BigInt& b) const noexcept { return is_equal(b); }
398 bool operator!=(const BigInt& b) const noexcept { return !is_equal(b); }
399 bool operator<=(const BigInt& b) const noexcept { return cmp(b) <= 0; }
400 bool operator>=(const BigInt& b) const noexcept { return cmp(b) >= 0; }
401 bool operator<(const BigInt& b) const noexcept { return is_less_than(b); }
402 bool operator>(const BigInt& b) const noexcept { return b.is_less_than(*this); }
403#if 0
404 std::strong_ordering operator<=>(const big_int_t& b) const noexcept {
405 const int r = cmp(b);
406 return 0 == r ? std::strong_ordering::equal : ( 0 > r ? std::strong_ordering::less : std::strong_ordering::greater);
407 }
408#endif
409
410 BigInt& operator++() noexcept { return *this += 1; }
411
412 BigInt& operator--() noexcept { return *this -= 1; }
413
414 BigInt operator++(int) noexcept { BigInt x = (*this); ++(*this); return x; }
415
416 BigInt operator--(int) noexcept { BigInt x = (*this); --(*this); return x; }
417
418 BigInt& operator+=(const BigInt& y ) noexcept {
419 return add(y.data(), y.sig_words(), y.sign());
420 }
421
422 BigInt& operator-=(const BigInt& y ) noexcept {
423 return add(y.data(), y.sig_words(), y.sign() == positive ? negative : positive);
424 }
425
426 BigInt operator+(const BigInt& y ) const noexcept {
427 return add2(*this, y.data(), y.sig_words(), y.sign());
428 }
429
430 BigInt operator-(const BigInt& y ) const noexcept {
431 return add2(*this, y.data(), y.sig_words(), y.reverse_sign());
432 }
433
434 BigInt& operator<<=(size_t shift) noexcept {
435 const size_t shift_words = shift / mp_word_bits;
436 const size_t shift_bits = shift % mp_word_bits;
437 const size_t size = sig_words();
438
439 const size_t bits_free = top_bits_free();
440
441 const size_t new_size = size + shift_words + (bits_free < shift_bits);
442
443 m_data.grow_to(new_size);
444
445 ops::bigint_shl1(m_data.mutable_data(), new_size, size, shift_words, shift_bits);
446
447 return *this;
448 }
449
450 BigInt& operator>>=(size_t shift) noexcept {
451 const size_t shift_words = shift / mp_word_bits;
452 const size_t shift_bits = shift % mp_word_bits;
453
454 ops::bigint_shr1(m_data.mutable_data(), m_data.size(), shift_words, shift_bits);
455
456 if(is_negative() && is_zero()) {
458 }
459 return *this;
460 }
461
462 BigInt operator<<(size_t shift) const {
463 const size_t shift_words = shift / mp_word_bits;
464 const size_t shift_bits = shift % mp_word_bits;
465 const size_t x_sw = sig_words();
466
467 BigInt y = BigInt::with_capacity(x_sw + shift_words + (shift_bits ? 1 : 0));
468 ops::bigint_shl2(y.mutable_data(), data(), x_sw, shift_words, shift_bits);
469 y.set_sign(sign());
470 return y;
471 }
472
473 BigInt operator>>(size_t shift) const {
474 const size_t shift_words = shift / mp_word_bits;
475 const size_t shift_bits = shift % mp_word_bits;
476 const size_t x_sw = sig_words();
477
478 if(shift_words >= x_sw) {
479 return BigInt::zero();
480 }
481
482 BigInt y = BigInt::with_capacity(x_sw - shift_words);
483 ops::bigint_shr2(y.mutable_data(), data(), x_sw, shift_words, shift_bits);
484
485 if(is_negative() && y.is_zero()) {
487 } else {
488 y.set_sign(sign());
489 }
490 return y;
491 }
492
493
494 BigInt& operator*=(const BigInt& y) noexcept {
495 std::vector<mp_word_t> ws;
496 return this->mul(y, ws);
497 }
498
499 BigInt operator*(const BigInt& y) noexcept
500 {
501 const size_t x_sw = sig_words();
502 const size_t y_sw = y.sig_words();
503
504 BigInt z;
505 z.resize(size() + y.size());
506
507 if(x_sw == 1 && y_sw) {
508 ops::bigint_linmul3(z.mutable_data(), y.data(), y_sw, word_at(0));
509 } else if(y_sw == 1 && x_sw) {
510 ops::bigint_linmul3(z.mutable_data(), data(), x_sw, y.word_at(0));
511 } else if(x_sw && y_sw) {
512 ops::basecase_mul(z.mutable_data(), z.size(), data(), x_sw, y.data(), y_sw);
513 }
514 z.cond_flip_sign(x_sw > 0 && y_sw > 0 && sign() != y.sign());
515 return z;
516 }
517
519 if(y.sig_words() == 1 && jau::is_power_of_2(y.word_at(0))) {
520 (*this) >>= (y.bits() - 1);
521 } else {
522 (*this) = (*this) / y;
523 }
524 return (*this);
525 }
526
527 BigInt operator/(const BigInt& y) const {
528 if(y.sig_words() == 1) {
529 return *this / y.word_at(0);
530 }
531 BigInt q, r;
532 vartime_divide(*this, y, q, r);
533 return q;
534 }
535 /**
536 * Modulo operator
537 * @param y the modulus to reduce this by
538 */
539 BigInt& operator%=(const BigInt& mod) {
540 return (*this = (*this) % mod);
541 }
542
543 BigInt operator%(const BigInt& mod) {
544 if(mod.is_zero()) {
546 }
547 if(mod.is_negative()) {
548 throw jau::math::MathDomainError("mod < 0", E_FILE_LINE);
549 }
550 if(is_positive() && mod.is_positive() && *this < mod) {
551 return *this;
552 }
553 if(mod.sig_words() == 1) {
554 return from_word(*this % mod.word_at(0));
555 }
556 BigInt q, r;
557 vartime_divide(*this, mod, q, r);
558 return r;
559 }
560
561 /**
562 * Returns (*this)^e, or pow(*this, e)
563 *
564 * Implementation is not optimized and naive, i.e. O(n)
565 *
566 * @param e the exponent
567 */
569 const BigInt& b = *this;
570 if( b.is_zero() ) {
571 return BigInt::zero();
572 }
573 const BigInt one_v = BigInt::one();
574 BigInt r = one_v;
575 bool is_negative;
576 if( e.is_negative() ) {
577 is_negative = true;
578 e.flip_sign();
579 } else {
580 is_negative = false;
581 }
582
583 while( e.is_nonzero() ) {
584 r *= b;
585 --e;
586 }
587
588 if( is_negative ) {
589 return one_v / r;
590 } else {
591 return r;
592 }
593 }
594
595 /**
596 * Returns (*this)^e % m, or pow(*this, e) % m
597 *
598 * Implementation is not optimized and naive, i.e. O(n)
599 *
600 * @param e the exponent
601 */
603 const BigInt& b = *this;
604 if( b.is_zero() ) {
605 return BigInt::zero();
606 }
607 const BigInt one_v = BigInt::one();
608 BigInt r = one_v;
609 bool is_negative;
610 if( e.is_negative() ) {
611 is_negative = true;
612 e.flip_sign();
613 } else {
614 is_negative = false;
615 }
616
617 while( e.is_nonzero() ) {
618 r *= b;
619 r %= m;
620 --e;
621 }
622
623 if( is_negative ) {
624 return one_v / r;
625 } else {
626 return r;
627 }
628 }
629
630 /**
631 * Square value of *this
632 * @param ws a temp workspace
633 */
634 BigInt& square(std::vector<mp_word_t>& ws); // TODO
635
636 /**
637 * Set *this to y - *this
638 * @param y the big_int_t to subtract from
639 * @param ws a temp workspace
640 */
641 BigInt& rev_sub(const BigInt& y, std::vector<mp_word_t>& ws); // TODO
642
643 /**
644 * Set *this to (*this + y) % mod
645 * This function assumes *this is >= 0 && < mod
646 * @param y the big_int_t to add - assumed y >= 0 and y < mod
647 * @param mod the positive modulus
648 * @param ws a temp workspace
649 */
650 BigInt& mod_add(const BigInt& y, const BigInt& mod, std::vector<mp_word_t>& ws); // TODO
651
652 /**
653 * Set *this to (*this - y) % mod
654 * This function assumes *this is >= 0 && < mod
655 * @param y the big_int_t to subtract - assumed y >= 0 and y < mod
656 * @param mod the positive modulus
657 * @param ws a temp workspace
658 */
659 BigInt& mod_sub(const BigInt& y, const BigInt& mod, std::vector<mp_word_t>& ws); // TODO
660
661 /**
662 * Set *this to (*this * y) % mod
663 * This function assumes *this is >= 0 && < mod
664 * y should be small, less than 16
665 * @param y the small integer to multiply by
666 * @param mod the positive modulus
667 * @param ws a temp workspace
668 */
669 BigInt& mod_mul(uint8_t y, const BigInt& mod, std::vector<mp_word_t>& ws); // TODO
670
671 /**
672 * @param rng a random number generator
673 * @param min the minimum value (must be non-negative)
674 * @param max the maximum value (must be non-negative and > min)
675 * @return random integer in [min,max)
676 static big_int_t random_integer(RandomNumberGenerator& rng,
677 const big_int_t& min,
678 const big_int_t& max); // TODO
679 */
680
681 std::string to_dec_string(bool add_details=false) const noexcept {
682 // Use the largest power of 10 that fits in a mp_word_t
683 mp_word_t conversion_radix, radix_digits;
684 if constexpr ( 64 == mp_word_bits ) {
685 conversion_radix = 10000000000000000000U;
686 radix_digits = 19;
687 } else {
688 conversion_radix = 1000000000U;
689 radix_digits = 9;
690 }
691
692 // (over-)estimate of the number of digits needed; log2(10) ~ 3.3219
693 const size_t digit_estimate = static_cast<size_t>(1 + (this->bits() / 3.32));
694
695 // (over-)estimate of db such that conversion_radix^db > *this
696 const size_t digit_blocks = (digit_estimate + radix_digits - 1) / radix_digits;
697
698 BigInt value = *this;
699 value.set_sign(positive);
700
701 // Extract groups of digits into words
702 std::vector<mp_word_t> digit_groups(digit_blocks);
703
704 for(size_t i = 0; i != digit_blocks; ++i) {
705 mp_word_t remainder = 0;
706 ct_divide_word(value, conversion_radix, value, remainder);
707 digit_groups[i] = remainder;
708 }
709 assert(value.is_zero());
710
711 // Extract digits from the groups
712 std::vector<uint8_t> digits(digit_blocks * radix_digits);
713
714 for(size_t i = 0; i != digit_blocks; ++i) {
715 mp_word_t remainder = digit_groups[i];
716 for(size_t j = 0; j != radix_digits; ++j)
717 {
718 // Compiler should convert div/mod by 10 into mul by magic constant
719 const mp_word_t digit = remainder % 10;
720 remainder /= 10;
721 digits[radix_digits*i + j] = static_cast<uint8_t>(digit);
722 }
723 }
724
725 // remove leading zeros
726 while(!digits.empty() && digits.back() == 0) {
727 digits.pop_back();
728 }
729
730 assert(digit_estimate >= digits.size());
731
732 // Reverse the digits to big-endian and format to text
733 std::string s;
734 s.reserve(1 + digits.size());
735
736 if(is_negative()) {
737 s += "-";
738 }
739
740 // Reverse and convert to textual digits
741 for(auto i = digits.rbegin(); i != digits.rend(); ++i) {
742 s.push_back(*i + '0'); // assumes ASCII
743 }
744
745 if(s.empty()) {
746 s += "0";
747 }
748 if( add_details ) {
749 append_detail(s);
750 }
751 return s;
752 }
753
754 std::string to_hex_string(bool add_details=false) const noexcept {
755 std::vector<uint8_t> bits;
756 const uint8_t* data;
757 size_t data_len;
758
759 if( is_zero() ) {
760 bits.push_back(0);
761 data = bits.data();
762 data_len = bits.size();
763 } else {
764 data = reinterpret_cast<const uint8_t*>(m_data.const_data());
765 data_len = bytes();
766 }
767
768 std::string s;
769 if(is_negative()) {
770 s += "-";
771 }
772 s.append( jau::bytesHexString(data, 0, data_len,
773 false /* lsbFirst */, true /* lowerCase */) );
774 if( add_details ) {
775 append_detail(s);
776 }
777 return s;
778 }
779
780 private:
781 class data_t
782 {
783 public:
784 data_t() noexcept
785 : m_reg(), m_sig_words(sig_words_npos) {}
786
787 data_t(const data_t& o) noexcept = default;
788
789 data_t(data_t&& o) noexcept {
790 swap(o);
791 }
792 data_t(std::vector<mp_word_t>&& reg) noexcept {
793 swap(reg);
794 }
795
796 ~data_t() noexcept = default;
797
798 data_t& operator=(const data_t& r) noexcept = default;
799
800 data_t& operator=(data_t&& other) noexcept {
801 if(this != &other) {
802 this->swap(other);
803 }
804 return *this;
805 }
806 data_t& operator=(std::vector<mp_word_t>&& other_reg) noexcept {
807 if(&m_reg != &other_reg) {
808 this->swap(other_reg);
809 }
810 return *this;
811 }
812
813 mp_word_t* mutable_data() noexcept {
814 invalidate_sig_words();
815 return m_reg.data();
816 }
817
818 const mp_word_t* const_data() const noexcept {
819 return m_reg.data();
820 }
821
822 std::vector<mp_word_t>& mutable_vector() noexcept {
823 invalidate_sig_words();
824 return m_reg;
825 }
826
827 const std::vector<mp_word_t>& const_vector() const noexcept {
828 return m_reg;
829 }
830
831 mp_word_t get_word_at(size_t n) const noexcept {
832 if(n < m_reg.size()) {
833 return m_reg[n];
834 }
835 return 0;
836 }
837
838 void set_word_at(size_t i, mp_word_t w) {
839 invalidate_sig_words();
840 if(i >= m_reg.size()) {
841 if(w == 0) {
842 return;
843 }
844 grow_to(i + 1);
845 }
846 m_reg[i] = w;
847 }
848
849 void set_words(const mp_word_t w[], size_t len) {
850 invalidate_sig_words();
851 m_reg.assign(w, w + len);
852 }
853
854 void set_to_zero() {
855 m_reg.resize(m_reg.capacity());
856 clear_mem(m_reg.data(), m_reg.size());
857 m_sig_words = 0;
858 }
859
860 void set_size(size_t s) {
861 invalidate_sig_words();
862 clear_mem(m_reg.data(), m_reg.size());
863 m_reg.resize(s + (8 - (s % 8)));
864 }
865
866 void mask_bits(size_t n) noexcept {
867 if(n == 0) { return set_to_zero(); }
868
869 const size_t top_word = n / mp_word_bits;
870
871 // if(top_word < sig_words()) ?
872 if(top_word < size())
873 {
874 const mp_word_t mask = (static_cast<mp_word_t>(1) << (n % mp_word_bits)) - 1;
875 const size_t len = size() - (top_word + 1);
876 if(len > 0)
877 {
878 clear_mem(&m_reg[top_word+1], len);
879 }
880 m_reg[top_word] &= mask;
881 invalidate_sig_words();
882 }
883 }
884
885 void grow_to(size_t n) const {
886 if(n > size()) {
887 if(n <= m_reg.capacity()) {
888 m_reg.resize(n);
889 } else {
890 m_reg.resize(n + (8 - (n % 8)));
891 }
892 }
893 }
894
895 size_t size() const noexcept { return m_reg.size(); }
896
897 void shrink_to_fit(size_t min_size = 0) {
898 const size_t words = std::max(min_size, sig_words());
899 m_reg.resize(words);
900 }
901
902 void resize(size_t s) {
903 m_reg.resize(s);
904 }
905
906 void swap(data_t& other) noexcept {
907 m_reg.swap(other.m_reg);
908 std::swap(m_sig_words, other.m_sig_words);
909 }
910
911 void swap(std::vector<mp_word_t>& reg) noexcept {
912 m_reg.swap(reg);
913 invalidate_sig_words();
914 }
915
916 void invalidate_sig_words() const noexcept {
917 m_sig_words = sig_words_npos;
918 }
919
920 size_t sig_words() const noexcept {
921 if(m_sig_words == sig_words_npos) {
922 m_sig_words = calc_sig_words();
923 } else {
924 assert(m_sig_words == calc_sig_words());
925 }
926 return m_sig_words;
927 }
928 private:
929 static const size_t sig_words_npos = static_cast<size_t>(-1);
930
931 size_t calc_sig_words() const noexcept {
932 const size_t sz = m_reg.size();
933 size_t sig = sz;
934
935 mp_word_t sub = 1;
936
937 for(size_t i = 0; i != sz; ++i)
938 {
939 const mp_word_t w = m_reg[sz - i - 1];
940 sub &= ct_is_zero(w);
941 sig -= sub;
942 }
943
944 /*
945 * This depends on the data so is poisoned, but unpoison it here as
946 * later conditionals are made on the size.
947 */
948 CT::unpoison(sig);
949
950 return sig;
951 }
952
953 mutable std::vector<mp_word_t> m_reg;
954 mutable size_t m_sig_words = sig_words_npos;
955 };
956 data_t m_data;
957 sign_t m_signedness = positive;
958
959 /**
960 * Byte extraction of big-endian value
961 * @param byte_num which byte to extract, 0 == highest byte
962 * @param input the value to extract from
963 * @return byte byte_num of input
964 */
965 template<typename T> static inline constexpr uint8_t get_byte_var_be(size_t byte_num, T input) noexcept {
966 return static_cast<uint8_t>( input >> (((~byte_num)&(sizeof(T)-1)) << 3) );
967 }
968 /**
969 * Byte extraction of little-endian value
970 * @param byte_num which byte to extract, 0 == lowest byte
971 * @param input the value to extract from
972 * @return byte byte_num of input
973 */
974 template<typename T> static inline constexpr uint8_t get_byte_var_le(size_t byte_num, T input) noexcept {
975 return static_cast<uint8_t>( input >> ( byte_num << 3 ) );
976 }
977
978 /**
979 * Zero out some bytes. Warning: use secure_scrub_memory instead if the
980 * memory is about to be freed or otherwise the compiler thinks it can
981 * elide the writes.
982 *
983 * @param ptr a pointer to memory to zero
984 * @param bytes the number of bytes to zero in ptr
985 */
986 static inline constexpr void clear_bytes(void* ptr, size_t bytes) noexcept {
987 if(bytes > 0) {
988 std::memset(ptr, 0, bytes);
989 }
990 }
991
992 /**
993 * Zero memory before use. This simply calls memset and should not be
994 * used in cases where the compiler cannot see the call as a
995 * side-effecting operation (for example, if calling clear_mem before
996 * deallocating memory, the compiler would be allowed to omit the call
997 * to memset entirely under the as-if rule.)
998 *
999 * @param ptr a pointer to an array of Ts to zero
1000 * @param n the number of Ts pointed to by ptr
1001 */
1002 template<typename T> static inline constexpr void clear_mem(T* ptr, size_t n) noexcept {
1003 clear_bytes(ptr, sizeof(T)*n);
1004 }
1005
1006 /**
1007 * Increase internal register buffer to at least n words
1008 * @param n new size of register
1009 */
1010 void grow_to(size_t n) const { m_data.grow_to(n); }
1011
1012 void resize(size_t s) { m_data.resize(s); }
1013
1014 void set_word_at(size_t i, mp_word_t w) {
1015 m_data.set_word_at(i, w);
1016 }
1017
1018 void set_words(const mp_word_t w[], size_t len) {
1019 m_data.set_words(w, len);
1020 }
1021 /**
1022 * Return a mutable pointer to the register
1023 * @result a pointer to the start of the internal register
1024 */
1025 mp_word_t* mutable_data() { return m_data.mutable_data(); }
1026
1027 /**
1028 * Set this number to the value in buf with given byte_len,
1029 * considering the given byte_order.
1030 *
1031 * The value is stored in the local storage format, see \ref bigint_storage_format
1032 */
1033 void binary_decode(const uint8_t buf[], size_t byte_len, const lb_endian_t byte_order) {
1034 const size_t full_words = byte_len / sizeof(mp_word_t);
1035 const size_t extra_bytes = byte_len % sizeof(mp_word_t);
1036
1037 // clear() + setting size
1038 m_signedness = positive;
1039 m_data.set_size( jau::round_up(full_words + (extra_bytes > 0 ? 1U : 0U), 8U) );
1040
1041 mp_word_t* sink = m_data.mutable_data();
1042 if( is_little_endian(byte_order) ) {
1043 // little-endian to local (words arranged as little-endian w/ word itself in native-endian)
1044 for(size_t i = 0; i < full_words; ++i) {
1045 sink[i] = jau::get_value<mp_word_t>(buf + sizeof(mp_word_t)*i, byte_order);
1046 }
1047 } else {
1048 // big-endian to local (words arranged as little-endian w/ word itself in native-endian)
1049 for(size_t i = 0; i < full_words; ++i) {
1050 sink[i] = jau::get_value<mp_word_t>(buf + byte_len - sizeof(mp_word_t)*(i+1), byte_order);
1051 }
1052 }
1053 mp_word_t le_w = 0;
1054 if( is_little_endian(byte_order) ) {
1055 for(size_t i = 0; i < extra_bytes; ++i) {
1056 le_w |= mp_word_t( buf[full_words*sizeof(mp_word_t) + i] ) << ( i * 8 ); // next lowest byte
1057 }
1058 } else {
1059 for(size_t i = 0; i < extra_bytes; ++i) {
1060 le_w = (le_w << 8) | mp_word_t( buf[i] ); // buf[0] highest byte
1061 }
1062 }
1063 sink[full_words] = jau::le_to_cpu( le_w );
1064 }
1065
1066 /**
1067 * Set this number to the decoded hex-string value of buf with given str_len,
1068 * considering the given byte_order.
1069 *
1070 * The value is stored in the local storage format, see \ref bigint_storage_format
1071 */
1072 static BigInt hex_decode(const uint8_t buf[], size_t str_len, const lb_endian_t byte_order) {
1073 BigInt r;
1074
1075 std::vector<uint8_t> bin_out;
1076 const size_t exp_blen = str_len / 2 + str_len % 2;
1077 const size_t blen = jau::hexStringBytes(bin_out, buf, str_len, is_little_endian(byte_order), false /* checkLeading0x */);
1078 if( exp_blen != blen ) {
1079 throw jau::math::MathDomainError("invalid hexadecimal char @ "+std::to_string(blen)+"/"+std::to_string(exp_blen)+" of '"+
1080 std::string(cast_uint8_ptr_to_char(buf), str_len)+"'", E_FILE_LINE);
1081 }
1082 r.binary_decode(bin_out.data(), bin_out.size(), lb_endian_t::little);
1083 return r;
1084 }
1085
1086 static BigInt dec_decode(const uint8_t buf[], size_t str_len) {
1087 BigInt r;
1088
1089 // This could be made faster using the same trick as to_dec_string
1090 for(size_t i = 0; i < str_len; ++i) {
1091 const char c = buf[i];
1092
1093 if(c < '0' || c > '9') {
1094 throw jau::math::MathDomainError("invalid decimal char", E_FILE_LINE);
1095 }
1096 const uint8_t x = c - '0';
1097 assert(x < 10);
1098
1099 r *= 10;
1100 r += x;
1101 }
1102 return r;
1103 }
1104
1105 public:
1106 /**
1107 * Stores this number to the value in buf with given byte_len,
1108 * considering the given byte_order.
1109 *
1110 * The value is read from the local storage in its format, see \ref bigint_storage_format
1111 *
1112 * If byte_len is less than the byt-esize of this integer, i.e. bytes(),
1113 * then it will be truncated.
1114 *
1115 * If byte_len is greater than the byte-size of this integer, i.e. bytes(), it will be zero-padded.
1116 *
1117 * @return actual number of bytes copied, i.e. min(byte_len, bytes());
1118 */
1119 size_t binary_encode(uint8_t output[], size_t byte_len, const lb_endian_t byte_order) const noexcept {
1120 const size_t full_words = byte_len / sizeof(mp_word_t);
1121 const size_t extra_bytes = byte_len % sizeof(mp_word_t);
1122
1123 if( is_little_endian( byte_order ) ) {
1124 // to little-endian from local (words arranged as little-endian w/ word itself in native-endian)
1125 for(size_t i = 0; i < full_words; ++i) {
1126 jau::put_value<mp_word_t>(output + i*sizeof(mp_word_t), word_at(i), byte_order);
1127 }
1128 } else {
1129 // to big-endian from local (words arranged as little-endian w/ word itself in native-endian)
1130 for(size_t i = 0; i < full_words; ++i) {
1131 jau::put_value<mp_word_t>(output + byte_len - (i+1)*sizeof(mp_word_t), word_at(i), byte_order);
1132 }
1133 }
1134 if(extra_bytes > 0) {
1135 const mp_word_t le_w = jau::cpu_to_le( word_at(full_words) );
1136 if( is_little_endian( byte_order ) ) {
1137 for(size_t i = 0; i < extra_bytes; ++i) {
1138 output[full_words*sizeof(mp_word_t) + i] = get_byte_var_le(i, le_w); // next lowest byte
1139 }
1140 } else {
1141 for(size_t i = 0; i < extra_bytes; ++i) {
1142 output[extra_bytes-1-i] = get_byte_var_le(i, le_w); // output[0] highest byte
1143 }
1144 }
1145 }
1146 return extra_bytes + full_words * sizeof(mp_word_t);
1147 }
1148
1149 private:
1150 size_t top_bits_free() const noexcept {
1151 const size_t words = sig_words();
1152
1153 const mp_word_t top_word = word_at(words - 1);
1154 const size_t bits_used = jau::high_bit(top_word);
1155 CT::unpoison(bits_used);
1156 return mp_word_bits - bits_used;
1157 }
1158
1159 /**
1160 * Compares this instance against other, considering sign depending on check_signs
1161 * Returns
1162 * - -1 if this < other
1163 * - 0 if this == other
1164 * - 1 if this > other
1165 */
1166 int cmp(const BigInt& other, bool check_signs = true) const noexcept {
1167 if(check_signs) {
1168 if(other.is_positive() && this->is_negative()) {
1169 return -1;
1170 }
1171 if(other.is_negative() && this->is_positive()) {
1172 return 1;
1173 }
1174 if(other.is_negative() && this->is_negative()) {
1175 return (-ops::bigint_cmp(this->data(), this->size(),
1176 other.data(), other.size()));
1177 }
1178 }
1179 return ops::bigint_cmp(this->data(), this->size(),
1180 other.data(), other.size());
1181 }
1182
1183 bool is_equal(const BigInt& other) const noexcept {
1184 if(this->sign() != other.sign()) {
1185 return false;
1186 }
1187 return ops::bigint_ct_is_eq(this->data(), this->sig_words(),
1188 other.data(), other.sig_words()).is_set();
1189 }
1190
1191 bool is_less_than(const BigInt& other) const noexcept {
1192 if(this->is_negative() && other.is_positive()) {
1193 return true;
1194 }
1195 if(this->is_positive() && other.is_negative()) {
1196 return false;
1197 }
1198 if(other.is_negative() && this->is_negative()) {
1199 return ops::bigint_ct_is_lt(other.data(), other.sig_words(),
1200 this->data(), this->sig_words()).is_set();
1201 }
1202 return ops::bigint_ct_is_lt(this->data(), this->sig_words(),
1203 other.data(), other.sig_words()).is_set();
1204 }
1205
1206 BigInt& add(const mp_word_t y[], size_t y_words, sign_t y_sign) {
1207 const size_t x_sw = sig_words();
1208 grow_to(std::max(x_sw, y_words) + 1);
1209
1210 if(sign() == y_sign)
1211 {
1212 ops::bigint_add2(mutable_data(), size() - 1, y, y_words);
1213 } else {
1214 const int32_t relative_size = ops::bigint_cmp(data(), x_sw, y, y_words);
1215
1216 if(relative_size >= 0)
1217 {
1218 // *this >= y
1219 ops::bigint_sub2(mutable_data(), x_sw, y, y_words);
1220 } else {
1221 // *this < y
1222 ops::bigint_sub2_rev(mutable_data(), y, y_words);
1223 }
1224 //this->sign_fixup(relative_size, y_sign);
1225 if(relative_size < 0) {
1226 set_sign(y_sign);
1227 } else if(relative_size == 0) {
1229 }
1230 }
1231 return (*this);
1232 }
1233
1234 BigInt& operator+=(mp_word_t y) noexcept {
1235 return add(&y, 1, sign_t::positive);
1236 }
1237 BigInt& operator-=(mp_word_t y) noexcept {
1238 return add(&y, 1, sign_t::negative);
1239 }
1240
1241 static BigInt add2(const BigInt& x, const mp_word_t y[], size_t y_words, sign_t y_sign) {
1242 const size_t x_sw = x.sig_words();
1243
1244 BigInt z = BigInt::with_capacity(std::max(x_sw, y_words) + 1);
1245
1246 if(x.sign() == y_sign) {
1247 ops::bigint_add3(z.mutable_data(), x.data(), x_sw, y, y_words);
1248 z.set_sign(x.sign());
1249 } else {
1250 const int32_t relative_size = ops::bigint_sub_abs(z.mutable_data(), x.data(), x_sw, y, y_words);
1251
1252 //z.sign_fixup(relative_size, y_sign);
1253 if(relative_size < 0) {
1254 z.set_sign(y_sign);
1255 } else if(relative_size == 0) {
1256 z.set_sign(positive);
1257 } else {
1258 z.set_sign(x.sign());
1259 }
1260 }
1261 return z;
1262 }
1263
1264 BigInt& mul(const BigInt& y, std::vector<mp_word_t>& ws) noexcept {
1265 const size_t x_sw = sig_words();
1266 const size_t y_sw = y.sig_words();
1267 set_sign((sign() == y.sign()) ? positive : negative);
1268
1269 if(x_sw == 0 || y_sw == 0)
1270 {
1271 clear();
1273 }
1274 else if(x_sw == 1 && y_sw)
1275 {
1276 grow_to(y_sw + 1);
1277 ops::bigint_linmul3(mutable_data(), y.data(), y_sw, word_at(0));
1278 }
1279 else if(y_sw == 1 && x_sw)
1280 {
1281 mp_word_t carry = ops::bigint_linmul2(mutable_data(), x_sw, y.word_at(0));
1282 set_word_at(x_sw, carry);
1283 }
1284 else
1285 {
1286 const size_t new_size = x_sw + y_sw + 1;
1287 // ws.resize(new_size);
1288 (void)ws;
1289 std::vector<mp_word_t> z_reg(new_size);
1290
1291 ops::basecase_mul(z_reg.data(), z_reg.size(), data(), x_sw, y.data(), y_sw);
1292
1293 this->swap_reg(z_reg);
1294 }
1295
1296 return (*this);
1297 }
1298
1300 const size_t x_sw = sig_words();
1301 BigInt z = BigInt::with_capacity(x_sw + 1);
1302
1303 if(x_sw && y) {
1304 ops::bigint_linmul3(z.mutable_data(), data(), x_sw, y);
1305 z.set_sign(sign());
1306 }
1307 return z;
1308 }
1309
1310 void cond_flip_sign(bool predicate) noexcept {
1311 // This code is assuming Negative == 0, Positive == 1
1312 const auto mask = CT::Mask<uint8_t>::expand(predicate);
1313 const uint8_t current_sign = static_cast<uint8_t>(sign());
1314 const uint8_t new_sign = mask.select(current_sign ^ 1, current_sign);
1315 set_sign(static_cast<sign_t>(new_sign));
1316 }
1317
1318 /**
1319 * Return *this % mod
1320 *
1321 * Assumes that *this is (if anything) only slightly larger than
1322 * mod and performs repeated subtractions. It should not be used if
1323 * *this is much larger than mod, instead use modulo operator.
1324 */
1325 inline size_t reduce_below(const BigInt& p, std::vector<mp_word_t>& ws) {
1326 if(p.is_negative() || this->is_negative()) {
1327 std::string msg;
1328 if( p.is_negative() ) {
1329 msg.append("p < 0");
1330 }
1331 if( this->is_negative() ) {
1332 if( msg.length() > 0 ) {
1333 msg.append(" and ");
1334 }
1335 msg.append("*this < 0");
1336 }
1338 }
1339 const size_t p_words = p.sig_words();
1340
1341 if(size() < p_words + 1) {
1342 grow_to(p_words + 1);
1343 }
1344 if(ws.size() < p_words + 1) {
1345 ws.resize(p_words + 1);
1346 }
1347 clear_mem(ws.data(), ws.size());
1348
1349 size_t reductions = 0;
1350
1351 for(;;)
1352 {
1353 mp_word_t borrow = ops::bigint_sub3(ws.data(), data(), p_words + 1, p.data(), p_words);
1354 if(borrow) {
1355 break;
1356 }
1357 ++reductions;
1358 swap_reg(ws);
1359 }
1360 return reductions;
1361 }
1362
1363 static void sign_fixup(const BigInt& x, const BigInt& y, BigInt& q, BigInt& r) {
1364 q.cond_flip_sign(x.sign() != y.sign());
1365
1366 if(x.is_negative() && r.is_nonzero())
1367 {
1368 q -= 1;
1369 r = y.abs() - r;
1370 }
1371 }
1372
1373 static bool division_check(mp_word_t q, mp_word_t y2, mp_word_t y1,
1374 mp_word_t x3, mp_word_t x2, mp_word_t x1) noexcept {
1375 /*
1376 Compute (y3,y2,y1) = (y2,y1) * q
1377 and return true if (y3,y2,y1) > (x3,x2,x1)
1378 */
1379
1380 mp_word_t y3 = 0;
1381 y1 = ops::word_madd2(q, y1, y3);
1382 y2 = ops::word_madd2(q, y2, y3);
1383
1384 const mp_word_t x[3] = { x1, x2, x3 };
1385 const mp_word_t y[3] = { y1, y2, y3 };
1386
1387 return ops::bigint_ct_is_lt(x, 3, y, 3).is_set();
1388 }
1389
1390 /*
1391 * Solve x = q * y + r
1392 *
1393 * See Handbook of Applied Cryptography section 14.2.5
1394 */
1395 static void vartime_divide(const BigInt& x, const BigInt& y_arg, BigInt& q_out, BigInt& r_out) {
1396 if( y_arg.is_zero() ) {
1397 throw jau::math::MathDivByZeroError("y_arg == 0", E_FILE_LINE);
1398 }
1399 const size_t y_words = y_arg.sig_words();
1400
1401 assert(y_words > 0);
1402
1403 BigInt y = y_arg;
1404
1405 BigInt r = x;
1406 BigInt q = BigInt::zero();
1407 std::vector<mp_word_t> ws;
1408
1409 r.set_sign(BigInt::positive);
1410 y.set_sign(BigInt::positive);
1411
1412 // Calculate shifts needed to normalize y with high bit set
1413 const size_t shifts = y.top_bits_free();
1414
1415 y <<= shifts;
1416 r <<= shifts;
1417
1418 // we know y has not changed size, since we only shifted up to set high bit
1419 const size_t t = y_words - 1;
1420 const size_t n = std::max(y_words, r.sig_words()) - 1; // r may have changed size however
1421 assert(n >= t);
1422
1423 q.grow_to(n - t + 1);
1424
1425 mp_word_t* q_words = q.mutable_data();
1426
1427 BigInt shifted_y = y << (mp_word_bits * (n-t));
1428
1429 // Set q_{n-t} to number of times r > shifted_y
1430 q_words[n-t] = r.reduce_below(shifted_y, ws);
1431
1432 const mp_word_t y_t0 = y.word_at(t);
1433 const mp_word_t y_t1 = y.word_at(t-1);
1434 assert((y_t0 >> (mp_word_bits-1)) == 1);
1435
1436 for(size_t j = n; j != t; --j) {
1437 const mp_word_t x_j0 = r.word_at(j);
1438 const mp_word_t x_j1 = r.word_at(j-1);
1439 const mp_word_t x_j2 = r.word_at(j-2);
1440
1441 mp_word_t qjt = ops::bigint_divop(x_j0, x_j1, y_t0);
1442
1443 qjt = CT::Mask<mp_word_t>::is_equal(x_j0, y_t0).select(mp_word_max, qjt);
1444
1445 // Per HAC 14.23, this operation is required at most twice
1446 qjt -= division_check(qjt, y_t0, y_t1, x_j0, x_j1, x_j2);
1447 qjt -= division_check(qjt, y_t0, y_t1, x_j0, x_j1, x_j2);
1448 assert(division_check(qjt, y_t0, y_t1, x_j0, x_j1, x_j2) == false);
1449
1450 shifted_y >>= mp_word_bits;
1451 // Now shifted_y == y << (BOTAN_MP_WORD_BITS * (j-t-1))
1452
1453 // TODO this sequence could be better
1454 r -= shifted_y * qjt;
1455 qjt -= r.is_negative();
1456 r += shifted_y * static_cast<mp_word_t>(r.is_negative());
1457
1458 q_words[j-t-1] = qjt;
1459 }
1460
1461 r >>= shifts;
1462
1463 sign_fixup(x, y_arg, q, r);
1464
1465 r_out = r;
1466 q_out = q;
1467 }
1468
1469 BigInt operator/(const mp_word_t& y) const {
1470 if(y == 0) {
1472 }
1473 BigInt q;
1474 mp_word_t r;
1475 ct_divide_word(*this, y, q, r);
1476 return q;
1477 }
1478
1479 static void ct_divide_word(const BigInt& x, mp_word_t y, BigInt& q_out, mp_word_t& r_out) {
1480 if(y == 0) {
1482 }
1483 const size_t x_words = x.sig_words();
1484 const size_t x_bits = x.bits();
1485
1486 BigInt q = BigInt::with_capacity(x_words);
1487 mp_word_t r = 0;
1488
1489 for(size_t i = 0; i != x_bits; ++i)
1490 {
1491 const size_t b = x_bits - 1 - i;
1492 const bool x_b = x.get_bit(b);
1493
1494 const auto r_carry = CT::Mask<mp_word_t>::expand(r >> (mp_word_bits - 1));
1495
1496 r *= 2;
1497 r += x_b;
1498
1499 const auto r_gte_y = CT::Mask<mp_word_t>::is_gte(r, y) | r_carry;
1500 q.conditionally_set_bit(b, r_gte_y.is_set());
1501 r = r_gte_y.select(r - y, r);
1502 }
1503
1504 if(x.is_negative()) {
1505 q.flip_sign();
1506 if(r != 0) {
1507 --q;
1508 r = y - r;
1509 }
1510 }
1511
1512 r_out = r;
1513 q_out = q;
1514 }
1515
1517 if(mod == 0) {
1518 throw jau::math::MathDivByZeroError("mod == 0", E_FILE_LINE);
1519 }
1520 if(mod == 1) {
1521 return 0;
1522 }
1523 mp_word_t remainder = 0;
1524
1525 if( jau::is_power_of_2(mod) ) {
1526 remainder = (word_at(0) & (mod - 1));
1527 } else {
1528 const size_t sw = sig_words();
1529 for(size_t i = sw; i > 0; --i) {
1530 remainder = ops::bigint_modop(remainder, word_at(i-1), mod);
1531 }
1532 }
1533 if(remainder && sign() == BigInt::negative) {
1534 return mod - remainder;
1535 }
1536 return remainder;
1537 }
1538
1539 void append_detail(std::string& s) const noexcept {
1540 s.append(", bits ").append(std::to_string(bits())).append(", ").append(std::to_string(sig_words())).append(" word(s): ");
1541 for(size_t i=0; i<sig_words(); ++i) {
1542 const mp_word_t w = word_at(i);
1543 s.append( jau::bytesHexString(&w, 0, mp_word_bits/CHAR_BIT, false /* lsbFirst */, true /* lowerCase */) )
1544 .append(", ");
1545 }
1546 }
1547 };
1548
1549 /**@}*/
1550}
1551
1552namespace jau {
1553
1554 /** \addtogroup Integer
1555 *
1556 * @{
1557 */
1558
1559 inline mp::BigInt abs(mp::BigInt x) noexcept { return x.abs(); }
1560 inline mp::BigInt pow(mp::BigInt b, mp::BigInt e) { return b.pow(e); }
1561
1562 inline const mp::BigInt& min(const mp::BigInt& x, const mp::BigInt& y) noexcept {
1563 return x < y ? x : y;
1564 }
1565 inline const mp::BigInt& max(const mp::BigInt& x, const mp::BigInt& y) noexcept {
1566 return x > y ? x : y;
1567 }
1568 inline const mp::BigInt& clamp(const mp::BigInt& x, const mp::BigInt& min_val, const mp::BigInt& max_val) noexcept {
1569 return min(max(x, min_val), max_val);
1570 }
1571
1572 inline mp::BigInt& min(mp::BigInt& x, mp::BigInt& y) noexcept {
1573 return x < y ? x : y;
1574 }
1575 inline mp::BigInt& max(mp::BigInt& x, mp::BigInt& y) noexcept {
1576 return x > y ? x : y;
1577 }
1578 inline mp::BigInt& clamp(mp::BigInt& x, mp::BigInt& min_val, mp::BigInt& max_val) noexcept {
1579 return min(max(x, min_val), max_val);
1580 }
1581
1582 inline mp::BigInt gcd(const mp::BigInt& a, const mp::BigInt& b) noexcept {
1583 mp::BigInt a_ = abs(a);
1584 mp::BigInt b_ = abs(b);
1585 while( b_.is_nonzero() ) {
1586 const mp::BigInt t = b_;
1587 b_ = a_ % b_;
1588 a_ = t;
1589 }
1590 return a_;
1591 }
1592
1593 /**@}*/
1594}
1595
1596namespace std {
1597 inline std::ostream& operator<<(std::ostream& out, const jau::mp::BigInt& v) {
1598 return out << v.to_dec_string();
1599 }
1600}
1601
1602#endif /** JAU_BIG_INT_HPP_ */
#define E_FILE_LINE
static Mask< T > is_gte(T x, T y) noexcept
Return a Mask<T> which is set if x >= y.
Definition: ct_utils.hpp:177
static Mask< T > is_equal(T x, T y) noexcept
Return a Mask<T> which is set if x == y.
Definition: ct_utils.hpp:145
static Mask< T > expand(T v) noexcept
Return a Mask<T> which is set if v is != 0.
Definition: ct_utils.hpp:119
math_error_t::div_by_zero, i.e.
Definition: math_error.hpp:78
math_error_t::invalid
Definition: math_error.hpp:72
Arbitrary precision integer type.
Definition: big_int.hpp:43
BigInt & set_sign(sign_t sign) noexcept
Set sign of the integer.
Definition: big_int.hpp:281
BigInt(uint64_t n)
Create big_int_t from an unsigned 64 bit integer.
Definition: big_int.hpp:129
bool is_even() const noexcept
Test if the integer has an even value.
Definition: big_int.hpp:341
size_t bits() const noexcept
Returns bit length of this integer.
Definition: big_int.hpp:310
bool operator==(const BigInt &b) const noexcept
Definition: big_int.hpp:397
bool is_odd() const noexcept
Test if the integer has an odd value.
Definition: big_int.hpp:347
bool get_bit(size_t n) const noexcept
Return bit value at specified position.
Definition: big_int.hpp:368
BigInt & operator-=(const BigInt &y) noexcept
Definition: big_int.hpp:422
BigInt abs() const noexcept
Returns absolute (positive) value of this instance.
Definition: big_int.hpp:290
BigInt & operator++() noexcept
Definition: big_int.hpp:410
BigInt operator-() const noexcept
Unary negation operator, returns new negative instance of this.
Definition: big_int.hpp:216
BigInt & mod_add(const BigInt &y, const BigInt &mod, std::vector< mp_word_t > &ws)
Set *this to (*this + y) % mod This function assumes *this is >= 0 && < mod.
BigInt pow(BigInt e)
Returns (*this)^e, or pow(*this, e)
Definition: big_int.hpp:568
BigInt & operator--() noexcept
Definition: big_int.hpp:412
size_t bytes() const
Returns byte length of this integer.
Definition: big_int.hpp:307
BigInt & operator<<=(size_t shift) noexcept
Definition: big_int.hpp:434
BigInt(const uint8_t buf[], size_t byte_len, const lb_endian_t byte_order)
Create a big_int_t from an integer in a byte array with given byte_len, considering the given byte_or...
Definition: big_int.hpp:178
std::string to_dec_string(bool add_details=false) const noexcept
Definition: big_int.hpp:681
BigInt & operator=(const BigInt &r)=default
BigInt operator>>(size_t shift) const
Definition: big_int.hpp:473
std::string to_hex_string(bool add_details=false) const noexcept
Definition: big_int.hpp:754
BigInt operator--(int) noexcept
Definition: big_int.hpp:416
sign_t reverse_sign() const noexcept
Definition: big_int.hpp:263
bool operator<=(const BigInt &b) const noexcept
Definition: big_int.hpp:399
BigInt operator+(const BigInt &y) const noexcept
Definition: big_int.hpp:426
sign_t
Sign symbol definitions for positive and negative numbers.
Definition: big_int.hpp:48
BigInt operator%(const BigInt &mod)
Definition: big_int.hpp:543
bool operator>(const BigInt &b) const noexcept
Definition: big_int.hpp:402
const mp_word_t * data() const
Return a const pointer to the register.
Definition: big_int.hpp:240
bool is_nonzero() const noexcept
Test if the integer is not zero.
Definition: big_int.hpp:353
size_t size() const noexcept
Give size of internal register.
Definition: big_int.hpp:298
static BigInt one()
Create a 1-value big_int.
Definition: big_int.hpp:65
sign_t sign() const noexcept
Return the sign of the integer.
Definition: big_int.hpp:258
bool is_negative() const noexcept
Tests if the sign of the integer is negative.
Definition: big_int.hpp:246
BigInt & operator>>=(size_t shift) noexcept
Definition: big_int.hpp:450
BigInt & rev_sub(const BigInt &y, std::vector< mp_word_t > &ws)
Set *this to y - *this.
BigInt & operator=(BigInt &&other) noexcept
Definition: big_int.hpp:192
static BigInt zero()
Create a 0-value big_int.
Definition: big_int.hpp:60
void conditionally_set_bit(size_t n, bool set_it) noexcept
Conditionally set bit at specified position.
Definition: big_int.hpp:388
size_t sig_words() const noexcept
Return how many words we need to hold this value.
Definition: big_int.hpp:304
void clear()
Zeroize the big_int.
Definition: big_int.hpp:324
BigInt(std::vector< mp_word_t > &&other_reg) noexcept
Definition: big_int.hpp:182
BigInt() noexcept=default
BigInt & operator%=(const BigInt &mod)
Modulo operator.
Definition: big_int.hpp:539
bool operator!() const noexcept
! operator, returns true iff this is zero, otherwise false.
Definition: big_int.hpp:395
BigInt & mod_sub(const BigInt &y, const BigInt &mod, std::vector< mp_word_t > &ws)
Set *this to (*this - y) % mod This function assumes *this is >= 0 && < mod.
bool is_positive() const noexcept
Tests if the sign of the integer is positive.
Definition: big_int.hpp:252
bool operator<(const BigInt &b) const noexcept
Definition: big_int.hpp:401
bool operator!=(const BigInt &b) const noexcept
Definition: big_int.hpp:398
bool operator>=(const BigInt &b) const noexcept
Definition: big_int.hpp:400
size_t binary_encode(uint8_t output[], size_t byte_len, const lb_endian_t byte_order) const noexcept
Stores this number to the value in buf with given byte_len, considering the given byte_order.
Definition: big_int.hpp:1119
void swap(BigInt &other) noexcept
Swap this value with another.
Definition: big_int.hpp:203
int compare(const BigInt &b) const noexcept
Compares this instance against other considering both sign() value Returns.
Definition: big_int.hpp:333
static BigInt from_s32(int32_t n)
Create big_int from a signed 32 bit integer.
Definition: big_int.hpp:96
BigInt & square(std::vector< mp_word_t > &ws)
Square value of *this.
BigInt mod_pow(BigInt e, BigInt m)
Returns (*this)^e % m, or pow(*this, e) % m.
Definition: big_int.hpp:602
BigInt & operator/=(const BigInt &y)
Definition: big_int.hpp:518
static BigInt from_word(mp_word_t n)
Create big_int from a mp_word_t (limb)
Definition: big_int.hpp:86
static BigInt with_capacity(size_t n)
Create big_int of specified size, all zeros.
Definition: big_int.hpp:108
bool is_zero() const noexcept
Test if the integer is zero.
Definition: big_int.hpp:359
mp_word_t word_at(size_t n) const noexcept
Return the mp_word_t at a specified position of the internal register.
Definition: big_int.hpp:232
BigInt operator/(const BigInt &y) const
Definition: big_int.hpp:527
BigInt operator-(const BigInt &y) const noexcept
Definition: big_int.hpp:430
BigInt & operator*=(const BigInt &y) noexcept
Definition: big_int.hpp:494
BigInt & mod_mul(uint8_t y, const BigInt &mod, std::vector< mp_word_t > &ws)
Set *this to (*this * y) % mod This function assumes *this is >= 0 && < mod y should be small,...
uint8_t byte_at(size_t n) const noexcept
Definition: big_int.hpp:222
BigInt & flip_sign() noexcept
Flip the sign of this big_int.
Definition: big_int.hpp:273
BigInt & operator+=(const BigInt &y) noexcept
Definition: big_int.hpp:418
void set_bit(size_t n) noexcept
Set bit at specified position.
Definition: big_int.hpp:376
BigInt operator*(const BigInt &y) noexcept
Definition: big_int.hpp:499
BigInt operator<<(size_t shift) const
Definition: big_int.hpp:462
BigInt(const std::string &str)
Construct a big_int_t from a string encoded as hexadecimal or decimal.
Definition: big_int.hpp:145
static BigInt from_u64(uint64_t n)
Create big_int from an unsigned 64 bit integer.
Definition: big_int.hpp:71
BigInt operator++(int) noexcept
Definition: big_int.hpp:414
BigInt(BigInt &&other) noexcept
Definition: big_int.hpp:186
static BigInt power_of_2(size_t n)
Create a power of two.
Definition: big_int.hpp:119
char * cast_uint8_ptr_to_char(uint8_t *b) noexcept
Definition: byte_util.hpp:148
constexpr uint16_t le_to_cpu(uint16_t const l) noexcept
Definition: byte_util.hpp:396
lb_endian_t
Simplified reduced endian type only covering little- and big-endian.
Definition: byte_util.hpp:227
constexpr bool is_little_endian(const endian_t byte_order) noexcept
Returns true if given byte_order equals endian::little, otherwise false.
Definition: byte_util.hpp:293
constexpr uint16_t cpu_to_le(uint16_t const h) noexcept
Definition: byte_util.hpp:405
const uint8_t * cast_char_ptr_to_uint8(const char *s) noexcept
Definition: byte_util.hpp:155
@ little
Identifier for little endian, equivalent to endian::little.
@ big
Identifier for big endian, equivalent to endian::big.
std::string to_string(const alphabet &v) noexcept
Definition: base_codec.hpp:97
constexpr T ct_is_zero(T x)
Returns ~0 (2-complement) if arg is zero, otherwise 0 (w/o branching) in O(1) and constant time (CT).
void swap(cow_darray< Value_type, Size_type, Alloc_type > &rhs, cow_darray< Value_type, Size_type, Alloc_type > &lhs) noexcept
constexpr T clamp(const T x, const T min_val, const T max_val) noexcept
Returns constrained integral value to lie between given min- and maximum value (w/ branching) in O(1)...
Definition: base_math.hpp:208
constexpr T round_up(const T n, const U align_to)
Round up w/ branching in O(1)
Definition: int_math.hpp:101
constexpr nsize_t high_bit(T x)
Return the index of the highest set bit w/ branching (loop) in O(n), actually O(n/2).
Definition: int_math.hpp:164
constexpr T min(const T x, const T y) noexcept
Returns the minimum of two integrals (w/ branching) in O(1)
Definition: base_math.hpp:177
mp::BigInt & max(mp::BigInt &x, mp::BigInt &y) noexcept
Definition: big_int.hpp:1575
constexpr T gcd(T a, T b) noexcept
Returns the greatest common divisor (GCD) of the two given integer values following Euclid's algorith...
Definition: int_math.hpp:298
constexpr T max(const T x, const T y) noexcept
Returns the maximum of two integrals (w/ branching) in O(1)
Definition: base_math.hpp:191
mp::BigInt pow(mp::BigInt b, mp::BigInt e)
Definition: big_int.hpp:1560
constexpr T abs(const T x) noexcept
Returns the absolute value of an arithmetic number (w/ branching) in O(1)
Definition: base_math.hpp:155
constexpr bool is_power_of_2(const T x) noexcept
Power of 2 test (w/o branching ?) in O(1)
Definition: int_math.hpp:140
std::string bytesHexString(const void *data, const nsize_t offset, const nsize_t length, const bool lsbFirst, const bool lowerCase=true) noexcept
Produce a hexadecimal string representation of the given byte values.
size_t hexStringBytes(std::vector< uint8_t > &out, const std::string &hexstr, const bool lsbFirst, const bool checkLeading0x) noexcept
Converts a given hexadecimal string representation into a byte vector.
void unpoison(const T *p, size_t n)
Definition: ct_utils.hpp:54
CT::Mask< mp_word_t > bigint_ct_is_lt(const mp_word_t x[], size_t x_size, const mp_word_t y[], size_t y_size, bool lt_or_equal=false) noexcept
Compare x and y Return ~0 if x[0:x_size] < y[0:y_size] or 0 otherwise If lt_or_equal is true,...
CT::Mask< mp_word_t > bigint_ct_is_eq(const mp_word_t x[], size_t x_size, const mp_word_t y[], size_t y_size) noexcept
int bigint_cmp(const mp_word_t x[], size_t x_size, const mp_word_t y[], size_t y_size) noexcept
Compare unsigned x and y mp_word_t.
mp_word_t bigint_divop(mp_word_t n1, mp_word_t n0, mp_word_t d)
Computes ((n1<<bits) + n0) / d.
int32_t bigint_sub_abs(mp_word_t z[], const mp_word_t x[], size_t x_size, const mp_word_t y[], size_t y_size) noexcept
Set z to abs(x-y), ie if x >= y, then compute z = x - y Otherwise compute z = y - x No borrow is poss...
mp_word_t bigint_add2(mp_word_t x[], size_t x_size, const mp_word_t y[], size_t y_size) noexcept
Two operand addition with carry out.
void bigint_shr1(mp_word_t x[], size_t x_size, size_t word_shift, size_t bit_shift) noexcept
mp_word_t bigint_linmul2(mp_word_t x[], size_t x_size, mp_word_t y) noexcept
void bigint_linmul3(mp_word_t z[], const mp_word_t x[], size_t x_size, mp_word_t y) noexcept
mp_word_t bigint_sub3(mp_word_t z[], const mp_word_t x[], size_t x_size, const mp_word_t y[], size_t y_size) noexcept
Three operand subtraction.
void bigint_shl2(mp_word_t y[], const mp_word_t x[], size_t x_size, size_t word_shift, size_t bit_shift) noexcept
void bigint_add3(mp_word_t z[], const mp_word_t x[], size_t x_size, const mp_word_t y[], size_t y_size) noexcept
Three operand addition.
mp_word_t bigint_sub2(mp_word_t x[], size_t x_size, const mp_word_t y[], size_t y_size) noexcept
Two operand subtraction.
mp_word_t bigint_modop(mp_word_t n1, mp_word_t n0, mp_word_t d)
Compute ((n1<<bits) + n0) % d.
mp_word_t word_madd2(mp_word_t a, mp_word_t b, mp_word_t &c) noexcept
Word Multiply/Add.
void bigint_shl1(mp_word_t x[], size_t x_size, size_t x_words, size_t word_shift, size_t bit_shift) noexcept
void basecase_mul(mp_word_t z[], size_t z_size, const mp_word_t x[], size_t x_size, const mp_word_t y[], size_t y_size) noexcept
void bigint_shr2(mp_word_t y[], const mp_word_t x[], size_t x_size, size_t word_shift, size_t bit_shift) noexcept
void bigint_sub2_rev(mp_word_t x[], const mp_word_t y[], size_t y_size) noexcept
Two operand subtraction, x = y - x; assumes y >= x.
big_int_t (this file) (c) 2024 Gothel Software e.K.
Definition: big_int.hpp:26
jau::uint_bytes< impl::best_word_byte_size()>::type mp_word_t
Definition: big_int_ops.hpp:42
constexpr const size_t mp_word_bits
Definition: big_int_ops.hpp:41
constexpr const mp_word_t mp_word_max
Definition: big_int_ops.hpp:60
__pack(...): Produces MSVC, clang and gcc compatible lead-in and -out macros.
Definition: backtrace.hpp:32
STL namespace.
std::ostream & operator<<(std::ostream &os, const jau::fraction_timespec &v) noexcept
Output stream operator for jau::fraction_timespec.