/* mpn_mul -- Multiply two natural numbers. Copyright (C) 1991 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 the GNU General Public License as published by the Free Software Foundation; either version 1, 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. You should have received a copy of the GNU General Public License along with the GNU MP Library; see the file COPYING. If not, write to the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA. */ #include "gmp.h" #include "gmp-impl.h" #include "longlong.h" /* Multiply the natural numbers u (pointed to by UP, with USIZE limbs) and v (pointed to by VP, with VSIZE limbs), and store the result at PRODP. USIZE + VSIZE limbs are always stored, but if the input operands are normalized, the return value will reflect the true result size (which is either USIZE + VSIZE, or USIZE + VSIZE -1). NOTE: The space pointed to by PRODP is overwritten before finished with U and V, so overlap is an error. Argument constraints: 1. USIZE >= VSIZE. 2. PRODP != UP and PRODP != VP, i.e. the destination must be distinct from the multiplier and the multiplicand. */ mp_size_t mpn_mul_classic (mp_ptr prodp, mp_srcptr up, mp_size_t usize, mp_srcptr vp, mp_size_t vsize) { mp_size_t n; mp_size_t prod_size; mp_limb cy; int i, j; mp_limb prod_low, prod_high; mp_limb cy_dig; mp_limb v_limb, c; if (vsize == 0) return 0; /* Offset UP and PRODP so that the inner loop can be faster. */ up += usize; prodp += usize; /* Multiply by the first limb in V separately, as the result can be stored (not added) to PROD. We also avoid a loop for zeroing. */ v_limb = vp[0]; cy_dig = 0; j = -usize; do { umul_ppmm (prod_high, prod_low, up[j], v_limb); add_ssaaaa (cy_dig, prodp[j], prod_high, prod_low, 0, cy_dig); j++; } while (j < 0); prodp[j] = cy_dig; prodp++; /* For each iteration in the outer loop, multiply one limb from U with one limb from V, and add it to PROD. */ for (i = 1; i < vsize; i++) { v_limb = vp[i]; cy_dig = 0; j = -usize; /* Inner loops. Simulate the carry flag by jumping between these loops. The first is used when there was no carry in the previois iteration; the second when there was carry. */ do { umul_ppmm (prod_high, prod_low, up[j], v_limb); add_ssaaaa (cy_dig, prod_low, prod_high, prod_low, 0, cy_dig); c = prodp[j]; prod_low += c; prodp[j] = prod_low; if (prod_low < c) goto cy_loop; ncy_loop: j++; } while (j < 0); prodp[j] = cy_dig; prodp++; continue; do { umul_ppmm (prod_high, prod_low, up[j], v_limb); add_ssaaaa (cy_dig, prod_low, prod_high, prod_low, 0, cy_dig); c = prodp[j]; prod_low += c + 1; prodp[j] = prod_low; if (prod_low > c) goto ncy_loop; cy_loop: j++; } while (j < 0); cy_dig += 1; prodp[j] = cy_dig; prodp++; } return usize + vsize - (cy_dig == 0); }