Island Life

2011/09/29

無限素数列、続き

「よりClojureらしい素数列」をGaucheで」のエントリでは、Gaucheのlazy sequenceの例を簡単に示すのに良いと思って単純な素数列生成のClojureコードをそのまま移植したのだけれど、何やら素数列生成が盛り上がっているようなのでもうちょい続けてみる。なおこのエントリのコードも、実行にはGaucheの開発版HEADが必要。

前エントリでは、無限の自然数列から、素数であるものだけを抽出して素数列を求めていたが、nが素数であるかの判定は極めてナイーブで、2から(+ 1 (/ n 2))までの自然数でnを割ってみる、というものだった。

(define (prime? n)
  (not (any zero? (lmap (pa$ mod n) (lrange 2 (+ 1 (/ n 2)))))))

でもこれはかなり無駄なことをしている。2と3で割れなければ6でも割れるわけないのだから、(ceiling √n) 以下の素数についてだけ試せばいい。

今、仮に、*primes* に無限素数列が既に束縛されているとしよう。 すると、prime? の定義は次のように書ける。

(define (prime? n)
  (not (any zero? (lmap (pa$ mod n) (take-while (^k (<= (* k k) n)) *primes*)))))

ある数k (kは3以上の奇数とする) が与えられた時、そこから先の素数列は次のように書ける。

(define (primes-from k)
  (if (prime? k)
    (lcons k (primes-from (+ k 2)))
    (primes-from (+ k 2))))

(+ k 2) してるのは、3以上の素数はすべて奇数なので偶数を試す必要はないからだ。

ここで、lcons はlazyなcons。Clojureなら (lazy-seq (cons n (recur (+ k 2)))) などと書くところだが、Gaucheの場合lconsというプリミティブでlazyなペアを作る。ちょっと注意が必要なのは、lconsはcdr部分についてだけ評価を遅延するということ。つまり (cons k (delay (loop (+ k 2)))) に近い。ただしdelayと違って、lconsで作られたlazyなペアのcdr部は必要になったら自動的に評価されるので、いちいちforce する必要はない。

primes-fromを使えば、*primes*を定義できる。primes-from の中でprime?を使ってて、その中で*primes*を参照しているから循環定義になりそうなものだが、 primes-from でk以上の素数を計算する時には、(ceiling √k) 以下の素数しか参照しない。つまりある程度kが大きければ、既に計算済みの部分しか使わないようにできる。

(define *primes* (list* 2 3 (lcons 5 (primes-from 7))))

(primes-from 7) を計算するには、*primes* の頭から2と3だけを試せば良いので、2と3をあらかじめリストの頭にくっつけて置けば良いというわけだ。全部まとめるとこうなる。

(define (prime? n)
  (not (any zero? (lmap (pa$ mod n) (take-while (^k (<= (* k k) n)) *primes*)))))

(define (primes-from k)
  (if (prime? k)
    (lcons k (primes-from (+ k 2)))
    (primes-from (+ k 2))))

(define *primes* (list* 2 3 (lcons 5 (primes-from 7))))

(lazyな計算に慣れてる人は、(cons 2 (lcons 3 (primes-from 5))) でも良いんじゃないかと思うかもしれない。 だが、前エントリ でちょっと触れたように、Gaucheのlazy sequenceは一つ先の要素を計算してしまう。(primes-from 5) を計算するには2と3を試す必要があり、lazy pairのcarである3を取り出そうとした瞬間、Gaucheはそのcdrの要素を一つ計算しようとしてしまうので、計算が終わらなくなる。「一つ先読み」が落とし穴になりそうだというのは、こういうケースがあるから。)

なお、Clojureのcontribにも無限素数列のコードがあるが(clojure.contrib.lazy-seqs/primes)、そちらはもうちょい賢いことをやっている。

