Botan 2.19.1
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
21namespace {
22
23void sub_abs(BigInt& z, const BigInt& x, const BigInt& y)
24 {
25 const size_t x_sw = x.sig_words();
26 const size_t y_sw = y.sig_words();
27 z.resize(std::max(x_sw, y_sw));
28
29 bigint_sub_abs(z.mutable_data(),
30 x.data(), x_sw,
31 y.data(), y_sw);
32 }
33
34}
35
36/*
37* Return the number of 0 bits at the end of n
38*/
39size_t low_zero_bits(const BigInt& n)
40 {
41 size_t low_zero = 0;
42
43 auto seen_nonempty_word = CT::Mask<word>::cleared();
44
45 for(size_t i = 0; i != n.size(); ++i)
46 {
47 const word x = n.word_at(i);
48
49 // ctz(0) will return sizeof(word)
50 const size_t tz_x = ctz(x);
51
52 // if x > 0 we want to count tz_x in total but not any
53 // further words, so set the mask after the addition
54 low_zero += seen_nonempty_word.if_not_set_return(tz_x);
55
56 seen_nonempty_word |= CT::Mask<word>::expand(x);
57 }
58
59 // if we saw no words with x > 0 then n == 0 and the value we have
60 // computed is meaningless. Instead return 0 in that case.
61 return seen_nonempty_word.if_set_return(low_zero);
62 }
63
64namespace {
65
66size_t safegcd_loop_bound(size_t f_bits, size_t g_bits)
67 {
68 const size_t d = std::max(f_bits, g_bits);
69
70 if(d < 46)
71 return (49*d + 80) / 17;
72 else
73 return (49*d + 57) / 17;
74 }
75
76}
77
78/*
79* Calculate the GCD
80*/
81BigInt gcd(const BigInt& a, const BigInt& b)
82 {
83 if(a.is_zero())
84 return abs(b);
85 if(b.is_zero())
86 return abs(a);
87 if(a == 1 || b == 1)
88 return 1;
89
90 // See https://gcd.cr.yp.to/safegcd-20190413.pdf fig 1.2
91
92 BigInt f = a;
93 BigInt g = b;
96
99
100 const size_t common2s = std::min(low_zero_bits(f), low_zero_bits(g));
101 CT::unpoison(common2s);
102
103 f >>= common2s;
104 g >>= common2s;
105
106 f.ct_cond_swap(f.is_even(), g);
107
108 int32_t delta = 1;
109
110 const size_t loop_cnt = safegcd_loop_bound(f.bits(), g.bits());
111
112 BigInt newg, t;
113 for(size_t i = 0; i != loop_cnt; ++i)
114 {
115 sub_abs(newg, f, g);
116
117 const bool need_swap = (g.is_odd() && delta > 0);
118
119 // if(need_swap) { delta *= -1 } else { delta *= 1 }
120 delta *= CT::Mask<uint8_t>::expand(need_swap).if_not_set_return(2) - 1;
121 f.ct_cond_swap(need_swap, g);
122 g.ct_cond_swap(need_swap, newg);
123
124 delta += 1;
125
126 g.ct_cond_add(g.is_odd(), f);
127 g >>= 1;
128 }
129
130 f <<= common2s;
131
134
136
137 return f;
138 }
139
140/*
141* Calculate the LCM
142*/
143BigInt lcm(const BigInt& a, const BigInt& b)
144 {
145 return ct_divide(a * b, gcd(a, b));
146 }
147
148/*
149* Modular Exponentiation
150*/
151BigInt power_mod(const BigInt& base, const BigInt& exp, const BigInt& mod)
152 {
153 if(mod.is_negative() || mod == 1)
154 {
155 return 0;
156 }
157
158 if(base.is_zero() || mod.is_zero())
159 {
160 if(exp.is_zero())
161 return 1;
162 return 0;
163 }
164
165 Modular_Reducer reduce_mod(mod);
166
167 const size_t exp_bits = exp.bits();
168
169 if(mod.is_odd())
170 {
171 const size_t powm_window = 4;
172
173 auto monty_mod = std::make_shared<Montgomery_Params>(mod, reduce_mod);
174 auto powm_base_mod = monty_precompute(monty_mod, reduce_mod.reduce(base), powm_window);
175 return monty_execute(*powm_base_mod, exp, exp_bits);
176 }
177
178 /*
179 Support for even modulus is just a convenience and not considered
180 cryptographically important, so this implementation is slow ...
181 */
182 BigInt accum = 1;
183 BigInt g = reduce_mod.reduce(base);
184 BigInt t;
185
186 for(size_t i = 0; i != exp_bits; ++i)
187 {
188 t = reduce_mod.multiply(g, accum);
189 g = reduce_mod.square(g);
190 accum.ct_cond_assign(exp.get_bit(i), t);
191 }
192 return accum;
193 }
194
195
197 {
198 if(C < 1)
199 throw Invalid_Argument("is_perfect_square requires C >= 1");
200 if(C == 1)
201 return 1;
202
203 const size_t n = C.bits();
204 const size_t m = (n + 1) / 2;
205 const BigInt B = C + BigInt::power_of_2(m);
206
207 BigInt X = BigInt::power_of_2(m) - 1;
208 BigInt X2 = (X*X);
209
210 for(;;)
211 {
212 X = (X2 + C) / (2*X);
213 X2 = (X*X);
214
215 if(X2 < B)
216 break;
217 }
218
219 if(X2 == C)
220 return X;
221 else
222 return 0;
223 }
224
225/*
226* Test for primality using Miller-Rabin
227*/
228bool is_prime(const BigInt& n,
230 size_t prob,
231 bool is_random)
232 {
233 if(n == 2)
234 return true;
235 if(n <= 1 || n.is_even())
236 return false;
237
238 const size_t n_bits = n.bits();
239
240 // Fast path testing for small numbers (<= 65521)
241 if(n_bits <= 16)
242 {
243 const uint16_t num = static_cast<uint16_t>(n.word_at(0));
244
245 return std::binary_search(PRIMES, PRIMES + PRIME_TABLE_SIZE, num);
246 }
247
248 Modular_Reducer mod_n(n);
249
250 if(rng.is_seeded())
251 {
252 const size_t t = miller_rabin_test_iterations(n_bits, prob, is_random);
253
254 if(is_miller_rabin_probable_prime(n, mod_n, rng, t) == false)
255 return false;
256
257 if(is_random)
258 return true;
259 else
260 return is_lucas_probable_prime(n, mod_n);
261 }
262 else
263 {
264 return is_bailie_psw_probable_prime(n, mod_n);
265 }
266 }
267
268}
#define BOTAN_ASSERT_NOMSG(expr)
Definition: assert.h:68
void ct_cond_add(bool predicate, const BigInt &value)
Definition: bigint.cpp:455
bool is_odd() const
Definition: bigint.h:409
void ct_cond_assign(bool predicate, const BigInt &other)
Definition: bigint.cpp:488
size_t size() const
Definition: bigint.h:580
word word_at(size_t n) const
Definition: bigint.h:508
size_t bits() const
Definition: bigint.cpp:296
bool is_even() const
Definition: bigint.h:403
static BigInt power_of_2(size_t n)
Definition: bigint.h:758
void const_time_poison() const
Definition: bigint.h:739
void ct_cond_swap(bool predicate, BigInt &other)
Definition: bigint.cpp:466
bool is_zero() const
Definition: bigint.h:421
bool is_negative() const
Definition: bigint.h:527
void const_time_unpoison() const
Definition: bigint.h:740
bool get_bit(size_t n) const
Definition: bigint.h:465
void set_sign(Sign sign)
Definition: bigint.h:563
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
void unpoison(const T *p, size_t n)
Definition: ct_utils.h:59
Definition: alg_id.cpp:13
BigInt power_mod(const BigInt &base, const BigInt &exp, const BigInt &mod)
Definition: numthry.cpp:151
BigInt lcm(const BigInt &a, const BigInt &b)
Definition: numthry.cpp:143
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:39
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:228
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:81
BigInt is_perfect_square(const BigInt &C)
Definition: numthry.cpp:196
size_t miller_rabin_test_iterations(size_t n_bits, size_t prob, bool random)
Definition: primality.cpp:165
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