sidef-scripts/Math/repunits_from_repunits.sf

#!/usr/bin/ruby

# Daniel "Trizen" Șuteu
# Date: 14 December 2018
# https://github.com/trizen

# a(n) = least positive `k`, such that `((n^k - 1)/(n-1)) / (Prod_{p^e | k} (n^(p^e) - 1)/(n-1))` is a prime number.

# For example:
#   a(19) = 34

# Let:
#   t = (19^34 - 1)/18

# Since 34 = 2*17, we have:
#   a = (19^2 - 1)/18
#   b = (19^17 - 1)/18

# By definition, `t / (a*b)` is a prime number, meaning that `t` can be factorized as:
#   t = a * b * (t / (a*b))

# Expanded form:
#   (19^34 - 1)/18 = (19^2 - 1)/18 * (19^17 - 1)/18 * (((19^34 - 1)/18) / ((19^2 - 1)/18 * (19^17 - 1)/18))
#   (19^34 - 1)/18 = 20 * 304465936543600121441 * 274019342889240109297

# See also:
#   https://oeis.org/A320914
#   https://en.wikipedia.org/wiki/Repunit

for n in (2..40) {

    for k in (2 .. Inf) {

        var t = (n**k - 1)/(n-1)

        t /= k.factor_prod {|p,e|
            (n**(p**e) - 1)/(n-1)
        }

        t.is_int || next

        if (t.is_prob_prime) {
            say "a(#{n}) = #{k} -> #{n**k - 1 / (n-1)}"
            break
        }
    }
}

__END__
a(2) = 6 -> 63
a(3) = 6 -> 364
a(4) = 6 -> 1365
a(5) = 10 -> 2441406
a(6) = 6 -> 9331
a(7) = 6 -> 19608
a(8) = 87 -> 529335265084874036222038788611144721614948501328642578466091812607602878353993
a(9) = 6 -> 66430
a(10) = 10 -> 1111111111
a(11) = 10 -> 2593742460
a(12) = 10 -> 5628851293
a(13) = 6 -> 402234
a(14) = 14 -> 854769755812155
a(15) = 6 -> 813616
a(16) = 6 -> 1118481
a(17) = 14 -> 10523614159962558
a(18) = 6 -> 2000719
a(19) = 34 -> 1668591117276687645552547004361146262739540
a(20) = 10 -> 538947368421