NumCpp  2.11.0
A Templatized Header Only C++ Implementation of the Python NumPy Library
cholesky.hpp
Go to the documentation of this file.
1 #pragma once
34 
37 #include "NumCpp/Core/Types.hpp"
38 #include "NumCpp/NdArray.hpp"
39 
40 namespace nc::linalg
41 {
42  //============================================================================
43  // Method Description:
53  template<typename dtype>
55  {
57 
58  const auto shape = inMatrix.shape();
59  if (!shape.issquare())
60  {
61  THROW_RUNTIME_ERROR("Input matrix should be square.");
62  }
63 
64  auto lMatrix = inMatrix.template astype<double>();
65 
66  for (uint32 row = 0; row < shape.rows; ++row)
67  {
68  for (uint32 col = row + 1; col < shape.cols; ++col)
69  {
70  lMatrix(row, col) = 0.;
71  }
72  }
73 
74  for (uint32 k = 0; k < shape.cols; ++k)
75  {
76  const double& a_kk = lMatrix(k, k);
77 
78  if (a_kk > 0.)
79  {
80  lMatrix(k, k) = std::sqrt(a_kk);
81 
82  for (uint32 i = k + 1; i < shape.rows; ++i)
83  {
84  lMatrix(i, k) /= lMatrix(k, k);
85 
86  for (uint32 j = k + 1; j <= i; ++j)
87  {
88  lMatrix(i, j) -= lMatrix(i, k) * lMatrix(j, k);
89  }
90  }
91  }
92  else
93  {
94  THROW_RUNTIME_ERROR("Matrix is not positive definite.");
95  }
96  }
97 
98  return lMatrix;
99  }
100 } // 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:4464
uint32 rows
Definition: Core/Shape.hpp:44
bool issquare() const noexcept
Definition: Core/Shape.hpp:125
uint32 cols
Definition: Core/Shape.hpp:45
constexpr auto j
Definition: Core/Constants.hpp:42
Definition: cholesky.hpp:41
NdArray< double > cholesky(const NdArray< dtype > &inMatrix)
Definition: cholesky.hpp:54
auto sqrt(dtype inValue) noexcept
Definition: sqrt.hpp:48
Shape shape(const NdArray< dtype > &inArray) noexcept
Definition: Functions/Shape.hpp:42
std::uint32_t uint32
Definition: Types.hpp:40