50 #ifndef MERCURY_SMALLMATRIX_H
51 #define MERCURY_SMALLMATRIX_H
66 template<
unsigned int numberOfRows,
unsigned int numberOfColumns>
80 logger.assert_debug(numberOfColumns == 1,
"Trying to construct a matrix with more than 1 columns from a vector");
81 std::copy(other.
data(), other.
data() + numberOfRows,
data_.begin());
94 logger.assert_debug(entries.size() == numberOfColumns,
"expected a matrix with % "
95 "columns, but got a matrix with % columns", numberOfColumns,
97 unsigned int column = 0;
100 for (
unsigned int i = 0;
i < numberOfRows; ++
i)
102 (*this)(
i, column) = entry[i];
119 for (
unsigned int i = 0;
i < numberOfRows; ++
i)
121 for (
unsigned int j = 0; j < numberOfColumns; ++j)
123 (*this)(
i, j) = entries[j][i];
137 logger.assert_debug(n < numberOfRows,
"Requested row number % for a matrix with only % rows", n, numberOfRows);
138 logger.assert_debug(m < numberOfColumns,
"Requested column number % for a matrix with only % columns", m,
140 return data_[n + m * numberOfRows];
146 logger.assert_debug(n < numberOfRows,
"Requested row number % for a matrix with only % rows", n, numberOfRows);
147 logger.assert_debug(m < numberOfColumns,
"Requested column number % for a matrix with only % columns", m,
149 return data_[n + m * numberOfRows];
155 logger.assert_debug(n < numberOfRows * numberOfColumns,
"Requested entry % for a matrix with only % entries", n,
156 numberOfRows * numberOfColumns);
162 logger.assert_debug(n < numberOfRows * numberOfColumns,
"Requested entry % for a matrix with only % entries", n,
163 numberOfRows * numberOfColumns);
177 std::bind(std::multiplies<Mdouble>(), std::placeholders::_1, right));
182 template<
unsigned int K>
185 template<
unsigned int K>
190 std::transform(
data_.begin(),
data_.end(), other.
data_.begin(),
data_.begin(), std::plus<Mdouble>());
196 std::transform(
data_.begin(),
data_.end(), other.
data_.begin(),
data_.begin(), std::minus<Mdouble>());
203 std::transform(
data_.begin(),
data_.end(), other.
data_.begin(), result.
data_.begin(), std::plus<Mdouble>());
210 std::transform(
data_.begin(),
data_.end(), other.
data_.begin(), result.
data_.begin(), std::minus<Mdouble>());
223 std::bind(std::multiplies<Mdouble>(), std::placeholders::_1, scalar));
235 std::bind(std::divides<Mdouble>(), std::placeholders::_1, scalar));
244 std::bind(std::divides<Mdouble>(), std::placeholders::_1, scalar));
258 std::move(right.data_.begin(), right.data_.end(),
data_.begin());
268 for (
unsigned int i = 0;
i < numberOfRows * numberOfColumns; ++
i)
277 return numberOfRows * numberOfColumns;
295 return numberOfColumns;
306 logger.assert_debug(j < numberOfColumns,
"Asked for column %, but there are only % columns", j, numberOfColumns);
313 logger.assert_debug(i < numberOfRows,
"Asked for row %, but there are only % rows", i, numberOfRows);
315 for (
unsigned int j = 0; j < numberOfColumns; ++j)
317 result[j] = (*this)(
i, j);
333 for (
unsigned int i = 0;
i < numberOfRows; ++
i)
335 for (
unsigned int j = 0; j < numberOfColumns; ++j)
337 result(j,
i) = (*this)(
i, j);
344 template<
unsigned int numberOfRightHandS
ideColumns>
363 std::array<Mdouble, numberOfRows * numberOfColumns>
data_;
367 template<
unsigned int numberOfRows,
unsigned int numberOfColumns>
368 std::ostream& operator<<(std::ostream& os, const SmallMatrix<numberOfRows, numberOfColumns>&
A)
370 for (
unsigned int i = 0;
i < numberOfRows; ++
i)
372 os <<
A.getRow(
i) << std::endl;
378 template<
unsigned int numberOfRows,
unsigned int numberOfColumns>
386 template<
unsigned int numberOfRows,
unsigned int numberOfColumns>
392 #endif //MERCURY_SMALLMATRIX_H
SmallMatrix(const SmallVector< numberOfRows > &other)
SmallMatrix operator*(const Mdouble &right) const
Does matrix A_ij=scalar*B_ij.
SmallMatrix & operator=(SmallMatrix &&right)
Assigns one matrix to another.
const Mdouble & operator()(unsigned int n, unsigned int m) const
defines the operator(n,m) to access the element on row n and column m
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
Definition of different loggers with certain modules. A user can define its own custom logger here...
SmallMatrix(SmallMatrix &&other)
Move Matrix from another Matrix.
const Mdouble * data() const
SmallMatrix(std::array< SmallVector< numberOfRows >, numberOfColumns > entries)
Glues one or more vectors with the same number of rows together.
SmallMatrix operator-(const SmallMatrix &other) const
SmallMatrix< numberOfColumns, numberOfRows > transpose() const
unsigned int size() const
Get total number of Matrix entries.
Mdouble & operator[](const unsigned int n)
Access the n linear element in the matrix.
const std::complex< Mdouble > i
Implementation of a 3D vector (by Vitaliy).
Mdouble & operator()(unsigned int n, unsigned int m)
defines the operator(n,m) to access the element on row n and column m
void axpy(Mdouble a, const SmallMatrix &x)
Applies the matrix y=ax + y, where x is another matrix and a is a scalar.
unsigned int getNumberOfColumns() const
Get the number of columns.
unsigned int getNCols() const
SmallVector< numberOfColumns > getRow(unsigned int i) const
get the i^th row
unsigned int getNumberOfRows() const
Get the number of rows.
std::array< Mdouble, numberOfRows *numberOfColumns > data_
The actually data of the matrix class.
SmallMatrix & operator/=(const Mdouble &scalar)
Does matrix A_ij=scalar*A_ij.
Mdouble determinant() const
const Mdouble & operator[](const unsigned int n) const
SmallMatrix & operator+=(const SmallMatrix &other)
const Mdouble * data() const
unsigned int getNRows() const
SmallVector< numberOfRows > getColumn(unsigned int j) const
get the j^th column
SmallMatrix< numberOfRows, numberOfColumns > operator*(const Mdouble d, const SmallMatrix< numberOfRows, numberOfColumns > &mat)
Multiplies a matrix with a Mdouble.
void solve(SmallMatrix< numberOfRows, numberOfRightHandSideColumns > &B) const
solves Ax=B where A is the current matrix and B is passed in. The result is returned in B...
SmallMatrix operator-() const
SmallMatrix & operator=(const SmallMatrix &right)
Assigns one matrix to another.
SmallMatrix inverse() const
return the inverse in the vector result. The size of result matches the matrix.
SmallMatrix()
Constructs a matrix of size n-rows by m-columns.
SmallMatrix & operator-=(const SmallMatrix &other)
SmallVector< numberOfRows > computeWedgeStuffVector() const
computeWedgeStuffVector.
SmallMatrix LUfactorisation() const
Return the LUfactorisation of the matrix.
SmallMatrix(const SmallMatrix &other)
Construct and copy Matrix from another Matrix i.e. B(A) where B and A are both matrices.
SmallMatrix & operator*=(const Mdouble &scalar)
Does matrix A_ij=scalar*A_ij.
SmallMatrix(const Mdouble &c)
Constructs a matrix of size n-rows by m-columns and initialises all entry to a constant.
Data type for small dense matrix.
SmallMatrix operator+(const SmallMatrix &other) const
SmallMatrix(const std::initializer_list< SmallVector< numberOfRows >> &entries)
SmallVector< numberOfRows > operator*(SmallVector< numberOfColumns > &right)
Defines Matrix A times vector B and return vector C i.e. C_,j= A_ij B_,j.
SmallMatrix operator/(const Mdouble &scalar) const
this does element by divided by a scalar