4 posts tagged “scheme”

1/21 の NHK スペシャルを見ていたら，ハイウェイ沿いに掲げられた Google の謎の求人広告の話題が挙がっていた．

{ first 10-digit prime found in consecutive digits of e }.com

あれー，確か昔に聞いたことがあったけど，当時はハナから「どうせ解けないだろうなー」と諦めていたような気がする．しかし，いま改めて向き合ってみると，e が分かっていれば大したことないと分かる．

10 桁というと 32bit で収まりきらないのでチト気を付ける必要があるが，gauche は多倍長整数を扱えるのでそちらに任せる．素数判定は，まぁその都度行っても問題ないだろう．ダメだったらそれはそれで．

年末は SICP を読んでいたので，素数判定はフェルマーテストで行ってみる．合格した場合は念のためシラミ潰す．

(define (square x) (* x x))

(define (power base exp m)
  (cond ((zero? exp) 1)
        ((even? exp) (modulo (square (power base (/ exp 2) m)) m))
        (else        (modulo (* base (power base (- exp 1) m)) m))))

(define (fermat-test x)
  (if (= 2 x) #t
      (= 2 (power 2 x x))))

(define (prime? x)
  (if (fermat-test x)
    (if (< x 2)
      #f
      (let loop ((n 2) (h (sqrt x)))
        (cond ((> n h) #t)
              ((zero? (remainder x n)) #f)
              (else (loop (if (< 2 n) (+ 2 n) 3) h)))))
    #f))

末尾再帰の if がどうも美しくないような気がするけどソコは触れない約束で．

で，問題の e だけど，ググれば何万桁オーダーで載っているサイトも見つかる．これをテキストファイルにして頭から読み込みながら計算してゆけばいずれ答えは出る．やや bruteforcely だが，他に (数学的にエレガントな?) 解法は思いつかない．

(define low (expt 10 9))
(define high (expt 10 10))

(display
 (let loop ((x 0))
   (if (and (>= x low) (prime? x))
       x
     (loop (remainder (+ (- (char->integer (read-char)) 48) (* 10 x)) high)))))

とりあえずコレを q.scm という名前で保存して

cat e.txt | gosh q.scm

と実行すると 7427466391 が求められる．小数点以下 99 桁目だ．意外と浅い場所にあったが，ググってみると "Google" というキーワードと共に沢山のページが引っかかるので，おそらくコレで正解なのだろう (http://7427466391.com 自体はもう存在しないようだ)．

しかし，これでいいのか? 既知の e を用いて求めるのはアリなのか!?

e 自体は lim(n → ∞)(1 + 1/n)^n で求められる (ググった)．しかし n が大きくなると結果が 2.71 にも届かないうちに NaN ってしまうので無理がある．ここはひとつ，マクローリン展開された Σ(n → ∞)(1/n!) を使ってジワジワ攻めてゆくのが良いはず (ググった)．というかマクローリン展開なつかしす(w

gauche の多倍長有理数を使えば，メモリの許す限りの精度で e を算出できる．

(define precision 120)

(define e
  (let loop ((n 1) (d 1) (x 1))
    (if (> (- precision) (/ (log d) (log 10)))
      x
      (loop (+ 1 n) (/ d n) (+ x (/ d n))))))

(let loop ((n 0) (e e))
  (if (> n precision)
    #t
    (let ((x (floor e)))
      (display x)
      (loop (+ 1 n) (* 10 (- e x))))))

さて，これで求められた結果を先ほどの q.scm に食わせて…って，おい，何か違うだろう．e.txt を作る手間が増えただけじゃないか．そうじゃなくて，もっと，こう，なんつーか，エレガンスが必要じゃないか? これだと予め決められた桁数の範囲内でしか試行できない．できれば，「何桁目で答えが出るのか分からない状態で試行を続ける」というアプローチでいきたい．

n = 0 から出発して，前回よりも非常に小さい値 (階乗の逆数) を加算し続けてゆくので e はどんどん収束してゆくはず．例えば，n = 3 において 0.16 を加算しても 1 の位への桁上がりは起きないので一桁目の "2" は確定する．また，n = 5 において 0.04 を加算しても 0.1 の位への桁上がりは起きないので二桁目の "7" は確定する．

n	e
0	1.00000000000000
1	2.00000000000000
2	2.50000000000000
3	2.66666666666667
4	2.70833333333333
5	2.71666666666667
6	2.71805555555556
7	2.71825396825397
8	2.71827876984127
9	2.71828152557319
10	2.71828180114638
11	2.71828182619849
12	2.71828182828617

こうやって，確定した値を次々と取り出してゆき，10 桁に達したら素数判定を行い，10 桁目を捨てる．この作業を素数が見つかるまで延々と繰り返す．コードはあまり綺麗ではないが，こんな感じ．

(let loop ((n 1) (d 1) (x 0) (y 0))
  (let ((d (/ d n)) (z (+ x d)))
    (cond ((= (floor (* 1000 z)) (floor (* 1000 y)))
           (let ((fz (floor z)))
             (display fz)
             (flush)
             (loop (+ n 1) (* 10 d) (* 10 (- z fz)) x)))
          (else (loop (+ n 1) d z x)))))

1,000 を掛けているのは，実は適当．頭 4 桁くらい一致してればまぁ一桁目は確定でいいかな～なんて．あはは．2 桁目までしか見ないと 960 桁目で実際の値 (e.txt) と食い違ったので，余裕を持ってみた次第 (こんなんでいいのか?

それはさて置き，上記ソースを e.scm という名前で保存して

gosh e.scm | gosh q.scm

とすると，見事 7427466391 が求められる．いやっほう，エレガント!! (微妙)．

本当は，e.scm の方をコルーチンにして「呼ばれるたびに次の e の桁を返す関数」みたいにすると格好いいんだろうけど，call/cc が未だに理解できていないのでまた後日．

図書館で SICP を借りた．「正月休みに読破するか～」とか考えてたら，甘すぎた (ぶはは)．まだ 2 章前半．Church 数かわいいよ Church 数!!

うむぅ，図書館だと 2 週間しか借りられないんだよなぁ．他区から取り寄せてもらった本だし，延滞するのは憚られる．しゃーないから買うかぁ…．

なんだろう，scheme 始めて SICP に転ぶというのは凄く分かりやすいお決まりパターンのような気がするけど，まぁそれはさて置き，読み進めるにつれ，ここ 10 年間で如何に無駄な時間を過ごしてきたかを身に染みて体感できて困る．こーゆーのは学生のうちに済ませておきたかったなぁ (てゆーか元々講義用っぽいし)．会社で読み出したら仕事にならんだろうし :-(．そうだ，会社サボって読めばいいんだ! (ばか

scheme の処理系 guile, gauche を試してみた．

どちらもインタプリタフロントエンドで readline が効かないようなので何とかする．

~/.guile を作る．

(use-modules (ice-9 readline))
(activate-readline)

Gauche-readline を入れて gosh-rl を起動する．

guile, gauche 両者とも，シンボルが case-sensitive なのは良いと思う．case-insensitive だと

(equal? "a" (symbol->string (string->symbol "a")))

が成り立たないかもしれなくて気持ち悪い．

ところで，guile に何故か nil という自己評価フォームが予め定義されている．

nil => nil
'nil => nil
(eqv? nil 'nil) => #t
(eqv? nil '()) => #f
(eqv? nil #f) => #f

なんだこれ…．何のためにあるんだろう．R4RS には nil という単語は見当たらないし，R5RS にも注釈として

Note: Programmers accustomed to other dialects of Lisp should be aware that Scheme distinguishes both #f and the empty list from the symbol nil.

という一文があるだけで，これでは nil の在り方を定めているようには読めない．てゆーか guile でも

(equal? nil '()) => #f

だし．

昔から .emacs や .xyzzy を「コピペじゃなくて，理解して書きたい!!」とは思ってたけど，そこら辺の入門サイトをチラ見しただけだと (+ 1 2) は理解できてもシンボルや特殊フォームが分からずに毎回いつの間にかフェードアウトしていたような気がする．

で，今回．
何がトリガだったか忘れたけど (/.j で珍しく LISP 系の話題が幾つか挙がったから?)，ここ一週間くらい失速せずに没頭して，いつもつまずいていた障害は知らないうちに乗り越えていた．

「評価させないためには quote する」って書いてあるのにココ quote してないじゃん!! と俺を憤慨させていたのは特殊フォームというやつだと分かった．「setq の第一引数が quote されていないことにお気づきでしょうか?」というお約束導入で特殊フォームを説明しているサイトをよく見かけたけど，setq の前に set を説明して貰えれば混乱も少なくて済むような気がする．

で，既存のコードなんかを読んでて (add-hook 'ed::*java-mode-hook* ...) なんてのに出くわしちゃうと，コロンやアスタリスクを識別子の一部だと思えずに「こんな syntax があるのか!? どう解釈されるんだ!?」とか脱線しちゃって「LISP 仕様」とかでググっちゃったりして．「さわりだけ知りたいだけなのに，入門時から言語仕様なんか参照しないだろふつー」とか思いつつも，結局は R5RS あたりを読んで「あぁ，ただの識別子なのか」と理解した気がする．

#'(lambda ...) とかもヤバい．未知との遭遇だ (未知なのは当たり前か…)．結局コレも function の省略形だと知ったが，なぜか function だろうが quote だろうが

(defvar f1 '(lambda (x) (+ x 1)))
(defvar f2 #'(lambda (x) (- x 1)))

みたいに区別せずに書けてしまう．なんだよー (違いが無いわけではないらしい)．

で，色々と惰性でウェブを彷徨っているうちに Common Lisp とか Emacs Lisp とか，延いては scheme なんてのも出てきて，define だの defun だの defvar だの set だの set! だの setq だの setf だのと混乱しているうちに，何となく

とか思ってしまい，今に至る．