comparison yetilab/matrix/matrix.yeti @ 271:c206de7c3018

Faster sparseProductLeft. Similar code would work for other sparseProducts
author Chris Cannam
date Thu, 23 May 2013 16:12:21 +0100
parents 53ff481f1a41
children 2ebda6646c40
comparison
equal deleted inserted replaced
270:d7dd391a90fd 271:c206de7c3018
450 newMatrix (typeOf m) (map vfilter (asColumns m)); 450 newMatrix (typeOf m) (map vfilter (asColumns m));
451 fi; 451 fi;
452 fi; 452 fi;
453 453
454 sparseProductLeft size m1 m2 = 454 sparseProductLeft size m1 m2 =
455 (e = enumerateSparse m1; 455 (d1 = case m1 of
456 SparseCSR d: d;
457 SparseCSC d: d;
458 _: failWith "sparseProductLeft called for non-sparse m1";
459 esac;
460 swap = not isRowMajor? m1;
456 data = array (map \(new double[size.rows]) [1..size.columns]); 461 data = array (map \(new double[size.rows]) [1..size.columns]);
462 vv = d1.values;
463 ii = d1.indices;
464 pp = d1.pointers;
457 for [0..size.columns - 1] do j': 465 for [0..size.columns - 1] do j':
458 c = getColumn j' m2; 466 c = getColumn j' m2;
459 for e do { v, i, j }: 467 var i = 0;
460 data[j'][i] := data[j'][i] + v * (vec.at c j); 468 for [0..length ii - 1] do ix:
469 j = ii[ix];
470 ix == pp[i+1] loop (i := i + 1);
471 if swap then
472 data[j'][j] := data[j'][j] + (vec.at vv ix) * (vec.at c i);
473 else
474 data[j'][i] := data[j'][i] + (vec.at vv ix) * (vec.at c j);
475 fi;
461 done; 476 done;
462 done; 477 done;
463 DenseCols (array (map vec.vector (list data)))); 478 DenseCols (array (map vec.vector (list data))));
464 479
465 sparseProductRight size m1 m2 = 480 sparseProductRight size m1 m2 =