(defvar primes
  (concat
   [2 3 5 7]
   (lazy-seq
    (let [primes-from
          (fn primes-from [n [f & r]]
              (if (some #(zero? (rem n %))
                        (take-while #(<= (* % %) n) primes))
                (recur (+ n f) r)
                (lazy-seq (cons n (primes-from (+ n f) r)))))
          wheel (cycle [2 4 2 4 6 2 6 4 2 4 6 6 2 6  4  2
                        6 4 6 8 4 2 4 2 4 8 6 4 6 2  4  6
                        2 6 6 4 2 4 6 2 6 4 2 4 2 10 2 10])]
      (primes-from 11 wheel))))
  "Lazy sequence of all the prime numbers.")

次の素数を試すのに、単純に2を足すのではなく、何やら妙な循環数列wheelを使っている。 これはざっくり言えば「ある素数から次の素数までの差は、この数の倍数になる」という列。 詳しくは Recursion on the cycle of gaps を参照。ここで使ってるのはG(7#)なので[2 3 5 7]をあらかじめ計算済みとしている。

さて最後に、「Gaucheらしい」無限素数列を。普通のlazyな言語で無限リストを作るには、余再帰と呼ばれる、(lazy-cons <要素> <再帰呼び出し>) というパターンが使われる。上に示したコードもこのパターンになっている。ただ、このパターンでは <再帰呼び出し> の部分を後で計算するために、その計算環境を保存しておかねばならない。計算環境を保存するのは、Gaucheではクロージャを作るのと同じコストがかかる (もともとlazyな言語なら処理系がもっと効率よく処理するだろうけど)。リストの1要素ごとに新しいクロージャを作るのは、性能上はあんまり嬉しくない。

そこで、Gaucheのlazy sequenceでは、最下層にジェネレータを置いた。ジェネレータは引数無しの手続きで、呼ばれる度に新しい値を返す。値が尽きたら#<eof>を返す。 例えばread-charは引数無しで呼べば次々と#<eof>に達するまで文字を返すので、ジェネレータとして使える。

プリミティブgenerator->lseqはジェネレータを取って、生成される値のlazy sequenceを返す。(generator->lseq read-char) とすれば、文字のlazy sequenceが返ってくるわけだ。ジェネレータの前にあらかじめわかっている要素を付け足すこともできる。(generator->lseq #\a read-char) とすると、最初に#\aが、その後に現在の入力ポートから読んだ文字がくるlazy sequenceとなる。

ジェネレータは作り方によってはnon-consing、つまり全くアロケーションを行わないように書けるので、一要素生成するごとにクロージャを作るのに比べるとかなり速い。

こちらがジェネレータベースの無限素数列。gen-primes-above がジェネレータを作って返す。

(define (gen-primes-above last)
  (^() (set! last (next-prime (+ last 2))) last))

(define *primes* (generator->lseq 2 3 5 (gen-primes-above 5)))

(define (next-prime last)
  (let loop ([ps *primes*])
    (let1 p (car ps)
      (cond [(> (* p p) last) last]
            [(zero? (modulo last p)) (next-prime (+ last 2))]
            [else (loop (cdr ps))]))))

前エントリのコード、lconsを使った余再帰によるコード、ジェネレータベースのコードの ベンチマークを示しておこう。無限素数列から最初の5000個の素数を取ってくるのに要した 時間 (Core2 quad 2.4GHz):

前エントリ 85.7s
余再帰 0.79s
ジェネレータ 0.046s

なお、エラトステネスの無限の篩 - 水底で思うことでは「無限リストに無限個の無限リストを適用する」という鮮やかな技を使ってる。時間ができたらこれもGaucheでやってみたい。

(追記2011/11/20 13:07:07 UTC: ジェネレータから遅延シーケンスを作るプリミティブは当初lseqという名前だったが、何を引数とするかわかりにくいのでgenerator->lseqと改名した。本文もそれに合わせて修正してある。lseqという名前は<sequence>から遅延シーケンスを作るコンストラクタにするかもしれん。)

Tags: Gauche, Clojure, Programming

2011/09/28

ピアノレッスン16回目/新エージェント

  • 基礎: スケール、アルペジオ MM=144
    • スケールは上昇で走り、下降で遅れがち
    • B♭ majorのアルペジオはまだ鬼門。成功率50%くらいだなあ。
  • Alkan: Op39-2 "En Rythme Molossique"
    • 半分の速度でメトロノームかけてとりあえず通した。
    • まあ、譜面どおりに弾いただけではひどく退屈な曲になるねぇ。
    • 一般的な注意として、タンタタタタの「タタタタ」の4連続スタッカートは同じようにひいてはいけない。基本クレッシェンド気味に。
    • 腕を使うスタッカートと手首を使うスタッカートがある。この曲の速度だと手首を使うべき。
    • 左に旋律が出てくるとこはせっかくだから左を歌う。今は右が強すぎ。
    • あとは、全体の組み立てを考えながら、歌うところ、息を抜くところ、など細かいニュアンスをつけてく。これについては来週から。セクションで切って考えてゆこう。

レッスンの後は、エージェントのインタビューへ。 面接してモノローグひとつ演じて、無事契約完了。 スモールビジネスの継続性でちょっと書いたけど、数年前から契約してたエージェントがマネジメント変更でちょっとごたごたしたので、 新しいエージェントを探さなくちゃ、と思いつつぐずぐず引き伸ばしていたのだった。 これでやっとタスクがひとつ片付いた。

デモリールもできれば作って欲しいと言われたんだけど、 映像になってるやつで台詞ががっつりあるシーンってあまり無いんだよな。

Tags: Piano, 芝居

2011/09/27

「よりClojureらしい素数列」をGaucheで

よりClojureらしい素数列 - OGINO Masanori@はてな

(defn prime?
  "Returns true if x is prime number, false otherwise."
  [x]
  (not-any? zero?
            (map #(rem x %)
                 (range 2 (inc (/ x 2))))))

(defn prime-numbers
  "Returns a lazy seq of prime numbers."
  []
  (filter prime? (drop 2 (range))))

効率的にベストではないけどとてもわかりやすいコードと解説だと思った。

で、無限素数列を生成して呼び出し側が必要なとこだけ取れるようにできるのが lazy sequenceの良いところなんだけど、現在のGaucheの開発版HEADでも これができるようになっている。

まだ試しにいろいろ書いてみて実装を調整してる段階なので、 APIとか揃ってないし確定もしてないんだけど。 上のコードをそのまま移植すると今のところはこんな感じで書ける:

(use srfi-1)
(use gauche.lazy)

(define (prime? n)
  (not (any zero? (lmap (pa$ mod n) (lrange 2 (+ 1 (/ n 2)))))))

(define (prime-numbers)
  (lfilter prime? (lrange 2)))

prime-numbersの戻り値は無限リストなのでREPLで表示させようとすると終わらなくなる。 *print-length*みたいのをそのうちつけるつもりだけど。 結果を表示するには、srfi-1のtakeで頭から必要なだけ取るとか (clojureのtakeと 引数の順序が違うことに注意):

gosh> (take (prime-numbers) 10)
(2 3 5 7 11 13 17 19 23 29)

100以下の素数だけ取るとか:

gosh> (take-while (pa$ > 100) (prime-numbers))
(2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97)

100番目の素数だけ取るとか:

gosh> (~ (prime-numbers) 99)
541

ストリームライブラリ (refj:util.stream) との違いは、 リストを期待している関数にlazy sequenceをそのまま渡せるところ。 なのでeagerに結果が欲しくなった時点で普通のリスト処理関数を使える。 carcdrももちろん使える。というより Schemeレベルではlazyなシーケンスか普通のペアかという区別はつかない (副作用や停止しないことを以って観測可能ではあるけれど)。

現在の実装はちょっと注意が必要で、 「必要とされる要素のひとつ先まで計算する」って性質がある。 これは性能上の要件としてわりとクリティカルで、今回の実装の肝でもあるんだけれど、 lazyなアルゴリズムを実装するのにshow stopperにならないかどうか、 まだ確信が持てない。 もし何らかのshow stopperが見つかったら、このlazy sequenceの機能を すぱっと諦めるかもしれない。

Tags: Gauche, Clojure

2011/09/26

It is the tale, not he who tells it

この話題は根強いなあ。

『作者の気持ちを読み取る』←そんなエスパーな!!

"「火垂るの墓」原作者、野坂昭如の娘の国語の授業で父の作品が扱われた。問題に「この時の著者の心境を答えよ」というものがあったので、娘は家に帰ってから父に訪ねた。「その時どんな気持ちだったの?」「締め切りに追われて必死だった」翌日のテストで答えにそう書いた娘は×をもらった。"

火垂るの墓のこのエピソードは何回か耳にしたことはあるけど原典を知らないので真偽はよくわからない。「作者の意図なんてわからんよ」という話なら、『ある物語の顛末』がおもしろい。シャーリー・ジャクスンが『くじ』の反響について語ったもの。

よくできた創作物は、歪んだ鏡のようなものかもしれない。人は目を凝らして何がそこにあるか見ようとするが、見えるものは変形された自分の一部なのだ。でもそれが、創作物の価値なのだろう。

つくる人を近くで見てると、優れた作品というのは作者がつくるのではなく、何か大きなものにつくらされるんだろうなという気がする。つくらせる力というのは巨大な潮流で、作者の意図なんてのは海に流される筏の上でぱちゃぱちゃやってる程度の、わずかな力にすぎない。でも、転覆せずに潮の流れ至るところまで筏を導くのは至難の技で、そこに創作者のすごさというのがあるわけだが。

関連するかもしれないエントリ: 登場人物の気持ちを述べよ

Tag: 表現

2011/09/22

ピアノレッスン15回目

  • 基礎: スケールとアルペジオ。MM=144。まだ半分位の調でよくスベるけど、 だんだん慣れてきた。レッスンの最初の頃はMM=96でさらってて、144なんて想像を越えてたけど、続けてればなんとかなるもんだな。
  • RavelのSonatine。第1楽章速く入りすぎてもつれたので弾き直し。急ぐよりは粒を揃えた方がベター。第2楽章、ペダルを濁らないように気をつけてたんだけどドライになりすぎた感じ。あと拍が強すぎて"bouncy"と言われた。第3楽章はok。

Tag: Piano

More entries ...