aboutsummaryrefslogtreecommitdiffstats
path: root/meowpp/math/Matrix.h
diff options
context:
space:
mode:
Diffstat (limited to 'meowpp/math/Matrix.h')
-rw-r--r--meowpp/math/Matrix.h536
1 files changed, 0 insertions, 536 deletions
diff --git a/meowpp/math/Matrix.h b/meowpp/math/Matrix.h
deleted file mode 100644
index a2f45aa..0000000
--- a/meowpp/math/Matrix.h
+++ /dev/null
@@ -1,536 +0,0 @@
-#ifndef math_Matrix_H__
-#define math_Matrix_H__
-
-#include "../Self.h"
-#include "../math/utility.h"
-
-#include <vector>
-#include <algorithm>
-
-#include <cstdlib>
-
-namespace meow {
-/*!
- * @brief \b matrix
- *
- * @author cat_leopard
- */
-template<class Entry>
-class Matrix {
-public:
- typedef typename std::vector<Entry>::reference EntryRef ;
- typedef typename std::vector<Entry>::const_reference EntryRefK;
-private:
- struct Myself {
- size_t rows_;
- size_t cols_;
- std::vector<Entry> entries_;
-
- Myself():
- rows_(0), cols_(0), entries_(0) {
- }
-
- Myself(Myself const& b):
- rows_(b.rows_), cols_(b.cols_), entries_(b.entries_) {
- }
-
- Myself(size_t r, size_t c, Entry const& e):
- rows_(r), cols_(c), entries_(r * c, e) {
- }
-
- ~Myself() {
- }
-
- size_t index(size_t r, size_t c) const {
- return r * cols_ + c;
- }
-
- void realSize() {
- std::vector<Entry> tmp(entries_);
- entries_.swap(tmp);
- }
- };
-
- Self<Myself> const self;
-public:
- /*!
- * @brief constructor
- *
- * Create an empty matrix with size \b 0x0.
- * In other world, create an \b invalid matrix
- */
- Matrix(): self() { }
-
- /*!
- * @brief constructor
- *
- * Copy data from another one
- *
- * @param [in] m another matrix
- */
- Matrix(Matrix const& m): self(m.self, Self<Myself>::COPY_FROM) {
- }
-
- /*!
- * @brief constructor
- *
- * Create an \a r x \a c matrix with all entry be \a e
- *
- * @param [in] r number of rows
- * @param [in] c number of columns
- * @param [in] e inital entry
- */
- Matrix(size_t r, size_t c, Entry const& e): self(Myself(r, c, e)) {
- }
-
- //! @brief destructor
- ~Matrix() { }
-
- /*!
- * @brief copy
- *
- * Copy data from another matrix
- *
- * @param [in] m matrix
- * @return *this
- */
- Matrix& copyFrom(Matrix const& m) {
- self().copyFrom(m.self);
- return *this;
- }
-
- /*!
- * @brief reference
- *
- * Reference itself to another matrix
- *
- * @param [in] m matrix
- * @return *this
- */
- Matrix& referenceFrom(Matrix const& m) {
- self().referenceFrom(m.self);
- return *this;
- }
-
- //! @brief reset the size of the matrix to \a r x \a c with entry all be \a e
- void reset(size_t r, size_t c, Entry const& e) {
- self()->rows_ = r;
- self()->cols_ = c;
- self()->entries_.clear();
- self()->entries_.resize(r * c, e);
- }
-
- //! @brief Return whether it is a \b valid matrix
- bool valid() const {
- return (rows() > 0 && cols() > 0);
- }
-
- //! @brief Return number of rows
- size_t rows() const {
- return self->rows_;
- }
-
- //! @brief Return number of cols
- size_t cols() const {
- return self->cols_;
- }
-
- //! @brief Return number of rows times number of cols
- size_t size() const {
- return rows() * cols();
- }
-
- /*!
- * @brief resize the matrix such that number of rows become \a r.
- *
- * New created entry will be \a e
- *
- * @param [in] r new number of rows
- * @param [in] e inital entry
- * @return new number of rows
- */
- size_t rows(size_t r, Entry const& e) {
- if (r != rows()) {
- self()->entries_.resize(r * cols(), e);
- self()->rows_ = r;
- }
- return rows();
- }
-
- /*!
- * @brief resize the matrix such that number of cols become \a c
- *
- * New created entry will be \a e
- *
- * @param [in] c new number of columns
- * @param [in] e inital entry
- * @return new number of columns
- */
- size_t cols(size_t c, Entry const& e) {
- if (c != cols()) {
- Self<Myself> const old(self, Self<Myself>::COPY_FROM);
- self()->entries_.resize(rows() * c);
- self()->cols_ = c;
- for (size_t i = 0, I = rows(); i < I; i++) {
- size_t j, J1 = std::min(old->cols_, cols()), J2 = cols();
- for (j = 0; j < J1; j++)
- self()->entries_[self->index(i, j)] = old->entries_[old->index(i, j)];
- for (j = J1; j < J2; j++)
- self()->entries_[self->index(i, j)] = e;
- }
- }
- return cols();
- }
-
- /*!
- * @brief resize
- *
- * Resize to \a r x \a c, with new created entry be \a e
- *
- * @param [in] r number of rows
- * @param [in] c number of rows
- * @param [in] e inital entry
- * @return \a r * \a c
- */
- size_t size(size_t r, size_t c, Entry const& e) {
- cols(c, e);
- rows(r, e);
- return rows() * cols();
- }
-
- /*!
- * @brief free the memory
- */
- void clear() {
- self()->rows_ = 0;
- self()->cols_ = 0;
- self()->entries_.clear();
- self()->realSize();
- }
-
- //! @brief Access the entry at \a r x \a c
- Entry entry(size_t r, size_t c) const {
- return self->entries_[self->index(r, c)];
- }
-
- //! @brief Change the entry at \a r x \a c
- Entry entry(size_t r, size_t c, Entry const& e) {
- self()->entries_[self->index(r, c)] = e;
- return entry(r, c);
- }
-
- //! @brief Get the entry at \a r x \a c
- EntryRef entryGet(size_t r, size_t c) {
- return self()->entries_[self->index(r, c)];
- }
-
- /*!
- * @brief Change the entries from \a rFirst x \a cFirst to \a rLast x \a cLast
- *
- * @param [in] rFirst
- * @param [in] rLast
- * @param [in] cFirst
- * @param [in] cLast
- * @param [in] e value
- * @return void
- */
- void entries(ssize_t rFirst, ssize_t rLast,
- ssize_t cFirst, ssize_t cLast,
- Entry const& e) {
- for (ssize_t r = rFirst; r <= rLast; r++) {
- for (ssize_t c = cFirst; c <=cFirst; c++) {
- entry(r, c, e);
- }
- }
- }
-
- /*!
- * @brief Return a \a rLast-rFirst+1 x \a cLast-cFirst+1 matrix
- *
- * With value be the entries from \a rFirst x \a cFirst to \a rLast x \a cLast
- *
- * @param [in] rFirst
- * @param [in] rLast
- * @param [in] cFirst
- * @param [in] cLast
- * @return a matrix
- */
- Matrix subMatrix(size_t rFirst, size_t rLast,
- size_t cFirst, size_t cLast) const {
- if (rFirst > rLast || cFirst > cLast) return Matrix();
- if (rFirst == 0 && cFirst == 0) {
- Matrix ret(*this);
- ret.size(rLast + 1, cLast + 1, Entry(0));
- return ret;
- }
- Matrix ret(rLast - rFirst + 1, cLast - cFirst + 1, entry(rFirst, cFirst));
- for (size_t r = rFirst; r <= rLast; r++)
- for (size_t c = cFirst; c <= cLast; c++)
- ret.entry(r - rFirst, c - cFirst, entry(r, c));
- return ret;
- }
-
- //! @brief Return the \a r -th row
- Matrix row(size_t r) const {
- return subMatrix(r, r, 0, cols() - 1);
- }
-
- //! @brief Return the \a c -th column
- Matrix col(size_t c) const {
- return subMatrix(0, rows() - 1, c, c);
- }
-
- //! @brief return +\a (*this)
- Matrix positive() const {
- return *this;
- }
-
- //! @brief return -\a (*this)
- Matrix negative() const {
- Matrix ret(*this);
- for (size_t r = 0, R = rows(); r < R; r++)
- for (size_t c = 0, C = cols(); c < C; c++)
- ret.entry(r, c, -ret.entry(r, c));
- return ret;
- }
-
- /*! @brief return \a (*this) + \a m.
- *
- * If the size not match, it will return an invalid matrix
- */
- Matrix add(Matrix const& m) const {
- if (rows() != m.rows() || cols() != m.cols()) return Matrix();
- Matrix ret(*this);
- for (size_t r = 0, R = rows(); r < R; r++)
- for (size_t c = 0, C = cols(); c < C; c++)
- ret.entry(r, c, ret.entry(r, c) + m.entry(r, c));
- return ret;
- }
-
- /*! @brief return \a (*this) - \a m.
- *
- * If the size not match, it will return an invalid matrix
- */
- Matrix sub(Matrix const& m) const {
- if (rows() != m.rows() || cols() != m.cols()) return Matrix();
- Matrix ret(*this);
- for (size_t r = 0, R = rows(); r < R; r++)
- for (size_t c = 0, C = cols(); c < C; c++)
- ret.entry(r, c, ret.entry(r, c) - m.entry(r, c));
- return ret;
- }
-
- /*! @brief return \a (*this) times \a m.
- *
- * If the size not match, it will return an invalid matrix
- */
- Matrix mul(Matrix const& m) const {
- if (cols() != m.rows()) return Matrix();
- Matrix ret(rows(), m.cols(), Entry(0));
- for (size_t r = 0, R = rows(); r < R; r++)
- for (size_t c = 0, C = m.cols(); c < C; c++)
- for (size_t k = 0, K = cols(); k < K; k++)
- ret.entry(r, c, ret.entry(r, c) + entry(r, k) * m.entry(k, c));
- return ret;
- }
-
- //! @brief return \a (*this) times \a s. \a s is a scalar
- Matrix mul(Entry const& s) const {
- Matrix ret(*this);
- for (size_t r = 0, R = rows(); r < R; r++)
- for (size_t c = 0, C = cols(); c < C; c++)
- ret.entry(r, c, ret.entry(r, c) * s);
- return ret;
- }
-
- //! @brief return \a (*this) / \a s. \a s is a scalar
- Matrix div(Entry const& s) const {
- Matrix ret(*this);
- for (size_t r = 0, R = rows(); r < R; r++)
- for (size_t c = 0, C = cols(); c < C; c++)
- ret.entry(r, c, ret.entry(r, c) / s);
- return ret;
- }
-
- //! @brief Return a identity matrix with size equal to itself
- Matrix identity() const {
- Matrix ret(*this);
- ret.identitied();
- return ret;
- }
-
- /*!
- * @brief Let itself be an identity matrix
- *
- * Our definition of Identity matrix is 1 for entry(i, i) and 0 otherwise.
- */
- Matrix& identitied() {
- for (size_t r = 0, R = rows(); r < R; r++)
- for (size_t c = 0, C = cols(); c < C; c++)
- entry(r, c, (r == c ? Entry(1) : Entry(0)));
- return *this;
- }
-
- /*!
- * @brief Let itself be an diagonal form of original itself
- */
- Matrix& diagonaled() {
- triangulared();
- for (size_t i = 0, I = rows(); i < I; ++i) {
- for (size_t j = i + 1, J = cols(); j < J; ++j) {
- entry(i, j, Entry(0));
- }
- }
- return *this;
- }
-
- /*!
- * @brief Return a matrix which is a diangonal form of me
- */
- Matrix diagonal() const {
- Matrix ret(*this);
- ret.diagonaled();
- return ret;
- }
-
- /*!
- * @brief Return a matrix which is an inverse matrix of \a (*this)
- *
- * If inverse matrix doesn't exist, it will return a invalid matrix
- */
- Matrix inverse() const {
- if (rows() != cols() || rows() == 0) return Matrix<Entry>();
- Matrix tmp(rows(), cols() * 2, Entry(0));
- for (size_t r = 0, R = rows(); r < R; r++) {
- for (size_t c = 0, C = cols(); c < C; c++) {
- tmp.entry(r, c, entry(r, c));
- tmp.entry(r, c + cols(), (r == c ? Entry(1) : Entry(0)));
- }
- }
- tmp.triangulared();
- for (ssize_t r = rows() - 1; r >= 0; r--) {
- if (tmp(r, r) == Entry(0)) return Matrix<Entry>();
- for (ssize_t r2 = r - 1; r2 >= 0; r2--) {
- Entry rat(-tmp.entry(r2, r) / tmp.entry(r, r));
- for (size_t c = r, C = tmp.cols(); c < C; c++) {
- tmp.entry(r2, c, tmp.entry(r2, c) + rat * tmp(r, c));
- }
- }
- Entry rat(tmp.entry(r, r));
- for (size_t c = cols(), C = tmp.cols(); c < C; c++) {
- tmp.entry(r, c - cols(), tmp.entry(r, c) / rat);
- }
- }
- tmp.size(cols(), rows(), Entry(0));
- return tmp;
- }
-
- //! @brief let itself become itself's inverse matrix
- Matrix& inversed() {
- copyFrom(inverse());
- return *this;
- }
-
- //! @brief return itself's transpose matrix
- Matrix transpose() const {
- Matrix ret(cols(), rows(), Entry(0));
- for (size_t r = 0, R = cols(); r < R; r++)
- for (size_t c = 0, C = rows(); c < C; c++)
- ret.entry(r, c, entry(c, r));
- return ret;
- }
-
- //! @brief Let itself become itself's transpose matrix
- Matrix& transposed() {
- copyFrom(transpose());
- return *this;
- }
-
- //! @brief return a matrix which is the triangular form of \a (*this)
- Matrix triangular() const {
- Matrix<Entry> ret(*this);
- ret.triangulared();
- return ret;
- }
-
- //! @brief triangluar itself
- Matrix& triangulared() {
- for (size_t r = 0, c = 0, R = rows(), C = cols(); r < R && c < C; r++) {
- ssize_t maxR;
- for ( ; c < C; c++) {
- maxR = -1;
- for (size_t r2 = r; r2 < R; r2++)
- if (maxR == -1 || tAbs(entry(r2, c)) > tAbs(entry(maxR, c)))
- maxR = r2;
- if (entry(maxR, c) != Entry(0)) break;
- }
- if (c >= C) break;
- if (maxR != (ssize_t)r) {
- for (size_t c2 = c; c2 < C; c2++)
- std::swap(self()->entries_[self->index( r, c2)],
- self()->entries_[self->index(maxR, c2)]);
- }
- for (size_t r2 = r + 1; r2 < R; r2++) {
- Entry rati = -entry(r2, c) / entry(r, c);
- entry(r2, c, Entry(0));
- for (size_t c2 = c + 1; c2 < C; c2++)
- entry(r2, c2, entry(r2, c2) + entry(r, c2) * rati);
- }
- }
- return *this;
- }
-
- //! @brief same as \a copyFrom
- Matrix& operator=(Matrix const& m) {
- return copyFrom(m);
- }
-
- //! @brief same as \a entry(r,c)
- Entry operator()(size_t r, size_t c) const {
- return entry(r, c);
- }
-
- //! @brief same as \a entry(r,c,e)
- Entry operator()(size_t r, size_t c, Entry const& e) {
- return entry(r, c, e);
- }
-
- //! @brief same as \a positive()
- Matrix operator+() const {
- return positive();
- }
-
- //! @brief same as \a negative()
- Matrix operator-() const {
- return negative();
- }
-
- //! @brief same as \a add(m)
- Matrix operator+(Matrix const& m) const {
- return add(m);
- }
-
- //! @brief same as \a sub(m)
- Matrix operator-(Matrix const& m) const {
- return sub(m);
- }
-
- //! @brief same as \a mul(m)
- Matrix operator*(Matrix const& m) const {
- return mul(m);
- }
-
- //! @brief same as \a mul(m)
- Matrix operator*(Entry const& s) const {
- return mul(s);
- }
-
- //! @brief same as \a div(s)
- Matrix operator/(Entry const& s) const {
- return div(s);
- }
-};
-
-} // meow
-
-#endif // math_Matrix_H__