From: Pol Lux on
On Jun 10, 1:57 am, Han de Bruijn <umum...(a)gmail.com> wrote:
> On Jun 10, 8:49 am, Han de Bruijn <umum...(a)gmail.com> wrote:
>
>
>
>
>
> > On Jun 9, 5:12 pm, Pollux <po....(a)gmail.com> wrote:
>
> > > (6/9/10 1:45 AM), Han de Bruijn wrote:
>
> > > > On Jun 7, 6:33 pm, Pollux<po....(a)gmail.com>  wrote:
> > > >> Hi -
>
> > > >> I have a non-symmetric 0/1 sparse matrix with uniform distribution of
> > > >> the non-zeros. I want to find a permutation of the columns that will
> > > >> result in a sparse matrix with the minimum distance between first and
> > > >> last non-zero on each row. (the hope is to implement a spmv algorithm
> > > >> that will take advantage of "semi-dense" blocks for speed, if there are
> > > >> enough semi-dense blocks, and if their density of non-zeros is high enough)
>
> > > >> Is there a known algorithm for that problem?
>
> > > >> Thanks!
>
> > > >> Pollux
>
> > > >> --- news://freenews.netfront.net/ - complaints: n...(a)netfront.net ---
>
> > > > Depends on what you are going to do with it. Solve a system of linear
> > > > equations and that's it? There do exist efficient _iterative_ solvers
> > > > as well ..
>
> > > > Han de Bruijn
>
> > > That's the thing. I am aware that there are algorithms to reduce fill-in
> > > while solving systems of linear equations, but that's not what I want to
> > > do. I just want to do matrix vector multiplication, but doing a
> > > "semi-dense" multiplication, instead of a sparse multiplication. The
> > > sparse multiplication has an indirection that doesn't play nice with the
> > > cache. The completely dense multiplication visits too many zeros in the
> > > case of a sparse matrix. So, I'm hoping I can find a happy middle ground
> > > where a permutation would locally increase the number of non-zeros in
> > > blocks. That way, I would use a dense multiplication in the blocks, even
> > > if I have to visit a few zeros (not too many), and I'll handle a few
> > > "stragglers" with the usual sparse indirection.
>
> > > Does that make sense?
>
> > > Pollux
>
> > > --- news://freenews.netfront.net/ - complaints: n...(a)netfront.net ---
>
> > Yes that makes sense. Can you predict (i.e. calculate) the coordinates
> > [i,j] of the non-zeros in the matrix beforehand ? I'm just asking some
> > questions because there are a hundred strategies possible in principle
> > and there's always a trade-off between space and time. You think there
> > is a "general" strategy, maybe, but I'm not so sure. The more details
> > you provide, the more likely that others can help.
>
> ----------------------------------------------------
>
> program Sparse;
>
> const
>   M : integer = 10; { number of columns }
>   N : integer = 15; { number of rows }
> type
>   store = record
>     i,j : integer;
>   end;
> var
>   matrix : array of store;
>   t,i,j,L : integer;
>   B,uit : array of double;
>   S : array of array of integer;
> begin
>   SetLength(matrix,M*N div 4); { Much less storage !! }
>   SetLength(B,M);
>   for i := 0 to M-1 do
>     B[i] := Random;
>   SetLength(uit,N);
>   for j := 0 to N-1 do
>     uit[j] := 0;
>
>   t := 0;
>   for i := 0 to M-1 do
>   begin
>     for j := 0 to N-1 do
>     begin
>       if Random > 0.9 then { sparse !! }
>       begin
>         matrix[t].i := i;
>         matrix[t].j := j;
>         t := t + 1;
>       end;
>     end;
>   end;
>   L := t;
>
>   SetLength(S,M,N); { Only for printout / SCRATCH }
>   Writeln('Matrix:');
>   for i := 0 to M-1 do
>   for j := 0 to N-1 do
>     S[i,j] := 0;
>   for t := 0 to L-1 do
>     S[matrix[t].i,matrix[t].j] := 1;
>   for j := 0 to N-1 do
>   begin
>     for i := 0 to M-1 do
>       Write(S[i,j]:2);
>     Writeln;
>   end;
>
>   Writeln('Vector:');
>   for i := 0 to M-1 do
>     Write(B[i]:7:3);
>   Writeln;
>   for i := 0 to 6 do
>     Write('1234567890'); { fit to page in Google groups }
>   Writeln;
>
> { !! Matrix-vector multiplication reduces to !! : }
>   for t := 0 to L-1 do
>     uit[matrix[t].j] := uit[matrix[t].j] + 1 * B[matrix[t].i];
>
>   Writeln('Result:');
>   for j := 0 to N-1 do
>     Writeln(uit[j]:7:3);
> end.
>
> ----------------------------------------------------
>
> Running the program gives the following.
>
> ----------------------------------------------------
>
> Matrix:
>  0 0 0 0 0 0 0 0 1 0
>  0 0 0 1 1 0 0 0 0 0
>  0 0 0 0 0 0 0 0 0 0
>  0 0 0 0 0 0 0 0 0 0
>  0 0 0 0 0 1 0 0 0 0
>  0 0 0 0 0 0 0 0 0 0
>  1 0 0 0 1 0 0 0 0 1
>  0 0 0 0 0 0 0 0 0 0
>  0 0 0 0 0 0 0 0 0 0
>  0 0 0 0 0 0 0 0 0 1
>  0 1 0 0 0 0 0 0 0 0
>  0 0 0 1 0 0 0 0 0 0
>  0 0 0 0 0 0 0 0 0 0
>  0 0 0 0 0 0 0 0 0 0
>  0 0 1 0 0 0 0 0 0 0
> Vector:
>   0.000  0.031  0.861  0.203  0.273  0.672  0.319  0.162  0.372  0.426
> 1234567890123456789012345678901234567890123456789012345678901234567890
> Result:
>   0.372
>   0.476
>   0.000
>   0.000
>   0.672
>   0.000
>   0.699
>   0.000
>   0.000
>   0.426
>   0.031
>   0.203
>   0.000
>   0.000
>   0.861
>
> ----------------------------------------------------
>
> Note that only the non-zeros of the matrix are stored.
> Disclaimer: from the top of my head / no optimization.
>
> Is that what you mean, approximately ?
>
> Han de Bruijn

