Coverage Report

Created: 2024-11-21 07:03

/src/libgmp/mpz/aorsmul.c
Line
Count
Source (jump to first uncovered line)
1
/* mpz_addmul, mpz_submul -- add or subtract multiple.
2
3
Copyright 2001, 2004, 2005, 2012, 2022 Free Software Foundation, Inc.
4
5
This file is part of the GNU MP Library.
6
7
The GNU MP Library is free software; you can redistribute it and/or modify
8
it under the terms of either:
9
10
  * the GNU Lesser General Public License as published by the Free
11
    Software Foundation; either version 3 of the License, or (at your
12
    option) any later version.
13
14
or
15
16
  * the GNU General Public License as published by the Free Software
17
    Foundation; either version 2 of the License, or (at your option) any
18
    later version.
19
20
or both in parallel, as here.
21
22
The GNU MP Library is distributed in the hope that it will be useful, but
23
WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
24
or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
25
for more details.
26
27
You should have received copies of the GNU General Public License and the
28
GNU Lesser General Public License along with the GNU MP Library.  If not,
29
see https://www.gnu.org/licenses/.  */
30
31
#include "gmp-impl.h"
32
33
34
/* expecting x and y both with non-zero high limbs */
35
#define mpn_cmp_twosizes_lt(xp,xsize, yp,ysize)                 \
36
30
  ((xsize) < (ysize)                                            \
37
30
   || ((xsize) == (ysize) && mpn_cmp (xp, yp, xsize) < 0))
38
39
40
/* sub>=0 means an addmul w += x*y, sub<0 means a submul w -= x*y.
41
42
   The signs of w, x and y are fully accounted for by each flipping "sub".
43
44
   The sign of w is retained for the result, unless the absolute value
45
   submul underflows, in which case it flips.  */
46
47
static void __gmpz_aorsmul (REGPARM_3_1 (mpz_ptr w, mpz_srcptr x, mpz_srcptr y, mp_size_t sub)) REGPARM_ATTR (1);
48
72
#define mpz_aorsmul(w,x,y,sub)  __gmpz_aorsmul (REGPARM_3_1 (w, x, y, sub))
49
50
REGPARM_ATTR (1) static void
51
mpz_aorsmul (mpz_ptr w, mpz_srcptr x, mpz_srcptr y, mp_size_t sub)
52
72
{
53
72
  mp_size_t  xsize, ysize, tsize, wsize, wsize_signed;
54
72
  mp_ptr     wp, tp;
55
72
  mp_limb_t  c;
56
72
  TMP_DECL;
57
58
  /* w unaffected if x==0 or y==0 */
59
72
  xsize = SIZ(x);
60
72
  ysize = SIZ(y);
61
72
  if (xsize == 0 || ysize == 0)
62
6
    return;
63
64
  /* make x the bigger of the two */
65
66
  if (ABS(ysize) > ABS(xsize))
66
37
    {
67
37
      MPZ_SRCPTR_SWAP (x, y);
68
37
      MP_SIZE_T_SWAP (xsize, ysize);
69
37
    }
70
71
66
  sub ^= ysize;
72
66
  ysize = ABS(ysize);
73
74
  /* use mpn_addmul_1/mpn_submul_1 if possible */
75
66
  if (ysize == 1)
76
9
    {
77
9
      mpz_aorsmul_1 (w, x, PTR(y)[0], sub);
78
9
      return;
79
9
    }
80
81
57
  sub ^= xsize;
82
57
  xsize = ABS(xsize);
83
84
57
  wsize_signed = SIZ(w);
85
57
  sub ^= wsize_signed;
86
57
  wsize = ABS(wsize_signed);
87
88
57
  tsize = xsize + ysize;
89
57
  wp = MPZ_REALLOC (w, MAX (wsize, tsize) + 1);
90
91
57
  if (wsize_signed == 0)
92
4
    {
93
4
      mp_limb_t  high;
94
      /* Nothing to add to, just set w=x*y.  No w==x or w==y overlap here,
95
   since we know x,y!=0 but w==0.  */
96
4
      if (x == y)
97
0
  {
98
0
    mpn_sqr (wp, PTR(x),xsize);
99
0
    high = wp[tsize-1];
100
0
  }
101
4
      else
102
4
  high = mpn_mul (wp, PTR(x),xsize, PTR(y),ysize);
103
4
      tsize -= (high == 0);
104
4
      SIZ(w) = (sub >= 0 ? tsize : -tsize);
105
4
      return;
106
4
    }
107
108
53
  TMP_MARK;
109
53
  tp = TMP_ALLOC_LIMBS (tsize);
110
111
53
  if (x == y)
112
0
    {
113
0
      mpn_sqr (tp, PTR(x),xsize);
114
0
      tsize -= (tp[tsize-1] == 0);
115
0
    }
116
53
  else
117
53
    {
118
53
      mp_limb_t high;
119
53
      high = mpn_mul (tp, PTR(x),xsize, PTR(y),ysize);
120
53
      tsize -= (high == 0);
121
53
    }
122
53
  ASSERT (tp[tsize-1] != 0);
123
53
  if (sub >= 0)
124
23
    {
125
23
      mp_srcptr up    = wp;
126
23
      mp_size_t usize = wsize;
127
128
23
      if (usize < tsize)
129
19
  {
130
19
    up  = tp;
131
19
    usize = tsize;
132
19
    tp  = wp;
133
19
    tsize = wsize;
134
135
19
    wsize = usize;
136
19
  }
137
138
23
      c = mpn_add (wp, up,usize, tp,tsize);
139
23
      wp[wsize] = c;
140
23
      wsize += (c != 0);
141
23
    }
142
30
  else
143
30
    {
144
30
      mp_srcptr up    = wp;
145
30
      mp_size_t usize = wsize;
146
147
30
      if (mpn_cmp_twosizes_lt (up,usize, tp,tsize))
148
24
  {
149
24
    up  = tp;
150
24
    usize = tsize;
151
24
    tp  = wp;
152
24
    tsize = wsize;
153
154
24
    wsize = usize;
155
24
    wsize_signed = -wsize_signed;
156
24
  }
157
158
30
      ASSERT_NOCARRY (mpn_sub (wp, up,usize, tp,tsize));
159
30
      wsize = usize;
160
30
      MPN_NORMALIZE (wp, wsize);
161
30
    }
162
163
53
  SIZ(w) = (wsize_signed >= 0 ? wsize : -wsize);
164
165
53
  TMP_FREE;
166
53
}
167
168
169
void
170
mpz_addmul (mpz_ptr w, mpz_srcptr u, mpz_srcptr v)
171
33
{
172
33
  mpz_aorsmul (w, u, v, (mp_size_t) 0);
173
33
}
174
175
void
176
mpz_submul (mpz_ptr w, mpz_srcptr u, mpz_srcptr v)
177
39
{
178
39
  mpz_aorsmul (w, u, v, (mp_size_t) -1);
179
39
}