1 /* mpn_toom3_sqr -- Square {ap,an}.
3 Contributed to the GNU project by Torbjorn Granlund.
4 Additional improvements by Marco Bodrato.
6 THE FUNCTION IN THIS FILE IS INTERNAL WITH A MUTABLE INTERFACE. IT IS ONLY
7 SAFE TO REACH IT THROUGH DOCUMENTED INTERFACES. IN FACT, IT IS ALMOST
8 GUARANTEED THAT IT WILL CHANGE OR DISAPPEAR IN A FUTURE GNU MP RELEASE.
10 Copyright 2006-2010, 2012 Free Software Foundation, Inc.
12 This file is part of the GNU MP Library.
14 The GNU MP Library is free software; you can redistribute it and/or modify
15 it under the terms of either:
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.
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
27 or both in parallel, as here.
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
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/. */
42 /* Evaluate in: -1, 0, +1, +2, +inf
49 v1 = (a0+ a1+ a2)^2 # A(1)^2 ah <= 2
50 vm1 = (a0- a1+ a2)^2 # A(-1)^2 |ah| <= 1
51 v2 = (a0+2a1+4a2)^2 # A(2)^2 ah <= 6
52 vinf= a2 ^2 # A(inf)^2
55 #if TUNE_PROGRAM_BUILD || WANT_FAT_BINARY
56 #define MAYBE_sqr_basecase 1
57 #define MAYBE_sqr_toom3 1
59 #define MAYBE_sqr_basecase \
60 (SQR_TOOM3_THRESHOLD < 3 * SQR_TOOM2_THRESHOLD)
61 #define MAYBE_sqr_toom3 \
62 (SQR_TOOM4_THRESHOLD >= 3 * SQR_TOOM3_THRESHOLD)
65 #define TOOM3_SQR_REC(p, a, n, ws) \
67 if (MAYBE_sqr_basecase \
68 && BELOW_THRESHOLD (n, SQR_TOOM2_THRESHOLD)) \
69 mpn_sqr_basecase (p, a, n); \
70 else if (! MAYBE_sqr_toom3 \
71 || BELOW_THRESHOLD (n, SQR_TOOM3_THRESHOLD)) \
72 mpn_toom2_sqr (p, a, n, ws); \
74 mpn_toom3_sqr (p, a, n, ws); \
78 mpn_toom3_sqr (mp_ptr pp,
79 mp_srcptr ap, mp_size_t an,
82 const int __gmpn_cpuvec_initialized = 1;
86 mp_ptr as1, asm1, as2;
92 n = (an + 2) / (size_t) 3;
96 ASSERT (0 < s && s <= n);
98 as1 = scratch + 4 * n + 4;
99 asm1 = scratch + 2 * n + 2;
104 /* Compute as1 and asm1. */
105 cy = mpn_add (gp, a0, n, a2, s);
106 #if HAVE_NATIVE_mpn_add_n_sub_n
107 if (cy == 0 && mpn_cmp (gp, a1, n) < 0)
109 cy = mpn_add_n_sub_n (as1, asm1, a1, gp, n);
116 cy2 = mpn_add_n_sub_n (as1, asm1, gp, a1, n);
117 as1[n] = cy + (cy2 >> 1);
118 asm1[n] = cy - (cy2 & 1);
121 as1[n] = cy + mpn_add_n (as1, gp, a1, n);
122 if (cy == 0 && mpn_cmp (gp, a1, n) < 0)
124 mpn_sub_n (asm1, a1, gp, n);
129 cy -= mpn_sub_n (asm1, gp, a1, n);
135 #if HAVE_NATIVE_mpn_rsblsh1_n
136 cy = mpn_add_n (as2, a2, as1, s);
138 cy = mpn_add_1 (as2 + s, as1 + s, n - s, cy);
140 cy = 2 * cy + mpn_rsblsh1_n (as2, a0, as2, n);
142 #if HAVE_NATIVE_mpn_addlsh1_n
143 cy = mpn_addlsh1_n (as2, a1, a2, s);
145 cy = mpn_add_1 (as2 + s, a1 + s, n - s, cy);
146 cy = 2 * cy + mpn_addlsh1_n (as2, a0, as2, n);
148 cy = mpn_add_n (as2, a2, as1, s);
150 cy = mpn_add_1 (as2 + s, as1 + s, n - s, cy);
152 cy = 2 * cy + mpn_lshift (as2, as2, n, 1);
153 cy -= mpn_sub_n (as2, as2, a0, n);
158 ASSERT (as1[n] <= 2);
159 ASSERT (asm1[n] <= 1);
161 #define v0 pp /* 2n */
162 #define v1 (pp + 2 * n) /* 2n+1 */
163 #define vinf (pp + 4 * n) /* s+s */
164 #define vm1 scratch /* 2n+1 */
165 #define v2 (scratch + 2 * n + 1) /* 2n+2 */
166 #define scratch_out (scratch + 5 * n + 5)
168 /* vm1, 2n+1 limbs */
169 #ifdef SMALLER_RECURSION
170 TOOM3_SQR_REC (vm1, asm1, n, scratch_out);
173 cy = asm1[n] + mpn_add_n (vm1 + n, vm1 + n, asm1, n);
175 cy += mpn_add_n (vm1 + n, vm1 + n, asm1, n);
178 TOOM3_SQR_REC (vm1, asm1, n + 1, scratch_out);
181 TOOM3_SQR_REC (v2, as2, n + 1, scratch_out); /* v2, 2n+1 limbs */
183 TOOM3_SQR_REC (vinf, a2, s, scratch_out); /* vinf, s+s limbs */
185 vinf0 = vinf[0]; /* v1 overlaps with this */
187 #ifdef SMALLER_RECURSION
189 TOOM3_SQR_REC (v1, as1, n, scratch_out);
192 cy = as1[n] + mpn_add_n (v1 + n, v1 + n, as1, n);
194 else if (as1[n] != 0)
196 #if HAVE_NATIVE_mpn_addlsh1_n
197 cy = 2 * as1[n] + mpn_addlsh1_n (v1 + n, v1 + n, as1, n);
199 cy = 2 * as1[n] + mpn_addmul_1 (v1 + n, as1, n, CNST_LIMB(2));
206 cy += mpn_add_n (v1 + n, v1 + n, as1, n);
208 else if (as1[n] != 0)
210 #if HAVE_NATIVE_mpn_addlsh1_n
211 cy += mpn_addlsh1_n (v1 + n, v1 + n, as1, n);
213 cy += mpn_addmul_1 (v1 + n, as1, n, CNST_LIMB(2));
219 TOOM3_SQR_REC (v1, as1, n + 1, scratch_out);
223 TOOM3_SQR_REC (v0, ap, n, scratch_out); /* v0, 2n limbs */
225 mpn_toom_interpolate_5pts (pp, v2, vm1, n, s + s, 0, vinf0);