Return to aorsmul.c CVS log

Up to [local] / OpenXM_contrib / gmp / mpz

Annotation of OpenXM_contrib/gmp/mpz/aorsmul.c, Revision 1.1

1.1     ! ohara       1: /* mpz_addmul, mpz_submul -- add or subtract multiple.
        !             2:
        !             3: Copyright 2001 Free Software Foundation, Inc.
        !             4:
        !             5: This file is part of the GNU MP Library.
        !             6:
        !             7: The GNU MP Library is free software; you can redistribute it and/or modify
        !             8: it under the terms of the GNU Lesser General Public License as published by
        !             9: the Free Software Foundation; either version 2.1 of the License, or (at your
        !            10: option) any later version.
        !            11:
        !            12: The GNU MP Library is distributed in the hope that it will be useful, but
        !            13: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
        !            14: or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
        !            15: License for more details.
        !            16:
        !            17: You should have received a copy of the GNU Lesser General Public License
        !            18: along with the GNU MP Library; see the file COPYING.LIB.  If not, write to
        !            19: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
        !            20: MA 02111-1307, USA. */
        !            21:
        !            22: #include "gmp.h"
        !            23: #include "gmp-impl.h"
        !            24:
        !            25:
        !            26: /* expecting x and y both with non-zero high limbs */
        !            27: #define mpn_cmp_twosizes_lt(xp,xsize, yp,ysize)                 \
        !            28:   ((xsize) < (ysize)                                            \
        !            29:    || ((xsize) == (ysize) && mpn_cmp (xp, yp, xsize) < 0))
        !            30:
        !            31:
        !            32: /* sub>=0 means an addmul w += x*y, sub<0 means a submul w -= x*y.
        !            33:
        !            34:    The signs of w, x and y are fully accounted for by each flipping "sub".
        !            35:
        !            36:    The sign of w is retained for the result, unless the absolute value
        !            37:    submul underflows, in which case it flips.  */
        !            38:
        !            39: static void __gmpz_aorsmul _PROTO ((REGPARM_3_1 (mpz_ptr w, mpz_srcptr x, mpz_srcptr y, mp_size_t sub))) REGPARM_ATTR (1);
        !            40: #define mpz_aorsmul(w,x,y,sub)  __gmpz_aorsmul (REGPARM_3_1 (w, x, y, sub))
        !            41:
        !            42: static void
        !            43: mpz_aorsmul (mpz_ptr w, mpz_srcptr x, mpz_srcptr y, mp_size_t sub)
        !            44: {
        !            45:   mp_size_t  xsize, ysize, tsize, wsize, wsize_signed;
        !            46:   mp_ptr     wp, tp;
        !            47:   mp_limb_t  c, high;
        !            48:   TMP_DECL (marker);
        !            49:
        !            50:   /* w unaffected if x==0 or y==0 */
        !            51:   xsize = SIZ(x);
        !            52:   ysize = SIZ(y);
        !            53:   if (xsize == 0 || ysize == 0)
        !            54:     return;
        !            55:
        !            56:   /* make x the bigger of the two */
        !            57:   if (ABS(ysize) > ABS(xsize))
        !            58:     {
        !            59:       MPZ_SRCPTR_SWAP (x, y);
        !            60:       MP_SIZE_T_SWAP (xsize, ysize);
        !            61:     }
        !            62:
        !            63:   sub ^= ysize;
        !            64:   ysize = ABS(ysize);
        !            65:
        !            66:   /* use mpn_addmul_1/mpn_submul_1 if possible */
        !            67:   if (ysize == 1)
        !            68:     {
        !            69:       mpz_aorsmul_1 (w, x, PTR(y)[0], sub);
        !            70:       return;
        !            71:     }
        !            72:
        !            73:   sub ^= xsize;
        !            74:   xsize = ABS(xsize);
        !            75:
        !            76:   wsize_signed = SIZ(w);
        !            77:   sub ^= wsize_signed;
        !            78:   wsize = ABS(wsize_signed);
        !            79:
        !            80:   tsize = xsize + ysize;
        !            81:   MPZ_REALLOC (w, MAX (wsize, tsize) + 1);
        !            82:   wp = PTR(w);
        !            83:
        !            84:   if (wsize_signed == 0)
        !            85:     {
        !            86:       /* Nothing to add to, just set w=x*y.  No w==x or w==y overlap here,
        !            87:          since we know x,y!=0 but w==0.  */
        !            88:       high = mpn_mul (wp, PTR(x),xsize, PTR(y),ysize);
        !            89:       tsize -= (high == 0);
        !            90:       SIZ(w) = (sub >= 0 ? tsize : -tsize);
        !            91:       return;
        !            92:     }
        !            93:
        !            94:   TMP_MARK (marker);
        !            95:   tp = TMP_ALLOC_LIMBS (tsize);
        !            96:
        !            97:   high = mpn_mul (tp, PTR(x),xsize, PTR(y),ysize);
        !            98:   tsize -= (high == 0);
        !            99:   ASSERT (tp[tsize-1] != 0);
        !           100:   if (sub >= 0)
        !           101:     {
        !           102:       mp_srcptr up    = wp;
        !           103:       mp_size_t usize = wsize;
        !           104:
        !           105:       if (usize < tsize)
        !           106:         {
        !           107:           up    = tp;
        !           108:           usize = tsize;
        !           109:           tp    = wp;
        !           110:           tsize = wsize;
        !           111:
        !           112:           wsize = usize;
        !           113:         }
        !           114:
        !           115:       c = mpn_add (wp, up,usize, tp,tsize);
        !           116:       wp[wsize] = c;
        !           117:       wsize += (c != 0);
        !           118:     }
        !           119:   else
        !           120:     {
        !           121:       mp_srcptr up    = wp;
        !           122:       mp_size_t usize = wsize;
        !           123:
        !           124:       if (mpn_cmp_twosizes_lt (up,usize, tp,tsize))
        !           125:         {
        !           126:           up    = tp;
        !           127:           usize = tsize;
        !           128:           tp    = wp;
        !           129:           tsize = wsize;
        !           130:
        !           131:           wsize = usize;
        !           132:           wsize_signed = -wsize_signed;
        !           133:         }
        !           134:
        !           135:       ASSERT_NOCARRY (mpn_sub (wp, up,usize, tp,tsize));
        !           136:       wsize = usize;
        !           137:       MPN_NORMALIZE (wp, wsize);
        !           138:     }
        !           139:
        !           140:   SIZ(w) = (wsize_signed >= 0 ? wsize : -wsize);
        !           141:
        !           142:   TMP_FREE (marker);
        !           143: }
        !           144:
        !           145:
        !           146: void
        !           147: mpz_addmul (mpz_ptr w, mpz_srcptr u, mpz_srcptr v)
        !           148: {
        !           149:   mpz_aorsmul (w, u, v, (mp_size_t) 0);
        !           150: }
        !           151:
        !           152: void
        !           153: mpz_submul (mpz_ptr w, mpz_srcptr u, mpz_srcptr v)
        !           154: {
        !           155:   mpz_aorsmul (w, u, v, (mp_size_t) -1);
        !           156: }