借学习 thunk 之机看看 delay-evaluation

SICP当然指的是大名鼎鼎的《计算机程序的构造与解释》啦。偶然间聊天的时候突然想到 SICP 在第三章中的流模型中举过一些例子，正好来作为本文的一部分。

（需要说明的是，要理解无穷流的表示需要一定的递归基础，而且我在文中并未对每个序列的定义都给出足够详细的解释。如果你对递归和一些数学递归定义不太熟悉的话，觉得某些例子难以理解是很正常的。实在看不下去的话可以去读读 SICP，或者直接看下一个例子。这部分只能算作我对于 SICP 知识的一个整理。）

在 SICP 中这样写道：

从抽象的观点来看，一个流就是一个序列，但是与列表不同的是，它将使我们能够用流去表示非常长（甚至是无穷的）序列。

流是一种非常巧妙的想法，使我们可以利用各种序列操作，但又不会带来将序列作为表操作而引起的代价。

利用流结构，我们能同时得到两个世界（指命令式和函数式）上最好的东西：如此形成的程序可以像序列操作那么优雅，同时又能得到递归计算的效率。这里的基本想法是做出一种安排，只是部分地构造出流的结构，并将这样的部分结构送给处理流的程序，如果使用者需要访问这个流的尚未构造出来的那个部分，那么这个流就会自动地继续构造下去，但是只做出满足需要的那部分。这一做法造成了一种假象，就好像整个流都存在着一样。

一言以蔽之，流是这样一种数据结构，它和表一样是一种序列，但是它仅仅在需要的时候才对它的一部分进行求值。表（或者说是 cons cell）有构造函数 cons 和析构函数 car 、 cdr （这里的构造和析构按字面意思理解即可）。对于流自然也应该有相似的函数，这里就叫它们 stream-car 、 stream-cdr 和 stream-cons 吧。下面给出一种实现：

;; stream-car
(defun stream-car (s) (car s))
;; stream-cdr
(defun stream-cdr (s) (thunk-force (cdr s)))
;; stream-cons
(defmacro stream-cons (a b) `(cons ,a (thunk-delay ,b)))
;; stream-nullp
(defun stream-nullp (s) (null s))
;; stream-null
(defvar stream-null nil)

上面的实现借助了 Elisp 原生数据类型 cons 来作为骨架，并且使用 nil 来作为流的结尾。如果你想玩的花一点也可以自己实现另一种 cons cell，这里就从简了。

从上面的定义中可以看出，一个流就是一个 cons cell，它的 car 部分是一个值，而 cdr 部分是一个待求值的 thunk，每当我们对一个流使用 cdr 操作时，流的 cdr 部分被求值得到一个新的 cons cell，或者是 null 值。

上面的定义看上去也没啥意思，但是就像 SICP 中说的，我们可以构建无穷流！使用上面的函数可以使用有限的代码来表示无穷这一概念，因为延时求值的原因，如果不需要取到某个值的话，是不需要对其进行计算的。使用以下代码我们就可以定义出全为 1 的无穷序列：

(setq ones (stream-cons 1 ones))

使用一个辅助函数，我们可以更加清楚地看出这个定义是不是真正起到了作用：

;; stream-ref, get the nth value of stream s, 0 based
;; with no tail-call optimism, use loop instead
(defun stream-ref (s n)
  (if (< n 0) (error "negative index" n)
    (cl-loop for i from 0 to (- n 1)
         if (stream-nullp s)
           return stream-null
         else
           do (setq s (stream-cdr s))
         finally return (stream-car s))))

(stream-ref ones 10000)
=> 1

上面的 stream-ref 写法有点不清真，但是没有尾调用优化还是悠着点吧，不过 cl-loop 的可读性很不错，就算没学过也看得懂。可以看到，就算到了 10000，结果也还是 1，再往上加也没问题。

下面我们再考虑另一个例子，使用流表示正整数序列。正整数序列指 1 2 3 4 5 ... ，若将其记为 N，很显然有 cdr (N - ones) = N（即 cdr(N) - ones = N）。根据这个关系式，我们借用一个辅助函数可以给出正整数序列的表示：

(defun stream-map (func &rest ss)
  (if (stream-nullp (car ss)) stream-null
    (stream-cons (apply func (mapcar 'stream-car ss))
         (apply 'stream-map func (mapcar 'stream-cdr ss)))))

(defun stream-add (sa sb)
  (stream-map '+ sa sb))

(setq pos-integer (stream-cons 1 (stream-add ones pos-integer)))

(stream-ref pos-integer 10000)
=> 10001

可惜的是，Elisp 的性能还是太贫弱了，我使用 114513 作为索引值会直接爆掉，连提示都不给。

看了上面的例子，如果感觉有意思的话可以去买一本 SICP 看看（笑），下面以一些数学运算来作为例子，以说明流的灵活与强大。

(defun fibgen (a b)
   (stream-cons a (fibgen b (+ a b))))

(setq fib (fibgen 1 1))
(stream-ref fib 5)
=> 8

;;another approach
(setq fib2 (stream-cons 0
            (stream-cons 1
              (stream-add (stream-cdr fib2)
                fib2))))
(stream-ref fib2 5)
=> 5
;; this seq is start with zero, so it's 5

(defun stream-scale (s n)
  (stream-map (lambda (x) (* x n)) s))

(setq fives (stream-scale ones 5))
=> 5

(setq double (stream-cons 1 (stream-scale double 2)))
(stream-ref double 10)
=> 1024

(setq three (stream-cons 1 (stream-scale three 3)))
(stream-ref three 3)
=> 27

(defun stream-filter (func s)
  (if (stream-nullp s) stream-null
    (if (funcall func (stream-car s))
    (stream-cons (stream-car s) (stream-filter func (stream-cdr s)))
      (stream-filter func (stream-cdr s)))))

(setq x5 (stream-filter (lambda (x) (zerop (mod x 5))) pos-integer))

(stream-ref x5 10)
=> 55


(setq start100 (stream-filter (lambda (x) (> x 100)) pos-integer))

(stream-ref start100 0)
=> 101

(defun stream-mul (s1 s2)
  (stream-cons (* (stream-car s1) (stream-car s2))
           (stream-mul (stream-cdr s1) (stream-cdr s2))))
(setq factor (stream-cons 1 (stream-mul factor pos-integer)))

(stream-ref factor 3)
=> 6
;;another approach
(defun stream-mul (s1 s2)
  (stream-map '* s1 s2))

(defun stream-const (n)
  (stream-cons n (stream-const n)))

(defun stream-partial-sums (s)
  (stream-cons (stream-car s) (stream-add (stream-const (stream-car s))
                      (stream-partial-sums (stream-cdr s)))))


(setq pos-int (stream-partial-sums pos-integer))

(stream-ref pos-int 3)
=> 10 ; 1 3 6 10 ...

这一部分是我觉得最有意思的部分，通过流来表示函数的泰勒级数展开，可以用流来表示各种各样的初等函数。

①级数的积分序列，对 a0 + a1x + a2x^2 + a3x^3 + .... 的积分为 c + a0x + 1/2*a1x^2 + 1/3*a2^3 ...

(defun stream-integrate (s)
  (stream-mul (stream-map (lambda (x) (/ 1.0 x)) pos-integer)
              s))

(setq one2 (stream-integrate pos-integer))
(stream-ref one2 200)
=> 1.0

② 三角函数级数 展开和 指数函数级数 展开，通过它们的定义和相互关系可以这样定义：

(setq exp-series
      (stream-cons 1 (stream-integrate exp-series)))

(stream-ref exp-series 3)
=> 0.16666667


(setq cosine-series
      (stream-cons 1
           (stream-scale
            (stream-integrate sine-series)
            -1)))
(setq sine-series
      (stream-cons 0
           (stream-integrate cosine-series)))

(stream-ref cosine-series 4)
=> 0.041666667 ;1 0 -1/2 0 1/24 ....

(stream-ref sine-series 3)
=> -0.1666666666 ;0 1 0 -1/6 ....

③ 1/S的级数 ，也就是说，若f(x)对应级数为 S，要求出 1/f(x)对应的级数。SICP上面给出了对应的解释，这里我直接给出实现：

(defun stream-mul-series (s1 s2)
  (stream-cons (* (stream-car s1) (stream-car s2))
           (stream-add (stream-mul-series s1 (stream-cdr s2))
               (stream-scale (stream-cdr s1)
                     (stream-car s2)))))

(defun stream-div0 (s)
  (stream-cons (/ 1.0 (stream-car s))
           (stream-mul-series
        (stream-scale (stream-cdr s) (/ -1.0 (stream-car s)))
        (stream-div0 s))))
(defun stream-div-series (s1 s2)
  (stream-mul-series s1 (stream-div0 s2)))

(setq tan-series (stream-div-series sine-series cosine-series))
(stream-ref tan-series 5)
=> 0.133333333333 ;0 1 0 1/3 0 2/15 ......

上面定义的 mul-series 和 div-series 实现了流的卷积和，由 sinx 和 cosx 的级数可以计算 tanx 的级数。

以上级数计算的正确性可以去网站【3】进行验证。

延时求值与动态规划看上去貌似很难联系到一起，实际上也应该如此。但是某天我在计算斐波那契数列时，偶然想到了一种表达方式，且听我细细道来。

已知 Fib(n) = Fib(n - 1) + Fib(n - 2) ，最简单的递归解法就是这样：

(defun fib (n)
  (cond
   ((= n 0) 0)
   ((= n 1) 1)
   (t (+ (fib (- n 1)) (fib (- n 2))))))

这种解法最为直观，但是由于过多的重复计算导致其效率不高，可以考虑把已经计算过的结果存储下来避免重复计算。这个时候就可以用到延时求值了：

(defun fib2 (n)
  (let ((table (make-vector (+ n 2) 0)))
    (cl-loop for i from 2 to n
             do (let ((i i))
                  (aset table i
                  (thunk-delay (+ (thunk-force (aref table (- i 1)))
                                  (thunk-force (aref table (- i 2))))))))
    (aset table 0 (thunk-delay 0))
    (aset table 1 (thunk-delay 1))
    (thunk-force (aref table n))))

我当然知道你知道迭代的方法怎么写，但是这里我使用了所谓的记忆表方法，把每个 fib 值以 thunk 的形式存储下来，需要时直接在数组中进行访问，以此避免了重复计算。

那么，这样做相对于没有使用 thunk 的记忆表方法有什么优点呢？如果没有 thunk 的话，我们就需要添加条件语句来判断某个 fib 值是否被求取，使用 thunk 可以让我们直接 force 而省略掉是否求值的判断逻辑。换言之，是否求值的逻辑被 转移到 thunk 里面去了 。

你也许会对上面代码中的 (let ((i i)) 感到迷惑，我一开始没有加上这一行，但是没有它的话不行。由于 thunk-delay 把表达式的求值延后了， (- i 1) 和 (- i 2) 的求值得等到 force 的时候再开始了，但是根据闭包作用域，它们所指向的 i 值并不是 aset 时的值，而是 cl-loop 结束后的 i 值，这个值并不等于它们的 index，由此会导致求值时进入死循环。所以需要在外面包一层 let 来保证 force 时使用的是 index。如果使用递归而不是迭代就不会有这样的问题了（还是 Scheme 好点）。

斐波那契数列毕竟只是一个非常简单的例子，这里介绍一个稍微复杂点的题目：

Projecteuler81 【4】。题目大意是：给出一个矩阵，要求你从左上角走到右下角，你只能向右或向下运动，要求求出最小路径和（即走过路径上的所有数的和）。

由于题目数据放在这里占的空间有点大，这里我直接使用题目示例中的数据并定义 thunk 矩阵：（加上一些辅助函数）

(setq mat0
      [[131 673 234 103 18]
       [201 96 342 965 150]
       [630 803 746 422 111]
       [537 699 497 121 956]
       [805 732 524 37 331]]
      )

(defun mat-ref (mat i j)
  (aref (aref mat i) j))
(defun mat-set (mat i j newval)
  (aset (aref mat i) j newval))
(defun make-matrix (i j)
  (let ((ret (make-vector i nil)))
    (cl-loop for k from 0 to (- j 1)
         do (aset ret k (make-vector j nil)))
    ret))

(setq thunk-mat (make-matrix 5 5))

thunk 矩阵的含义是：它的 (i, j) 元素的值是从 (0, 0) 到 (i, j) 的最小路径和。

对于 (i, j) 处的最小路径和，我们只需要知道 thunk 矩阵在 (i, j - 1) 和 (i - 1, j) 处的值就行了，因为 (i, j) 处的最小路径和就是原矩阵在 (i, j) 的值加上它上面和左边的最小路径和中的较小值。因此，我们可以这样定义 thunk-mat 矩阵：

(cl-loop
 for i from 0 to 4
 do (cl-loop
     for j from 0 to 4
     do (mat-set thunk-mat i j
         (let ((i i)
               (j j))
           (thunk-delay
            (+ (mat-ref mat0 i j)
               (min (thunk-force (mat-ref thunk-mat i (- j 1)))
                    (thunk-force (mat-ref thunk-mat (- i 1) j)))))))))

同时，要考虑到矩阵的第一列左边没有元素，第一列上面没有元素，所以需要对第一行和第一列做特殊处理。因为从 (0, 0) 到 (0, i) 只能一直向右，从 (0, 0) 到 (i, 0) 只能一直向下，所以边界条件的确定是很容易的：

(cl-loop
 for i from 0 to 4
 do (mat-set thunk-mat i 0
         (let ((i i))
           (thunk-delay (+ (thunk-force (mat-ref thunk-mat (- i 1) 0))
                           (mat-ref mat0 i 0)))))
 do (mat-set thunk-mat 0 i
         (let ((i i))
           (thunk-delay (+ (thunk-force (mat-ref thunk-mat 0 (- i 1)))
                           (mat-ref mat0 0 i))))))
(mat-set thunk-mat 0 0
     (thunk-delay (mat-ref mat0 0 0)))

接下来，我们只需要一条简单的表达式就可以了：

(thunk-force (mat-ref thunk-mat 4 4))
=> 2427

结果与问题页面的答案相同。

要说的话， thunk 的作用也就是减少一层 if 逻辑，其他的也没啥。用 thunk 来解决这类问题的步骤就是先确定一般的递推式子，然后把握好边界条件就好了。