Band too slow when constructing a SparseArray
- To: mathgroup at smc.vnet.net
- Subject: [mg106980] Band too slow when constructing a SparseArray
- From: "Norbert P." <bertapozar at gmail.com>
- Date: Sat, 30 Jan 2010 07:12:50 -0500 (EST)
Hi everyone,
this has been driving me crazy for a while. I'm using Mathematica
7.0.1, Win XP. I had the same problem on 6.0.2.
I've implemented a nonlinear PDE solver that uses an implicit method,
i.e. I need to solve a system of linear equations at every step, say
~40,000. I'm surprised how fast Mathematica can be when everything is
implemented using SparseArray and packed arrays.
Now my problem: I need to construct a new matrix at each iteration
(nonlinear system). It is a sparse matrix with 5 bands. But using Band
to do this is prohibitively slow. The matrix construction using the
simple SparseArray[Band[{1,1}]->vector,...] takes much more than the
actual LinearSolve!
I found a way how to speed this up by at least 100x. It turns out that
the construction of a diagonal matrix using SparseArray[{i_,i_}:>vector
[[i]], {n,n}] is incredibly fast (this is the only pattern that is
fast). I can't use DiagonalMatrix here since that requires a
SparseArray as an argument to produce SparseArray, slowing the process
down (unfortunate new feature in v7). Then matrix operations like
addition or PadLeft, PadRight are blazing fast. So I construct 5
diagonal matrices using my vectors and combine them together this way.
Works nice.
My question: Why doesn't Band work at least as fast? SparseArray[Band
[{1,1}]->v] looks like it should work better then using pattern
SparseArray[{i_,i_}:>vector[[i]], {n,n}] . Is there a good reason why
I have to construct a banded matrices this way? Why isn't this
mentioned in the documentation?
A test case:
In[1]:=
v=RandomReal[1.,1000000];
(m1=SparseArray[{i_,i_}:>v[[i]],{Length[v],Length[v]}])//Timing
Tr[m1]-Tr[v]//Timing
SparseArray[{i_,i_}->1.,{Length[v],Length[v]}]//Timing
(m2=SparseArray[Band[{1,1}]->v,{Length[v],Length[v]}])//Timing
m1==m2
Out[2]= {0.016,SparseArray[<1000000>,{1000000,1000000}]}
Out[3]= {0.031,-9.31323*10^-9}
Out[4]= {0.078,SparseArray[<1000000>,{1000000,1000000}]}
Out[5]= {12.047,SparseArray[<1000000>,{1000000,1000000}]}
Out[6]= True
As you can see, using the pattern is 750x faster than using Band!! And
the construction of a diagonal matrix with values from a list is 3
times faster than doing the same with constants...
Best,
Norbert