1 /* Written by Lenka Fibikova <fibikova (at) exp-math.uni-essen.de> 2 * and Bodo Moeller for the OpenSSL project. */ 3 /* ==================================================================== 4 * Copyright (c) 1998-2000 The OpenSSL Project. All rights reserved. 5 * 6 * Redistribution and use in source and binary forms, with or without 7 * modification, are permitted provided that the following conditions 8 * are met: 9 * 10 * 1. Redistributions of source code must retain the above copyright 11 * notice, this list of conditions and the following disclaimer. 12 * 13 * 2. Redistributions in binary form must reproduce the above copyright 14 * notice, this list of conditions and the following disclaimer in 15 * the documentation and/or other materials provided with the 16 * distribution. 17 * 18 * 3. All advertising materials mentioning features or use of this 19 * software must display the following acknowledgment: 20 * "This product includes software developed by the OpenSSL Project 21 * for use in the OpenSSL Toolkit. (http://www.openssl.org/)" 22 * 23 * 4. The names "OpenSSL Toolkit" and "OpenSSL Project" must not be used to 24 * endorse or promote products derived from this software without 25 * prior written permission. For written permission, please contact 26 * openssl-core (at) openssl.org. 27 * 28 * 5. Products derived from this software may not be called "OpenSSL" 29 * nor may "OpenSSL" appear in their names without prior written 30 * permission of the OpenSSL Project. 31 * 32 * 6. Redistributions of any form whatsoever must retain the following 33 * acknowledgment: 34 * "This product includes software developed by the OpenSSL Project 35 * for use in the OpenSSL Toolkit (http://www.openssl.org/)" 36 * 37 * THIS SOFTWARE IS PROVIDED BY THE OpenSSL PROJECT ``AS IS'' AND ANY 38 * EXPRESSED OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 39 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR 40 * PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE OpenSSL PROJECT OR 41 * ITS CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, 42 * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT 43 * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; 44 * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 45 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, 46 * STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) 47 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED 48 * OF THE POSSIBILITY OF SUCH DAMAGE. 49 * ==================================================================== 50 * 51 * This product includes cryptographic software written by Eric Young 52 * (eay (at) cryptsoft.com). This product includes software written by Tim 53 * Hudson (tjh (at) cryptsoft.com). */ 54 55 #include <openssl/bn.h> 56 57 #include <openssl/err.h> 58 59 #include "internal.h" 60 61 62 BIGNUM *BN_mod_sqrt(BIGNUM *in, const BIGNUM *a, const BIGNUM *p, BN_CTX *ctx) { 63 /* Compute a square root of |a| mod |p| using the Tonelli/Shanks algorithm 64 * (cf. Henri Cohen, "A Course in Algebraic Computational Number Theory", 65 * algorithm 1.5.1). |p| is assumed to be a prime. */ 66 67 BIGNUM *ret = in; 68 int err = 1; 69 int r; 70 BIGNUM *A, *b, *q, *t, *x, *y; 71 int e, i, j; 72 73 if (!BN_is_odd(p) || BN_abs_is_word(p, 1)) { 74 if (BN_abs_is_word(p, 2)) { 75 if (ret == NULL) { 76 ret = BN_new(); 77 } 78 if (ret == NULL) { 79 goto end; 80 } 81 if (!BN_set_word(ret, BN_is_bit_set(a, 0))) { 82 if (ret != in) { 83 BN_free(ret); 84 } 85 return NULL; 86 } 87 return ret; 88 } 89 90 OPENSSL_PUT_ERROR(BN, BN_R_P_IS_NOT_PRIME); 91 return (NULL); 92 } 93 94 if (BN_is_zero(a) || BN_is_one(a)) { 95 if (ret == NULL) { 96 ret = BN_new(); 97 } 98 if (ret == NULL) { 99 goto end; 100 } 101 if (!BN_set_word(ret, BN_is_one(a))) { 102 if (ret != in) { 103 BN_free(ret); 104 } 105 return NULL; 106 } 107 return ret; 108 } 109 110 BN_CTX_start(ctx); 111 A = BN_CTX_get(ctx); 112 b = BN_CTX_get(ctx); 113 q = BN_CTX_get(ctx); 114 t = BN_CTX_get(ctx); 115 x = BN_CTX_get(ctx); 116 y = BN_CTX_get(ctx); 117 if (y == NULL) { 118 goto end; 119 } 120 121 if (ret == NULL) { 122 ret = BN_new(); 123 } 124 if (ret == NULL) { 125 goto end; 126 } 127 128 /* A = a mod p */ 129 if (!BN_nnmod(A, a, p, ctx)) { 130 goto end; 131 } 132 133 /* now write |p| - 1 as 2^e*q where q is odd */ 134 e = 1; 135 while (!BN_is_bit_set(p, e)) { 136 e++; 137 } 138 /* we'll set q later (if needed) */ 139 140 if (e == 1) { 141 /* The easy case: (|p|-1)/2 is odd, so 2 has an inverse 142 * modulo (|p|-1)/2, and square roots can be computed 143 * directly by modular exponentiation. 144 * We have 145 * 2 * (|p|+1)/4 == 1 (mod (|p|-1)/2), 146 * so we can use exponent (|p|+1)/4, i.e. (|p|-3)/4 + 1. 147 */ 148 if (!BN_rshift(q, p, 2)) { 149 goto end; 150 } 151 q->neg = 0; 152 if (!BN_add_word(q, 1) || 153 !BN_mod_exp_mont(ret, A, q, p, ctx, NULL)) { 154 goto end; 155 } 156 err = 0; 157 goto vrfy; 158 } 159 160 if (e == 2) { 161 /* |p| == 5 (mod 8) 162 * 163 * In this case 2 is always a non-square since 164 * Legendre(2,p) = (-1)^((p^2-1)/8) for any odd prime. 165 * So if a really is a square, then 2*a is a non-square. 166 * Thus for 167 * b := (2*a)^((|p|-5)/8), 168 * i := (2*a)*b^2 169 * we have 170 * i^2 = (2*a)^((1 + (|p|-5)/4)*2) 171 * = (2*a)^((p-1)/2) 172 * = -1; 173 * so if we set 174 * x := a*b*(i-1), 175 * then 176 * x^2 = a^2 * b^2 * (i^2 - 2*i + 1) 177 * = a^2 * b^2 * (-2*i) 178 * = a*(-i)*(2*a*b^2) 179 * = a*(-i)*i 180 * = a. 181 * 182 * (This is due to A.O.L. Atkin, 183 * <URL: 184 *http://listserv.nodak.edu/scripts/wa.exe?A2=ind9211&L=nmbrthry&O=T&P=562>, 185 * November 1992.) 186 */ 187 188 /* t := 2*a */ 189 if (!BN_mod_lshift1_quick(t, A, p)) { 190 goto end; 191 } 192 193 /* b := (2*a)^((|p|-5)/8) */ 194 if (!BN_rshift(q, p, 3)) { 195 goto end; 196 } 197 q->neg = 0; 198 if (!BN_mod_exp_mont(b, t, q, p, ctx, NULL)) { 199 goto end; 200 } 201 202 /* y := b^2 */ 203 if (!BN_mod_sqr(y, b, p, ctx)) { 204 goto end; 205 } 206 207 /* t := (2*a)*b^2 - 1*/ 208 if (!BN_mod_mul(t, t, y, p, ctx) || 209 !BN_sub_word(t, 1)) { 210 goto end; 211 } 212 213 /* x = a*b*t */ 214 if (!BN_mod_mul(x, A, b, p, ctx) || 215 !BN_mod_mul(x, x, t, p, ctx)) { 216 goto end; 217 } 218 219 if (!BN_copy(ret, x)) { 220 goto end; 221 } 222 err = 0; 223 goto vrfy; 224 } 225 226 /* e > 2, so we really have to use the Tonelli/Shanks algorithm. 227 * First, find some y that is not a square. */ 228 if (!BN_copy(q, p)) { 229 goto end; /* use 'q' as temp */ 230 } 231 q->neg = 0; 232 i = 2; 233 do { 234 /* For efficiency, try small numbers first; 235 * if this fails, try random numbers. 236 */ 237 if (i < 22) { 238 if (!BN_set_word(y, i)) { 239 goto end; 240 } 241 } else { 242 if (!BN_pseudo_rand(y, BN_num_bits(p), 0, 0)) { 243 goto end; 244 } 245 if (BN_ucmp(y, p) >= 0) { 246 if (!(p->neg ? BN_add : BN_sub)(y, y, p)) { 247 goto end; 248 } 249 } 250 /* now 0 <= y < |p| */ 251 if (BN_is_zero(y)) { 252 if (!BN_set_word(y, i)) { 253 goto end; 254 } 255 } 256 } 257 258 r = bn_jacobi(y, q, ctx); /* here 'q' is |p| */ 259 if (r < -1) { 260 goto end; 261 } 262 if (r == 0) { 263 /* m divides p */ 264 OPENSSL_PUT_ERROR(BN, BN_R_P_IS_NOT_PRIME); 265 goto end; 266 } 267 } while (r == 1 && ++i < 82); 268 269 if (r != -1) { 270 /* Many rounds and still no non-square -- this is more likely 271 * a bug than just bad luck. 272 * Even if p is not prime, we should have found some y 273 * such that r == -1. 274 */ 275 OPENSSL_PUT_ERROR(BN, BN_R_TOO_MANY_ITERATIONS); 276 goto end; 277 } 278 279 /* Here's our actual 'q': */ 280 if (!BN_rshift(q, q, e)) { 281 goto end; 282 } 283 284 /* Now that we have some non-square, we can find an element 285 * of order 2^e by computing its q'th power. */ 286 if (!BN_mod_exp_mont(y, y, q, p, ctx, NULL)) { 287 goto end; 288 } 289 if (BN_is_one(y)) { 290 OPENSSL_PUT_ERROR(BN, BN_R_P_IS_NOT_PRIME); 291 goto end; 292 } 293 294 /* Now we know that (if p is indeed prime) there is an integer 295 * k, 0 <= k < 2^e, such that 296 * 297 * a^q * y^k == 1 (mod p). 298 * 299 * As a^q is a square and y is not, k must be even. 300 * q+1 is even, too, so there is an element 301 * 302 * X := a^((q+1)/2) * y^(k/2), 303 * 304 * and it satisfies 305 * 306 * X^2 = a^q * a * y^k 307 * = a, 308 * 309 * so it is the square root that we are looking for. 310 */ 311 312 /* t := (q-1)/2 (note that q is odd) */ 313 if (!BN_rshift1(t, q)) { 314 goto end; 315 } 316 317 /* x := a^((q-1)/2) */ 318 if (BN_is_zero(t)) /* special case: p = 2^e + 1 */ 319 { 320 if (!BN_nnmod(t, A, p, ctx)) { 321 goto end; 322 } 323 if (BN_is_zero(t)) { 324 /* special case: a == 0 (mod p) */ 325 BN_zero(ret); 326 err = 0; 327 goto end; 328 } else if (!BN_one(x)) { 329 goto end; 330 } 331 } else { 332 if (!BN_mod_exp_mont(x, A, t, p, ctx, NULL)) { 333 goto end; 334 } 335 if (BN_is_zero(x)) { 336 /* special case: a == 0 (mod p) */ 337 BN_zero(ret); 338 err = 0; 339 goto end; 340 } 341 } 342 343 /* b := a*x^2 (= a^q) */ 344 if (!BN_mod_sqr(b, x, p, ctx) || 345 !BN_mod_mul(b, b, A, p, ctx)) { 346 goto end; 347 } 348 349 /* x := a*x (= a^((q+1)/2)) */ 350 if (!BN_mod_mul(x, x, A, p, ctx)) { 351 goto end; 352 } 353 354 while (1) { 355 /* Now b is a^q * y^k for some even k (0 <= k < 2^E 356 * where E refers to the original value of e, which we 357 * don't keep in a variable), and x is a^((q+1)/2) * y^(k/2). 358 * 359 * We have a*b = x^2, 360 * y^2^(e-1) = -1, 361 * b^2^(e-1) = 1. 362 */ 363 364 if (BN_is_one(b)) { 365 if (!BN_copy(ret, x)) { 366 goto end; 367 } 368 err = 0; 369 goto vrfy; 370 } 371 372 373 /* find smallest i such that b^(2^i) = 1 */ 374 i = 1; 375 if (!BN_mod_sqr(t, b, p, ctx)) { 376 goto end; 377 } 378 while (!BN_is_one(t)) { 379 i++; 380 if (i == e) { 381 OPENSSL_PUT_ERROR(BN, BN_R_NOT_A_SQUARE); 382 goto end; 383 } 384 if (!BN_mod_mul(t, t, t, p, ctx)) { 385 goto end; 386 } 387 } 388 389 390 /* t := y^2^(e - i - 1) */ 391 if (!BN_copy(t, y)) { 392 goto end; 393 } 394 for (j = e - i - 1; j > 0; j--) { 395 if (!BN_mod_sqr(t, t, p, ctx)) { 396 goto end; 397 } 398 } 399 if (!BN_mod_mul(y, t, t, p, ctx) || 400 !BN_mod_mul(x, x, t, p, ctx) || 401 !BN_mod_mul(b, b, y, p, ctx)) { 402 goto end; 403 } 404 e = i; 405 } 406 407 vrfy: 408 if (!err) { 409 /* verify the result -- the input might have been not a square 410 * (test added in 0.9.8) */ 411 412 if (!BN_mod_sqr(x, ret, p, ctx)) { 413 err = 1; 414 } 415 416 if (!err && 0 != BN_cmp(x, A)) { 417 OPENSSL_PUT_ERROR(BN, BN_R_NOT_A_SQUARE); 418 err = 1; 419 } 420 } 421 422 end: 423 if (err) { 424 if (ret != in) { 425 BN_clear_free(ret); 426 } 427 ret = NULL; 428 } 429 BN_CTX_end(ctx); 430 return ret; 431 } 432 433 int BN_sqrt(BIGNUM *out_sqrt, const BIGNUM *in, BN_CTX *ctx) { 434 BIGNUM *estimate, *tmp, *delta, *last_delta, *tmp2; 435 int ok = 0, last_delta_valid = 0; 436 437 if (in->neg) { 438 OPENSSL_PUT_ERROR(BN, BN_R_NEGATIVE_NUMBER); 439 return 0; 440 } 441 if (BN_is_zero(in)) { 442 BN_zero(out_sqrt); 443 return 1; 444 } 445 446 BN_CTX_start(ctx); 447 if (out_sqrt == in) { 448 estimate = BN_CTX_get(ctx); 449 } else { 450 estimate = out_sqrt; 451 } 452 tmp = BN_CTX_get(ctx); 453 last_delta = BN_CTX_get(ctx); 454 delta = BN_CTX_get(ctx); 455 if (estimate == NULL || tmp == NULL || last_delta == NULL || delta == NULL) { 456 OPENSSL_PUT_ERROR(BN, ERR_R_MALLOC_FAILURE); 457 goto err; 458 } 459 460 /* We estimate that the square root of an n-bit number is 2^{n/2}. */ 461 if (!BN_lshift(estimate, BN_value_one(), BN_num_bits(in)/2)) { 462 goto err; 463 } 464 465 /* This is Newton's method for finding a root of the equation |estimate|^2 - 466 * |in| = 0. */ 467 for (;;) { 468 /* |estimate| = 1/2 * (|estimate| + |in|/|estimate|) */ 469 if (!BN_div(tmp, NULL, in, estimate, ctx) || 470 !BN_add(tmp, tmp, estimate) || 471 !BN_rshift1(estimate, tmp) || 472 /* |tmp| = |estimate|^2 */ 473 !BN_sqr(tmp, estimate, ctx) || 474 /* |delta| = |in| - |tmp| */ 475 !BN_sub(delta, in, tmp)) { 476 OPENSSL_PUT_ERROR(BN, ERR_R_BN_LIB); 477 goto err; 478 } 479 480 delta->neg = 0; 481 /* The difference between |in| and |estimate| squared is required to always 482 * decrease. This ensures that the loop always terminates, but I don't have 483 * a proof that it always finds the square root for a given square. */ 484 if (last_delta_valid && BN_cmp(delta, last_delta) >= 0) { 485 break; 486 } 487 488 last_delta_valid = 1; 489 490 tmp2 = last_delta; 491 last_delta = delta; 492 delta = tmp2; 493 } 494 495 if (BN_cmp(tmp, in) != 0) { 496 OPENSSL_PUT_ERROR(BN, BN_R_NOT_A_SQUARE); 497 goto err; 498 } 499 500 ok = 1; 501 502 err: 503 if (ok && out_sqrt == in && !BN_copy(out_sqrt, estimate)) { 504 ok = 0; 505 } 506 BN_CTX_end(ctx); 507 return ok; 508 } 509
