program GeneratePrimes(output, primes, table) ;
{ Tom Verhoeff, Eindhoven University of Technology }

const
  ND = 5 ; { generate primes with ND digits }

var
  primes, table: text ;

procedure Init ;
  begin
  writeln('Generate ', ND:1, '-digit primes to file and tabulate counts') ;
  assign(primes, 'primes.dat') ; rewrite(primes) ;
  assign(table, 'counts.dat') ; rewrite(table)
  end { Init } ;

procedure Fini ;
  begin
  close(primes) ; close(table)
  end { Fini } ;

var 
  count: array[2..ND*9, 1..9] of integer ;
    { count[s, f] = # ND-digit primes with digit sum s and first digit f }
 
function IsPrime(n: integer): boolean ; { pre: n >= 4 }
  var i, j, h: integer ; 
  begin  
  h := trunc(sqrt(n)) + 1 ;
  i := 2 ; j := h ; { bounded linear search for smallest divisor }
  while i <> j do 
    if n mod i = 0 then j := i { i divides n, n is not prime } 
    else inc(i) ; 
  IsPrime := (i = h) 
  end { IsPrime } ;

procedure Generate ; 
  var first, i, s, f: integer ; 
  begin 
  first := 1 ;
  for i := 1 to pred(ND) do first := 10*first ; { first = 10^(ND-1) }
  for i := first to pred(10*first) do
    if IsPrime(i) then begin 
      f := i ; s := i mod 10 ;
      while f >= 10 do begin f := f div 10 ; s := s + f mod 10 end ;   
      { s = digit sum ; f = first digit }    
      writeln(primes, i:ND, ' ', s:2) ;
      inc(count[s, f])
      end { if }  
  end { Generate } ;

procedure Tabulate ;
  var s, f, t, m, mt: integer ; ct: array [1..9] of integer ;
  begin
  writeln('Writing table') ;
  writeln(table, 'Counts for ', ND:1, '-digit primes') ;
  writeln(table, 'Vertical: digit sum (multiples of 3 omitted)') ;
  writeln(table, 'Horizontal: first digit') ;
  writeln(table) ;
  writeln(table, '    |   1   2   3   4   5   6   7   8   9 | total') ;
  writeln(table, '----+-------------------------------------+------') ;
  mt := 0 ; { mt = max count }
  for f := 1 to 9 do ct[f] := 0 ; { ct[f] = total count for column f }
  for s := 2 to ND*9 do if s mod 3 <> 0 then begin
      write(table, s:3, ' |') ;
      t := 0 ;
      for f := 1 to 9 do begin
        if count[s, f] = 0 then write(table, ' ':4)
        else write(table, count[s, f]:4) ;
        inc(t, count[s, f]) ; inc(ct[f], count[s, f])
        end { for f } ;
      writeln(table, ' |', t:5) ;
      if t > mt then begin mt := t ; m := s end
      end { for s if } ;
  writeln(table, '----+-------------------------------------+------') ;
  write(table, 'tot |') ;
  t := 0 ;
  for f := 1 to 9 do begin write(table, ct[f]:4) ; inc(t, ct[f]) end ;
  writeln(table, ' |', t:5) ;
  writeln(table) ;
  writeln(table, 'Max count is ', mt:1, ' for digit sum ', m:1)
  end { Tabulate } ;

begin
Init ;
Generate ;
Tabulate ;
Fini
end.