Coverage Report

Created: 2024-11-21 07:03

/src/libgmp/mpz/prodlimbs.c
Line
Count
Source
1
/* mpz_prodlimbs(RESULT, V, LEN) -- Set RESULT to V[0]*V[1]*...*V[LEN-1].
2
3
Contributed to the GNU project by Marco Bodrato.
4
5
THE FUNCTION IN THIS FILE IS INTERNAL WITH A MUTABLE INTERFACE.
6
IT IS ONLY SAFE TO REACH IT THROUGH DOCUMENTED INTERFACES.
7
IN FACT, IT IS ALMOST GUARANTEED THAT IT WILL CHANGE OR
8
DISAPPEAR IN A FUTURE GNU MP RELEASE.
9
10
Copyright 2010-2012 Free Software Foundation, Inc.
11
12
This file is part of the GNU MP Library.
13
14
The GNU MP Library is free software; you can redistribute it and/or modify
15
it under the terms of either:
16
17
  * the GNU Lesser General Public License as published by the Free
18
    Software Foundation; either version 3 of the License, or (at your
19
    option) any later version.
20
21
or
22
23
  * the GNU General Public License as published by the Free Software
24
    Foundation; either version 2 of the License, or (at your option) any
25
    later version.
26
27
or both in parallel, as here.
28
29
The GNU MP Library is distributed in the hope that it will be useful, but
30
WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
31
or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
32
for more details.
33
34
You should have received copies of the GNU General Public License and the
35
GNU Lesser General Public License along with the GNU MP Library.  If not,
36
see https://www.gnu.org/licenses/.  */
37
38
#include "gmp-impl.h"
39
40
/*********************************************************/
41
/* Section list-prod: product of a list -> mpz_t         */
42
/*********************************************************/
43
44
/* FIXME: should be tuned */
45
#ifndef RECURSIVE_PROD_THRESHOLD
46
#define RECURSIVE_PROD_THRESHOLD (MUL_TOOM22_THRESHOLD)
47
#endif
48
49
/* Computes the product of the j>1 limbs pointed by factors, puts the
50
 * result in x. It assumes that all limbs are non-zero. Above
51
 * Karatsuba's threshold it uses a binary splitting strategy, to gain
52
 * speed by the asymptotically fast multiplication algorithms.
53
 *
54
 * The list in  {factors, j} is overwritten.
55
 * Returns the size of the result
56
 */
57
58
mp_size_t
59
mpz_prodlimbs (mpz_ptr x, mp_ptr factors, mp_size_t j)
60
1.88k
{
61
1.88k
  mp_limb_t cy;
62
1.88k
  mp_size_t size, i;
63
1.88k
  mp_ptr    prod;
64
65
1.88k
  ASSERT (j > 1);
66
1.88k
  ASSERT (RECURSIVE_PROD_THRESHOLD > 3);
67
68
1.88k
  if (BELOW_THRESHOLD (j, RECURSIVE_PROD_THRESHOLD)) {
69
986
    j--;
70
986
    size = 1;
71
72
16.8k
    for (i = 1; i < j; i++)
73
15.8k
      {
74
15.8k
  cy = mpn_mul_1 (factors, factors, size, factors[i]);
75
15.8k
  factors[size] = cy;
76
15.8k
  size += cy != 0;
77
15.8k
      };
78
79
986
    prod = MPZ_NEWALLOC (x, size + 1);
80
81
986
    cy = mpn_mul_1 (prod, factors, size, factors[i]);
82
986
    prod[size] = cy;
83
986
    return SIZ (x) = size + (cy != 0);
84
986
  } else {
85
899
    mpz_t x1, x2;
86
899
    TMP_DECL;
87
88
899
    i = j >> 1;
89
899
    j -= i;
90
899
    TMP_MARK;
91
92
899
    MPZ_TMP_INIT (x2, j);
93
94
899
    PTR (x1) = factors + i;
95
899
    ALLOC (x1) = j;
96
899
    j = mpz_prodlimbs (x2, factors + i, j);
97
899
    i = mpz_prodlimbs (x1, factors, i);
98
899
    size = i + j;
99
899
    prod = MPZ_NEWALLOC (x, size);
100
899
    if (i >= j)
101
420
      cy = mpn_mul (prod, PTR(x1), i, PTR(x2), j);
102
479
    else
103
479
      cy = mpn_mul (prod, PTR(x2), j, PTR(x1), i);
104
899
    TMP_FREE;
105
106
899
    return SIZ (x) = size - (cy == 0);
107
899
  }
108
1.88k
}