Also, I think it is more standard to use m for number of rows, and n
for number of columns, with i from 0 to m, and j from 0 to n.

Pollux
From: Pol Lux on
On Jun 7, 3:11 pm, Nicolas Bonneel <nbonn...(a)cs.ubc.ca> wrote:
> Pol Lux wrote:
> > On Jun 7, 12:55 pm, Nicolas Bonneel <nbonn...(a)cs.ubc.ca> wrote:
> >> Pollux wrote:
> >>> Hi -
> >>> I have a non-symmetric 0/1 sparse matrix with uniform distribution of
> >>> the non-zeros. I want to find a permutation of the columns that will
> >>> result in a sparse matrix with the minimum distance between first and
> >>> last non-zero on each row. (the hope is to implement a spmv algorithm
> >>> that will take advantage of "semi-dense" blocks for speed, if there are
> >>> enough semi-dense blocks, and if their density of non-zeros is high enough)
> >>> Is there a known algorithm for that problem?
> >> I would say a Cuthill MacKee algorithm ? It allows to minimize the
> >> profile of the matrix.
>
> > Hmmm.... I am aware of this algorithm, but it seems to be for
> > symmetric matrices only. Is it worth giving it a shot on a non-
> > symmetric matrix?
>
> Hi,
> indeed sorry, I read too quickly.
> The package "Metis" has a reordering routine for non symmetric matrices.
> I think the geometry library CGAL also has some reodering routines but I
> don't know if it's for symmetric matrices only or not.
>
> Cheers

Thanks. Metis seems to be the most promising.

Pollux