Botan 2.17.3
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
64/*
65* Calculate the GCD
66*/
67BigInt gcd(const BigInt& a, const BigInt& b)
68 {
69 if(a.is_zero() || b.is_zero())
70 return 0;
71 if(a == 1 || b == 1)
72 return 1;
73
74 // See https://gcd.cr.yp.to/safegcd-20190413.pdf fig 1.2
75
76 BigInt f = a;
77 BigInt g = b;
80
83
84 const size_t common2s = std::min(low_zero_bits(f), low_zero_bits(g));
85 CT::unpoison(common2s);
86
87 f >>= common2s;
88 g >>= common2s;
89
90 f.ct_cond_swap(f.is_even(), g);
91
92 int32_t delta = 1;
93
94 const size_t loop_cnt = 4 + 3*std::max(f.bits(), g.bits());
95
96 BigInt newg, t;
97 for(size_t i = 0; i != loop_cnt; ++i)
98 {
99 sub_abs(newg, f, g);
100
101 const bool need_swap = (g.is_odd() && delta > 0);
102
103 // if(need_swap) delta *= -1
104 delta *= CT::Mask<uint8_t>::expand(need_swap).select(0, 2) - 1;
105 f.ct_cond_swap(need_swap, g);
106 g.ct_cond_swap(need_swap, newg);
107
108 delta += 1;
109
110 g.ct_cond_add(g.is_odd(), f);
111 g >>= 1;
112 }
113
114 f <<= common2s;
115
118
119 return f;
120 }
121
122/*
123* Calculate the LCM
124*/
125BigInt lcm(const BigInt& a, const BigInt& b)
126 {
127 return ct_divide(a * b, gcd(a, b));
128 }
129
130/*
131* Modular Exponentiation
132*/
133BigInt power_mod(const BigInt& base, const BigInt& exp, const BigInt& mod)
134 {
135 if(mod.is_negative() || mod == 1)
136 {
137 return 0;
138 }
139
140 if(base.is_zero() || mod.is_zero())
141 {
142 if(exp.is_zero())
143 return 1;
144 return 0;
145 }
146
147 Modular_Reducer reduce_mod(mod);
148
149 const size_t exp_bits = exp.bits();
150
151 if(mod.is_odd())
152 {
153 const size_t powm_window = 4;
154
155 auto monty_mod = std::make_shared<Montgomery_Params>(mod, reduce_mod);
156 auto powm_base_mod = monty_precompute(monty_mod, reduce_mod.reduce(base), powm_window);
157 return monty_execute(*powm_base_mod, exp, exp_bits);
158 }
159
160 /*
161 Support for even modulus is just a convenience and not considered
162 cryptographically important, so this implementation is slow ...
163 */
164 BigInt accum = 1;
165 BigInt g = reduce_mod.reduce(base);
166 BigInt t;
167
168 for(size_t i = 0; i != exp_bits; ++i)
169 {
170 t = reduce_mod.multiply(g, accum);
171 g = reduce_mod.square(g);
172 accum.ct_cond_assign(exp.get_bit(i), t);
173 }
174 return accum;
175 }
176
177
179 {
180 if(C < 1)
181 throw Invalid_Argument("is_perfect_square requires C >= 1");
182 if(C == 1)
183 return 1;
184
185 const size_t n = C.bits();
186 const size_t m = (n + 1) / 2;
187 const BigInt B = C + BigInt::power_of_2(m);
188
189 BigInt X = BigInt::power_of_2(m) - 1;
190 BigInt X2 = (X*X);
191
192 for(;;)
193 {
194 X = (X2 + C) / (2*X);
195 X2 = (X*X);
196
197 if(X2 < B)
198 break;
199 }
200
201 if(X2 == C)
202 return X;
203 else
204 return 0;
205 }
206
207/*
208* Test for primality using Miller-Rabin
209*/
210bool is_prime(const BigInt& n,
212 size_t prob,
213 bool is_random)
214 {
215 if(n == 2)
216 return true;
217 if(n <= 1 || n.is_even())
218 return false;
219
220 const size_t n_bits = n.bits();
221
222 // Fast path testing for small numbers (<= 65521)
223 if(n_bits <= 16)
224 {
225 const uint16_t num = static_cast<uint16_t>(n.word_at(0));
226
227 return std::binary_search(PRIMES, PRIMES + PRIME_TABLE_SIZE, num);
228 }
229
230 Modular_Reducer mod_n(n);
231
232 if(rng.is_seeded())
233 {
234 const size_t t = miller_rabin_test_iterations(n_bits, prob, is_random);
235
236 if(is_miller_rabin_probable_prime(n, mod_n, rng, t) == false)
237 return false;
238
239 if(is_random)
240 return true;
241 else
242 return is_lucas_probable_prime(n, mod_n);
243 }
244 else
245 {
246 return is_bailie_psw_probable_prime(n, mod_n);
247 }
248 }
249
250}
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:133
BigInt lcm(const BigInt &a, const BigInt &b)
Definition: numthry.cpp:125
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
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:210
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:67
BigInt is_perfect_square(const BigInt &C)
Definition: numthry.cpp:178
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