ROOT logo
ROOT » MATH » MATRIX » TDecompBase

class TDecompBase: public TObject


Decomposition Base class

This class forms the base for all the decompositions methods in the
linear algebra package .
It or its derived classes have installed the methods to solve
equations,invert matrices and calculate determinants while monitoring
the accuracy.

Each derived class has always the following methods available:

Condition() :
In an iterative scheme the condition number for matrix inversion is
calculated . This number is of interest for estimating the accuracy
of x in the equation Ax=b
For example:
A is a (10x10) Hilbert matrix which looks deceivingly innocent
and simple, A(i,j) = 1/(i+j+1)
b(i) = Sum_j A(i,j), so a sum of a row in A

the solution is x(i) = 1. i=0,.,9

However,
TMatrixD m....; TVectorD b.....
TDecompLU lu(m); lu.SetTol(1.0e-12); lu.Solve(b); b.Print()
gives,

{1.000,1.000,1.000,1.000,0.998,1.000,0.993,1.001,0.996,1.000}

Looking at the condition number, this is in line with expected the
accuracy . The condition number is 3.957e+12 . As a simple rule of
thumb, a condition number of 1.0e+n means that you lose up to n
digits of accuracy in a solution . Since doubles are stored with 15
digits, we can expect the accuracy to be as small as 3 digits .

Det(Double_t &d1,Double_t &d2)
The determinant is d1*TMath::Power(2.,d2)
Expressing the determinant this way makes under/over-flow very
unlikely .

Decompose()
Here the actually decomposition is performed . One can change the
matrix A after the decomposition constructor has been called
without effecting the decomposition result

Solve(TVectorD &b)
Solve A x = b . x is supplied through the argument and replaced with
the solution .

TransSolve(TVectorD &b)
Solve A^T x = b . x is supplied through the argument and replaced
with the solution .

MultiSolve(TMatrixD    &B)
Solve A X = B . where X and are now matrices . X is supplied through
the argument and replaced with the solution .

Invert(TMatrixD &inv)
This is of course just a call to MultiSolve with as input argument
the unit matrix . Note that for a matrix a(m,n) with m > n  a
pseudo-inverse is calculated .

Tolerances and Scaling

The tolerance parameter (which is a member of this base class) plays
a crucial role in all operations of the decomposition classes . It
gives the user a powerful tool to monitor and steer the operations
Its default value is sqrt(epsilon) where 1+epsilon = 1

If you do not want to be bothered by the following considerations,
like in most other linear algebra packages, just set the tolerance
with SetTol to an arbitrary small number .

The tolerance number is used by each decomposition method to decide
whether the matrix is near singular, except of course SVD which can
handle singular matrices .
For each decomposition this will be checked in a different way; in LU
the matrix is considered singular when, at some point in the
 decomposition, a diagonal element < fTol . Therefore, we had to set in
the example above of the (10x10) Hilbert, which is near singular, the
tolerance on 10e-12 . (The fact that we have to set the tolerance <
sqrt(epsilon) is a clear indication that we are losing precision .)

If the matrix is flagged as being singular, operations with the
decomposition will fail and will return matrices/vectors that are
invalid .

The observant reader will notice that by scaling the complete matrix
 by some small number the decomposition will detect a singular matrix .
 In this case the user will have to reduce the tolerance number by this
 factor . (For CPU time saving we decided not to make this an automatic
procedure) .

Code for this could look as follows:
const Double_t max_abs = Abs(a).Max();
const Double_t scale = TMath::Min(max_abs,1.);
a.SetTol(a.GetTol()*scale);

For usage examples see $ROOTSYS/test/stressLinear.cxx

Function Members (Methods)

 
    This is an abstract class, constructors will not be documented.
    Look at the header to check for available constructors.

