2007/07/06

Some Folds & Unfolds

吃晚餐時心血來潮,想用 fold/unfold 亂湊一陣,寫個矩陣乘法出來,然後看有沒有機會做點有趣的 derivations(最好可以玩到 fold fusion XD)。我把矩陣定義為 list of lists,採 row-major order。以下是目前成果(in OCaml):

(* some helper functions *)

let null = function [] -> true | x :: xs -> false

let first = function (x, y) -> x

let second = function (x, y) -> y

(* foldr and unfoldr on lists *)

let rec foldr f e xs =
  match xs with
    [] -> e
  | x :: xs -> f x (foldr f e xs)

let rec unfoldr p f s =
  if p s then [] else let (x, ns) = f s in x :: unfoldr p f ns

(* now folds and unfolds follow *)

let sum = foldr (fun x y -> x + y) 0

let map f = foldr (fun x xs -> f x :: xs) []

let couple xs ys = unfoldr
  (function ([], ys) -> true | (xs, []) -> true | (xs, ys) -> false)
  (function (x :: xs, y :: ys) -> ((x, y), (xs, ys)))
  (xs, ys)

let allpairs xs ys =
    let rec aps = function
        ([], y :: ys, zs) -> aps (zs, ys, zs)
      | (x :: xs, y :: ys, zs) -> ((x, y), (xs, y :: ys, zs))
    in
      unfoldr (function ([], y :: [], zs) -> true
               | (xs, ys, zs) -> (null ys) || (null zs)) aps (xs, ys, xs)

let innerproduct xs ys =
  foldr (fun (x, y) n -> x * y + n) 0 (couple xs ys)

let transpose =
  let rec t xs yss =
    match (xs, yss) with
      ([], yss) -> []
    | (x :: xs, []) -> [x] :: t xs []
    | (x :: xs, ys :: yss) -> (x :: ys) :: t xs yss
  in
    foldr t []

其中 innerproduct 已經是 hylomorphism 了!可惜最後我寫昏了頭,把矩陣乘法寫成

let multiply a b =
  sum (map (fun (x, y) -> innerproduct x y) (allpairs a (transpose b)))

竟然把乘得的矩陣元素全加起來了 XD。再想想,因為要產生 list of lists,矩陣乘法大約是雙重 unfold,恐怕不會太漂亮,而且寫出來應該會複雜到讓我很難做 derivation,所以先就此打住 XD。

--
當作寫 fold/unfold 的練習 XD。


發現 allpairstranspose 寫得太複雜,其實用 doubly-nested fold 就可以變得很簡潔:

let allpairs xs ys =
  let a x xs =
    foldr (fun y ys -> (x, y) :: ys) [] ys @ xs
  in
    foldr a [] xs

let transpose =
  let t = fun xs -> function
                      [] -> foldr (fun x xss -> [x] :: xss) [] xs
                   | yss -> map (fun (x, xs) -> x :: xs) (couple xs yss)
  in foldr t []

兩個函式的新版本看起來都有改善空間,來導一下試試看吧 XD。


transpose 的改良很容易。t xsyss case 看得出來是 hylomorphism,目前在 unfold/fold 之間會產生一個 temporary list。解決辦法當然是明確引進 hylomorphism:

(* hylor f e p g s = foldr f e (unfoldr p g s) *)
let rec hylor f e p g s =
  if p s then e else let (x, ns) = g s in f x (hylor f e p g ns)

另外我把 couple 的定義改寫為

let hasnull =
  (function ([], ys) -> true | (xs, []) -> true | (xs, ys) -> false)

let couple xs ys = unfoldr
  hasnull
  (fun (x :: xs, y :: ys) -> ((x, y), (xs, ys)))
  (xs, ys)

讓最後的 hylomorphism 短一點。把 map 寫成 fold 化簡、把 couple 寫成 unfold,求式就長成

foldr (fun (x, xs) xss -> (x :: xs) :: xss) []
 (unfoldr nullpair (function (x :: xs, y :: ys) -> ((x, y), (xs, ys))) (xs, yss))

