Alex Pirciu
a
include/Linalg/linalg.h@1:ceee5a608e7c, 2019-03-28 (annotated)
- Committer:
- alexpirciu
- Date:
- Thu Mar 28 07:44:42 2019 +0000
- Revision:
- 1:ceee5a608e7c
assa
Who changed what in which revision?
| User | Revision | Line number | New contents of line |
|---|---|---|---|
| alexpirciu | 1:ceee5a608e7c | 1 | #ifndef LINALG_H |
| alexpirciu | 1:ceee5a608e7c | 2 | #define LINALG_H |
| alexpirciu | 1:ceee5a608e7c | 3 | |
| alexpirciu | 1:ceee5a608e7c | 4 | // #include <mbed.h> |
| alexpirciu | 1:ceee5a608e7c | 5 | #include <stdint.h> |
| alexpirciu | 1:ceee5a608e7c | 6 | #include <array> |
| alexpirciu | 1:ceee5a608e7c | 7 | |
| alexpirciu | 1:ceee5a608e7c | 8 | namespace linalg |
| alexpirciu | 1:ceee5a608e7c | 9 | { |
| alexpirciu | 1:ceee5a608e7c | 10 | template <class T, uint32_t M, uint32_t N> |
| alexpirciu | 1:ceee5a608e7c | 11 | class CMatrix |
| alexpirciu | 1:ceee5a608e7c | 12 | { |
| alexpirciu | 1:ceee5a608e7c | 13 | public: |
| alexpirciu | 1:ceee5a608e7c | 14 | using CThisType = CMatrix<T,M,N>; |
| alexpirciu | 1:ceee5a608e7c | 15 | using CTransposeType = CMatrix<T,N,M>; |
| alexpirciu | 1:ceee5a608e7c | 16 | using CContainerType = std::array<std::array<T,N>,M>; |
| alexpirciu | 1:ceee5a608e7c | 17 | using CDataType =T; |
| alexpirciu | 1:ceee5a608e7c | 18 | |
| alexpirciu | 1:ceee5a608e7c | 19 | template <uint32_t P> |
| alexpirciu | 1:ceee5a608e7c | 20 | using CLeftMultipliableType = CMatrix<T,P,M>; |
| alexpirciu | 1:ceee5a608e7c | 21 | |
| alexpirciu | 1:ceee5a608e7c | 22 | template <uint32_t P> |
| alexpirciu | 1:ceee5a608e7c | 23 | using CRightMultipliableType = CMatrix<T,N,P>; |
| alexpirciu | 1:ceee5a608e7c | 24 | |
| alexpirciu | 1:ceee5a608e7c | 25 | template <uint32_t P> |
| alexpirciu | 1:ceee5a608e7c | 26 | using CLeftMultiplicationResultType = CMatrix<T,P,N>; |
| alexpirciu | 1:ceee5a608e7c | 27 | |
| alexpirciu | 1:ceee5a608e7c | 28 | template <uint32_t P> |
| alexpirciu | 1:ceee5a608e7c | 29 | using CRightMultiplicationResultType = CMatrix<T,M,P>; |
| alexpirciu | 1:ceee5a608e7c | 30 | |
| alexpirciu | 1:ceee5a608e7c | 31 | friend class CMatrix<T,N,M>; |
| alexpirciu | 1:ceee5a608e7c | 32 | |
| alexpirciu | 1:ceee5a608e7c | 33 | CMatrix() : m_data() {} |
| alexpirciu | 1:ceee5a608e7c | 34 | CMatrix(const CThisType& f_matrix) : m_data(f_matrix.m_data) {} |
| alexpirciu | 1:ceee5a608e7c | 35 | CMatrix(const CThisType&& f_matrix) : m_data(f_matrix.m_data) {} |
| alexpirciu | 1:ceee5a608e7c | 36 | CMatrix(const CContainerType& f_data) : m_data(f_data) {} |
| alexpirciu | 1:ceee5a608e7c | 37 | CMatrix(const CContainerType&& f_data) : m_data(f_data) {} |
| alexpirciu | 1:ceee5a608e7c | 38 | |
| alexpirciu | 1:ceee5a608e7c | 39 | CThisType& operator=(const CThisType& f_matrix) |
| alexpirciu | 1:ceee5a608e7c | 40 | { |
| alexpirciu | 1:ceee5a608e7c | 41 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 42 | { |
| alexpirciu | 1:ceee5a608e7c | 43 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 44 | { |
| alexpirciu | 1:ceee5a608e7c | 45 | this->m_data[l_row][l_col] = f_matrix.m_data[l_row][l_col]; |
| alexpirciu | 1:ceee5a608e7c | 46 | } |
| alexpirciu | 1:ceee5a608e7c | 47 | } |
| alexpirciu | 1:ceee5a608e7c | 48 | return *this; |
| alexpirciu | 1:ceee5a608e7c | 49 | } |
| alexpirciu | 1:ceee5a608e7c | 50 | CThisType& operator=(const CThisType&& f_matrix) |
| alexpirciu | 1:ceee5a608e7c | 51 | { |
| alexpirciu | 1:ceee5a608e7c | 52 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 53 | { |
| alexpirciu | 1:ceee5a608e7c | 54 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 55 | { |
| alexpirciu | 1:ceee5a608e7c | 56 | this->m_data[l_row][l_col] = f_matrix.m_data[l_row][l_col]; |
| alexpirciu | 1:ceee5a608e7c | 57 | } |
| alexpirciu | 1:ceee5a608e7c | 58 | } |
| alexpirciu | 1:ceee5a608e7c | 59 | return *this; |
| alexpirciu | 1:ceee5a608e7c | 60 | } |
| alexpirciu | 1:ceee5a608e7c | 61 | |
| alexpirciu | 1:ceee5a608e7c | 62 | std::array<T,N>& operator[](uint32_t f_row) |
| alexpirciu | 1:ceee5a608e7c | 63 | { |
| alexpirciu | 1:ceee5a608e7c | 64 | return m_data[f_row]; |
| alexpirciu | 1:ceee5a608e7c | 65 | } |
| alexpirciu | 1:ceee5a608e7c | 66 | |
| alexpirciu | 1:ceee5a608e7c | 67 | CDataType& operator()(uint32_t f_row, uint32_t f_col) |
| alexpirciu | 1:ceee5a608e7c | 68 | { |
| alexpirciu | 1:ceee5a608e7c | 69 | return m_data[f_row][f_col]; |
| alexpirciu | 1:ceee5a608e7c | 70 | } |
| alexpirciu | 1:ceee5a608e7c | 71 | |
| alexpirciu | 1:ceee5a608e7c | 72 | const std::array<T,N>& operator[](uint32_t f_row) const |
| alexpirciu | 1:ceee5a608e7c | 73 | { |
| alexpirciu | 1:ceee5a608e7c | 74 | return m_data[f_row]; |
| alexpirciu | 1:ceee5a608e7c | 75 | } |
| alexpirciu | 1:ceee5a608e7c | 76 | |
| alexpirciu | 1:ceee5a608e7c | 77 | const CDataType& operator()(uint32_t f_row, uint32_t f_col) const |
| alexpirciu | 1:ceee5a608e7c | 78 | { |
| alexpirciu | 1:ceee5a608e7c | 79 | return m_data[f_row][f_col]; |
| alexpirciu | 1:ceee5a608e7c | 80 | } |
| alexpirciu | 1:ceee5a608e7c | 81 | |
| alexpirciu | 1:ceee5a608e7c | 82 | // template<> |
| alexpirciu | 1:ceee5a608e7c | 83 | // CMatrix<T,1,1>::operator CDataType() |
| alexpirciu | 1:ceee5a608e7c | 84 | // { |
| alexpirciu | 1:ceee5a608e7c | 85 | // return m_data[0][0]; |
| alexpirciu | 1:ceee5a608e7c | 86 | // } |
| alexpirciu | 1:ceee5a608e7c | 87 | |
| alexpirciu | 1:ceee5a608e7c | 88 | CThisType& operator+() {return *this;} |
| alexpirciu | 1:ceee5a608e7c | 89 | CThisType operator-() |
| alexpirciu | 1:ceee5a608e7c | 90 | { |
| alexpirciu | 1:ceee5a608e7c | 91 | CThisType l_matrix; |
| alexpirciu | 1:ceee5a608e7c | 92 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 93 | { |
| alexpirciu | 1:ceee5a608e7c | 94 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 95 | { |
| alexpirciu | 1:ceee5a608e7c | 96 | l_matrix.m_data[l_row][l_col] = -this->m_data[l_row][l_col]; |
| alexpirciu | 1:ceee5a608e7c | 97 | } |
| alexpirciu | 1:ceee5a608e7c | 98 | } |
| alexpirciu | 1:ceee5a608e7c | 99 | return l_matrix; |
| alexpirciu | 1:ceee5a608e7c | 100 | } |
| alexpirciu | 1:ceee5a608e7c | 101 | CThisType operator+(const CThisType& f_matrix) |
| alexpirciu | 1:ceee5a608e7c | 102 | { |
| alexpirciu | 1:ceee5a608e7c | 103 | CThisType l_matrix; |
| alexpirciu | 1:ceee5a608e7c | 104 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 105 | { |
| alexpirciu | 1:ceee5a608e7c | 106 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 107 | { |
| alexpirciu | 1:ceee5a608e7c | 108 | l_matrix.m_data[l_row][l_col] = this->m_data[l_row][l_col] + f_matrix.m_data[l_row][l_col]; |
| alexpirciu | 1:ceee5a608e7c | 109 | } |
| alexpirciu | 1:ceee5a608e7c | 110 | } |
| alexpirciu | 1:ceee5a608e7c | 111 | return l_matrix; |
| alexpirciu | 1:ceee5a608e7c | 112 | } |
| alexpirciu | 1:ceee5a608e7c | 113 | CThisType operator-(const CThisType& f_matrix) |
| alexpirciu | 1:ceee5a608e7c | 114 | { |
| alexpirciu | 1:ceee5a608e7c | 115 | CThisType l_matrix; |
| alexpirciu | 1:ceee5a608e7c | 116 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 117 | { |
| alexpirciu | 1:ceee5a608e7c | 118 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 119 | { |
| alexpirciu | 1:ceee5a608e7c | 120 | l_matrix.m_data[l_row][l_col] = this->m_data[l_row][l_col] - f_matrix.m_data[l_row][l_col]; |
| alexpirciu | 1:ceee5a608e7c | 121 | } |
| alexpirciu | 1:ceee5a608e7c | 122 | } |
| alexpirciu | 1:ceee5a608e7c | 123 | return l_matrix; |
| alexpirciu | 1:ceee5a608e7c | 124 | } |
| alexpirciu | 1:ceee5a608e7c | 125 | CThisType& operator+=(const CThisType& f_matrix) |
| alexpirciu | 1:ceee5a608e7c | 126 | { |
| alexpirciu | 1:ceee5a608e7c | 127 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 128 | { |
| alexpirciu | 1:ceee5a608e7c | 129 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 130 | { |
| alexpirciu | 1:ceee5a608e7c | 131 | this->m_data[l_row][l_col] += f_matrix.m_data[l_row][l_col]; |
| alexpirciu | 1:ceee5a608e7c | 132 | } |
| alexpirciu | 1:ceee5a608e7c | 133 | } |
| alexpirciu | 1:ceee5a608e7c | 134 | return *this; |
| alexpirciu | 1:ceee5a608e7c | 135 | } |
| alexpirciu | 1:ceee5a608e7c | 136 | CThisType& operator-=(const CThisType& f_matrix) |
| alexpirciu | 1:ceee5a608e7c | 137 | { |
| alexpirciu | 1:ceee5a608e7c | 138 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 139 | { |
| alexpirciu | 1:ceee5a608e7c | 140 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 141 | { |
| alexpirciu | 1:ceee5a608e7c | 142 | this->m_data[l_row][l_col] -= f_matrix.m_data[l_row][l_col]; |
| alexpirciu | 1:ceee5a608e7c | 143 | } |
| alexpirciu | 1:ceee5a608e7c | 144 | } |
| alexpirciu | 1:ceee5a608e7c | 145 | return *this; |
| alexpirciu | 1:ceee5a608e7c | 146 | } |
| alexpirciu | 1:ceee5a608e7c | 147 | |
| alexpirciu | 1:ceee5a608e7c | 148 | CThisType operator+(const CDataType& f_val) |
| alexpirciu | 1:ceee5a608e7c | 149 | { |
| alexpirciu | 1:ceee5a608e7c | 150 | CThisType l_matrix; |
| alexpirciu | 1:ceee5a608e7c | 151 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 152 | { |
| alexpirciu | 1:ceee5a608e7c | 153 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 154 | { |
| alexpirciu | 1:ceee5a608e7c | 155 | l_matrix.m_data[l_row][l_col] = this->m_data[l_row][l_col] + f_val; |
| alexpirciu | 1:ceee5a608e7c | 156 | } |
| alexpirciu | 1:ceee5a608e7c | 157 | } |
| alexpirciu | 1:ceee5a608e7c | 158 | return l_matrix; |
| alexpirciu | 1:ceee5a608e7c | 159 | } |
| alexpirciu | 1:ceee5a608e7c | 160 | CThisType operator-(const CDataType& f_val) |
| alexpirciu | 1:ceee5a608e7c | 161 | { |
| alexpirciu | 1:ceee5a608e7c | 162 | CThisType l_matrix; |
| alexpirciu | 1:ceee5a608e7c | 163 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 164 | { |
| alexpirciu | 1:ceee5a608e7c | 165 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 166 | { |
| alexpirciu | 1:ceee5a608e7c | 167 | l_matrix.m_data[l_row][l_col] = this->m_data[l_row][l_col] - f_val; |
| alexpirciu | 1:ceee5a608e7c | 168 | } |
| alexpirciu | 1:ceee5a608e7c | 169 | } |
| alexpirciu | 1:ceee5a608e7c | 170 | return l_matrix; |
| alexpirciu | 1:ceee5a608e7c | 171 | } |
| alexpirciu | 1:ceee5a608e7c | 172 | CThisType& operator+=(const CDataType& f_val) |
| alexpirciu | 1:ceee5a608e7c | 173 | { |
| alexpirciu | 1:ceee5a608e7c | 174 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 175 | { |
| alexpirciu | 1:ceee5a608e7c | 176 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 177 | { |
| alexpirciu | 1:ceee5a608e7c | 178 | this->m_data[l_row][l_col] += f_val; |
| alexpirciu | 1:ceee5a608e7c | 179 | } |
| alexpirciu | 1:ceee5a608e7c | 180 | } |
| alexpirciu | 1:ceee5a608e7c | 181 | return *this; |
| alexpirciu | 1:ceee5a608e7c | 182 | } |
| alexpirciu | 1:ceee5a608e7c | 183 | CThisType& operator-=(const CDataType& f_val) |
| alexpirciu | 1:ceee5a608e7c | 184 | { |
| alexpirciu | 1:ceee5a608e7c | 185 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 186 | { |
| alexpirciu | 1:ceee5a608e7c | 187 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 188 | { |
| alexpirciu | 1:ceee5a608e7c | 189 | this->m_data[l_row][l_col] -= f_val; |
| alexpirciu | 1:ceee5a608e7c | 190 | } |
| alexpirciu | 1:ceee5a608e7c | 191 | } |
| alexpirciu | 1:ceee5a608e7c | 192 | return *this; |
| alexpirciu | 1:ceee5a608e7c | 193 | } |
| alexpirciu | 1:ceee5a608e7c | 194 | CThisType& operator*=(const CDataType& f_val) |
| alexpirciu | 1:ceee5a608e7c | 195 | { |
| alexpirciu | 1:ceee5a608e7c | 196 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 197 | { |
| alexpirciu | 1:ceee5a608e7c | 198 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 199 | { |
| alexpirciu | 1:ceee5a608e7c | 200 | this->m_data[l_row][l_col] *= f_val; |
| alexpirciu | 1:ceee5a608e7c | 201 | } |
| alexpirciu | 1:ceee5a608e7c | 202 | } |
| alexpirciu | 1:ceee5a608e7c | 203 | return *this; |
| alexpirciu | 1:ceee5a608e7c | 204 | } |
| alexpirciu | 1:ceee5a608e7c | 205 | CThisType& operator*=(const CThisType& f_val) |
| alexpirciu | 1:ceee5a608e7c | 206 | { |
| alexpirciu | 1:ceee5a608e7c | 207 | CThisType& l_thisRef(*this); |
| alexpirciu | 1:ceee5a608e7c | 208 | l_thisRef = l_thisRef * f_val; |
| alexpirciu | 1:ceee5a608e7c | 209 | return l_thisRef; |
| alexpirciu | 1:ceee5a608e7c | 210 | } |
| alexpirciu | 1:ceee5a608e7c | 211 | CThisType& operator/=(const CDataType& f_val) |
| alexpirciu | 1:ceee5a608e7c | 212 | { |
| alexpirciu | 1:ceee5a608e7c | 213 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 214 | { |
| alexpirciu | 1:ceee5a608e7c | 215 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 216 | { |
| alexpirciu | 1:ceee5a608e7c | 217 | this->m_data[l_row][l_col] /= f_val; |
| alexpirciu | 1:ceee5a608e7c | 218 | } |
| alexpirciu | 1:ceee5a608e7c | 219 | } |
| alexpirciu | 1:ceee5a608e7c | 220 | return *this; |
| alexpirciu | 1:ceee5a608e7c | 221 | } |
| alexpirciu | 1:ceee5a608e7c | 222 | CThisType operator*(const CDataType& f_val) |
| alexpirciu | 1:ceee5a608e7c | 223 | { |
| alexpirciu | 1:ceee5a608e7c | 224 | CThisType l_matrix; |
| alexpirciu | 1:ceee5a608e7c | 225 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 226 | { |
| alexpirciu | 1:ceee5a608e7c | 227 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 228 | { |
| alexpirciu | 1:ceee5a608e7c | 229 | l_matrix.m_data[l_row][l_col] = this->m_data[l_row][l_col] * f_val; |
| alexpirciu | 1:ceee5a608e7c | 230 | } |
| alexpirciu | 1:ceee5a608e7c | 231 | } |
| alexpirciu | 1:ceee5a608e7c | 232 | return l_matrix; |
| alexpirciu | 1:ceee5a608e7c | 233 | } |
| alexpirciu | 1:ceee5a608e7c | 234 | template <uint32_t P> |
| alexpirciu | 1:ceee5a608e7c | 235 | CRightMultiplicationResultType<P> operator*(const CRightMultipliableType<P>& f_matrix) |
| alexpirciu | 1:ceee5a608e7c | 236 | { |
| alexpirciu | 1:ceee5a608e7c | 237 | CRightMultiplicationResultType<P> l_matrix; |
| alexpirciu | 1:ceee5a608e7c | 238 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 239 | { |
| alexpirciu | 1:ceee5a608e7c | 240 | for (uint32_t l_col = 0; l_col < P; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 241 | { |
| alexpirciu | 1:ceee5a608e7c | 242 | l_matrix[l_row][l_col] = 0; |
| alexpirciu | 1:ceee5a608e7c | 243 | for (uint32_t l_idx = 0; l_idx < N; ++l_idx) |
| alexpirciu | 1:ceee5a608e7c | 244 | { |
| alexpirciu | 1:ceee5a608e7c | 245 | l_matrix[l_row][l_col] += this->m_data[l_row][l_idx] * f_matrix[l_idx][l_col]; |
| alexpirciu | 1:ceee5a608e7c | 246 | } |
| alexpirciu | 1:ceee5a608e7c | 247 | } |
| alexpirciu | 1:ceee5a608e7c | 248 | } |
| alexpirciu | 1:ceee5a608e7c | 249 | return l_matrix; |
| alexpirciu | 1:ceee5a608e7c | 250 | } |
| alexpirciu | 1:ceee5a608e7c | 251 | CThisType inv(); |
| alexpirciu | 1:ceee5a608e7c | 252 | |
| alexpirciu | 1:ceee5a608e7c | 253 | CTransposeType transpose() |
| alexpirciu | 1:ceee5a608e7c | 254 | { |
| alexpirciu | 1:ceee5a608e7c | 255 | CTransposeType l_trsp; |
| alexpirciu | 1:ceee5a608e7c | 256 | |
| alexpirciu | 1:ceee5a608e7c | 257 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 258 | { |
| alexpirciu | 1:ceee5a608e7c | 259 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 260 | { |
| alexpirciu | 1:ceee5a608e7c | 261 | l_trsp.m_data[l_col][l_row] = this->m_data[l_row][l_col]; |
| alexpirciu | 1:ceee5a608e7c | 262 | } |
| alexpirciu | 1:ceee5a608e7c | 263 | } |
| alexpirciu | 1:ceee5a608e7c | 264 | |
| alexpirciu | 1:ceee5a608e7c | 265 | return l_trsp; |
| alexpirciu | 1:ceee5a608e7c | 266 | } |
| alexpirciu | 1:ceee5a608e7c | 267 | |
| alexpirciu | 1:ceee5a608e7c | 268 | template <uint32_t P> |
| alexpirciu | 1:ceee5a608e7c | 269 | CRightMultiplicationResultType<P> solve(const CRightMultipliableType<P>& f_B) |
| alexpirciu | 1:ceee5a608e7c | 270 | { |
| alexpirciu | 1:ceee5a608e7c | 271 | return CLUDecomposition(*this).solve(); |
| alexpirciu | 1:ceee5a608e7c | 272 | } |
| alexpirciu | 1:ceee5a608e7c | 273 | |
| alexpirciu | 1:ceee5a608e7c | 274 | static CThisType zeros() |
| alexpirciu | 1:ceee5a608e7c | 275 | { |
| alexpirciu | 1:ceee5a608e7c | 276 | CThisType l_res; |
| alexpirciu | 1:ceee5a608e7c | 277 | |
| alexpirciu | 1:ceee5a608e7c | 278 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 279 | { |
| alexpirciu | 1:ceee5a608e7c | 280 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 281 | { |
| alexpirciu | 1:ceee5a608e7c | 282 | l_res[l_row][l_col] = 0; |
| alexpirciu | 1:ceee5a608e7c | 283 | } |
| alexpirciu | 1:ceee5a608e7c | 284 | } |
| alexpirciu | 1:ceee5a608e7c | 285 | |
| alexpirciu | 1:ceee5a608e7c | 286 | return l_res; |
| alexpirciu | 1:ceee5a608e7c | 287 | } |
| alexpirciu | 1:ceee5a608e7c | 288 | |
| alexpirciu | 1:ceee5a608e7c | 289 | static CThisType ones() |
| alexpirciu | 1:ceee5a608e7c | 290 | { |
| alexpirciu | 1:ceee5a608e7c | 291 | CThisType l_res; |
| alexpirciu | 1:ceee5a608e7c | 292 | |
| alexpirciu | 1:ceee5a608e7c | 293 | for (uint32_t l_row = 0; l_row < M; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 294 | { |
| alexpirciu | 1:ceee5a608e7c | 295 | for (uint32_t l_col = 0; l_col < N; ++l_col) |
| alexpirciu | 1:ceee5a608e7c | 296 | { |
| alexpirciu | 1:ceee5a608e7c | 297 | l_res[l_row][l_col] = 1; |
| alexpirciu | 1:ceee5a608e7c | 298 | } |
| alexpirciu | 1:ceee5a608e7c | 299 | } |
| alexpirciu | 1:ceee5a608e7c | 300 | |
| alexpirciu | 1:ceee5a608e7c | 301 | return l_res; |
| alexpirciu | 1:ceee5a608e7c | 302 | } |
| alexpirciu | 1:ceee5a608e7c | 303 | |
| alexpirciu | 1:ceee5a608e7c | 304 | static CThisType eye() |
| alexpirciu | 1:ceee5a608e7c | 305 | { |
| alexpirciu | 1:ceee5a608e7c | 306 | CThisType l_res(CThisType::zeros()); |
| alexpirciu | 1:ceee5a608e7c | 307 | uint32_t l_minDim = std::min(M,N); |
| alexpirciu | 1:ceee5a608e7c | 308 | |
| alexpirciu | 1:ceee5a608e7c | 309 | for (uint32_t l_row = 0; l_row < l_minDim; ++l_row) |
| alexpirciu | 1:ceee5a608e7c | 310 | { |
| alexpirciu | 1:ceee5a608e7c | 311 | l_res[l_row][l_row] = 1; |
| alexpirciu | 1:ceee5a608e7c | 312 | } |
| alexpirciu | 1:ceee5a608e7c | 313 | |
| alexpirciu | 1:ceee5a608e7c | 314 | return l_res; |
| alexpirciu | 1:ceee5a608e7c | 315 | } |
| alexpirciu | 1:ceee5a608e7c | 316 | |
| alexpirciu | 1:ceee5a608e7c | 317 | protected: |
| alexpirciu | 1:ceee5a608e7c | 318 | std::array<std::array<T,N>,M> m_data; |
| alexpirciu | 1:ceee5a608e7c | 319 | }; |
| alexpirciu | 1:ceee5a608e7c | 320 | |
| alexpirciu | 1:ceee5a608e7c | 321 | template<class T, uint32_t N> |
| alexpirciu | 1:ceee5a608e7c | 322 | class CLUDecomposition |
| alexpirciu | 1:ceee5a608e7c | 323 | { |
| alexpirciu | 1:ceee5a608e7c | 324 | public: |
| alexpirciu | 1:ceee5a608e7c | 325 | using CThisType = CLUDecomposition<T,N>; |
| alexpirciu | 1:ceee5a608e7c | 326 | using COriginalType = CMatrix<T,N,N>; |
| alexpirciu | 1:ceee5a608e7c | 327 | using CDataType =T; |
| alexpirciu | 1:ceee5a608e7c | 328 | |
| alexpirciu | 1:ceee5a608e7c | 329 | template <uint32_t P> |
| alexpirciu | 1:ceee5a608e7c | 330 | using CLeftMultipliableType = CMatrix<T,P,N>; |
| alexpirciu | 1:ceee5a608e7c | 331 | |
| alexpirciu | 1:ceee5a608e7c | 332 | template <uint32_t P> |
| alexpirciu | 1:ceee5a608e7c | 333 | using CRightMultipliableType = CMatrix<T,N,P>; |
| alexpirciu | 1:ceee5a608e7c | 334 | |
| alexpirciu | 1:ceee5a608e7c | 335 | template <uint32_t P> |
| alexpirciu | 1:ceee5a608e7c | 336 | using CLeftMultiplicationResultType = CMatrix<T,P,N>; |
| alexpirciu | 1:ceee5a608e7c | 337 | |
| alexpirciu | 1:ceee5a608e7c | 338 | template <uint32_t P> |
| alexpirciu | 1:ceee5a608e7c | 339 | using CRightMultiplicationResultType = CMatrix<T,N,P>; |
| alexpirciu | 1:ceee5a608e7c | 340 | |
| alexpirciu | 1:ceee5a608e7c | 341 | CLUDecomposition(const CThisType& f_decomposition) : m_L(f_decomposition.m_L), m_U(f_decomposition.m_U), m_P(f_decomposition.m_P) {} |
| alexpirciu | 1:ceee5a608e7c | 342 | CLUDecomposition(const CThisType&& f_decomposition) : m_L(f_decomposition.m_L), m_U(f_decomposition.m_U), m_P(f_decomposition.m_P) {} |
| alexpirciu | 1:ceee5a608e7c | 343 | CLUDecomposition(const COriginalType& f_matrix) : m_L(), m_U(), m_P() {decompose(f_matrix);} |
| alexpirciu | 1:ceee5a608e7c | 344 | CLUDecomposition(const COriginalType&& f_matrix) : m_L(), m_U(), m_P() {decompose(f_matrix);} |
| alexpirciu | 1:ceee5a608e7c | 345 | |
| alexpirciu | 1:ceee5a608e7c | 346 | operator COriginalType() |
| alexpirciu | 1:ceee5a608e7c | 347 | { |
| alexpirciu | 1:ceee5a608e7c | 348 | return m_L * m_U; |
| alexpirciu | 1:ceee5a608e7c | 349 | } |
| alexpirciu | 1:ceee5a608e7c | 350 | |
| alexpirciu | 1:ceee5a608e7c | 351 | COriginalType inv() |
| alexpirciu | 1:ceee5a608e7c | 352 | { |
| alexpirciu | 1:ceee5a608e7c | 353 | return triUInv() * triLInv(); |
| alexpirciu | 1:ceee5a608e7c | 354 | } |
| alexpirciu | 1:ceee5a608e7c | 355 | |
| alexpirciu | 1:ceee5a608e7c | 356 | COriginalType triLInv() |
| alexpirciu | 1:ceee5a608e7c | 357 | { |
| alexpirciu | 1:ceee5a608e7c | 358 | COriginalType l_invL(m_L); |
| alexpirciu | 1:ceee5a608e7c | 359 | |
| alexpirciu | 1:ceee5a608e7c | 360 | for (uint32_t l_jdx = 0; l_jdx < N; ++l_jdx) |
| alexpirciu | 1:ceee5a608e7c | 361 | { |
| alexpirciu | 1:ceee5a608e7c | 362 | l_invL[l_jdx][l_jdx] = 1/l_invL[l_jdx][l_jdx]; |
| alexpirciu | 1:ceee5a608e7c | 363 | |
| alexpirciu | 1:ceee5a608e7c | 364 | for (uint32_t l_idx = l_jdx+1; l_idx < N; ++l_idx) |
| alexpirciu | 1:ceee5a608e7c | 365 | { |
| alexpirciu | 1:ceee5a608e7c | 366 | l_invL[l_idx][l_jdx] = 0; |
| alexpirciu | 1:ceee5a608e7c | 367 | |
| alexpirciu | 1:ceee5a608e7c | 368 | for (uint32_t l_kdx = l_jdx; l_kdx < l_idx; ++l_kdx) |
| alexpirciu | 1:ceee5a608e7c | 369 | { |
| alexpirciu | 1:ceee5a608e7c | 370 | l_invL[l_idx][l_jdx] -= m_L[l_idx][l_kdx]*l_invL[l_kdx][l_jdx]; |
| alexpirciu | 1:ceee5a608e7c | 371 | } |
| alexpirciu | 1:ceee5a608e7c | 372 | |
| alexpirciu | 1:ceee5a608e7c | 373 | l_invL[l_idx][l_jdx] /= l_invL[l_idx][l_idx]; |
| alexpirciu | 1:ceee5a608e7c | 374 | } |
| alexpirciu | 1:ceee5a608e7c | 375 | } |
| alexpirciu | 1:ceee5a608e7c | 376 | |
| alexpirciu | 1:ceee5a608e7c | 377 | return l_invL; |
| alexpirciu | 1:ceee5a608e7c | 378 | } |
| alexpirciu | 1:ceee5a608e7c | 379 | |
| alexpirciu | 1:ceee5a608e7c | 380 | COriginalType triUInv() |
| alexpirciu | 1:ceee5a608e7c | 381 | { |
| alexpirciu | 1:ceee5a608e7c | 382 | COriginalType l_invU(m_U); |
| alexpirciu | 1:ceee5a608e7c | 383 | |
| alexpirciu | 1:ceee5a608e7c | 384 | for (int32_t l_jdx = (N-1); l_jdx >= 0; --l_jdx) |
| alexpirciu | 1:ceee5a608e7c | 385 | { |
| alexpirciu | 1:ceee5a608e7c | 386 | l_invU[l_jdx][l_jdx] = 1/l_invU[l_jdx][l_jdx]; |
| alexpirciu | 1:ceee5a608e7c | 387 | |
| alexpirciu | 1:ceee5a608e7c | 388 | for (int32_t l_idx = (l_jdx-1); l_idx >= 0; --l_idx) |
| alexpirciu | 1:ceee5a608e7c | 389 | { |
| alexpirciu | 1:ceee5a608e7c | 390 | l_invU[l_idx][l_jdx] = 0; |
| alexpirciu | 1:ceee5a608e7c | 391 | |
| alexpirciu | 1:ceee5a608e7c | 392 | for (int32_t l_kdx = (l_idx+1); l_kdx < (l_jdx+1); ++l_kdx) |
| alexpirciu | 1:ceee5a608e7c | 393 | { |
| alexpirciu | 1:ceee5a608e7c | 394 | l_invU[l_idx][l_jdx] -= m_U[l_idx][l_kdx]*l_invU[l_kdx][l_jdx]; |
| alexpirciu | 1:ceee5a608e7c | 395 | } |
| alexpirciu | 1:ceee5a608e7c | 396 | |
| alexpirciu | 1:ceee5a608e7c | 397 | l_invU[l_idx][l_jdx] /= l_invU[l_idx][l_idx]; |
| alexpirciu | 1:ceee5a608e7c | 398 | } |
| alexpirciu | 1:ceee5a608e7c | 399 | } |
| alexpirciu | 1:ceee5a608e7c | 400 | |
| alexpirciu | 1:ceee5a608e7c | 401 | return l_invU; |
| alexpirciu | 1:ceee5a608e7c | 402 | } |
| alexpirciu | 1:ceee5a608e7c | 403 | |
| alexpirciu | 1:ceee5a608e7c | 404 | template <uint32_t P> |
| alexpirciu | 1:ceee5a608e7c | 405 | CRightMultiplicationResultType<P> solve(const CRightMultipliableType<P>& f_B) |
| alexpirciu | 1:ceee5a608e7c | 406 | { |
| alexpirciu | 1:ceee5a608e7c | 407 | CRightMultipliableType<P> l_B; |
| alexpirciu | 1:ceee5a608e7c | 408 | |
| alexpirciu | 1:ceee5a608e7c | 409 | for (uint32_t l_idx = 0; l_idx < N; ++l_idx) |
| alexpirciu | 1:ceee5a608e7c | 410 | { |
| alexpirciu | 1:ceee5a608e7c | 411 | for (uint32_t l_jdx = 0; l_jdx < P; ++l_jdx) |
| alexpirciu | 1:ceee5a608e7c | 412 | { |
| alexpirciu | 1:ceee5a608e7c | 413 | l_B[l_idx][l_jdx] = f_B[m_P[l_idx]][l_jdx]; |
| alexpirciu | 1:ceee5a608e7c | 414 | } |
| alexpirciu | 1:ceee5a608e7c | 415 | } |
| alexpirciu | 1:ceee5a608e7c | 416 | |
| alexpirciu | 1:ceee5a608e7c | 417 | return sTriL(sTriU(l_B)); |
| alexpirciu | 1:ceee5a608e7c | 418 | } |
| alexpirciu | 1:ceee5a608e7c | 419 | |
| alexpirciu | 1:ceee5a608e7c | 420 | template <uint32_t P> |
| alexpirciu | 1:ceee5a608e7c | 421 | CRightMultiplicationResultType<P> sTriU(const CRightMultipliableType<P>& f_B) |
| alexpirciu | 1:ceee5a608e7c | 422 | { |
| alexpirciu | 1:ceee5a608e7c | 423 | return triUInv()*f_B; |
| alexpirciu | 1:ceee5a608e7c | 424 | } |
| alexpirciu | 1:ceee5a608e7c | 425 | |
| alexpirciu | 1:ceee5a608e7c | 426 | template <uint32_t P> |
| alexpirciu | 1:ceee5a608e7c | 427 | CRightMultiplicationResultType<P> sTriL(const CRightMultipliableType<P>& f_B) |
| alexpirciu | 1:ceee5a608e7c | 428 | { |
| alexpirciu | 1:ceee5a608e7c | 429 | // return triLInv()*f_B; |
| alexpirciu | 1:ceee5a608e7c | 430 | } |
| alexpirciu | 1:ceee5a608e7c | 431 | |
| alexpirciu | 1:ceee5a608e7c | 432 | |
| alexpirciu | 1:ceee5a608e7c | 433 | // private: |
| alexpirciu | 1:ceee5a608e7c | 434 | void decompose(const CMatrix<T,N,N>& f_matrix) |
| alexpirciu | 1:ceee5a608e7c | 435 | { |
| alexpirciu | 1:ceee5a608e7c | 436 | m_U = f_matrix; |
| alexpirciu | 1:ceee5a608e7c | 437 | for (uint32_t l_kdx = 0; l_kdx < (N-1); ++l_kdx) |
| alexpirciu | 1:ceee5a608e7c | 438 | { |
| alexpirciu | 1:ceee5a608e7c | 439 | m_L[l_kdx][l_kdx] = 1; |
| alexpirciu | 1:ceee5a608e7c | 440 | m_P[l_kdx] = l_kdx; |
| alexpirciu | 1:ceee5a608e7c | 441 | |
| alexpirciu | 1:ceee5a608e7c | 442 | for(uint32_t l_idx = l_kdx+1; l_idx < N; ++l_idx) |
| alexpirciu | 1:ceee5a608e7c | 443 | { |
| alexpirciu | 1:ceee5a608e7c | 444 | m_L[l_idx][l_kdx] = m_U[l_idx][l_kdx]/m_U[l_kdx][l_kdx]; |
| alexpirciu | 1:ceee5a608e7c | 445 | m_U[l_idx][l_kdx] = 0; |
| alexpirciu | 1:ceee5a608e7c | 446 | for(uint32_t l_jdx = l_kdx+1; l_jdx < N; ++l_jdx) |
| alexpirciu | 1:ceee5a608e7c | 447 | { |
| alexpirciu | 1:ceee5a608e7c | 448 | m_U[l_idx][l_jdx] = m_U[l_idx][l_jdx] - m_L[l_idx][l_kdx]*m_U[l_kdx][l_jdx]; |
| alexpirciu | 1:ceee5a608e7c | 449 | } |
| alexpirciu | 1:ceee5a608e7c | 450 | } |
| alexpirciu | 1:ceee5a608e7c | 451 | } |
| alexpirciu | 1:ceee5a608e7c | 452 | m_L[N-1][N-1] = 1; |
| alexpirciu | 1:ceee5a608e7c | 453 | m_P[N-1] = N-1; |
| alexpirciu | 1:ceee5a608e7c | 454 | } |
| alexpirciu | 1:ceee5a608e7c | 455 | CMatrix<T,N,N> m_L; |
| alexpirciu | 1:ceee5a608e7c | 456 | CMatrix<T,N,N> m_U; |
| alexpirciu | 1:ceee5a608e7c | 457 | std::array<uint32_t,N> m_P; |
| alexpirciu | 1:ceee5a608e7c | 458 | }; |
| alexpirciu | 1:ceee5a608e7c | 459 | |
| alexpirciu | 1:ceee5a608e7c | 460 | template <class T, uint32_t N> |
| alexpirciu | 1:ceee5a608e7c | 461 | using CColVector = CMatrix<T,N,1>; |
| alexpirciu | 1:ceee5a608e7c | 462 | |
| alexpirciu | 1:ceee5a608e7c | 463 | template <class T, uint32_t N> |
| alexpirciu | 1:ceee5a608e7c | 464 | using CVector = CColVector<T,N>; |
| alexpirciu | 1:ceee5a608e7c | 465 | |
| alexpirciu | 1:ceee5a608e7c | 466 | template <class T, uint32_t N> |
| alexpirciu | 1:ceee5a608e7c | 467 | using CRowVector = CMatrix<T,1,N>; |
| alexpirciu | 1:ceee5a608e7c | 468 | }; |
| alexpirciu | 1:ceee5a608e7c | 469 | |
| alexpirciu | 1:ceee5a608e7c | 470 | #include "linalg.inl" |
| alexpirciu | 1:ceee5a608e7c | 471 | |
| alexpirciu | 1:ceee5a608e7c | 472 | #endif // LINALG_H |
| alexpirciu | 1:ceee5a608e7c | 473 |