1 #ifndef _theplu_yat_utility_matrix_ 2 #define _theplu_yat_utility_matrix_ 34 #include "BLAS_level2.h" 35 #include "BLAS_level3.h" 37 #include "Exception.h" 38 #include "MatrixExpression.h" 39 #include "StrideIterator.h" 41 #include "VectorConstView.h" 42 #include "VectorExpression.h" 43 #include "VectorView.h" 44 #include "yat_assert.h" 46 #include <gsl/gsl_matrix.h> 144 Matrix(
const size_t& r,
const size_t& c,
double init_value=0);
175 #ifdef YAT_HAVE_RVALUE 200 : blas_result_(NULL), m_(NULL)
237 explicit Matrix(std::istream &,
char sep=
'\0')
248 void all(const
double value);
256 iterator
begin(
void);
264 const_iterator
begin(
void) const;
292 const_iterator
begin_row(
size_t i) const;
326 const_iterator
end(
void) const;
346 const_iterator
end_row(
size_t i) const;
357 bool equal(const
Matrix&, const
double precision=0) const;
376 void mul(const Matrix& b);
388 void resize(
size_t r,
size_t c,
double init_value=0);
393 size_t rows(
void) const;
427 void swap_rows(const
size_t i, const
size_t j);
449 double& operator()(
size_t row,
size_t column);
461 const
double& operator()(
size_t row,
size_t column) const;
474 bool operator==(const Matrix& other) const;
487 bool operator!=(const Matrix& other) const;
497 const Matrix& operator=(const Matrix& other);
520 #ifdef YAT_HAVE_RVALUE 659 gsl_matrix* create_gsl_matrix_copy(
void)
const;
666 void delete_allocated_memory(
void);
672 gsl_matrix* blas_result_;
683 void inverse_svd(
const Matrix& input, Matrix& result);
693 bool isnull(
const Matrix&);
702 double max(
const Matrix&);
711 double min(
const Matrix&);
724 std::pair<size_t,size_t>& min,
725 std::pair<size_t,size_t>& max);
742 bool nan(
const Matrix& templat, Matrix& flag);
757 void swap(Matrix&, Matrix&);
764 std::ostream&
operator<< (std::ostream& s,
const Matrix&);
775 double trace(
const Matrix&);
780 : blas_result_(NULL), m_(NULL)
792 detail::deallocate(m_);
802 rhs.
get(blas_result_);
803 std::swap(m_, blas_result_);
const VectorConstView column_const_view(size_t) const
double trace(const Matrix &)
Trace of matrix.
void mul(const Matrix &b)
const double & const_reference
Definition: Matrix.h:96
void inverse_svd(const Matrix &input, Matrix &result)
iterator end_row(size_t i)
Class for errors reported from underlying GSL calls.
Definition: Exception.h:87
The Department of Theoretical Physics namespace as we define it.
StrideIterator< const double * > const_iterator
Definition: Matrix.h:106
This is the yat interface to gsl_vector_view.
Definition: VectorView.h:79
Definition: stl_utility.h:64
void resize(size_t r, size_t c, double init_value=0)
Resize Matrix.
iterator begin_row(size_t i)
const gsl_matrix * gsl_matrix_p(void) const
Definition: MatrixExpression.h:140
StrideIterator< const double * > const_row_iterator
Definition: Matrix.h:126
bool nan(const Matrix &templat, Matrix &flag)
Create a Matrix flag indicating NaN's in another Matrix templat.
StrideIterator< double * > row_iterator
Definition: Matrix.h:121
double & reference
Definition: Matrix.h:89
void minmax_index(const Matrix &, std::pair< size_t, size_t > &min, std::pair< size_t, size_t > &max)
Locate the maximum and minimum element in the Matrix.
Matrix(void)
The default constructor.
void all(const double value)
Matrix & operator-=(const MatrixExpression< T > &rhs)
Subtraction and assign operator.
Definition: Matrix.h:611
Matrix & operator+=(const MatrixExpression< T > &rhs)
Addition and assign operator.
Definition: Matrix.h:578
Read-only view.
Definition: VectorConstView.h:56
iterator begin_column(size_t i)
double min(const Matrix &)
Get the minimum value of the Matrix.
double value_type
Definition: Matrix.h:82
bool equal(const Matrix &, const double precision=0) const
Check whether matrices are equal within a user defined precision, set by precision.
void swap(Matrix &, Matrix &)
Exchange all elements between the matrices by copying.
StrideIterator< const double * > const_column_iterator
Definition: Matrix.h:116
const Matrix & operator-=(const Matrix &)
Subtract and assign operator.
std::ostream & operator<<(std::ostream &s, const Matrix &)
The output operator for the Matrix class.
VectorView column_view(size_t i)
An expression that can be converted to a Matrix.
Definition: MatrixExpression.h:46
void get(gsl_matrix *&m) const
Definition: MatrixExpression.h:112
void div(const Matrix &b)
Class to report errors associated with IO operations.
Definition: Exception.h:109
VectorView row_view(size_t)
const VectorConstView row_const_view(size_t) const
Interface to GSL matrix.
Definition: Matrix.h:74
bool isnull(const Matrix &)
Check if all elements of the Matrix are zero.
Definition: BasicMatrix.h:38
const Matrix & operator=(const Matrix &other)
The assignment operator.
const gsl_matrix * gsl_matrix_p(void) const
StrideIterator< double * > iterator
Definition: Matrix.h:101
StrideIterator< double * > column_iterator
Definition: Matrix.h:111
const Matrix & operator+=(const Matrix &b)
Add and assign operator.
void swap_rows(const size_t i, const size_t j)
Swap rows i and j.
const Matrix & operator*=(const Matrix &)
Multiply and assignment operator.
iterator end_column(size_t i)
void swap_columns(const size_t i, const size_t j)
Swap columns i and j.
void transpose(void)
Transpose the matrix.
Adaptor using a stride on underlying iterator.
Definition: StrideIterator.h:50
double max(const Matrix &)
Get the maximum value of the Matrix.
void swap_rowcol(const size_t i, const size_t j)
Swap row i and column j.
size_t columns(void) const