# WilsonPrimes.py
'''
   Generate the Wilson isPrimes for orders n = 1 to 11 and 
   isPrimes < 11,000. 

   A Wilson prime of order n is a prime number p such that:
       p^2 divides (n - 1)! * (p - n)! - (- 1)^n

   If n is 1, the formula reduces to (p - n)! + 1
   where the only known p's are 5, 13, and 563.

   https://en.wikipedia.org/wiki/Wilson_prime
   and
   https://rosettacode.org/wiki/Wilson_primes_of_order_n#Java
'''

import math
from prime import *


def isWilson(facts, n, p, powMinus1):
  # Does p^2 divide (n - 1)! * (p - n)! - (- 1)^n  ?
  wilson = (facts[n-1]*(facts[p - n]) - powMinus1) % (p * p)
  return (wilson == 0)



LIMIT = 11000


# ------ main ----------

facts = [0]*LIMIT
facts[0] = 1
for i in range(1, LIMIT):
  facts[i] = i*facts[i-1]

prims = primes(LIMIT)

print(" n | Wilson Primes")
powMinus1 = -1
for n in range(1,12):
  print(f"{n:2d} | ", end='')
  for p in prims:
    if (p >= n) and isWilson(facts, n, p, powMinus1):
      print(p, end=' ')
  powMinus1 = -powMinus1
  print()