NumCpp  2.7.0
A Templatized Header Only C++ Implementation of the Python NumPy Library
cholesky.hpp
Go to the documentation of this file.
1 #pragma once
33 
36 #include "NumCpp/Core/Types.hpp"
37 #include "NumCpp/NdArray.hpp"
38 
39 namespace nc
40 {
41  namespace linalg
42  {
43  //============================================================================
44  // 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.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.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 linalg
101 } // namespace nc
#define THROW_RUNTIME_ERROR(msg)
Definition: Error.hpp:37
#define STATIC_ASSERT_ARITHMETIC(dtype)
Definition: StaticAsserts.hpp:37
Shape shape() const noexcept
Definition: NdArrayCore.hpp:4283
uint32 rows
Definition: Core/Shape.hpp:44
bool issquare() const noexcept
Definition: Core/Shape.hpp:123
uint32 cols
Definition: Core/Shape.hpp:45
constexpr auto j
Definition: Constants.hpp:45
NdArray< double > cholesky(const NdArray< dtype > &inMatrix)
Definition: cholesky.hpp:54
Definition: Coordinate.hpp:45
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