1 /* mpn_sec_invert 2 3 Contributed to the GNU project by Niels Möller 4 5 Copyright 2013 Free Software Foundation, Inc. 6 7 This file is part of the GNU MP Library. 8 9 The GNU MP Library is free software; you can redistribute it and/or modify 10 it under the terms of either: 11 12 * the GNU Lesser General Public License as published by the Free 13 Software Foundation; either version 3 of the License, or (at your 14 option) any later version. 15 16 or 17 18 * the GNU General Public License as published by the Free Software 19 Foundation; either version 2 of the License, or (at your option) any 20 later version. 21 22 or both in parallel, as here. 23 24 The GNU MP Library is distributed in the hope that it will be useful, but 25 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY 26 or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License 27 for more details. 28 29 You should have received copies of the GNU General Public License and the 30 GNU Lesser General Public License along with the GNU MP Library. If not, 31 see https://www.gnu.org/licenses/. */ 32 33 #include "gmp.h" 34 #include "gmp-impl.h" 35 36 #if 0 37 /* Currently unused. Should be resurrected once mpn_cnd_neg is 38 advertised. */ 39 static mp_size_t 40 mpn_cnd_neg_itch (mp_size_t n) 41 { 42 return n; 43 } 44 #endif 45 46 /* FIXME: Ought to return carry */ 47 static void 48 mpn_cnd_neg (int cnd, mp_limb_t *rp, const mp_limb_t *ap, mp_size_t n, 49 mp_ptr scratch) 50 { 51 mpn_lshift (scratch, ap, n, 1); 52 mpn_cnd_sub_n (cnd, rp, ap, scratch, n); 53 } 54 55 static int 56 mpn_sec_eq_ui (mp_srcptr ap, mp_size_t n, mp_limb_t b) 57 { 58 mp_limb_t d; 59 ASSERT (n > 0); 60 61 d = ap[0] ^ b; 62 63 while (--n > 0) 64 d |= ap[n]; 65 66 return d == 0; 67 } 68 69 mp_size_t 70 mpn_sec_invert_itch (mp_size_t n) 71 { 72 return 4*n; 73 } 74 75 /* Compute V <-- A^{-1} (mod M), in data-independent time. M must be 76 odd. Returns 1 on success, and 0 on failure (i.e., if gcd (A, m) != 77 1). Inputs and outputs of size n, and no overlap allowed. The {ap, 78 n} area is destroyed. For arbitrary inputs, bit_size should be 79 2*n*GMP_NUMB_BITS, but if A or M are known to be smaller, e.g., if 80 M = 2^521 - 1 and A < M, bit_size can be any bound on the sum of 81 the bit sizes of A and M. */ 82 int 83 mpn_sec_invert (mp_ptr vp, mp_ptr ap, mp_srcptr mp, 84 mp_size_t n, mp_bitcnt_t bit_size, 85 mp_ptr scratch) 86 { 87 ASSERT (n > 0); 88 ASSERT (bit_size > 0); 89 ASSERT (mp[0] & 1); 90 ASSERT (! MPN_OVERLAP_P (ap, n, vp, n)); 91 #define bp (scratch + n) 92 #define up (scratch + 2*n) 93 #define m1hp (scratch + 3*n) 94 95 /* Maintain 96 97 a = u * orig_a (mod m) 98 b = v * orig_a (mod m) 99 100 and b odd at all times. Initially, 101 102 a = a_orig, u = 1 103 b = m, v = 0 104 */ 105 106 107 up[0] = 1; 108 mpn_zero (up+1, n - 1); 109 mpn_copyi (bp, mp, n); 110 mpn_zero (vp, n); 111 112 ASSERT_CARRY (mpn_rshift (m1hp, mp, n, 1)); 113 ASSERT_NOCARRY (mpn_sec_add_1 (m1hp, m1hp, n, 1, scratch)); 114 115 while (bit_size-- > 0) 116 { 117 mp_limb_t odd, swap, cy; 118 119 /* Always maintain b odd. The logic of the iteration is as 120 follows. For a, b: 121 122 odd = a & 1 123 a -= odd * b 124 if (underflow from a-b) 125 { 126 b += a, assigns old a 127 a = B^n-a 128 } 129 130 a /= 2 131 132 For u, v: 133 134 if (underflow from a - b) 135 swap u, v 136 u -= odd * v 137 if (underflow from u - v) 138 u += m 139 140 u /= 2 141 if (a one bit was shifted out) 142 u += (m+1)/2 143 144 As long as a > 0, the quantity 145 146 (bitsize of a) + (bitsize of b) 147 148 is reduced by at least one bit per iteration, hence after (bit_size of 149 orig_a) + (bit_size of m) - 1 iterations we surely have a = 0. Then b 150 = gcd(orig_a, m) and if b = 1 then also v = orig_a^{-1} (mod m). 151 */ 152 153 ASSERT (bp[0] & 1); 154 odd = ap[0] & 1; 155 156 swap = mpn_cnd_sub_n (odd, ap, ap, bp, n); 157 mpn_cnd_add_n (swap, bp, bp, ap, n); 158 mpn_cnd_neg (swap, ap, ap, n, scratch); 159 160 mpn_cnd_swap (swap, up, vp, n); 161 cy = mpn_cnd_sub_n (odd, up, up, vp, n); 162 cy -= mpn_cnd_add_n (cy, up, up, mp, n); 163 ASSERT (cy == 0); 164 165 cy = mpn_rshift (ap, ap, n, 1); 166 ASSERT (cy == 0); 167 cy = mpn_rshift (up, up, n, 1); 168 cy = mpn_cnd_add_n (cy, up, up, m1hp, n); 169 ASSERT (cy == 0); 170 } 171 /* Should be all zeros, but check only extreme limbs */ 172 ASSERT ( (ap[0] | ap[n-1]) == 0); 173 /* Check if indeed gcd == 1. */ 174 return mpn_sec_eq_ui (bp, n, 1); 175 #undef bp 176 #undef up 177 #undef m1hp 178 } 179