Project Euler Problem 19

問題19: 1900/1/1は月曜日である。では20世紀で月の初めが日曜日になるのは何回あるか。[=>問題文]

問題19の解答

;; Emacs Lisp
(defun problem019 ()
  (let* ((table [31 31 28 31 30 31 30 31 31 30 31 30])
         (day (- (1+ 365) (aref table 0)))
         (count 0))
    (do ((y 1901 (1+ y)))
        ((< 2000 y))
      (let ((is-leap (zerop (% y 4))))
        (do ((m 1 (1+ m)))
            ((< 12 m))
          (incf day (aref table (1- m)))
          (if (and is-leap (= m 3))
              (incf day))
          (if (zerop (% day 7))
                (incf count)))))
    count))

　毎月の１日が基準日(1900/1/1)から何日目かを求め、7で割った余りで曜日を判定しています。「前月の日数を加算する」という処理の都合から、table の第一要素は12月の日数で、1月〜11月の日数がそのあとに続くようにしてあります。
　1901年から1999年の間は 「100で割れる場合云々」「400で割れる場合云々」は無関係で、4年おきのうるう年さえ考えればよく、2000年も「例外の例外年」にあたり通常通りうるう年となります。

既存機能を利用する

　カレンダーライブラリを持つ言語なら面倒な計算をすることなく日付から曜日を取得できます。Emacs Lisp も日時を計算するライブラリ関数 (current-time decode-time encode-time等)がありますが、これらは1970年以降しか扱うことができません。ただこれとは別に M-x calendar で起動する万年カレンダー（日記/スケジューラ）機能が Emacs Lisp には用意されています。これをうまく利用できないでしょうか？

(defun problem019-2 ()
  (let ((buffer))
    (labels
        (;; 準備
         (start ()
                (calendar)
                (setq buffer (current-buffer))
                (buffer-disable-undo buffer))  ;; undo を記録しない
         ;; 後始末
         (finish ()
                 (buffer-enable-undo buffer)
                 (kill-buffer buffer))
         ;; 指定日のバッファ内表示位置
         (get-point (year month day)
                    (calendar-other-month month year)
                    (calendar-goto-date (list month day year))
                    (point)))
      (start)
      (let* ((count 0)
             (monday-point (get-point 1900 1 1)))
        (do ((y 1901 (1+ y)))
            ((< 2000 y))
          (do ((m 1 (1+ m)))
              ((< 12 m))
            (if (= monday-point (get-point y m 2))
                (incf count))))
        (finish)
        count))))

　まず基準日(1900/1/1)から第一週の月曜日のバッファ内表示位置を得ます。そして対象期間のカレンダーを順に表示して「月の２日目が月曜日の位置に表示される」ことを判定します。
　指定日を表示する calendar-goto-date の前に calendar-other-month を評価しているのは対象月をバッファの中央に表示するためです。カレンダーは３か月分並んで表示されるので、これをしないと位置が固定できません。
　buffer-disable-undo はバッファ変更情報の記録を無効化する関数です。頻繁なバッファ表示更新で溜まった undo 記録が limit を越えてしまいワニングが表示されてしまうのでその回避策です。
　バッファ内でカレンダーがパラパラと変わっていくのは見ていてちょっと面白いですが、処理には1分以上かかります。

Project Euler Problem 18, 67

問題18: 数字で構成される三角形を頂点から最下段まで移動する。経路の合計がもっとも大きくなる時の値を求めよ。[=>問題文]

問題18の解答

「経路の総当りをする」という発想がむしろ出ませんでした。

