open Jog open Zone;;
let top = [[0; 1; 2; 3; 4]; [0; 1; 2; 4; 5]; [0; 1; 2; 5; 6]; [0; 1; 2; 6; 3]]
and ini = ref (fun _ -> raise (OffTrack 3))
and ds2 i j = if i=j then 0. else if i>2 && j>2 then 4. /. 3. else 1.
and fac = [0; 1; 2; 3]
and bid = 3 in
let fa = Array.of_list fac in
let n = len fa in
let gd = aini n (fun i -> aini n (fun j -> 0.5 *. (
  (ds2 fa.(n-1) fa.(i)) +. (ds2 fa.(n-1) fa.(j)) -. (ds2 fa.(i) fa.(j))))) in
  for j = 0 to n-2 do gd.(j).(n-1) <- 0.; gd.(n-1).(j) <- 0. done; gd.(n-1).(n-1) <- 1.;
let ip a b = let s = ref 0. in for i=0 to n-1 do for j=0 to n-1 do
   s := !s +. a.(i) *. b.(j) *. gd.(i).(j) done done; !s in
let pON = Linear.gs ip (let fON = Linear.iden n in
     (let t = fON.(bid) in fON.(bid) <- fON.(n-2); fON.(n-2) <- t); fON) in
let bndry = Glue.morphgen (nz ds2) top
    [(fac, ((fun  x -> x),
      (fun t -> ini := t.curl; {t with curl = (fun (k, ball) -> ball)})))] in
    (pON, Linear.mul (Linear.inv pON) (!ini (bid, pON)));;

(* should yield (_, iden 4) *)