1 /* mpn_toom_eval_pm2 -- Evaluate a polynomial in +2 and -2
3 Contributed to the GNU project by Niels Möller and Marco Bodrato
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.
9 Copyright 2009 Free Software Foundation, Inc.
11 This file is part of the GNU MP Library.
13 The GNU MP Library is free software; you can redistribute it and/or modify
14 it under the terms of the GNU Lesser General Public License as published by
15 the Free Software Foundation; either version 3 of the License, or (at your
16 option) any later version.
18 The GNU MP Library is distributed in the hope that it will be useful, but
19 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
20 or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
21 License for more details.
23 You should have received a copy of the GNU Lesser General Public License
24 along with the GNU MP Library. If not, see http://www.gnu.org/licenses/. */
29 /* DO_addlsh2(d,a,b,n,cy) computes cy,{d,n} <- {a,n} + 4*(cy,{b,n}), it
30 can be used as DO_addlsh2(d,a,d,n,d[n]), for accumulation on {d,n+1}. */
31 #if HAVE_NATIVE_mpn_addlsh2_n
32 #define DO_addlsh2(d, a, b, n, cy) \
35 (cy) += mpn_addlsh2_n(d, a, b, n); \
38 #if HAVE_NATIVE_mpn_addlsh_n
39 #define DO_addlsh2(d, a, b, n, cy) \
42 (cy) += mpn_addlsh_n(d, a, b, n, 2); \
45 /* The following is not a general substitute for addlsh2.
46 It is correct if d == b, but it is not if d == a. */
47 #define DO_addlsh2(d, a, b, n, cy) \
50 (cy) += mpn_lshift(d, b, n, 2); \
51 (cy) += mpn_add_n(d, d, a, n); \
56 /* Evaluates a polynomial of degree 2 < k < GMP_NUMB_BITS, in the
59 mpn_toom_eval_pm2 (mp_ptr xp2, mp_ptr xm2, unsigned k,
60 mp_srcptr xp, mp_size_t n, mp_size_t hn, mp_ptr tp)
67 ASSERT (k < GMP_NUMB_BITS);
72 /* The degree k is also the number of full-size coefficients, so
73 * that last coefficient, of size hn, starts at xp + k*n. */
76 DO_addlsh2 (xp2, xp + (k-2) * n, xp + k * n, hn, cy);
78 cy = mpn_add_1 (xp2 + hn, xp + (k-2) * n + hn, n - hn, cy);
79 for (i = k - 4; i >= 0; i -= 2)
80 DO_addlsh2 (xp2, xp + i * n, xp2, n, cy);
86 DO_addlsh2 (tp, xp + (k-2) * n, xp + k * n, n, cy);
87 for (i = k - 4; i >= 0; i -= 2)
88 DO_addlsh2 (tp, xp + i * n, tp, n, cy);
92 ASSERT_NOCARRY(mpn_lshift (tp , tp , n + 1, 1));
94 ASSERT_NOCARRY(mpn_lshift (xp2, xp2, n + 1, 1));
96 neg = (mpn_cmp (xp2, tp, n + 1) < 0) ? ~0 : 0;
98 #if HAVE_NATIVE_mpn_add_n_sub_n
100 mpn_add_n_sub_n (xp2, xm2, tp, xp2, n + 1);
102 mpn_add_n_sub_n (xp2, xm2, xp2, tp, n + 1);
103 #else /* !HAVE_NATIVE_mpn_add_n_sub_n */
105 mpn_sub_n (xm2, tp, xp2, n + 1);
107 mpn_sub_n (xm2, xp2, tp, n + 1);
109 mpn_add_n (xp2, xp2, tp, n + 1);
110 #endif /* !HAVE_NATIVE_mpn_add_n_sub_n */
112 ASSERT (xp2[n] < (1<<(k+2))-1);
113 ASSERT (xm2[n] < ((1<<(k+3))-1 - (1^k&1))/3);
115 neg ^= ((k & 1) - 1);