library: libMatrix
#include "TDecompLU.h"

TDecompLU


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

class TDecompLU : public TDecompBase

Inheritance Chart:
TObject
<-
TDecompBase
<-
TDecompLU

    protected:
static Bool_t DecomposeLUCrout(TMatrixD& lu, Int_t* index, Double_t& sign, Double_t tol, Int_t& nrZeros) static Bool_t DecomposeLUGauss(TMatrixD& lu, Int_t* index, Double_t& sign, Double_t tol, Int_t& nrZeros) virtual const TMatrixDBase& GetDecompMatrix() const public:
TDecompLU() TDecompLU(Int_t nrows) TDecompLU(Int_t row_lwb, Int_t row_upb) TDecompLU(const TMatrixD& m, Double_t tol = 0.0, Int_t implicit = 1) TDecompLU(const TDecompLU& another) virtual ~TDecompLU() static TClass* Class() virtual Bool_t Decompose() virtual void Det(Double_t& d1, Double_t& d2) const TMatrixD& GetLU() const const TMatrixD GetMatrix() const virtual Int_t GetNcols() const virtual Int_t GetNrows() const void Invert(TMatrixD& inv) TMatrixD Invert() static Bool_t InvertLU(TMatrixD& a, Double_t tol, Double_t* det = 0) virtual TClass* IsA() const TDecompLU& operator=(const TDecompLU& source) virtual void Print(Option_t* opt = "") const virtual void SetMatrix(const TMatrixD& a) 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:
Int_t fImplicitPivot control to determine implicit row scale before Int_t fNIndex size of row permutation index Int_t* fIndex [fNIndex] row permutation index Double_t fSign = +/- 1 reflecting even/odd row permutations, resp. TMatrixD fLU decomposed matrix so that a = l u where

Class Description

                                                                       
 LU Decomposition class                                                
                                                                       
 Decompose  a general n x n matrix A into P A = L U                    
                                                                       
 where P is a permutation matrix, L is unit lower triangular and U     
 is upper triangular.                                                  
 L is stored in the strict lower triangular part of the matrix fLU.    
 The diagonal elements of L are unity and are not stored.              
 U is stored in the diagonal and upper triangular part of the matrix   
 fU.                                                                   
 P is stored in the index array fIndex : j = fIndex[i] indicates that  
 row j and row i should be swapped .                                   
                                                                       
 fSign gives the sign of the permutation, (-1)^n, where n is the       
 number of interchanges in the permutation.                            
                                                                       
 fLU has the same indexing range as matrix A .                         
                                                                       
 The decomposition fails if a diagonal element of abs(fLU) is == 0,    
 The matrix fUL is made invalid .                                      
                                                                       


TDecompLU()

TDecompLU(Int_t nrows)

TDecompLU(Int_t row_lwb,Int_t row_upb)

TDecompLU(const TMatrixD &a,Double_t tol,Int_t implicit)

TDecompLU(const TDecompLU &another) : TDecompBase(another)

Bool_t Decompose()

const TMatrixD GetMatrix()
 Reconstruct the original matrix using the decomposition parts

void SetMatrix(const TMatrixD &a)

Bool_t Solve(TVectorD &b)
 Solve Ax=b assuming the LU form of A is stored in fLU, but assume b has *not*
 been transformed.  Solution returned in b.

Bool_t Solve(TMatrixDColumn &cb)
 Solve Ax=b assuming the LU form of A is stored in fLU, but assume b has *not*
 been transformed.  Solution returned in b.

Bool_t TransSolve(TVectorD &b)
 Solve A^T x=b assuming the LU form of A^T is stored in fLU, but assume b has *not*
 been transformed.  Solution returned in b.

Bool_t TransSolve(TMatrixDColumn &cb)
 Solve A^T x=b assuming the LU form of A^T is stored in fLU, but assume b has *not*
 been transformed.  Solution returned in b.

void Det(Double_t &d1,Double_t &d2)

void Invert(TMatrixD &inv)
 For a matrix A(m,m), its inverse A_inv is defined as A * A_inv = A_inv * A = unit
 (m x m) Ainv is returned .

TMatrixD Invert()
 For a matrix A(m,n), its inverse A_inv is defined as A * A_inv = A_inv * A = unit
 (n x m) Ainv is returned .

void Print(Option_t *opt) const

Bool_t DecomposeLUCrout(TMatrixD &lu,Int_t *index,Double_t &sign, Double_t tol,Int_t &nrZeros)
 Crout/Doolittle algorithm of LU decomposing a square matrix, with implicit partial
 pivoting.  The decomposition is stored in fLU: U is explicit in the upper triag
 and L is in multiplier form in the subdiagionals .
 Row permutations are mapped out in fIndex. fSign, used for calculating the
 determinant, is +/- 1 for even/odd row permutations. .

Bool_t DecomposeLUGauss(TMatrixD &lu,Int_t *index,Double_t &sign, Double_t tol,Int_t &nrZeros)
 LU decomposition using Gaussain Elimination with partial pivoting (See Golub &
 Van Loan, Matrix Computations, Algorithm 3.4.1) of a square matrix .
 The decomposition is stored in fLU: U is explicit in the upper triag and L is in
 multiplier form in the subdiagionals . Row permutations are mapped out in fIndex.
 fSign, used for calculating the determinant, is +/- 1 for even/odd row permutations.
 Since this algorithm uses partial pivoting without scaling like in Crout/Doolitle.
 it is somewhat faster but less precise .

Bool_t InvertLU(TMatrixD &lu,Double_t tol,Double_t *det)
 Calculate matrix inversion through in place forward/backward substitution



Inline Functions


                       void ~TDecompLU()
        const TMatrixDBase& GetDecompMatrix() const
                      Int_t GetNrows() const
                      Int_t GetNcols() const
            const TMatrixD& GetLU() const
                     Bool_t Solve(TMatrixDColumn& b)
                     Bool_t TransSolve(TMatrixDColumn& b)
                 TDecompLU& operator=(const TDecompLU& 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: TDecompLU.cxx,v 1.23 2005/09/02 11:04:45 brun Exp $
Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *


ROOT page - Class index - Class Hierarchy - Top of the page

This page has been automatically generated. If you have any comments or suggestions about the page layout send a mail to ROOT support, or contact the developers with any questions or problems regarding ROOT.