Moltiplicazione matrice in Common Lisp

Question

Sto scrivendo il programma in CL (con SBCL 1.2.15) che utilizza algebra lineare. Durante il corso dell'esecuzione, moltiplica spesso una matrice per un vettore.

Profiler ha mostrato che la maggior parte delle volte (80%) il programma spende esattamente ciò, moltiplicando la matrice per vettore. Mostra anche che questa funzione fa molto consing (80.000.000 per ~ 100 chiamate per matrici 100x100), che innesca anche GC. Avendo fatto qualcosa di simile in F # (che ha un vantaggio nel controllo dei tipi statici, ma di solito non supera in modo efficace SBCL), il programma F # viene eseguito 10 volte più velocemente.

Sto sbagliando?

(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. Ho anche provato a usare DOTIMES anziché LOOP interno - nessuna differenza

PPS. L'opzione successiva può utilizzare BLAS da CL, ma (i) non vedo l'ora di farlo funzionare su Windows, (ii) avranno bisogno di matrici/vettori di marshalling avanti e indietro.

Answer 1

Il solito problema è che l'aritmetica fluttuante, qui con doppio float (indipendente dal codice circostante come moltiplicazione della matrice) può essere attivata.

In genere la prima cosa da fare con SBCL in questi casi:

mettere il codice in un file e compilarlo

il compilatore poi uscita un sacco di problemi di ottimizzazione, dato che uno compilazioni per la velocità. Dovresti quindi esaminare le note e vedere cosa puoi fare.

Qui ad esempio la somma LOOP manca di informazioni sul tipo.

In realtà esiste una sintassi LOOP per dichiarare il tipo della variabile somma. Non so se SBCL ne approfitti:

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

SBCL 1.3.0 su ARM a 32 bit per il codice:

* (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)

;  (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 SETQ THE 
; ==> 
; (+ #:LOOP-SUM-8 (THE DOUBLE-FLOAT (* (AREF MATRIX I J) (AREF VECTOR J)))) 
; 
; note: unable to 
; optimize 
; due to type uncertainty: 
; The first argument is a NUMBER, not a (COMPLEX SINGLE-FLOAT). 
; 
; note: unable to 
; optimize 
; due to type uncertainty: 
; The first argument is a NUMBER, not a (COMPLEX DOUBLE-FLOAT). 

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