library: libMatrix #include "TDecompBase.h"

# TDecompBase

class description - source file - inheritance tree (.pdf)

## class TDecompBase : public TObject

| Inheritance Chart:
`TObject`
| <- | `TDecompBase` |
<- |
`TDecompBK` |
`TDecompChol` |
`TDecompLU` |
`TDecompQRH` |
`TDecompSparse` |
`TDecompSVD` |

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

        protected:

                          static void DiagProd(const TVectorD& diag, Double_t tol, Double_t& d1, Double_t& d2)
          virtual const TMatrixDBase& GetDecompMatrix() const
                                Int_t Hager(Double_t& est, Int_t iter = 5)
                                 void ResetStatus()

        public:

                   virtual ~TDecompBase()
            static TClass* Class()
          virtual Double_t Condition()
            virtual Bool_t Decompose()
              virtual void Det(Double_t& d1, Double_t& d2)
                     Int_t GetColLwb() const
                  Double_t GetCondition() const
                  Double_t GetDet1() const
                  Double_t GetDet2() const
             virtual Int_t GetNcols() const
             virtual Int_t GetNrows() const
                     Int_t GetRowLwb() const
                  Double_t GetTol() const
           virtual TClass* IsA() const
            virtual Bool_t MultiSolve(TMatrixD& B)
              TDecompBase& operator=(const TDecompBase& source)
              virtual void Print(Option_t* opt = "") const
                  Double_t SetTol(Double_t tol)
              virtual void ShowMembers(TMemberInspector& insp, char* parent)
            virtual Bool_t Solve(TVectorD& b)
          virtual TVectorD Solve(const TVectorD& b, Bool_t& ok)
            virtual Bool_t Solve(TMatrixDColumn& b)
              virtual void Streamer(TBuffer& b)
                      void StreamerNVirtual(TBuffer& b)
            virtual Bool_t TransSolve(TVectorD& b)
          virtual TVectorD TransSolve(const TVectorD& b, Bool_t& ok)
            virtual Bool_t TransSolve(TMatrixDColumn& b)

### Data Members

        protected:

          Double_t fTol        sqrt(epsilon); epsilon is smallest number number so that  1+epsilon > 1
          Double_t fDet1       determinant mantissa
          Double_t fDet2       determinant exponent for powers of 2
          Double_t fCondition  matrix condition number
             Int_t fRowLwb     Row    lower bound of decomposed matrix
             Int_t fColLwb     Column lower bound of decomposed matrix

        public:

          static const TDecompBase::EMatrixDecompStat kInit
          static const TDecompBase::EMatrixDecompStat kPatternSet
          static const TDecompBase::EMatrixDecompStat kValuesSet
          static const TDecompBase::EMatrixDecompStat kMatrixSet
          static const TDecompBase::EMatrixDecompStat kDecomposed
          static const TDecompBase::EMatrixDecompStat kDetermined
          static const TDecompBase::EMatrixDecompStat kCondition
          static const TDecompBase::EMatrixDecompStat kSingular
                      static const enum TDecompBase:: kWorkMax

## Class Description

     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

Int_t Hager(Double_t &est,Int_t iter)

void DiagProd(const TVectorD &diag,Double_t tol,Double_t &d1,Double_t &d2)

Double_t Condition()

Bool_t MultiSolve(TMatrixD &B)

     Solve set of equations with RHS in columns of B

void Det(Double_t &d1,Double_t &d2)

void Print(Option_t * /*opt*/) const

### Inline Functions

                           void ~TDecompBase()
                           void ResetStatus()
            const TMatrixDBase& GetDecompMatrix() const
                       Double_t GetTol() const
                       Double_t GetDet1() const
                       Double_t GetDet2() const
                       Double_t GetCondition() const
                          Int_t GetNrows() const
                          Int_t GetNcols() const
                          Int_t GetRowLwb() const
                          Int_t GetColLwb() const
                       Double_t SetTol(Double_t tol)
                         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)
                   TDecompBase& operator=(const TDecompBase& source)
                        TClass* Class()
                        TClass* IsA() const
                           void ShowMembers(TMemberInspector& insp, char* parent)
                           void Streamer(TBuffer& b)
                           void StreamerNVirtual(TBuffer& b)

Last update: root/matrix:$Name: $:$Id: TDecompBase.cxx,v 1.16 2004/10/16 18:09:16 brun Exp $
Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *