Balanced binary search trees

The type of "association tables" (binary search trees).

type (α, β) t =
| Empty
| Node of (α , β) t * α * β * (α, β) t * int

There are two cases : a tree that is empty or, a node consisting of a left sub-tree, a key, the value associated with that key, a right sub-tree and, an integer representing the "height" of the tree (the number of nodes to traverse before reaching the most distant leaf).

The binary search tree invariant will be made to apply in that for any non empty tree $n$, every node in the left sub-tree is ordered less than $n$ and every node in the right sub-tree of $n$ is ordered greater than $n$ (in this program, ordering of keys is performed using the Pervasives.compare function).

This function, height, given a tree, extracts its height.

let height : (α, β) t -> int = function
  | Empty -> 0
  | Node (_, _, _, _, h) -> h

The value empty, is a constant, the empty tree.

let empty : (α, β) t = Empty

create l x d r creates a new non-empty tree with left sub-tree l, right sub-tree r and the binding of key x to the data d. The height of the tree created is computed from the heights of the two sub-trees.

let create (l : (α, β) t) (x : α) (d : β) (r : (α, β) t) : (α, β) t =
  let hl = height l and hr = height r in
  Node (l, x, d, r, (max hl hr) + 1)

This next function, balance is where all the action is at. Like the preceding function create, it is a factory function for interior nodes and so takes the same argument list as create. It has an additional duty though in that the tree that it produces takes balancing into consideration.

let balance (l : (α, β) t) (x : α) (d : β) (r : (α, β) t) : (α, β) t =
  let hl = height l and hr = height r in
  if hl > hr + 1 then
    match l with

In this branch of the program, it has determined that production of a node with the given left and right sub-trees (denoted $l$ and $r$ respectively) would be unbalanced because $h(l) > hr(1) + 1$ (where $h$ denotes the height function).

There are two possible reasons to account for this. They are considered in turn.

    (*Case 1*)
    | Node (ll, lv, ld, lr, _) when height ll >= height lr ->
      create ll lv ld (create lr x d r)

So here, we find that $h(l) > h(r) + 1$, because of the height of the left sub-tree of $l$.

    (*Case 2*)
    | Node (ll, lv, ld, Node (lrl, lrv, lrd, lrr, _), _) ->
      create (create ll lv ld lrl) lrv lrd (create lrr x d r)

In this case, $h(l) > h(r) + 1$ because of the height of the right sub-tree of $l$.

    | _ -> assert false

We assert false for all other patterns as we aim to admit by construction no further possibilities.

We now consider the case $h(r) > h(l) + 1$, that is the right sub-tree being "too long".

  else if hr > hl + 1 then
    match r with

There are two possible reasons.

    (*Case 3*)
    | Node (rl, rv, rd, rr, _) when height rr >= height rl ->
      create (create l x d rl) rv rd rr

Here $h(r) > h(l) + 1$ because of the right sub-tree of $r$.

    (*Case 4*)
    | Node (Node (rll, rlv, rld, rlr, _), rv, rd, rr, _) ->
      create (create l x d rll) rlv rld (create rlr rv rd rr)

Lastly, $h(r) > h(l) + 1$ because of the left sub-tree of $r$.

    | _ -> assert false

Again, all other patterns are (if we write this program correctly according to our intentions,) impossible and so, assert false as there are no further possibilities.

In the last case, neither $h(l) > h(r) + 1$ or $h(r) > h(l) + 1$ so no rotation is required.

  else
    create l x d r

add x data t computes a new tree from t containing a binding of x to data. It resembles standard insertion into a binary search tree except that it propagates rotations through the tree to maintain balance after the insertion.

let rec add (x : α) (data : β) : (α, β) t -> (α, β) t = function
    | Empty -> Node (Empty, x, data, Empty, 1)
    | Node (l, v, d, r, h) ->
      let c = compare x v in
      if c = 0 then
        Node (l, x, data, r, h)
      else if c < 0 then
        balance (add x data l) v d r
      else 
        balance l v d (add x data r)

To implement removal of nodes from a tree, we'll find ourselves needing a function to "merge" two binary searchtrees $l$ and $r$ say where we can assume that all the elements of $l$ are ordered before the elements of $r$.

