Coverage Report

Created: 2024-06-28 06:19

/src/gmp-6.2.1/mpn/toom2_sqr.c
Line
Count
Source
1
/* mpn_toom2_sqr -- Square {ap,an}.
2
3
   Contributed to the GNU project by Torbjorn Granlund.
4
5
   THE FUNCTION IN THIS FILE IS INTERNAL WITH A MUTABLE INTERFACE.  IT IS ONLY
6
   SAFE TO REACH IT THROUGH DOCUMENTED INTERFACES.  IN FACT, IT IS ALMOST
7
   GUARANTEED THAT IT WILL CHANGE OR DISAPPEAR IN A FUTURE GNU MP RELEASE.
8
9
Copyright 2006-2010, 2012, 2014, 2018 Free Software Foundation, Inc.
10
11
This file is part of the GNU MP Library.
12
13
The GNU MP Library is free software; you can redistribute it and/or modify
14
it under the terms of either:
15
16
  * the GNU Lesser General Public License as published by the Free
17
    Software Foundation; either version 3 of the License, or (at your
18
    option) any later version.
19
20
or
21
22
  * the GNU General Public License as published by the Free Software
23
    Foundation; either version 2 of the License, or (at your option) any
24
    later version.
25
26
or both in parallel, as here.
27
28
The GNU MP Library is distributed in the hope that it will be useful, but
29
WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
30
or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
31
for more details.
32
33
You should have received copies of the GNU General Public License and the
34
GNU Lesser General Public License along with the GNU MP Library.  If not,
35
see https://www.gnu.org/licenses/.  */
36
37
38
#include "gmp-impl.h"
39
40
/* Evaluate in: -1, 0, +inf
41
42
  <-s--><--n-->
43
   ____ ______
44
  |_a1_|___a0_|
45
46
  v0  =  a0     ^2  #   A(0)^2
47
  vm1 = (a0- a1)^2  #  A(-1)^2
48
  vinf=      a1 ^2  # A(inf)^2
49
*/
50
51
#if TUNE_PROGRAM_BUILD || WANT_FAT_BINARY
52
#define MAYBE_sqr_toom2   1
53
#else
54
#define MAYBE_sqr_toom2             \
55
81.5M
  (SQR_TOOM3_THRESHOLD >= 2 * SQR_TOOM2_THRESHOLD)
56
#endif
57
58
#define TOOM2_SQR_REC(p, a, n, ws)          \
59
40.7M
  do {                 \
60
40.7M
    if (! MAYBE_sqr_toom2            \
61
40.7M
  || BELOW_THRESHOLD (n, SQR_TOOM2_THRESHOLD))      \
62
40.7M
      mpn_sqr_basecase (p, a, n);         \
63
40.7M
    else                \
64
40.7M
      mpn_toom2_sqr (p, a, n, ws);         \
65
40.7M
  } while (0)
66
67
void
68
mpn_toom2_sqr (mp_ptr pp,
69
         mp_srcptr ap, mp_size_t an,
70
         mp_ptr scratch)
71
13.5M
{
72
13.5M
  const int __gmpn_cpuvec_initialized = 1;
73
13.5M
  mp_size_t n, s;
74
13.5M
  mp_limb_t cy, cy2;
75
13.5M
  mp_ptr asm1;
76
77
32.8M
#define a0  ap
78
21.6M
#define a1  (ap + n)
79
80
13.5M
  s = an >> 1;
81
13.5M
  n = an - s;
82
83
13.5M
  ASSERT (0 < s && s <= n && s >= n - 1);
84
85
13.5M
  asm1 = pp;
86
87
  /* Compute asm1.  */
88
13.5M
  if (s == n)
89
7.95M
    {
90
7.95M
      if (mpn_cmp (a0, a1, n) < 0)
91
701k
  {
92
701k
    mpn_sub_n (asm1, a1, a0, n);
93
701k
  }
94
7.25M
      else
95
7.25M
  {
96
7.25M
    mpn_sub_n (asm1, a0, a1, n);
97
7.25M
  }
98
7.95M
    }
99
5.64M
  else /* n - s == 1 */
100
5.64M
    {
101
5.64M
      if (a0[s] == 0 && mpn_cmp (a0, a1, s) < 0)
102
18.5k
  {
103
18.5k
    mpn_sub_n (asm1, a1, a0, s);
104
18.5k
    asm1[s] = 0;
105
18.5k
  }
106
5.62M
      else
107
5.62M
  {
108
5.62M
    asm1[s] = a0[s] - mpn_sub_n (asm1, a0, a1, s);
109
5.62M
  }
110
5.64M
    }
111
112
27.1M
#define v0  pp        /* 2n */
113
27.1M
#define vinf  (pp + 2 * n)      /* s+s */
114
13.5M
#define vm1 scratch        /* 2n */
115
13.5M
#define scratch_out scratch + 2 * n
116
117
  /* vm1, 2n limbs */
118
13.5M
  TOOM2_SQR_REC (vm1, asm1, n, scratch_out);
119
120
  /* vinf, s+s limbs */
121
13.5M
  TOOM2_SQR_REC (vinf, a1, s, scratch_out);
122
123
  /* v0, 2n limbs */
124
13.5M
  TOOM2_SQR_REC (v0, ap, n, scratch_out);
125
126
  /* H(v0) + L(vinf) */
127
13.5M
  cy = mpn_add_n (pp + 2 * n, v0 + n, vinf, n);
128
129
  /* L(v0) + H(v0) */
130
13.5M
  cy2 = cy + mpn_add_n (pp + n, pp + 2 * n, v0, n);
131
132
  /* L(vinf) + H(vinf) */
133
13.5M
  cy += mpn_add (pp + 2 * n, pp + 2 * n, n, vinf + n, s + s - n);
134
135
13.5M
  cy -= mpn_sub_n (pp + n, pp + n, vm1, 2 * n);
136
137
13.5M
  ASSERT (cy + 1 <= 3);
138
13.5M
  ASSERT (cy2 <= 2);
139
140
13.5M
  if (LIKELY (cy <= 2)) {
141
13.5M
    MPN_INCR_U (pp + 2 * n, s + s, cy2);
142
13.5M
    MPN_INCR_U (pp + 3 * n, s + s - n, cy);
143
13.5M
  } else { /* cy is negative */
144
    /* The total contribution of v0+vinf-vm1 can not be negative. */
145
64.4k
#if WANT_ASSERT
146
    /* The borrow in cy stops the propagation of the carry cy2, */
147
64.4k
    ASSERT (cy2 == 1);
148
64.4k
    cy += mpn_add_1 (pp + 2 * n, pp + 2 * n, n, cy2);
149
64.4k
    ASSERT (cy == 0);
150
#else
151
    /* we simply fill the area with zeros. */
152
    MPN_FILL (pp + 2 * n, n, 0);
153
#endif
154
64.4k
  }
155
13.5M
}