Botan 2.19.4
Crypto and TLS for C&
numthry.cpp
Go to the documentation of this file.
1/*
2* Number Theory Functions
3* (C) 1999-2011,2016,2018,2019 Jack Lloyd
4*
5* Botan is released under the Simplified BSD License (see license.txt)
6*/
7
8#include <botan/numthry.h>
9#include <botan/reducer.h>
10#include <botan/monty.h>
11#include <botan/divide.h>
12#include <botan/rng.h>
13#include <botan/internal/ct_utils.h>
14#include <botan/internal/mp_core.h>
15#include <botan/internal/monty_exp.h>
16#include <botan/internal/primality.h>
17#include <algorithm>
18
19namespace Botan {
20
21/*
22* Return the number of 0 bits at the end of n
23*/
24size_t low_zero_bits(const BigInt& n)
25 {
26 size_t low_zero = 0;
27
28 auto seen_nonempty_word = CT::Mask<word>::cleared();
29
30 for(size_t i = 0; i != n.size(); ++i)
31 {
32 const word x = n.word_at(i);
33
34 // ctz(0) will return sizeof(word)
35 const size_t tz_x = ctz(x);
36
37 // if x > 0 we want to count tz_x in total but not any
38 // further words, so set the mask after the addition
39 low_zero += seen_nonempty_word.if_not_set_return(tz_x);
40
41 seen_nonempty_word |= CT::Mask<word>::expand(x);
42 }
43
44 // if we saw no words with x > 0 then n == 0 and the value we have
45 // computed is meaningless. Instead return 0 in that case.
46 return seen_nonempty_word.if_set_return(low_zero);
47 }
48
49/*
50* Calculate the GCD in constant time
51*/
52BigInt gcd(const BigInt& a, const BigInt& b)
53 {
54 if(a.is_zero())
55 return abs(b);
56 if(b.is_zero())
57 return abs(a);
58
59 if(a == 1 || b == 1)
60 return 1;
61
62 const size_t sz = std::max(a.sig_words(), b.sig_words());
63 auto u = BigInt::with_capacity(sz);
64 auto v = BigInt::with_capacity(sz);
65 u += a;
66 v += b;
67
69 v.const_time_poison();
70
71 u.set_sign(BigInt::Positive);
72 v.set_sign(BigInt::Positive);
73
74 // In the worst case we have two fully populated big ints. After right
75 // shifting so many times, we'll have reached the result for sure.
76 const size_t loop_cnt = u.bits() + v.bits();
77
78 using WordMask = CT::Mask<word>;
79
80 // This temporary is big enough to hold all intermediate results of the
81 // algorithm. No reallocation will happen during the loop.
82 // Note however, that `ct_cond_assign()` will invalidate the 'sig_words'
83 // cache, which _does not_ shrink the capacity of the underlying buffer.
84 auto tmp = BigInt::with_capacity(sz);
85 size_t factors_of_two = 0;
86 for(size_t i = 0; i != loop_cnt; ++i) {
87 auto both_odd = WordMask::expand(u.is_odd()) & WordMask::expand(v.is_odd());
88
89 // Subtract the smaller from the larger if both are odd
90 auto u_gt_v = WordMask::expand(bigint_cmp(u.data(), u.size(), v.data(), v.size()) > 0);
91 bigint_sub_abs(tmp.mutable_data(), u.data(), sz, v.data(), sz);
92 u.ct_cond_assign((u_gt_v & both_odd).is_set(), tmp);
93 v.ct_cond_assign((~u_gt_v & both_odd).is_set(), tmp);
94
95 const auto u_is_even = WordMask::expand(u.is_even());
96 const auto v_is_even = WordMask::expand(v.is_even());
97 BOTAN_DEBUG_ASSERT((u_is_even | v_is_even).is_set());
98
99 // When both are even, we're going to eliminate a factor of 2.
100 // We have to reapply this factor to the final result.
101 factors_of_two += (u_is_even & v_is_even).if_set_return(1);
102
103 // remove one factor of 2, if u is even
104 bigint_shr2(tmp.mutable_data(), u.data(), sz, 0, 1);
105 u.ct_cond_assign(u_is_even.is_set(), tmp);
106
107 // remove one factor of 2, if v is even
108 bigint_shr2(tmp.mutable_data(), v.data(), sz, 0, 1);
109 v.ct_cond_assign(v_is_even.is_set(), tmp);
110 }
111
112 // The GCD (without factors of two) is either in u or v, the other one is
113 // zero. The non-zero variable _must_ be odd, because all factors of two were
114 // removed in the loop iterations above.
115 BOTAN_DEBUG_ASSERT(u.is_zero() || v.is_zero());
116 BOTAN_DEBUG_ASSERT(u.is_odd() || v.is_odd());
117
118 // make sure that the GCD (without factors of two) is in u
119 u.ct_cond_assign(u.is_even() /* .is_zero() would not be constant time */, v);
120
121 // re-apply the factors of two
122 u.ct_shift_left(factors_of_two);
123
124 u.const_time_unpoison();
125 v.const_time_unpoison();
126
127 return u;
128}
129
130/*
131* Calculate the LCM
132*/
133BigInt lcm(const BigInt& a, const BigInt& b)
134 {
135 return ct_divide(a * b, gcd(a, b));
136 }
137
138/*
139* Modular Exponentiation
140*/
141BigInt power_mod(const BigInt& base, const BigInt& exp, const BigInt& mod)
142 {
143 if(mod.is_negative() || mod == 1)
144 {
145 return 0;
146 }
147
148 if(base.is_zero() || mod.is_zero())
149 {
150 if(exp.is_zero())
151 return 1;
152 return 0;
153 }
154
155 Modular_Reducer reduce_mod(mod);
156
157 const size_t exp_bits = exp.bits();
158
159 if(mod.is_odd())
160 {
161 const size_t powm_window = 4;
162
163 auto monty_mod = std::make_shared<Montgomery_Params>(mod, reduce_mod);
164 auto powm_base_mod = monty_precompute(monty_mod, reduce_mod.reduce(base), powm_window);
165 return monty_execute(*powm_base_mod, exp, exp_bits);
166 }
167
168 /*
169 Support for even modulus is just a convenience and not considered
170 cryptographically important, so this implementation is slow ...
171 */
172 BigInt accum = 1;
173 BigInt g = reduce_mod.reduce(base);
174 BigInt t;
175
176 for(size_t i = 0; i != exp_bits; ++i)
177 {
178 t = reduce_mod.multiply(g, accum);
179 g = reduce_mod.square(g);
180 accum.ct_cond_assign(exp.get_bit(i), t);
181 }
182 return accum;
183 }
184
185
187 {
188 if(C < 1)
189 throw Invalid_Argument("is_perfect_square requires C >= 1");
190 if(C == 1)
191 return 1;
192
193 const size_t n = C.bits();
194 const size_t m = (n + 1) / 2;
195 const BigInt B = C + BigInt::power_of_2(m);
196
197 BigInt X = BigInt::power_of_2(m) - 1;
198 BigInt X2 = (X*X);
199
200 for(;;)
201 {
202 X = (X2 + C) / (2*X);
203 X2 = (X*X);
204
205 if(X2 < B)
206 break;
207 }
208
209 if(X2 == C)
210 return X;
211 else
212 return 0;
213 }
214
215/*
216* Test for primality using Miller-Rabin
217*/
218bool is_prime(const BigInt& n,
220 size_t prob,
221 bool is_random)
222 {
223 if(n == 2)
224 return true;
225 if(n <= 1 || n.is_even())
226 return false;
227
228 const size_t n_bits = n.bits();
229
230 // Fast path testing for small numbers (<= 65521)
231 if(n_bits <= 16)
232 {
233 const uint16_t num = static_cast<uint16_t>(n.word_at(0));
234
235 return std::binary_search(PRIMES, PRIMES + PRIME_TABLE_SIZE, num);
236 }
237
238 Modular_Reducer mod_n(n);
239
240 if(rng.is_seeded())
241 {
242 const size_t t = miller_rabin_test_iterations(n_bits, prob, is_random);
243
244 if(is_miller_rabin_probable_prime(n, mod_n, rng, t) == false)
245 return false;
246
247 if(is_random)
248 return true;
249 else
250 return is_lucas_probable_prime(n, mod_n);
251 }
252 else
253 {
254 return is_bailie_psw_probable_prime(n, mod_n);
255 }
256 }
257
258}
#define BOTAN_DEBUG_ASSERT(expr)
Definition: assert.h:123
size_t sig_words() const
Definition: bigint.h:592
bool is_odd() const
Definition: bigint.h:415
void ct_cond_assign(bool predicate, const BigInt &other)
Definition: bigint.cpp:530
size_t size() const
Definition: bigint.h:586
word word_at(size_t n) const
Definition: bigint.h:514
size_t bits() const
Definition: bigint.cpp:304
bool is_even() const
Definition: bigint.h:409
static BigInt power_of_2(size_t n)
Definition: bigint.h:770
void const_time_poison() const
Definition: bigint.h:751
bool is_zero() const
Definition: bigint.h:427
bool is_negative() const
Definition: bigint.h:533
bool get_bit(size_t n) const
Definition: bigint.h:471
static BigInt with_capacity(size_t n)
Definition: bigint.cpp:49
static Mask< T > expand(T v)
Definition: ct_utils.h:123
static Mask< T > cleared()
Definition: ct_utils.h:115
BigInt square(const BigInt &x) const
Definition: reducer.h:39
BigInt multiply(const BigInt &x, const BigInt &y) const
Definition: reducer.h:31
BigInt reduce(const BigInt &x) const
Definition: reducer.cpp:37
virtual bool is_seeded() const =0
fe X
Definition: ge.cpp:27
Definition: alg_id.cpp:13
BigInt power_mod(const BigInt &base, const BigInt &exp, const BigInt &mod)
Definition: numthry.cpp:141
BigInt lcm(const BigInt &a, const BigInt &b)
Definition: numthry.cpp:133
size_t ctz(T n)
Definition: bit_ops.h:99
const uint16_t PRIMES[]
Definition: primes.cpp:12
size_t low_zero_bits(const BigInt &n)
Definition: numthry.cpp:24
BigInt abs(const BigInt &n)
Definition: numthry.h:58
CT::Mask< word > bigint_sub_abs(word z[], const word x[], const word y[], size_t N, word ws[])
Definition: mp_core.h:377
const size_t PRIME_TABLE_SIZE
Definition: numthry.h:287
bool is_prime(const BigInt &n, RandomNumberGenerator &rng, size_t prob, bool is_random)
Definition: numthry.cpp:218
bool is_miller_rabin_probable_prime(const BigInt &n, const Modular_Reducer &mod_n, RandomNumberGenerator &rng, size_t test_iterations)
Definition: primality.cpp:143
void ct_divide(const BigInt &x, const BigInt &y, BigInt &q_out, BigInt &r_out)
Definition: divide.cpp:52
std::shared_ptr< const Montgomery_Exponentation_State > monty_precompute(std::shared_ptr< const Montgomery_Params > params, const BigInt &g, size_t window_bits, bool const_time)
Definition: monty_exp.cpp:157
bool is_bailie_psw_probable_prime(const BigInt &n, const Modular_Reducer &mod_n)
Definition: primality.cpp:92
BigInt gcd(const BigInt &a, const BigInt &b)
Definition: numthry.cpp:52
void bigint_shr2(word y[], const word x[], size_t x_size, size_t word_shift, size_t bit_shift)
Definition: mp_core.h:466
BigInt is_perfect_square(const BigInt &C)
Definition: numthry.cpp:186
size_t miller_rabin_test_iterations(size_t n_bits, size_t prob, bool random)
Definition: primality.cpp:165
int32_t bigint_cmp(const word x[], size_t x_size, const word y[], size_t y_size)
Definition: mp_core.h:525
bool is_lucas_probable_prime(const BigInt &C, const Modular_Reducer &mod_C)
Definition: primality.cpp:17
BigInt monty_execute(const Montgomery_Exponentation_State &precomputed_state, const BigInt &k, size_t max_k_bits)
Definition: monty_exp.cpp:165