diff options
Diffstat (limited to 'gmp/mpz/aorsmul_i.c')
| -rw-r--r-- | gmp/mpz/aorsmul_i.c | 145 |
1 files changed, 68 insertions, 77 deletions
diff --git a/gmp/mpz/aorsmul_i.c b/gmp/mpz/aorsmul_i.c index eafd5f36de..c23495c405 100644 --- a/gmp/mpz/aorsmul_i.c +++ b/gmp/mpz/aorsmul_i.c @@ -4,33 +4,22 @@ ALMOST CERTAIN TO BE SUBJECT TO INCOMPATIBLE CHANGES OR DISAPPEAR COMPLETELY IN FUTURE GNU MP RELEASES. -Copyright 2001, 2002, 2004, 2005, 2012 Free Software Foundation, Inc. +Copyright 2001, 2002, 2004, 2005 Free Software Foundation, Inc. This file is part of the GNU MP Library. The GNU MP Library is free software; you can redistribute it and/or modify -it under the terms of either: - - * the GNU Lesser General Public License as published by the Free - Software Foundation; either version 3 of the License, or (at your - option) any later version. - -or - - * the GNU General Public License as published by the Free Software - Foundation; either version 2 of the License, or (at your option) any - later version. - -or both in parallel, as here. +it under the terms of the GNU Lesser General Public License as published by +the Free Software Foundation; either version 3 of the License, or (at your +option) any later version. The GNU MP Library is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY -or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License -for more details. +or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public +License for more details. -You should have received copies of the GNU General Public License and the -GNU Lesser General Public License along with the GNU MP Library. If not, -see https://www.gnu.org/licenses/. */ +You should have received a copy of the GNU Lesser General Public License +along with the GNU MP Library. If not, see http://www.gnu.org/licenses/. */ #include "gmp.h" #include "gmp-impl.h" @@ -58,13 +47,13 @@ see https://www.gnu.org/licenses/. */ The final w will retain its sign, unless an underflow occurs in a submul of absolute values, in which case it's flipped. - If x has more limbs than w, then mpn_submul_1 followed by mpn_com is + If x has more limbs than w, then mpn_submul_1 followed by mpn_com_n is used. The alternative would be mpn_mul_1 into temporary space followed by mpn_sub_n. Avoiding temporary space seem good, and submul+com stands a chance of being faster since it involves only one set of carry propagations, not two. Note that doing an addmul_1 with a twos-complement negative y doesn't work, because it effectively adds an - extra x * 2^GMP_LIMB_BITS. */ + extra x * 2^BITS_PER_MP_LIMB. */ REGPARM_ATTR(1) void mpz_aorsmul_1 (mpz_ptr w, mpz_srcptr x, mp_limb_t y, mp_size_t sub) @@ -86,7 +75,8 @@ mpz_aorsmul_1 (mpz_ptr w, mpz_srcptr x, mp_limb_t y, mp_size_t sub) if (wsize_signed == 0) { /* nothing to add to, just set x*y, "sub" gives the sign */ - wp = MPZ_REALLOC (w, xsize+1); + MPZ_REALLOC (w, xsize+1); + wp = PTR (w); cy = mpn_mul_1 (wp, PTR(x), xsize, y); wp[xsize] = cy; xsize += (cy != 0); @@ -98,7 +88,8 @@ mpz_aorsmul_1 (mpz_ptr w, mpz_srcptr x, mp_limb_t y, mp_size_t sub) wsize = ABS (wsize_signed); new_wsize = MAX (wsize, xsize); - wp = MPZ_REALLOC (w, new_wsize+1); + MPZ_REALLOC (w, new_wsize+1); + wp = PTR (w); xp = PTR (x); min_size = MIN (wsize, xsize); @@ -113,25 +104,25 @@ mpz_aorsmul_1 (mpz_ptr w, mpz_srcptr x, mp_limb_t y, mp_size_t sub) dsize = xsize - wsize; #if HAVE_NATIVE_mpn_mul_1c if (dsize > 0) - cy = mpn_mul_1c (wp, xp, dsize, y, cy); + cy = mpn_mul_1c (wp, xp, dsize, y, cy); else if (dsize < 0) - { - dsize = -dsize; - cy = mpn_add_1 (wp, wp, dsize, cy); - } + { + dsize = -dsize; + cy = mpn_add_1 (wp, wp, dsize, cy); + } #else if (dsize != 0) - { - mp_limb_t cy2; - if (dsize > 0) - cy2 = mpn_mul_1 (wp, xp, dsize, y); - else - { - dsize = -dsize; - cy2 = 0; - } - cy = cy2 + mpn_add_1 (wp, wp, dsize, cy); - } + { + mp_limb_t cy2; + if (dsize > 0) + cy2 = mpn_mul_1 (wp, xp, dsize, y); + else + { + dsize = -dsize; + cy2 = 0; + } + cy = cy2 + mpn_add_1 (wp, wp, dsize, cy); + } #endif wp[dsize] = cy; @@ -143,52 +134,52 @@ mpz_aorsmul_1 (mpz_ptr w, mpz_srcptr x, mp_limb_t y, mp_size_t sub) cy = mpn_submul_1 (wp, xp, min_size, y); if (wsize >= xsize) - { - /* if w bigger than x, then propagate borrow through it */ - if (wsize != xsize) - cy = mpn_sub_1 (wp+xsize, wp+xsize, wsize-xsize, cy); - - if (cy != 0) - { - /* Borrow out of w, take twos complement negative to get - absolute value, flip sign of w. */ - wp[new_wsize] = ~-cy; /* extra limb is 0-cy */ - mpn_com (wp, wp, new_wsize); - new_wsize++; - MPN_INCR_U (wp, new_wsize, CNST_LIMB(1)); - wsize_signed = -wsize_signed; - } - } + { + /* if w bigger than x, then propagate borrow through it */ + if (wsize != xsize) + cy = mpn_sub_1 (wp+xsize, wp+xsize, wsize-xsize, cy); + + if (cy != 0) + { + /* Borrow out of w, take twos complement negative to get + absolute value, flip sign of w. */ + wp[new_wsize] = ~-cy; /* extra limb is 0-cy */ + mpn_com_n (wp, wp, new_wsize); + new_wsize++; + MPN_INCR_U (wp, new_wsize, CNST_LIMB(1)); + wsize_signed = -wsize_signed; + } + } else /* wsize < xsize */ - { - /* x bigger than w, so want x*y-w. Submul has given w-x*y, so - take twos complement and use an mpn_mul_1 for the rest. */ + { + /* x bigger than w, so want x*y-w. Submul has given w-x*y, so + take twos complement and use an mpn_mul_1 for the rest. */ - mp_limb_t cy2; + mp_limb_t cy2; - /* -(-cy*b^n + w-x*y) = (cy-1)*b^n + ~(w-x*y) + 1 */ - mpn_com (wp, wp, wsize); - cy += mpn_add_1 (wp, wp, wsize, CNST_LIMB(1)); - cy -= 1; + /* -(-cy*b^n + w-x*y) = (cy-1)*b^n + ~(w-x*y) + 1 */ + mpn_com_n (wp, wp, wsize); + cy += mpn_add_1 (wp, wp, wsize, CNST_LIMB(1)); + cy -= 1; - /* If cy-1 == -1 then hold that -1 for latter. mpn_submul_1 never - returns cy==MP_LIMB_T_MAX so that value always indicates a -1. */ - cy2 = (cy == MP_LIMB_T_MAX); - cy += cy2; - MPN_MUL_1C (cy, wp+wsize, xp+wsize, xsize-wsize, y, cy); - wp[new_wsize] = cy; - new_wsize += (cy != 0); + /* If cy-1 == -1 then hold that -1 for latter. mpn_submul_1 never + returns cy==MP_LIMB_T_MAX so that value always indicates a -1. */ + cy2 = (cy == MP_LIMB_T_MAX); + cy += cy2; + MPN_MUL_1C (cy, wp+wsize, xp+wsize, xsize-wsize, y, cy); + wp[new_wsize] = cy; + new_wsize += (cy != 0); - /* Apply any -1 from above. The value at wp+wsize is non-zero - because y!=0 and the high limb of x will be non-zero. */ - if (cy2) - MPN_DECR_U (wp+wsize, new_wsize-wsize, CNST_LIMB(1)); + /* Apply any -1 from above. The value at wp+wsize is non-zero + because y!=0 and the high limb of x will be non-zero. */ + if (cy2) + MPN_DECR_U (wp+wsize, new_wsize-wsize, CNST_LIMB(1)); - wsize_signed = -wsize_signed; - } + wsize_signed = -wsize_signed; + } /* submul can produce high zero limbs due to cancellation, both when w - has more limbs or x has more */ + has more limbs or x has more */ MPN_NORMALIZE (wp, new_wsize); } |