sidef-scripts/Math/miller-rabin_primality_test.sf

#!/usr/bin/ruby

# A simple implementation of the Miller-Rabin primality test.

# See aslo:
#   https://en.wikipedia.org/wiki/Miller%E2%80%93Rabin_primality_test
#   https://rosettacode.org/wiki/Miller%E2%80%93Rabin_primality_test

func miller_rabin(n, k=10) {

    return false if (n <= 1)
    return true  if (n == 2)
    return false if (n.is_even)

    var t = n-1
    var s = t.valuation(2)
    var d = t>>s

    k.times {
        var a = irand(2, t)
        var x = powmod(a, d, n)

        next if (x ~~ [1, t])

        (s-1).times {
            x.powmod!(2, n)
            return false if (x == 1)
            break if (x == t)
        }

        return false if (x != t)
    }

    return true
}

var numbers = [
    61_794_479,
    2867561004669023153611,
    803_086_491,
    171_659_461_843,
    902_802_468_593,
    3_539_679_283_117,
    12_905_496_217_051,
    103_497_586_783_721,
]

numbers.each { |n|
    var p = miller_rabin(n, 12)
    say ("#{n} is" + (p ? ' ' : ' NOT ') + 'prime')
    assert_eq(p, n.is_prime)
}