Pure Haskell 10x-100x più veloce di HMatrix per piccole matrici?

Question

Passiamo la maggior parte dei nostri cicli di CPU su operazioni che coinvolgono piccole matrici, quindi mi chiedevo se fosse possibile ottimizzare per quel caso. Si consideri il seguente codice:

module Main where 

import Numeric.LinearAlgebra.HMatrix 
import Criterion.Main 


data Matrix2x2 = Matrix2x2 {-# UNPACK #-} !Double !Double !Double !Double 

mul2x2p :: Matrix2x2 -> Matrix2x2 -> Matrix2x2 
mul2x2p (Matrix2x2 a1 b1 c1 d1) (Matrix2x2 a2 b2 c2 d2) = 
    Matrix2x2 (a1*a2 + b1*c2) (a1*b2 + b1*d2) (c1*a2 + d1*c2) (c1*b2 + d1*d2) 

inv2x2 :: Matrix2x2 -> Matrix2x2 
inv2x2 (Matrix2x2 a b c d) = 
    let detInv = a * d - b * c 
    in Matrix2x2 (d/detInv) (-b/detInv) (-c/detInv) (a/detInv) 

add2x2 (Matrix2x2 a1 b1 c1 d1) (Matrix2x2 a2 b2 c2 d2) = 
    Matrix2x2 (a1+a2) (b1+b2) (c1+c2) (d1+d2) 

hm1 = matrix 2 [1, 2, 3, 4] 
hm2 = matrix 2 [5, 6, 7, 8] 

pm1 = Matrix2x2 1 2 3 4 
pm2 = Matrix2x2 5 6 7 8 

main = defaultMain [ 
    bgroup "matrix tests" [ bench "pure mult"  $ whnf (mul2x2p pm1) pm2 
         , bench "hmatrix mult" $ whnf (hm1 <>) hm2 
         , bench "pure add"  $ whnf (add2x2 pm1) pm2 
         , bench "hmatrix add" $ whnf (hm1 +) hm2 
         , bench "pure inv"  $ whnf inv2x2 pm1 
         , bench "hmatrix inv" $ whnf inv hm1 
         ]] 

I risultati:

benchmarking matrix tests/pure mult 
time     6.461 ns (6.368 ns .. 6.553 ns) 
        0.999 R² (0.998 R² .. 0.999 R²) 
mean     6.482 ns (6.394 ns .. 6.594 ns) 
std dev    345.1 ps (271.4 ps .. 477.3 ps) 
variance introduced by outliers: 77% (severely inflated) 

benchmarking matrix tests/hmatrix mult 
time     180.6 ns (178.2 ns .. 183.1 ns) 
        0.999 R² (0.998 R² .. 0.999 R²) 
mean     183.0 ns (180.6 ns .. 186.3 ns) 
std dev    9.363 ns (7.405 ns .. 12.73 ns) 
variance introduced by outliers: 71% (severely inflated) 

benchmarking matrix tests/pure add 
time     6.262 ns (6.223 ns .. 6.297 ns) 
        0.999 R² (0.999 R² .. 1.000 R²) 
mean     6.281 ns (6.220 ns .. 6.355 ns) 
std dev    235.0 ps (183.3 ps .. 321.0 ps) 
variance introduced by outliers: 62% (severely inflated) 

benchmarking matrix tests/hmatrix add 
time     116.4 ns (115.0 ns .. 117.9 ns) 
        0.999 R² (0.998 R² .. 0.999 R²) 
mean     116.3 ns (115.2 ns .. 117.7 ns) 
std dev    4.176 ns (3.447 ns .. 5.150 ns) 
variance introduced by outliers: 55% (severely inflated) 

benchmarking matrix tests/pure inv 
time     7.811 ns (7.718 ns .. 7.931 ns) 
        0.999 R² (0.998 R² .. 0.999 R²) 
mean     7.895 ns (7.808 ns .. 7.988 ns) 
std dev    296.4 ps (247.2 ps .. 358.3 ps) 
variance introduced by outliers: 62% (severely inflated) 

benchmarking matrix tests/hmatrix inv 
time     908.5 ns (901.3 ns .. 916.6 ns) 
        0.999 R² (0.998 R² .. 0.999 R²) 
mean     934.0 ns (917.6 ns .. 961.3 ns) 
std dev    73.92 ns (50.53 ns .. 108.6 ns) 
variance introduced by outliers: 84% (severely inflated) 

Le mie domande sono:

1) è la velocità fino reale o a causa di un artefatto con il processo di benchmarking?

2) Se la velocità è reale, esiste una libreria esistente che gestirà, ad esempio, matrici 1x1, 2x2, 3x3, 4x4 come casi speciali?

3) In caso contrario, qual è il modo migliore per avvolgere HMatrix in modo da poter prendere il percorso veloce quando la matrice è piccola? ghc può solo decomprimere i record con un costruttore. C'è un modo per generare automaticamente diverse versioni del nostro codice, ecc

esempio-test.cabal:

name:    example-test 
version:    0.1.0.0 
build-type:   Simple 
cabal-version:  >=1.10 
executable example-test 
    main-is: 
    Main.hs 
    build-depends: 
    base >=4.7 && <4.8, 
    criterion, 
    hmatrix 
    default-language: 
    Haskell2010 
    ghc-options: 
    -H12G -O3 -optc-O3 -fllvm -rtsopts -threaded -fexcess-precision -j6 +RTS -N6 -RTS -fno-ignore-asserts -fcontext-stack=150 
    -- -fforce-recomp 

Answer 1

HMatrix fornisce un'interfaccia funzionale all'algebra standard di densa lineare (SVD, autovalori, sistemi lineari , ecc.) basato sulle eccellenti e ottimizzate librerie BLAS/LAPACK BLAS/LAPACK.

C'è un overhead causato dalla FFI chiamate, allocazione di memoria, controlli di errore, codice per evitare alcuni matrice traspone, ecc Questo è trascurabile rispetto al il tempo richiesto da tipici calcoli di matrice di dimensioni contenute e grandi, ma per gli array molto piccoli e le operazioni semplici il sovraccarico può essere grande.

Non ci sono piani immediati per reimplementare calcoli per piccoli array. È non facile fare cose migliori di BLAS/LAPACK, e non sono sicuro che il guadagno della velocità per una particolare applicazione sia più importante della semplicità del codice e della generalità .

Se l'applicazione richiede un numero elevato di operazioni semplici su matrici 2x2, o 3x3, allora altre librerie sono probabilmente più appropriate per questa attività .

Tuttavia, si consiglia di eseguire il programma in modalità di profilo per individuare i veri colli di bottiglia . Un benchmark che utilizza dati costanti potrebbe non essere completamente rappresentativo dei tempi trascorsi in un programma "normale".

Anche se si trova che la maggior parte del tempo viene impiegata in questi calcoli con piccole matrici , è possibile riscrivere l'algoritmo utilizzando un calcolo equivalente di dimensioni maggiori della matrice .

Answer 2

Per i piccoli vettori e matrici, il pacchetto lineare è spesso usato: https://hackage.haskell.org/package/linear