// WilsonPrimes.java
// Andrew Davison, ad@coe.psu.ac.th, April 2023

/* Generate the Wilson primes for orders n = 1 to 11 and 
   primes < 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

   Usage:
     > compile WilsonPrimes.java
     > run WilsonPrimes
*/

import java.math.*;
import java.util.*;

import de.tilman_neumann.jml.combinatorics.Factorial;


public class WilsonPrimes
{
  private static final int LIMIT = 11000;


  public static void main(String[] args)
  {
    System.out.println("Generating factorials for 0-" + LIMIT + "...");
    BigInteger[] facts = new BigInteger[LIMIT];
    for (int i = 0; i < LIMIT; i++)
      facts[i] = Factorial.factorial(i);

    System.out.println("Generating primes <= " + LIMIT + "...");
    int[] primes = Nut.getPrimes(LIMIT);
    // int[] primes = Nut.getPrimes(LIMIT, Nut.loadPrimes("primes1M.txt"));


    System.out.println(" n | Wilson primes\n--------------------");
    BigInteger powMinus1 = BigInteger.valueOf(-1);
    for (int n = 1; n <= 11; ++n) {
      System.out.printf("%2d |", n);
      for (int p : primes)
        if ((p >= n) && isWilson(facts, n, p, powMinus1))
          System.out.printf(" %d", p);
      powMinus1 = powMinus1.negate();
      System.out.println();
    }
  } /* main */

 
  private static boolean isWilson(BigInteger[] facts, int n, int p,
                                        BigInteger powMinus1)
  // Does p^2 divide (n - 1)! * (p - n)! - (- 1)^n  ?
  {
    BigInteger wilson = facts[n-1].multiply(facts[p - n]).
                         subtract(powMinus1).
                           mod(BigInteger.valueOf(p * p));
    return wilson.equals(BigInteger.ZERO);
  }  // end of isWilson()


}  // end of WilsonPrimes class