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