9#include <kernel/lafem/tuple_vector.hpp>
10#include <kernel/lafem/sparse_layout.hpp>
11#include <kernel/lafem/meta_element.hpp>
37 template<
typename,
typename...>
55 template <
typename DT2_ = DataType,
typename IT2_ = IndexType>
61 template <
typename DT2_,
typename IT2_>
79 _first(std::forward<First_>(the_first)),
98 _first(std::forward<First_>(the_first)),
99 _rest(std::forward<Rest_>(the_rest)...)
139 std::getline(file, line);
141 }
while (line.find(
"%%") == 0 || line ==
"");
143 _first = First_(mode, directory + line);
156 auto found = filename.rfind(
"/");
157 if (found != std::string::npos)
159 directory = filename.substr(0, found + 1);
162 std::ifstream file(filename.c_str(), std::ifstream::in);
163 if (! file.is_open())
164 XABORTM(
"Unable to open Matrix file " + filename);
167 std::getline(file, line);
168 if (line.find(
"%%MatrixMarket tuplediagmatrix coordinate real general") == String::npos)
169 XABORTM(
"Input-file is not a complatible file");
173 _first = std::move(other._first);
174 _rest = std::move(other._rest);
184 _first = std::forward<First_>(other._first);
185 _rest = std::forward<RestClass>(other._rest);
208 std::ofstream file(filename.c_str(), std::ofstream::out);
209 if (! file.is_open())
210 XABORTM(
"Unable to open Matrix file " + filename);
213 auto found = filename.rfind(
".");
214 if (found != std::string::npos)
216 suffix = filename.substr(found);
217 filename.erase(found);
219 found = filename.rfind(
"/");
220 if (found != std::string::npos)
222 directory = filename.substr(0, found + 1);
223 filename.erase(0, found + 1);
226 file <<
"%%MatrixMarket tuplediagmatrix coordinate real general\n";
229 file << filename <<
"_td" << i << suffix <<
"\n";
285 template<
int i_,
int j_>
288 static_assert(i_ == j_,
"invalid sub-matrix index");
293 template<
int i_,
int j_>
296 static_assert(i_ == j_,
"invalid sub-matrix index");
306 const First_& first()
const
321 int row_blocks()
const
326 int col_blocks()
const
335 return first().rows() + rest().rows();
341 return first().columns() + rest().columns();
347 return first().used_elements() + rest().used_elements();
369 first().format(
value);
370 rest().format(
value);
385 first().extract_diag(diag.first());
386 rest().extract_diag(diag.rest());
403 first().apply(r.first(), x.first());
404 rest().apply(r.rest(), x.rest());
409 XASSERTM(r.
size() == this->rows(),
"Vector size of r does not match!");
410 XASSERTM(x.
size() == this->columns(),
"Vector size of x does not match!");
418 first().apply(r_first, x_first);
419 rest().apply(r_rest, x_rest);
436 first().apply_transposed(r.first(), x.first());
437 rest().apply_transposed(r.rest(), x.rest());
442 XASSERTM(r.
size() == this->columns(),
"Vector size of r does not match!");
443 XASSERTM(x.
size() == this->rows(),
"Vector size of x does not match!");
451 first().apply_transposed(r_first, x_first);
452 rest().apply_transposed(r_rest, x_rest);
473 first().apply(r.first(), x.first(), y.first(), alpha);
474 rest().apply(r.rest(), x.rest(), y.rest(), alpha);
480 XASSERTM(r.
size() == this->rows(),
"Vector size of r does not match!");
481 XASSERTM(x.
size() == this->columns(),
"Vector size of x does not match!");
482 XASSERTM(y.
size() == this->rows(),
"Vector size of y does not match!");
493 first().apply(r_first, x_first, y_first, alpha);
494 rest().apply(r_rest, x_rest, y_rest, alpha);
515 first().apply_transposed(r.first(), x.first(), y.first(), alpha);
516 rest().apply_transposed(r.rest(), x.rest(), y.rest(), alpha);
522 XASSERTM(r.
size() == this->columns(),
"Vector size of r does not match!");
523 XASSERTM(x.
size() == this->rows(),
"Vector size of x does not match!");
524 XASSERTM(y.
size() == this->columns(),
"Vector size of y does not match!");
535 first().apply_transposed(r_first, x_first, y_first, alpha);
536 rest().apply_transposed(r_rest, x_rest, y_rest, alpha);
553 first().scale_rows(a.first(), w.first());
554 rest().scale_rows(a.rest(), w.rest());
557 void scale_cols(
const TupleDiagMatrix& a,
const VectorTypeR& w)
559 first().scale_cols(a.first(), w.first());
560 rest().scale_cols(a.rest(), w.rest());
570 return this->first().get_length_of_line(row);
574 return this->rest().get_length_of_line(row - brows);
581 const Index col_start,
const Index stride = 1)
const
588 this->first().set_line(row, pval_set, pcol_set, col_start, stride);
592 this->rest().set_line(row - brows, pval_set, pcol_set, col_start + bcolumns, stride);
596 void set_line_reverse(
const Index row,
DataType *
const pval_set,
const Index stride = 1)
602 this->first().set_line_reverse(row, pval_set, stride);
606 this->
rest().set_line_reverse(row - brows, pval_set, stride);
612 const Index first_rows = first().template rows<Perspective::pod>();
614 return first().row_degree(row);
616 return rest().row_degree(row - first_rows);
619 template<
typename IT2_>
620 Index get_row_col_indices(
const Index row, IT2_*
const pcol_idx,
const IT2_ col_offset)
const
622 const Index first_rows = first().template rows<Perspective::pod>();
623 const Index first_cols = first().template columns<Perspective::pod>();
625 return first().get_row_col_indices(row, pcol_idx, col_offset);
627 return rest().get_row_col_indices(row - first_rows, pcol_idx, col_offset + IT2_(first_cols));
630 template<
typename DT2_>
631 Index get_row_values(
const Index row, DT2_ *
const pvals)
const
633 const Index first_rows = first().template rows<Perspective::pod>();
635 return first().get_row_values(row, pvals);
637 return rest().get_row_values(row - first_rows, pvals);
640 template<
typename DT2_>
641 Index set_row_values(
const Index row,
const DT2_ *
const pvals)
643 const Index first_rows = first().template rows<Perspective::pod>();
645 return first().set_row_values(row, pvals);
647 return rest().set_row_values(row - first_rows, pvals);
661 std::uint64_t isize = *(std::uint64_t*) data.data();
662 std::vector<char>::iterator start = std::begin(data) +
sizeof(std::uint64_t);
663 std::vector<char>::iterator last_of_first = std::begin(data) +
sizeof(std::uint64_t) + (
int) isize;
664 std::vector<char> buffer_first(start, last_of_first);
665 _first.restore_from_checkpoint_data(buffer_first);
667 data.erase(std::begin(data), last_of_first);
674 std::size_t old_size = data.size();
675 data.insert(std::end(data),
sizeof(std::uint64_t), 0);
676 std::uint64_t ireal_size =
_first.set_checkpoint_data(data, config);
677 char* csize =
reinterpret_cast<char*
>(&ireal_size);
678 for(std::size_t i(0) ; i <
sizeof(std::uint64_t) ; ++i)
680 data[old_size +i] = csize[i];
693 template <
typename First2_,
typename... Rest2_>
696 this->first().convert(other.first());
697 this->rest().convert(other.rest());
700 template <
typename First2_,
typename... Rest2_>
703 this->first().convert_reverse(other.first());
704 this->rest().convert_reverse(other.rest());
717 return (this->
name() == x.
name()) && (this->first().same_layout(x.first())) && (this->rest().same_layout(x.rest()));
722 template<
typename First_>
723 class TupleDiagMatrix<First_>
725 template<
typename,
typename...>
726 friend class TupleDiagMatrix;
729 typedef typename First_::DataType
DataType;
730 typedef typename First_::IndexType
IndexType;
732 typedef TupleVector<typename First_::VectorTypeL>
VectorTypeL;
734 typedef TupleVector<typename First_::VectorTypeR>
VectorTypeR;
736 template <
typename DT2_ = DataType,
typename IT2_ = IndexType>
737 using ContainerType = TupleDiagMatrix<typename First_::template ContainerType<DT2_, IT2_> >;
747 return First_::name();
756 explicit TupleDiagMatrix(First_&& the_first) :
757 _first(std::forward<First_>(the_first))
762 TupleDiagMatrix(TupleDiagMatrix&& other) :
768 explicit TupleDiagMatrix(
FileMode mode,
const String& filename)
774 explicit TupleDiagMatrix(
FileMode mode, std::istream& file,
const String& directory =
"")
780 std::getline(file, line);
782 }
while (line.find(
"%%") == 0 || line ==
"");
784 _first = First_(mode, directory + line);
790 auto found = filename.rfind(
"/");
791 if (found != std::string::npos)
793 directory = filename.substr(0, found + 1);
796 std::ifstream file(filename.c_str(), std::ifstream::in);
797 if (! file.is_open())
798 XABORTM(
"Unable to open Matrix file " + filename);
800 TupleDiagMatrix
other(mode, file, directory);
808 TupleDiagMatrix&
operator=(TupleDiagMatrix&& other)
818 TupleDiagMatrix(
const TupleDiagMatrix&) =
delete;
820 TupleDiagMatrix&
operator=(
const TupleDiagMatrix&) =
delete;
829 std::ofstream file(filename.c_str(), std::ofstream::out);
830 if (! file.is_open())
831 XABORTM(
"Unable to open Matrix file " + filename);
833 String suffix, directory;
834 auto found = filename.rfind(
".");
835 if (found != std::string::npos)
837 suffix = filename.substr(found);
838 filename.erase(found);
840 found = filename.rfind(
"/");
841 if (found != std::string::npos)
843 directory = filename.substr(0, found + 1);
844 filename.erase(0, found + 1);
847 file <<
"%%MatrixMarket tuplediagmatrix coordinate real general\n";
848 file << filename <<
"_td" << 1 << suffix <<
"\n";
857 _first.write_out(mode, directory + prefix +
"_td" +
stringify(length) + suffix);
862 return TupleDiagMatrix(
_first.clone(mode));
866 std::size_t
bytes()
const
873 return TupleDiagMatrix(
_first.transpose());
876 template<
int i,
int j>
877 typename TupleElement<i, First_>::Type&
at()
879 static_assert(i == 0,
"invalid sub-matrix index");
880 static_assert(j == 0,
"invalid sub-matrix index");
884 template<
int i,
int j>
885 const typename TupleElement<i, First_>::Type&
at()
const
887 static_assert(i == 0,
"invalid sub-matrix index");
888 static_assert(j == 0,
"invalid sub-matrix index");
897 const First_& first()
const
902 int row_blocks()
const
907 int col_blocks()
const
914 return first().rows();
919 return first().columns();
924 return first().used_elements();
929 first().format(
value);
950 first().extract_diag(diag.first());
955 first().apply(r.first(), x.first());
960 first().apply(r.first(), x.first(), y.first(), alpha);
965 first().apply_transposed(r.first(), x.first());
970 first().apply_transposed(r.first(), x.first(), y.first(), alpha);
985 void scale_rows(
const TupleDiagMatrix& a,
const VectorTypeL& w)
987 first().scale_rows(a.first(), w.first());
990 void scale_cols(
const TupleDiagMatrix& a,
const VectorTypeR& w)
992 first().scale_cols(a.first(), w.first());
998 return this->first().get_length_of_line(row);
1003 const Index col_start,
const Index stride = 1)
const
1005 this->first().set_line(row, pval_set, pcol_set, col_start, stride);
1008 void set_line_reverse(
const Index row,
DataType *
const pval_set,
const Index stride = 1)
const
1010 this->first().set_line_reverse(row, pval_set, stride);
1015 return first().row_degree(row);
1018 template<
typename IT2_>
1019 Index get_row_col_indices(
const Index row, IT2_*
const pcol_idx,
const IT2_ col_offset)
const
1021 return first().get_row_col_indices(row, pcol_idx, col_offset);
1024 template<
typename DT2_>
1025 Index get_row_values(
const Index row, DT2_ *
const pvals)
const
1027 return first().get_row_values(row, pvals);
1030 template<
typename DT2_>
1031 Index set_row_values(
const Index row,
const DT2_ *
const pvals)
1033 return first().set_row_values(row, pvals);
1039 return _first.get_checkpoint_size(config);
1045 _first.restore_from_checkpoint_data(data);
1051 return _first.set_checkpoint_data(data, config);
1061 template <
typename First2_>
1062 void convert(
const TupleDiagMatrix<First2_>& other)
1064 this->first().convert(
other.first());
1067 template <
typename First2_>
1068 void convert_reverse(TupleDiagMatrix<First2_>& other)
const
1070 this->first().convert_reverse(
other.first());
1083 return (this->
name() == x.name()) && (this->first().same_layout(x.first()));
#define XABORTM(msg)
Abortion macro definition with custom message.
#define XASSERTM(expr, msg)
Assertion macro definition with custom message.
Index size() const
Returns the containers size.
Dense data vector class template.
Config class for serialize parameter.
Tuple-Diag-Matrix meta class template.
void write_out(FileMode mode, String filename) const
Write out matrix to file.
TupleDiagMatrix< Rest_... > RestClass
rest-class typedef
TupleDiagMatrix(FileMode mode, std::istream &file, const String &directory="")
Constructor.
TupleDiagMatrix & operator=(TupleDiagMatrix &&other)
move-assign operator
TupleDiagMatrix()
default ctor
void convert(const TupleDiagMatrix< First2_, Rest2_... > &other)
Conversion method.
TupleDiagMatrix(const TupleDiagMatrix &)=delete
deleted copy-ctor
TupleDiagMatrix clone(LAFEM::CloneMode mode=LAFEM::CloneMode::Weak) const
Creates and returns a deep copy of this matrix.
TupleDiagMatrix(First_ &&the_first, RestClass &&the_rest)
base-class constructor; this one is protected for a reason
void restore_from_checkpoint_data(std::vector< char > &data)
Extract object from checkpoint.
Index columns() const
Returns the total number of columns in this matrix.
void format(DataType value=DataType(0))
Reset all elements of the container to a given value or zero if missing.
Index get_length_of_line(const Index row) const
Returns the number of NNZ-elements of the selected row.
TupleDiagMatrix< typename First_::template ContainerType< DT2_, IT2_ >, typename Rest_::template ContainerType< DT2_, IT2_ >... > ContainerType
Our 'base' class type.
TupleDiagMatrix & operator=(const TupleDiagMatrix &)=delete
deleted copy-assign operator
First_ _first
the first sub-matrix
std::uint64_t get_checkpoint_size(SerialConfig &config)
Calculate size.
VectorTypeR create_vector_r() const
Returns a new compatible R-Vector.
const TupleElement< i_, First_, Rest_... >::Type & at() const
Returns a sub-matrix block.
void apply(VectorTypeL &r, const VectorTypeR &x, const VectorTypeL &y, DataType alpha=DataType(1)) const
Applies this matrix onto a vector.
void apply(VectorTypeL &r, const VectorTypeR &x) const
Applies this matrix onto a vector.
void read_from(FileMode mode, const String &filename)
Read in matrix from file.
std::uint64_t set_checkpoint_data(std::vector< char > &data, SerialConfig &config)
TupleDiagMatrix transpose() const
Creates and returns the transpose of this matrix.
First_::DataType DataType
sub-matrix data type
bool same_layout(const TupleDiagMatrix &x) const
Checks if the structural layout of this matrix matches that of another matrix. This excludes comparis...
void apply_transposed(VectorTypeR &r, const VectorTypeL &x, const VectorTypeR &y, DataType alpha=DataType(1)) const
Applies this matrix onto a vector.
RestClass _rest
the remaining part
Index rows() const
Returns the total number of rows in this matrix.
TupleElement< i_, First_, Rest_... >::Type & at()
Returns a sub-matrix block.
void clear()
Free all allocated arrays.
TupleDiagMatrix(FileMode mode, const String &filename)
Constructor.
static String name()
Returns a descriptive string for this container.
void extract_diag(VectorTypeL &diag) const
extract main diagonal vector from matrix
virtual ~TupleDiagMatrix()
virtual destructor
std::size_t bytes() const
Returns the total amount of bytes allocated.
TupleVector< typename First_::VectorTypeL, typename Rest_::VectorTypeL... > VectorTypeL
Compatible L-vector type.
TupleVector< typename First_::VectorTypeR, typename Rest_::VectorTypeR... > VectorTypeR
Compatible R-vector type.
void write_out_submatrices(FileMode mode, const String &directory, const String &prefix, const String &suffix, Index length=num_row_blocks) const
Write out submatrices to file.
static String sub_name_list()
Returns a list of all sub-matrix type names.
void apply_transposed(VectorTypeR &r, const VectorTypeL &x) const
Applies this matrix onto a vector.
static constexpr int num_row_blocks
number of row blocks (vertical size)
TupleDiagMatrix(TupleDiagMatrix &&other)
move ctor
First_::IndexType IndexType
sub-matrix index type
static constexpr int num_col_blocks
number of column blocks (horizontal size)
Index used_elements() const
Returns the total number of non-zeros in this matrix.
TupleDiagMatrix(First_ &&the_first, Rest_ &&... the_rest)
Sub-Matrix emplacement constructor.
VectorTypeL create_vector_l() const
Returns a new compatible L-Vector.
Variadic TupleVector class template.
String class implementation.
String & trim_me(const String &charset)
Trims this string.
@ other
generic/other permutation strategy
@ rest
restriction (multigrid)
String stringify(const T_ &item)
Converts an item into a String.
@ value
specifies whether the space should supply basis function values
std::uint64_t Index
Index data type.
Tuple container element helper class template.