sidef-scripts/Math/difference_of_powers_factorization_method.sf

#!/usr/bin/ruby

# Daniel "Trizen" Șuteu
# Date: 04 July 2019
# Edit: 28 July 2019
# https://github.com/trizen

# A simple factorization method for numbers that can be expressed as a difference of powers.

# Very effective for numbers of the form:
#
#   n^k - 1
#
# where k has many divisors.

# See also:
#   https://trizenx.blogspot.com/2019/08/special-purpose-factorization-algorithms.html

define MIN_FACTOR = 1e6  # ignore small factors

func diff_power_factorization (n) {

    var orig = n
    var params = []

    func f(r1, e1, r2, e2) {
        var factors = Set()

        var divs1 = divisors(e1)
        var divs2 = divisors(e2)

        for d1 in (divs1), d2 in (divs2), j in (1, -1) {
            var t = (powmod(r1, d1, n) - (j * r2**d2))
            var g = gcd(t, n)

            if (g>MIN_FACTOR && g<n) {
                factors << g
            }
        }

        factors.to_a
    }

    func check_power (r1, e1) {

        var u = r1**e1
        var dx = abs(u - n)

        if (dx.is_power) {

            var e2 = dx.perfect_power
            var r2 = dx.perfect_root

            params << [r1, e1, r2, e2]
        }
    }

    for r1 in (2 .. n.ilog2) {
        var t = n.ilog(r1)
        check_power(r1, t)
        check_power(r1, t+1)
    }

    for e1 in (2 .. n.ilog2) {
        var t = n.iroot(e1)
        check_power(t,   e1)
        check_power(t+1, e1)
    }

    var factors = []
    var divisors = params.uniq.map { f(_...) }.flat.sort.uniq

    divisors.each {|d|
        var g = gcd(n, d)

        if (g>MIN_FACTOR && g<n) {
            while (g `divides` n) {
                n /= g
                factors << g
            }
        }
    }

    factors << (orig / factors.prod)
    factors.sort
}

if (ARGV) {
    say diff_power_factorization(Num(ARGV[0]))
    return 1
}

say diff_power_factorization(1009**24 + 29**12)
say diff_power_factorization(1009**24 - 29**12)
say diff_power_factorization(2**256 - 1)
say diff_power_factorization(10**120 + 1)
say diff_power_factorization(10**120 - 1)
say diff_power_factorization(10**120 - 25)
say diff_power_factorization(10**105 - 1)
say diff_power_factorization(10**105 + 1)
say diff_power_factorization(10**120 - 2134**2)

__END__
[1074309286591662655506002, 1154140443257087164049583013000044736320575461201]
[1018052, 1018110, 1036459399053, 1036488923402, 1036518447751, 1074309285719975471632201]
[4294967295, 4294967297, 18446744073709551617, 340282366920938463463374607431768211457]
[100000001, 9999999900000001, 99999999000000009999999900000001, 10000000099999999999999989999999899999999000000000000000100000001]
[1000001, 99009901, 9999900001, 9999999999, 10000100001, 1000000000001, 9999000099990001, 10099989899000101, 100009999999899989999000000010001]
[999999999999999999999999999999999999999999999999999999999995, 1000000000000000000000000000000000000000000000000000000000005]
[9999999, 900900990991, 111111111111111, 900009090090909909099991, 1109988789001111109989898989900111110998878900111]
[10000001, 1098900989011, 90909090909091, 1099988890111109888900011, 910009191000909089989898989899909091000919100091]
[999999999999999999999999999999999999999999999999999999997866, 1000000000000000000000000000000000000000000000000000000002134]