1 /* mpn_powm -- Compute R = U^E mod M. 2 3 Contributed to the GNU project by Torbjorn Granlund. 4 5 THE FUNCTIONS IN THIS FILE ARE INTERNAL WITH MUTABLE INTERFACES. IT IS ONLY 6 SAFE TO REACH THEM THROUGH DOCUMENTED INTERFACES. IN FACT, IT IS ALMOST 7 GUARANTEED THAT THEY WILL CHANGE OR DISAPPEAR IN A FUTURE GNU MP RELEASE. 8 9 Copyright 2007-2012 Free Software Foundation, Inc. 10 11 This file is part of the GNU MP Library. 12 13 The GNU MP Library is free software; you can redistribute it and/or modify 14 it under the terms of either: 15 16 * the GNU Lesser General Public License as published by the Free 17 Software Foundation; either version 3 of the License, or (at your 18 option) any later version. 19 20 or 21 22 * the GNU General Public License as published by the Free Software 23 Foundation; either version 2 of the License, or (at your option) any 24 later version. 25 26 or both in parallel, as here. 27 28 The GNU MP Library is distributed in the hope that it will be useful, but 29 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY 30 or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License 31 for more details. 32 33 You should have received copies of the GNU General Public License and the 34 GNU Lesser General Public License along with the GNU MP Library. If not, 35 see https://www.gnu.org/licenses/. */ 36 37 38 /* 39 BASIC ALGORITHM, Compute U^E mod M, where M < B^n is odd. 40 41 1. W <- U 42 43 2. T <- (B^n * U) mod M Convert to REDC form 44 45 3. Compute table U^1, U^3, U^5... of E-dependent size 46 47 4. While there are more bits in E 48 W <- power left-to-right base-k 49 50 51 TODO: 52 53 * Make getbits a macro, thereby allowing it to update the index operand. 54 That will simplify the code using getbits. (Perhaps make getbits' sibling 55 getbit then have similar form, for symmetry.) 56 57 * Write an itch function. Or perhaps get rid of tp parameter since the huge 58 pp area is allocated locally anyway? 59 60 * Choose window size without looping. (Superoptimize or think(tm).) 61 62 * Handle small bases with initial, reduction-free exponentiation. 63 64 * Call new division functions, not mpn_tdiv_qr. 65 66 * Consider special code for one-limb M. 67 68 * How should we handle the redc1/redc2/redc_n choice? 69 - redc1: T(binvert_1limb) + e * (n) * (T(mullo-1x1) + n*T(addmul_1)) 70 - redc2: T(binvert_2limbs) + e * (n/2) * (T(mullo-2x2) + n*T(addmul_2)) 71 - redc_n: T(binvert_nlimbs) + e * (T(mullo-nxn) + T(M(n))) 72 This disregards the addmul_N constant term, but we could think of 73 that as part of the respective mullo. 74 75 * When U (the base) is small, we should start the exponentiation with plain 76 operations, then convert that partial result to REDC form. 77 78 * When U is just one limb, should it be handled without the k-ary tricks? 79 We could keep a factor of B^n in W, but use U' = BU as base. After 80 multiplying by this (pseudo two-limb) number, we need to multiply by 1/B 81 mod M. 82 */ 83 84 #include "gmp.h" 85 #include "gmp-impl.h" 86 #include "longlong.h" 87 88 #undef MPN_REDC_1 89 #define MPN_REDC_1(rp, up, mp, n, invm) \ 90 do { \ 91 mp_limb_t cy; \ 92 cy = mpn_redc_1 (rp, up, mp, n, invm); \ 93 if (cy != 0) \ 94 mpn_sub_n (rp, rp, mp, n); \ 95 } while (0) 96 97 #undef MPN_REDC_2 98 #define MPN_REDC_2(rp, up, mp, n, mip) \ 99 do { \ 100 mp_limb_t cy; \ 101 cy = mpn_redc_2 (rp, up, mp, n, mip); \ 102 if (cy != 0) \ 103 mpn_sub_n (rp, rp, mp, n); \ 104 } while (0) 105 106 #if HAVE_NATIVE_mpn_addmul_2 || HAVE_NATIVE_mpn_redc_2 107 #define WANT_REDC_2 1 108 #endif 109 110 #define getbit(p,bi) \ 111 ((p[(bi - 1) / GMP_LIMB_BITS] >> (bi - 1) % GMP_LIMB_BITS) & 1) 112 113 static inline mp_limb_t 114 getbits (const mp_limb_t *p, mp_bitcnt_t bi, int nbits) 115 { 116 int nbits_in_r; 117 mp_limb_t r; 118 mp_size_t i; 119 120 if (bi < nbits) 121 { 122 return p[0] & (((mp_limb_t) 1 << bi) - 1); 123 } 124 else 125 { 126 bi -= nbits; /* bit index of low bit to extract */ 127 i = bi / GMP_NUMB_BITS; /* word index of low bit to extract */ 128 bi %= GMP_NUMB_BITS; /* bit index in low word */ 129 r = p[i] >> bi; /* extract (low) bits */ 130 nbits_in_r = GMP_NUMB_BITS - bi; /* number of bits now in r */ 131 if (nbits_in_r < nbits) /* did we get enough bits? */ 132 r += p[i + 1] << nbits_in_r; /* prepend bits from higher word */ 133 return r & (((mp_limb_t ) 1 << nbits) - 1); 134 } 135 } 136 137 static inline int 138 win_size (mp_bitcnt_t eb) 139 { 140 int k; 141 static mp_bitcnt_t x[] = {0,7,25,81,241,673,1793,4609,11521,28161,~(mp_bitcnt_t)0}; 142 for (k = 1; eb > x[k]; k++) 143 ; 144 return k; 145 } 146 147 /* Convert U to REDC form, U_r = B^n * U mod M */ 148 static void 149 redcify (mp_ptr rp, mp_srcptr up, mp_size_t un, mp_srcptr mp, mp_size_t n) 150 { 151 mp_ptr tp, qp; 152 TMP_DECL; 153 TMP_MARK; 154 155 TMP_ALLOC_LIMBS_2 (tp, un + n, qp, un + 1); 156 157 MPN_ZERO (tp, n); 158 MPN_COPY (tp + n, up, un); 159 mpn_tdiv_qr (qp, rp, 0L, tp, un + n, mp, n); 160 TMP_FREE; 161 } 162 163 /* rp[n-1..0] = bp[bn-1..0] ^ ep[en-1..0] mod mp[n-1..0] 164 Requires that mp[n-1..0] is odd. 165 Requires that ep[en-1..0] is > 1. 166 Uses scratch space at tp of MAX(mpn_binvert_itch(n),2n) limbs. */ 167 void 168 mpn_powm (mp_ptr rp, mp_srcptr bp, mp_size_t bn, 169 mp_srcptr ep, mp_size_t en, 170 mp_srcptr mp, mp_size_t n, mp_ptr tp) 171 { 172 mp_limb_t ip[2], *mip; 173 int cnt; 174 mp_bitcnt_t ebi; 175 int windowsize, this_windowsize; 176 mp_limb_t expbits; 177 mp_ptr pp, this_pp; 178 long i; 179 TMP_DECL; 180 181 ASSERT (en > 1 || (en == 1 && ep[0] > 1)); 182 ASSERT (n >= 1 && ((mp[0] & 1) != 0)); 183 184 TMP_MARK; 185 186 MPN_SIZEINBASE_2EXP(ebi, ep, en, 1); 187 188 #if 0 189 if (bn < n) 190 { 191 /* Do the first few exponent bits without mod reductions, 192 until the result is greater than the mod argument. */ 193 for (;;) 194 { 195 mpn_sqr (tp, this_pp, tn); 196 tn = tn * 2 - 1, tn += tp[tn] != 0; 197 if (getbit (ep, ebi) != 0) 198 mpn_mul (..., tp, tn, bp, bn); 199 ebi--; 200 } 201 } 202 #endif 203 204 windowsize = win_size (ebi); 205 206 #if WANT_REDC_2 207 if (BELOW_THRESHOLD (n, REDC_1_TO_REDC_2_THRESHOLD)) 208 { 209 mip = ip; 210 binvert_limb (mip[0], mp[0]); 211 mip[0] = -mip[0]; 212 } 213 else if (BELOW_THRESHOLD (n, REDC_2_TO_REDC_N_THRESHOLD)) 214 { 215 mip = ip; 216 mpn_binvert (mip, mp, 2, tp); 217 mip[0] = -mip[0]; mip[1] = ~mip[1]; 218 } 219 #else 220 if (BELOW_THRESHOLD (n, REDC_1_TO_REDC_N_THRESHOLD)) 221 { 222 mip = ip; 223 binvert_limb (mip[0], mp[0]); 224 mip[0] = -mip[0]; 225 } 226 #endif 227 else 228 { 229 mip = TMP_ALLOC_LIMBS (n); 230 mpn_binvert (mip, mp, n, tp); 231 } 232 233 pp = TMP_ALLOC_LIMBS (n << (windowsize - 1)); 234 235 this_pp = pp; 236 redcify (this_pp, bp, bn, mp, n); 237 238 /* Store b^2 at rp. */ 239 mpn_sqr (tp, this_pp, n); 240 #if WANT_REDC_2 241 if (BELOW_THRESHOLD (n, REDC_1_TO_REDC_2_THRESHOLD)) 242 MPN_REDC_1 (rp, tp, mp, n, mip[0]); 243 else if (BELOW_THRESHOLD (n, REDC_2_TO_REDC_N_THRESHOLD)) 244 MPN_REDC_2 (rp, tp, mp, n, mip); 245 #else 246 if (BELOW_THRESHOLD (n, REDC_1_TO_REDC_N_THRESHOLD)) 247 MPN_REDC_1 (rp, tp, mp, n, mip[0]); 248 #endif 249 else 250 mpn_redc_n (rp, tp, mp, n, mip); 251 252 /* Precompute odd powers of b and put them in the temporary area at pp. */ 253 for (i = (1 << (windowsize - 1)) - 1; i > 0; i--) 254 { 255 mpn_mul_n (tp, this_pp, rp, n); 256 this_pp += n; 257 #if WANT_REDC_2 258 if (BELOW_THRESHOLD (n, REDC_1_TO_REDC_2_THRESHOLD)) 259 MPN_REDC_1 (this_pp, tp, mp, n, mip[0]); 260 else if (BELOW_THRESHOLD (n, REDC_2_TO_REDC_N_THRESHOLD)) 261 MPN_REDC_2 (this_pp, tp, mp, n, mip); 262 #else 263 if (BELOW_THRESHOLD (n, REDC_1_TO_REDC_N_THRESHOLD)) 264 MPN_REDC_1 (this_pp, tp, mp, n, mip[0]); 265 #endif 266 else 267 mpn_redc_n (this_pp, tp, mp, n, mip); 268 } 269 270 expbits = getbits (ep, ebi, windowsize); 271 if (ebi < windowsize) 272 ebi = 0; 273 else 274 ebi -= windowsize; 275 276 count_trailing_zeros (cnt, expbits); 277 ebi += cnt; 278 expbits >>= cnt; 279 280 MPN_COPY (rp, pp + n * (expbits >> 1), n); 281 282 #define INNERLOOP \ 283 while (ebi != 0) \ 284 { \ 285 while (getbit (ep, ebi) == 0) \ 286 { \ 287 MPN_SQR (tp, rp, n); \ 288 MPN_REDUCE (rp, tp, mp, n, mip); \ 289 ebi--; \ 290 if (ebi == 0) \ 291 goto done; \ 292 } \ 293 \ 294 /* The next bit of the exponent is 1. Now extract the largest \ 295 block of bits <= windowsize, and such that the least \ 296 significant bit is 1. */ \ 297 \ 298 expbits = getbits (ep, ebi, windowsize); \ 299 this_windowsize = windowsize; \ 300 if (ebi < windowsize) \ 301 { \ 302 this_windowsize -= windowsize - ebi; \ 303 ebi = 0; \ 304 } \ 305 else \ 306 ebi -= windowsize; \ 307 \ 308 count_trailing_zeros (cnt, expbits); \ 309 this_windowsize -= cnt; \ 310 ebi += cnt; \ 311 expbits >>= cnt; \ 312 \ 313 do \ 314 { \ 315 MPN_SQR (tp, rp, n); \ 316 MPN_REDUCE (rp, tp, mp, n, mip); \ 317 this_windowsize--; \ 318 } \ 319 while (this_windowsize != 0); \ 320 \ 321 MPN_MUL_N (tp, rp, pp + n * (expbits >> 1), n); \ 322 MPN_REDUCE (rp, tp, mp, n, mip); \ 323 } 324 325 326 #if WANT_REDC_2 327 if (REDC_1_TO_REDC_2_THRESHOLD < MUL_TOOM22_THRESHOLD) 328 { 329 if (BELOW_THRESHOLD (n, REDC_1_TO_REDC_2_THRESHOLD)) 330 { 331 if (REDC_1_TO_REDC_2_THRESHOLD < SQR_BASECASE_THRESHOLD 332 || BELOW_THRESHOLD (n, SQR_BASECASE_THRESHOLD)) 333 { 334 #undef MPN_MUL_N 335 #undef MPN_SQR 336 #undef MPN_REDUCE 337 #define MPN_MUL_N(r,a,b,n) mpn_mul_basecase (r,a,n,b,n) 338 #define MPN_SQR(r,a,n) mpn_mul_basecase (r,a,n,a,n) 339 #define MPN_REDUCE(rp,tp,mp,n,mip) MPN_REDC_1 (rp, tp, mp, n, mip[0]) 340 INNERLOOP; 341 } 342 else 343 { 344 #undef MPN_MUL_N 345 #undef MPN_SQR 346 #undef MPN_REDUCE 347 #define MPN_MUL_N(r,a,b,n) mpn_mul_basecase (r,a,n,b,n) 348 #define MPN_SQR(r,a,n) mpn_sqr_basecase (r,a,n) 349 #define MPN_REDUCE(rp,tp,mp,n,mip) MPN_REDC_1 (rp, tp, mp, n, mip[0]) 350 INNERLOOP; 351 } 352 } 353 else if (BELOW_THRESHOLD (n, MUL_TOOM22_THRESHOLD)) 354 { 355 if (MUL_TOOM22_THRESHOLD < SQR_BASECASE_THRESHOLD 356 || BELOW_THRESHOLD (n, SQR_BASECASE_THRESHOLD)) 357 { 358 #undef MPN_MUL_N 359 #undef MPN_SQR 360 #undef MPN_REDUCE 361 #define MPN_MUL_N(r,a,b,n) mpn_mul_basecase (r,a,n,b,n) 362 #define MPN_SQR(r,a,n) mpn_mul_basecase (r,a,n,a,n) 363 #define MPN_REDUCE(rp,tp,mp,n,mip) MPN_REDC_2 (rp, tp, mp, n, mip) 364 INNERLOOP; 365 } 366 else 367 { 368 #undef MPN_MUL_N 369 #undef MPN_SQR 370 #undef MPN_REDUCE 371 #define MPN_MUL_N(r,a,b,n) mpn_mul_basecase (r,a,n,b,n) 372 #define MPN_SQR(r,a,n) mpn_sqr_basecase (r,a,n) 373 #define MPN_REDUCE(rp,tp,mp,n,mip) MPN_REDC_2 (rp, tp, mp, n, mip) 374 INNERLOOP; 375 } 376 } 377 else if (BELOW_THRESHOLD (n, REDC_2_TO_REDC_N_THRESHOLD)) 378 { 379 #undef MPN_MUL_N 380 #undef MPN_SQR 381 #undef MPN_REDUCE 382 #define MPN_MUL_N(r,a,b,n) mpn_mul_n (r,a,b,n) 383 #define MPN_SQR(r,a,n) mpn_sqr (r,a,n) 384 #define MPN_REDUCE(rp,tp,mp,n,mip) MPN_REDC_2 (rp, tp, mp, n, mip) 385 INNERLOOP; 386 } 387 else 388 { 389 #undef MPN_MUL_N 390 #undef MPN_SQR 391 #undef MPN_REDUCE 392 #define MPN_MUL_N(r,a,b,n) mpn_mul_n (r,a,b,n) 393 #define MPN_SQR(r,a,n) mpn_sqr (r,a,n) 394 #define MPN_REDUCE(rp,tp,mp,n,mip) mpn_redc_n (rp, tp, mp, n, mip) 395 INNERLOOP; 396 } 397 } 398 else 399 { 400 if (BELOW_THRESHOLD (n, MUL_TOOM22_THRESHOLD)) 401 { 402 if (MUL_TOOM22_THRESHOLD < SQR_BASECASE_THRESHOLD 403 || BELOW_THRESHOLD (n, SQR_BASECASE_THRESHOLD)) 404 { 405 #undef MPN_MUL_N 406 #undef MPN_SQR 407 #undef MPN_REDUCE 408 #define MPN_MUL_N(r,a,b,n) mpn_mul_basecase (r,a,n,b,n) 409 #define MPN_SQR(r,a,n) mpn_mul_basecase (r,a,n,a,n) 410 #define MPN_REDUCE(rp,tp,mp,n,mip) MPN_REDC_1 (rp, tp, mp, n, mip[0]) 411 INNERLOOP; 412 } 413 else 414 { 415 #undef MPN_MUL_N 416 #undef MPN_SQR 417 #undef MPN_REDUCE 418 #define MPN_MUL_N(r,a,b,n) mpn_mul_basecase (r,a,n,b,n) 419 #define MPN_SQR(r,a,n) mpn_sqr_basecase (r,a,n) 420 #define MPN_REDUCE(rp,tp,mp,n,mip) MPN_REDC_1 (rp, tp, mp, n, mip[0]) 421 INNERLOOP; 422 } 423 } 424 else if (BELOW_THRESHOLD (n, REDC_1_TO_REDC_2_THRESHOLD)) 425 { 426 #undef MPN_MUL_N 427 #undef MPN_SQR 428 #undef MPN_REDUCE 429 #define MPN_MUL_N(r,a,b,n) mpn_mul_n (r,a,b,n) 430 #define MPN_SQR(r,a,n) mpn_sqr (r,a,n) 431 #define MPN_REDUCE(rp,tp,mp,n,mip) MPN_REDC_1 (rp, tp, mp, n, mip[0]) 432 INNERLOOP; 433 } 434 else if (BELOW_THRESHOLD (n, REDC_2_TO_REDC_N_THRESHOLD)) 435 { 436 #undef MPN_MUL_N 437 #undef MPN_SQR 438 #undef MPN_REDUCE 439 #define MPN_MUL_N(r,a,b,n) mpn_mul_n (r,a,b,n) 440 #define MPN_SQR(r,a,n) mpn_sqr (r,a,n) 441 #define MPN_REDUCE(rp,tp,mp,n,mip) MPN_REDC_2 (rp, tp, mp, n, mip) 442 INNERLOOP; 443 } 444 else 445 { 446 #undef MPN_MUL_N 447 #undef MPN_SQR 448 #undef MPN_REDUCE 449 #define MPN_MUL_N(r,a,b,n) mpn_mul_n (r,a,b,n) 450 #define MPN_SQR(r,a,n) mpn_sqr (r,a,n) 451 #define MPN_REDUCE(rp,tp,mp,n,mip) mpn_redc_n (rp, tp, mp, n, mip) 452 INNERLOOP; 453 } 454 } 455 456 #else /* WANT_REDC_2 */ 457 458 if (REDC_1_TO_REDC_N_THRESHOLD < MUL_TOOM22_THRESHOLD) 459 { 460 if (BELOW_THRESHOLD (n, REDC_1_TO_REDC_N_THRESHOLD)) 461 { 462 if (REDC_1_TO_REDC_N_THRESHOLD < SQR_BASECASE_THRESHOLD 463 || BELOW_THRESHOLD (n, SQR_BASECASE_THRESHOLD)) 464 { 465 #undef MPN_MUL_N 466 #undef MPN_SQR 467 #undef MPN_REDUCE 468 #define MPN_MUL_N(r,a,b,n) mpn_mul_basecase (r,a,n,b,n) 469 #define MPN_SQR(r,a,n) mpn_mul_basecase (r,a,n,a,n) 470 #define MPN_REDUCE(rp,tp,mp,n,mip) MPN_REDC_1 (rp, tp, mp, n, mip[0]) 471 INNERLOOP; 472 } 473 else 474 { 475 #undef MPN_MUL_N 476 #undef MPN_SQR 477 #undef MPN_REDUCE 478 #define MPN_MUL_N(r,a,b,n) mpn_mul_basecase (r,a,n,b,n) 479 #define MPN_SQR(r,a,n) mpn_sqr_basecase (r,a,n) 480 #define MPN_REDUCE(rp,tp,mp,n,mip) MPN_REDC_1 (rp, tp, mp, n, mip[0]) 481 INNERLOOP; 482 } 483 } 484 else if (BELOW_THRESHOLD (n, MUL_TOOM22_THRESHOLD)) 485 { 486 if (MUL_TOOM22_THRESHOLD < SQR_BASECASE_THRESHOLD 487 || BELOW_THRESHOLD (n, SQR_BASECASE_THRESHOLD)) 488 { 489 #undef MPN_MUL_N 490 #undef MPN_SQR 491 #undef MPN_REDUCE 492 #define MPN_MUL_N(r,a,b,n) mpn_mul_basecase (r,a,n,b,n) 493 #define MPN_SQR(r,a,n) mpn_mul_basecase (r,a,n,a,n) 494 #define MPN_REDUCE(rp,tp,mp,n,mip) mpn_redc_n (rp, tp, mp, n, mip) 495 INNERLOOP; 496 } 497 else 498 { 499 #undef MPN_MUL_N 500 #undef MPN_SQR 501 #undef MPN_REDUCE 502 #define MPN_MUL_N(r,a,b,n) mpn_mul_basecase (r,a,n,b,n) 503 #define MPN_SQR(r,a,n) mpn_sqr_basecase (r,a,n) 504 #define MPN_REDUCE(rp,tp,mp,n,mip) mpn_redc_n (rp, tp, mp, n, mip) 505 INNERLOOP; 506 } 507 } 508 else 509 { 510 #undef MPN_MUL_N 511 #undef MPN_SQR 512 #undef MPN_REDUCE 513 #define MPN_MUL_N(r,a,b,n) mpn_mul_n (r,a,b,n) 514 #define MPN_SQR(r,a,n) mpn_sqr (r,a,n) 515 #define MPN_REDUCE(rp,tp,mp,n,mip) mpn_redc_n (rp, tp, mp, n, mip) 516 INNERLOOP; 517 } 518 } 519 else 520 { 521 if (BELOW_THRESHOLD (n, MUL_TOOM22_THRESHOLD)) 522 { 523 if (MUL_TOOM22_THRESHOLD < SQR_BASECASE_THRESHOLD 524 || BELOW_THRESHOLD (n, SQR_BASECASE_THRESHOLD)) 525 { 526 #undef MPN_MUL_N 527 #undef MPN_SQR 528 #undef MPN_REDUCE 529 #define MPN_MUL_N(r,a,b,n) mpn_mul_basecase (r,a,n,b,n) 530 #define MPN_SQR(r,a,n) mpn_mul_basecase (r,a,n,a,n) 531 #define MPN_REDUCE(rp,tp,mp,n,mip) MPN_REDC_1 (rp, tp, mp, n, mip[0]) 532 INNERLOOP; 533 } 534 else 535 { 536 #undef MPN_MUL_N 537 #undef MPN_SQR 538 #undef MPN_REDUCE 539 #define MPN_MUL_N(r,a,b,n) mpn_mul_basecase (r,a,n,b,n) 540 #define MPN_SQR(r,a,n) mpn_sqr_basecase (r,a,n) 541 #define MPN_REDUCE(rp,tp,mp,n,mip) MPN_REDC_1 (rp, tp, mp, n, mip[0]) 542 INNERLOOP; 543 } 544 } 545 else if (BELOW_THRESHOLD (n, REDC_1_TO_REDC_N_THRESHOLD)) 546 { 547 #undef MPN_MUL_N 548 #undef MPN_SQR 549 #undef MPN_REDUCE 550 #define MPN_MUL_N(r,a,b,n) mpn_mul_n (r,a,b,n) 551 #define MPN_SQR(r,a,n) mpn_sqr (r,a,n) 552 #define MPN_REDUCE(rp,tp,mp,n,mip) MPN_REDC_1 (rp, tp, mp, n, mip[0]) 553 INNERLOOP; 554 } 555 else 556 { 557 #undef MPN_MUL_N 558 #undef MPN_SQR 559 #undef MPN_REDUCE 560 #define MPN_MUL_N(r,a,b,n) mpn_mul_n (r,a,b,n) 561 #define MPN_SQR(r,a,n) mpn_sqr (r,a,n) 562 #define MPN_REDUCE(rp,tp,mp,n,mip) mpn_redc_n (rp, tp, mp, n, mip) 563 INNERLOOP; 564 } 565 } 566 #endif /* WANT_REDC_2 */ 567 568 done: 569 570 MPN_COPY (tp, rp, n); 571 MPN_ZERO (tp + n, n); 572 573 #if WANT_REDC_2 574 if (BELOW_THRESHOLD (n, REDC_1_TO_REDC_2_THRESHOLD)) 575 MPN_REDC_1 (rp, tp, mp, n, mip[0]); 576 else if (BELOW_THRESHOLD (n, REDC_2_TO_REDC_N_THRESHOLD)) 577 MPN_REDC_2 (rp, tp, mp, n, mip); 578 #else 579 if (BELOW_THRESHOLD (n, REDC_1_TO_REDC_N_THRESHOLD)) 580 MPN_REDC_1 (rp, tp, mp, n, mip[0]); 581 #endif 582 else 583 mpn_redc_n (rp, tp, mp, n, mip); 584 585 if (mpn_cmp (rp, mp, n) >= 0) 586 mpn_sub_n (rp, rp, mp, n); 587 588 TMP_FREE; 589 } 590