[BACK]Return to pow_1.c CVS log [TXT][DIR] Up to [local] / OpenXM_contrib / gmp / mpn / generic

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>