PSW_primality_test.sf

#!/usr/bin/ruby

# The PSW primality test, named after Carl Pomerance, John Selfridge, and Samuel Wagstaff.

# No counter-examples are known to this test.

# Algorithm: given an odd integer n, that is not a perfect power:
#   1. Perform a (strong) base-2 Fermat test.
#   2. Find the first P > 0 such that kronecker(P^2 + 4, n) = -1.
#   3. If the Lucas U function: U(P, -1, n+1) = 0 (mod n), then n is probably prime.

# See also:
#   https://en.wikipedia.org/wiki/Lucas_pseudoprime
#   https://en.wikipedia.org/wiki/Baillie%E2%80%93PSW_primality_test

func PSW_primality_test(n) {

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

    return false if n.is_even
    return false if n.is_power

    2.powmod(n-1, n) == 1 || return false

    var P = (1..Inf -> first_by {|k|
        kronecker(k*k + 4, n) == -1
    })

    lucasUmod(P, -1, n+1, n) == 0
}

var count = 0
var from  = 1
var upto  = 1e4

for n in (from .. upto) {

    if (PSW_primality_test(n)) {
        if (!n.is_prime) {
            say "Counter-example: #{n}"
        }
        ++count
    }
    elsif (n.is_prime) {
        say "Missed a prime: #{n}"
    }
}

say "There are #{count} primes between #{from} and #{upto}."