diff --git a/ChangeLog b/ChangeLog index 4f78bb9..e69de29 100644 --- a/ChangeLog +++ b/ChangeLog @@ -1,52 +0,0 @@ -2016-08-19 14:52 - - * COPYING, INSTALL, ax_cxx_compile_stdcxx_11.m4, - ax_init_standard_project.m4, bootstrap.sh, makefile_test.inc.am, - src/matrix.hxx, test/basic.cxx: more operators more tests - -2016-08-18 22:03 - - * COPYING, ChangeLog, INSTALL, src/matrix.hxx: more operators - -2016-08-17 07:26 - - * configure.ac: only requires c++11 - -2016-08-16 14:41 - - * COPYING, INSTALL, ax_init_standard_project.m4, configure.ac, - src/matrix.hxx, test/basic.cxx: remove redundancy, collect common - functionality in base class - -2016-08-08 20:03 - - * src/matrix.hxx, test/basic.cxx: more operator, more checks passed - -2016-08-03 18:43 - - * configure.ac, test/makefile.am: all tests passed - -2016-08-03 18:39 - - * COPYING, ChangeLog, INSTALL, ax_cxx_compile_stdcxx_11.m4, - ax_init_standard_project.m4, configure.ac, examples/makefile.am, - examples/matrix-sample.cxx, src/makefile.am, src/matrix.hxx, - test/basic.cxx, test/makefile.am: first approach including first - tests - -2016-07-30 08:50 - - * ., AUTHORS, ChangeLog, NEWS, README, autogen.sh, ax_check_qt.m4, - ax_cxx_compile_stdcxx_11.m4, ax_init_standard_project.m4, - bootstrap.sh, build-in-docker.conf, build-in-docker.sh, - build-resource-file.sh, configure.ac, debian, - debian/changelog.in, debian/compat, debian/control.in, - debian/docs, debian/libmatricxx-dev.install, - debian/libmatricxx.install, debian/rules, doc, doc/doxyfile.in, - doc/makefile.am, examples, examples/makefile.am, - libmatricxx.desktop.in, libmatricxx.spec.in, - mac-create-app-bundle.sh, makefile.am, resolve-debbuilddeps.sh, - resolve-rpmbuilddeps.sh, sql-to-dot.sed, src, - src/libmatricxx.pc.in, src/makefile.am, src/version.cxx, - src/version.hxx, test, test/makefile.am: initial project - diff --git a/src/matrix.hxx b/src/matrix.hxx index c331739..6fa2e34 100644 --- a/src/matrix.hxx +++ b/src/matrix.hxx @@ -177,26 +177,32 @@ template class MatrixBase { ///{ TYPE det() { - // rule of sarrus - TYPE res(0); - assert(ROWS==COLUMNS); // not really necessary - // 1) terms to add - for (size_t i(0); i=(TYPE*)(_c); p-=COLUMNS) res *= *p; return res; } + TYPE gauss() { + /// calculate using gauss algorithmus + /// @see http://www.mathebibel.de/determinante-berechnen-nach-gauss + /// 1. normalize first line to first value + TYPE lambda(at(0, 0)); + at(0, 0) = 1; + for (TYPE *p((TYPE*)(_c)+COLUMNS); p>(TYPE*)(_c)+1;) *--p/=lambda; + /// 2. nullify lower triangle + for (size_t column(0); column class Matrix: return res; } + static Matrix i() { + Matrix res; + for (size_t row(0); rowROWS; ++row) + if (rowsCOLUMNS) { + /// 2. normalize first line to first value + TYPE pivot(o(row, row)); + if (pivot!=1) { + o(row, row) = 1; + for (size_t column(row+1); columnCOLUMNS; ++column) o(row, column)/=pivot; + for (size_t column(row); columnCOLUMNS; ++column) this->at(row, column)/=pivot; + } + /// 3. nullify lower triangle + for (size_t row(column+1); rowROWS; ++row) { + TYPE pivot(o(row, column)); + if (pivot!=0) { + o(row, column) = 0; + for (size_t pos(column+1); posCOLUMNS; ++pos) o(row, pos) -= pivot*o(0, pos); + for (size_t pos(column); posCOLUMNS; ++pos) this->at(row, pos) -= pivot*this->at(0, pos); + } + } + } + /// 4. nullify the upper triangle + + return *this; + } + ///} }; @@ -306,6 +348,8 @@ template class Matrix: public MatrixBase { Matrix(size_t rows, size_t columns): Parent(rows, columns) { + assert(rows>0); + assert(columns>0); Parent::_c = new TYPE[rows*columns]; memset(Parent::_c, 0, Parent::MEM_SIZE); } @@ -347,14 +391,6 @@ template class Matrix: public MatrixBase { return Parent::operator=(o); } - Matrix t() const { - Matrix res(Parent::COLUMNS, Parent::ROWS); - for (size_t row(0); rowat(row, column); - return res; - } - Matrix operator-() const { Matrix res(Parent::COLUMNS, Parent::ROWS); for (TYPE *to((TYPE*)(res._c)+Parent::SIZE); to>(TYPE*)(Parent::_c); --to) @@ -374,16 +410,36 @@ template class Matrix: public MatrixBase { ///} + ///@name special operations + ///{ + + Matrix t() const { + Matrix res(this->COLUMNS, this->ROWS); + for (size_t row(0); rowROWS; ++row) + for (size_t column(0); columnCOLUMNS; ++column) + res(column, row) = this->at(row, column); + return res; + } + + Matrix i() const { + Matrix res(this->ROWS, this->COLUMNS); + for (size_t row(0); rowROWS&&rowCOLUMNS; ++row) + res(row, row) = 1; + return res; + } + + ///} + //................................................................methods protected: virtual void assert_check(const Matrix& o) const { - assert(o.ROWS==Parent::ROWS); - assert(o.COLUMNS==Parent::COLUMNS); + assert(o.ROWS==this->ROWS); + assert(o.COLUMNS==this->COLUMNS); } virtual bool check(const Matrix& o) const { - return o.ROWS==Parent::ROWS && o.COLUMNS==Parent::COLUMNS; + return o.ROWS==this->ROWS && o.COLUMNS==this->COLUMNS; } void copy_args(TYPE*) {} diff --git a/test/basic.cxx b/test/basic.cxx index a678871..3664281 100644 --- a/test/basic.cxx +++ b/test/basic.cxx @@ -169,6 +169,53 @@ class TemplateMatrixTest: public CppUnit::TestFixture { 3, 6); CPPUNIT_ASSERT_EQUAL(res, m.t()); } + template + void gauss() { + Matrix m(2, -2, 4, + -2, 1, -6, + 1, 0, -2); + const Matrix res(1, -1, 2, + 0, -1, -2, + 0, 0, -6); + T lambda(m.gauss()); + CPPUNIT_ASSERT_EQUAL(res, m); + CPPUNIT_ASSERT_EQUAL((T)2, lambda); + } + template + void det() { + Matrix m(2, -2, 4, + -2, 1, -6, + 1, 0, -2); + CPPUNIT_ASSERT_EQUAL((T)12, m.det()); + } + template + void i() { + const Matrix m1(1, 0, 0, + 0, 1, 0, + 0, 0, 1); + const Matrix m2(1, 0, + 0, 1, + 0, 0); + const Matrix m3(1, 0, 0, + 0, 1, 0); + CPPUNIT_ASSERT_EQUAL(m1, m1.i()); + CPPUNIT_ASSERT_EQUAL(m2, m2.i()); + CPPUNIT_ASSERT_EQUAL(m3, m3.i()); + } + template + void inv() { + Matrix m(2, -1, 0, + 1, 2, -2, + 0, -1, 1); + const Matrix res(0.5, 0, 0, + -0.2, 0.4, 0, + -1, 2, 5); + // const Matrix res(0, 1, 2, + // -1, 2, 4, + // -1, 2, 5); + m.inv(); + CPPUNIT_ASSERT_EQUAL(res, m); + } CPPUNIT_TEST_SUITE(TemplateMatrixTest); CPPUNIT_TEST(initFromArray1); CPPUNIT_TEST(initFromArray1); @@ -236,6 +283,22 @@ class TemplateMatrixTest: public CppUnit::TestFixture { CPPUNIT_TEST(transpose); CPPUNIT_TEST(transpose); CPPUNIT_TEST(transpose); + CPPUNIT_TEST(gauss); + CPPUNIT_TEST(gauss); + CPPUNIT_TEST(gauss); + CPPUNIT_TEST(gauss); + CPPUNIT_TEST(det); + CPPUNIT_TEST(det); + CPPUNIT_TEST(det); + CPPUNIT_TEST(det); + CPPUNIT_TEST(i); + CPPUNIT_TEST(i); + CPPUNIT_TEST(i); + CPPUNIT_TEST(i); + CPPUNIT_TEST(i); + CPPUNIT_TEST(i); + CPPUNIT_TEST(inv); + CPPUNIT_TEST(inv); CPPUNIT_TEST_SUITE_END(); }; CPPUNIT_TEST_SUITE_REGISTRATION(TemplateMatrixTest); @@ -400,6 +463,45 @@ class VariableMatrixTest: public CppUnit::TestFixture { 3, 6); CPPUNIT_ASSERT_EQUAL(res, m.t()); } + template + void gauss() { + Matrix m(3, 3, + 2, -2, 4, + -2, 1, -6, + 1, 0, -2); + const Matrix res(3,3, + 1, -1, 2, + 0, -1, -2, + 0, 0, -6); + T lambda(m.gauss()); + CPPUNIT_ASSERT_EQUAL(res, m); + CPPUNIT_ASSERT_EQUAL((T)2, lambda); + } + template + void det() { + Matrix m(3, 3, + 2, -2, 4, + -2, 1, -6, + 1, 0, -2); + CPPUNIT_ASSERT_EQUAL((T)12, m.det()); + } + template + void i() { + const Matrix m1(3, 3, + 1, 0, 0, + 0, 1, 0, + 0, 0, 1); + const Matrix m2(3, 2, + 1, 0, + 0, 1, + 0, 0); + const Matrix m3(2, 3, + 1, 0, 0, + 0, 1, 0); + CPPUNIT_ASSERT_EQUAL(m1, m1.i()); + CPPUNIT_ASSERT_EQUAL(m2, m2.i()); + CPPUNIT_ASSERT_EQUAL(m3, m3.i()); + } CPPUNIT_TEST_SUITE(VariableMatrixTest); CPPUNIT_TEST(initFromArray1); CPPUNIT_TEST(initFromArray1); @@ -467,6 +569,20 @@ class VariableMatrixTest: public CppUnit::TestFixture { CPPUNIT_TEST(transpose); CPPUNIT_TEST(transpose); CPPUNIT_TEST(transpose); + CPPUNIT_TEST(gauss); + CPPUNIT_TEST(gauss); + CPPUNIT_TEST(gauss); + CPPUNIT_TEST(gauss); + CPPUNIT_TEST(det); + CPPUNIT_TEST(det); + CPPUNIT_TEST(det); + CPPUNIT_TEST(det); + CPPUNIT_TEST(i); + CPPUNIT_TEST(i); + CPPUNIT_TEST(i); + CPPUNIT_TEST(i); + CPPUNIT_TEST(i); + CPPUNIT_TEST(i); CPPUNIT_TEST_SUITE_END(); }; CPPUNIT_TEST_SUITE_REGISTRATION(VariableMatrixTest);