;; Emacs Lisp
(require 'cl)

(defun problem018 (data)
  (let ((prev '(0 0)) ; 一段上のリスト
        (next))       ; 次に処理する段のリスト
    (dolist (next data)
      (let ((tmp '(0)))
        (while (consp next)
          (push
           (+ (pop next)
              (max (first prev) (second prev)))
           tmp)
          (pop prev))
        (setq prev (reverse (cons 0 tmp)))))
    (apply #'max prev)))

(let ((data
       '((75)
         (95 64)
         (17 47 82)
         (18 35 87 10)
         (20 04 82 47 65)
         (19 01 23 75 03 34)
         (88 02 77 73 07 63 67)
         (99 65 04 28 06 16 70 92)
         (41 41 26 56 83 40 80 70 33)
         (41 48 72 33 47 32 37 16 94 29)
         (53 71 44 65 25 43 91 52 97 51 14)
         (70 11 33 28 77 73 17 78 39 68 17 57)
         (91 71 52 38 17 14 91 43 58 50 27 29 48)
         (63 66 04 68 89 53 67 30 73 16 69 87 40 31)
         (04 62 98 27 23 09 70 98 73 93 38 53 60 04 23))))
  (problem018 data))

　上の段の値を次の段の値に加算する。それを上から順に繰り返して雪崩のように三角形を更新してゆきます。三角形内部の任意の数字は、加算候補が右上・左上の 2 つあるので、値の大きいほうを選んで加えて行けば可能な限り大きな数字が最後の段にリストアップされます。最終段の max をとればそれが求める答えとなります。
　各段の端の数字も他と同じロジックで扱えるように　リストの前後に 0 を追加しながら処理をしています。

別解・逆さバージョン

　投稿してから気づきましたが、三角形を逆さまにして処理すればリストに 0 を追加する必要もないし、最後は頂点に収束するから max をとる必要もありませんね。

(defun problem018-2 (data)
  (let* ((rdata (reverse data))
         (prev (car rdata))
         (next))
    (dolist (next (cdr rdata) (car prev))
      (let ((tmp))
        (while (consp next)
          (push
           (+ (pop next)
              (max (first prev) (second prev)))
           tmp)
          (pop prev))
        (setq prev (reverse tmp))))))

問題67の解答

　やることは一緒なので同じ関数を使ってついでにやってしまいます。こちらはファイルで用意されている100段の数値で構成される三角形が対象となります。[=>問題文]

;;problem 67 (Emacs Lisp)
(let* ((str (let ((str)  ;;ファイルからデータ取得
                  (buf (find-file "triangle.txt")))
              (setq str (buffer-string))
              (kill-buffer buf)
              str))
       (lis (mapcar      ;;読み込んだデータを数値のリストへ変換
             #'(lambda (line)
                 (mapcar #'string-to-number
                         (split-string line)))
             (split-string str "\n"))))

  ;; 末尾のゴミを除去
  (setq lis (reverse (cdr (reverse lis))))

  (problem018-2 lis))

　Emacs Lisp ではファイルから直接データを取得することはできず、バッファにデータを読み込んで、そこからデータを切り出すしかないようです。 read 関数を使えば バッファ内の数字列を数値として読み込むことができますがそれだと改行位置がわからなくなるため文字列として取得 (buffer-string) して自前でパースしています。ここに来て初めて Emacs Lisp ならではといえるコードになりました。
　ファイル末尾に改行だけの行があり、これを読み込んで三角形の最後に空リストが入ってしまうため処理する前にゴミ取りが必要です。

Project Euler Problem 17

問題17: 1 から 1000 までの数字をすべて英単語で書くと全部で何文字になるか。[=>問題文]

問題17の解答

　数詞の法則が奇妙な国って意外に多いです。ゼロの発明以前から人類が「数」を使っていた証でしょうか。

;; Emacs Lisp
(require 'cl)
(defun problem017 ()
  (let ((name-table (make-hash-table)))
    (puthash 1 "one" name-table)
    (puthash 2 "two" name-table)
    (puthash 3 "three" name-table)
    (puthash 4 "four" name-table)
    (puthash 5 "five" name-table)
    (puthash 6 "six" name-table)
    (puthash 7 "seven" name-table)
    (puthash 8 "eight" name-table)
    (puthash 9 "nine" name-table)
    (puthash 10 "ten" name-table)
    (puthash 11 "eleven" name-table)
    (puthash 12 "twelve" name-table)
    (puthash 13 "thirteen" name-table)
    (puthash 14 "fourteen" name-table)
    (puthash 15 "fifteen" name-table)
    (puthash 16 "sixteen" name-table)
    (puthash 17 "seventeen" name-table)
    (puthash 18 "eighteen" name-table)
    (puthash 19 "nineteen" name-table)
    (puthash 20 "twenty" name-table)
    (puthash 30 "thirty" name-table)
    (puthash 40 "forty" name-table)
    (puthash 50 "fifty" name-table)
    (puthash 60 "sixty" name-table)
    (puthash 70 "seventy" name-table)
    (puthash 80 "eighty" name-table)
    (puthash 90 "ninety" name-table)
    (puthash 100 "hundred" name-table)
    (puthash 1000 "thousand" name-table)
    (let ((work-vect (make-vector 1001 nil)))
      (labels
          ((set-vect (i lis) (aset work-vect i lis))
           (ref-vect (i) (aref work-vect i))
           (ref-table (i) (list (gethash i name-table))))
        ;; 1 - 20
        (do ((i 1 (1+ i)))
            ((< 20 i))
          (set-vect i (ref-table i)))
        ;; 21 - 99
        (do ((i 20 (+ i 10)))
            ((< 90 i))
          (let ((lis (ref-table i)))
            (set-vect i lis)
            (do ((j 0 (1+ j)))
                ((< 9 j))
              (set-vect (+ i j) (append lis (ref-vect j))))))
        ;; 100 - 999
        (do ((i 1 (1+ i)) (name100 (ref-table 100)))
            ((< 9 i))
          (let* ((j (* i 100))
                 (lis (append (ref-vect i) name100)))
            (set-vect j lis)
            (do ((k 1 (1+ k)))
                ((< 99 k))
              (set-vect (+ j k)
                        (append lis (cons "and" (ref-vect k)))))))
        ;; 1000
        (set-vect 1000 (append (ref-vect 1) (ref-table 1000)))

        ;; 全文字列を連結して長さを返す
        (length (apply #'concat
                       (mapcar #'(lambda (lis)
                                   (apply #'concat lis))
                               work-vect)))))))
(problem017)

　たとえば 256 なら ("two" "hundred" "and" "fifty" "six") というリストを作って work-vect のインデックス 256 の位置に保存します。そうやって 1000 までリストを作成して、最後に全文字列を連結して長さを取得しています。必要なのは文字列長だけなので、実際に文字列を作る必要は無く 文字列の長さだけ記録していけばいいんですが、デバッグすることを考えるとそれはちょっと避けたいですね。
　100 以降の下二桁は既に作成済みの 99 までの数詞を参照することで少しだけ処理を省略できます。
　あ、そうそう。今回ループに do を使ってます。ちょっと括弧の対応がややこしくて今まで敬遠してきたけど一念発起。ループ式の中で let を使わずレキシカル変数が作れるってのは便利ですね。

Project Euler Problem 16

問題16: 2^1000 の各桁の数字を合計した値を求めよ。[=>問題文]

bigint ライブラリの改良

　かなりヘビーな計算量になりそうなので、問題13で作った多倍長整数ライブラリの改良と、機能追加をします。改良点は２つ。まず、多倍長数をあらわすリストの要素を0〜9の1桁から、0〜9999の4桁に変更。計算速度が4倍になります。もう一つの改良は加算関数 bigint+ の効率化。以前のものは桁数の大きく異なる2数の和計算でとても非効率でした。
　追加機能は、乗算と累乗計算。ここまで作れば問題は解けたも同然です。

(defun bigint-expt (n m)
  "累乗計算。n:bigint, m:0以上の整数"
  (let ((result '(1)) (i 0) (cnt 0))
    (while (zerop (/ m 2))
      (incf cnt)
      (setq m (/ m 2)))
    (while (<= 0 (decf m))
      (setq result (bigint* result n)))
    (while (<= 0 (decf cnt))
      (setq result (bigint* result result)))
    result))

　累乗計算は、最初に指数 m を割れなくなるまで 2 で割ってその回数をカウントしてます。残りの m の回数だけ n を繰り返し積算したその結果を、最初に2で割った回数だけ二乗すれば答えがでます。これは掛け算の回数を減らす工夫ですが、パフォーマンスアップはさほどでもありませんでした。
　※多倍長演算ライブラリの全コードはこちらへ移動しました。

問題16の解答

(defun problem016 ()
  (let ((str (bigint-to-string
              (bigint-expt
               (string-to-bigint "2")
               1000))))
    (apply #'+
           (mapcar #'(lambda (c) (- c ?0)) str))))

　処理にかかった時間は 0.035 秒。なかなかのもんじゃないでしょうか。

Project Euler Problem 15

問題15: 20x20のマトリクスにおいて左上端から右下端へ格子を通って到達する経路は何通りあるか。ただし、左から右、上から下へしか移動できないものとする。[=>問題文]

図を描いてみる

　簡単な図を書くだけですぐに仕組みが分かります。

　右から左、下から上へは戻れないため上辺に位置する格子点は常に左からの１通りしか到達手段がありません。同様に左辺の点もみんな１。内部の点に関しては到達経路が上からと左から２通りになります。したがってある点への到達経路の数は、その上の点への経路数と左の点への経路数の和で求められます。いくつかの格子点に経路の数を記入していくと、これは「パスカルの三角形」と同じことをやってることに気づきます。図は左に45度傾いたパスカルの三角形そのものです。
　ゴールは対角線上にあるので割り当てられる値は、三角形の行の要素数が奇数となる段の中央の値となります。つまりこの問題は次のように読みかえることが出来ます。
　パスカルの三角形の 2n 段目の中央の値を答えよ（ただし三角形の頂点は 0段目とする）

問題15の解答

(defun problem015 (N)
  (let ((prev '(1.0))
        (next nil)
        (i 0))
    (dotimes (i (* 2 N))
      (setq next '(1.0))
      (while (< 1 (length prev))
        (push (+ (pop prev) (car prev)) next))
      (setq prev (cons 1.0 next)))
    (nth (/ (length prev) 2) prev)))

(problem015 20)

　これまた 28bitでは収まらないので float を使っています。

Project Euler Problem 14

問題14: 100万未満の数の中でコラッツ数列がもっとも長くなる数を求める。[=>問題文]
　最初の10個を書いてみると。

  1 (1) : 1
  2 (2) : 2 - 1 
  3 (8) : 3 - 10 - 5 - 16 - 8 - 4 - 2 - 1
  4 (3) : 4 - 2 - 1
  5 (6) : 5 - 16 - 8 - 4 - 2 - 1
  6 (9) : 6 - 3 - 10 - 5 - 16 - 8 - 4 - 2 - 1
  7 (17): 7 - 22 - 11 - 34 - 17 - 52 - 26 - 13 - 40 - 20 - 10 - 5 - 16 - 8 - 4 - 2 - 1
  8 (4) : 8 - 4 - 2 - 1
  9 (20): 9 - 28 - 14 - 7 - 22 - 11 - 34 - 17 - 52 - 26 - 13 - 40 - 20 - 10 - 5 - 16 - 8 - 4 - 2 - 1
 10 (7) : 10 - 5 - 16 - 8 - 4 - 2 - 1

　「9」のコラッツ列は、4番目で「7」が現れます。数列の長さを「1」から順にチェックするならば、「9」の時点では「7」は調べ終わっているはずです。したがって「9」から「7」へ到達するまでの長さ 3 と既にわかっている「7」のコラッツ列長 17 を足した 20 が「9」の長さとわかります。17回の計算が省略できたことになります。このように過去のデータを記録しておき計算量を減らすことを「メモ化」というそうです。

問題14の解答

vectoer によるメモ化

(defun problem014 (N)
  (let ((memo (make-vector (1+ N) 0))
        (ans 1) (max-cnt 0) (n 1))
    (aset memo 1 1)
    (while (< (incf n) N)
      (if (evenp n)
          (aset memo n (1+ (aref memo (/ n 2))))
        (let ((cnt 0) (m (float n)))
          (while (>= m n)
            (incf cnt)
            (setq m (if (zerop (mod m 2))
                        (/ m 2)
                      (1+ (* 3 m)))))
          (incf cnt (aref memo (truncate m)))
          (aset memo n cnt)
          (when (> cnt max-cnt)
            (setq max-cnt cnt)
            (setq ans n)))))
    ans))

(problem014 1000000)

　小さいほうから順に100万まで n を増加させコラッツ列を調べます。数列を計算して n より小さな値が現れたら、memo を参照し、そこまでの cnt と加算し列の長さを求めます。
　偶数は最初にいきなり半分にされ、必ず自身より小さな値になるので memo の インデックス (/ n 2)に保存されている値に 1 を加えるだけで計算がすみます。
　100万までのコラッツ列に現れる数の最大は 56991483520 で、これは 28bitどころか32bitで表現できる整数の範囲も超えてしまいます。そのためfloat型で数列作成を行っています。ただし、vector の インデックスには float 型が使えないため、memo を参照するさいに truncate をかける必要があります。

hash table 版

　ところで「3」のコラッツ列(3-10-5-16-8-4-2-1)を調べると「10」のコラッツ列の長さは 7,「5」の長さは 6、「16」の長さは 5 と、自分(n)より大きなコラッツ数列も調べていることになります。はじめに書いたプログラムでは、n より大きい数についての情報は保存しないので、これらの情報は捨てられ、「5」や「10」を調べる時、もう一度計算しなおしています。これらも記録したらもっと早くなるでしょうか？

(defun problem014-2 (N)
  (let ((memo (make-hash-table))
        (ans 1) (max-cnt 0) (n 1))
    (puthash 1.0 1 memo)
    (while (< (incf n) N)
      (let ((m (float n)))
        (when (null (gethash m memo))
          (if (evenp n)
              (puthash m
                       (1+ (gethash (/ m 2) memo))
                       memo))
          (let ((pool nil) (cnt 0))
            (while (null (setq cnt (gethash m memo)))
              (push m pool)
              (setq m (if (zerop (mod m 2))
                          (/ m 2)
                        (1+ (* 3 m)))))
            (let ((x 0))
              (dolist (x pool)
                (incf cnt)
                (puthash x cnt memo)))
            (when (> cnt max-cnt)
              (setq max-cnt cnt)
              (setq ans n))))))
    ans))

　n より大きい数といってもどれだけ大きい値が出るのか事前にはわかりません。固定サイズのvectorでは扱えません。そこで vector の代わりに　ハッシュテーブルを使うことにします。ハッシュに記録が無い数値が続く間は、リスト pool に値を push してゆきます。ハッシュに記録が見つかったらそのカウント値を取得し、while を終了。次の dolist で pool に保存されている値を順に取り出し、cnt を増やしながら ハッシュに保存します。これで n より大きな数のコラッツ長が保存されてゆきます。
　実行して処理にかかる時間を計測すると・・・

Function       Avg time per call
=============  =================
problem014     5.328500   (vector)
problem014-2   9.031000   (hash table)

　むしろ、大きく速度低下してしまいました。ループ回数は85%ほどに減っているけど、それ以上に ハッシュデータの作成と参照にかかるコストが上回ってしまっているのでしょうね。

vector + hash table 版

　では、100万以下の数は vector に保存し、100万を越えたらhash tableに保存するという方法ではどうでしょうか？ 両者のいいとこ取りで速くなるかも？

(defun problem014-3 (N)
  (let ((memo-vect (make-vector (1+ N) 0))
        (memo-hash (make-hash-table))
        (ans 1) (max-cnt 0) (n 1))
    (aset memo-vect 1 1)
    (while (< (incf n) N)
      (if (zerop (aref memo-vect n))
          (if (evenp n)
              (aset memo-vect n
                    (1+ (aref memo-vect (/ n 2))))
            (let ((m (float n))
                  (pool nil) (cnt-from-hash nil))
              (while (and
                      (>= m n)
                      (null (setq cnt-from-hash
                                  (gethash m memo-hash))))
                (push m pool)
                (setq m (if (zerop (mod m 2))
                            (/ m 2)
                          (1+ (* 3 m)))))
              (let ((cnt (if (null cnt-from-hash)
                             (aref memo-vect (truncate m))
                           (gethash m memo-hash)))
                    (x 0))
                (dolist (x pool)
                  (incf cnt)
                  (if (< N x)
                      (puthash x cnt memo-hash)
                    (aset memo-vect (truncate x) cnt)))
                (when (> cnt max-cnt)
                  (setq max-cnt cnt)
                  (setq ans n)))))))
    ans))

　プロファイルすると・・・

Function       Avg time per call
=============  =================
problem014     5.328500   (vector)
problem014-3   5.680000   (vector + hash table)
problem014-2   9.031000   (hash table)

　だいぶ肉薄しましたが結局　vector だけの方式が一番速いようです。
　hashtable は生成時にキャパシティを指定して、メモリ再確保のコストをなくすオプションがありますが　:size 200万 を指定しても高速化は微々たる物でした。

Project Euler Problem 13

問題13: 提示された50桁の数字100個の総和の上位10桁を求める。[=>問題文]
　とうとうこの手の来てしまいましたね。多倍長演算ってライブラリや言語サポートがないとめんどくさいんだよな・・・。足し算だけなので最低限機能でサクっといきましょう。

問題13の解答

(require 'cl)

(defun string-to-bigint (str)
  (reverse (mapcar #'(lambda (c) (- c ?0))
                   (string-to-list str))))

(defun bigint-to-string (n)
  (apply #'concat
         (mapcar #'number-to-string (reverse n))))

(defun bigint+ (n1 n2)
  (let* ((len-n1 (length n1))
         (len-n2 (length n2))
         (m1 (if (< len-n1 len-n2) n2 n1))
         (m2 (append
              (if (eq m1 n1) n2 n1)
              (make-list (abs (- len-n1 len-n2)) 0)))
         (carry 0)
         (lis (mapcar*
               #'(lambda (a b)
                   (let ((n (+ a b carry)))
                     (setq carry (/ n 10))
                     (% n 10)))
               m1 m2)))
    (if (zerop carry)
        lis
      (append lis (list carry)))))

(defun problem013 (numbers)
  (let ((n nil)
        (sum (string-to-bigint (car numbers))))
    (dolist (n (cdr numbers))
      (setq sum (bigint+ sum (string-to-bigint n))))
    (subseq (bigint-to-string sum) 0 10)))

(defvar *numbers*)
(setq *numbers*
      '("37107287533902102798797998220837590246510135740250"
        "46376937677490009712648124896970078050417018260538"
        "74324986199524741059474233309513058123726617309629"
        "91942213363574161572522430563301811072406154908250"
        ;; 中略
        "77158542502016545090413245809786882778948721859617"
        "72107838435069186155435662884062257473692284509516"
        "20849603980134001723930671666823555245252804609722"
        "53503534226472524250874054075591789781264330331690"))

(problem013 *numbers*)

　多倍長整数と言っても、桁を一桁づつばらして逆順にした整数リストです。足し算でポイントとなるのは繰り上がりをどう処理するか。今回作った加算関数 bigint+ では、mapcar で 2つのリストを回しつつ、外側に置いた carry 変数で繰り上がりを受け渡しています。
　Emacs Lisp の mapcar はリストを一つしかとれない劣化版ですが、cl パッケージに mapcar* という Common Lisp 仕様の高階関数が用意されていました。
　mapcar* に渡すリストは同じ長さでなければ切り捨てられるため、事前に桁の少ないほうのリストを 0 埋めしてしています。

再帰版

　作ってみましたが相変わらず再帰が深いと怒られます。

(setq max-lisp-eval-depth 3000)
(setq max-specpdl-size 3000)

(defun bigint+ (n1 n2)
  (let* ((len-n1 (length n1))
         (len-n2 (length n2))
         (m1 (if (< len-n1 len-n2) n2 n1))
         (m2 (append
              (if (eq m1 n1) n2 n1)
              (make-list (abs (- len-n1 len-n2)) 0))))
    (labels
        ((repeat (m1 m2 carry)
                 (if (null m1)
                     (if (zerop carry)
                         nil
                       (list carry))
                   (let ((n (+ (car m1) (car m2) carry)))
                     (cons (% n 10)
                           (repeat (cdr m1) (cdr m2)
                                   (/ n 10)))))))
      (repeat m1 m2 0))))

プログラム修正

　文字列をリストにする時 reverse をかけ逆順になるようにしました。加算では下の桁から処理するのでこのほうが無駄がありません（修正前は加算前と後に２回 reverseをかけていた）。減算、積算もおそらくこの形式のほうがやりやすいんじゃないと思います。除算は･･････あんまり考えたくない･･･。