9#include <kernel/lafem/power_vector.hpp>
10#include <kernel/lafem/sparse_layout.hpp>
11#include <kernel/lafem/meta_element.hpp>
12#include <kernel/lafem/container.hpp>
41 static_assert(blocks_ > 1,
"invalid block size");
44 template<
typename,
int>
54 typedef typename SubMatrixType::DataType
DataType;
56 typedef typename SubMatrixType::IndexType
IndexType;
64 template <
typename DT2_ = DataType,
typename IT2_ = IndexType>
80 _first(std::move(the_first)),
81 _rest(std::move(the_rest))
134 std::getline(file, line);
136 }
while (line.find(
"%%") == 0 || line ==
"");
139 _first = std::move(tmp_first);
141 RestClass tmp_rest(mode, file, directory);
142 _rest = std::move(tmp_rest);
154 auto found = filename.rfind(
"/");
155 if (found != std::string::npos)
157 directory = filename.substr(0, found + 1);
160 std::ifstream file(filename.c_str(), std::ifstream::in);
161 if (! file.is_open())
162 XABORTM(
"Unable to open Matrix file " + filename);
165 std::getline(file, line);
166 if (line.find(
"%%MatrixMarket powerrowmatrix coordinate real general") == String::npos)
167 XABORTM(
"Input-file is not a complatible file");
171 _first = std::move(other._first);
172 _rest = std::move(other._rest);
182 _first = std::move(other._first);
183 _rest = std::move(other._rest);
206 std::ofstream file(filename.c_str(), std::ofstream::out);
207 if (! file.is_open())
208 XABORTM(
"Unable to open Matrix file " + filename);
211 auto found = filename.rfind(
".");
212 if (found != std::string::npos)
214 suffix = filename.substr(found);
215 filename.erase(found);
217 found = filename.rfind(
"/");
218 if (found != std::string::npos)
220 directory = filename.substr(0, found + 1);
221 filename.erase(0, found + 1);
224 file <<
"%%MatrixMarket powerrowmatrix coordinate real general\n";
225 for (
Index i(1); i <= blocks_; ++i)
227 file << filename <<
"_pr" << i << suffix <<
"\n";
245 _first.write_out(mode, directory + prefix +
"_pr" +
stringify(length + 1 - blocks_) + suffix);
275 template<
int i_,
int j_>
278 static_assert(i_ == 0,
"invalid sub-matrix index");
279 static_assert((0 <= j_) && (j_ < blocks_),
"invalid sub-matrix index");
284 template<
int i_,
int j_>
287 static_assert(i_ == 0,
"invalid sub-matrix index");
288 static_assert((0 <= j_) && (j_ < blocks_),
"invalid sub-matrix index");
303 XASSERTM((i == 0) && (0 <= j) && (j < blocks_),
"invalid sub-matrix index");
310 XASSERTM((i == 0) && (0 <= j) && (j < blocks_),
"invalid sub-matrix index");
335 int row_blocks()
const
340 int col_blocks()
const
352 template <Perspective perspective_ = Perspective::native>
355 return first().template rows<perspective_>();
364 template <Perspective perspective_ = Perspective::native>
367 return first().template columns<perspective_>() + rest().template columns<perspective_>();
376 template <Perspective perspective_ = Perspective::native>
379 return first().template used_elements<perspective_>() + rest().template used_elements<perspective_>();
385 return String(
"PowerRowMatrix<") + SubMatrixType::name() +
"," +
stringify(blocks_) +
">";
388 template <Perspective perspective_ = Perspective::native>
402 first().format(
value);
403 rest().format(
value);
429 first().apply(r, x.first());
430 rest().apply(r, x.rest(), r,
DataType(1));
435 XASSERTM(r.
size() == this->rows(),
"Vector size of r does not match!");
436 XASSERTM(x.
size() == this->columns(),
"Vector size of x does not match!");
441 first().apply(r, x_first);
442 rest().apply(r, x_rest, r,
DataType(1));
459 first().apply_transposed(r.first(), x);
460 rest().apply_transposed(r.rest(), x);
465 XASSERTM(r.
size() == this->columns(),
"Vector size of r does not match!");
466 XASSERTM(x.
size() == this->rows(),
"Vector size of x does not match!");
471 first().apply_transposed(r_first, x);
472 rest().apply_transposed(r_rest, x);
493 first().apply(r, x.first(), y, alpha);
494 rest().apply(r, x.rest(), r, alpha);
500 XASSERTM(r.
size() == this->rows(),
"Vector size of r does not match!");
501 XASSERTM(x.
size() == this->columns(),
"Vector size of x does not match!");
502 XASSERTM(y.
size() == this->rows(),
"Vector size of y does not match!");
507 first().apply(r, x_first, y, alpha);
508 rest().apply(r, x_rest, r, alpha);
529 first().apply_transposed(r.first(), x, y.first(), alpha);
530 rest().apply_transposed(r.rest(), x, y.rest(), alpha);
536 XASSERTM(r.
size() == this->columns(),
"Vector size of r does not match!");
537 XASSERTM(x.
size() == this->rows(),
"Vector size of x does not match!");
538 XASSERTM(y.
size() == this->columns(),
"Vector size of y does not match!");
546 first().apply(r_first, x, y_first, alpha);
547 rest().apply(r_rest, x, y_rest, alpha);
572 return (this->
name() == x.
name()) && (this->first().same_layout(x.first())) && (this->rest().same_layout(x.rest()));
578 return first().create_vector_l();
590 return this->first().get_length_of_line(row) + this->rest().get_length_of_line(row);
596 const Index col_start,
const Index stride = 1)
const
600 this->first().set_line(row, pval_set, pcol_set, col_start, stride);
601 this->rest().set_line(row, pval_set + stride * length_of_base, pcol_set + stride * length_of_base, col_start + this->first().
template columns<Perspective::pod>(), stride);
604 void set_line_reverse(
const Index row,
const DataType *
const pval_set,
const Index stride = 1)
608 this->first().set_line_reverse(row, pval_set, stride);
609 this->rest().set_line_reverse(row, pval_set + stride * length_of_base, stride);
614 return first().row_degree(row) + rest().row_degree(row);
617 template<
typename IT2_>
618 Index get_row_col_indices(
const Index row, IT2_*
const pcol_idx,
const IT2_ col_offset)
const
620 const Index first_cols = first().template columns<Perspective::pod>();
621 const Index first_degree = first().get_row_col_indices(row, pcol_idx, col_offset);
622 return first_degree +
rest().get_row_col_indices(row, pcol_idx + first_degree, col_offset + IT2_(first_cols));
625 template<
typename DT2_>
626 Index get_row_values(
const Index row, DT2_ *
const pvals)
const
628 const Index first_degree = first().get_row_values(row, pvals);
629 return first_degree +
rest().get_row_values(row, pvals + first_degree);
632 template<
typename DT2_>
633 Index set_row_values(
const Index row,
const DT2_ *
const pvals)
635 const Index first_degree = first().set_row_values(row, pvals);
636 return first_degree +
rest().set_row_values(row, pvals + first_degree);
650 std::uint64_t isize = *(std::uint64_t*) data.data();
651 std::vector<char>::iterator start = std::begin(data) +
sizeof(std::uint64_t);
652 std::vector<char>::iterator last_of_first = std::begin(data) +
sizeof(std::uint64_t) + (
int) isize;
653 std::vector<char> buffer_first(start, last_of_first);
654 _first.restore_from_checkpoint_data(buffer_first);
656 data.erase(std::begin(data), last_of_first);
663 std::size_t old_size = data.size();
664 data.insert(std::end(data),
sizeof(std::uint64_t), 0);
665 std::uint64_t ireal_size =
_first.set_checkpoint_data(data, config);
666 char* csize =
reinterpret_cast<char*
>(&ireal_size);
667 for(std::size_t i(0) ; i <
sizeof(std::uint64_t) ; ++i)
669 data[old_size +i] = csize[i];
682 template <
typename SubType2_>
685 this->first().convert(other.first());
686 this->rest().convert(other.rest());
689 template <
typename SubType2_>
692 this->first().convert_reverse(other.first());
693 this->rest().convert_reverse(other.rest());
698 template<
typename SubType_>
699 class PowerRowMatrix<SubType_, 1>
701 template<
typename,
int>
702 friend class PowerRowMatrix;
706 typedef typename SubMatrixType::DataType
DataType;
707 typedef typename SubMatrixType::IndexType
IndexType;
711 typedef typename SubMatrixType::VectorTypeL
VectorTypeL;
713 typedef PowerVector<typename SubMatrixType::VectorTypeR, 1>
VectorTypeR;
715 template <
typename DT2_ = DataType,
typename IT2_ = IndexType>
716 using ContainerType = PowerRowMatrix<typename SubType_::template ContainerType<DT2_, IT2_>, 1>;
726 _first(std::move(the_first))
737 explicit PowerRowMatrix(
const SparseLayout<IndexType, layout_id>& layout) :
743 PowerRowMatrix(PowerRowMatrix&& other) :
749 explicit PowerRowMatrix(
FileMode mode,
const String& filename)
755 explicit PowerRowMatrix(
FileMode mode, std::istream& file,
const String& directory =
"")
761 std::getline(file, line);
763 }
while (line.find(
"%%") == 0 || line ==
"");
766 _first = std::move(tmp_first);
772 auto found = filename.rfind(
"/");
773 if (found != std::string::npos)
775 directory = filename.substr(0, found + 1);
778 std::ifstream file(filename.c_str(), std::ifstream::in);
779 if (! file.is_open())
780 XABORTM(
"Unable to open Matrix file " + filename);
782 PowerRowMatrix
other(mode, file, directory);
790 PowerRowMatrix&
operator=(PowerRowMatrix&& other)
800 PowerRowMatrix(
const PowerRowMatrix&) =
delete;
802 PowerRowMatrix&
operator=(
const PowerRowMatrix&) =
delete;
811 std::ofstream file(filename.c_str(), std::ofstream::out);
812 if (! file.is_open())
813 XABORTM(
"Unable to open Matrix file " + filename);
815 String suffix, directory;
816 auto found = filename.rfind(
".");
817 if (found != std::string::npos)
819 suffix = filename.substr(found);
820 filename.erase(found);
822 found = filename.rfind(
"/");
823 if (found != std::string::npos)
825 directory = filename.substr(0, found + 1);
826 filename.erase(0, found + 1);
829 file <<
"%%MatrixMarket powerrowmatrix coordinate real general\n";
830 file << filename <<
"_pr" << 1 << suffix <<
"\n";
839 _first.write_out(mode, directory + prefix +
"_pr" +
stringify(length) + suffix);
844 return PowerRowMatrix(
_first.clone(mode));
848 std::size_t
bytes()
const
854 template<
int i,
int j>
857 static_assert(i == 0,
"invalid sub-matrix index");
858 static_assert(j == 0,
"invalid sub-matrix index");
862 template<
int i,
int j>
865 static_assert(i == 0,
"invalid sub-matrix index");
866 static_assert(j == 0,
"invalid sub-matrix index");
872 XASSERTM((i == 0) && (j == 0),
"invalid sub-matrix index");
878 XASSERTM((i == 0) && (j == 0),
"invalid sub-matrix index");
892 int row_blocks()
const
897 int col_blocks()
const
902 template <Perspective perspective_ = Perspective::native>
905 return first().template rows<perspective_>();
908 template <Perspective perspective_ = Perspective::native>
911 return first().template columns<perspective_>();
914 template <Perspective perspective_ = Perspective::native>
917 return first().template used_elements<perspective_>();
922 return String(
"PowerRowMatrix<") + SubMatrixType::name() +
"," +
stringify(1) +
">";
925 template <Perspective perspective_ = Perspective::native>
933 first().format(
value);
943 first().apply(r, x.first());
946 void apply(DenseVector<DataType, IndexType>& r,
const DenseVector<DataType, IndexType>& x)
const
953 first().apply(r, x.first(), y, alpha);
956 void apply(DenseVector<DataType, IndexType>& r,
const DenseVector<DataType, IndexType>& x,
957 const DenseVector<DataType, IndexType>& y,
DataType alpha =
DataType(1))
const
959 first().apply(r, x, y, alpha);
964 first().apply_transposed(r.first(), x);
967 void apply_transposed(DenseVector<DataType, IndexType>& r,
const DenseVector<DataType, IndexType>& x)
const
969 first().apply_transposed(r, x);
974 first().apply_transposed(r.first(), x, y.first(), alpha);
977 void apply_transposed(DenseVector<DataType, IndexType>& r,
const DenseVector<DataType, IndexType>& x,
978 const DenseVector<DataType, IndexType>& y,
DataType alpha =
DataType(1))
const
980 first().apply_transposed(r, x, y, alpha);
991 return this->first().max_rel_diff(x.first());
1004 return (this->
name() == x.name()) && (this->first().same_layout(x.first()));
1022 return this->first().get_length_of_line(row);
1027 const Index col_start,
const Index stride = 1)
const
1029 this->first().set_line(row, pval_set, pcol_set, col_start, stride);
1032 void set_line_reverse(
const Index row,
const DataType *
const pval_set,
const Index stride = 1)
1034 this->first().set_line_reverse(row, pval_set, stride);
1039 return first().row_degree(row);
1042 template<
typename IT2_>
1043 Index get_row_col_indices(
const Index row, IT2_*
const pcol_idx,
const IT2_ col_offset)
const
1045 return first().get_row_col_indices(row, pcol_idx, col_offset);
1048 template<
typename DT2_>
1049 Index get_row_values(
const Index row, DT2_ *
const pvals)
const
1051 return first().get_row_values(row, pvals);
1054 template<
typename DT2_>
1055 Index set_row_values(
const Index row,
const DT2_ *
const pvals)
1057 return first().set_row_values(row, pvals);
1063 return _first.get_checkpoint_size(config);
1069 _first.restore_from_checkpoint_data(data);
1075 return _first.set_checkpoint_data(data, config);
1085 template <
typename SubType2_>
1086 void convert(
const PowerRowMatrix<SubType2_, 1> & other)
1088 this->first().convert(
other.first());
1091 template <
typename SubType2_>
1092 void convert_reverse(PowerRowMatrix<SubType2_, 1> & other)
const
1094 this->first().convert_reverse(
other.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.
Power-Row-Matrix meta class template.
Index rows() const
Returns the total number of rows in this matrix.
void convert(const PowerRowMatrix< SubType2_, blocks_ > &other)
Conversion method.
const SubMatrixType & get(int i, int j) const
Returns a sub-matrix block.
SubMatrixType & get(int i, int j)
Returns a sub-matrix block.
PowerRowMatrix< SubType_, blocks_-1 > RestClass
rest-class typedef
PowerRowMatrix(const PowerRowMatrix &)=delete
deleted copy-ctor
void write_out_submatrices(FileMode mode, String directory, String prefix, String suffix, Index length=blocks_) const
Write out submatrices to file.
void apply(VectorTypeL &r, const VectorTypeR &x) const
Applies this matrix onto a vector.
SubMatrixType & at()
Returns a sub-matrix block.
bool same_layout(const PowerRowMatrix &x) const
Checks if the structural layout of this matrix matches that of another matrix. This excludes comparis...
Index used_elements() const
Returns the total number of non-zeros in this matrix.
PowerRowMatrix(FileMode mode, String filename)
Constructor.
void read_from(FileMode mode, String filename)
Read in matrix from file.
void write_out(FileMode mode, String filename) const
Write out matrix to file.
PowerRowMatrix & operator=(const PowerRowMatrix &)=delete
deleted copy-assign operator
VectorTypeR create_vector_r() const
Returns a new compatible R-Vector.
PowerRowMatrix()
default ctor
SubMatrixType::IndexType IndexType
sub-matrix index type
PowerVector< typename SubMatrixType::VectorTypeR, blocks_ > VectorTypeR
Compatible R-vector type.
Index columns() const
Returns the total number of columns in this matrix.
void apply(VectorTypeL &r, const VectorTypeR &x, const VectorTypeL &y, DataType alpha=DataType(1)) const
Applies this matrix onto a vector.
VectorTypeL create_vector_l() const
Returns a new compatible L-Vector.
std::size_t bytes() const
Returns the total amount of bytes allocated.
PowerRowMatrix(PowerRowMatrix &&other)
move ctor
void restore_from_checkpoint_data(std::vector< char > &data)
Extract object from checkpoint.
PowerRowMatrix(SubMatrixType &&the_first, RestClass &&the_rest)
base-class constructor; this one is protected for a reason
SubMatrixType::VectorTypeL VectorTypeL
Compatible L-vector type.
void format(DataType value=DataType(0))
Reset all elements of the container to a given value or zero if missing.
static constexpr int num_col_blocks
number of column blocks (horizontal size)
static constexpr int num_row_blocks
number of row blocks (vertical size)
void apply_transposed(VectorTypeR &r, const VectorTypeL &x, const VectorTypeR &y, DataType alpha=DataType(1)) const
Applies this matrix onto a vector.
std::uint64_t get_checkpoint_size(SerialConfig &config)
Calculate size.
SubType_ SubMatrixType
sub-matrix type
std::uint64_t set_checkpoint_data(std::vector< char > &data, SerialConfig &config)
DataType max_rel_diff(const PowerRowMatrix &x) const
Retrieve the maximum relative difference of this matrix and another one y.max_rel_diff(x) returns .
PowerRowMatrix(const SparseLayout< IndexType, layout_id > &layout)
sub-matrix layout ctor
static String name()
Returns a descriptive string for this container.
SubMatrixType::DataType DataType
sub-matrix data type
static constexpr SparseLayoutId layout_id
sub-matrix layout type
void clear()
Free all allocated arrays.
SubMatrixType _first
the first sub-matrix
virtual ~PowerRowMatrix()
virtual destructor
Index get_length_of_line(const Index row) const
Returns the number of NNZ-elements of the selected row.
PowerRowMatrix< typename SubType_::template ContainerType< DT2_, IT2_ >, blocks_ > ContainerType
Our 'base' class type.
PowerRowMatrix & operator=(PowerRowMatrix &&other)
move-assign operator
RestClass _rest
the remaining part
PowerRowMatrix(FileMode mode, std::istream &file, String directory="")
Constructor.
PowerRowMatrix clone(LAFEM::CloneMode mode=LAFEM::CloneMode::Weak) const
Creates and returns a deep copy of this matrix.
const SubMatrixType & at() const
Returns a sub-matrix block.
void apply_transposed(VectorTypeR &r, const VectorTypeL &x) const
Applies this matrix onto a vector.
Power-Vector meta class template.
Config class for serialize parameter.
Layout scheme for sparse matrix containers.
String class implementation.
String & trim_me(const String &charset)
Trims this string.
@ other
generic/other permutation strategy
T_ max(T_ a, T_ b)
Returns the maximum of two values.
@ 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.
Power container element helper class template.