Home | History | Annotate | Line # | Download | only in mpz
mul.c revision 1.1.1.1
      1  1.1  mrg /* mpz_mul -- Multiply two integers.
      2  1.1  mrg 
      3  1.1  mrg Copyright 1991, 1993, 1994, 1996, 2000, 2001, 2005, 2009, 2011 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 #ifdef BERKELEY_MP
     25  1.1  mrg #include "mp.h"
     26  1.1  mrg #endif
     27  1.1  mrg 
     28  1.1  mrg 
     29  1.1  mrg void
     30  1.1  mrg #ifndef BERKELEY_MP
     31  1.1  mrg mpz_mul (mpz_ptr w, mpz_srcptr u, mpz_srcptr v)
     32  1.1  mrg #else /* BERKELEY_MP */
     33  1.1  mrg mult (mpz_srcptr u, mpz_srcptr v, mpz_ptr w)
     34  1.1  mrg #endif /* BERKELEY_MP */
     35  1.1  mrg {
     36  1.1  mrg   mp_size_t usize;
     37  1.1  mrg   mp_size_t vsize;
     38  1.1  mrg   mp_size_t wsize;
     39  1.1  mrg   mp_size_t sign_product;
     40  1.1  mrg   mp_ptr up, vp;
     41  1.1  mrg   mp_ptr wp;
     42  1.1  mrg   mp_ptr free_me;
     43  1.1  mrg   size_t free_me_size;
     44  1.1  mrg   mp_limb_t cy_limb;
     45  1.1  mrg   TMP_DECL;
     46  1.1  mrg 
     47  1.1  mrg   usize = SIZ (u);
     48  1.1  mrg   vsize = SIZ (v);
     49  1.1  mrg   sign_product = usize ^ vsize;
     50  1.1  mrg   usize = ABS (usize);
     51  1.1  mrg   vsize = ABS (vsize);
     52  1.1  mrg 
     53  1.1  mrg   if (usize < vsize)
     54  1.1  mrg     {
     55  1.1  mrg       MPZ_SRCPTR_SWAP (u, v);
     56  1.1  mrg       MP_SIZE_T_SWAP (usize, vsize);
     57  1.1  mrg     }
     58  1.1  mrg 
     59  1.1  mrg   if (vsize == 0)
     60  1.1  mrg     {
     61  1.1  mrg       SIZ(w) = 0;
     62  1.1  mrg       return;
     63  1.1  mrg     }
     64  1.1  mrg 
     65  1.1  mrg #if HAVE_NATIVE_mpn_mul_2
     66  1.1  mrg   if (vsize <= 2)
     67  1.1  mrg     {
     68  1.1  mrg       MPZ_REALLOC (w, usize+vsize);
     69  1.1  mrg       wp = PTR(w);
     70  1.1  mrg       if (vsize == 1)
     71  1.1  mrg         cy_limb = mpn_mul_1 (wp, PTR(u), usize, PTR(v)[0]);
     72  1.1  mrg       else
     73  1.1  mrg         {
     74  1.1  mrg           cy_limb = mpn_mul_2 (wp, PTR(u), usize, PTR(v));
     75  1.1  mrg           usize++;
     76  1.1  mrg         }
     77  1.1  mrg       wp[usize] = cy_limb;
     78  1.1  mrg       usize += (cy_limb != 0);
     79  1.1  mrg       SIZ(w) = (sign_product >= 0 ? usize : -usize);
     80  1.1  mrg       return;
     81  1.1  mrg     }
     82  1.1  mrg #else
     83  1.1  mrg   if (vsize == 1)
     84  1.1  mrg     {
     85  1.1  mrg       MPZ_REALLOC (w, usize+1);
     86  1.1  mrg       wp = PTR(w);
     87  1.1  mrg       cy_limb = mpn_mul_1 (wp, PTR(u), usize, PTR(v)[0]);
     88  1.1  mrg       wp[usize] = cy_limb;
     89  1.1  mrg       usize += (cy_limb != 0);
     90  1.1  mrg       SIZ(w) = (sign_product >= 0 ? usize : -usize);
     91  1.1  mrg       return;
     92  1.1  mrg     }
     93  1.1  mrg #endif
     94  1.1  mrg 
     95  1.1  mrg   TMP_MARK;
     96  1.1  mrg   free_me = NULL;
     97  1.1  mrg   up = PTR(u);
     98  1.1  mrg   vp = PTR(v);
     99  1.1  mrg   wp = PTR(w);
    100  1.1  mrg 
    101  1.1  mrg   /* Ensure W has space enough to store the result.  */
    102  1.1  mrg   wsize = usize + vsize;
    103  1.1  mrg   if (ALLOC(w) < wsize)
    104  1.1  mrg     {
    105  1.1  mrg       if (wp == up || wp == vp)
    106  1.1  mrg 	{
    107  1.1  mrg 	  free_me = wp;
    108  1.1  mrg 	  free_me_size = ALLOC(w);
    109  1.1  mrg 	}
    110  1.1  mrg       else
    111  1.1  mrg 	(*__gmp_free_func) (wp, ALLOC(w) * BYTES_PER_MP_LIMB);
    112  1.1  mrg 
    113  1.1  mrg       ALLOC(w) = wsize;
    114  1.1  mrg       wp = (mp_ptr) (*__gmp_allocate_func) (wsize * BYTES_PER_MP_LIMB);
    115  1.1  mrg       PTR(w) = wp;
    116  1.1  mrg     }
    117  1.1  mrg   else
    118  1.1  mrg     {
    119  1.1  mrg       /* Make U and V not overlap with W.  */
    120  1.1  mrg       if (wp == up)
    121  1.1  mrg 	{
    122  1.1  mrg 	  /* W and U are identical.  Allocate temporary space for U.  */
    123  1.1  mrg 	  up = TMP_ALLOC_LIMBS (usize);
    124  1.1  mrg 	  /* Is V identical too?  Keep it identical with U.  */
    125  1.1  mrg 	  if (wp == vp)
    126  1.1  mrg 	    vp = up;
    127  1.1  mrg 	  /* Copy to the temporary space.  */
    128  1.1  mrg 	  MPN_COPY (up, wp, usize);
    129  1.1  mrg 	}
    130  1.1  mrg       else if (wp == vp)
    131  1.1  mrg 	{
    132  1.1  mrg 	  /* W and V are identical.  Allocate temporary space for V.  */
    133  1.1  mrg 	  vp = TMP_ALLOC_LIMBS (vsize);
    134  1.1  mrg 	  /* Copy to the temporary space.  */
    135  1.1  mrg 	  MPN_COPY (vp, wp, vsize);
    136  1.1  mrg 	}
    137  1.1  mrg     }
    138  1.1  mrg 
    139  1.1  mrg   if (up == vp)
    140  1.1  mrg     {
    141  1.1  mrg       mpn_sqr (wp, up, usize);
    142  1.1  mrg       cy_limb = wp[wsize - 1];
    143  1.1  mrg     }
    144  1.1  mrg   else
    145  1.1  mrg     {
    146  1.1  mrg       cy_limb = mpn_mul (wp, up, usize, vp, vsize);
    147  1.1  mrg     }
    148  1.1  mrg 
    149  1.1  mrg   wsize -= cy_limb == 0;
    150  1.1  mrg 
    151  1.1  mrg   SIZ(w) = sign_product < 0 ? -wsize : wsize;
    152  1.1  mrg   if (free_me != NULL)
    153  1.1  mrg     (*__gmp_free_func) (free_me, free_me_size * BYTES_PER_MP_LIMB);
    154  1.1  mrg   TMP_FREE;
    155  1.1  mrg }
    156