Botan  2.19.1
Crypto and TLS for C++11
primality.cpp
Go to the documentation of this file.
1 /*
2 * (C) 2016,2018 Jack Lloyd
3 *
4 * Botan is released under the Simplified BSD License (see license.txt)
5 */
6 
7 #include <botan/internal/primality.h>
8 #include <botan/internal/monty_exp.h>
9 #include <botan/bigint.h>
10 #include <botan/monty.h>
11 #include <botan/reducer.h>
12 #include <botan/rng.h>
13 #include <algorithm>
14 
15 namespace Botan {
16 
17 bool is_lucas_probable_prime(const BigInt& C, const Modular_Reducer& mod_C)
18  {
19  if(C <= 1)
20  return false;
21  else if(C == 2)
22  return true;
23  else if(C.is_even())
24  return false;
25  else if(C == 3 || C == 5 || C == 7 || C == 11 || C == 13)
26  return true;
27 
28  BigInt D = 5;
29 
30  for(;;)
31  {
32  int32_t j = jacobi(D, C);
33  if(j == 0)
34  return false;
35 
36  if(j == -1)
37  break;
38 
39  // Check 5, -7, 9, -11, 13, -15, 17, ...
40  if(D.is_negative())
41  {
42  D.flip_sign();
43  D += 2;
44  }
45  else
46  {
47  D += 2;
48  D.flip_sign();
49  }
50 
51  if(D == 17 && is_perfect_square(C).is_nonzero())
52  return false;
53  }
54 
55  const BigInt K = C + 1;
56  const size_t K_bits = K.bits() - 1;
57 
58  BigInt U = 1;
59  BigInt V = 1;
60 
61  BigInt Ut, Vt, U2, V2;
62 
63  for(size_t i = 0; i != K_bits; ++i)
64  {
65  const bool k_bit = K.get_bit(K_bits - 1 - i);
66 
67  Ut = mod_C.multiply(U, V);
68 
69  Vt = mod_C.reduce(mod_C.square(V) + mod_C.multiply(D, mod_C.square(U)));
70  Vt.ct_cond_add(Vt.is_odd(), C);
71  Vt >>= 1;
72  Vt = mod_C.reduce(Vt);
73 
74  U = Ut;
75  V = Vt;
76 
77  U2 = mod_C.reduce(Ut + Vt);
78  U2.ct_cond_add(U2.is_odd(), C);
79  U2 >>= 1;
80 
81  V2 = mod_C.reduce(Vt + Ut*D);
82  V2.ct_cond_add(V2.is_odd(), C);
83  V2 >>= 1;
84 
85  U.ct_cond_assign(k_bit, U2);
86  V.ct_cond_assign(k_bit, V2);
87  }
88 
89  return (U == 0);
90  }
91 
93  {
94  auto monty_n = std::make_shared<Montgomery_Params>(n, mod_n);
95  return passes_miller_rabin_test(n, mod_n, monty_n, 2) && is_lucas_probable_prime(n, mod_n);
96  }
97 
99  {
100  Modular_Reducer mod_n(n);
101  return is_bailie_psw_probable_prime(n, mod_n);
102  }
103 
105  const Modular_Reducer& mod_n,
106  const std::shared_ptr<Montgomery_Params>& monty_n,
107  const BigInt& a)
108  {
109  BOTAN_ASSERT_NOMSG(n > 1);
110 
111  const BigInt n_minus_1 = n - 1;
112  const size_t s = low_zero_bits(n_minus_1);
113  const BigInt nm1_s = n_minus_1 >> s;
114  const size_t n_bits = n.bits();
115 
116  const size_t powm_window = 4;
117 
118  auto powm_a_n = monty_precompute(monty_n, a, powm_window);
119 
120  BigInt y = monty_execute(*powm_a_n, nm1_s, n_bits);
121 
122  if(y == 1 || y == n_minus_1)
123  return true;
124 
125  for(size_t i = 1; i != s; ++i)
126  {
127  y = mod_n.square(y);
128 
129  if(y == 1) // found a non-trivial square root
130  return false;
131 
132  /*
133  -1 is the trivial square root of unity, so ``a`` is not a
134  witness for this number - give up
135  */
136  if(y == n_minus_1)
137  return true;
138  }
139 
140  return false;
141  }
142 
144  const Modular_Reducer& mod_n,
146  size_t test_iterations)
147  {
148  BOTAN_ASSERT_NOMSG(n > 1);
149 
150  auto monty_n = std::make_shared<Montgomery_Params>(n, mod_n);
151 
152  for(size_t i = 0; i != test_iterations; ++i)
153  {
154  const BigInt a = BigInt::random_integer(rng, 2, n);
155 
156  if(!passes_miller_rabin_test(n, mod_n, monty_n, a))
157  return false;
158  }
159 
160  // Failed to find a counterexample
161  return true;
162  }
163 
164 
165 size_t miller_rabin_test_iterations(size_t n_bits, size_t prob, bool random)
166  {
167  const size_t base = (prob + 2) / 2; // worst case 4^-t error rate
168 
169  /*
170  * If the candidate prime was maliciously constructed, we can't rely
171  * on arguments based on p being random.
172  */
173  if(random == false)
174  return base;
175 
176  /*
177  * For randomly chosen numbers we can use the estimates from
178  * http://www.math.dartmouth.edu/~carlp/PDF/paper88.pdf
179  *
180  * These values are derived from the inequality for p(k,t) given on
181  * the second page.
182  */
183  if(prob <= 128)
184  {
185  if(n_bits >= 1536)
186  return 4; // < 2^-133
187  if(n_bits >= 1024)
188  return 6; // < 2^-133
189  if(n_bits >= 512)
190  return 12; // < 2^-129
191  if(n_bits >= 256)
192  return 29; // < 2^-128
193  }
194 
195  /*
196  If the user desires a smaller error probability than we have
197  precomputed error estimates for, just fall back to using the worst
198  case error rate.
199  */
200  return base;
201  }
202 
203 }
bool is_lucas_probable_prime(const BigInt &C, const Modular_Reducer &mod_C)
Definition: primality.cpp:17
bool is_odd() const
Definition: bigint.h:409
bool is_even() const
Definition: bigint.h:403
size_t low_zero_bits(const BigInt &n)
Definition: numthry.cpp:39
size_t miller_rabin_test_iterations(size_t n_bits, size_t prob, bool random)
Definition: primality.cpp:165
SIMD_8x32 D
#define BOTAN_ASSERT_NOMSG(expr)
Definition: assert.h:68
bool is_negative() const
Definition: bigint.h:527
static BigInt random_integer(RandomNumberGenerator &rng, const BigInt &min, const BigInt &max)
Definition: big_rand.cpp:45
bool is_bailie_psw_probable_prime(const BigInt &n, const Modular_Reducer &mod_n)
Definition: primality.cpp:92
size_t bits() const
Definition: bigint.cpp:296
BigInt multiply(const BigInt &x, const BigInt &y) const
Definition: reducer.h:31
BigInt is_perfect_square(const BigInt &C)
Definition: numthry.cpp:196
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
Definition: alg_id.cpp:13
BigInt reduce(const BigInt &x) const
Definition: reducer.cpp:37
bool passes_miller_rabin_test(const BigInt &n, const Modular_Reducer &mod_n, const std::shared_ptr< Montgomery_Params > &monty_n, const BigInt &a)
Definition: primality.cpp:104
bool is_miller_rabin_probable_prime(const BigInt &n, const Modular_Reducer &mod_n, RandomNumberGenerator &rng, size_t test_iterations)
Definition: primality.cpp:143
bool get_bit(size_t n) const
Definition: bigint.h:465
void ct_cond_add(bool predicate, const BigInt &value)
Definition: bigint.cpp:455
SIMD_8x32 C
void flip_sign()
Definition: bigint.h:554
BigInt monty_execute(const Montgomery_Exponentation_State &precomputed_state, const BigInt &k, size_t max_k_bits)
Definition: monty_exp.cpp:165
BigInt square(const BigInt &x) const
Definition: reducer.h:39
int32_t jacobi(const BigInt &a, const BigInt &n)
Definition: jacobi.cpp:15
void ct_cond_assign(bool predicate, const BigInt &other)
Definition: bigint.cpp:488