(Usually sparse matrices are stored a bit different: for every row and every column there's a linked list of non-zero items together with their indexes, so that walking all non-zero elements is fast, ie following the links).

I'd like to quickly jump in to advocate Math::MatrixSparse, a pure-Perl sparse matrix library. It does just what its name promises. Note, however, that the original code has some annoying bugs, therefore, I maintain a separate version here.

rg0now

Maybe you should ask the author for permissions to maintain the official module on CPAN? If you can't reach the author, have a look at the official PAUSE FAQ on this and mail the modules mailing list.

Cheers,
Steffen

I'd suggest PDL::CCS (collapsed-column storage).

I usually reach for PDL in situations like this.

(Wow, it's like deja vu all over again...)

So, you want to multiply your sparse array by a not-sparse vector @p. Seems straightforward enough:

  my @r = (0) x scalar(@p) ;

  for my $i (0..$#row_index) {
    $r[$row_index[$i]-1] += $realval[$i] * $p[$col_index[$i]-1] ;
  } ;
[download]

Dear Oshalla,
Just to let you know. I've ported your ingenious code to C++. They are so indispensable! Thanks a million again.

#include <iostream>
#include <vector>
#include <boost/assign/std/vector.hpp>
#include <vector>
#include <stdio.h>
#include <stdlib.h>
using namespace boost::assign;
using namespace std;


int main  ( int arg_count, char *arg_vec[] ) {
    vector <double> realValue;
    realValue +=
       0.994704478,
       0.989459074,
       0.994717023,
       1.000000000,
       1.000000000,
       0.002647761,
       0.005282977,
       0.000882587,
       0.005270463,
       0.000882587,
       0.002635231,
       0.000882587,
       0.002635231;

    vector<int>rowIndex;
    rowIndex +=  1, 2, 3, 4, 5, 1, 3, 1, 2, 1, 2, 1, 2;

    vector<int>colIndex;
    colIndex  +=   1, 2, 3, 4 ,5, 5, 3, 2, 1, 3, 3, 4, 4;

    vector<int>matDim;
    matDim  += 5,5; // M always equal to N
    
    vector <double> theP;
    theP += 0.4, 0.2, 0.2, 0.2, 0.2;
    
  
    vector <double> Result; 
    Result.assign(theP.size(), 0);

     for (int i= 0; i < rowIndex.size(); i++) {
        Result[rowIndex[i]-1] += realValue[i] * theP[colIndex[i]-1];
     }    

    // print it
    for (int k =0; k < Result.size(); k++){
        cout << Result[k] << endl;
    }
    
 
    
    
    return 0;
}
[download]




---
neversaint and everlastingly indebted.......

Dear oshalla,
Thank you so much for your reply. It is truly useful.

I also tried to modify your code to do multiplication of transpose of the sparse matrix. With the following code:

use Data::Dumper;
my @realval = (
    0.994704478,
    0.989459074,
    0.994717023,
    1.000000000, 
    1.000000000,
    0.002647761,
    0.005282977, 
    0.000882587, 
    0.005270463, 
    0.000882587, 
    0.002635231, 
    0.000882587,
    0.002635231,);

my @row_index = qw( 1 2 3 4 5 1 3 1 2 1 2 1 2);
my @col_index =  qw(1 2 3 4 5 5 3 2 1 3 3 4 4);
my @mat_dim  = (5,5); # M always equal to N

my @trans_r = (0) x scalar(@p);

for my $j (0..$#row_index) {
 $trans_r[$col_index[$i]-1] += $realval[$i] * $p[$row_index[$i]-1] ;
}

print Dumper \@trans_r;
[download]

Not sure why it doesn't give this desired result:

@result = (0.3989359, 0.1982448, 0.2008801, 0.2008801, 0.2010591);
[download]

Truly need your expert advice on this. Thanks and hope to hear from you again.


---
neversaint and everlastingly indebted.......

You really should make friends with use strict and its cousin use warnings... With the code as given:

• use strict complains about $i, which is used in the loop, but doesn't seem to be declared anywhere

  -- that's what it means by 'Global symbol "$i" requires explicit package name'... we just have to accept that when looked at from the correct angle, this makes perfect sense.

• use warnings would warn you that each time $i is used it is uninitialised (though with use strict you wouldn't get that far !).</c>

please accept my sincere apologies if you find the following question trivial.

Could anyone tell me why a self-made sparse matrix multiplication without PDL is bound to fail, as moritz said? In my opinion, anyone can use twice-nested hashs in order to represent a sparse 2d-matrix:

$matrix{$row}{$column}=$value
[download]

I don't see why this solution should be slower than using piddles.