27 #ifndef CDPL_MATH_LUDECOMPOSITION_HPP
28 #define CDPL_MATH_LUDECOMPOSITION_HPP
47 typedef typename E::SizeType SizeType;
48 typedef typename E::ValueType ValueType;
53 SizeType size1 = e().getSize1();
54 SizeType size2 = e().getSize2();
55 SizeType size = std::min(size1, size2);
60 SizeType singular = 0;
62 for (SizeType i = 0; i < size; i++) {
66 if (e()(i, i) != ValueType(0)) {
67 ValueType m_inv = ValueType(1) / e()(i, i);
68 range(col_i, ColumnRangeType(i + 1, size1)) *= m_inv;
70 }
else if (singular == 0)
73 range(e, MatrixRangeType(i + 1, size1), MatrixRangeType(i + 1, size2))
74 .minusAssign(
outerProd(
range(col_i, ColumnRangeType(i + 1, size1)),
75 range(row_i, RowRangeType(i + 1, size2))));
81 template <
typename E,
typename PV,
typename T>
85 typedef typename E::SizeType SizeType;
86 typedef typename E::ValueType ValueType;
91 SizeType size1 = e().getSize1();
92 SizeType size2 = e().getSize2();
93 SizeType size = std::min(size1, size2);
100 SizeType singular = 0;
102 for (SizeType i = 0; i < size; i++) {
105 SizeType norm_inf_idx = i +
normInfIndex(
range(col_i, ColumnRangeType(i, size1)));
107 if (e()(norm_inf_idx, i) != ValueType(0)) {
108 pv[i] = norm_inf_idx;
110 if (norm_inf_idx != i) {
111 row(e, norm_inf_idx).swap(row_i);
115 ValueType m_inv = ValueType(1) / e()(i, i);
116 range(col_i, ColumnRangeType(i + 1, size1)) *= m_inv;
118 }
else if (singular == 0)
121 range(e, MatrixRangeType(i + 1, size1), MatrixRangeType(i + 1, size2))
122 .minusAssign(
outerProd(
range(col_i, ColumnRangeType(i + 1, size1)),
123 range(row_i, RowRangeType(i + 1, size2))));
129 template <
typename E,
typename PV>
133 typedef typename E::SizeType SizeType;
135 for (SizeType i = 0, size = e().getSize(); i < size; i++) {
136 if (i != SizeType(pv[i]))
137 std::swap(e()(i), e()(pv[i]));
141 template <
typename E,
typename PV>
145 typedef typename E::SizeType SizeType;
147 for (SizeType i = 0, size = e().getSize1(); i < size; i++) {
148 if (i != SizeType(pv[i])) {
150 row(e, i).swap(other_row);
155 template <
typename E1,
typename E2>
165 template <
typename E1,
typename E2,
typename PV>
174 template <
typename E1,
typename E2>
184 template <
typename E1,
typename E2,
typename PV>
195 #endif // CDPL_MATH_LUDECOMPOSITION_HPP