CIS 301. Probabilistic program termination and correctness.

Consider the program below:

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

public class Prime_Gen {

  private static KeyboardReader keyboard = new KeyboardReader();

  public static void main(String[] args) {

    Random seed;
    BigInteger p,q,n,x;

    int l = keyboard.readInt("number of bits of prime to be
    generated: ");

    do {
        seed = new Random();
        p = new BigInteger(l, seed);
    } while (!p.isProbablePrime(50));

    System.out.println("prime p is " + p);
  }
}

• The method isProbablePrime can render incorrect replies, resulting in Prime_Gen being incorrect, although the probability of that happening is very small.
• Even if isProbablePrime were correct for all its inputs, we still could not prove termination of Prime_Gen with the tools of chapter 4. Observe that the program uses a random seed for generating a potential prime p. Thus, the while-statement above will simply not terminate if none of the randomly generated p turn out to be primes (according to isProbablePrime). However, one can estimate the number of primes between 1 and 2**l and that estimate is good enough to conclude that the program above terminates with probability 1. The pseudo-random number generator used in constructing the seed will "hit" a prime number within a practically acceptable expected number of trials. See"How many primes are there?" for details.

nunk% java Prime_Gen
number of bits of prime to be generated: 10
prime p is 521

number of bits of prime to be generated: 50
prime p is 779086272649297

number of bits of prime to be generated: 80
prime p is 294755374770727725802903

number of bits of prime to be generated: 512
prime p is 1582669745767556434598467840683850400995794483300909462432520859046879044255382564136703926567700854183918516779360732250925951294970703995515688910221737

We could implement the method isProbablePrime with the Miller-Rabin algorithm:

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

public class Miller_Rabin {

  private static KeyboardReader keyboard = new KeyboardReader();
  private static final BigInteger ZERO = new BigInteger("0");
  private static final BigInteger ONE = new BigInteger("1");
  private static final BigInteger TWO = new BigInteger("2");

  public static void main(String[] args) {
    String ns = keyboard.readString("number to be tested for primality: ");
    BigInteger n = new BigInteger(ns);
    int s = keyboard.readInt("number of tests: ");
    if (is_prime(n,s)) { System.out.println(n + " is prime with probability at least " + (1 - Math.pow(2,-s))); }
    else { System.out.println(n + " is definitely not prime"); }
  }

  public static boolean is_prime(BigInteger n, int s) {
    BigInteger a;
    for (int i = 1; i <= s; ++i) {
      Random seed = new Random();
      do { a = new BigInteger(n.bitLength()-1, seed);
      } while (a.compareTo(ZERO) <= 0 || a.compareTo(n) >= 0);
      if (witness(a,n)) { return false; }
    }
    return true;
  }

  public static boolean witness(BigInteger a, BigInteger n) {

    BigInteger b = n.subtract(ONE);
    int k = b.bitLength() - 1;
    BigInteger t = new BigInteger("1");
    BigInteger x;

    for (int i = k; i >= 0; --i) {

      x = t;
      t = t.modPow(TWO,n);

      if (t.equals(ONE) && !x.equals(ONE) && !x.equals(b)) return true;
      if (b.testBit(i)) t = t.multiply(a).mod(n);
    }

    if (!t.equals(ONE)) { return true; }
    return false;
  }

}

• The main method of this program reads in a number, n, to be tested for primality, and a number, s, of tests to be done; it then calls the method is_prime(n,s) which returns "false" if n is not prime, and has the property that a "true" return value for a non-prime n occurs with probability less than 2**(-s).
• The method is_prime(n,s) simple engages in s many calls to the method witness(a,n), where, for each call, a is randomly generated to be between 1 and n -1. The witness method returns "true" if a is conclusive evidence for n not being prime; otherwise, it returns "false".
• Among other things, the witness method computes a**(n-1) modulo n. Since n is really big, this cannot possibly be computed in a naive iterative way. The program uses a technique called iterative squaring. If we only focus on the aspect of computing "a**(n-1) modulo n", and if we present the core code without the syntax for the arbitrary precision integer class, the program looks like this:

          // b[1..k] is the binary representation of n - 1 and b[k] stores the
  	// most significant bit
  	t = 1;
  	i = k;
  	while (i >= 1) {
  	  t = (t * t) mod n;
  	  if (b[i] == 1) t = (t * a) mod n;
  	  i = i - 1;
          }
  

  Let us first consider an example simulation. To compute 3**22 modulo 23, note that the binary representation of 22 is 10110. So 3**22 modulo 23 may be written as ((((1**2 * 3)**2)**2 * 3)**2 * 3)**2 modulo 23. Note that this expression suggests the manner in which the program will compute this: the "net effect" exponent is 10110 in binary. Note that there are only two squaring that are not immediately multiplied by 3; they correspond to the two 0 bits in 10110. This computation is an instance of a**(p-1) modulo p, where p is prime and a != 0 modulo p. By Fermat's theorem, all such expressions evaluate to 1. This is one of the things the witness method checks!
• In showing the correctness of this program, we seek an informative invariant. Intuitively, we expect that t equals a**j modulo n, where j is b[k..i+1] in binary. However, this won't be true right before the while-statement! Thus, our invariant candidate has to be an implication: 0 <= i < k --> t = a**b[k..i+1] modulo n.
  • If we place this invariant candidate right before the while-statement and push it upwards, we obtain a true formula, namely "0 <= k < k --> ...", this being true as its premise is false.
  • If we place our invariant candidate right after the while-statement, together with the negation of the boolean guard, it does not imply the desired statement, that t equal a**b[k..1] modulo n. The reason is that we cannot infer that i+1 equals 1. This situation is similar to what we discovered in the program DOWNFAC on assignment #11. Thus, our invariant candidate is the conjunction i >= 0 and 0 <= i < k --> t = a**b[k..i+1] modulo n.
  • It remains to show that this is indeed an invariant for the while-statement in question. This is an example, where the truth of the boolean guard in the "invariant hypothesis" plays an active role in securing a crucial implication. Note that this while-body contains an if-statement, making this more complicated than the exercises in the textbook. Try it!