1 /* 2 * Copyright 2000-2019 The OpenSSL Project Authors. All Rights Reserved. 3 * 4 * Licensed under the OpenSSL license (the "License"). You may not use 5 * this file except in compliance with the License. You can obtain a copy 6 * in the file LICENSE in the source distribution or at 7 * https://www.openssl.org/source/license.html 8 */ 9 10 #include "internal/cryptlib.h" 11 #include "bn_local.h" 12 13 BIGNUM *BN_mod_sqrt(BIGNUM *in, const BIGNUM *a, const BIGNUM *p, BN_CTX *ctx) 14 /* 15 * Returns 'ret' such that ret^2 == a (mod p), using the Tonelli/Shanks 16 * algorithm (cf. Henri Cohen, "A Course in Algebraic Computational Number 17 * Theory", algorithm 1.5.1). 'p' must be prime! 18 */ 19 { 20 BIGNUM *ret = in; 21 int err = 1; 22 int r; 23 BIGNUM *A, *b, *q, *t, *x, *y; 24 int e, i, j; 25 26 if (!BN_is_odd(p) || BN_abs_is_word(p, 1)) { 27 if (BN_abs_is_word(p, 2)) { 28 if (ret == NULL) 29 ret = BN_new(); 30 if (ret == NULL) 31 goto end; 32 if (!BN_set_word(ret, BN_is_bit_set(a, 0))) { 33 if (ret != in) 34 BN_free(ret); 35 return NULL; 36 } 37 bn_check_top(ret); 38 return ret; 39 } 40 41 BNerr(BN_F_BN_MOD_SQRT, BN_R_P_IS_NOT_PRIME); 42 return NULL; 43 } 44 45 if (BN_is_zero(a) || BN_is_one(a)) { 46 if (ret == NULL) 47 ret = BN_new(); 48 if (ret == NULL) 49 goto end; 50 if (!BN_set_word(ret, BN_is_one(a))) { 51 if (ret != in) 52 BN_free(ret); 53 return NULL; 54 } 55 bn_check_top(ret); 56 return ret; 57 } 58 59 BN_CTX_start(ctx); 60 A = BN_CTX_get(ctx); 61 b = BN_CTX_get(ctx); 62 q = BN_CTX_get(ctx); 63 t = BN_CTX_get(ctx); 64 x = BN_CTX_get(ctx); 65 y = BN_CTX_get(ctx); 66 if (y == NULL) 67 goto end; 68 69 if (ret == NULL) 70 ret = BN_new(); 71 if (ret == NULL) 72 goto end; 73 74 /* A = a mod p */ 75 if (!BN_nnmod(A, a, p, ctx)) 76 goto end; 77 78 /* now write |p| - 1 as 2^e*q where q is odd */ 79 e = 1; 80 while (!BN_is_bit_set(p, e)) 81 e++; 82 /* we'll set q later (if needed) */ 83 84 if (e == 1) { 85 /*- 86 * The easy case: (|p|-1)/2 is odd, so 2 has an inverse 87 * modulo (|p|-1)/2, and square roots can be computed 88 * directly by modular exponentiation. 89 * We have 90 * 2 * (|p|+1)/4 == 1 (mod (|p|-1)/2), 91 * so we can use exponent (|p|+1)/4, i.e. (|p|-3)/4 + 1. 92 */ 93 if (!BN_rshift(q, p, 2)) 94 goto end; 95 q->neg = 0; 96 if (!BN_add_word(q, 1)) 97 goto end; 98 if (!BN_mod_exp(ret, A, q, p, ctx)) 99 goto end; 100 err = 0; 101 goto vrfy; 102 } 103 104 if (e == 2) { 105 /*- 106 * |p| == 5 (mod 8) 107 * 108 * In this case 2 is always a non-square since 109 * Legendre(2,p) = (-1)^((p^2-1)/8) for any odd prime. 110 * So if a really is a square, then 2*a is a non-square. 111 * Thus for 112 * b := (2*a)^((|p|-5)/8), 113 * i := (2*a)*b^2 114 * we have 115 * i^2 = (2*a)^((1 + (|p|-5)/4)*2) 116 * = (2*a)^((p-1)/2) 117 * = -1; 118 * so if we set 119 * x := a*b*(i-1), 120 * then 121 * x^2 = a^2 * b^2 * (i^2 - 2*i + 1) 122 * = a^2 * b^2 * (-2*i) 123 * = a*(-i)*(2*a*b^2) 124 * = a*(-i)*i 125 * = a. 126 * 127 * (This is due to A.O.L. Atkin, 128 * Subject: Square Roots and Cognate Matters modulo p=8n+5. 129 * URL: https://listserv.nodak.edu/cgi-bin/wa.exe?A2=ind9211&L=NMBRTHRY&P=4026 130 * November 1992.) 131 */ 132 133 /* t := 2*a */ 134 if (!BN_mod_lshift1_quick(t, A, p)) 135 goto end; 136 137 /* b := (2*a)^((|p|-5)/8) */ 138 if (!BN_rshift(q, p, 3)) 139 goto end; 140 q->neg = 0; 141 if (!BN_mod_exp(b, t, q, p, ctx)) 142 goto end; 143 144 /* y := b^2 */ 145 if (!BN_mod_sqr(y, b, p, ctx)) 146 goto end; 147 148 /* t := (2*a)*b^2 - 1 */ 149 if (!BN_mod_mul(t, t, y, p, ctx)) 150 goto end; 151 if (!BN_sub_word(t, 1)) 152 goto end; 153 154 /* x = a*b*t */ 155 if (!BN_mod_mul(x, A, b, p, ctx)) 156 goto end; 157 if (!BN_mod_mul(x, x, t, p, ctx)) 158 goto end; 159 160 if (!BN_copy(ret, x)) 161 goto end; 162 err = 0; 163 goto vrfy; 164 } 165 166 /* 167 * e > 2, so we really have to use the Tonelli/Shanks algorithm. First, 168 * find some y that is not a square. 169 */ 170 if (!BN_copy(q, p)) 171 goto end; /* use 'q' as temp */ 172 q->neg = 0; 173 i = 2; 174 do { 175 /* 176 * For efficiency, try small numbers first; if this fails, try random 177 * numbers. 178 */ 179 if (i < 22) { 180 if (!BN_set_word(y, i)) 181 goto end; 182 } else { 183 if (!BN_priv_rand(y, BN_num_bits(p), 0, 0)) 184 goto end; 185 if (BN_ucmp(y, p) >= 0) { 186 if (!(p->neg ? BN_add : BN_sub) (y, y, p)) 187 goto end; 188 } 189 /* now 0 <= y < |p| */ 190 if (BN_is_zero(y)) 191 if (!BN_set_word(y, i)) 192 goto end; 193 } 194 195 r = BN_kronecker(y, q, ctx); /* here 'q' is |p| */ 196 if (r < -1) 197 goto end; 198 if (r == 0) { 199 /* m divides p */ 200 BNerr(BN_F_BN_MOD_SQRT, BN_R_P_IS_NOT_PRIME); 201 goto end; 202 } 203 } 204 while (r == 1 && ++i < 82); 205 206 if (r != -1) { 207 /* 208 * Many rounds and still no non-square -- this is more likely a bug 209 * than just bad luck. Even if p is not prime, we should have found 210 * some y such that r == -1. 211 */ 212 BNerr(BN_F_BN_MOD_SQRT, BN_R_TOO_MANY_ITERATIONS); 213 goto end; 214 } 215 216 /* Here's our actual 'q': */ 217 if (!BN_rshift(q, q, e)) 218 goto end; 219 220 /* 221 * Now that we have some non-square, we can find an element of order 2^e 222 * by computing its q'th power. 223 */ 224 if (!BN_mod_exp(y, y, q, p, ctx)) 225 goto end; 226 if (BN_is_one(y)) { 227 BNerr(BN_F_BN_MOD_SQRT, BN_R_P_IS_NOT_PRIME); 228 goto end; 229 } 230 231 /*- 232 * Now we know that (if p is indeed prime) there is an integer 233 * k, 0 <= k < 2^e, such that 234 * 235 * a^q * y^k == 1 (mod p). 236 * 237 * As a^q is a square and y is not, k must be even. 238 * q+1 is even, too, so there is an element 239 * 240 * X := a^((q+1)/2) * y^(k/2), 241 * 242 * and it satisfies 243 * 244 * X^2 = a^q * a * y^k 245 * = a, 246 * 247 * so it is the square root that we are looking for. 248 */ 249 250 /* t := (q-1)/2 (note that q is odd) */ 251 if (!BN_rshift1(t, q)) 252 goto end; 253 254 /* x := a^((q-1)/2) */ 255 if (BN_is_zero(t)) { /* special case: p = 2^e + 1 */ 256 if (!BN_nnmod(t, A, p, ctx)) 257 goto end; 258 if (BN_is_zero(t)) { 259 /* special case: a == 0 (mod p) */ 260 BN_zero(ret); 261 err = 0; 262 goto end; 263 } else if (!BN_one(x)) 264 goto end; 265 } else { 266 if (!BN_mod_exp(x, A, t, p, ctx)) 267 goto end; 268 if (BN_is_zero(x)) { 269 /* special case: a == 0 (mod p) */ 270 BN_zero(ret); 271 err = 0; 272 goto end; 273 } 274 } 275 276 /* b := a*x^2 (= a^q) */ 277 if (!BN_mod_sqr(b, x, p, ctx)) 278 goto end; 279 if (!BN_mod_mul(b, b, A, p, ctx)) 280 goto end; 281 282 /* x := a*x (= a^((q+1)/2)) */ 283 if (!BN_mod_mul(x, x, A, p, ctx)) 284 goto end; 285 286 while (1) { 287 /*- 288 * Now b is a^q * y^k for some even k (0 <= k < 2^E 289 * where E refers to the original value of e, which we 290 * don't keep in a variable), and x is a^((q+1)/2) * y^(k/2). 291 * 292 * We have a*b = x^2, 293 * y^2^(e-1) = -1, 294 * b^2^(e-1) = 1. 295 */ 296 297 if (BN_is_one(b)) { 298 if (!BN_copy(ret, x)) 299 goto end; 300 err = 0; 301 goto vrfy; 302 } 303 304 /* find smallest i such that b^(2^i) = 1 */ 305 i = 1; 306 if (!BN_mod_sqr(t, b, p, ctx)) 307 goto end; 308 while (!BN_is_one(t)) { 309 i++; 310 if (i == e) { 311 BNerr(BN_F_BN_MOD_SQRT, BN_R_NOT_A_SQUARE); 312 goto end; 313 } 314 if (!BN_mod_mul(t, t, t, p, ctx)) 315 goto end; 316 } 317 318 /* t := y^2^(e - i - 1) */ 319 if (!BN_copy(t, y)) 320 goto end; 321 for (j = e - i - 1; j > 0; j--) { 322 if (!BN_mod_sqr(t, t, p, ctx)) 323 goto end; 324 } 325 if (!BN_mod_mul(y, t, t, p, ctx)) 326 goto end; 327 if (!BN_mod_mul(x, x, t, p, ctx)) 328 goto end; 329 if (!BN_mod_mul(b, b, y, p, ctx)) 330 goto end; 331 e = i; 332 } 333 334 vrfy: 335 if (!err) { 336 /* 337 * verify the result -- the input might have been not a square (test 338 * added in 0.9.8) 339 */ 340 341 if (!BN_mod_sqr(x, ret, p, ctx)) 342 err = 1; 343 344 if (!err && 0 != BN_cmp(x, A)) { 345 BNerr(BN_F_BN_MOD_SQRT, BN_R_NOT_A_SQUARE); 346 err = 1; 347 } 348 } 349 350 end: 351 if (err) { 352 if (ret != in) 353 BN_clear_free(ret); 354 ret = NULL; 355 } 356 BN_CTX_end(ctx); 357 bn_check_top(ret); 358 return ret; 359 } 360