Common Lispでの行列乗算

Question

私は線形代数を使うCL（SBCL 1.2.15）でプログラムを書いています。実行の過程で、行列にベクトルを掛けることがよくあります。

プロファイラでは、ほとんどの時間（80％）はプログラムが正確にそれを行い、行列をベクトルで乗算するのに費やしたことがわかりました。また、この関数はGCをトリガーする多くのコンサイン（100x100行列に対して〜100コールに対して80,000,000）を実行することも示しています。 F＃（静的型チェックの利点がありますが、そうでなければSBCLを上回ることはありません）で同様のことをしていれば、F＃プログラムは10倍高速に実行されます。

私は間違っていますか？

(defun matrix-mul (matrix vector dest) 
    "Multiply MATRIX by VECTOR putting the result into DEST. 
Optimized for DOUBLE-FLOAT vectors and matrices" 
    (declare (type (array double-float (* *)) matrix) 
      (type (vector double-float *) vector dest) 
      (optimize (speed 3) (debug 0) (safety 0))) 
    (destructuring-bind (rows cols) (array-dimensions matrix) 
    (dotimes (i rows) 
    (setf (aref dest i) 
      (loop for j below cols sum (the double-float 
              (* (aref matrix i j) (aref vector j)))))))) 

PS。私も内部ループの代わりにDOTIMESを使用しようとしました - 違いはありません

PPS次のオプションは、CLからのBLASを使用することができますが、（i）Windows上で動作させることを期待していません。（ii）行列/ベクトルを前後にマーシャリングする必要があります。

Answer 1

通常の問題は、浮動小数点演算（行列乗算のような周囲のコードとは独立したdouble-floats）がconsesであることです。一般

このような場合にはSBCLで最初に行うべきこと：

をファイルにコードを入れてコンパイルし、それ

コンパイラは、出力その1つのコンパイル与えられた最適化問題の多くは、速度のために。あなたはメモを調べ、あなたができることを見る必要があります。

たとえば、LOOPの合計にはタイプ情報がありません。

sum変数の型を宣言する実際にはLOOPの構文があります。 SBCLがそれを利用するかどうかは分かりません。

(loop repeat 10 
     sum 1.0d0 of-type double-float) 

SBCL 1.3。

* (compile-file "/tmp/test.lisp") 

; compiling file "/tmp/test.lisp" (written 13 DEC 2015 11:34:26 AM): 
; compiling (DEFUN MATRIX-MUL ...) 
; file: /tmp/test.lisp 

1）

; in: DEFUN MATRIX-MUL 
;  (SETF (AREF DEST I) 
;    (LOOP FOR J BELOW COLS 
;     SUM (THE DOUBLE-FLOAT (* # #)))) 
; --> LET* FUNCALL SB-C::%FUNCALL (SETF AREF) 
; ==> 
; (SB-KERNEL:HAIRY-DATA-VECTOR-SET ARRAY SB-INT:INDEX SB-C::NEW-VALUE) 
; 
; note: unable to 
; avoid runtime dispatch on array element type 
; due to type uncertainty: 
; The first argument is a (VECTOR DOUBLE-FLOAT), not a SIMPLE-ARRAY. 

2）

;  (AREF MATRIX I J) 
; --> LET* 
; ==> 
; (SB-KERNEL:HAIRY-DATA-VECTOR-REF ARRAY SB-INT:INDEX) 
; 
; note: unable to 
; avoid runtime dispatch on array element type 
; due to type uncertainty: 
; The first argument is a (ARRAY DOUBLE-FLOAT (* *)), not a SIMPLE-ARRAY. 

3）

;  (AREF VECTOR J) 
; ==> 
; (SB-KERNEL:HAIRY-DATA-VECTOR-REF ARRAY SB-INT:INDEX) 
; 
; note: unable to 
; avoid runtime dispatch on array element type 
; due to type uncertainty: 
; The first argument is a (VECTOR DOUBLE-FLOAT), not a SIMPLE-ARRAY. 

4）

012：コードのための32ビットARMに

5）

; --> BLOCK LET SB-LOOP::WITH-SUM-COUNT LET SB-LOOP::LOOP-BODY TAGBODY WHEN IF 
; --> >= OR LET IF OR THE = IF 
; ==> 
; (= SB-C::X SB-C::Y) 
; 
; note: unable to open code because: The operands might not be the same type. 

6）

;  (DOTIMES (I ROWS) 
;  (SETF (AREF DEST I) 
;    (LOOP FOR J BELOW COLS 
;      SUM (THE DOUBLE-FLOAT #)))) 
; --> DO BLOCK LET TAGBODY UNLESS IF >= IF 
; ==> 
; (< SB-C::X SB-C::Y) 
; 
; note: forced to do static-fun Two-arg-< (cost 53) 
;  unable to do inline fixnum comparison (cost 4) because: 
;  The second argument is a INTEGER, not a FIXNUM. 
;  unable to do inline (signed-byte 32) comparison (cost 6) because: 
;  The second argument is a INTEGER, not a (SIGNED-BYTE 32). 
;  etc. 

7）

;  (LOOP FOR J BELOW COLS 
;   SUM (THE DOUBLE-FLOAT (* (AREF MATRIX I J) (AREF VECTOR J)))) 
; --> BLOCK LET SB-LOOP::WITH-SUM-COUNT LET SB-LOOP::LOOP-BODY TAGBODY WHEN IF 
; --> >= OR LET > IF 
; ==> 
; (> SB-C::X SB-C::Y) 
; 
; note: forced to do static-fun Two-arg-> (cost 53) 
;  unable to do inline fixnum comparison (cost 4) because: 
;  The second argument is a REAL, not a FIXNUM. 
;  unable to do inline (signed-byte 32) comparison (cost 6) because: 
;  The second argument is a REAL, not a (SIGNED-BYTE 32). 
;  etc. 

8）

; --> BLOCK LET SB-LOOP::WITH-SUM-COUNT LET SB-LOOP::LOOP-BODY TAGBODY SETQ THE 
; ==> 
; (+ #:LOOP-SUM-8 (THE DOUBLE-FLOAT (* (AREF MATRIX I J) (AREF VECTOR J)))) 
; 
; note: forced to do static-fun Two-arg-+ (cost 53) 
;  unable to do inline float arithmetic (cost 2) because: 
;  The first argument is a NUMBER, not a DOUBLE-FLOAT. 
;  The result is a (VALUES NUMBER &OPTIONAL), not a (VALUES DOUBLE-FLOAT 
;                &REST T). 
; 
; note: doing float to pointer coercion (cost 13), for: 
;  the second argument of static-fun Two-arg-+ 
; 
; compilation unit finished 
; printed 10 notes