1 /* mpn_toom22_mul -- Multiply {ap,an} and {bp,bn} where an >= bn. Or more
2 accurately, bn <= an < 2bn.
4 Contributed to the GNU project by Torbjorn Granlund.
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, 2007, 2008 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 the GNU Lesser General Public License as published by
16 the Free Software Foundation; either version 3 of the License, or (at your
17 option) any later version.
19 The GNU MP Library is distributed in the hope that it will be useful, but
20 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
21 or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
22 License for more details.
24 You should have received a copy of the GNU Lesser General Public License
25 along with the GNU MP Library. If not, see http://www.gnu.org/licenses/. */
31 /* Evaluate in: -1, 0, +inf
39 v0 = a0 * b0 # A(0)*B(0)
40 vm1 = (a0- a1)*(b0- b1) # A(-1)*B(-1)
41 vinf= a1 * b1 # A(inf)*B(inf)
44 #if TUNE_PROGRAM_BUILD
45 #define MAYBE_mul_toom22 1
47 #define MAYBE_mul_toom22 \
48 (MUL_TOOM33_THRESHOLD >= 2 * MUL_TOOM22_THRESHOLD)
51 #define TOOM22_MUL_N_REC(p, a, b, n, ws) \
53 if (! MAYBE_mul_toom22 \
54 || BELOW_THRESHOLD (n, MUL_KARATSUBA_THRESHOLD)) \
55 mpn_mul_basecase (p, a, n, b, n); \
57 mpn_toom22_mul (p, a, n, b, n, ws); \
61 mpn_toom22_mul (mp_ptr pp,
62 mp_srcptr ap, mp_size_t an,
63 mp_srcptr bp, mp_size_t bn,
83 ASSERT (0 < s && s <= n);
84 ASSERT (0 < t && t <= s);
94 if (mpn_cmp (a0, a1, n) < 0)
96 mpn_sub_n (asm1, a1, a0, n);
101 mpn_sub_n (asm1, a0, a1, n);
106 if (mpn_zero_p (a0 + s, n - s) && mpn_cmp (a0, a1, s) < 0)
108 mpn_sub_n (asm1, a1, a0, s);
109 MPN_ZERO (asm1 + s, n - s);
114 mpn_sub (asm1, a0, n, a1, s);
121 if (mpn_cmp (b0, b1, n) < 0)
123 mpn_sub_n (bsm1, b1, b0, n);
128 mpn_sub_n (bsm1, b0, b1, n);
133 if (mpn_zero_p (b0 + t, n - t) && mpn_cmp (b0, b1, t) < 0)
135 mpn_sub_n (bsm1, b1, b0, t);
136 MPN_ZERO (bsm1 + t, n - t);
141 mpn_sub (bsm1, b0, n, b1, t);
145 #define v0 pp /* 2n */
146 #define vinf (pp + 2 * n) /* s+t */
147 #define vm1 scratch /* 2n */
148 #define scratch_out scratch + 2 * n
151 TOOM22_MUL_N_REC (vm1, asm1, bsm1, n, scratch_out);
153 /* vinf, s+t limbs */
154 mpn_mul (vinf, a1, s, b1, t);
157 TOOM22_MUL_N_REC (v0, ap, bp, n, scratch_out);
159 /* H(v0) + L(vinf) */
160 cy = mpn_add_n (pp + 2 * n, v0 + n, vinf, n);
163 cy2 = cy + mpn_add_n (pp + n, pp + 2 * n, v0, n);
165 /* L(vinf) + H(vinf) */
166 cy += mpn_add (pp + 2 * n, pp + 2 * n, n, vinf + n, s + t - n);
169 cy += mpn_add_n (pp + n, pp + n, vm1, 2 * n);
171 cy -= mpn_sub_n (pp + n, pp + n, vm1, 2 * n);
173 ASSERT (cy + 1 <= 3);
176 mpn_incr_u (pp + 2 * n, cy2);
177 if (LIKELY (cy <= 2))
178 mpn_incr_u (pp + 3 * n, cy);
180 mpn_decr_u (pp + 3 * n, 1);