演習として、 Darren Wilkinsonによるブログ投稿Gibbssamplerのサンプルプログラムをさまざまな言語で書き直しました(再検討) 。
コードは以下のとおりです。このコードは、SBCL 1.0.56を使用して、私の(5歳の)マシンで約53秒で実行され、buildappを使用してコアイメージを作成し、
time ./gibbs > gibbs.dat
これが投稿の他の言語のタイミングの計算方法だったので、私はそれに匹敵する何かをするだろうと思いました。投稿のCコードは約25秒で実行されます。可能であれば、Lispコードを高速化してみたいと思います。
##############################
gibbs.lisp
##############################
(eval-when (:compile-toplevel :load-toplevel :execute)
(require :cl-rmath) (setf *read-default-float-format* 'double-float))
(defun gibbs (N thin)
(declare (fixnum N thin))
(declare (optimize (speed 3) (safety 1)))
(let ((x 0.0) (y 0.0))
(declare (double-float x y))
(print "Iter x y")
(dotimes (i N)
(dotimes (j thin)
(declare (fixnum i j))
(setf x (cl-rmath::rgamma 3.0 (/ 1.0 (+ (* y y) 4))))
(setf y (cl-rmath::rnorm (/ 1.0 (+ x 1.0)) (/ 1.0 (sqrt (+ (* 2 x) 2))))))
(format t "~a ~a ~a~%" i x y))))
(defun main (argv)
(declare (ignore argv))
(gibbs 50000 1000))
次に、 asgibbs
を呼び出しsh gibbs.sh
て実行可能ファイルをビルドしましたgibbs.sh
##################
gibbs.sh
##################
buildapp --output gibbs --asdf-tree /usr/share/common-lisp/source/ --asdf-tree /usr/local/share/common-lisp/source/ --load-system cl-rmath --load gibbs.lisp --entry main
SBCL 1.0.56でコンパイルすると、6つのコンパイラノートが表示されます。これは以下で再現されます。それらについてどうしたらよいかわかりませんが、ヒントをいただければ幸いです。
; compiling file "/home/faheem/lisp/gibbs.lisp" (written 30 MAY 2012 02:00:55 PM):
; file: /home/faheem/lisp/gibbs.lisp
; in: DEFUN GIBBS
; (SQRT (+ (* 2 X) 2))
;
; note: unable to
; optimize
; due to type uncertainty:
; The result is a (VALUES (OR (DOUBLE-FLOAT 0.0) (COMPLEX DOUBLE-FLOAT))
; &OPTIONAL), not a (VALUES FLOAT &REST T).
; (/ 1.0d0 (SQRT (+ (* 2 X) 2)))
;
; note: unable to
; optimize
; due to type uncertainty:
; The second argument is a (OR (DOUBLE-FLOAT 0.0)
; (COMPLEX DOUBLE-FLOAT)), not a (COMPLEX
; DOUBLE-FLOAT).
;
; note: forced to do static-fun Two-arg-/ (cost 53)
; unable to do inline float arithmetic (cost 12) because:
; The second argument is a (OR (DOUBLE-FLOAT 0.0) (COMPLEX DOUBLE-FLOAT)), not a DOUBLE-FLOAT.
; The result is a (VALUES (OR (COMPLEX DOUBLE-FLOAT) (DOUBLE-FLOAT 0.0))
; &OPTIONAL), not a (VALUES DOUBLE-FLOAT &REST T).
; (CL-RMATH:RGAMMA 3.0d0 (/ 1.0d0 (+ (* Y Y) 4)))
;
; note: doing float to pointer coercion (cost 13)
; (SQRT (+ (* 2 X) 2))
;
; note: doing float to pointer coercion (cost 13)
; (CL-RMATH:RNORM (/ 1.0d0 (+ X 1.0d0)) (/ 1.0d0 (SQRT (+ (* 2 X) 2))))
;
; note: doing float to pointer coercion (cost 13)
;
; compilation unit finished
; printed 6 notes
; /home/faheem/lisp/gibbs.fasl written
; compilation finished in 0:00:00.073
更新1:Rainer Joswigの回答は、SQRTの議論が否定的である可能性があることを指摘しました。これは、私が見たあいまいなコンパイラノートのソースでした。
The result is a (VALUES (OR (DOUBLE-FLOAT 0.0) (COMPLEX DOUBLE-FLOAT))
; &OPTIONAL), not a (VALUES FLOAT &REST T).
コンパイラは、引数の値が正であるかどうかがわからないため、結果が複素数になる可能性があると不満を漏らしていました。この例では、値x
はガンマ分布からのサンプル変量であるため、常に0より大きくなります。SBCLユーザーメーリングリストのStephanによって有益に指摘されました(スレッド「Optimizing a simple Common Lisp」の2番目のメッセージを参照してください)。ギブスサンプラープログラム」、これは次のようにxがゼロより大きいかゼロであると宣言することで解決できることを示しています。
(declare (type (double-float 0.0 *) x))
FLOATタイプ と 間隔指定子に関する関連ドキュメントについては、CommonLispHyperspecを参照してください 。
これはコードを少しスピードアップするようです。現在は確実に52秒未満ですが、それでもそれほど大きなメリットはありません。これはまたについてのメモを残します
note: doing float to pointer coercion (cost 13)
なんらかの理由で修正できない場合は、その理由を教えてください。また、メモの意味の説明は、関係なく興味深いでしょう。具体的には、pointer
ここでその言葉はどういう意味ですか?これは、C関数が呼び出されているという事実に関連していますか?また、コスト13はあまり役に立たないようです。何が測定されていますか?
また、Rainerは、実行時間を短縮する可能性のあるconsingを削減できる可能性があることを示唆しました。コンシングリダクションが可能か、それともランタイムを短縮できるかはわかりませんが、意見やアプローチに興味があります。全体として、この機能のパフォーマンスを改善するためにできることはあまりないようです。多分それは小さすぎて単純です。