NumCpp  2.12.1
A Templatized Header Only C++ Implementation of the Python NumPy Library
lu_decomposition.hpp
Go to the documentation of this file.
1
33#pragma once
34
35#include <cmath>
36#include <utility>
37
40#include "NumCpp/Core/Types.hpp"
42#include "NumCpp/NdArray.hpp"
44
45namespace nc::linalg
46{
47 //============================================================================
48 // Method Description:
55 template<typename dtype>
56 std::pair<NdArray<double>, NdArray<double>> lu_decomposition(const NdArray<dtype>& inMatrix)
57 {
59
60 const auto shape = inMatrix.shape();
61 if (!shape.issquare())
62 {
63 THROW_RUNTIME_ERROR("Input matrix should be square.");
64 }
65
66 NdArray<double> lMatrix = zeros_like<double>(inMatrix);
67 NdArray<double> uMatrix = inMatrix.template astype<double>();
68
69 for (uint32 col = 0; col < shape.cols; ++col)
70 {
71 lMatrix(col, col) = 1;
72
73 for (uint32 row = col + 1; row < shape.rows; ++row)
74 {
75 const double& divisor = uMatrix(col, col);
76 if (utils::essentiallyEqual(divisor, double{ 0. }))
77 {
78 THROW_RUNTIME_ERROR("Division by 0.");
79 }
80
81 lMatrix(row, col) = uMatrix(row, col) / divisor;
82
83 for (uint32 col2 = col; col2 < shape.cols; ++col2)
84 {
85 uMatrix(row, col2) -= lMatrix(row, col) * uMatrix(col, col2);
86 }
87 }
88 }
89
90 return std::make_pair(lMatrix, uMatrix);
91 }
92} // namespace nc::linalg
#define THROW_RUNTIME_ERROR(msg)
Definition: Error.hpp:40
#define STATIC_ASSERT_ARITHMETIC(dtype)
Definition: StaticAsserts.hpp:39
const Shape & shape() const noexcept
Definition: NdArrayCore.hpp:4511
uint32 rows
Definition: Core/Shape.hpp:44
bool issquare() const noexcept
Definition: Core/Shape.hpp:125
uint32 cols
Definition: Core/Shape.hpp:45
Definition: cholesky.hpp:41
std::pair< NdArray< double >, NdArray< double > > lu_decomposition(const NdArray< dtype > &inMatrix)
Definition: lu_decomposition.hpp:56
bool essentiallyEqual(dtype inValue1, dtype inValue2) noexcept
Definition: essentiallyEqual.hpp:49
Shape shape(const NdArray< dtype > &inArray) noexcept
Definition: Functions/Shape.hpp:42
std::uint32_t uint32
Definition: Types.hpp:40