Introduction to Programming in ATS:
Prev	Chapter 12. Programming with Theorem-Proving	Next

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 = half(n)
    val i2 = n-(2*n2)
  in
    if i2 > 0 then ifastpow (x*x, n2) else x * ifastpow (x*x, n2)
  end else 1 // end of [if]
// 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 ELT as follows:

sortdef elt = int // [elt] is just an alias for [int]
abst@ype ELT(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, ELT(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. Let us 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: ELT(a, x), y: ELT(a, y)): [xy:elt] (MUL(x, y, xy) | ELT(a, xy))
// end of [mul_elt_elt]
//

Please do not confuse MUL with the one of the same name that is declared in arith_prf.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}
(
  MUL(x, y, xy), MUL(xy, z, xy_z), MUL(y, z, yz), MUL(x, yz, x_yz)
) : [xy_z==x_yz] void // end of [mul_assoc]

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: ELT(a, x), n: int n): [p:elt] (POW(x, n, p) | ELT(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}(*tmp*)
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)
//
overload * with mul_elt_elt // [*] loaded with mul_elt_elt
//
in
//
if
n = 0
then let
  val res = mulunit<a> () in (POWbas () | res) // res = 1
end // end of [then]
else let
  val n2 = half n
  val (pfxx | 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 n=2*n2
    then (pfpow | res)
    else let
      val (pfmul | xres) = x * res in (POWind(pfpow, pfmul) | xres)
    end // end of [else]
end // end of [else]
//
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 (): ELT(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 to implement lemma before taking a look at the given implementation. Note that the following axioms are needed to implement lemma:

//
praxi
mul_istot // MUL is total
  {x,y:elt} ((*void*)): [xy:elt] MUL(x, y, xy)
//
praxi
mul_isfun // MUL is functional
  {x,y:elt}{z1,z2:elt}(MUL(x, y, z1), 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 it be used? 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 the Fibonacci numbers can be computed in a highly efficient manner.

Prev	Home	Next
Example: Another Verified Factorial	Up	Programming with Views and Viewtypes