let rec merge (l : (α, β) t) (r : (α, β) t) : (α, β) t = 
  match (l, r) with
  | Empty, t -> t
  | t, Empty -> t
  | Node (l1, v1, d1, r1, h1), Node (l2, v2, d2, r2, h2) ->
    balance l1 v1 d1 (balance (merge r1 l2) v2 d2 r2)

Again, rotations are propagated through the tree to ensure the result of the merge results in a balanced tree.

With merge available, implementing remove becomes tractable.

let remove (id : α) (t : (α, β) t) : (α, β) t = 
  let rec remove_rec = function
    | Empty -> Empty
    | Node (l, k, d, r, _) ->
      let c = compare id k in
      if c = 0 then merge l r else
        if c < 0 then balance (remove_rec l) k d r
        else balance l k d (remove_rec r) in
  remove_rec t

The remaining algorithms below are "stock" algorithms for binary search trees with no particular consideration of balancing necessary and so we won't dwell on them here.

let rec find (x : α) : (α, β) t -> β = function
  | Empty ->  raise Not_found
  | Node (l, v, d, r, _) ->
    let c = compare x v in
    if c = 0 then d
    else find x (if c < 0 then l else r)

let rec mem (x : α) : (α, β) t -> bool = function
  | Empty -> false
  | Node (l, v, d, r, _) ->
    let c = compare x v in
    c = 0 || mem x (if c < 0 then l else r)
    
let rec iter (f : α -> β -> unit) : (α, β) t -> unit = function
  | Empty -> ()
  | Node (l, v, d, r, _) ->
    iter f l; f v d; iter f r

let rec map (f : α -> β -> γ) : (α, β) t -> (α, γ) t = function
  | Empty -> Empty
  | Node (l, k, d, r, h) -> 
    Node (map f l, k, f k d, map f r, h)

let rec fold (f : α -> β -> γ -> γ) (m : (α, β) t) (acc : γ) : γ =
  match m with
  | Empty -> acc
  | Node (l, k, d, r, _) -> fold f r (f k d (fold f l acc))

open Format

let print 
    (print_key : formatter -> α -> unit)
    (print_data : formatter -> β -> unit)
    (ppf : formatter)
    (tbl : (α, β) t) : unit =
  let print_tbl ppf tbl =
    iter (fun k d -> 
           fprintf ppf "@[<2>%a ->@ %a;@]@ " print_key k print_data d)
      tbl in
  fprintf ppf "@[[[%a]]@]" print_tbl tbl

The source code for this post can be found in the file 'ocaml/misc/tbl.ml' in the OCaml source distribution. More information on balanced binary search trees including similar but different implementation techniques and complexity analyses can be found in this Cornell lecture and this one.

Rotate

Rotate

This post is inspired by one of those classic "99 problems in Prolog".What we are looking for here are two functions that satisfy these signatures.

val rotate_left : int -> α list -> α list
val rotate_right : int -> α list -> α list 

rotate_left n rotates a list $n$ places to the left, rotate_right n rotates a list $n$ places to the right. Examples:

# rotate_left 3 ['a';'b';'c';'d';'e';'f';'g';'h'] ;;
- : char list = ['d'; 'e'; 'f'; 'g'; 'h'; 'a'; 'b'; 'c']

# rotate_left (-2) ['a';'b';'c';'d';'e';'f';'g';'h'] ;;
- : char list = ['g'; 'h'; 'a'; 'b'; 'c'; 'd'; 'e'; 'f']

Of course, rotate_left and rotate_right are inverse functions of each other so we expect, for any int $x$ and list $l$, rotate_right x @@ rotate_left x l $=$ rotate_left x @@ rotate_right x l $=$ l.

Well, there are a variety of solutions to this problem with differing degrees of verbosity, complexity and efficiency. My own attempt at a solution resulted in this.

let rec drop (k : int) (l : α list) : α list =
  match k, l with
  | i, _ when i <= 0 -> l
  | _, [] -> []
  | _, (_ :: xs) -> drop (k - 1) xs

let rec take (k : int) (l : α list) : α list =
  match k, l with
  | i, _ when i <= 0 -> []
  | _, [] -> []
  | _, (x :: xs)  -> x :: take (k - 1) xs

let split_at (n : int) (l : α list) : α list * α list = 
  (take n l), (drop n l)

