Coverage Report

Created: 2025-03-18 06:55

/src/gmp/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
0
  ((xsize) < (ysize)                                            \
37
0
   || ((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
0
#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
0
{
53
0
  mp_size_t  xsize, ysize, tsize, wsize, wsize_signed;
54
0
  mp_ptr     wp, tp;
55
0
  mp_limb_t  c;
56
0
  TMP_DECL;
57
58
  /* w unaffected if x==0 or y==0 */
59
0
  xsize = SIZ(x);
60
0
  ysize = SIZ(y);
61
0
  if (xsize == 0 || ysize == 0)
62
0
    return;
63
64
  /* make x the bigger of the two */
65
0
  if (ABS(ysize) > ABS(xsize))
66
0
    {
67
0
      MPZ_SRCPTR_SWAP (x, y);
68
0
      MP_SIZE_T_SWAP (xsize, ysize);
69
0
    }
70
71
0
  sub ^= ysize;
72
0
  ysize = ABS(ysize);
73
74
  /* use mpn_addmul_1/mpn_submul_1 if possible */
75
0
  if (ysize == 1)
76
0
    {
77
0
      mpz_aorsmul_1 (w, x, PTR(y)[0], sub);
78
0
      return;
79
0
    }
80
81
0
  sub ^= xsize;
82
0
  xsize = ABS(xsize);
83
84
0
  wsize_signed = SIZ(w);
85
0
  sub ^= wsize_signed;
86
0
  wsize = ABS(wsize_signed);
87
88
0
  tsize = xsize + ysize;
89
0
  wp = MPZ_REALLOC (w, MAX (wsize, tsize) + 1);
90
91
0
  if (wsize_signed == 0)
92
0
    {
93
0
      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
0
      if (x == y)
97
0
  {
98
0
    mpn_sqr (wp, PTR(x),xsize);
99
0
    high = wp[tsize-1];
100
0
  }
101
0
      else
102
0
  high = mpn_mul (wp, PTR(x),xsize, PTR(y),ysize);
103
0
      tsize -= (high == 0);
104
0
      SIZ(w) = (sub >= 0 ? tsize : -tsize);
105
0
      return;
106
0
    }
107
108
0
  TMP_MARK;
109
0
  tp = TMP_ALLOC_LIMBS (tsize);
110
111
0
  if (x == y)
112
0
    {
113
0
      mpn_sqr (tp, PTR(x),xsize);
114
0
      tsize -= (tp[tsize-1] == 0);
115
0
    }
116
0
  else
117
0
    {
118
0
      mp_limb_t high;
119
0
      high = mpn_mul (tp, PTR(x),xsize, PTR(y),ysize);
120
0
      tsize -= (high == 0);
121
0
    }
122
0
  ASSERT (tp[tsize-1] != 0);
123
0
  if (sub >= 0)
124
0
    {
125
0
      mp_srcptr up    = wp;
126
0
      mp_size_t usize = wsize;
127
128
0
      if (usize < tsize)
129
0
  {
130
0
    up  = tp;
131
0
    usize = tsize;
132
0
    tp  = wp;
133
0
    tsize = wsize;
134
135
0
    wsize = usize;
136
0
  }
137
138
0
      c = mpn_add (wp, up,usize, tp,tsize);
139
0
      wp[wsize] = c;
140
0
      wsize += (c != 0);
141
0
    }
142
0
  else
143
0
    {
144
0
      mp_srcptr up    = wp;
145
0
      mp_size_t usize = wsize;
146
147
0
      if (mpn_cmp_twosizes_lt (up,usize, tp,tsize))
148
0
  {
149
0
    up  = tp;
150
0
    usize = tsize;
151
0
    tp  = wp;
152
0
    tsize = wsize;
153
154
0
    wsize = usize;
155
0
    wsize_signed = -wsize_signed;
156
0
  }
157
158
0
      ASSERT_NOCARRY (mpn_sub (wp, up,usize, tp,tsize));
159
0
      wsize = usize;
160
0
      MPN_NORMALIZE (wp, wsize);
161
0
    }
162
163
0
  SIZ(w) = (wsize_signed >= 0 ? wsize : -wsize);
164
165
0
  TMP_FREE;
166
0
}
167
168
169
void
170
mpz_addmul (mpz_ptr w, mpz_srcptr u, mpz_srcptr v)
171
0
{
172
0
  mpz_aorsmul (w, u, v, (mp_size_t) 0);
173
0
}
174
175
void
176
mpz_submul (mpz_ptr w, mpz_srcptr u, mpz_srcptr v)
177
0
{
178
0
  mpz_aorsmul (w, u, v, (mp_size_t) -1);
179
0
}