public:
virtual~TDecompBase()
voidTObject::AbstractMethod(const char* method) const
virtual voidTObject::AppendPad(Option_t* option = "")
virtual voidTObject::Browse(TBrowser* b)
static TClass*Class()
static TClass*TObject::Class()
virtual const char*TObject::ClassName() const
virtual voidTObject::Clear(Option_t* = "")
virtual TObject*TObject::Clone(const char* newname = "") const
virtual Int_tTObject::Compare(const TObject* obj) const
virtual Double_tCondition()
virtual voidTObject::Copy(TObject& object) const
virtual Bool_tDecompose()
virtual voidTObject::Delete(Option_t* option = "")MENU
virtual voidDet(Double_t& d1, Double_t& d2)
virtual Int_tTObject::DistancetoPrimitive(Int_t px, Int_t py)
virtual voidTObject::Draw(Option_t* option = "")
virtual voidTObject::DrawClass() constMENU
virtual TObject*TObject::DrawClone(Option_t* option = "") constMENU
virtual voidTObject::Dump() constMENU
virtual voidTObject::Error(const char* method, const char* msgfmt) const
virtual voidTObject::Execute(const char* method, const char* params, Int_t* error = 0)
virtual voidTObject::Execute(TMethod* method, TObjArray* params, Int_t* error = 0)
virtual voidTObject::ExecuteEvent(Int_t event, Int_t px, Int_t py)
virtual voidTObject::Fatal(const char* method, const char* msgfmt) const
virtual TObject*TObject::FindObject(const char* name) const
virtual TObject*TObject::FindObject(const TObject* obj) const
Int_tGetColLwb() const
Double_tGetCondition() const
Double_tGetDet1() const
Double_tGetDet2() const
virtual Option_t*TObject::GetDrawOption() const
static Long_tTObject::GetDtorOnly()
virtual const char*TObject::GetIconName() const
virtual const char*TObject::GetName() const
virtual Int_tGetNcols() const
virtual Int_tGetNrows() const
virtual char*TObject::GetObjectInfo(Int_t px, Int_t py) const
static Bool_tTObject::GetObjectStat()
virtual Option_t*TObject::GetOption() const
Int_tGetRowLwb() const
virtual const char*TObject::GetTitle() const
Double_tGetTol() const
virtual UInt_tTObject::GetUniqueID() const
virtual Bool_tTObject::HandleTimer(TTimer* timer)
virtual ULong_tTObject::Hash() const
virtual voidTObject::Info(const char* method, const char* msgfmt) const
virtual Bool_tTObject::InheritsFrom(const char* classname) const
virtual Bool_tTObject::InheritsFrom(const TClass* cl) const
virtual voidTObject::Inspect() constMENU
voidTObject::InvertBit(UInt_t f)
virtual TClass*IsA() const
virtual TClass*TObject::IsA() const
virtual Bool_tTObject::IsEqual(const TObject* obj) const
virtual Bool_tTObject::IsFolder() const
Bool_tTObject::IsOnHeap() const
virtual Bool_tTObject::IsSortable() const
Bool_tTObject::IsZombie() const
virtual voidTObject::ls(Option_t* option = "") const
voidTObject::MayNotUse(const char* method) const
virtual Bool_tMultiSolve(TMatrixD& B)
virtual Bool_tTObject::Notify()
static voidTObject::operator delete(void* ptr)
static voidTObject::operator delete(void* ptr, void* vp)
static voidTObject::operator delete[](void* ptr)
static voidTObject::operator delete[](void* ptr, void* vp)
void*TObject::operator new(size_t sz)
void*TObject::operator new(size_t sz, void* vp)
void*TObject::operator new[](size_t sz)
void*TObject::operator new[](size_t sz, void* vp)
TDecompBase&operator=(const TDecompBase& source)
TObject&TObject::operator=(const TObject& rhs)
virtual voidTObject::Paint(Option_t* option = "")
virtual voidTObject::Pop()
virtual voidPrint(Option_t* opt = "") const
virtual voidTObject::Print(Option_t* option = "") const
virtual Int_tTObject::Read(const char* name)
virtual voidTObject::RecursiveRemove(TObject* obj)
voidTObject::ResetBit(UInt_t f)
virtual voidTObject::SaveAs(const char* filename = "", Option_t* option = "") constMENU
virtual voidTObject::SavePrimitive(basic_ostream<char,char_traits<char> >& out, Option_t* option = "")
voidTObject::SetBit(UInt_t f)
voidTObject::SetBit(UInt_t f, Bool_t set)
virtual voidTObject::SetDrawOption(Option_t* option = "")MENU
static voidTObject::SetDtorOnly(void* obj)
static voidTObject::SetObjectStat(Bool_t stat)
Double_tSetTol(Double_t newTol)
virtual voidTObject::SetUniqueID(UInt_t uid)
virtual voidShowMembers(TMemberInspector& insp, char* parent)
virtual voidTObject::ShowMembers(TMemberInspector& insp, char* parent)
virtual Bool_tSolve(TVectorD& b)
virtual Bool_tSolve(TMatrixDColumn& b)
virtual TVectorDSolve(const TVectorD& b, Bool_t& ok)
virtual voidStreamer(TBuffer& b)
virtual voidTObject::Streamer(TBuffer& b)
voidStreamerNVirtual(TBuffer& b)
voidTObject::StreamerNVirtual(TBuffer& b)
virtual voidTObject::SysError(const char* method, const char* msgfmt) const
Bool_tTObject::TestBit(UInt_t f) const
Int_tTObject::TestBits(UInt_t f) const
virtual Bool_tTransSolve(TVectorD& b)
virtual Bool_tTransSolve(TMatrixDColumn& b)
virtual TVectorDTransSolve(const TVectorD& b, Bool_t& ok)
virtual voidTObject::UseCurrentStyle()
virtual voidTObject::Warning(const char* method, const char* msgfmt) const
virtual Int_tTObject::Write(const char* name = 0, Int_t option = 0, Int_t bufsize = 0)
virtual Int_tTObject::Write(const char* name = 0, Int_t option = 0, Int_t bufsize = 0) const
protected:
static voidDiagProd(const TVectorD& diag, Double_t tol, Double_t& d1, Double_t& d2)
virtual voidTObject::DoError(int level, const char* location, const char* fmt, va_list va) const
virtual const TMatrixDBase&GetDecompMatrix() const
Int_tHager(Double_t& est, Int_t iter = 5)
voidTObject::MakeZombie()
voidResetStatus()

