CPP_utils/source/matrix.cpp at 45bed6e6c33744b842f6c2003565b177f3a360b4 · sci-code-711/CPP_utils · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
#include "matrix.hpp"

namespace cpp_utils {

    Matrix::Matrix(
        double xx, double xy, double xz,
        double yx, double yy, double yz,
        double zx, double zy, double zz
    ): xx(xx), xy(xy), xz(xz), yx(yx), yy(yy), yz(yz), zx(zx), zy(zy), zz(zz) {};

    Matrix::Matrix(double xx, double yy, double zz): xx(xx), xy(0), xz(0), yx(0), yy(yy), yz(0), zx(0), zy(0), zz(zz) {
        Matrix diagonalMatrix(xx, 0, 0, 0, yy, 0, 0, 0, zz);
    };

    double Matrix::det() const {
        return xx*(yy*zz-yz*zy)-xy*(yx*zz-yz*zx)+xz*(yx*zy-yy*zx);
    };

    bool Matrix::operator==(const Matrix &that) const {
        double testPrecision = std::fmax(precision, that.precision);
        return !(std::abs(xx - that.xx) > testPrecision ||
                 std::abs(xy - that.xy) > testPrecision ||
                 std::abs(xz - that.xz) > testPrecision ||
                 std::abs(yx - that.yx) > testPrecision ||
                 std::abs(yy - that.yy) > testPrecision ||
                 std::abs(yz - that.yz) > testPrecision ||
                 std::abs(zx - that.zx) > testPrecision ||
                 std::abs(zy - that.zy) > testPrecision ||
                 std::abs(zz - that.zz) > testPrecision);
    };

    Matrix Matrix::operator+(const Matrix &that) const {
        return Matrix(this->xx + that.xx, this->xy + that.xy, this->xz + that.xz,
                      this->yx + that.yx, this->yy + that.yy, this->yz + that.yz,
                      this->zx + that.zx, this->zy + that.zy, this->zz + that.zz);
    };

    Matrix Matrix::operator-(const Matrix &that) const {
        return Matrix(this->xx - that.xx, this->xy - that.xy, this->xz - that.xz,
                      this->yx - that.yx, this->yy - that.yy, this->yz - that.yz,
                      this->zx - that.zx, this->zy - that.zy, this->zz - that.zz);
    };

    Matrix Matrix::operator*(const double &that) const {
        return Matrix(xx * that, xy * that, xz * that, yx * that, yy * that, yz * that, zx * that, zy * that, zz * that);
    };

    Matrix operator*(const double that, const Matrix &those) {
        return those * that;
    };

    Matrix Matrix::operator*(const Matrix &that) const {
        return Matrix((xx * that.xx) + (xy * that.yx) + (xz * that.zx),
                      (xx * that.xy) + (xy * that.yy) + (xz * that.zy),
                      (xx * that.xz) + (xy * that.yz) + (xz * that.zz),
                      (yx * that.xx) + (yy * that.yx) + (yz * that.zx),
                      (yx * that.xy) + (yy * that.yy) + (yz * that.zy),
                      (yx * that.xz) + (yy * that.yz) + (yz * that.zz),
                      (zx * that.xx) + (zy * that.yx) + (zz * that.zx),
                      (zx * that.xy) + (zy * that.yy) + (zz * that.zy),
                      (zx * that.xz) + (zy * that.yz) + (zz * that.zz));
    };

    Vector Matrix::operator*(const Vector &that) const {
        return Vector((xx * that.x) + (xy * that.y) + (xz * that.z),
                      (yx * that.x) + (yy * that.y) + (yz * that.z),
                      (zx * that.x) + (zy * that.y) + (zz * that.z));
    };

    Vector operator*(const Vector that, const Matrix &those) {
        return Vector((that.x * those.xx) + (that.y * those.yx) + (that.z * those.zx),
                      (that.x * those.xy) + (that.y * those.yy) + (that.z * those.zy),
                      (that.x * those.xz) + (that.y * those.yz) + (that.z * those.zz));
    };

    Matrix Matrix::inverse() const {
        if (this->det() == 0) {
            throw value_error("No inverse due to zero value of determinant.");
        } else {
            return (1/det()) * Matrix( ((yy * zz) - (zy * yz)), -((xy * zz) - (zy * xz)),  ((xy * yz) - (yy * xz)),
                                      -((yx * zz) - (zx * yz)),  ((xx * zz) - (zx * xz)), -((xx * yz) - (yx * xz)),
                                       ((yx * zy) - (zx * yy)), -((xx * zy) - (zx * xy)),  ((xx * yy) - (yx * xy)));
        };
    };

}