Loading...
Searching...
No Matches
base_matrix_with_column_compression.h
Go to the documentation of this file.
1/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
2 * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
3 * Author(s): Hannah Schreiber
4 *
5 * Copyright (C) 2022-24 Inria
6 *
7 * Modification(s):
8 * - YYYY/MM Author: Description of the modification
9 */
10
17#ifndef PM_BASE_MATRIX_COMPRESSION_H
18#define PM_BASE_MATRIX_COMPRESSION_H
19
20#include <iostream> //print() only
21#include <vector>
22#include <utility> //std::swap, std::move & std::exchange
23
24#include <boost/intrusive/set.hpp>
25#include <boost/pending/disjoint_sets.hpp>
26
27#include <gudhi/Simple_object_pool.h>
28
29namespace Gudhi {
30namespace persistence_matrix {
31
44template <class Master_matrix>
45class Base_matrix_with_column_compression : protected Master_matrix::Matrix_row_access_option
46{
47 public:
48 using index = typename Master_matrix::index;
49 using dimension_type = typename Master_matrix::dimension_type;
53 using Field_operators = typename Master_matrix::Field_operators;
54 using Field_element_type = typename Master_matrix::element_type;
55 using Row_type = typename Master_matrix::Row_type;
57 using Cell_constructor = typename Master_matrix::Cell_constructor;
58 using Column_settings = typename Master_matrix::Column_settings;
65 : public Master_matrix::Column_type,
66 public boost::intrusive::set_base_hook<boost::intrusive::link_mode<boost::intrusive::normal_link> >
67 {
68 public:
69 using Base = typename Master_matrix::Column_type;
70
71 Column_type(Column_settings* colSettings = nullptr)
72 : Base(colSettings) {}
73 template <class Container_type>
74 Column_type(const Container_type& nonZeroRowIndices, Column_settings* colSettings)
75 : Base(nonZeroRowIndices, colSettings) {}
76 template <class Container_type, class Row_container_type>
77 Column_type(index columnIndex, const Container_type& nonZeroRowIndices, Row_container_type& rowContainer,
78 Column_settings* colSettings)
79 : Base(columnIndex, nonZeroRowIndices, rowContainer, colSettings) {}
80 template <class Container_type>
81 Column_type(const Container_type& nonZeroRowIndices, dimension_type dimension, Column_settings* colSettings)
82 : Base(nonZeroRowIndices, dimension, colSettings) {}
83 template <class Container_type, class Row_container_type>
84 Column_type(index columnIndex, const Container_type& nonZeroRowIndices, dimension_type dimension,
85 Row_container_type& rowContainer, Column_settings* colSettings)
86 : Base(columnIndex, nonZeroRowIndices, dimension, rowContainer, colSettings) {}
87 Column_type(const Column_type& column, Column_settings* colSettings = nullptr)
88 : Base(static_cast<const Base&>(column), colSettings) {}
89 template <class Row_container_type>
90 Column_type(const Column_type& column, index columnIndex, Row_container_type& rowContainer,
91 Column_settings* colSettings = nullptr)
92 : Base(static_cast<const Base&>(column), columnIndex, rowContainer, colSettings) {}
93 Column_type(Column_type&& column) noexcept : Base(std::move(static_cast<Base&>(column))) {}
94
95 //TODO: is it possible to make this work?
96 // template <class... U>
97 // Column_type(U&&... u) : Base(std::forward<U>(u)...) {}
98
99 index get_rep() const { return rep_; }
100 void set_rep(const index& rep) { rep_ = rep; }
101
102 struct Hasher {
103 size_t operator()(const Column_type& c) const { return std::hash<Base>()(static_cast<Base>(c)); }
104 };
105
106 private:
107 index rep_;
108 };
109
129 template <class Container_type>
130 Base_matrix_with_column_compression(const std::vector<Container_type>& columns,
131 Column_settings* colSettings);
139 Base_matrix_with_column_compression(unsigned int numberOfColumns,
140 Column_settings* colSettings);
151 Column_settings* colSettings = nullptr);
162
171 template <class Container_type>
172 void insert_column(const Container_type& column);
181 template <class Boundary_type>
182 void insert_boundary(const Boundary_type& boundary, dimension_type dim = -1);
193 const Column_type& get_column(index columnIndex);
203 const Row_type& get_row(index rowIndex) const;
217 void erase_empty_row(index rowIndex);
218
225
237 template <class Cell_range_or_column_index>
238 void add_to(const Cell_range_or_column_index& sourceColumn, index targetColumnIndex);
252 template <class Cell_range_or_column_index>
253 void multiply_target_and_add_to(const Cell_range_or_column_index& sourceColumn,
254 const Field_element_type& coefficient,
255 index targetColumnIndex);
269 template <class Cell_range_or_column_index>
270 void multiply_source_and_add_to(const Field_element_type& coefficient,
271 const Cell_range_or_column_index& sourceColumn,
272 index targetColumnIndex);
273
282 bool is_zero_cell(index columnIndex, index rowIndex);
290 bool is_zero_column(index columnIndex);
291
298 void reset(Column_settings* colSettings) {
299 columnToRep_.clear_and_dispose(delete_disposer(this));
300 columnClasses_ = boost::disjoint_sets_with_storage<>();
301 repToColumn_.clear();
302 nextColumnIndex_ = 0;
303 colSettings_ = colSettings;
304 }
305
314 matrix1.columnToRep_.swap(matrix2.columnToRep_);
315 std::swap(matrix1.columnClasses_, matrix2.columnClasses_);
316 matrix1.repToColumn_.swap(matrix2.repToColumn_);
317 std::swap(matrix1.nextColumnIndex_, matrix2.nextColumnIndex_);
318 std::swap(matrix1.colSettings_, matrix2.colSettings_);
319 std::swap(matrix1.columnPool_, matrix2.columnPool_);
320
321 if constexpr (Master_matrix::Option_list::has_row_access) {
322 swap(static_cast<typename Master_matrix::Matrix_row_access_option&>(matrix1),
323 static_cast<typename Master_matrix::Matrix_row_access_option&>(matrix2));
324 }
325 }
326
327 void print(); // for debug
328
329 private:
333 struct delete_disposer {
334 delete_disposer(Base_matrix_with_column_compression* matrix) : matrix_(matrix) {}
335
336 void operator()(Column_type* delete_this) { matrix_->columnPool_->destroy(delete_this); }
337
338 Base_matrix_with_column_compression* matrix_;
339 };
340
341 using ra_opt = typename Master_matrix::Matrix_row_access_option;
342 using col_dict_type = boost::intrusive::set<Column_type, boost::intrusive::constant_time_size<false> >;
343
344 col_dict_type columnToRep_;
345 boost::disjoint_sets_with_storage<> columnClasses_;
347 std::vector<Column_type*> repToColumn_;
349 index nextColumnIndex_;
350 Column_settings* colSettings_;
355 std::unique_ptr<Simple_object_pool<Column_type> > columnPool_;
356 inline static const Column_type empty_column_;
358 void _insert_column(index columnIndex);
359 void _insert_double_column(index columnIndex, typename col_dict_type::iterator& doubleIt);
360};
361
362template <class Master_matrix>
364 Column_settings* colSettings)
365 : ra_opt(), nextColumnIndex_(0), colSettings_(colSettings), columnPool_(new Simple_object_pool<Column_type>)
366{}
367
368template <class Master_matrix>
369template <class Container_type>
371 const std::vector<Container_type>& columns, Column_settings* colSettings)
372 : ra_opt(columns.size()),
373 columnClasses_(columns.size()),
374 repToColumn_(columns.size(), nullptr),
375 nextColumnIndex_(0),
376 colSettings_(colSettings),
377 columnPool_(new Simple_object_pool<Column_type>)
378{
379 for (const Container_type& c : columns) {
380 insert_column(c);
381 }
382}
383
384template <class Master_matrix>
386 unsigned int numberOfColumns, Column_settings* colSettings)
387 : ra_opt(numberOfColumns),
388 columnClasses_(numberOfColumns),
389 repToColumn_(numberOfColumns, nullptr),
390 nextColumnIndex_(0),
391 colSettings_(colSettings),
392 columnPool_(new Simple_object_pool<Column_type>)
393{}
394
395template <class Master_matrix>
397 const Base_matrix_with_column_compression& matrixToCopy, Column_settings* colSettings)
398 : ra_opt(static_cast<const ra_opt&>(matrixToCopy)),
399 columnClasses_(matrixToCopy.columnClasses_),
400 repToColumn_(matrixToCopy.repToColumn_.size(), nullptr),
401 nextColumnIndex_(0),
402 colSettings_(colSettings == nullptr ? matrixToCopy.colSettings_ : colSettings),
403 columnPool_(new Simple_object_pool<Column_type>)
404{
405 for (const Column_type* col : matrixToCopy.repToColumn_) {
406 if (col != nullptr) {
407 if constexpr (Master_matrix::Option_list::has_row_access) {
408 repToColumn_[nextColumnIndex_] =
409 columnPool_->construct(*col, col->get_column_index(), ra_opt::rows_, colSettings_);
410 } else {
411 repToColumn_[nextColumnIndex_] = columnPool_->construct(*col, colSettings_);
412 }
413 columnToRep_.insert(columnToRep_.end(), *repToColumn_[nextColumnIndex_]);
414 repToColumn_[nextColumnIndex_]->set_rep(nextColumnIndex_);
415 }
416 nextColumnIndex_++;
417 }
418}
419
420template <class Master_matrix>
423 : ra_opt(std::move(static_cast<ra_opt&>(other))),
424 columnToRep_(std::move(other.columnToRep_)),
425 columnClasses_(std::move(other.columnClasses_)),
426 repToColumn_(std::move(other.repToColumn_)),
427 nextColumnIndex_(std::exchange(other.nextColumnIndex_, 0)),
428 colSettings_(std::exchange(other.colSettings_, nullptr)),
429 columnPool_(std::exchange(other.columnPool_, nullptr))
430{}
431
432template <class Master_matrix>
434{
435 columnToRep_.clear_and_dispose(delete_disposer(this));
436}
437
438template <class Master_matrix>
439template <class Container_type>
441{
442 insert_boundary(column);
443}
444
445template <class Master_matrix>
446template <class Boundary_type>
448 dimension_type dim)
449{
450 // handles a dimension which is not actually stored.
451 // TODO: verify if this is not a problem for the cohomology algorithm and if yes,
452 // change how Column_dimension_option is choosen.
453 if (dim == -1) dim = boundary.size() == 0 ? 0 : boundary.size() - 1;
454
455 if constexpr (Master_matrix::Option_list::has_row_access && !Master_matrix::Option_list::has_removable_rows) {
456 if (boundary.begin() != boundary.end()) {
457 index pivot;
458 if constexpr (Master_matrix::Option_list::is_z2) {
459 pivot = *std::prev(boundary.end());
460 } else {
461 pivot = std::prev(boundary.end())->first;
462 }
463 if (ra_opt::rows_->size() <= pivot) ra_opt::rows_->resize(pivot + 1);
464 }
465 }
466
467 if (repToColumn_.size() == nextColumnIndex_) {
468 // could perhaps be avoided, if find_set returns something special when it does not find
469 columnClasses_.link(nextColumnIndex_, nextColumnIndex_);
470 if constexpr (Master_matrix::Option_list::has_row_access) {
471 repToColumn_.push_back(
472 columnPool_->construct(nextColumnIndex_, boundary, dim, ra_opt::rows_, colSettings_));
473 } else {
474 repToColumn_.push_back(columnPool_->construct(boundary, dim, colSettings_));
475 }
476 } else {
477 if constexpr (Master_matrix::Option_list::has_row_access) {
478 repToColumn_[nextColumnIndex_] =
479 columnPool_->construct(nextColumnIndex_, boundary, dim, ra_opt::rows_, colSettings_);
480 } else {
481 repToColumn_[nextColumnIndex_] = columnPool_->construct(boundary, dim, colSettings_);
482 }
483 }
484 _insert_column(nextColumnIndex_);
485
486 nextColumnIndex_++;
487}
488
489template <class Master_matrix>
492{
493 auto col = repToColumn_[columnClasses_.find_set(columnIndex)];
494 if (col == nullptr) return empty_column_;
495 return *col;
496}
497
498template <class Master_matrix>
501{
502 static_assert(Master_matrix::Option_list::has_row_access, "Row access has to be enabled for this method.");
503
504 return ra_opt::get_row(rowIndex);
505}
506
507template <class Master_matrix>
509{
510 if constexpr (Master_matrix::Option_list::has_row_access && Master_matrix::Option_list::has_removable_rows) {
511 ra_opt::erase_empty_row(rowIndex);
512 }
513}
514
515template <class Master_matrix>
521
522template <class Master_matrix>
523template <class Cell_range_or_column_index>
524inline void Base_matrix_with_column_compression<Master_matrix>::add_to(const Cell_range_or_column_index& sourceColumn,
525 index targetColumnIndex)
526{
527 // handle case where targetRep == sourceRep?
528 index targetRep = columnClasses_.find_set(targetColumnIndex);
529 Column_type& target = *repToColumn_[targetRep];
530 columnToRep_.erase(target);
531 if constexpr (std::is_integral_v<Cell_range_or_column_index>) {
532 target += get_column(sourceColumn);
533 } else {
534 target += sourceColumn;
535 }
536 _insert_column(targetRep);
537}
538
539template <class Master_matrix>
540template <class Cell_range_or_column_index>
542 const Cell_range_or_column_index& sourceColumn, const Field_element_type& coefficient, index targetColumnIndex)
543{
544 // handle case where targetRep == sourceRep?
545 index targetRep = columnClasses_.find_set(targetColumnIndex);
546 Column_type& target = *repToColumn_[targetRep];
547 columnToRep_.erase(target);
548 if constexpr (std::is_integral_v<Cell_range_or_column_index>) {
549 target.multiply_target_and_add(coefficient, get_column(sourceColumn));
550 } else {
551 target.multiply_target_and_add(coefficient, sourceColumn);
552 }
553 _insert_column(targetRep);
554}
555
556template <class Master_matrix>
557template <class Cell_range_or_column_index>
559 const Field_element_type& coefficient, const Cell_range_or_column_index& sourceColumn, index targetColumnIndex)
560{
561 // handle case where targetRep == sourceRep?
562 index targetRep = columnClasses_.find_set(targetColumnIndex);
563 Column_type& target = *repToColumn_[targetRep];
564 columnToRep_.erase(target);
565 if constexpr (std::is_integral_v<Cell_range_or_column_index>) {
566 target.multiply_source_and_add(get_column(sourceColumn), coefficient);
567 } else {
568 target.multiply_source_and_add(sourceColumn, coefficient);
569 }
570 _insert_column(targetRep);
571}
572
573template <class Master_matrix>
575{
576 auto col = repToColumn_[columnClasses_.find_set(columnIndex)];
577 if (col == nullptr) return true;
578 return !col->is_non_zero(rowIndex);
579}
580
581template <class Master_matrix>
583{
584 auto col = repToColumn_[columnClasses_.find_set(columnIndex)];
585 if (col == nullptr) return true;
586 return col->is_empty();
587}
588
589template <class Master_matrix>
592{
593 for (auto col : repToColumn_) {
594 if (col != nullptr) {
595 columnPool_->destroy(col);
596 col = nullptr;
597 }
598 }
599 ra_opt::operator=(other);
600 columnClasses_ = other.columnClasses_;
601 columnToRep_.reserve(other.columnToRep_.size());
602 repToColumn_.resize(other.repToColumn_.size(), nullptr);
603 nextColumnIndex_ = 0;
604 colSettings_ = other.colSettings_;
605 for (const Column_type* col : other.repToColumn_) {
606 if constexpr (Master_matrix::Option_list::has_row_access) {
607 repToColumn_[nextColumnIndex_] =
608 columnPool_->construct(*col, col->get_column_index(), ra_opt::rows_, colSettings_);
609 } else {
610 repToColumn_[nextColumnIndex_] = columnPool_->construct(*col, colSettings_);
611 }
612 columnToRep_.insert(columnToRep_.end(), *repToColumn_[nextColumnIndex_]);
613 repToColumn_[nextColumnIndex_]->set_rep(nextColumnIndex_);
614 nextColumnIndex_++;
615 }
616 return *this;
617}
618
619template <class Master_matrix>
621{
622 std::cout << "Compressed_matrix:\n";
623 for (Column_type& col : columnToRep_) {
624 for (auto e : col->get_content(nextColumnIndex_)) {
625 if (e == 0u)
626 std::cout << "- ";
627 else
628 std::cout << e << " ";
629 }
630 std::cout << "(";
631 for (index i = 0; i < nextColumnIndex_; ++i) {
632 if (columnClasses_.find_set(i) == col.get_rep()) std::cout << i << " ";
633 }
634 std::cout << ")\n";
635 }
636 std::cout << "\n";
637 std::cout << "Row Matrix:\n";
638 for (index i = 0; i < ra_opt::rows_->size(); ++i) {
639 const Row_type& row = ra_opt::rows_[i];
640 for (const auto& cell : row) {
641 std::cout << cell.get_column_index() << " ";
642 }
643 std::cout << "(" << i << ")\n";
644 }
645 std::cout << "\n";
646}
647
648template <class Master_matrix>
649inline void Base_matrix_with_column_compression<Master_matrix>::_insert_column(index columnIndex)
650{
651 Column_type& col = *repToColumn_[columnIndex];
652
653 if (col.is_empty()) {
654 columnPool_->destroy(&col); // delete curr_col;
655 repToColumn_[columnIndex] = nullptr;
656 return;
657 }
658
659 col.set_rep(columnIndex);
660 auto res = columnToRep_.insert(col);
661 if (res.first->get_rep() != columnIndex) { //if true, then redundant column
662 _insert_double_column(columnIndex, res.first);
663 }
664}
665
666template <class Master_matrix>
667inline void Base_matrix_with_column_compression<Master_matrix>::_insert_double_column(
668 index columnIndex, typename col_dict_type::iterator& doubleIt)
669{
670 index doubleRep = doubleIt->get_rep();
671 columnClasses_.link(columnIndex, doubleRep); // both should be representatives
672 index newRep = columnClasses_.find_set(columnIndex);
673
674 columnPool_->destroy(repToColumn_[columnIndex]);
675 repToColumn_[columnIndex] = nullptr;
676
677 if (newRep == columnIndex) {
678 std::swap(repToColumn_[doubleRep], repToColumn_[columnIndex]);
679 doubleIt->set_rep(columnIndex);
680 }
681}
682
683} // namespace persistence_matrix
684} // namespace Gudhi
685
686#endif // PM_BASE_MATRIX_COMPRESSION_H
Type for columns. Only one for each "column class" is explicitely constructed.
Definition base_matrix_with_column_compression.h:67
A base matrix (also see Base_matrix), but with column compression. That is, all identical columns in ...
Definition base_matrix_with_column_compression.h:46
bool is_zero_cell(index columnIndex, index rowIndex)
Indicates if the cell at given coordinates has value zero.
Definition base_matrix_with_column_compression.h:574
void reset(Column_settings *colSettings)
Resets the matrix to an empty matrix.
Definition base_matrix_with_column_compression.h:298
typename Master_matrix::Column_settings Column_settings
Definition base_matrix_with_column_compression.h:59
bool is_zero_column(index columnIndex)
Indicates if the column at given index has value zero.
Definition base_matrix_with_column_compression.h:582
const Column_type & get_column(index columnIndex)
Returns the column at the given MatIdx index. The type of the column depends on the choosen options,...
Definition base_matrix_with_column_compression.h:491
void erase_empty_row(index rowIndex)
If PersistenceMatrixOptions::has_row_access and PersistenceMatrixOptions::has_removable_rows are true...
Definition base_matrix_with_column_compression.h:508
typename Master_matrix::dimension_type dimension_type
Definition base_matrix_with_column_compression.h:49
void multiply_source_and_add_to(const Field_element_type &coefficient, const Cell_range_or_column_index &sourceColumn, index targetColumnIndex)
Multiplies the source column with the coefficiant before adding it to the target column....
Definition base_matrix_with_column_compression.h:558
void add_to(const Cell_range_or_column_index &sourceColumn, index targetColumnIndex)
Adds column represented by sourceColumn onto the column at targetColumnIndex in the matrix.
Definition base_matrix_with_column_compression.h:524
Base_matrix_with_column_compression(Column_settings *colSettings)
Constructs an empty matrix.
Definition base_matrix_with_column_compression.h:363
typename Master_matrix::Field_operators Field_operators
Field operators class. Necessary only if PersistenceMatrixOptions::is_z2 is false.
Definition base_matrix_with_column_compression.h:53
void insert_column(const Container_type &column)
Inserts a new ordered column at the end of the matrix by copying the given range of Matrix::cell_rep_...
Definition base_matrix_with_column_compression.h:440
void insert_boundary(const Boundary_type &boundary, dimension_type dim=-1)
Same as insert_column, only for interface purposes. The given dimension is ignored and not stored.
Definition base_matrix_with_column_compression.h:447
typename Master_matrix::Cell_constructor Cell_constructor
Definition base_matrix_with_column_compression.h:57
void multiply_target_and_add_to(const Cell_range_or_column_index &sourceColumn, const Field_element_type &coefficient, index targetColumnIndex)
Multiplies the target column with the coefficiant and then adds the source column to it....
Definition base_matrix_with_column_compression.h:541
typename Master_matrix::Row_type Row_type
Definition base_matrix_with_column_compression.h:56
friend void swap(Base_matrix_with_column_compression &matrix1, Base_matrix_with_column_compression &matrix2)
Swap operator.
Definition base_matrix_with_column_compression.h:313
const Row_type & get_row(index rowIndex) const
Only available if PersistenceMatrixOptions::has_row_access is true. Returns the row at the given row ...
Definition base_matrix_with_column_compression.h:500
typename Master_matrix::index index
Definition base_matrix_with_column_compression.h:48
~Base_matrix_with_column_compression()
Destructor.
Definition base_matrix_with_column_compression.h:433
typename Master_matrix::element_type Field_element_type
Definition base_matrix_with_column_compression.h:54
index get_number_of_columns() const
Returns the current number of columns in the matrix, counting also the redundant columns.
Definition base_matrix_with_column_compression.h:517
Base_matrix_with_column_compression & operator=(const Base_matrix_with_column_compression &other)
Assign operator.
Definition base_matrix_with_column_compression.h:591
Gudhi namespace.
Definition SimplicialComplexForAlpha.h:14