Annotation of OpenXM_contrib/gmp/mpn/generic/pow_1.c, Revision 1.1
1.1 ! ohara 1: /* mpn_pow_1 -- Compute powers R = U^exp.
! 2:
! 3: THE FUNCTIONS IN THIS FILE ARE FOR INTERNAL USE ONLY. THEY'RE ALMOST
! 4: CERTAIN TO BE SUBJECT TO INCOMPATIBLE CHANGES OR DISAPPEAR COMPLETELY IN
! 5: FUTURE GNU MP RELEASES.
! 6:
! 7: Copyright 2002 Free Software Foundation, Inc.
! 8:
! 9: This file is part of the GNU MP Library.
! 10:
! 11: The GNU MP Library is free software; you can redistribute it and/or modify it
! 12: under the terms of the GNU Lesser General Public License as published by the
! 13: Free Software Foundation; either version 2.1 of the License, or (at your
! 14: option) any later version.
! 15:
! 16: The GNU MP Library is distributed in the hope that it will be useful, but
! 17: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
! 18: FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License
! 19: for more details.
! 20:
! 21: You should have received a copy of the GNU Lesser General Public License along
! 22: with the GNU MP Library; see the file COPYING.LIB. If not, write to the Free
! 23: Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
! 24: USA. */
! 25:
! 26:
! 27: #include "gmp.h"
! 28: #include "gmp-impl.h"
! 29: #include "longlong.h"
! 30:
! 31: mp_size_t
! 32: mpn_pow_1 (mp_ptr rp, mp_srcptr bp, mp_size_t bn, mp_limb_t exp, mp_ptr tp)
! 33: {
! 34: mp_limb_t x;
! 35: int cnt, i;
! 36: mp_size_t rn;
! 37: int par;
! 38:
! 39: if (exp <= 1)
! 40: {
! 41: if (exp == 0)
! 42: {
! 43: rp[0] = 1;
! 44: return 1;
! 45: }
! 46: else
! 47: {
! 48: MPN_COPY (rp, bp, bn);
! 49: return bn;
! 50: }
! 51: }
! 52:
! 53: /* Count number of bits in exp, and compute where to put initial square in
! 54: order to magically get results in the entry rp. Use simple code,
! 55: optimized for small exp. For large exp, the bignum operations will take
! 56: so much time that the slowness of this code will be negligible. */
! 57: par = 0;
! 58: cnt = GMP_LIMB_BITS;
! 59: for (x = exp; x != 0; x >>= 1)
! 60: {
! 61: par ^= x & 1;
! 62: cnt--;
! 63: }
! 64: exp <<= cnt;
! 65:
! 66: if (bn == 1)
! 67: {
! 68: mp_limb_t bl = bp[0];
! 69:
! 70: if ((cnt & 1) != 0)
! 71: MP_PTR_SWAP (rp, tp);
! 72:
! 73: mpn_sqr_n (rp, bp, bn);
! 74: rn = 2 * bn; rn -= rp[rn - 1] == 0;
! 75:
! 76: for (i = GMP_LIMB_BITS - cnt - 1;;)
! 77: {
! 78: exp <<= 1;
! 79: if ((exp & GMP_LIMB_HIGHBIT) != 0)
! 80: {
! 81: rp[rn] = mpn_mul_1 (rp, rp, rn, bl);
! 82: rn += rp[rn] != 0;
! 83: }
! 84:
! 85: if (--i == 0)
! 86: break;
! 87:
! 88: mpn_sqr_n (tp, rp, rn);
! 89: rn = 2 * rn; rn -= tp[rn - 1] == 0;
! 90: MP_PTR_SWAP (rp, tp);
! 91: }
! 92: }
! 93: else
! 94: {
! 95: if (((par ^ cnt) & 1) == 0)
! 96: MP_PTR_SWAP (rp, tp);
! 97:
! 98: mpn_sqr_n (rp, bp, bn);
! 99: rn = 2 * bn; rn -= rp[rn - 1] == 0;
! 100:
! 101: for (i = GMP_LIMB_BITS - cnt - 1;;)
! 102: {
! 103: exp <<= 1;
! 104: if ((exp & GMP_LIMB_HIGHBIT) != 0)
! 105: {
! 106: rn = rn + bn - (mpn_mul (tp, rp, rn, bp, bn) == 0);
! 107: MP_PTR_SWAP (rp, tp);
! 108: }
! 109:
! 110: if (--i == 0)
! 111: break;
! 112:
! 113: mpn_sqr_n (tp, rp, rn);
! 114: rn = 2 * rn; rn -= tp[rn - 1] == 0;
! 115: MP_PTR_SWAP (rp, tp);
! 116: }
! 117: }
! 118:
! 119: return rn;
! 120: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>