123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107 |
- #!/usr/bin/ruby
- # Daniel "Trizen" Șuteu
- # Date: 07 July 2022
- # https://github.com/trizen
- # Compute the n-th prime power, using binary search and the prime-power counting function.
- # See also:
- # https://oeis.org/A143039
- func prime_power_count_lower(n) {
- sum(1..n.ilog2, {|k|
- prime_count_lower(n.iroot(k))
- })
- }
- func prime_power_count_upper(n) {
- sum(1..n.ilog2, {|k|
- prime_count_upper(n.iroot(k))
- })
- }
- func nth_prime_power_lower(n) {
- bsearch_min(n, n.nth_prime_upper, {|k|
- prime_power_count_upper(k) <=> n
- })
- }
- func nth_prime_power_upper(n) {
- bsearch_max(n, n.nth_prime_upper, {|k|
- prime_power_count_lower(k) <=> n
- })
- }
- func nth_prime_power(n) {
- n == 0 && return 1 # not a prime power, but...
- n <= 0 && return NaN
- n == 1 && return 2
- var min = nth_prime_power_lower(n)
- var max = nth_prime_power_upper(n)
- var k = 0
- var c = 0
- loop {
- k = (min + max)>>1
- c = prime_power_count(k)
- if (abs(c - n) <= n.iroot(3)) {
- break
- }
- given (c <=> n) {
- when (+1) { max = k-1 }
- when (-1) { min = k+1 }
- else { break }
- }
- }
- while (!k.is_prime_power) {
- --k
- }
- while (c != n) {
- var j = (n <=> c)
- k += j
- c += j
- k += j while !k.is_prime_power
- }
- return k
- }
- for n in (1..10) {
- var pp = nth_prime_power(10**n)
- assert(pp.is_prime_power)
- assert_eq(pp.prime_power_count, 10**n)
- assert_eq(10**n -> nth_prime_power, pp)
- say "PP(10^#{n}) = #{pp}"
- }
- assert_eq(
- nth_prime_power.map(1..100),
- 100.by { .is_prime_power }
- )
- __END__
- PP(10^1) = 16
- PP(10^2) = 419
- PP(10^3) = 7517
- PP(10^4) = 103511
- PP(10^5) = 1295953
- PP(10^6) = 15474787
- PP(10^7) = 179390821
- PP(10^8) = 2037968761
- PP(10^9) = 22801415981
- PP(10^10) = 252096675073
- PP(10^11) = 2760723662941
- PP(10^12) = 29996212395727
- PP(10^13) = 323780470283789
- PP(10^14) = 3475385632514321
- PP(10^15) = 37124507635789309
|