行列（scipy sparse） - 行列（高密度配列）乗算効率

Question

私は地球物理学の逆転に取り組んでいる研究者です。どの線形システムを解決する必要があります：Au = rhs。ここでは希薄な行列ですが、rhsとuは密行列またはベクトルです。グラジエントベースの逆変換を実行するには、感度計算が必要であり、多数の行列行列と行列ベクトルの乗算が必要です。最近、私は行列（スパース）で奇妙な行動を発見した - マトリックス（密）の乗算を、以下の例である：私は密行列のサイズを大きくしていたときに

import numpy as np 
import scipy.sparse as sp 
n = int(1e6) 
m = int(100) 
e = np.ones(n) 
A = sp.spdiags(np.vstack((e, e, e)), np.array([-1, 0, 1]), n, n) 
A = A.tocsr() 
u = np.random.randn(n,m) 

%timeit rhs = A*u[:,0] 
#10 loops, best of 3: 22 ms per loop  
%timeit rhs = A*u[:,:10] 
#10 loops, best of 3: 98.4 ms per loop 
%timeit rhs = A*u 
#1 loop, best of 3: 570 ms per loop​ 

私はcompution時間のほぼ直線的な増加を期待していましたuに疎行列Aを掛けたものです（たとえば、2番目のものは私には220m秒、最後の1つはA*u[:,:10] 2.2sです）。しかし、それは私が予想したよりはるかに高速です。逆に、行列 - ベクトルの乗算は、行列 - 行列の乗算よりもずっと遅い。なぜ誰かが説明できますか？さらに、行列 - 行列乗算と同様の効率レベルを行列 - 行列乗算に高める効果的な方法がありますか？

Answer 1

あなたはsource code見れば（行列 - 行列乗算を実装）csr_matvecsを呼び出しとして実装されている間、あなたは、（行列 - ベクトル乗算を実装）csr_matvecは、Cコードで簡単な和ループとして実装されていることがわかりますaxpy BLASルーチンに転送します。インストールされているBLASライブラリの種類によっては、farが、行列 - ベクトル乗算に使用される単純なC実装よりも効率的です。それはおそらく、あなたが行列 - ベクトル乗算が非常に遅いのを見ている理由です。

マトリクスベクトルの場合にBLASを呼び出すようにscipyを変更すると、パッケージに役立つ効果があります。