From 406f2abbc835a9dc69796a8067ebfc172670d381 Mon Sep 17 00:00:00 2001 From: nojhan Date: Wed, 14 Dec 2011 15:53:57 +0100 Subject: [PATCH] beautifying code for Cholesky classes --- edo/src/utils/edoCholesky.h | 17 +++++++++++------ 1 file changed, 11 insertions(+), 6 deletions(-) diff --git a/edo/src/utils/edoCholesky.h b/edo/src/utils/edoCholesky.h index ea0c2bca1..be17f5fcd 100644 --- a/edo/src/utils/edoCholesky.h +++ b/edo/src/utils/edoCholesky.h @@ -129,6 +129,7 @@ protected: template< typename T > class CholeskyLLT : public CholeskyBase { +public: virtual void operator()( const CovarMat& V ) { unsigned int N = assert_properties( V ); @@ -164,6 +165,8 @@ class CholeskyLLT : public CholeskyBase } // for i in [1,N[ } + + /** The step of the standard LLT algorithm where round off errors may appear */ inline virtual L_i_i( const CovarMat& V, const unsigned int& i, const double& sum ) const { // round-off errors may appear here @@ -173,7 +176,6 @@ class CholeskyLLT : public CholeskyBase }; - /** This standard algorithm makes use of square root but do not fail * if the covariance matrix is very ill-conditioned. * Here, we propagate the error by using the absolute value before @@ -185,6 +187,7 @@ class CholeskyLLT : public CholeskyBase template< typename T > class CholeskyLLTabs : public CholeskyLLT { +public: inline virtual L_i_i( const CovarMat& V, const unsigned int& i, const double& sum ) const { /***** ugly hack *****/ @@ -203,6 +206,7 @@ class CholeskyLLTabs : public CholeskyLLT template< typename T > class CholeskyLLTzero : public CholeskyLLT { +public: inline virtual L_i_i( const CovarMat& V, const unsigned int& i, const double& sum ) const { T Lii; @@ -226,6 +230,7 @@ class CholeskyLLTzero : public CholeskyLLT template< typename T > class CholeskyLDLT : public CholeskyBase { +public: virtual void operator()( const CovarMat& V ) { // use "int" everywhere, because of the "j-1" operation @@ -256,24 +261,24 @@ class CholeskyLDLT : public CholeskyBase } // for i in rows } // for j in columns - _L = compute_L( L, D ); + _L = root( L, D ); } - inline FactorMat compute_L( FactorMat& L, FactorMat& D ) + inline FactorMat root( FactorMat& L, FactorMat& D ) { // now compute the final symetric matrix: _L = L D^1/2 // remember that V = ( L D^1/2) ( L D^1/2)^T // fortunately, the square root of a diagonal matrix is the square // root of all its elements - FactorMat D12 = D; + FactorMat sqrt_D = D; for(int i=0; i