A triangle counting assignment for A.U.TH Parallel and distributed systems class.
Du kannst nicht mehr als 25 Themen auswählen Themen müssen entweder mit einem Buchstaben oder einer Ziffer beginnen. Sie können Bindestriche („-“) enthalten und bis zu 35 Zeichen lang sein.

vor 4 Jahren
vor 4 Jahren
vor 4 Jahren
vor 4 Jahren
vor 4 Jahren
vor 4 Jahren
vor 4 Jahren
vor 4 Jahren
vor 4 Jahren
123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722
  1. /**
  2. * \file impl.hpp
  3. * \brief Implementation common header file
  4. *
  5. * \author
  6. * Christos Choutouridis AEM:8997
  7. * <cchoutou@ece.auth.gr>
  8. */
  9. #ifndef IMPL_HPP_
  10. #define IMPL_HPP_
  11. #include <type_traits>
  12. #include <utility>
  13. #include <algorithm>
  14. #include <vector>
  15. #include <tuple>
  16. #include <cstddef>
  17. #include <iostream>
  18. #include <fstream>
  19. using std::size_t;
  20. /*
  21. * Small helper to strip types
  22. */
  23. template<typename T>
  24. struct remove_cvref {
  25. typedef std::remove_cv_t<std::remove_reference_t<T>> type;
  26. };
  27. template<typename T>
  28. using remove_cvref_t = typename remove_cvref<T>::type;
  29. /*!
  30. * Enumerator to denote the storage type of the array to use.
  31. */
  32. enum class MatrixType {
  33. FULL, /*!< Matrix is asymmetric */
  34. SYMMETRIC, /*!< Matrix is symmetric */
  35. };
  36. /*
  37. * Forward type declerations
  38. */
  39. template<typename DataType, typename IndexType, MatrixType Type = MatrixType::SYMMETRIC> struct Matrix;
  40. template<typename DataType, typename IndexType, MatrixType Type = MatrixType::SYMMETRIC> struct SpMat;
  41. template<typename DataType, typename IndexType> struct SpMatCol;
  42. template<typename DataType, typename IndexType> struct SpMatRow;
  43. template<typename DataType, typename IndexType> struct SpMatVal;
  44. /*!
  45. * 2D-array wrapper for v1 and v2 part of the exercise used as RAII
  46. * and copy-prevention.
  47. *
  48. * This is a very thin abstraction layer over a native array.
  49. * This is tested using compiler explorer and our template produce
  50. * almost identical assembly.
  51. *
  52. * The penalty hit we have is due to the fact that we use a one dimension array
  53. * and we have to calculate the actual position from an (i,j) pair.
  54. * The use of 1D array was our intention from the beginning, so the penalty
  55. * was pretty much unavoidable.
  56. *
  57. * \tparam DataType The underling data type of the array
  58. * \tparam Type The storage type of the array
  59. * \arg FULL For full matrix
  60. * \arg SYMMETRIC For symmetric matrix (we use only the lower part)
  61. */
  62. template<typename DataType, typename IndexType, MatrixType Type>
  63. struct Matrix {
  64. using dataType = DataType; //!< meta:export of underling data type
  65. using indexType = IndexType; //!< meta:export of underling index type
  66. static constexpr MatrixType matrixType = Type; //!< export of array type
  67. /*!
  68. * \name Obj lifetime
  69. */
  70. //! @{
  71. //! Constructor using data type and size
  72. Matrix(DataType* t, IndexType s) noexcept : m_(t), size_(s) { }
  73. //! RAII deleter
  74. ~Matrix() { delete m_; }
  75. //! move ctor
  76. Matrix(Matrix&& a) noexcept { a.swap(*this); }
  77. //! move
  78. Matrix& operator=(Matrix&& a) noexcept { a.swap(*this); return *this; }
  79. Matrix(const Matrix& a) = delete; //!< No copy ctor
  80. Matrix& operator=(const Matrix& a) = delete; //!< No copy
  81. //! @}
  82. //! \name Data exposure
  83. //! @{
  84. //! memory capacity of the matrix with diagonal data
  85. template<MatrixType AT=Type> std::enable_if_t<AT==MatrixType::SYMMETRIC, IndexType>
  86. static constexpr capacity(IndexType N) noexcept { return (N+1)*N/2; }
  87. //! memory capacity for full matrix
  88. template<MatrixType AT=Type> std::enable_if_t<AT==MatrixType::FULL, IndexType>
  89. static constexpr capacity(IndexType N) noexcept { return N*N; }
  90. //! Get the size of each dimension
  91. IndexType size() noexcept { return size_; }
  92. /*
  93. * virtual 2D accessors
  94. */
  95. template<MatrixType AT=Type>
  96. std::enable_if_t <AT==MatrixType::SYMMETRIC, DataType>
  97. get (IndexType i, IndexType j) {
  98. if (i<j) return m_[j*(j+1)/2 + i]; // Upper, use opposite index
  99. else return m_[i*(i+1)/2 + j]; // Lower, use our notation
  100. }
  101. template<MatrixType AT=Type>
  102. std::enable_if_t <AT==MatrixType::FULL, DataType>
  103. get(IndexType i, IndexType j) {
  104. return m_[i*size_ + j];
  105. }
  106. template<MatrixType AT=Type>
  107. std::enable_if_t <AT==MatrixType::SYMMETRIC, DataType>
  108. set(DataType v, IndexType i, IndexType j) {
  109. if (i<j) return m_[j*(j+1)/2 + i] =v; // Upper, use opposite index
  110. else return m_[i*(i+1)/2 + j] =v; // Lower, use our notation
  111. }
  112. template<MatrixType AT=Type>
  113. std::enable_if_t <AT==MatrixType::FULL, DataType>
  114. set(DataType v, IndexType i, IndexType j) {
  115. return m_[i*size_ + j] =v;
  116. }
  117. DataType operator()(IndexType i, IndexType j) { return get(i, j); }
  118. // a basic serial iterator support
  119. DataType* begin() noexcept { return m_; }
  120. DataType* end() noexcept { return m_ + Matrix::capacity(size_); }
  121. //! @}
  122. /*!
  123. * \name Safe iteration API
  124. *
  125. * This api automates the iteration over the array based on
  126. * MatrixType
  127. */
  128. //! @{
  129. template<typename F, typename... Args>
  130. void for_each_in (IndexType begin, IndexType end, F&& lambda, Args&&... args) {
  131. for (IndexType it=begin ; it<end ; ++it) {
  132. std::forward<F>(lambda)(std::forward<Args>(args)..., it);
  133. }
  134. }
  135. //! @}
  136. // move helper
  137. void swap(Matrix& src) noexcept {
  138. std::swap(m_, src.m_);
  139. std::swap(size_, src.size_);
  140. }
  141. private:
  142. DataType* m_; //!< Pointer to actual data.
  143. IndexType size_; //!< the virtual size of each dimension.
  144. };
  145. /*!
  146. * RAII allocation helper, the smart_ptr way.
  147. */
  148. template<typename DataType, typename IndexType, MatrixType Type = MatrixType::SYMMETRIC>
  149. Matrix<DataType, IndexType, Type> make_Matrix(IndexType s) {
  150. return Matrix<DataType, IndexType, Type>(new DataType[Matrix<DataType, IndexType, Type>::capacity(s)], s);
  151. }
  152. /**
  153. * A simple sparse matrix implementation.
  154. *
  155. * We use CSC format and provide get/set functionalities for each (i,j) item
  156. * on the matrix. We also provide a () overload using a proxy SpMatVal object.
  157. * This way the user can:
  158. * \code
  159. * auto v = A(3,4);
  160. * A(3, 4) = 7;
  161. * \endcode
  162. *
  163. * We also provide getCol() and getRow() functions witch return a viewer/iterator to rows and
  164. * columns of the matrix. In the case of a symmetric matrix instead of a row we return the
  165. * equivalent column. This way we gain speed due to CSC format nature.
  166. *
  167. * @tparam DataType The type for values
  168. * @tparam IndexType The type for indexes
  169. * @tparam Type The Matrix type (FULL or SYMMETRIC)
  170. */
  171. template<typename DataType, typename IndexType, MatrixType Type>
  172. struct SpMat {
  173. using dataType = DataType; //!< meta:export of underling data type
  174. using indexType = IndexType; //!< meta:export of underling index type
  175. static constexpr MatrixType matrixType = Type; //!< export of array type
  176. friend struct SpMatCol<DataType, IndexType>;
  177. friend struct SpMatRow<DataType, IndexType>;
  178. friend struct SpMatVal<DataType, IndexType>;
  179. /*!
  180. * \name Obj lifetime
  181. */
  182. //! @{
  183. //! Default ctor with optional memory allocations
  184. SpMat(IndexType n=IndexType{}, IndexType nnz=IndexType{}) :
  185. values(nnz, DataType{}),
  186. rows(nnz, IndexType{}),
  187. col_ptr((n)? n+1:2, IndexType{}),
  188. N(n),
  189. NNZ(nnz) { }
  190. //! A ctor using csc array data
  191. SpMat(IndexType n, IndexType nnz, const IndexType* row, const IndexType* col) :
  192. values(nnz, 1),
  193. rows(row, row+nnz),
  194. col_ptr(col, col+n+1),
  195. N(n),
  196. NNZ(nnz) { }
  197. //! ctor using csc array data with value array
  198. SpMat(IndexType n, IndexType nnz, const DataType* v, const IndexType* row, const IndexType* col) :
  199. values(v, v+nnz),
  200. rows(row, row+nnz),
  201. col_ptr(col, col+n+1),
  202. N(n),
  203. NNZ(nnz) { }
  204. //! ctor vectors of row/col and default value for values array
  205. SpMat(IndexType n, IndexType nnz, const DataType v, const std::vector<IndexType>& row, const std::vector<IndexType>& col) :
  206. values(nnz, v),
  207. rows (row),
  208. col_ptr(col),
  209. N(n),
  210. NNZ(nnz) { }
  211. //! move ctor
  212. SpMat(SpMat&& a) noexcept { moves(std::move(a)); }
  213. //! move
  214. SpMat& operator=(SpMat&& a) noexcept { moves(std::move(a)); return *this; }
  215. SpMat(const SpMat& a) = delete; //!< make sure there are no copies
  216. SpMat& operator=(const SpMat& a) = delete; //!< make sure there are no copies
  217. //! @}
  218. //! \name Data exposure
  219. //! @{
  220. //! \return the dimension of the matrix
  221. IndexType size() noexcept { return N; }
  222. //! After construction size configuration tool
  223. IndexType size(IndexType n) {
  224. col_ptr.resize(n+1);
  225. return N = n;
  226. }
  227. //! \return the NNZ of the matrix
  228. IndexType capacity() noexcept { return NNZ; }
  229. //! After construction NNZ size configuration tool
  230. IndexType capacity(IndexType nnz) {
  231. values.reserve(nnz);
  232. rows.reserve(nnz);
  233. return NNZ;
  234. }
  235. // getters for row arrays of the struct (unused)
  236. std::vector<DataType>& getValues() noexcept { return values; }
  237. std::vector<IndexType>& getRows() noexcept { return rows; }
  238. std::vector<IndexType>& getCols() noexcept { return col_ptr; }
  239. /*!
  240. * Return a proxy SpMatVal object with read and write capabilities.
  241. * @param i The row number
  242. * @param j The column number
  243. * @return tHE SpMatVal object
  244. */
  245. SpMatVal<DataType, IndexType> operator()(IndexType i, IndexType j) {
  246. return SpMatVal<DataType, IndexType>(this, get(i, j), i, j);
  247. }
  248. /*!
  249. * A read item functionality using binary search to find the correct row
  250. *
  251. * @param i The row number
  252. * @param j The column number
  253. * @return The value of the item or DataType{} if is not present.
  254. */
  255. DataType get(IndexType i, IndexType j) {
  256. IndexType idx; bool found;
  257. std::tie(idx, found) =find_idx(rows, col_ptr[j], col_ptr[j+1], i);
  258. return (found) ? values[idx] : 0;
  259. }
  260. /*!
  261. * A read item functionality using linear search to find the correct row
  262. *
  263. * @param i The row number
  264. * @param j The column number
  265. * @return The value of the item or DataType{} if is not present.
  266. */
  267. DataType get_lin(IndexType i, IndexType j) {
  268. IndexType idx; bool found;
  269. std::tie(idx, found) =find_lin_idx(rows, col_ptr[j], col_ptr[j+1], i);
  270. return (found) ? values[idx] : 0;
  271. }
  272. /*!
  273. * A write item functionality.
  274. *
  275. * First we search if the matrix has already a value in (i, j) position.
  276. * If so we just change it to a new value. If not we add the item on the matrix.
  277. *
  278. * @note
  279. * We don't increase the NNZ value of the struct. We expect the user has already
  280. * change the NNZ value to the right one using @see capacity() function.
  281. *
  282. * @param i The row number
  283. * @param j The column number
  284. * @return The new value of the item .
  285. */
  286. DataType set(DataType v, IndexType i, IndexType j) {
  287. IndexType idx; bool found;
  288. std::tie(idx, found) = find_lin_idx(rows, col_ptr[j], col_ptr[j+1], i);
  289. if (found)
  290. return values[idx] = v; // we don't change NNZ even if we write "0"
  291. else {
  292. values.insert(values.begin()+idx, v);
  293. rows.insert(rows.begin()+idx, i);
  294. std::transform(col_ptr.begin()+j+1, col_ptr.end(), col_ptr.begin()+j+1, [](IndexType it) {
  295. return ++it;
  296. });
  297. ++NNZ; // we increase the NNZ even if we write "0"
  298. return v;
  299. }
  300. }
  301. /*!
  302. * Get a view of a CSC column
  303. * @param j The column to get
  304. * @return The SpMatCol object @see SpMatCol
  305. */
  306. SpMatCol<DataType, IndexType> getCol(IndexType j) {
  307. return SpMatCol<DataType, IndexType>(this, col_ptr[j], col_ptr[j+1]);
  308. }
  309. /*!
  310. * Get a view of a CSC row
  311. *
  312. * In case of a SYMMETRIC matrix we can return a column instead.
  313. *
  314. * @param j The row to get
  315. * @return The SpMatCol object @see SpMatCol
  316. */
  317. template<MatrixType AT= Type>
  318. std::enable_if_t<AT==MatrixType::SYMMETRIC, SpMatCol<DataType, IndexType>>
  319. getRow(IndexType i) {
  320. return getCol(i);
  321. }
  322. /*!
  323. * Get a view of a CSC row
  324. *
  325. * @param j The row to get
  326. * @return The SpMatRow object @see SpMatRow
  327. */
  328. template<MatrixType AT= Type>
  329. std::enable_if_t<AT==MatrixType::FULL, SpMatCol<DataType, IndexType>>
  330. getRow(IndexType i) {
  331. return SpMatRow<DataType, IndexType>(this, i);
  332. }
  333. // values only iterator support
  334. DataType* begin() noexcept { return values.begin(); }
  335. DataType* end() noexcept { return values.end(); }
  336. //! @}
  337. //! A small iteration helper
  338. template<typename F, typename... Args>
  339. void for_each_in (IndexType begin, IndexType end, F&& lambda, Args&&... args) {
  340. for (IndexType it=begin ; it<end ; ++it) {
  341. std::forward<F>(lambda)(std::forward<Args>(args)..., it);
  342. }
  343. }
  344. // friend operations for printing
  345. template<typename D, typename I, MatrixType T> friend void print(SpMat<D, I, T>& mat);
  346. template<typename D, typename I, MatrixType T> friend void print_dense(SpMat<D, I, T>& mat);
  347. private:
  348. /*!
  349. * A small binary search implementation using index for begin-end instead of iterators.
  350. *
  351. * \param v Reference to vector to search
  352. * \param begin The vector's index to begin
  353. * \param end The vector's index to end
  354. * \param match What to search
  355. * @return The index of the item or end on failure.
  356. */
  357. std::pair<IndexType, bool> find_idx(const std::vector<IndexType>& v, IndexType begin, IndexType end, IndexType match) {
  358. IndexType b = begin, e = end-1;
  359. while (true) {
  360. IndexType m = (b+e)/2;
  361. if (v[m] == match) return std::make_pair(m, true);
  362. else if (b >= e) return std::make_pair(end, false);
  363. else {
  364. if (v[m] < match) b = m +1;
  365. else e = m -1;
  366. }
  367. }
  368. return std::make_pair(end, false);
  369. }
  370. /*!
  371. * find helper for set using index for begin-end instead of iterators.
  372. *
  373. * We search for the item or a place to add the item.
  374. * So we return the index if we find it. Otherwise we return the place to add it.
  375. *
  376. * \param v Reference to vector to search
  377. * \param begin The vector's index to begin
  378. * \param end The vector's index to end
  379. * \param match What to search
  380. */
  381. std::pair<IndexType, bool> find_lin_idx(const std::vector<IndexType>& v, IndexType begin, IndexType end, IndexType match) {
  382. for ( ; begin < end ; ++begin) {
  383. if (match == v[begin]) return std::make_pair(begin, true);
  384. else if (match < v[begin]) return std::make_pair(begin, false);
  385. }
  386. return std::make_pair(end, false);
  387. }
  388. // move helper
  389. void moves(SpMat&& src) noexcept {
  390. values = std::move(src.values);
  391. rows = std::move(src.rows);
  392. col_ptr = std::move(src.col_ptr);
  393. N = std::move(src.N); // redundant for primitives
  394. NNZ = std::move(src.NNZ); //
  395. }
  396. //! \name Data
  397. //! @{
  398. std::vector<DataType> values {}; //!< vector to store the values of the matrix
  399. std::vector<IndexType> rows{}; //!< vector to store the row information
  400. std::vector<IndexType> col_ptr{1,0}; //!< vector to stor the column pointers
  401. IndexType N{0}; //!< The dimension of the matrix (square)
  402. IndexType NNZ{0}; //!< The NNZ (capacity of the matrix)
  403. //! @}
  404. };
  405. /*!
  406. * A view/iterator hybrid object for SpMat columns.
  407. *
  408. * This object provides access to a column of a SpMat. The public functionalities
  409. * allow data access using indexes instead of iterators. We prefer indexes over iterators
  410. * because we can apply the same index to different inner vector of SpMat without conversion.
  411. *
  412. * @tparam DataType
  413. * @tparam IndexType
  414. */
  415. template<typename DataType, typename IndexType>
  416. struct SpMatCol {
  417. using owner_t = SpMat<DataType, IndexType>;
  418. /*!
  419. * ctor using column pointers for begin-end. own is pointer to SpMat.
  420. */
  421. SpMatCol(owner_t* own, const IndexType begin, const IndexType end) noexcept :
  422. owner_(own), index_(begin), begin_(begin), end_(end) {
  423. vindex_ = vIndexCalc(index_);
  424. }
  425. SpMatCol() = default;
  426. SpMatCol(const SpMatCol&) = delete; //!< make sure there are no copies
  427. SpMatCol& operator=(const SpMatCol&)= delete; //!< make sure there are no copies
  428. SpMatCol(SpMatCol&&) = default;
  429. SpMatCol& operator=(SpMatCol&&) = default;
  430. //! a simple dereference operator, like an iterator
  431. DataType operator* () {
  432. return get();
  433. }
  434. //! Increment operator acts on index(), like an iterator
  435. SpMatCol& operator++ () { advance(); return *this; }
  436. SpMatCol& operator++ (int) { SpMatCol& p = *this; advance(); return p; }
  437. //! () operator acts as member access (like a view)
  438. DataType operator()(IndexType x) {
  439. return (x == index())? get() : DataType{};
  440. }
  441. //! = operator acts as member assignment (like a view)
  442. DataType operator= (DataType v) { return owner_->values[index_] = v; }
  443. // iterator like handlers
  444. // these return a virtual index value based on the items position on the full matrix
  445. // but the move of the index is just a ++ away.
  446. IndexType index() noexcept { return vindex_; }
  447. const IndexType index() const noexcept { return vindex_; }
  448. IndexType begin() noexcept { return vIndexCalc(begin_); }
  449. const IndexType begin() const noexcept { return vIndexCalc(begin_); }
  450. IndexType end() noexcept { return owner_->N; }
  451. const IndexType end() const noexcept { return owner_->N; }
  452. /*!
  453. * Multiplication operator
  454. * @tparam C Universal reference for the type right half site column
  455. *
  456. * @param c The right hand site matrix
  457. * @return The value of the inner product of two vectors
  458. */
  459. template <typename C>
  460. DataType operator* (C&& c) {
  461. static_assert(std::is_same<remove_cvref_t<C>, SpMatCol<DataType, IndexType>>(), "");
  462. DataType v{};
  463. while (index() != end() && c.index() != c.end()) {
  464. if (index() < c.index()) advance();
  465. else if (index() > c.index()) ++c;
  466. else { //index() == c.index()
  467. v += get() * *c;
  468. ++c;
  469. advance();
  470. }
  471. }
  472. return v;
  473. }
  474. private:
  475. //! small tool to increase the index pointers to SpMat matrix
  476. void advance() noexcept {
  477. ++index_;
  478. vindex_ = vIndexCalc(index_);
  479. }
  480. //! tool to translate between col_ptr indexes and SpMat "virtual" full matrix indexes
  481. IndexType vIndexCalc(IndexType idx) {
  482. return (idx < end_) ? owner_->rows[idx] : end();
  483. }
  484. //! small get tool
  485. DataType get() { return owner_->values[index_]; }
  486. owner_t* owner_ {nullptr}; //!< Pointer to owner SpMat. SpMatCol is just a view
  487. IndexType vindex_ {IndexType{}}; //!< Virtual index of full matrix
  488. IndexType index_ {IndexType{}}; //!< index to SpMat::rows
  489. IndexType begin_ {IndexType{}}; //!< beginning index of the column in SpMat::rows
  490. IndexType end_ {IndexType{}}; //!< ending index of the column in SpMat::rows
  491. };
  492. /*!
  493. * A view/iterator hybrid object for SpMat rows.
  494. *
  495. * This object provides access to a column of a SpMat. The public functionalities
  496. * allow data access using indexes instead of iterators. We prefer indexes over iterators
  497. * because we can apply the same index to different inner vector of SpMat without conversion.
  498. *
  499. * @tparam DataType
  500. * @tparam IndexType
  501. */
  502. template<typename DataType, typename IndexType>
  503. struct SpMatRow {
  504. using owner_t = SpMat<DataType, IndexType>;
  505. /*!
  506. * ctor using virtual full matrix row index. own is pointer to SpMat.
  507. */
  508. SpMatRow(owner_t* own, const IndexType row) noexcept :
  509. owner_(own), vindex_(IndexType{}), row_(row), index_(IndexType{}),
  510. begin_(IndexType{}), end_(owner_->NNZ) {
  511. // place begin
  512. while(begin_ != end_ && owner_->rows[begin_] != row_)
  513. ++begin_;
  514. // place index_ and vindex_
  515. if (owner_->rows[index_] != row_)
  516. advance();
  517. }
  518. SpMatRow() = default;
  519. SpMatRow(const SpMatRow&) = delete; //!< make sure there are no copies
  520. SpMatRow& operator=(const SpMatRow&)= delete; //!< make sure there are no copies
  521. SpMatRow(SpMatRow&&) = default;
  522. SpMatRow& operator=(SpMatRow&&) = default;
  523. //! a simple dereference operator, like an iterator
  524. DataType operator* () {
  525. return get();
  526. }
  527. //! Increment operator acts on index(), like an iterator
  528. //! here the increment is a O(N) process.
  529. SpMatRow& operator++ () { advance(); return *this; }
  530. SpMatRow& operator++ (int) { SpMatRow& p = *this; advance(); return p; }
  531. //! () operator acts as member access (like a view)
  532. DataType operator()(IndexType x) {
  533. return (x == index())? get() : DataType{};
  534. }
  535. //! = operator acts as member assignment (like a view)
  536. DataType operator= (DataType v) { return owner_->values[index_] = v; }
  537. // iterator like handlers
  538. // these return a virtual index value based on the items position on the full matrix
  539. // but the move of the index is just a ++ away.
  540. IndexType index() noexcept { return vindex_; }
  541. const IndexType index() const noexcept { return vindex_; }
  542. IndexType begin() noexcept { return vIndexCalc(begin_); }
  543. const IndexType begin() const noexcept { return vIndexCalc(begin_); }
  544. IndexType end() noexcept { return owner_->N; }
  545. const IndexType end() const noexcept { return owner_->N; }
  546. /*!
  547. * Multiplication operator
  548. * @tparam C Universal reference for the type right half site column
  549. *
  550. * @param c The right hand site matrix
  551. * @return The value of the inner product of two vectors
  552. */
  553. template <typename C>
  554. DataType operator* (C&& c) {
  555. static_assert(std::is_same<remove_cvref_t<C>, SpMatCol<DataType, IndexType>>(), "");
  556. DataType v{};
  557. while (index() != end() && c.index() != c.end()) {
  558. if (index() < c.index()) advance();
  559. else if (index() > c.index()) ++c;
  560. else { //index() == c.index()
  561. v += get()* *c;
  562. ++c;
  563. advance();
  564. }
  565. }
  566. return v;
  567. }
  568. private:
  569. //! small tool to increase the index pointers to SpMat matrix
  570. //! We have to search the entire rows vector in SpMat to find the next
  571. //! virtual row position.
  572. //! time complexity O(N)
  573. void advance() noexcept {
  574. do
  575. ++index_;
  576. while(index_ != end_ && owner_->rows[index_] != row_);
  577. vindex_ = vIndexCalc(index_);
  578. }
  579. //! tool to translate between col_ptr indexes and SpMat "virtual" full matrix indexes
  580. IndexType vIndexCalc(IndexType idx) {
  581. for(IndexType i =0 ; i<(owner_->N+1) ; ++i)
  582. if (idx < owner_->col_ptr[i])
  583. return i-1;
  584. return end();
  585. }
  586. //! small get tool
  587. DataType get() { return owner_->values[index_]; }
  588. owner_t* owner_ {nullptr}; //!< Pointer to owner SpMat. SpMatCol is just a view
  589. IndexType vindex_ {IndexType{}}; //!< Virtual index of full matrix
  590. IndexType row_ {IndexType{}}; //!< The virtual full matrix row of the object
  591. IndexType index_ {IndexType{}}; //!< index to SpMat::rows
  592. IndexType begin_ {IndexType{}}; //!< beginning index of the column in SpMat::rows
  593. IndexType end_ {IndexType{}}; //!< ending index of the column in SpMat::rows
  594. };
  595. /*!
  596. * A proxy SpMat value object/view.
  597. *
  598. * This object acts as proxy to provide read/write access to an SpMat item.
  599. *
  600. * @tparam DataType The type of the values of the SpMat matrix
  601. * @tparam IndexType The type of the indexes of the SpMat matrix
  602. */
  603. template<typename DataType, typename IndexType>
  604. struct SpMatVal {
  605. using owner_t = SpMat<DataType, IndexType>;
  606. //!< ctor using all value-row-column data, plus a pointer to owner SpMat object
  607. SpMatVal(owner_t* own, DataType v, IndexType i, IndexType j) :
  608. owner_(own), v_(v), i_(i), j_(j) { }
  609. SpMatVal() = default;
  610. SpMatVal(const SpMatVal&) = delete; //!< make sure there are no copies
  611. SpMatVal& operator=(const SpMatVal&) = delete; //!< make sure there are no copies
  612. SpMatVal(SpMatVal&&) = default;
  613. SpMatVal& operator=(SpMatVal&&) = default;
  614. //! Operator to return the DataType value implicitly
  615. operator DataType() { return v_; }
  616. //! Operator to write back to owner the assigned value
  617. //! for ex: A(2,3) = 5;
  618. SpMatVal& operator=(DataType v) {
  619. v_ = v;
  620. owner_->set(v_, i_, j_);
  621. return *this;
  622. }
  623. private:
  624. owner_t* owner_{nullptr}; //!< Pointer to owner SpMat. SpMatVal is just a view.
  625. DataType v_{DataType{}}; //!< The value of the row-column pair (for speed)
  626. IndexType i_{IndexType{}}; //!< The row
  627. IndexType j_{IndexType{}}; //!< the column
  628. };
  629. //! enumerator for input matrix type.
  630. enum class InputMatrix{ UNSPECIFIED, GENERATE, MTX };
  631. //! enumerator for output handling
  632. enum class OutputMode{ STD, FILE };
  633. /*!
  634. * Session option for each invocation of the executable
  635. */
  636. struct session_t {
  637. InputMatrix inputMatrix {InputMatrix::UNSPECIFIED}; //!< Source of the matrix
  638. std::ifstream mtxFile {}; //!< matrix file in MatrixMarket format
  639. std::size_t gen_size {}; //!< size of the matrix if we select to generate a random one
  640. double gen_prob {}; //!< probability of the binomial distribution for the matrix
  641. //!< if we generate one
  642. OutputMode outputMode {OutputMode::STD}; //!< Type of the output file
  643. std::ofstream outFile {}; //!< File to use for output
  644. std::size_t max_threads {}; //!< Maximum threads to use
  645. std::size_t repeat {1}; //!< How many times we execute the calculations part
  646. bool timing {false}; //!< Enable timing prints of the program
  647. bool verbose {false}; //!< Flag to enable verbose output to stdout
  648. #if CODE_VERSION == 3
  649. bool makeSymmetric {false}; //!< symmetric matrix creation flag (true by default)
  650. #else
  651. bool makeSymmetric {true}; //!< symmetric matrix creation flag (true by default)
  652. #endif
  653. bool validate_mtx {false}; //!< Flag to request mtx input data triangular validation.
  654. bool print_count {false}; //!< Flag to request total count printing
  655. bool mtx_print {false}; //!< matrix print flag
  656. std::size_t mtx_print_size {}; //!< matrix print size
  657. bool dynamic {false}; //!< Selects dynamic scheduling for OpenMP and pthreads.
  658. };
  659. extern session_t session;
  660. #endif /* IMPL_HPP_ */