Annotation of OpenXM_contrib/gmp/mpz/addmul_ui.c, Revision 1.1
1.1 ! maekawa 1: /* mpz_addmul_ui(prodsum, multiplier, small_multiplicand) --
! 2: Add MULTIPLICATOR times SMALL_MULTIPLICAND to PRODSUM.
! 3:
! 4: Copyright (C) 1997, 2000 Free Software Foundation, Inc.
! 5:
! 6: This file is part of the GNU MP Library.
! 7:
! 8: The GNU MP Library is free software; you can redistribute it and/or modify
! 9: it under the terms of the GNU Lesser General Public License as published by
! 10: the Free Software Foundation; either version 2.1 of the License, or (at your
! 11: option) any later version.
! 12:
! 13: The GNU MP Library is distributed in the hope that it will be useful, but
! 14: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
! 15: or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
! 16: License for more details.
! 17:
! 18: You should have received a copy of the GNU Lesser General Public License
! 19: along with the GNU MP Library; see the file COPYING.LIB. If not, write to
! 20: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
! 21: MA 02111-1307, USA. */
! 22:
! 23: #include "gmp.h"
! 24: #include "gmp-impl.h"
! 25:
! 26: static mp_limb_t mpn_neg1 _PROTO ((mp_ptr, mp_size_t));
! 27:
! 28: #if 0
! 29: #undef MPN_NORMALIZE
! 30: #define MPN_NORMALIZE(DST, NLIMBS) \
! 31: do { \
! 32: while (--(NLIMBS) >= 0 && (DST)[NLIMBS] == 0) \
! 33: ; \
! 34: (NLIMBS)++; \
! 35: } while (0)
! 36: #undef MPN_NORMALIZE_NOT_ZERO
! 37: #define MPN_NORMALIZE_NOT_ZERO(DST, NLIMBS) \
! 38: do { \
! 39: while ((DST)[--(NLIMBS)] == 0) \
! 40: ; \
! 41: (NLIMBS)++; \
! 42: } while (0)
! 43: #endif
! 44:
! 45: void
! 46: #if __STDC__
! 47: mpz_addmul_ui (mpz_ptr rz, mpz_srcptr az, unsigned long int bu)
! 48: #else
! 49: mpz_addmul_ui (rz, az, bu)
! 50: mpz_ptr rz;
! 51: mpz_srcptr az;
! 52: unsigned long int bu;
! 53: #endif
! 54: {
! 55: mp_size_t rn, an;
! 56: mp_ptr rp, ap;
! 57:
! 58: an = SIZ (az);
! 59:
! 60: /* If either multiplier is zero, result is unaffected. */
! 61: if (bu == 0 || an == 0)
! 62: return;
! 63:
! 64: rn = SIZ (rz);
! 65:
! 66: if (rn == 0)
! 67: {
! 68: mp_limb_t cy;
! 69:
! 70: an = ABS (an);
! 71: if (ALLOC (rz) <= an)
! 72: _mpz_realloc (rz, an + 1);
! 73: rp = PTR (rz);
! 74: ap = PTR (az);
! 75: cy = mpn_mul_1 (rp, ap, an, (mp_limb_t) bu);
! 76: rp[an] = cy;
! 77: an += cy != 0;
! 78: SIZ (rz) = SIZ (az) >= 0 ? an : -an;
! 79: return;
! 80: }
! 81:
! 82: if ((an ^ rn) >= 0)
! 83: {
! 84: /* Sign of operands are the same--really add. */
! 85: an = ABS (an);
! 86: rn = ABS (rn);
! 87: if (rn > an)
! 88: {
! 89: mp_limb_t cy;
! 90: if (ALLOC (rz) <= rn)
! 91: _mpz_realloc (rz, rn + 1);
! 92: rp = PTR (rz);
! 93: ap = PTR (az);
! 94: cy = mpn_addmul_1 (rp, ap, an, (mp_limb_t) bu);
! 95: cy = mpn_add_1 (rp + an, rp + an, rn - an, cy);
! 96: rp[rn] = cy;
! 97: rn += cy != 0;
! 98: SIZ (rz) = SIZ (rz) >= 0 ? rn : -rn;
! 99: return;
! 100: }
! 101: else
! 102: {
! 103: mp_limb_t cy;
! 104: if (ALLOC (rz) <= an)
! 105: _mpz_realloc (rz, an + 1);
! 106: rp = PTR (rz);
! 107: ap = PTR (az);
! 108: cy = mpn_addmul_1 (rp, ap, rn, (mp_limb_t) bu);
! 109: if (an != rn)
! 110: {
! 111: mp_limb_t cy2;
! 112: cy2 = mpn_mul_1 (rp + rn, ap + rn, an - rn, (mp_limb_t) bu);
! 113: cy = cy2 + mpn_add_1 (rp + rn, rp + rn, an - rn, cy);
! 114: }
! 115: rn = an;
! 116: rp[rn] = cy;
! 117: rn += cy != 0;
! 118: SIZ (rz) = SIZ (rz) >= 0 ? rn : -rn;
! 119: return;
! 120: }
! 121: }
! 122: else
! 123: {
! 124: /* Sign of operands are different--actually subtract. */
! 125: an = ABS (an);
! 126: rn = ABS (rn);
! 127: if (rn > an)
! 128: {
! 129: mp_limb_t cy;
! 130: rp = PTR (rz);
! 131: ap = PTR (az);
! 132: cy = mpn_submul_1 (rp, ap, an, (mp_limb_t) bu);
! 133: cy = mpn_sub_1 (rp + an, rp + an, rn - an, cy);
! 134: if (cy != 0)
! 135: {
! 136: mpn_neg1 (rp, rn);
! 137: MPN_NORMALIZE_NOT_ZERO (rp, rn);
! 138: }
! 139: else
! 140: {
! 141: MPN_NORMALIZE (rp, rn);
! 142: rn = -rn;
! 143: }
! 144:
! 145: SIZ (rz) = SIZ (rz) >= 0 ? -rn : rn;
! 146: return;
! 147: }
! 148: else
! 149: {
! 150: /* Tricky case. We need to subtract an operand that might be larger
! 151: than the minuend. To avoid allocating temporary space, we compute
! 152: a*b-r instead of r-a*b and then negate. */
! 153: mp_limb_t cy;
! 154: if (ALLOC (rz) <= an)
! 155: _mpz_realloc (rz, an + 1);
! 156: rp = PTR (rz);
! 157: ap = PTR (az);
! 158: cy = mpn_submul_1 (rp, ap, rn, (mp_limb_t) bu);
! 159: if (an != rn)
! 160: {
! 161: mp_limb_t cy2;
! 162: cy -= mpn_neg1 (rp, rn);
! 163: cy2 = mpn_mul_1 (rp + rn, ap + rn, an - rn, (mp_limb_t) bu);
! 164: if (cy == ~(mp_limb_t) 0)
! 165: cy = cy2 - mpn_sub_1 (rp + rn, rp + rn, an - rn, (mp_limb_t) 1);
! 166: else
! 167: cy = cy2 + mpn_add_1 (rp + rn, rp + rn, an - rn, cy);
! 168: rp[an] = cy;
! 169: rn = an + (cy != 0);
! 170: rn -= rp[rn - 1] == 0;
! 171: }
! 172: else if (cy != 0)
! 173: {
! 174: cy -= mpn_neg1 (rp, rn);
! 175: rp[an] = cy;
! 176: rn = an + 1;
! 177: MPN_NORMALIZE_NOT_ZERO (rp, rn);
! 178: }
! 179: else
! 180: {
! 181: rn = an;
! 182: MPN_NORMALIZE (rp, rn);
! 183: rn = -rn;
! 184: }
! 185:
! 186: SIZ (rz) = SIZ (rz) >= 0 ? -rn : rn;
! 187: return;
! 188: }
! 189: }
! 190: }
! 191:
! 192: static mp_limb_t
! 193: #if __STDC__
! 194: mpn_neg1 (mp_ptr rp, mp_size_t rn)
! 195: #else
! 196: mpn_neg1 (rp, rn)
! 197: mp_ptr rp;
! 198: mp_size_t rn;
! 199: #endif
! 200: {
! 201: mp_size_t i;
! 202:
! 203: while (rn != 0 && rp[0] == 0)
! 204: rp++, rn--;
! 205:
! 206: if (rn != 0)
! 207: {
! 208: rp[0] = -rp[0];
! 209: for (i = 1; i < rn; i++)
! 210: rp[i] = ~rp[i];
! 211: return 1;
! 212: }
! 213: return 0;
! 214: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to addmul_ui.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / gmp / mpz