yann@1250
|
1 |
Original patch from: perfpow.c.diff
|
yann@1250
|
2 |
|
yann@1250
|
3 |
-= BEGIN original header =-
|
yann@1250
|
4 |
Copyright 1998, 1999, 2000, 2001, 2005, 2008 Free Software Foundation, Inc.
|
yann@1250
|
5 |
|
yann@1250
|
6 |
This file is part of the GNU MP Library.
|
yann@1250
|
7 |
|
yann@1250
|
8 |
The GNU MP Library is free software; you can redistribute it and/or modify
|
yann@1250
|
9 |
it under the terms of the GNU Lesser General Public License as published by
|
yann@1250
|
10 |
the Free Software Foundation; either version 3 of the License, or (at your
|
yann@1250
|
11 |
option) any later version.
|
yann@1250
|
12 |
|
yann@1250
|
13 |
The GNU MP Library is distributed in the hope that it will be useful, but
|
yann@1250
|
14 |
WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
|
yann@1250
|
15 |
or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
|
yann@1250
|
16 |
License for more details.
|
yann@1250
|
17 |
|
yann@1250
|
18 |
You should have received a copy of the GNU Lesser General Public License
|
yann@1250
|
19 |
along with the GNU MP Library. If not, see http://www.gnu.org/licenses/.
|
yann@1250
|
20 |
|
yann@1250
|
21 |
-= END original header =-
|
yann@1250
|
22 |
|
yann@1250
|
23 |
diff -durN gmp-4.2.4.orig/mpz/perfpow.c gmp-4.2.4/mpz/perfpow.c
|
yann@1250
|
24 |
--- gmp-4.2.4.orig/mpz/perfpow.c 2007-08-30 20:31:41.000000000 +0200
|
yann@1250
|
25 |
+++ gmp-4.2.4/mpz/perfpow.c 2009-03-08 18:36:16.000000000 +0100
|
yann@1250
|
26 |
@@ -1,7 +1,7 @@
|
yann@1250
|
27 |
/* mpz_perfect_power_p(arg) -- Return non-zero if ARG is a perfect power,
|
yann@1250
|
28 |
zero otherwise.
|
yann@1250
|
29 |
|
yann@1250
|
30 |
-Copyright 1998, 1999, 2000, 2001, 2005 Free Software Foundation, Inc.
|
yann@1250
|
31 |
+Copyright 1998, 1999, 2000, 2001, 2005, 2008 Free Software Foundation, Inc.
|
yann@1250
|
32 |
|
yann@1250
|
33 |
This file is part of the GNU MP Library.
|
yann@1250
|
34 |
|
yann@1250
|
35 |
@@ -59,6 +59,8 @@
|
yann@1250
|
36 |
#define SMALLEST_OMITTED_PRIME 1009
|
yann@1250
|
37 |
|
yann@1250
|
38 |
|
yann@1250
|
39 |
+#define POW2P(a) (((a) & ((a) - 1)) == 0)
|
yann@1250
|
40 |
+
|
yann@1250
|
41 |
int
|
yann@1250
|
42 |
mpz_perfect_power_p (mpz_srcptr u)
|
yann@1250
|
43 |
{
|
yann@1250
|
44 |
@@ -72,16 +74,13 @@
|
yann@1250
|
45 |
mp_size_t usize = SIZ (u);
|
yann@1250
|
46 |
TMP_DECL;
|
yann@1250
|
47 |
|
yann@1250
|
48 |
- if (usize == 0)
|
yann@1250
|
49 |
- return 1; /* consider 0 a perfect power */
|
yann@1250
|
50 |
+ if (mpz_cmpabs_ui (u, 1) <= 0)
|
yann@1250
|
51 |
+ return 1; /* -1, 0, and +1 are perfect powers */
|
yann@1250
|
52 |
|
yann@1250
|
53 |
n2 = mpz_scan1 (u, 0);
|
yann@1250
|
54 |
if (n2 == 1)
|
yann@1250
|
55 |
return 0; /* 2 divides exactly once. */
|
yann@1250
|
56 |
|
yann@1250
|
57 |
- if (n2 != 0 && (n2 & 1) == 0 && usize < 0)
|
yann@1250
|
58 |
- return 0; /* 2 has even multiplicity with negative U */
|
yann@1250
|
59 |
-
|
yann@1250
|
60 |
TMP_MARK;
|
yann@1250
|
61 |
|
yann@1250
|
62 |
uns = ABS (usize) - n2 / BITS_PER_MP_LIMB;
|
yann@1250
|
63 |
@@ -89,6 +88,14 @@
|
yann@1250
|
64 |
MPZ_TMP_INIT (u2, uns);
|
yann@1250
|
65 |
|
yann@1250
|
66 |
mpz_tdiv_q_2exp (u2, u, n2);
|
yann@1250
|
67 |
+ mpz_abs (u2, u2);
|
yann@1250
|
68 |
+
|
yann@1250
|
69 |
+ if (mpz_cmp_ui (u2, 1) == 0)
|
yann@1250
|
70 |
+ {
|
yann@1250
|
71 |
+ TMP_FREE;
|
yann@1250
|
72 |
+ /* factoring completed; consistent power */
|
yann@1250
|
73 |
+ return ! (usize < 0 && POW2P(n2));
|
yann@1250
|
74 |
+ }
|
yann@1250
|
75 |
|
yann@1250
|
76 |
if (isprime (n2))
|
yann@1250
|
77 |
goto n2prime;
|
yann@1250
|
78 |
@@ -97,6 +104,9 @@
|
yann@1250
|
79 |
{
|
yann@1250
|
80 |
prime = primes[i];
|
yann@1250
|
81 |
|
yann@1250
|
82 |
+ if (mpz_cmp_ui (u2, prime) < 0)
|
yann@1250
|
83 |
+ break;
|
yann@1250
|
84 |
+
|
yann@1250
|
85 |
if (mpz_divisible_ui_p (u2, prime)) /* divisible by this prime? */
|
yann@1250
|
86 |
{
|
yann@1250
|
87 |
rem = mpz_tdiv_q_ui (q, u2, prime * prime);
|
yann@1250
|
88 |
@@ -115,12 +125,6 @@
|
yann@1250
|
89 |
n++;
|
yann@1250
|
90 |
}
|
yann@1250
|
91 |
|
yann@1250
|
92 |
- if ((n & 1) == 0 && usize < 0)
|
yann@1250
|
93 |
- {
|
yann@1250
|
94 |
- TMP_FREE;
|
yann@1250
|
95 |
- return 0; /* even multiplicity with negative U, reject */
|
yann@1250
|
96 |
- }
|
yann@1250
|
97 |
-
|
yann@1250
|
98 |
n2 = gcd (n2, n);
|
yann@1250
|
99 |
if (n2 == 1)
|
yann@1250
|
100 |
{
|
yann@1250
|
101 |
@@ -128,10 +132,11 @@
|
yann@1250
|
102 |
return 0; /* we have multiplicity 1 of some factor */
|
yann@1250
|
103 |
}
|
yann@1250
|
104 |
|
yann@1250
|
105 |
- if (mpz_cmpabs_ui (u2, 1) == 0)
|
yann@1250
|
106 |
+ if (mpz_cmp_ui (u2, 1) == 0)
|
yann@1250
|
107 |
{
|
yann@1250
|
108 |
TMP_FREE;
|
yann@1250
|
109 |
- return 1; /* factoring completed; consistent power */
|
yann@1250
|
110 |
+ /* factoring completed; consistent power */
|
yann@1250
|
111 |
+ return ! (usize < 0 && POW2P(n2));
|
yann@1250
|
112 |
}
|
yann@1250
|
113 |
|
yann@1250
|
114 |
/* As soon as n2 becomes a prime number, stop factoring.
|
yann@1250
|
115 |
@@ -169,6 +174,10 @@
|
yann@1250
|
116 |
else
|
yann@1250
|
117 |
{
|
yann@1250
|
118 |
unsigned long int nth;
|
yann@1250
|
119 |
+
|
yann@1250
|
120 |
+ if (usize < 0 && POW2P(n2))
|
yann@1250
|
121 |
+ return 0;
|
yann@1250
|
122 |
+
|
yann@1250
|
123 |
/* We found some factors above. We just need to consider values of n
|
yann@1250
|
124 |
that divides n2. */
|
yann@1250
|
125 |
for (nth = 2; nth <= n2; nth++)
|
yann@1250
|
126 |
@@ -184,8 +193,11 @@
|
yann@1250
|
127 |
exact = mpz_root (q, u2, nth);
|
yann@1250
|
128 |
if (exact)
|
yann@1250
|
129 |
{
|
yann@1250
|
130 |
- TMP_FREE;
|
yann@1250
|
131 |
- return 1;
|
yann@1250
|
132 |
+ if (! (usize < 0 && POW2P(nth)))
|
yann@1250
|
133 |
+ {
|
yann@1250
|
134 |
+ TMP_FREE;
|
yann@1250
|
135 |
+ return 1;
|
yann@1250
|
136 |
+ }
|
yann@1250
|
137 |
}
|
yann@1250
|
138 |
if (mpz_cmp_ui (q, SMALLEST_OMITTED_PRIME) < 0)
|
yann@1250
|
139 |
{
|
yann@1250
|
140 |
@@ -199,6 +211,9 @@
|
yann@1250
|
141 |
}
|
yann@1250
|
142 |
|
yann@1250
|
143 |
n2prime:
|
yann@1250
|
144 |
+ if (usize < 0 && POW2P(n2))
|
yann@1250
|
145 |
+ return 0;
|
yann@1250
|
146 |
+
|
yann@1250
|
147 |
exact = mpz_root (NULL, u2, n2);
|
yann@1250
|
148 |
TMP_FREE;
|
yann@1250
|
149 |
return exact;
|