| 1552 | |
| 1553 | template <typename TM_Res, typename TM1, typename TM2> |
| 1554 | shared_ptr<SparseMatrixTM<TM_Res>> |
| 1555 | MatMult (const SparseMatrixTM<TM1> & mata, const SparseMatrixTM<TM2> & matb, bool sort_output) |
| 1556 | { |
| 1557 | static Timer t ("sparse matrix multiplication"); |
| 1558 | static Timer t1a ("sparse matrix multiplication - setup a"); |
| 1559 | static Timer t1b ("sparse matrix multiplication - setup b"); |
| 1560 | static Timer t1b1 ("sparse matrix multiplication - setup b1"); |
| 1561 | static Timer t2 ("sparse matrix multiplication - mult"); |
| 1562 | |
| 1563 | static Timer tsort ("sparse matrix multiplication - sort"); |
| 1564 | RegionTimer reg(t); |
| 1565 | |
| 1566 | t1a.Start(); |
| 1567 | |
| 1568 | // find graph of product |
| 1569 | Array<int> cnt(mata.Height()); |
| 1570 | cnt = 0; |
| 1571 | |
| 1572 | /* |
| 1573 | ParallelForRange |
| 1574 | (mata.Height(), [&] (IntRange r) |
| 1575 | { |
| 1576 | Array<BaseSparseMatrix::ColIdx*> ptrs; |
| 1577 | Array<int> sizes; |
| 1578 | for (int i : r) |
| 1579 | { |
| 1580 | auto mata_ci = mata.GetRowIndices(i); |
| 1581 | ptrs.SetSize(mata_ci.Size()); |
| 1582 | sizes.SetSize(mata_ci.Size()); |
| 1583 | for (int j : Range(mata_ci)) |
| 1584 | { |
| 1585 | ptrs[j] = matb.GetRowIndices(mata_ci[j]).Addr(0); |
| 1586 | sizes[j] = matb.GetRowIndices(mata_ci[j]).Size(); |
| 1587 | } |
| 1588 | int cnti = 0; |
| 1589 | MergeArrays(ptrs, sizes, [&cnti] (int col) { cnti++; } ); |
| 1590 | cnt[i] = cnti; |
| 1591 | } |
| 1592 | }, |
| 1593 | TasksPerThread(10)); |
| 1594 | */ |
| 1595 | |
| 1596 | ParallelForRange |
| 1597 | (mata.Height(), [&] (IntRange r) |
| 1598 | { |
| 1599 | Array<int8_t> flags(matb.Width()); |
| 1600 | flags = 0; |
| 1601 | Array<int> list; |
| 1602 | |
| 1603 | for (int i : r) |
| 1604 | { |
| 1605 | for (int cola : mata.GetRowIndices(i)) |
| 1606 | for (int col : matb.GetRowIndices(cola)) |
| 1607 | { |
| 1608 | if (!flags[col]) |
| 1609 | { |
| 1610 | flags[col]=1; |
| 1611 | list += col; |