Home | History | Annotate | Line # | Download | only in mpz
mul.c revision 1.1.1.2
      1      1.1  mrg /* mpz_mul -- Multiply two integers.
      2      1.1  mrg 
      3  1.1.1.2  mrg Copyright 1991, 1993, 1994, 1996, 2000, 2001, 2005, 2009, 2011, 2012 Free
      4      1.1  mrg Software Foundation, Inc.
      5      1.1  mrg 
      6      1.1  mrg This file is part of the GNU MP Library.
      7      1.1  mrg 
      8      1.1  mrg The GNU MP Library is free software; you can redistribute it and/or modify
      9      1.1  mrg it under the terms of the GNU Lesser General Public License as published by
     10      1.1  mrg the Free Software Foundation; either version 3 of the License, or (at your
     11      1.1  mrg option) any later version.
     12      1.1  mrg 
     13      1.1  mrg The GNU MP Library is distributed in the hope that it will be useful, but
     14      1.1  mrg WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
     15      1.1  mrg or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
     16      1.1  mrg License for more details.
     17      1.1  mrg 
     18      1.1  mrg You should have received a copy of the GNU Lesser General Public License
     19      1.1  mrg along with the GNU MP Library.  If not, see http://www.gnu.org/licenses/.  */
     20      1.1  mrg 
     21      1.1  mrg #include <stdio.h> /* for NULL */
     22      1.1  mrg #include "gmp.h"
     23      1.1  mrg #include "gmp-impl.h"
     24      1.1  mrg 
     25      1.1  mrg 
     26      1.1  mrg void
     27      1.1  mrg mpz_mul (mpz_ptr w, mpz_srcptr u, mpz_srcptr v)
     28      1.1  mrg {
     29      1.1  mrg   mp_size_t usize;
     30      1.1  mrg   mp_size_t vsize;
     31      1.1  mrg   mp_size_t wsize;
     32      1.1  mrg   mp_size_t sign_product;
     33      1.1  mrg   mp_ptr up, vp;
     34      1.1  mrg   mp_ptr wp;
     35      1.1  mrg   mp_ptr free_me;
     36      1.1  mrg   size_t free_me_size;
     37      1.1  mrg   mp_limb_t cy_limb;
     38      1.1  mrg   TMP_DECL;
     39      1.1  mrg 
     40      1.1  mrg   usize = SIZ (u);
     41      1.1  mrg   vsize = SIZ (v);
     42      1.1  mrg   sign_product = usize ^ vsize;
     43      1.1  mrg   usize = ABS (usize);
     44      1.1  mrg   vsize = ABS (vsize);
     45      1.1  mrg 
     46      1.1  mrg   if (usize < vsize)
     47      1.1  mrg     {
     48      1.1  mrg       MPZ_SRCPTR_SWAP (u, v);
     49      1.1  mrg       MP_SIZE_T_SWAP (usize, vsize);
     50      1.1  mrg     }
     51      1.1  mrg 
     52      1.1  mrg   if (vsize == 0)
     53      1.1  mrg     {
     54  1.1.1.2  mrg       SIZ (w) = 0;
     55      1.1  mrg       return;
     56      1.1  mrg     }
     57      1.1  mrg 
     58      1.1  mrg #if HAVE_NATIVE_mpn_mul_2
     59      1.1  mrg   if (vsize <= 2)
     60      1.1  mrg     {
     61  1.1.1.2  mrg       wp = MPZ_REALLOC (w, usize+vsize);
     62      1.1  mrg       if (vsize == 1)
     63  1.1.1.2  mrg 	cy_limb = mpn_mul_1 (wp, PTR (u), usize, PTR (v)[0]);
     64      1.1  mrg       else
     65  1.1.1.2  mrg 	{
     66  1.1.1.2  mrg 	  cy_limb = mpn_mul_2 (wp, PTR (u), usize, PTR (v));
     67  1.1.1.2  mrg 	  usize++;
     68  1.1.1.2  mrg 	}
     69      1.1  mrg       wp[usize] = cy_limb;
     70      1.1  mrg       usize += (cy_limb != 0);
     71  1.1.1.2  mrg       SIZ (w) = (sign_product >= 0 ? usize : -usize);
     72      1.1  mrg       return;
     73      1.1  mrg     }
     74      1.1  mrg #else
     75      1.1  mrg   if (vsize == 1)
     76      1.1  mrg     {
     77  1.1.1.2  mrg       wp = MPZ_REALLOC (w, usize+1);
     78  1.1.1.2  mrg       cy_limb = mpn_mul_1 (wp, PTR (u), usize, PTR (v)[0]);
     79      1.1  mrg       wp[usize] = cy_limb;
     80      1.1  mrg       usize += (cy_limb != 0);
     81  1.1.1.2  mrg       SIZ (w) = (sign_product >= 0 ? usize : -usize);
     82      1.1  mrg       return;
     83      1.1  mrg     }
     84      1.1  mrg #endif
     85      1.1  mrg 
     86      1.1  mrg   TMP_MARK;
     87      1.1  mrg   free_me = NULL;
     88  1.1.1.2  mrg   up = PTR (u);
     89  1.1.1.2  mrg   vp = PTR (v);
     90  1.1.1.2  mrg   wp = PTR (w);
     91      1.1  mrg 
     92      1.1  mrg   /* Ensure W has space enough to store the result.  */
     93      1.1  mrg   wsize = usize + vsize;
     94  1.1.1.2  mrg   if (ALLOC (w) < wsize)
     95      1.1  mrg     {
     96      1.1  mrg       if (wp == up || wp == vp)
     97      1.1  mrg 	{
     98      1.1  mrg 	  free_me = wp;
     99  1.1.1.2  mrg 	  free_me_size = ALLOC (w);
    100      1.1  mrg 	}
    101      1.1  mrg       else
    102  1.1.1.2  mrg 	(*__gmp_free_func) (wp, ALLOC (w) * BYTES_PER_MP_LIMB);
    103      1.1  mrg 
    104  1.1.1.2  mrg       ALLOC (w) = wsize;
    105      1.1  mrg       wp = (mp_ptr) (*__gmp_allocate_func) (wsize * BYTES_PER_MP_LIMB);
    106  1.1.1.2  mrg       PTR (w) = wp;
    107      1.1  mrg     }
    108      1.1  mrg   else
    109      1.1  mrg     {
    110      1.1  mrg       /* Make U and V not overlap with W.  */
    111      1.1  mrg       if (wp == up)
    112      1.1  mrg 	{
    113      1.1  mrg 	  /* W and U are identical.  Allocate temporary space for U.  */
    114      1.1  mrg 	  up = TMP_ALLOC_LIMBS (usize);
    115      1.1  mrg 	  /* Is V identical too?  Keep it identical with U.  */
    116      1.1  mrg 	  if (wp == vp)
    117      1.1  mrg 	    vp = up;
    118      1.1  mrg 	  /* Copy to the temporary space.  */
    119      1.1  mrg 	  MPN_COPY (up, wp, usize);
    120      1.1  mrg 	}
    121      1.1  mrg       else if (wp == vp)
    122      1.1  mrg 	{
    123      1.1  mrg 	  /* W and V are identical.  Allocate temporary space for V.  */
    124      1.1  mrg 	  vp = TMP_ALLOC_LIMBS (vsize);
    125      1.1  mrg 	  /* Copy to the temporary space.  */
    126      1.1  mrg 	  MPN_COPY (vp, wp, vsize);
    127      1.1  mrg 	}
    128      1.1  mrg     }
    129      1.1  mrg 
    130      1.1  mrg   if (up == vp)
    131      1.1  mrg     {
    132      1.1  mrg       mpn_sqr (wp, up, usize);
    133      1.1  mrg       cy_limb = wp[wsize - 1];
    134      1.1  mrg     }
    135      1.1  mrg   else
    136      1.1  mrg     {
    137      1.1  mrg       cy_limb = mpn_mul (wp, up, usize, vp, vsize);
    138      1.1  mrg     }
    139      1.1  mrg 
    140      1.1  mrg   wsize -= cy_limb == 0;
    141      1.1  mrg 
    142  1.1.1.2  mrg   SIZ (w) = sign_product < 0 ? -wsize : wsize;
    143      1.1  mrg   if (free_me != NULL)
    144      1.1  mrg     (*__gmp_free_func) (free_me, free_me_size * BYTES_PER_MP_LIMB);
    145      1.1  mrg   TMP_FREE;
    146      1.1  mrg }
    147