Coverage Report

Created: 2021-05-04 09:02

/src/botan/src/lib/math/numbertheory/mod_inv.cpp
Line
Count
Source (jump to first uncovered line)
1
/*
2
* (C) 1999-2011,2016,2018,2019,2020 Jack Lloyd
3
*
4
* Botan is released under the Simplified BSD License (see license.txt)
5
*/
6
7
#include <botan/numthry.h>
8
#include <botan/internal/divide.h>
9
#include <botan/internal/ct_utils.h>
10
#include <botan/internal/mp_core.h>
11
#include <botan/internal/rounding.h>
12
13
namespace Botan {
14
15
namespace {
16
17
BigInt inverse_mod_odd_modulus(const BigInt& n, const BigInt& mod)
18
58.9k
   {
19
   // Caller should assure these preconditions:
20
58.9k
   BOTAN_DEBUG_ASSERT(n.is_positive());
21
58.9k
   BOTAN_DEBUG_ASSERT(mod.is_positive());
22
58.9k
   BOTAN_DEBUG_ASSERT(n < mod);
23
58.9k
   BOTAN_DEBUG_ASSERT(mod >= 3 && mod.is_odd());
24
25
   /*
26
   This uses a modular inversion algorithm designed by Niels Möller
27
   and implemented in Nettle. The same algorithm was later also
28
   adapted to GMP in mpn_sec_invert.
29
30
   It can be easily implemented in a way that does not depend on
31
   secret branches or memory lookups, providing resistance against
32
   some forms of side channel attack.
33
34
   There is also a description of the algorithm in Appendix 5 of "Fast
35
   Software Polynomial Multiplication on ARM Processors using the NEON Engine"
36
   by Danilo Câmara, Conrado P. L. Gouvêa, Julio López, and Ricardo
37
   Dahab in LNCS 8182
38
      https://conradoplg.cryptoland.net/files/2010/12/mocrysen13.pdf
39
40
   Thanks to Niels for creating the algorithm, explaining some things
41
   about it, and the reference to the paper.
42
   */
43
44
58.9k
   const size_t mod_words = mod.sig_words();
45
58.9k
   BOTAN_ASSERT(mod_words > 0, "Not empty");
46
47
58.9k
   secure_vector<word> tmp_mem(5*mod_words);
48
49
58.9k
   word* v_w = &tmp_mem[0];
50
58.9k
   word* u_w = &tmp_mem[1*mod_words];
51
58.9k
   word* b_w = &tmp_mem[2*mod_words];
52
58.9k
   word* a_w = &tmp_mem[3*mod_words];
53
58.9k
   word* mp1o2 = &tmp_mem[4*mod_words];
54
55
58.9k
   CT::poison(tmp_mem.data(), tmp_mem.size());
56
57
58.9k
   copy_mem(a_w, n.data(), std::min(n.size(), mod_words));
58
58.9k
   copy_mem(b_w, mod.data(), std::min(mod.size(), mod_words));
59
58.9k
   u_w[0] = 1;
60
   // v_w = 0
61
62
   // compute (mod + 1) / 2 which [because mod is odd] is equal to
63
   // (mod / 2) + 1
64
58.9k
   copy_mem(mp1o2, mod.data(), std::min(mod.size(), mod_words));
65
58.9k
   bigint_shr1(mp1o2, mod_words, 0, 1);
66
58.9k
   word carry = bigint_add2_nc(mp1o2, mod_words, u_w, 1);
67
58.9k
   BOTAN_ASSERT_NOMSG(carry == 0);
68
69
   // Only n.bits() + mod.bits() iterations are required, but avoid leaking the size of n
70
58.9k
   const size_t execs = 2 * mod.bits();
71
72
43.2M
   for(size_t i = 0; i != execs; ++i)
73
43.2M
      {
74
43.2M
      const word odd_a = a_w[0] & 1;
75
76
      //if(odd_a) a -= b
77
43.2M
      word underflow = bigint_cnd_sub(odd_a, a_w, b_w, mod_words);
78
79
      //if(underflow) { b -= a; a = abs(a); swap(u, v); }
80
43.2M
      bigint_cnd_add(underflow, b_w, a_w, mod_words);
81
43.2M
      bigint_cnd_abs(underflow, a_w, mod_words);
82
43.2M
      bigint_cnd_swap(underflow, u_w, v_w, mod_words);
83
84
      // a >>= 1
85
43.2M
      bigint_shr1(a_w, mod_words, 0, 1);
86
87
      //if(odd_a) u -= v;
88
43.2M
      word borrow = bigint_cnd_sub(odd_a, u_w, v_w, mod_words);
89
90
      // if(borrow) u += p
91
43.2M
      bigint_cnd_add(borrow, u_w, mod.data(), mod_words);
92
93
43.2M
      const word odd_u = u_w[0] & 1;
94
95
      // u >>= 1
96
43.2M
      bigint_shr1(u_w, mod_words, 0, 1);
97
98
      //if(odd_u) u += mp1o2;
99
43.2M
      bigint_cnd_add(odd_u, u_w, mp1o2, mod_words);
100
43.2M
      }
101
102
58.9k
   auto a_is_0 = CT::Mask<word>::set();
103
401k
   for(size_t i = 0; i != mod_words; ++i)
104
342k
      a_is_0 &= CT::Mask<word>::is_zero(a_w[i]);
105
106
58.9k
   auto b_is_1 = CT::Mask<word>::is_equal(b_w[0], 1);
107
342k
   for(size_t i = 1; i != mod_words; ++i)
108
283k
      b_is_1 &= CT::Mask<word>::is_zero(b_w[i]);
109
110
58.9k
   BOTAN_ASSERT(a_is_0.is_set(), "A is zero");
111
112
   // if b != 1 then gcd(n,mod) > 1 and inverse does not exist
113
   // in which case zero out the result to indicate this
114
58.9k
   (~b_is_1).if_set_zero_out(v_w, mod_words);
115
116
   /*
117
   * We've placed the result in the lowest words of the temp buffer.
118
   * So just clear out the other values and then give that buffer to a
119
   * BigInt.
120
   */
121
58.9k
   clear_mem(&tmp_mem[mod_words], 4*mod_words);
122
123
58.9k
   CT::unpoison(tmp_mem.data(), tmp_mem.size());
124
125
58.9k
   BigInt r;
126
58.9k
   r.swap_reg(tmp_mem);
127
58.9k
   return r;
128
58.9k
   }
129
130
BigInt inverse_mod_pow2(const BigInt& a1, size_t k)
131
1.15k
   {
132
   /*
133
   * From "A New Algorithm for Inversion mod p^k" by Çetin Kaya Koç
134
   * https://eprint.iacr.org/2017/411.pdf sections 5 and 7.
135
   */
136
137
1.15k
   if(a1.is_even() || k == 0)
138
0
      return BigInt::zero();
139
1.15k
   if(k == 1)
140
7
      return BigInt::one();
141
142
1.14k
   BigInt a = a1;
143
1.14k
   a.mask_bits(k);
144
145
1.14k
   BigInt b = BigInt::one();
146
1.14k
   BigInt X = BigInt::zero();
147
1.14k
   BigInt newb;
148
149
1.14k
   const size_t a_words = a.sig_words();
150
151
1.14k
   X.grow_to(round_up(k, BOTAN_MP_WORD_BITS) / BOTAN_MP_WORD_BITS);
152
1.14k
   b.grow_to(a_words);
153
154
   /*
155
   Hide the exact value of k. k is anyway known to word length
156
   granularity because of the length of a, so no point in doing more
157
   than this.
158
   */
159
1.14k
   const size_t iter = round_up(k, BOTAN_MP_WORD_BITS);
160
161
676k
   for(size_t i = 0; i != iter; ++i)
162
675k
      {
163
675k
      const bool b0 = b.get_bit(0);
164
675k
      X.conditionally_set_bit(i, b0);
165
675k
      newb = b - a;
166
675k
      b.ct_cond_assign(b0, newb);
167
675k
      b >>= 1;
168
675k
      }
169
170
1.14k
   X.mask_bits(k);
171
1.14k
   X.const_time_unpoison();
172
1.14k
   return X;
173
1.14k
   }
174
175
}
176
177
BigInt inverse_mod(const BigInt& n, const BigInt& mod)
178
59.0k
   {
179
59.0k
   if(mod.is_zero())
180
0
      throw Invalid_Argument("inverse_mod modulus cannot be zero");
181
59.0k
   if(mod.is_negative() || n.is_negative())
182
0
      throw Invalid_Argument("inverse_mod: arguments must be non-negative");
183
59.0k
   if(n.is_zero() || (n.is_even() && mod.is_even()))
184
9
      return BigInt::zero();
185
186
59.0k
   if(mod.is_odd())
187
58.2k
      {
188
      /*
189
      Fastpath for common case. This leaks if n is greater than mod or
190
      not, but we don't guarantee const time behavior in that case.
191
      */
192
58.2k
      if(n < mod)
193
58.1k
         return inverse_mod_odd_modulus(n, mod);
194
78
      else
195
78
         return inverse_mod_odd_modulus(ct_modulo(n, mod), mod);
196
760
      }
197
198
   // If n is even and mod is even we already returned 0
199
   // If n is even and mod is odd we jumped directly to odd-modulus algo
200
760
   BOTAN_DEBUG_ASSERT(n.is_odd());
201
202
760
   const size_t mod_lz = low_zero_bits(mod);
203
760
   BOTAN_ASSERT_NOMSG(mod_lz > 0);
204
760
   const size_t mod_bits = mod.bits();
205
760
   BOTAN_ASSERT_NOMSG(mod_bits > mod_lz);
206
207
760
   if(mod_lz == mod_bits - 1)
208
64
      {
209
      // In this case we are performing an inversion modulo 2^k
210
64
      return inverse_mod_pow2(n, mod_lz);
211
64
      }
212
213
696
   if(mod_lz == 1)
214
103
      {
215
      /*
216
      Inversion modulo 2*o is an easier special case of CRT
217
218
      This is exactly the main CRT flow below but taking advantage of
219
      the fact that any odd number ^-1 modulo 2 is 1. As a result both
220
      inv_2k and c can be taken to be 1, m2k is 2, and h is always
221
      either 0 or 1, and its value depends only on the low bit of inv_o.
222
223
      This is worth special casing because we generate RSA primes such
224
      that phi(n) is of this form. However this only works for keys
225
      that we generated in this way; pre-existing keys will typically
226
      fall back to the general algorithm below.
227
      */
228
229
103
      const BigInt o = mod >> 1;
230
103
      const BigInt n_redc = ct_modulo(n, o);
231
103
      const BigInt inv_o = inverse_mod_odd_modulus(n_redc, o);
232
233
      // No modular inverse in this case:
234
103
      if(inv_o == 0)
235
11
         return BigInt::zero();
236
237
92
      BigInt h = inv_o;
238
92
      h.ct_cond_add(!inv_o.get_bit(0), o);
239
92
      return h;
240
92
      }
241
242
   /*
243
   * In this case we are performing an inversion modulo 2^k*o for
244
   * some k >= 2 and some odd (not necessarily prime) integer.
245
   * Compute the inversions modulo 2^k and modulo o, then combine them
246
   * using CRT, which is possible because 2^k and o are relatively prime.
247
   */
248
249
593
   const BigInt o = mod >> mod_lz;
250
593
   const BigInt n_redc = ct_modulo(n, o);
251
593
   const BigInt inv_o = inverse_mod_odd_modulus(n_redc, o);
252
593
   const BigInt inv_2k = inverse_mod_pow2(n, mod_lz);
253
254
   // No modular inverse in this case:
255
593
   if(inv_o == 0 || inv_2k == 0)
256
96
      return BigInt::zero();
257
258
497
   const BigInt m2k = BigInt::power_of_2(mod_lz);
259
   // Compute the CRT parameter
260
497
   const BigInt c = inverse_mod_pow2(o, mod_lz);
261
262
   // Compute h = c*(inv_2k-inv_o) mod 2^k
263
497
   BigInt h = c * (inv_2k - inv_o);
264
497
   const bool h_neg = h.is_negative();
265
497
   h.set_sign(BigInt::Positive);
266
497
   h.mask_bits(mod_lz);
267
497
   const bool h_nonzero = h.is_nonzero();
268
497
   h.ct_cond_assign(h_nonzero && h_neg, m2k - h);
269
270
   // Return result inv_o + h * o
271
497
   h *= o;
272
497
   h += inv_o;
273
497
   return h;
274
497
   }
275
276
}