(* primes.sml *)

(* implementation of generation of primes *)

structure Primes :> PRIMES =
struct

(* val divisible : int * int list -> bool

   if ms is an ascending sequence of primes with no gaps, and, if ms
   <> [], then l is not divisible by any primes < hd ms, and l is >
   the last element of ms, then divisible(l, ms) tests whether l is
   divisible by an element of ms

   it returns true as soon as the next element divides l; otherwise it
   returns false if the square of that next element is > l

   tail recursive, using tail of second argument *)

fun divisible(_, nil)     = false
  | divisible(l, m :: ms) =
      l mod m = 0 orelse
      (m * m < l andalso divisible(l, ms))

(* val next :
         int list * int * int list * int ->
         int list * int * int list * int

   if ms @ rev ls are the first length ms + length ls primes, listed
   in ascending order, ms is nonempty and p is the square of the last
   element of ms, k is > the last element of ms @ rev ls, and there
   are no primes that are > the last element of ms @ rev ls but < k,
   then next(ms, p, ls, k) returns (ms', p', ls', k') where k' is the
   smallest prime that is >= k (and so the smallest prime that is >
   the elements of ms @ rev ls), ms' @ rev ls' are the first length
   ms' + length ls' primes, listed in ascending order, ms' is nonempty
   and p' is the square of the last element of ms', k' is > the last
   element of ms' @ rev ls', and there are no primes that are > the
   last element of ms' @ rev ls' but < k'

   tail recursive; distance between fourth argument and smallest
   prime >= fourth argument goes down *)

fun next(ms, p, ls, k) =
      let val (ms, p, ls) =
                if null ls orelse p >= k
                then (ms, p, ls)
                else (ms @ rev ls, hd ls * hd ls, nil)
      in if divisible(k, ms)
         then next(ms, p, ls, k + 1)
         else (ms, p, ls, k)
      end

(* val prs :
         int * int * int list * int * int list * int ->
         int list

   if 1 <= i <= n, ms @ rev ls are the first i primes, listed in
   ascending order, ms is nonempty and p is the square of the last
   element of ms, and k is the last element of ms @ rev ls, then
   prs(n, i, ms, p, ls, k) returns the first n primes, listed in
   ascending order

   tail recursive *)

fun prs(n, i, ms, p, ls, k) =
      if i = n
      then ms @ rev ls
      else let val (ms, p, ls, k) = next(ms, p, ls, k + 1)
           in prs(n, i + 1, ms, p, k :: ls, k) end

(* val primes : int -> int list

   if n >= 0, then primes n returns the first n prime numbers,
   listed in ascending order *)

fun primes n =
      if n = 0 then nil else prs(n, 1, [2], 4, [], 2)

end;