70 CSparseMatrixElement::operator
const C_FLOAT64 & ()
const
80 mThreshold(std::numeric_limits<
C_FLOAT64 >::epsilon()),
89 mElement(*this, this->mSearchRow, this->mSearchCol, 0.0)
93 mThreshold(std::numeric_limits<
C_FLOAT64 >::epsilon()),
102 mElement(*this, this->mSearchRow, this->mSearchCol, 0.0)
110 std::vector< std::vector< CSparseMatrixElement * > >::iterator itRow;
111 std::vector< std::vector< CSparseMatrixElement * > >::iterator endRow;
113 std::vector< CSparseMatrixElement * >::iterator itElement;
114 std::vector< CSparseMatrixElement * >::iterator endElement;
116 for (itRow =
mRows.begin(), endRow =
mRows.end(); itRow != endRow; ++itRow)
117 for (itElement = itRow->begin(), endElement = itRow->end(); itElement != endElement; ++itElement)
137 mRows.resize(mNumRows);
138 mCols.resize(mNumCols);
153 std::vector< std::vector< CSparseMatrixElement * > >::const_iterator itRow;
154 std::vector< std::vector< CSparseMatrixElement * > >::const_iterator endRow;
156 for (itRow =
mRows.begin(), endRow =
mRows.end(); itRow != endRow; ++itRow)
157 NonZeros += itRow->size();
185 for (j = 0, i = 0; j <
mNumCols; j++, pColumnStart++)
187 mCols[j].resize(*pColumnStart - * (pColumnStart - 1));
189 for (k = 0; i < *pColumnStart; i++, k++, pRowIndex++, pValue++)
192 mCols[j][k] = pElement;
193 mRows[*pRowIndex].push_back(pElement);
212 for (j = 0; j <
mNumCols; j++, pTmp++)
216 mRows[i].push_back(pElement);
217 mCols[j].push_back(pElement);
229 M.
resize(mNumRows, mNumCols);
234 std::vector< std::vector< CSparseMatrixElement * > >::const_iterator itRow;
235 std::vector< std::vector< CSparseMatrixElement * > >::const_iterator endRow;
237 std::vector< CSparseMatrixElement * >::const_iterator itElement;
238 std::vector< CSparseMatrixElement * >::const_iterator endElement;
240 for (itRow = mRows.begin(), endRow = mRows.end(); itRow != endRow; ++itRow)
242 for (i = 0, itElement = itRow->begin(), endElement = itRow->end();
243 itElement != endElement; ++itElement, i++, pTmp++)
245 for (; i < (*itElement)->col(); i++, pTmp++)
251 for (; i < mNumCols; i++, pTmp++)
262 std::vector< std::vector< CSparseMatrixElement * > >::iterator itRow;
263 std::vector< std::vector< CSparseMatrixElement * > >::iterator endRow;
265 std::vector< CSparseMatrixElement * >::reverse_iterator itElement;
266 std::vector< CSparseMatrixElement * >::reverse_iterator endElement;
268 for (itRow =
mRows.begin(), endRow =
mRows.end(); itRow != endRow; ++itRow)
269 for (itElement = itRow->rbegin(), endElement = itRow->rend(); itElement != endElement; ++itElement)
270 if (fabs(**itElement) < threshold)
271 remove((*itElement)->row(), (*itElement)->col());
288 std::vector< CSparseMatrixElement * >::iterator found
291 if (found !=
mCols[col].end() &&
292 (*found)->row() == row)
return **found;
298 const size_t & col)
const
303 std::vector< CSparseMatrixElement * >::const_iterator found
306 if (found !=
mCols[col].end() &&
307 (*found)->row() == row)
return **found;
319 std::vector< CSparseMatrixElement * >::iterator found
322 if (found !=
mCols[col].end() &&
323 (*found)->row() == row)
return false;
326 mCols[col].insert(found, pElement);
329 mRows[row].insert(found, pElement);
340 std::vector< CSparseMatrixElement * >::iterator found
343 if (found ==
mCols[col].end())
return false;
345 mCols[col].erase(found);
348 mRows[row].erase(found);
361 std::vector< std::vector< CSparseMatrixElement * > >::const_iterator itRow;
362 std::vector< std::vector< CSparseMatrixElement * > >::const_iterator endRow;
364 std::vector< CSparseMatrixElement * >::const_iterator itElement;
365 std::vector< CSparseMatrixElement * >::const_iterator endElement;
367 for (itRow = A.
mRows.begin(), endRow = A.
mRows.end(); itRow != endRow; ++itRow)
369 for (i = 0, itElement = itRow->begin(), endElement = itRow->end();
370 itElement != endElement; ++itElement, i++)
372 for (; i < (*itElement)->col(); i++)
375 os <<
"\t" << **itElement;
390 const size_t & rowIndex):
400 mRowIndex < mpMatrix->
numRows())
410 mpMatrix(src.mpMatrix),
411 mRowIndex(src.mRowIndex),
413 mColumnIndex(src.mColumnIndex),
414 mpColumnIndex(src.mpColumnIndex),
415 mpCurrent(src.mpCurrent)
444 const size_t * pRowIndexEnd = mpMatrix->getRowIndex() + mpMatrix->numNonZeros();
452 size_t index =
mpRowIndex - mpMatrix->getRowIndex();
453 mpCurrent = mpMatrix->getValues() + index;
455 while (*mpColumnIndex <= index) mpColumnIndex++;
457 mColumnIndex = mpColumnIndex - mpMatrix->getColumnStart() - 1;
466 {
return mColumnIndex;}
477 const size_t & columns,
478 const size_t & nonZeros):
481 mpColumnStart(new size_t[columns + 1]),
482 mpRowIndex(nonZeros ? new size_t[nonZeros] : NULL),
483 mpValue(nonZeros ? new
C_FLOAT64[nonZeros] : NULL)
554 std::vector< std::vector< CSparseMatrixElement * > >::const_iterator itCol;
555 std::vector< std::vector< CSparseMatrixElement * > >::const_iterator endCol;
557 std::vector< CSparseMatrixElement * >::const_iterator itElement;
558 std::vector< CSparseMatrixElement * >::const_iterator endElement;
564 itCol != endCol; ++itCol, i++)
568 for (itElement = itCol->begin(), endElement = itCol->end();
569 itElement != endElement; ++itElement, k++, pRowIndex++, pValue++)
571 *pRowIndex = (*itElement)->row();
572 *pValue = **itElement;
590 bool SparseMatrixTest(
const size_t & size,
594 const bool & dgemmFlag,
598 size_t i, j, l, loop = 1;
605 if (Sparseness == 0.0) Sparseness = 4.0 / size;
613 for (i = 0; i < size - 3; i++)
614 for (j = 0; j < size; j++)
620 for (i = 0; i < size; i++)
621 for (j = 0; j < size + 3; j++)
633 std::cout <<
"Memory requirements for sparseness:\t" << Sparseness << std::endl;
636 std::cout <<
"Matrix(" << size <<
"x" << size <<
"):\t" << tmp << std::endl;
639 + 2 * size *
sizeof(std::vector<CSparseMatrixElement *>)
642 std::cout <<
"Sparse(" << size <<
"x" << size <<
"):\t" << tmp2 << std::endl;
643 std::cout <<
"Sparse/Matrix:\t" << tmp2 / tmp << std::endl;
646 + 2 * C.numNonZeros() *
sizeof(
C_FLOAT64)
648 std::cout <<
"CompressedColumnFormat(" << size <<
"x" << size <<
"):\t" << tmp2 << std::endl;
649 std::cout <<
"CompressedColumnFormat/Matrix:\t" << tmp2 / tmp << std::endl << std::endl;
660 for (l = 0; l < loop; l++)
663 const C_FLOAT64 *pTmp1, *pTmp2, *pTmp4, *pTmp5;
667 size_t LDA = M.numCols();
668 size_t LDB = MM.numCols();
671 pEnd1 = pTmp1 + M.
numRows() * LDA;
673 pEnd2 = MM.array() + LDB;
676 for (; pTmp1 < pEnd1; pTmp1 += LDA)
677 for (pTmp2 = MM.
array(); pTmp2 < pEnd2; pTmp2++, pTmp3++)
681 for (pTmp4 = pTmp1, pTmp5 = pTmp2, pEnd4 = pTmp4 + LDA;
682 pTmp4 < pEnd4; pTmp4++, pTmp5 += LDB)
683 * pTmp3 += *pTmp4 * *pTmp5;
689 std::cout <<
"Matrix * Matrix:\t";
690 CPU.print(&std::cout);
692 WALL.print(&std::cout);
693 std::cout << std::endl;
701 for (l = 0; l < loop; l++)
714 M.array(), &k, &Beta, dgemmR.array(), &m);
725 std::cout <<
"dgemm(Matrix, Matrix):\t";
726 CPU.print(&std::cout);
728 WALL.print(&std::cout);
729 std::cout << std::endl;
738 for (l = 0; l < loop; l++)
742 std::vector< std::vector< CSparseMatrixElement * > >::const_iterator itRow;
743 std::vector< std::vector< CSparseMatrixElement * > >::const_iterator endRow;
744 std::vector< CSparseMatrixElement * >::const_iterator itRowElement;
745 std::vector< CSparseMatrixElement * >::const_iterator endRowElement;
746 std::vector< std::vector< CSparseMatrixElement * > >::const_iterator itCol;
747 std::vector< std::vector< CSparseMatrixElement * > >::const_iterator endCol;
748 std::vector< CSparseMatrixElement * >::const_iterator itColElement;
749 std::vector< CSparseMatrixElement * >::const_iterator endColElement;
751 for (itRow = S.getRows().begin(), endRow = S.getRows().end(); itRow != endRow; ++itRow)
753 endRowElement = itRow->end();
755 for (itCol = Ss.getColumns().begin(), endCol = Ss.getColumns().end(); itCol != endCol; ++itCol)
758 itRowElement = itRow->begin();
759 itColElement = itCol->begin();
760 endColElement = itCol->end();
762 while (itRowElement != endRowElement &&
763 itColElement != endColElement)
765 while (itRowElement != endRowElement &&
766 (*itRowElement)->col() < (*itColElement)->row()) ++itRowElement;
768 if (itRowElement == endRowElement)
break;
770 while (itColElement != endColElement &&
771 (*itColElement)->row() < (*itRowElement)->col()) ++itColElement;
773 if (itColElement == endColElement)
break;
775 if ((*itRowElement)->col() != (*itColElement)->row())
continue;
777 Tmp += **itRowElement * **itColElement;
782 if (fabs(Tmp) < SR.getTreshold())
continue;
784 SR.insert((*itRow->begin())->row(), (*itCol->begin())->col(), Tmp);
791 std::cout <<
"Sparse * Sparse:\t";
792 CPU.print(&std::cout);
794 WALL.print(&std::cout);
795 std::cout << std::endl;
810 for (l = 0; l < loop; l++)
815 size_t imax = CR.numRows();
816 size_t jmax = CR.numCols();
817 C_FLOAT64 * pColElement, * pEndColElement;
818 size_t * pColElementRow, * pEndColElementRow;
823 for (j = 0, pColStart = CC.getColumnStart(); j < jmax; j++, pColStart++)
825 for (i = 0; i < imax; i++)
829 itRowElement = C.beginRow(i);
830 pColElement = CC.getValues() + *pColStart;
831 pEndColElement = CC.getValues() + *(pColStart + 1);
832 pColElementRow = CC.getRowIndex() + *pColStart;
833 pEndColElementRow = CC.getRowIndex() + *(pColStart + 1);
835 while (itRowElement != endRowElement &&
836 pColElement != pEndColElement)
838 while (itRowElement != endRowElement &&
841 if (!(itRowElement != endRowElement))
break;
843 while (pColElement != pEndColElement &&
850 if (pColElement == pEndColElement)
break;
854 Tmp += *itRowElement * *pColElement;
860 if (fabs(Tmp) < TmpR.getTreshold())
continue;
862 TmpR.insert(i, j, Tmp);
871 std::cout <<
"Compressed * Compressed:\t";
872 CPU.print(&std::cout);
874 WALL.print(&std::cout);
875 std::cout << std::endl;
884 std::cout << std::endl;
885 std::cout << std::endl;
const size_t & col() const
const C_FLOAT64 & operator=(const C_FLOAT64 &value)
bool resize(const size_t &rows, const size_t &columns)
CSparseMatrix & operator=(const CCompressedColumnFormat &ccf)
const C_FLOAT64 & getTreshold() const
const size_t & numCols() const
static bool compareRow(const CSparseMatrixElement *pLhs, const CSparseMatrixElement *pRhs)
virtual size_t numRows() const
const size_t & numRows() const
std::vector< std::vector< CSparseMatrixElement * > > mRows
bool remove(const size_t &row, const size_t &col)
static CRandom * createGenerator(CRandom::Type type=CRandom::mt19937, unsigned C_INT32 seed=0)
virtual C_FLOAT64 getRandomCC()
static bool compareCol(const CSparseMatrixElement *pLhs, const CSparseMatrixElement *pRhs)
const std::vector< std::vector< CSparseMatrixElement * > > & getColumns() const
std::ostream & operator<<(std::ostream &os, const CSparseMatrix &A)
const std::vector< std::vector< CSparseMatrixElement * > > & getRows() const
int dgemm_(char *transa, char *transb, integer *m, integer *n, integer *k, doublereal *alpha, doublereal *a, integer *lda, doublereal *b, integer *ldb, doublereal *beta, doublereal *c, integer *ldc)
virtual void resize(size_t rows, size_t cols, const bool ©=false)
CSparseMatrix(const size_t &rows=0, const size_t &cols=0)
const size_t & row() const
std::vector< size_t > mColIndex
std::vector< size_t > mRowIndex
bool insert(const size_t &row, const size_t &col, const C_FLOAT64 &value)
std::vector< std::vector< CSparseMatrixElement * > > mCols
bool setTreshold(const C_FLOAT64 &threshold)
CSparseMatrixElement & operator()(const size_t &row, const size_t &col)
virtual size_t numCols() const
CSparseMatrixElement mElement
size_t numNonZeros() const