This is an automated email from the git hooks/post-receive script. ppm-guest pushed a commit to annotated tag v0.38 in repository libmath-prime-util-perl.
commit f174ddec326be9d142b6db4e86d7f1ad4df73598 Author: Dana Jacobsen <d...@acm.org> Date: Thu Jan 30 09:16:07 2014 -0800 Change to is_power, still needs work --- TODO | 2 ++ XS.xs | 25 +++++++++---------------- lib/Math/Prime/Util.pm | 19 ++++++++++++++++--- util.c | 32 ++++++++++++++++++++------------ util.h | 3 +-- 5 files changed, 48 insertions(+), 33 deletions(-) diff --git a/TODO b/TODO index 31f4abc..21f8915 100644 --- a/TODO +++ b/TODO @@ -78,3 +78,5 @@ - consider following Pari and doing is_power($n) / is_power($n, $k) instead of the three: is_perfect_power, is_perfect_square, is_perfect_cube. Document power. Tests for power. + From Sage: ispower(3089265681159475043336839581081873360674602365963130114355701114591322241990483812812582393906477998611814245513881) should return 14. + Pari's ispower had some problems in 2.4 diff --git a/XS.xs b/XS.xs index eba9916..666e28a 100644 --- a/XS.xs +++ b/XS.xs @@ -13,7 +13,6 @@ #include "cache.h" #include "sieve.h" #define FUNC_gcd_ui 1 -#define FUNC_is_perfect_square 1 #include "util.h" #include "primality.h" #include "factor.h" @@ -542,11 +541,9 @@ is_prime(IN SV* svn, ...) is_extra_strong_lucas_pseudoprime = 5 is_frobenius_underwood_pseudoprime = 6 is_aks_prime = 7 - is_perfect_square = 8 - is_perfect_cube = 9 - is_perfect_power = 10 - is_pseudoprime = 11 - is_almost_extra_strong_lucas_pseudoprime = 12 + is_power = 8 + is_pseudoprime = 9 + is_almost_extra_strong_lucas_pseudoprime = 10 PREINIT: int status; PPCODE: @@ -565,11 +562,9 @@ is_prime(IN SV* svn, ...) case 5: ret = _XS_is_lucas_pseudoprime(n, 2); break; case 6: ret = _XS_is_frobenius_underwood_pseudoprime(n); break; case 7: ret = _XS_is_aks_prime(n); break; - case 8: ret = is_perfect_square(n); break; - case 9: ret = is_perfect_cube(n); break; - case 10: ret = is_perfect_power(n); break; - case 11: ret = _XS_is_pseudoprime(n, (items == 1) ? 2 : a); break; - case 12: + case 8: ret = (a == 0) ? is_power(n) : !(is_power(n) % a); break; + case 9: ret = _XS_is_pseudoprime(n, (items == 1) ? 2 : a); break; + case 10: default: ret = _XS_is_almost_extra_strong_lucas_pseudoprime (n, (items == 1) ? 1 : a); break; } @@ -585,11 +580,9 @@ is_prime(IN SV* svn, ...) case 5: _vcallsub_with_gmp("is_extra_strong_lucas_pseudoprime"); break; case 6: _vcallsub_with_gmp("is_frobenius_underwood_pseudoprime"); break; case 7: _vcallsub_with_gmp("is_aks_prime"); break; - case 8: _vcallsub_with_gmp("is_perfect_square"); break; - case 9: _vcallsub_with_gmp("is_perfect_cube"); break; - case 10:_vcallsub_with_gmp("is_perfect_power"); break; - case 11:_vcallsub_with_gmp("is_pseudoprime"); break; - case 12: + case 8: _vcallsub_with_gmp("is_power"); break; + case 9:_vcallsub_with_gmp("is_pseudoprime"); break; + case 10: default:_vcallsub_with_gmp("is_almost_extra_strong_lucas_pseudoprime"); break; } return; /* skip implicit PUTBACK */ diff --git a/lib/Math/Prime/Util.pm b/lib/Math/Prime/Util.pm index 661daf0..51519b7 100644 --- a/lib/Math/Prime/Util.pm +++ b/lib/Math/Prime/Util.pm @@ -23,7 +23,7 @@ our @EXPORT_OK = is_almost_extra_strong_lucas_pseudoprime is_frobenius_underwood_pseudoprime is_aks_prime - is_perfect_power + is_power miller_rabin miller_rabin_random lucas_sequence primes @@ -1761,9 +1761,22 @@ Brent, 2010: "AKS is not a practical algorithm. ECPP is much faster." We have ECPP, and indeed it is much faster. -=head2 is_perfect_power +=head2 is_power -Given a positive integer input n, returns 1 if C<n = p^k> for C<k E<gt> 1>. + say "$n is a perfect square" if is_power($n, 2); + say "$n is a perfect cube" if is_power($n, 3); + say "$n is a ", is_power($n), "-th power"; + +Given a single positive integer input C<n>, returns k if C<n = p^k> for +some integer C<p E<gt> 1, k E<gt> 1>, and 0 otherwise. The k returned is +the largest possible. This can be used in a boolean statement to +determine if C<n> is a perfect power. + +If given two arguments C<n> and C<k>, returns 1 if C<n> is a C<k-th> power, +and 0 otherwise. + +This corresponds to Pari/GP's C<ispower> function, with the limitations of +only integer arguments and no third argument may be given returning the root. =head2 lucas_sequence diff --git a/util.c b/util.c index a31ac02..9ad1a76 100644 --- a/util.c +++ b/util.c @@ -1015,11 +1015,19 @@ static const UV perfect_powers[] = }; #define NPOWERS (sizeof(perfect_powers)/sizeof(perfect_powers[0])) -int is_perfect_power(UV n) { +/* TODO: This has to be redone to properly return the highest power */ +int is_power(UV n) { + int next; if ((n <= 3) || (n == UV_MAX)) return 0; - if ((n & (n-1)) == 0) return 1; /* powers of 2 */ - if (is_perfect_square(n)) return 1; /* perfect square */ - if (is_perfect_cube(n)) return 1; /* perfect cube */ + if ((n & (n-1)) == 0) return ctz(n); /* powers of 2 */ + if (is_perfect_square(n)) { + next = is_power(isqrt(n)); + return (next == 0) ? 2 : 2*next; + } + if (is_perfect_cube(n)) { + next = is_power(icbrt(n)); + return (next == 0) ? 3 : 3*next; + } { UV lo = 0; UV hi = NPOWERS-1; @@ -1033,17 +1041,17 @@ int is_perfect_power(UV n) { } #if BITS_PER_WORD == 64 { /* n > 2**32. If n = p^k, then p in (3 .. 7131) and k in (5,7,11,13) */ - int ib, k; + int ib; for (ib = 3; ib <= 6; ib++) { /* prime exponents from 5 to 13 */ - UV b = primes_small[ib]; + UV k, b = primes_small[ib]; UV root = (UV) ( powl(n, 1.0L / b ) + 0.01 ); - UV pk = 1; - while (b) { - if (b & 1) pk *= root; - b >>= 1; - if (b) root *= root; + UV pk = root * root * root * root * root; + for (k = 5; k < b; k++) + pk *= root; + if (n == pk) { + next = is_power(root); + return (next == 0) ? b : b*next; } - if (n == pk) return 1; } } #endif diff --git a/util.h b/util.h index 61627d7..1243fa0 100644 --- a/util.h +++ b/util.h @@ -21,8 +21,7 @@ extern UV prime_count_lower(UV x); extern UV prime_count_upper(UV x); extern UV prime_count_approx(UV x); -extern int is_perfect_cube(UV n); -extern int is_perfect_power(UV n); +extern int is_power(UV n); extern signed char* _moebius_range(UV low, UV high); extern UV* _totient_range(UV low, UV high); -- Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/pkg-perl/packages/libmath-prime-util-perl.git _______________________________________________ Pkg-perl-cvs-commits mailing list Pkg-perl-cvs-commits@lists.alioth.debian.org http://lists.alioth.debian.org/cgi-bin/mailman/listinfo/pkg-perl-cvs-commits