Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TDecompChol.h
Go to the documentation of this file.
1// @(#)root/matrix:$Id$
2// Authors: Fons Rademakers, Eddy Offermann Dec 2003
3
4/*************************************************************************
5 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12#ifndef ROOT_TDecompChol
13#define ROOT_TDecompChol
14
15///////////////////////////////////////////////////////////////////////////
16// //
17// Cholesky Decomposition class //
18// //
19///////////////////////////////////////////////////////////////////////////
20
21#include "TDecompBase.h"
22#include "TMatrixDSym.h"
23
25{
26protected :
27
28 TMatrixD fU; // decomposed matrix fU so that a = fU^T fU
29
30 const TMatrixDBase &GetDecompMatrix() const override { return fU; }
31
32public :
33
34 TDecompChol() : fU() {}
35 explicit TDecompChol(Int_t nrows);
36 TDecompChol(Int_t row_lwb,Int_t row_upb);
37 TDecompChol(const TMatrixDSym &a,Double_t tol = 0.0);
38 TDecompChol(const TMatrixD &a,Double_t tol = 0.0);
39 TDecompChol(const TDecompChol &another);
40 ~TDecompChol() override {}
41
42 const TMatrixDSym GetMatrix ();
43 Int_t GetNrows () const override { return fU.GetNrows(); }
44 Int_t GetNcols () const override { return fU.GetNcols(); }
45 const TMatrixD &GetU () const { return fU; }
46
47 virtual void SetMatrix (const TMatrixDSym &a);
48
49 Bool_t Decompose () override;
50 Bool_t Solve ( TVectorD &b) override;
51 TVectorD Solve (const TVectorD& b,Bool_t &ok) override { TVectorD x = b; ok = Solve(x); return x; }
52 Bool_t Solve ( TMatrixDColumn &b) override;
53 Bool_t TransSolve ( TVectorD &b) override { return Solve(b); }
54 TVectorD TransSolve (const TVectorD& b,Bool_t &ok) override { TVectorD x = b; ok = Solve(x); return x; }
55 Bool_t TransSolve ( TMatrixDColumn &b) override { return Solve(b); }
56 void Det (Double_t &d1,Double_t &d2) override;
57
59 TMatrixDSym Invert (Bool_t &status);
60 TMatrixDSym Invert () { Bool_t status; return Invert(status); }
61
62 void Print(Option_t *opt ="") const override; // *MENU*
63
64 TDecompChol &operator= (const TDecompChol &source);
65
66 ClassDefOverride(TDecompChol,2) // Matrix Decompositition Cholesky
67};
68
69TVectorD NormalEqn(const TMatrixD &A,const TVectorD &b);
70TVectorD NormalEqn(const TMatrixD &A,const TVectorD &b,const TVectorD &std);
71TMatrixD NormalEqn(const TMatrixD &A,const TMatrixD &b);
72TMatrixD NormalEqn(const TMatrixD &A,const TMatrixD &B,const TVectorD &std);
73
74#endif
#define b(i)
Definition RSha256.hxx:100
#define a(i)
Definition RSha256.hxx:99
bool Bool_t
Definition RtypesCore.h:63
int Int_t
Definition RtypesCore.h:45
double Double_t
Definition RtypesCore.h:59
const char Option_t
Definition RtypesCore.h:66
#define ClassDefOverride(name, id)
Definition Rtypes.h:341
TVectorD NormalEqn(const TMatrixD &A, const TVectorD &b)
Solve min {(A .
Decomposition Base class.
Definition TDecompBase.h:34
Cholesky Decomposition class.
Definition TDecompChol.h:25
Int_t GetNcols() const override
Definition TDecompChol.h:44
void Print(Option_t *opt="") const override
Print class members .
Bool_t Decompose() override
Matrix A is decomposed in component U so that A = U^T * U If the decomposition succeeds,...
TMatrixD fU
Definition TDecompChol.h:28
Int_t GetNrows() const override
Definition TDecompChol.h:43
Bool_t TransSolve(TMatrixDColumn &b) override
Definition TDecompChol.h:55
virtual void SetMatrix(const TMatrixDSym &a)
Set the matrix to be decomposed, decomposition status is reset.
Bool_t Solve(TVectorD &b) override
Solve equations Ax=b assuming A has been factored by Cholesky.
const TMatrixD & GetU() const
Definition TDecompChol.h:45
TDecompChol & operator=(const TDecompChol &source)
Assignment operator.
TVectorD Solve(const TVectorD &b, Bool_t &ok) override
Definition TDecompChol.h:51
TVectorD TransSolve(const TVectorD &b, Bool_t &ok) override
Definition TDecompChol.h:54
~TDecompChol() override
Definition TDecompChol.h:40
const TMatrixDSym GetMatrix()
Reconstruct the original matrix using the decomposition parts.
const TMatrixDBase & GetDecompMatrix() const override
Definition TDecompChol.h:30
void Det(Double_t &d1, Double_t &d2) override
Matrix determinant det = d1*TMath::Power(2.,d2) is square of diagProd of cholesky factor.
TMatrixDSym Invert()
Definition TDecompChol.h:60
Bool_t TransSolve(TVectorD &b) override
Definition TDecompChol.h:53
Int_t GetNrows() const
Int_t GetNcols() const
Double_t x[n]
Definition legend1.C:17
void inv(rsa_NUMBER *, rsa_NUMBER *, rsa_NUMBER *)
Definition rsaaux.cxx:949