let rec rotate_left (n : int) (l : α list) : α list =
  match n with
  | _ when n = 0 -> l
  | _ when n < 0 ->  rotate_right (-n) l
  | _ -> 
    let m : int = List.length l in
    let k : int = n mod m in
    let (l : α list), (r : α list) = split_at k l in 
    r @ l

and rotate_right (n : int) (l : α list) : α list =
  match n with
  | _ when n = 0 -> l
  | _ when n < 0 ->  rotate_left (-n) l
  | _ -> 
    let m : int = List.length l in
    let k : int = m - n mod m in
    let (l : α list), (r : α list) = split_at k l in 
    r @ l

So far so good, but then I was shown the following solution in Haskell.

rotateLeft n xs 
  | n >= 0     = take (length xs) $ drop n $ concat $ repeat xs
  | otherwise  = rotateLeft (length xs + n) xs

rotateRight n = rotateLeft (-n)

I found that pretty nifty! See, in the function rotateLeft, repeat xs creates an infinite list of lists, (each a copy of xs), "joins" that infinite list of lists into one infinite list, then the first $n$ elements are dropped from that the list and we take the next length xs which gets us the original list rotated left $n$ places.

I felt compelled to attempt to emulate the program above in OCaml.

The phrasing "works" in Haskell due to the feature of lazy evaluation. OCaml on the other hand is eagerly evaluated. Lazy evaluation is possible in OCaml however, you just need to be explicit about it. Here's a type for "lazy lists" aka "streams".

type α stream =  Nil | Cons of α * α stream Lazy.t

A value of type α Lazy.t is a deferred computation, called a suspension that has the result type α. The syntax lazy$(expr)$ makes a suspension of $expr$, without yet evaluating $expr$. "Forcing" the suspension (using Lazy.force) evaluates $expr$ and returns its result.

Next up, functions to get the head and tail of a stream.

let hd = function | Nil -> failwith "hd" | Cons (h, _) -> h
let tl = function | Nil -> failwith "tl" | Cons (_, t) -> Lazy.force t

Also useful, a function to lift an α list to an α stream.

let from_list (l : α list) : α stream =
  List.fold_right (fun x s -> Cons (x, lazy s)) l Nil

Those are the basic building blocks. Now we turn attention to implementing repeat x to create infinite lists of the repeated value $x$.

let rec repeat (x : α) : α stream = Cons (x, lazy (repeat x))

Now to implement concat (I prefer to call this function by its alternative name flatten).

The characteristic operation of flatten is the joining together of two lists. For eager lists, we can write a function join that appends two lists like this.

let rec join l m =
  match l with
  | [] -> m
  | h :: t -> h :: (join t m)

This generalizes naturally to streams.

let rec join (l : α stream) (m : α stream) =
  match l with
  | Nil -> m
  | Cons (h, t) -> Cons (h, lazy (join (Lazy.force t) m))

For eager lists, we can write flatten in terms of join.

let rec flatten : α list list -> α list = function
   | [] -> []
   | (h :: tl) -> join h (flatten tl)

Emboldened by our earlier success we might try to generalize it to streams like this.

let rec flatten (l : α stream stream) : α stream =
   match l with
   | Nil -> lazy Nil
   | Cons (l, r) ->  join l (flatten (Lazy.force r))

Sadly, no. This definition is going to result in stack overflow. There is an alternative phrasing of flatten we might try.

let rec flatten = function
  | [] -> []
  | [] :: t -> flatten t
  | (x :: xs) :: t -> x :: (flatten (xs :: t))

Happy to say, this one generalizes and gets around the eager evaluation problem that causes the unbounded recursion.

let rec flatten : α stream stream -> α stream = function
  | Nil -> Nil
  | Cons (Nil, t) -> flatten (Lazy.force t)
  | Cons (Cons (x, xs), t) ->
      Cons (x, lazy (flatten (Cons (Lazy.force xs, t))))

take and drop are straight forward generalizations of their eager counterparts.

let rec drop (n : int) (lst : α stream ) : α stream = 
  match (n, lst) with
  | (n, _) when n < 0 -> invalid_arg "negative index in drop"
  | (n, xs) when n = 0 -> xs
  | (_, Nil) -> Nil
  | (n, Cons (_, t)) -> drop (n - 1) (Lazy.force t)

let rec take (n : int) (lst : α stream) : α list = 
  match (n, lst) with
  | (n, _) when n < 0 -> invalid_arg "negative index in take"
  | (n, _) when n = 0 -> []
  | (_, Nil) -> []
  | (n, Cons (h, t)) -> h :: (take (n - 1) (Lazy.force t))

Which brings us to the lazy version of rotate expressed in about the same number of lines of code!

let rec rotate_left (k : int) (l : α list) : α list =
  let n = List.length l in
  if k >= 0 then
    l |> from_list |> repeat |> flatten |> drop k |> take n
  else rotate_left (n + k) l

let rotate_right (n : int) : α list -> α list = rotate_left (-n)

Heap sort

Given the existence of a priority queue data structure, the heap sort algorithm is trivially implemented by loading the unsorted sequence into a queue then successively pulling of the minimum element from the queue until the queue is exhausted.

There are many ways to implement a priority queue and so we seek an expression for a function for heap sorting that is polymorphic over those choices.

To begin, a module type for priority queues.

(**Priority queues over ordered types*)
module type PRIORITY_QUEUE = sig

    (**Output signature of the functor [Make]*)
    module type S = sig
        exception Empty

        type element (*Abstract type of elements of the queue*)
        type t (*Abstract type of a queue*)

        val empty : t (*The empty queue*)
        val is_empty : t -> bool (*Check if queue is empty*)
        val insert : t -> element -> t (*Insert item into queue*)
        val delete_min : t -> t (*Delete the minimum element*)
        val find_min : t -> element (*Return the minimum element*)
        val of_list : element list -> t
      end

    (**Input signature of the functor [Make]*)
    module type Ordered_type = sig
        type t
        val compare : t -> t -> int
      end

   (**Functor building an implementation of the priority queue structure
      given a totally ordered type*)
    module Make : 
       functor (Ord : Ordered_type) -> S with type element = Ord.t
  end

An implementation of this signature using "leftist heaps" is described for the interested in this Caltech lab but such details are omitted here.

module Priority_queue : PRIORITY_QUEUE = struct 
  module type S = sig .. end
  module type Ordered_type = sig .. end
  module Make (Elt : Ordered_type) : (S with type element = Elt.t) = struct .. end
end

What I really want to show you is this. We start with the following module type abbreviation.

type 'a queue_impl = (module Priority_queue.S with type element = 'a)

Then, the heap_sort function can be written such that it takes a module as a first class value and uses a locally abstract type to connect it with the element type of the list to be sorted.

let heap_sort (type a) (queue : a queue_impl) (l : a list) : a list =
  let module Queue = 
     (val queue : Priority_queue.S with type element = a) in
  let rec loop acc h =
    if Queue.is_empty h then acc
    else
      let p = Queue.find_min h in
      loop (p :: acc) (Queue.delete_min h) in
  List.rev (loop [] (Queue.of_list l))

There we have it. The objective has been achieved : we have written a heap sorting function that is polymorphic in the implementation of the priority queue with which it is implemented.

Usage (testing) proceeds as in this example.

(*Prepare an [Priority_queue.Ordered_type] module to pass as argument
  to [Priority_queue.Make]*)
module Int : Priority_queue.Ordered_type with type t = int = struct 
  type t = int let compare = Pervasives.compare 
end

(*Make a priority queue module*)
module Int_prioqueue : (Priority_queue.S with type element = int) = Priority_queue.Make (Int)

(*Make a first class value of the module by packing it*)
let queue = (module Int_prioqueue : Priority_queue.S with type element = int)

(*Now, pass the module to [heap_sort]*)
let sorted = heap_sort queue [-1; -2; 2] (*Produces the list [-2; -1; 2]*)

Addendum :

These ideas can be pushed a little further yielding a simpler syntax for the parametric heapsort algorithm.

(*Type abbreviations*)

type 'a order_impl = (module Priority_queue.Ordered_type with type t = 'a)
type 'a queue_impl = (module Priority_queue.S with type element = 'a)

(*Module factory functions*)

let mk_ord : 'a. unit -> 'a order_impl =
  fun (type s) () ->
    (module 
     struct 
       type t = s 
       let compare = Pervasives.compare 
     end : Priority_queue.Ordered_type with type t = s
    )

let mk_queue : 'a. unit -> 'a queue_impl =
  fun (type s) ord ->
    let module Ord = 
      (val mk_ord () : Priority_queue.Ordered_type with type t = s) in
    (module Priority_queue.Make (Ord) :
                               Priority_queue.S with type element = s)

For example, now we can write

# heap_sort (mk_queue ()) [-3; 1; 5] ;;
  - : int list = [-3; 1; 5]

Recursive list

Recursive list

Way, way back we looked at a function to estimate the value of the sin function utilizing its representation as a Taylor polynomial about x = 0. When one looks closely at that Taylor polynomial, one can observe a pattern in the coefficients involving the infinite sequence 0, 1, 0, -1, 0, 1, 0, -1, .... My buddy Marcelo in his formulation of that Taylor polynomial wrote this function that for a given n produces the first n values of this sequence as a list:

let rec sin_coeff n =
  match n with
  | 0 -> []
  | 1 -> [0.]
  | 2 -> [0.; 1.]
  | 3 -> [0.; 1.; 0. ]
  | _  -> 0. :: 1. :: 0. :: -1. :: sin_coeff (n-4)

He mused to me on the day that "surely there is a more elegant way to produce this list?". Oh, by golly Marcelo there certainly is!

let rec sin_coeff = 0 :: 1 :: 0 :: -1 :: sin_coeff

This little construction builds the list value sin_coeff recursively (that is, in terms of itself). Now, of course it's not impossible to define structures like this in C or C++ but good luck matching the brevity and elegance afforded to us in OCaml for this sort of thing! We still need a little function of course that can pull off this endless list a list containing just the first n elements. The ubiquitous take function will do for this purpose and I show it here just for completeness.

let rec list_take k l =
  if k <= 0 then []
  else
    match l with
    | [] -> []
    | (x :: xs)  -> x :: (list_take (k-1) xs)

Terms With Variables (C++)

Terms with Variables (C++)

In this earlier post I showed a nifty OCaml type for modeling terms with variables for problems involving substitutions. I got interested in what it would take to implement the type in C++(03) (doing 'sum' types in C++ elegantly is a perennial problem). It ain't nowhere near as succinct but we got there nonetheless.

#include <list>

#include <boost/variant.hpp>

/*
  type ('a, 'b) term =
    | Term of 'a * ('a, 'b) term list
    | Var of 'b
*/

template <class A, class B> struct term;
template <class B> struct var;

template <class A, class B>
struct make_tree
{
  typedef boost::variant <
    boost::recursive_wrapper<term <A, B> >, var<B> >   type;
};

template <class A, class B>
struct term
{
  typedef typename make_tree <A, B>::type tree;
  A a;
  std::list <tree> children;
  term (A a
    , std::list<tree> const& children)
    :       a (a), children (children)
  {}
};

template <class A, class B> 
inline term <A, B> make_term (
  A a, std::list<typename make_tree<A, B>::type> const& c, B const*)
{
  return term<A, B> (a, c);
}

template <class B> 
struct var
{
  B tag;
  var (B tag) : tag (tag) {}
};

template <class B>
inline var<B> make_var (B tag) { return var<B> (tag); }

For example, this little program builds the term represented by the concrete syntax "a(b(), c)".

int main ()
{
  typedef make_tree<std::string, std::string>::type tree;
  typedef std::list<tree> term_list;
  std::string const* tag_str=(std::string const*)0L;

  // a(b(), c)
  term_list l;
  l.push_back (make_term(std::string("b"), term_list (), tag_str));
  l.push_back (make_var(std::string("c")));
  tree t = make_term(std::string("a"), l, tag_str);
  
  return 0;
}

Terms with Variables

Cosineau and Mauny present a type suitable for modeling terms with variables convenient for applications involving variables and substitutions. It is defined like this:

(*term_types.ml*)
type (α, β) term = 
| Term of α * (α, β) term list
| Var of β

Analyzing syntax

To proceed, let us restrict our attention to the type (string, string) term. The concrete syntax we use for such expressions is parenthesized.

• a (t₁, ..., t_n) is the syntax of Term (a, [t1; ...; tn])

We disambiguate variables from constants by representing constants as functions applied to an empty list of arguments

• a (b (), c) is the syntax of Term ("a", [Term ("b", []); Var "c"])

We start towards a function to analyze syntax. First the lexer.

(*lexer.mll*)
{
  open Parser  (*The type token is defined in parser.mli.*)
}

let alpha=['a'-'z' 'A'-'Z']

rule token = parse
    [' ' '\t' '\n']         (* skip blanks *) { token lexbuf }
  | '('                                             { LPAREN }
  | ')'                                             { RPAREN }
  | ','                                              { COMMA }
  | ((alpha)(alpha)*) as s                         { STRING s}
  | eof                                                { EOI }

Now the parser.

/*parser.mly*/

%token  STRING
%token LPAREN RPAREN COMMA EOI

%start main
%type <(string, string) Term_types.term> main
%%
main:
 | exp EOI { $1 }
 ;
exp:
 | tag exp_more           
     { 
       match $2 with
       | None -> Term_types.Var $1
       | Some l -> Term_types.Term ($1, l)
     }
 ;
exp_more:
 | /*nothing*/ { None }
 | LPAREN exp_list RPAREN { Some $2 }
;
exp_list:
 | /*nothing*/ { [] }
 | exp exp_list_more { $1::$2 }
 ;
exp_list_more:
 | /*nothing*/ {[]}
 | COMMA exp exp_list_more { $2::$3 }
 ;
tag:
 | STRING { $1 }
 ;

With these pieces in place we can write the analyzer term_of_string like so

(*term.ml*)

let term_of_string s =
  let parse_buf lexbuf =
    try Parser.main Lexer.token lexbuf
    with 
    | Parsing.Parse_error ->
      begin
        let curr = lexbuf.Lexing.lex_curr_p in
        let line = curr.Lexing.pos_lnum in
        let cnum = curr.Lexing.pos_cnum - curr.Lexing.pos_bol in
        let tok = Lexing.lexeme lexbuf in
        raise 
          (Failure
             (Printf.sprintf 
                "file \"\", line %d, character %d\n\Error : Syntax error \"%s\"" line cnum tok))
      end
  in parse_buf (Lexing.from_string s)

Reverse parsing is done with string_of_term:

let rec string_of_term (t:(string, string)term) : string = 
  let string_of_list (f:α -> string) (l:α list) : string = 
    "(" ^ String.concat "," (List.map f l) ^ ")"  in
  match t with
  | Term (s, tl) -> Printf.sprintf "%s" (s^(string_of_list string_of_term tl))
  | Var s -> s

A simple function for printing a term then goes like this:

let print_term (t:(string, string) term) : unit = Printf.printf "%s" (string_of_term t)

For example, this program let tt = term_of_string "a(b(), c)" in print_term tt prints "a(b(), c)" recovering the input as we'd hope.

Traversal

The function term_trav is a workhorse for doing computations over terms. It generalizes all functions defined by recursion on the structure of a term.

(* [term_trav] is a generalization of all functions defined by
   recusion on the structure of an [(α, β) term] - a homomorphism.

   The definition of [term_trav] below gives the following signature.

   {[
   val term_trav : (α * β -> γ) ->
   (γ -> β -> β) -> β -> (δ -> γ) -> (α, δ) term -> γ
   ]}

   Proof:

   Assume [term_trav] to operate on arguments of type [(α, δ)
   term]. It "returns" either via [f] or [v], and so [f], [v] and
   [term_trav] must share the same return type, [γ] say.

   Since the intermediate list results from recursive calls to
   [term_trav] then it must have type [γ list]. This means [g] must
   be of type [γ -> β -> β] for some type [β] which fixes [x] to
   [β] whilst completing [f] as [α * β -> γ].

   Lastly [v] takes arguments of type [δ] meaning it types to [δ ->
   γ].

*)
let rec term_trav f g x v = function
  | Term (a, tl) -> 
    let l = List.map (term_trav f g x v) tl in
    let res = List.fold_right g l x in
    f (a, res)
  | Var b -> v b

For a trivial example, we can express a function to compute the "size" of a term this way:

let term_size (t:(α, β) term) : int = 
  (term_trav (fun (_, s) -> s + 1) (fun x y -> x + y) 0 (fun _ -> 1))
    t

I appreciate there's an efficiency gain to be had by eliminating the intermediate list in term_trav but the resulting function definition is more challenging to reason about (well, at least for me).