Monday, January 20, 2014

Computing a US holiday calendar

Computing a US holiday calendar

Assume the existence of a type date and a function nth_kday_of_month with semantics described by Boost.Datetime documentation. We also assume sums type day_of_week=|Sunday|Monday|... and type ord=|First|second|....

These are utility functions I find myself writing time and time again.

let range l u =
    let rec range_aux acc x y =
      if x = y then acc
      else range_aux (acc@[x]) (x+1) y
    in
    range_aux [] l u

let string_of_list f l = "[" ^ String.concat ";" (List.map f l) ^ "]"

The approach taken is to use date generators to convert abstract specification into a concrete set of dates.

let test () =
  let us_holiday_calendar (s, e) = 
    let generators = [
      nth_kday_of_month First Monday 09;  (*US labor day*)
      nth_kday_of_month Third Monday 01;  (*Martin Luther King day*)
      nth_kday_of_month Second Tuesday 02;  (*President's day*)
      nth_kday_of_month Fourth Thursday 11;  (*Thanksgiving*) 
    ] in
    let years = range s e in
    let apply_gen acc g = (acc@List.map g years) in
    let holidays = List.fold_left apply_gen [] generators in  
    List.sort compare holidays in
  assert_equal "[2014-01-20;2014-02-11;2014-09-01;2014-11-27;2015-01-19;2015-02-10;2015-09-07;2015-11-26;2016-01-18;2016-02-09;2016-09-05;2016-11-24;2017-01-16;2017-02-14;2017-09-04;2017-11-23;2018-01-15;2018-02-13;2018-09-03;2018-11-22;2019-01-21;2019-02-12;2019-09-02;2019-11-28;2020-01-20;2020-02-11;2020-09-07;2020-11-26;2021-01-18;2021-02-09;2021-09-06;2021-11-25;2022-01-17;2022-02-08;2022-09-05;2022-11-24;2023-01-16;2023-02-14;2023-09-04;2023-11-23]"
  @@ (us_holiday_calendar (2014, 2024) |> string_of_list string_of_date)

This code calculates a 10y calendar of US holidays between 2014 and 2024. If the operators @@ and |> are unfamiliar you can read about them here.

Saturday, December 7, 2013

Pipelining with the |> operator in OCaml

Pipelining with |> in OCaml

Reading through Real World OCaml (which I'm really enjoying and can heartily recommend by the way), I came across this clever definition for pipe notation.
let ( |> ) x f = f x
This little operator allows for (maybe misusing a term I think may have been coined by my friend MS,) "reverse application". To get a sense of the impact on readability that it may have, I applied it to the program of this earlier blog entry. At the heart of the program is the following little fragment.
match (filesin !root) with
| Some names ->
  let n=String.length !suffix in
  let pred e =
    let i = (String.length e - n) in
    (i >= 0) && (String.sub e i n) = !suffix
  in process_files (List.filter pred (Array.to_list names) )
Now, rewritten to use pipelining by virtue of |> we get this.
match (filesin !root) with
| Some names ->
  let n=String.length !suffix in
  let pred e =
    let i = (String.length e - n) in
    (i >= 0) && (String.sub e i n) = !suffix
  in names |> Array.to_list |> List.filter pred |> process_files 
 
I might be easily amused but... wow!

Coincidentally, reading from the OCaml PRO blog I learned of this related operator (also shown independently to me by my buddy -- thanks Stefano!).
let ( @@ ) f x = f x
This one can be used to good effect reducing the syntactic noise of lots of parentheses as per their example.
List.iter print_int @@ List.map (fun x -> x + 1 ) @@ [1; 2; 3]
(I assume it's obvious what this does and also, I added spaces in accordance with this advice).

Of course, if I'd had my wits about me, I'd have realized that I've known this operator all this time as $ from the Felix programming language... No excuse man, I mean, it's right there in the introductory tutorial!
println$ "Hello " + Env::getenv "USER";
Sorry 'bout being so slow to catch on Felix!

Now, I must admit I prefer $ to @@ but careful examination of the syntax reference seems to me that it won't be appropriate for OCaml (in that $ is left associative). As a final note, |> and @@ are "builtin" operators in OCaml compilers today.

Sunday, November 3, 2013

Matrix multiplication (functional approach in Python)

Matrix multiplication

The last entry looked at Gaussian elimination modeling matrices and vectors as tuples and using some idioms of functional programming. The linear algebra theme is continued here with programs to compute matrix/vector and matrix/matrix multiplications.

First some utilities to warm up on. Functions to project the rows or columns of a matrix.
def proj_r (A, i): return A[i]
def proj_c (A, j): return tuple (map (lambda t: (t[j],), A))
The function that scales a row vector by a constant is here again along with its symmetric partner for scaling a column vector.
def scale_r (u, lam): 
  return tuple (map (lambda ui: lam*ui, u))

def scale_c (c, k):
  def f (e):
    (xi,) = e
    return k*xi,
  return tuple (map (f, c))
For example, if c=((1,),(2,)), then scale_c (c, 2) is ((2,), (4,)).

The utility flatten reinterprets a sequence of column vectors as a matrix of row vectors.
def flatten (t):
  def row(i):
    def f (acc, y): 
        return acc + y[i]
    return functools.reduce (f, t, ())
  return tuple (map (row, range (0, len (t[0]))))
Here is an example. If t1=(1,),(2,) and t2=(3,),(4,) then s=(t1, t2) is a column vector sequence and flatten (s) is the matrix ((1, 3), (2, 4)).

Matrix vector multiplication comes next.
def mat_vec_mul (A, x):
  def f (i):
      (xi,) = x[i]
      return scale_c (proj_c (A, i), xi)
  M = tuple (map (f, range (0, len (x))))
  return tuple (map (lambda r : (sum (r), ), flatten (M)))
This computation utilizes the interpretation of the product as the sum of the column vectors of A weighted by the components of x.

Matrix multiplication follows immediately by seeing it as the flattening of a sequence of matrix/vector multiplications (that is of A on the individual columns of B).
def mat_mat_mul (A, B):
  f = lambda j : mat_vec_mul (A, proj_c (B, j))
  return flatten (tuple (map (f, range (0, len (B[0])))))
The transcript below shows the function at work in the Python top-level.
>>> A = ((2, 3), (4, 0) )
>>> B = ((1,  2, 0), (5, -1, 0) )
>>> print (str (mat_mat_mul (A, B)))
((17, 1, 0), (4, 8, 0))

Saturday, November 2, 2013

Gaussian elimination (functional approach in Python)

Gaussian elimination

The goal here is to implement simple Gaussian elimination in Python, in a functional style just using tuples.

We view (a, b, c) a row vector and interpret ((a,),(b,),(c,)) as a column vector.

These first two elementary operations (scaling a row by a scalar and subtracting one row from another) come easily.
def scale_row (u, lam): 
    return tuple (
      map (lambda ui: lam*ui, u))

def subtract_rows (u, v) :
  return tuple (
    map (lambda i : u[i] - v[i], range (0, len (u))))
These alone enable us to proceed directly to forward elimination of an augmented coefficient matrix.
def elimination_phase (t):
  def elimination_step (tp, k): 
    pr = tp[k]
    def g (e):
      (i, r) = e
      if i <= k : 
        return r
      else:
        lam = r[k]/pr[k]
        return subtract_rows (r, scale_row (pr, lam))
    return tuple (
             map (g, tuple (
                      zip (range (0, len (tp)), tp))))
  return functools.reduce (
              elimination_step, range (0, len (t)), t)
Vector dot products are a straight forward fold.
def dot (u, v):
  def f (acc, i):
    (a, (b,)) = (u[i], v[i])
    return acc + a*b
  return functools.reduce (f, range (0, len (u)), 0.0)
With this, the back substitution phase can be written like this.
def back_substitution_phase (t):
  n = len (t)
  def back_substitution_step (x, k):
    bk = (t[k])[n]
    akk = (t[k])[k]
    xk = bk/akk if k == n-1 \
      else (bk - dot ((t[k][k+1:])[0:n-k-1], x))/akk
    return ((xk,),) + x
  return functools.reduce (
    back_substitution_step, range (len (t)-1, -1, -1), ())
Can we test it? Yes we can!
#Solve Ax=b with 
#
#        6   -4     1
#  A =  -4    6    -1
#        1   -4     6
#
#and b = (-14, 36, 6)^T.

A = (
    ( 6.,  -4.,  1.,  -14.)
  , (-4.,   6., -4.,   36.)
  , ( 1.,  -4.,  6.,    6.))

print (str (back_substitution_phase (elimination_phase (A))))

#Should give (10, 22, 14)^T.

Term matching

Term matching

Extending one substitution with another can only work if the two substitutions are consistent with each other.
exception Term_match_exc

let extend_subst s1 s2 =
  let f acc (x, t)= 
    try 
      let u = List.assoc x acc in
      if t = u then acc
      else raise Term_match_exc
    with Not_found -> (x, t)::acc
in List.fold_left f s1 s2
A pattern match of a term u by a term t is a substitution that transforms t into u. If t is a variable s then the substitution (s, u) is an obvious pattern match for example. More generally though, if t is of the form f (t1, ..., tn), then
  • u must be of the form f (u1, ..., un);
  • for each i there must exist a match of ui with ti;
and all the matches must be mutually compatible for a solution to exist.
let term_match (t1, t2) =
  let rec term_match' subst = function
    | (Var v, t) -> extend_subst [v, t] subst
    | (t, Var v) -> raise Term_match_exc
    | (Term (f, l), Term (g, r)) ->
      if f = g then 
        List.fold_left term_match' subst (List.combine l r)
      else raise Term_match_exc
  in term_match' [] (t1, t2)
For example, a pattern match of h(y(z, k(x, s))) and h(y(l(m, n,o()), k(t(), u))) is {(x,t()),(z, l(m,n,o()), (s, u)}.

Saturday, October 19, 2013

Substitutions on terms with variables

Substitutions on terms with variables

This picks up on the earlier terms with variables blog entry. Remember the term_trav function? I'll recall it here but this time sharpen it up a bit with named parameters.
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

Application of substitutions

Substitutions are modeled as values of type ('b, ('a, 'b) term) list. Application of a substitution is about the finding of some of the variables of a term and replacing them with different terms to find the term's image under the substitution.
let apply_substitution subst t = 
  term_trav 
    ~f:(fun (f, l) -> Term (f, l)) 
    ~g:(fun x acc -> x::acc) 
    ~x:[] 
    ~v:(fun s -> try List.assoc s subst with _ -> Var s) 
    t

Of course, the definition of apply_substitution is equivalent to this one
let rec apply_substituion' subst = function
  | Term (f, ns) -> Term (f, List.map (apply_subst subst) ns)
  | Var x as v -> try List.assoc x subst with _ -> v
and I guess is more efficient too (in that a fold_right has been eliminated). 

For example, this program
let s = ["x", term_of_string "g(z)"] in
print_term (apply_substitution s (term_of_string "f (x, y, x)"))
shows that f (x, y, x) -> f (g(z), y, g(z)) under the substitution {(x, g(z))}.  

Composition of substitutions

If sig1 and sig2 are two substitutions, we can compute their composite sig1 o sig2 like this
let compose_subst sig1 sig2 =
  (List.map (fun (v, t) -> (v, apply_substitution sig1 t)) sig2)
  @(let vs = List.map fst sig2 in 
    List.filter (fun (x, t) -> not (List.mem x vs)) sig1)

That's quite subtle! Informally this is saying : first apply to the terms in sig2, the substitutions of sig1. Then, filter out those substitutions in sig1 that are in variables modified by sig2. The result is the concatenation of those lists. 

Here's an example. If sig1={(x, g(x,y))} and sig2={(y, h(x,z)), (x, k(x))} then we'd expect sig1 o sig2 = {(y, h(g(x,y), z)), (x, k(g(x, y)))}. The following program can be used to test our expectations.

 let sig1 = ["x", term_of_string "g (x, y)"] in
let sig2 = ["y", term_of_string "h (x, z)"; "x", term_of_string "k(x)"] in
let subst=compose_subst sig1 sig2 in
print_string "y -> " ; print_term (List.assoc "y" subst) ; print_newline ();
print_string "x -> " ; print_term (List.assoc "x" subst) ; print_newline ()

Saturday, October 12, 2013

Tuple matching

Tuple matching

Here is an algorithm that can be used to match tuples. For example, on matching (a, b, (c, d)) against (1, 2, (3, 4)) we'd hope to get the set (a = 1, b = 2, c = 3, d = 4).
module type S = sig

  type expr =
  | E_var of string
  | E_const of int
  | E_tuple of expr list

  val const : int -> expr
  val var : string -> expr
  val tuple : expr list -> expr

  val tuple_match : (string * expr) list -> expr -> expr -> (string * expr) list
  val tuple_of_expr : expr -> expr

end

module Expr : S = 
struct

  type expr =
  | E_var of string
  | E_const of int
  | E_tuple of expr list

  let var s = E_var s
  let const i = E_const i
  let tuple t = E_tuple t

  let as_tuple = function
    | E_tuple l as t -> t
    | _ ->  failwith "Tuple expected"

  let rec tuple_match acc x y =
    match x with
    | E_var (s) -> (s, y)::acc
    | E_tuple ((h::t)) ->
      let (E_tuple l) = as_tuple y in
      let acc = tuple_match acc (tuple t) (tuple (List.tl l)) in
      tuple_match acc h (List.hd l)
    | E_tuple [] -> acc
    | _ as unk -> failwith "Match failure"

end
For example,
(* (a, b, t) |- (1, 2, (3, 4)) *)
let x = Expr.tuple 
  [Expr.var "a"; Expr.var "b"; Expr.var "t" ] in
let y = Expr.tuple [Expr.const 1; Expr.const 2; 
   Expr.tuple [Expr.const 3; Expr. const 4] ] in 
Expr.tuple_match [] x y
;;
in the top-level prints
- : (string * Expr.expr) list =
[("a", Expr.E_const 1); ("b", Expr.E_const 2);
 ("t", Expr.E_tuple [Expr.E_const 3; Expr.E_const 4])]
whereas
(* (a, b, (c, d)) |- (1, 2, (3, 4))*)
let x = Expr.tuple 
  [Expr.var "a"; Expr.var "b"; 
   Expr.tuple [Expr.var "c"; Expr.var "d"]] in
let y = Expr.tuple 
  [Expr.const 1; Expr.const 2; 
   Expr.tuple [Expr.const 3; Expr. const 4]] in 
Expr.tuple_match [] x y
;;
results in this
- : (string * Expr.expr) list =
[("a", Expr.E_const 1); ("b", Expr.E_const 2); ("c", Expr.E_const 3);
 ("d", Expr.E_const 4)]