Introduction to Programming in ATS:
<<< PreviousProgramming with Theorem-ProvingNext >>>

Example: Verified Fast Exponentiation

Given an integer x, pow(x, n), the nth power of x, can be defined inductively as follows: 

pow (x, 0) = 1
pow (x, n) = x * pow (x, n-1) (for all n > 0)

A direct implementation of this definition is given as follows: 

fun ipow {n:nat} .<n>.
  (x: int, n: int n): int = if n > 0 then x * ipow (x, n-1) else 1
// end of [ipow]

which is of time-complexity O(n) (assuming multiplication is O(1)). A more efficient implmentation can be given as follows: 

fun ifastpow {n:nat} .<n>.
  (x: int, n: int n): int =
  if n > 0 then let
    val n2 = n/2; i = n-(2*n2)
  in
    if i > 0 then pow (x*x, n2) else x * pow (x*x, n2)
  end else 1
// end of [ifastpow]

which makes use of the property that pow(x, n) equals pow(x*x, n/2) if n is even or x * pow(x*x, n/2) if n is odd. This is referred to as fast exponentiation. Note that ifastpow is of time-complexity O(log(n)).

Clearly, what is done above is not restricted to exponentiation on integers. As long as the underlying multiplication is associative, fast exponentiation can be employed to compute powers of any given element. In particular, powers of square matrices can be computed in this way. I now present as follows a verified generic implementation of fast exponentiation.

Handling generic data properly in a verified implementation often requires some finesse with the type system of ATS. Let us first introduce an abstract type constructor E as follows: 

sortdef elt = int // [elt] is just an alias for [int]
abst@ype E (a:t@ype, x:elt) = a // [x] is an imaginary stamp

This is often referred to as stamping. For each type T and stamp x, E(T, x) is just T as far as data representation is concerned. The stamps are imaginary and they are solely used for the purpose of specification. We next introduce an abstract prop-type MUL and a function template mul_elt_elt: 

absprop MUL (elt, elt, elt) // abstract mul relation

fun{a:t@ype}
mul_elt_elt {x,y:elt}
  (x: E (a, x), y: E (a, y)): [xy:elt] (MUL (x, y, xy) | E (a, xy))
// end of [mul_elt_elt]

Please do not confuse MUL with the one of the same name that is declared in prelude/SATS/arith.sats. To state that the encoded multiplication is associative, we can introduce the following proof function: 

praxi mul_assoc
  {x,y,z:elt} {xy,yz:elt} {xy_z,x_yz:elt} (
  pf1: MUL (x, y, xy), pf2: MUL (xy, z, xy_z)
, pf3: MUL (y, z, yz), pf4: MUL (x, yz, x_yz)
) : [xy_z==x_yz] void

The keyword praxi indicates that mul_assoc is treated as a form of axiom, which is not expected to be implemented.

The abstract power function can be readily specified in terms of the abstract prop-type MUL: 

dataprop POW (
  elt(*base*), int(*exp*), elt(*res*)
) = // res = base^exp
  | {x:elt} POWbas (x, 0, 1(*unit*))
  | {x:elt} {n:nat} {p,p1:elt}
    POWind (x, n+1, p1) of (POW (x, n, p), MUL (x, p, p1))
// end of [POW]

As can be expected, generic fast exponentiation is given the following interface: 

fun{a:t@ype}
fastpow_elt_int {x:elt} {n:nat}
  (x: E (a, x), n: int n): [p:elt] (POW (x, n, p) | E (a, p))
// end of [fastpow_elt_int]

With the preparation done above, a straightforward implementation of fastpow_elt_int can now be presented as follows: 

implement{a}
fastpow_elt_int (x, n) = let
//
// lemma: (x*x)^n = x^(2n)
//
extern prfun
lemma {x:elt} {xx:elt} {n:nat} {y:elt}
  (pfxx: MUL (x, x, xx), pfpow: POW (xx, n, y)): POW (x, 2*n, y)
//
in
  if n > 0 then let
    val n2 = n / 2; val i = n - (n2+n2) // i = 0 or 1
    val (pfxx | xx) = mul_elt_elt (x, x) // xx = x*x
    val (pfpow2 | res) = fastpow_elt_int<a> (xx, n2) // xx^n2 = res
    prval pfpow = lemma (pfxx, pfpow2) // pfpow: x^(2*n2) = res
  in
    if i > 0 then let
      val (pfmul | xres) = mul_elt_elt<a> (x, res) // xres = x*res
    in
      (POWind (pfpow, pfmul) | xres)
    end else (pfpow | res)
  end else let
    val res = mulunit<a> () in (POWbas () | res) // res = 1
  end (* end of [if] *)
end // end of [fastpow_elt_int]

Note that this implementation of fastpow_elt_int is not tail-recursive. The function template mulunit, which is called to produce a unit for the underlying multiplication, is assigned the following interface: 

fun{a:t@ype} mulunit (): E (a, 1(*stamp*))

The proof function lemma simply establishes that pow(x, 2*n)=pow(x*x, n) for each natural number n. I have made an implementation of lemma available on-line but I suggest that the interested reader give it a try first before taking a look. Note that the following axioms are needed to implement lemma: 

praxi mul_istot // [MUL] is total
  {x,y:elt} (): [xy:elt] MUL (x, y, xy)
praxi mul_isfun {x,y:elt} {z1,z2:elt} // MUL is functional
  (pf1: MUL (x, y, z1), pf2: MUL (x, y, z2)): [z1==z2] void

Another interesting (and possibly a bit challenging) exercise is to implement fastpow_elt_int in a tail-recursive fashion.

Please find on-line the two files fastexp.sats and fastexp.dats that contain the entirety of the above presented code.

Now we have implemented fastpow_elt_int. How can we use it? Please find on-line an example in which fastpow_elt_int is called to implement fast exponentiation on a 2-by-2 matrix so that Fibonacci numbers can be computed in a highly efficient manner.

<<< PreviousHomeNext >>>
Example: Another Verified Factorial ImplementationUpSummary