Climb To A Prime
June 13, 2017
The central action of the program is to bring down powers and calculate the next number in the chain; we use flatten and uniq-c from the Standard Prelude:
(define (bring-down-powers fs)
(string->number
(apply string-append
(map number->string
(filter (lambda (n) (not (= n 1)))
(flatten
(uniq-c = fs)))))))
Then the program to climb to a prime just calls bring-down-powers repeatedly until it reaches a prime:
(define (climb-to-a-prime n)
(let loop ((n n))
(let ((fs (factors n)))
(display fs) (newline)
(when (pair? (cdr fs))
(loop (bring-down-powers fs))))))
The hard part of the program is factoring the numbers, which increase quickly. For demonstration, the code at ideone uses a 2,3,5-wheel; in practice, you will need something better. Here are some examples:
> (climb-to-a-prime 90) (2 3 3 5) (3 5 5 31) (7 7 719) (72719) > (climb-to-a-prime 284) (2 2 71) (3 757) (13 17 17) (2 2 37 89) (31 7219) (7 45317) (3 3 82813) (3 3 23 15859) (23 14013733) (3 3 29 8865953) (3 2963 3633577) (3 34327 3200917) (13 181 1420855589) (131811420855589)
And here’s the infinite loop:
> (climb-to-a-prime 13532385396179) (13 53 53 3853 96179) (13 53 53 3853 96179) (13 53 53 3853 96179) (13 53 53 3853 96179) ...
You can run the program at http://ideone.com/woHTqi.
Advertisements
Pages: 1 2
That was a very nice video indeed! Here is my take on the problem, using Julia. I could have written my own factoring method, but life is too short for replicating stuff that’s already on the base package…
function BringDown(x::Int64)
F = factor(x)
factors = collect(keys(F))
exponents = collect(values(F))
F = sortrows(hcat(factors, exponents))
n = length(factors)
Z = Array(AbstractString, n)
for i = 1:n
if F[i, 2] == 1
Z[i] = string(F[i, 1])
else
Z[i] = join(F[i,:])
end
end
return parse(Int64, join(Z))
end
function main(x::Int64)
y = string(x)
if !isprime(x)
if x != 13532385396179
while !isprime(x)
x = BringDown(x)
y = string(y, “, “, x)
end
else
return “This number cannot be processed!”
end
end
if y[end] == ‘ ‘
y = y[1:(end-2)]
end
return y
end
In Python. For n=20, my program runs rather long and I don’t know whether it will return (well, the conjecture is disproven in general, but still a curious on n=20 already). Perhaps a faster factorization algorithm might help.
from collections import Counter def factorize(number): # gives a list of all prime factors of number factors = [] while not (number % 2): number = number // 2 factors.append(2) i = 3 while i <= number**0.5 and number-1: if not (number % i): factors.append(i) number = number // i else: i += 2 if number != 1 and i >= number**0.5: factors.append(number) return factors def bring_down(factors): result = "" c = Counter(factors) for x in sorted(c): result += str(x) if c[x] > 1: result += str(c[x]) return int(result) for i in range(1,100): n = i factors = factorize(n) while len(factors) > 1: n_new = bring_down(factors) if n_new == n: print('-- recursive case') break n = n_new factors = factorize(n) print(i, ' result: ', n) #output # 1 result: 1 # 2 result: 2 # 3 result: 3 # 4 result: 211 # 5 result: 5 # 6 result: 23 # 7 result: 7 # 8 result: 23 # 9 result: 2213 # 10 result: 2213 # 11 result: 11 # 12 result: 223 # 13 result: 13 # 14 result: 311 # 15 result: 1129 # 16 result: 233 # 17 result: 17 # 18 result: 17137 # 19 result: 19In Python. I added function to find the counter example. It prints the counterexample after a few seconds.
def f(n): if not is_prime(n): nums = [] for k, g in itertools.groupby(brent_factors(n)): nums.append(k) L = len(list(g)) if L > 1: nums.append(L) n = int("".join(str(i) for i in nums)) return n def conway(n): while not is_prime(n): oldn, n = n, f(n) if oldn == n: raise ValueError("Failed for {}".format(n)) return n def find_counter_example(): """ try n = x*p = f(x) * 10^y + p with x = m * 10^y + 1 then p = fx / m see: http://aperiodical.com/ """ for y in range(2, 10): print("y=", y) ty = 10 ** y for m in range(2, 10**4): x = m * ty + 1 fx = f(x) p = fx // m cand = p * x if fx == m * p and is_prime(p) and f(cand) == cand: return cand ans = find_counter_example() print(ans) print(f(ans))@Rutger: Your program at 20 will run for a while: A195265 gives the first 110 values in the sequence for 20, and ends in a 178-digit composite that has not yet been factored. At A195264 you will find everything that is known about all numbers up to 10000.
Here’s a Haskell version.
import Data.Bool (bool) import Data.List (foldl', unfoldr) import Data.Tuple (swap) import Math.NumberTheory.Primes (factorise, isPrime) import Text.Printf (printf) import System.Environment (getArgs) -- Return the "climb to a prime" list for a given number. climb :: Integer -> [Integer] climb n = let ns = iterate bringDown n in n : map snd (takeWhile (uncurry (/=)) $ zip ns (tail ns)) -- Perform one step of "climb to a prime" by converting a number's sequence of -- factors and their multiplicities into a new number. bringDown :: Integer -> Integer bringDown = fromDigits . concatMap toDigits . concatMap toNums . factorise -- Convert a factor and its multiplicity to a list of two numbers. toNums :: (Integer, Int) -> [Integer] toNums (x, 1) = [x] toNums (x, y) = [x, fromIntegral y] -- Return the sequence of digits corresponding to a number. -- E.g. toDigits 123 = [1,2,3] toDigits :: Integral a => a -> [a] toDigits = reverse . unfoldr step where step 0 = Nothing step n = Just $ swap $ n `quotRem` 10 -- Return the number corresponding to a sequence of digits. -- E.g. fromDigits [1,2,3] = 123 fromDigits :: Integral a => [a] -> a fromDigits [] = 0 fromDigits ns = foldl' (\p d -> 10*p + d) 0 ns printClimb :: Integer -> IO () printClimb n = let ms = climb n m = last ms in printf "%d => %s: climbed to a %sprime\n" n (show ms) (bool "non-" "" (isPrime m)) main :: IO () main = getArgs >>= mapM_ (printClimb . read)Hmm, I just noticed that my first case in fromDigits is redundant. The function could be written more succinctly as: