M.Hiroi's Home Page

Common Lisp Programming

お気楽 CLOS プログラミング入門

[ PrevPage | CLOS | NextPage ]

平衡木

木構造のお話です。二分木は左右の部分木のバランスが崩れると、性能が劣化する欠点があります。極端な例ですが、ソートされたデータを二分木に挿入していくと、データは右側の木にしか挿入されず、連結リストと同じ線形探索になってしまいます。

これを補うため、木のバランスを一定の範囲に収める「平衡木 (balanced tree)」が考案されています。有名なところでは、AVL 木、赤黒木、AA 木、2-3 木、2-3-4 木、B 木、B* 木などがあります。AVL 木、赤黒木、AA 木は二分木ですが、2-3 木、2-3-4 木、B 木、B* 木は多分木になります。

このほかに、自己調整二分木と呼ばれるスプレー木 (splay tree) があります。スプレー木は、木のバランスが一定の範囲内に収まる保障はありませんが、複数回アクセスしたときの平均実行時間が log2 N に比例するという面白い性質があります。ようするに、一回あたり長い時間がかかる処理があったとしても、全体で平均してみると O(log2 N) になります。

今回は CLOS の例題として、赤黒木、AA 木、スプレー木を使った連想配列 treemap を作ってみましょう。アルゴリズムの詳しい説明は以下に示す拙作のページをお読みくださいませ。

●仕様

●プログラムリスト

●実行例

* (require :aatree "aatree.lsp")
("AATREE")
* (use-package :aatree)
T
* (defvar a (make-treemap #'string= #'string<))
A
* (tree-emptyp a)
T
* (tree-count a)
0
* (defvar ks '("foo" "bar" "baz" "oops"))
KS
* (defvar vs '(10 20 30 40))
VS
* (mapc (lambda (k v) (tree-set a k v)) ks vs)
("foo" "bar" "baz" "oops")
* (tree-emptyp a)
NIL
* (tree-count a)
4

* (mapc (lambda (k) (print (tree-exists a k))) '("foo" "FOO"))

T
NIL
("foo" "FOO")
* (mapc (lambda (k) (print (tree-exists a k))) ks)

T
T
T
T
("foo" "bar" "baz" "oops")
* (mapc (lambda (k) (print (tree-get a k))) '("foo" "FOO"))

10
NIL
("foo" "FOO")
* (mapc (lambda (k) (print (tree-get a k))) ks)

10
20
30
40
("foo" "bar" "baz" "oops")
* (tree-max a)
("oops" . 40)
* (tree-min a)
("bar" . 20)

* (tree-mapc a (lambda (k v) (print (cons k v))))

("bar" . 20)
("baz" . 30)
("foo" . 10)
("oops" . 40)
NIL
* (defun to-alist (tm) (tree-fold-right tm nil (lambda (a k v) (cons (cons k v) a))))
TO-ALIST
* (to-alist a)
(("bar" . 20) ("baz" . 30) ("foo" . 10) ("oops" . 40))

* (tree-max-del a)
T
* (to-alist a)
(("bar" . 20) ("baz" . 30) ("foo" . 10))
* (tree-count a)
3
* (tree-min-del a)
T
* (to-alist a)
(("baz" . 30) ("foo" . 10))
* (tree-count a)
2
* (tree-del a "foo")
T
* (to-alist a)
(("baz" . 30))
* (tree-count a)
1
* (tree-del a "baz")
T
* (to-alist a)
NIL
* (tree-count a)
0
* (tree-emptyp a)
T
* (tree-del a "foo")
NIL

* (mapc (lambda (k v) (tree-set a k v)) ks vs)
("foo" "bar" "baz" "oops")
* (to-alist a)
(("bar" . 20) ("baz" . 30) ("foo" . 10) ("oops" . 40))
* (mapc (lambda (k v) (tree-set a k (* v 10))) ks vs)
("foo" "bar" "baz" "oops")
* (to-alist a)
(("bar" . 200) ("baz" . 300) ("foo" . 100) ("oops" . 400))
* (tree-push a "FOO" 0)
T
* (tree-get a "FOO")
(0)
* (tree-push a "FOO" 1)
NIL
* (tree-get a "FOO")
(1 0)
* (tree-pop a "FOO")
1
* (tree-get a "FOO")
(0)
* (tree-pop a "FOO")
0
* (tree-get a "FOO")
NIL
* (tree-exists a "FOO")
T
* (tree-pop a "FOO")
NIL
* (tree-del a "FOO")
T
* (tree-exists a "FOO")
NIL
* (tree-pop a "FOO")
NIL

赤黒木とスプレー木の実行結果は同じになるので省略します。


●簡単なテスト

●AA 木

リスト : AA 木の簡単なテスト

(require :aatree "aatree.lsp")
(use-package :aatree)

;; 数列の生成
(defun iota (s e &optional (a nil))
  (if (> s e)
      a
    (iota s (1- e) (cons e a))))

(defun tabulate (fn s e &optional (a nil))
  (if (> s e)
      a
    (tabulate fn s (1- e) (cons (funcall fn e) a))))

;;; AA 木のチェック
(defun check-aatree (node)
  (cond
   ((not (aatree::node-emptyp node))
    (if (not (= (aatree::node-height node) (+ (aatree::node-height (aatree::node-left node)) 1)))
        (error "aa tree error1"))
    (if (and (not (= (aatree::node-height node) (aatree::node-height (aatree::node-right node))))
             (not (= (aatree::node-height node) (+ (aatree::node-height (aatree::node-right node)) 1))))
        (error "aa tree error2"))
    (if (= (aatree::node-height node) (aatree::node-height (aatree::node-right (aatree::node-right node))))
        (error "aa tree error3"))
    (check-aatree (aatree::node-left node))
    (check-aatree (aatree::node-right node)))))

(defun test-sub (xs)
  (let ((a (make-treemap #'= #'<)))
    (print "insert")
    (dolist (x xs)
      (tree-set a x t)
      (check-aatree (aatree::tree-root a)))
    (print (tree-emptyp a))
    (print (tree-count a))
    (print "search")
    (dolist (x xs)
      (if (not (tree-get a x))
          (error "test search error")))
    (print "delete")
    (dolist (x xs)
         (tree-del a x)
         (check-aatree (aatree::tree-root a)))
    (print (tree-emptyp a))
    (print (tree-count a))))

(defun test (n)
  (print "昇順")
  (test-sub (iota 1 n))
  (print "逆順")
  (test-sub (nreverse (iota 1 n)))
  (print "ランダム")
  (test-sub (tabulate (lambda (x) (random 1d0)) 1 n)))

;;;
;;; 表示
;;;
(defun print-aatree (n node)
  (unless (aatree::node-emptyp node)
    (print-aatree n (aatree::node-left node))
    (dotimes (x (- n (aatree::node-height node))) (princ "  "))
    (format t "(~A, ~A)~%" (aatree::node-key node) (aatree::node-height node))
    (print-aatree n (aatree::node-right node))))

(defun print-tree (tm)
  (format t "--------~%")
  (print-aatree (aatree::node-height (aatree::tree-root tm)) (aatree::tree-root tm)))

(defun test1 (xs)
  (let ((a (make-treemap #'= #'<)))
    (dolist (x xs)
      (tree-set a x t)
      (print-tree a))
    (dolist (x xs)
      (format t "~A~%" (tree-exists a x)))
    (dolist (x xs)
      (tree-del a x)
      (print-tree a))))
(load "testaa.lsp")

; ・・・略・・・
T
* (test 10000)

"昇順"
"insert"
NIL
10000
"search"
"delete"
T
0
"逆順"
"insert"
NIL
10000
"search"
"delete"
T
0
"ランダム"
"insert"
NIL
10000
"search"
"delete"
T
0
0
* (test1 (iota 1 10))
--------
(1, 1)
--------
(1, 1)
(2, 1)
--------
  (1, 1)
(2, 2)
  (3, 1)
--------
  (1, 1)
(2, 2)
  (3, 1)
  (4, 1)
--------
  (1, 1)
(2, 2)
  (3, 1)
(4, 2)
  (5, 1)
--------
  (1, 1)
(2, 2)
  (3, 1)
(4, 2)
  (5, 1)
  (6, 1)
--------
    (1, 1)
  (2, 2)
    (3, 1)
(4, 3)
    (5, 1)
  (6, 2)
    (7, 1)
--------
    (1, 1)
  (2, 2)
    (3, 1)
(4, 3)
    (5, 1)
  (6, 2)
    (7, 1)
    (8, 1)
--------
    (1, 1)
  (2, 2)
    (3, 1)
(4, 3)
    (5, 1)
  (6, 2)
    (7, 1)
  (8, 2)
    (9, 1)
--------
    (1, 1)
  (2, 2)
    (3, 1)
(4, 3)
    (5, 1)
  (6, 2)
    (7, 1)
  (8, 2)
    (9, 1)
    (10, 1)
T
T
T
T
T
T
T
T
T
T
--------
    (2, 1)
    (3, 1)
  (4, 2)
    (5, 1)
(6, 3)
    (7, 1)
  (8, 2)
    (9, 1)
    (10, 1)
--------
    (3, 1)
  (4, 2)
    (5, 1)
(6, 3)
    (7, 1)
  (8, 2)
    (9, 1)
    (10, 1)
--------
  (4, 1)
  (5, 1)
(6, 2)
  (7, 1)
(8, 2)
  (9, 1)
  (10, 1)
--------
  (5, 1)
(6, 2)
  (7, 1)
(8, 2)
  (9, 1)
  (10, 1)
--------
  (6, 1)
(7, 2)
  (8, 1)
(9, 2)
  (10, 1)
--------
  (7, 1)
(8, 2)
  (9, 1)
  (10, 1)
--------
  (8, 1)
(9, 2)
  (10, 1)
--------
(9, 1)
(10, 1)
--------
(10, 1)
--------
NIL

●赤黒木

リスト : 赤黒木の簡単なテスト

(require :rbtree "rbtree.lsp")
(use-package :rbtree)

;; 数列の生成
(defun iota (s e &optional (a nil))
  (if (> s e)
      a
    (iota s (1- e) (cons e a))))

(defun tabulate (fn s e &optional (a nil))
  (if (> s e)
      a
    (tabulate fn s (1- e) (cons (funcall fn e) a))))

;;; 赤黒木のチェック
(defun check-rbtree (node)
  (cond
   ((rbtree::node-emptyp node) 0)
   (t
    (when (rbtree::node-redp node)
      (when (or (rbtree::node-redp (rbtree::node-left node))
                (rbtree::node-redp (rbtree::node-right node)))
        (error "rbtree error1")))
    (let ((a (check-rbtree (rbtree::node-left node)))
          (b (check-rbtree (rbtree::node-right node))))
      (when (/= a b)
        (error "rbtree error2"))
      (when (rbtree::node-blackp node)
        (incf a))
      a))))

(defun test-sub (xs)
  (let ((a (make-treemap #'= #'<)))
    (print "insert")
    (dolist (x xs)
      (tree-set a x t)
      (check-rbtree (rbtree::tree-root a)))
    (print (tree-emptyp a))
    (print (tree-count a))
    (print "search")
    (dolist (x xs)
      (if (not (tree-get a x))
          (error "test search error")))
    (print "delete")
    (dolist (x xs)
         (tree-del a x)
         (check-rbtree (rbtree::tree-root a)))
    (print (tree-emptyp a))
    (print (tree-count a))))

(defun test (n)
  (print "昇順")
  (test-sub (iota 1 n))
  (print "逆順")
  (test-sub (nreverse (iota 1 n)))
  (print "ランダム")
  (test-sub (tabulate (lambda (x) (random 1d0)) 1 n)))

;;;
;;; 表示
;;;
(defun print-rb-tree (n node)
  (unless (rbtree::node-emptyp node)
    (print-rb-tree (1+ n) (rbtree::node-left node))
    (dotimes (x n) (princ "  "))
    (format t "(~A, ~A)~%" (rbtree::node-key node) (rbtree::node-color node))
    (print-rb-tree (1+ n) (rbtree::node-right node))))

(defun print-tree (tm)
  (format t "--------~%")
  (print-rb-tree 0 (rbtree::tree-root tm)))

(defun test1 (xs)
  (let ((a (make-treemap #'= #'<)))
    (dolist (x xs)
      (tree-set a x t)
      (print-tree a))
    (dolist (x xs)
      (format t "~A~%" (tree-exists a x)))
    (dolist (x xs)
      (tree-del a x)
      (print-tree a))))
(load "testrb.lsp")

; ・・・略・・・
T
(test 10000)

"昇順"
"insert"
NIL
10000
"search"
"delete"
T
0
"逆順"
"insert"
NIL
10000
"search"
"delete"
T
0
"ランダム"
"insert"
NIL
10000
"search"
"delete"
T
0
0
(test1 (iota 1 10))
--------
(1, BLACK)
--------
(1, BLACK)
  (2, RED)
--------
  (1, RED)
(2, BLACK)
  (3, RED)
--------
  (1, BLACK)
(2, BLACK)
  (3, BLACK)
    (4, RED)
--------
  (1, BLACK)
(2, BLACK)
    (3, RED)
  (4, BLACK)
    (5, RED)
--------
  (1, BLACK)
(2, BLACK)
    (3, BLACK)
  (4, RED)
    (5, BLACK)
      (6, RED)
--------
  (1, BLACK)
(2, BLACK)
    (3, BLACK)
  (4, RED)
      (5, RED)
    (6, BLACK)
      (7, RED)
--------
    (1, BLACK)
  (2, RED)
    (3, BLACK)
(4, BLACK)
    (5, BLACK)
  (6, RED)
    (7, BLACK)
      (8, RED)
--------
    (1, BLACK)
  (2, RED)
    (3, BLACK)
(4, BLACK)
    (5, BLACK)
  (6, RED)
      (7, RED)
    (8, BLACK)
      (9, RED)
--------
    (1, BLACK)
  (2, BLACK)
    (3, BLACK)
(4, BLACK)
    (5, BLACK)
  (6, BLACK)
      (7, BLACK)
    (8, RED)
      (9, BLACK)
        (10, RED)
T
T
T
T
T
T
T
T
T
T
--------
    (2, BLACK)
      (3, RED)
  (4, BLACK)
    (5, BLACK)
(6, BLACK)
    (7, BLACK)
  (8, BLACK)
    (9, BLACK)
      (10, RED)
--------
    (3, BLACK)
  (4, BLACK)
    (5, BLACK)
(6, BLACK)
    (7, BLACK)
  (8, BLACK)
    (9, BLACK)
      (10, RED)
--------
  (4, BLACK)
    (5, RED)
(6, BLACK)
    (7, BLACK)
  (8, RED)
    (9, BLACK)
      (10, RED)
--------
  (5, BLACK)
(6, BLACK)
    (7, BLACK)
  (8, RED)
    (9, BLACK)
      (10, RED)
--------
  (6, BLACK)
    (7, RED)
(8, BLACK)
  (9, BLACK)
    (10, RED)
--------
  (7, BLACK)
(8, BLACK)
  (9, BLACK)
    (10, RED)
--------
  (8, BLACK)
(9, BLACK)
  (10, BLACK)
--------
(9, BLACK)
  (10, RED)
--------
(10, BLACK)
--------
NIL

●スプレー木

リスト : スプレー木の簡単なテスト

(require :splay "splay.lsp")
(use-package :splay)

;; 数列の生成
(defun iota (s e &optional (a nil))
  (if (> s e)
      a
    (iota s (1- e) (cons e a))))

(defun tabulate (fn s e &optional (a nil))
  (if (> s e)
      a
    (tabulate fn s (1- e) (cons (funcall fn e) a))))

(defun test-sub (xs)
  (let ((a (make-treemap #'= #'<)))
    (print "insert")
    (dolist (x xs)
      (tree-set a x t))
    (print (tree-emptyp a))
    (print (tree-count a))
    (print "search")
    (dolist (x xs)
      (if (not (tree-get a x))
          (error "test search error")))
    (print "delete")
    (dolist (x xs)
         (tree-del a x))
    (print (tree-emptyp a))
    (print (tree-count a))))

(defun test (n)
  (print "昇順")
  (test-sub (iota 1 n))
  (print "逆順")
  (test-sub (nreverse (iota 1 n)))
  (print "ランダム")
  (test-sub (tabulate (lambda (x) (random 1d0)) 1 n)))


;;;
;;; 表示
;;;
(defun print-splay-tree (n node)
  (unless (null node)
    (print-splay-tree (1+ n) (splay::node-left node))
    (dotimes (x n) (princ "  "))
    (format t "(~A, ~A)~%" (splay::node-key node) (splay::node-value node))
    (print-splay-tree (1+ n) (splay::node-right node))))

(defun print-tree (tm)
  (format t "--------~%")
  (print-splay-tree 0 (splay::tree-root tm)))

(defun test1 (xs)
  (let ((a (make-treemap #'= #'<)))
    (dolist (x xs)
      (tree-set a x t)
      (print-tree a))
    (dolist (x xs)
      (format t "~A~%" (tree-exists a x))
      (print-tree a))
    (dolist (x xs)
      (tree-del a x)
      (print-tree a))))
(load "testsplay.lsp")

; ・・・略・・・
T
* (test 10000)

"昇順"
"insert"
NIL
10000
"search"
"delete"
T
0
"逆順"
"insert"
NIL
10000
"search"
"delete"
T
0
"ランダム"
"insert"
NIL
10000
"search"
"delete"
T
0
0
(test1 (iota 1 10))
--------
(1, T)
--------
  (1, T)
(2, T)
--------
    (1, T)
  (2, T)
(3, T)
--------
      (1, T)
    (2, T)
  (3, T)
(4, T)
--------
        (1, T)
      (2, T)
    (3, T)
  (4, T)
(5, T)
--------
          (1, T)
        (2, T)
      (3, T)
    (4, T)
  (5, T)
(6, T)
--------
            (1, T)
          (2, T)
        (3, T)
      (4, T)
    (5, T)
  (6, T)
(7, T)
--------
              (1, T)
            (2, T)
          (3, T)
        (4, T)
      (5, T)
    (6, T)
  (7, T)
(8, T)
--------
                (1, T)
              (2, T)
            (3, T)
          (4, T)
        (5, T)
      (6, T)
    (7, T)
  (8, T)
(9, T)
--------
                  (1, T)
                (2, T)
              (3, T)
            (4, T)
          (5, T)
        (6, T)
      (7, T)
    (8, T)
  (9, T)
(10, T)
T
--------
(1, T)
          (2, T)
        (3, T)
          (4, T)
      (5, T)
        (6, T)
    (7, T)
      (8, T)
  (9, T)
    (10, T)
T
--------
  (1, T)
(2, T)
    (3, T)
        (4, T)
      (5, T)
        (6, T)
  (7, T)
      (8, T)
    (9, T)
      (10, T)
T
--------
    (1, T)
  (2, T)
(3, T)
      (4, T)
    (5, T)
      (6, T)
  (7, T)
      (8, T)
    (9, T)
      (10, T)
T
--------
      (1, T)
    (2, T)
  (3, T)
(4, T)
  (5, T)
      (6, T)
    (7, T)
        (8, T)
      (9, T)
        (10, T)
T
--------
        (1, T)
      (2, T)
    (3, T)
  (4, T)
(5, T)
    (6, T)
  (7, T)
      (8, T)
    (9, T)
      (10, T)
T
--------
          (1, T)
        (2, T)
      (3, T)
    (4, T)
  (5, T)
(6, T)
  (7, T)
      (8, T)
    (9, T)
      (10, T)
T
--------
            (1, T)
          (2, T)
        (3, T)
      (4, T)
    (5, T)
  (6, T)
(7, T)
    (8, T)
  (9, T)
    (10, T)
T
--------
              (1, T)
            (2, T)
          (3, T)
        (4, T)
      (5, T)
    (6, T)
  (7, T)
(8, T)
  (9, T)
    (10, T)
T
--------
                (1, T)
              (2, T)
            (3, T)
          (4, T)
        (5, T)
      (6, T)
    (7, T)
  (8, T)
(9, T)
  (10, T)
T
--------
                  (1, T)
                (2, T)
              (3, T)
            (4, T)
          (5, T)
        (6, T)
      (7, T)
    (8, T)
  (9, T)
(10, T)
--------
        (2, T)
      (3, T)
        (4, T)
    (5, T)
      (6, T)
  (7, T)
    (8, T)
(9, T)
  (10, T)
--------
  (3, T)
      (4, T)
    (5, T)
      (6, T)
(7, T)
    (8, T)
  (9, T)
    (10, T)
--------
    (4, T)
  (5, T)
    (6, T)
(7, T)
    (8, T)
  (9, T)
    (10, T)
--------
(5, T)
    (6, T)
  (7, T)
      (8, T)
    (9, T)
      (10, T)
--------
  (6, T)
(7, T)
    (8, T)
  (9, T)
    (10, T)
--------
(7, T)
    (8, T)
  (9, T)
    (10, T)
--------
  (8, T)
(9, T)
  (10, T)
--------
(9, T)
  (10, T)
--------
(10, T)
--------
NIL

●タクシー数

それでは連想配列の簡単な例題として、「タクシー数」を取り上げます。不定方程式 X3 + Y3 = Z には整数解が存在します。たとえば、Z = 2 のときの解は 13 + 13 の一通りしかありませんが、Z の値によっては複数の解が存在します。そして、n 通りの解が存在する Z の中で、最小の値をタクシー数 (taxicab number) といい、n 番目のタクシー数を Ta(n) と書きます。Ta(1) と Ta(2) を示します。

Ta(1) = 2 = 13 + 13
Ta(2) = 1729 = 13 + 123
             = 93 + 103

参考 URL 1 によると、Ta(n) はすべての整数 n に対して存在することが証明されていて、今までに Ta(1) から Ta(6) までのタクシー数が知られているそうです。

ここでは問題を簡単にして、次式を満たす整数 a, b, c, d を列挙することにします。

a3 + b3 = c3 + d3

●プログラムの作成

最近のパソコンは高性能なので、単純な四重ループでも簡単に解くことができます。次のリストを見てください。

リスト : タクシー数

(defun taxi (n)
  (do ((a 1 (+ a 1)))
      ((< n a))
      (do ((b a (+ b 1)))
          ((< n b))
          (do ((c (+ a 1) (+ c 1)))
              ((< n c))
              (do ((d c (+ d 1)))
                  ((<= b d))
                  (let ((e (+ (* a a a) (* b b b))))
                    (when (= (+ (* c c c) (* d d d)) e)
                      (format t "~D: (~D, ~D), (~D,~D)~%" e a b c d))))))))

引数 n が整数の上限値を表します。重複解を削除するため、以下の条件を設定してます。

1. a <= b     # a と b を交換した式を削除
2. c <= d     # c と d を交換した式を削除
3. a < c      # (a, b) と (c, d) を交換した式を削除
4. d < b      # 不要な探索を削除

条件 3 で c は a よりも大きくなるので、d は必ず b よりも小さくなります。これを条件 4 で表しています。あとはとくに難しいところはないと思います。

それでは実行してみましょう。

* (time (taxi 50))
1729: (1, 12), (9,10)
4104: (2, 16), (9,15)
13832: (2, 24), (18,20)
39312: (2, 34), (15,33)
46683: (3, 36), (27,30)
32832: (4, 32), (18,30)
110656: (4, 48), (36,40)
110808: (6, 48), (27,45)
40033: (9, 34), (16,33)
20683: (10, 27), (19,24)
65728: (12, 40), (31,33)
64232: (17, 39), (26,36)
Evaluation took:
  0.000 seconds of real time
  0.009082 seconds of total run time (0.008717 user, 0.000365 system)
  100.00% CPU
  21,717,388 processor cycles
  0 bytes consed

NIL

* (time (taxi 100))
1729: (1, 12), (9,10)
4104: (2, 16), (9,15)
13832: (2, 24), (18,20)

・・・省略・・・

1016496: (47, 97), (66,90)
1009736: (50, 96), (59,93)
684019: (51, 82), (64,75)
Evaluation took:
  0.100 seconds of real time
  0.107375 seconds of total run time (0.107375 user, 0.000000 system)
  107.00% CPU
  257,681,154 processor cycles
  0 bytes consed

NIL

SBCL の場合、引数 n の値を 50 から 100 に増やしても高速に解くことができます。

●3 通りの解を求める

次は Ta(3) を含む 3 通りの解を求めるプログラムを作りましょう。つまり、次式を満たす整数 a, b, c, d, e, f を求めます。

a3 + b3 = c3 + d3 = e3 + f3

ループを六重にするとプログラムがちょっと面倒になるので、連想配列 (hash-table, treemap) を使うことにしましょう。次のリストを見てください。

リスト : タクシー数 (2)

(require :aatree "aatree.lsp")
(require :rbtree "rbtree.lsp")
(require :splay  "splay.lsp")

;;; hash-table
(defun taxi-fast (n m)
  (let ((ht (make-hash-table)))
    (do ((a 1 (+ a 1)))
        ((< n a))
        (do ((b a (+ b 1)))
            ((< n b))
            (let* ((k (+ (* a a a) (* b b b)))
                   (v (gethash k ht)))
              (setf (gethash k ht) (cons (list a b) v)))))
    (maphash
     (lambda (k v)
       (if (>= (length v) m)
           (format t "~D: ~A~%" k v)))
     ht)))

;;; 平衡木
(defun taxi-aatree (n m)
  (let ((tm (aatree:make-treemap #'= #'<)))
    (do ((a 1 (+ a 1)))
        ((< n a))
        (do ((b a (+ b 1)))
            ((< n b))
            (let ((k (+ (* a a a) (* b b b))))
              (aatree:tree-push tm k (list a b)))))
    (aatree:tree-mapc
     tm
     (lambda (k v)
       (when (>= (length v) m)
         (format t "~D: ~A~%" k v))))))

(defun taxi-rbtree (n m)
  (let ((tm (rbtree:make-treemap #'= #'<)))
    (do ((a 1 (+ a 1)))
        ((< n a))
        (do ((b a (+ b 1)))
            ((< n b))
            (let ((k (+ (* a a a) (* b b b))))
              (rbtree:tree-push tm k (list a b)))))
    (rbtree:tree-mapc
     tm
     (lambda (k v)
       (when (>= (length v) m)
         (format t "~D: ~A~%" k v))))))

(defun taxi-splay (n m)
  (let ((tm (splay:make-treemap #'= #'<)))
    (do ((a 1 (+ a 1)))
        ((< n a))
        (do ((b a (+ b 1)))
            ((< n b))
            (let ((k (+ (* a a a) (* b b b))))
              (splay:tree-push tm k (list a b)))))
    (splay:tree-mapc
     tm
     (lambda (k v)
       (when (>= (length v) m)
         (format t "~D: ~A~%" k v))))))

関数 taxi-fast は hash-table を、taxi-aatree, taxi-rbtree, taxi-splay は treemap (AA 木、赤黒木、スプレー木) を使います。どの関数も二重ループで 2 つの整数 a と b を選択します。それから、その 3 乗和を計算して変数 k にセットします。これが連想配列 (ht または tm) のキーになります。そして、k の値に (list a b) を追加していきます。treemap の場合、メソッド tree-push を使うと簡単です。最後に、連想配列を巡回して、整数の組が m 個以上ある値を format で表示します。

それでは実行してみましょう。

* (load "taxi.lsp")
T
* (time (taxi-fast 1000 3))
119824488: ((346 428) (90 492) (11 493))
804360375: ((295 920) (198 927) (15 930))
958595904: ((692 856) (180 984) (22 986))
175959000: ((315 525) (198 552) (70 560))
143604279: ((408 423) (359 460) (111 522))
87539319: ((255 414) (228 423) (167 436))
327763000: ((510 580) (339 661) (300 670))
700314552: ((510 828) (456 846) (334 872))
Evaluation took:
  0.310 seconds of real time
  0.307457 seconds of total run time (0.158936 user, 0.148521 system)
  [ Run times consist of 0.127 seconds GC time, and 0.181 seconds non-GC time. ]
  99.03% CPU
  737,688,225 processor cycles
  81,171,776 bytes consed

NIL
* (time (taxi-aatree 1000 3))
87539319: ((255 414) (228 423) (167 436))
119824488: ((346 428) (90 492) (11 493))
143604279: ((408 423) (359 460) (111 522))
175959000: ((315 525) (198 552) (70 560))
327763000: ((510 580) (339 661) (300 670))
700314552: ((510 828) (456 846) (334 872))
804360375: ((295 920) (198 927) (15 930))
958595904: ((692 856) (180 984) (22 986))
Evaluation took:
  1.490 seconds of real time
  1.493818 seconds of total run time (1.314148 user, 0.179670 system)
  [ Run times consist of 0.158 seconds GC time, and 1.336 seconds non-GC time. ]
  100.27% CPU
  112 lambdas converted
  3,585,347,510 processor cycles
  67,664,400 bytes consed

NIL
* (time (taxi-rbtree 1000 3))
87539319: ((255 414) (228 423) (167 436))
119824488: ((346 428) (90 492) (11 493))
143604279: ((408 423) (359 460) (111 522))
175959000: ((315 525) (198 552) (70 560))
327763000: ((510 580) (339 661) (300 670))
700314552: ((510 828) (456 846) (334 872))
804360375: ((295 920) (198 927) (15 930))
958595904: ((692 856) (180 984) (22 986))
Evaluation took:
  1.029 seconds of real time
  1.037785 seconds of total run time (1.026287 user, 0.011498 system)
  [ Run times consist of 0.116 seconds GC time, and 0.922 seconds non-GC time. ]
  100.87% CPU
  90 lambdas converted
  2,490,725,626 processor cycles
  66,575,296 bytes consed

NIL
* (time (taxi-splay 1000 3))
87539319: ((255 414) (228 423) (167 436))
119824488: ((346 428) (90 492) (11 493))
143604279: ((408 423) (359 460) (111 522))
175959000: ((315 525) (198 552) (70 560))
327763000: ((510 580) (339 661) (300 670))
700314552: ((510 828) (456 846) (334 872))
804360375: ((295 920) (198 927) (15 930))
958595904: ((692 856) (180 984) (22 986))
Evaluation took:
  1.239 seconds of real time
  1.241743 seconds of total run time (1.161554 user, 0.080189 system)
  [ Run times consist of 0.175 seconds GC time, and 1.067 seconds non-GC time. ]
  100.24% CPU
  78 lambdas converted
  2,980,279,258 processor cycles
  90,336,080 bytes consed

NIL

87539319 が Ta(3) になります。実行時間は hash-table の方が速いですね。treemap はちょっと時間がかかりますが、平衡木を使っているので、解を昇順に表示することができます。平衡木の中では赤黒木が一番速く、次がスプレー木、最後が AA 木でした。思っていたよりもスプレー木が速くて驚きました。スプレー木は優れたデータ構造だと思います。

ところで、引数 n の値を増やすと、実行時間だけではなくメモリの使用量も多くなります。treemap は効率化などの工夫はしていないので、メモリ消費量は hash-table よりも多くなると思います。Ta(4) は大きな値になるので、このような力任せのプログラムだと M.Hiroi の実行環境では Ta(3) くらいが限界のように思います。タクシー数を求めるのは難しい問題だと実感しました。

●参考 URL

  1. タクシー数 - Wikipedia
  2. 天才数学者ラマヌジャンのタクシー数の研究, (猫野さん)

Copyright (C) 2023 Makoto Hiroi
All rights reserved.

[ PrevPage | CLOS | NextPage ]