最後依 hylomorphism 的定義把上式收到 hylortranspose 就變成:

let transpose =
  let t = fun xs -> function
                      [] -> foldr (fun x xss -> [x] :: xss) [] xs
                   | yss -> hylor (fun (x, xs) xss -> (x :: xs) :: xss)
                                  []
                                  hasnull
                                  (fun (x :: xs, y :: ys) -> ((x, y), (xs, ys)))
                                  (xs, yss)
  in foldr t []

看起來還不錯 XD。innerproduct 也可以循此途徑改善一點效率,結果是

let innerproduct xs ys = hylor
  (fun (x, y) n -> x * y + n)
  0
  hasnull
  fun (x :: xs, y :: ys) -> ((x, y), (xs, ys)))
  (xs, ys)

--
顯然 type system 得等明天了 XD。


allpairs 比較礙眼的地方是那個 @,雖然對時間複雜度沒有影響,但常數能消還是盡量消。其實用直覺很快,不過我照規矩來做一次 XD。首先,@ 也是一個 fold,let cons x xs = x :: xs,從定義容易推得

a @ b = foldr cons b a

為了方便,let f = fun x y ys -> (x, y) :: ys,於是欲化簡的算式就變成 foldr (f x) [] ys @ xs。再把「@ 是 fold」的事實代入,就得到 foldr cons xs (foldr (f x) [] ys)。現在進行 fold fusion:

  foldr cons xs (f x y ys)
= { definition of f }
  foldr cons xs ((x, y) :: ys)
= { definition of foldr }
  cons (x, y) (foldr cons xs ys)
= { definition of cons }
  (x, y) :: foldr cons xs ys
= { definition of f (!!) }
  f x y (foldr cons xs ys)

因此依 fold-fusion theorem,我們得到

  foldr (f x) [] ys @ xs
= { @ is a fold }
  foldr cons xs (foldr (f x) [] ys)
= { fold fusion }
  foldr (f x) (foldr cons xs []) ys
= { definition of foldr }
  foldr (f x) xs ys

直觀上當然很簡單,本來 fold 先累積在 empty list 內,再和 xs 串接,現在則直接把 fold 的結果累積在 xs 裡面。不過如穆老師所說,"it's how we get there that matters" XD。(其實這道題不算結果簡單,而是過程簡單…不管啦 XD。)最後的 allpairs 就變成

let allpairs xs ys =
  let a x xs =
    foldr (fun y ys -> (x, y) :: ys) xs ys
  in
    foldr a [] xs

--
無論如何,總算是玩到 fold fusion 了 XD。


a @ b = foldr cons b a 其實也可以用(很簡單的)fold fusion 得到。除了 id 以外的所有(簡單的)folds 都可以這麼做,即利用 "id is a fold" 的性質。邏輯上可以先確立 id 是 fold,其他所有的(簡單的)folds 就可以用 fold fusion 導出其 fold 形式。把 @ b 視為一個 function of type 'a list -> 'a list,如下推導:

  a @ b
= { identity function }
  (id a) @ b
= { id is a fold }
  (foldr cons [] a) @ b
= { fold fusion, see below }
  foldr cons ([] @ b) a
= { definition of @ }
  foldr cons b a

(* condition check for fold fusion *)

  (cons x xs) @ b
= { definition of cons }
  (x :: xs) @ b
= { definition of @ }
  x :: (xs @ b)
= { definition of cons }
  cons x (xs @ b)

另外穆老師上課補充一個定理:

    h = foldr f e    for some f, e
iff h xs = h ys ==> h (x :: xs) = h (x :: ys)

我把這個定理套到 couple xs 上,它告訴我這是一個 fold,可是我卻造不太出來。現在只想到借用 tupling 之類的手法,可是又不太對勁 XD。

Labels: ,

Anonymous Anonymous7/07/2007 12:58 am 說:

Weijin改用blogger了...為啥我都不知道

 

<< 回到主頁