Data Members

private:
enum EMatrixDecompStat { kInit
kPatternSet
kValuesSet
kMatrixSet
kDecomposed
kDetermined
kCondition
kSingular
};
enum { kWorkMax
};
enum TObject::EStatusBits { kCanDelete
kMustCleanup
kObjInCanvas
kIsReferenced
kHasUUID
kCannotPick
kNoContextMenu
kInvalidObject
};
enum TObject::[unnamed] { kIsOnHeap
kNotDeleted
kZombie
kBitMask
kSingleKey
kOverwrite
kWriteDelete
};
protected:
Int_tfColLwbColumn lower bound of decomposed matrix
Double_tfConditionmatrix condition number
Double_tfDet1determinant mantissa
Double_tfDet2determinant exponent for powers of 2
Int_tfRowLwbRow lower bound of decomposed matrix
Double_tfTolsqrt(epsilon); epsilon is smallest number number so that 1+epsilon > 1

Class Charts

Inheritance Inherited Members Includes Libraries
Class Charts

Function documentation

Int_t Hager(Double_t& est, Int_t iter = 5)
void DiagProd(const TVectorD& diag, Double_t tol, Double_t& d1, Double_t& d2)
Double_t Condition()
 Matrix condition number
Bool_t MultiSolve(TMatrixD& B)
 Solve set of equations with RHS in columns of B
void Det(Double_t& d1, Double_t& d2)
 Matrix determinant det = d1*TMath::Power(2.,d2)
void Print(Option_t* opt = "") const
 Print class members
TDecompBase & operator=(const TDecompBase& source)
 Assignment operator
Double_t SetTol(Double_t newTol)
void ResetStatus()
{ for (Int_t i = 14; i < 22; i++) ResetBit(BIT(i)); }
const TMatrixDBase & GetDecompMatrix() const
virtual ~TDecompBase()
{}
Double_t GetTol() const
{ return fTol; }
Double_t GetDet1() const
{ return fDet1; }
Double_t GetDet2() const
{ return fDet2; }
Double_t GetCondition() const
{ return fCondition; }
Int_t GetNrows() const
Int_t GetNcols() const
Int_t GetRowLwb() const
{ return fRowLwb; }
Int_t GetColLwb() const
{ return fColLwb; }
Bool_t Decompose()
Bool_t Solve( TVectorD &b)
TVectorD Solve(const TVectorD& b, Bool_t& ok)
Bool_t Solve( TMatrixDColumn& b)
Bool_t TransSolve( TVectorD &b)
TVectorD TransSolve(const TVectorD& b, Bool_t& ok)
Bool_t TransSolve( TMatrixDColumn& b)