Logo ROOT   6.18/05
Reference Guide
TDecompBK.h
Go to the documentation of this file.
1// @(#)root/matrix:$Id$
2// Authors: Fons Rademakers, Eddy Offermann Sep 2004
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_TDecompBK
13#define ROOT_TDecompBK
14
15///////////////////////////////////////////////////////////////////////////
16// //
17// Bunch-Kaufman Decomposition class //
18// //
19///////////////////////////////////////////////////////////////////////////
20
21#include "Rtypes.h"
22#include "TDecompBase.h"
23#include "TMatrixDSym.h"
24#include "TVectorD.h"
25
26class TDecompBK : public TDecompBase
27{
28protected :
29
30 Int_t fNIpiv; // size of row permutation index
31 Int_t *fIpiv; //[fNIpiv] row permutation index
32 TMatrixD fU; // decomposed matrix so that a = u d u^T
33
34 virtual const TMatrixDBase &GetDecompMatrix() const { return fU; }
35
36public :
37
38 TDecompBK();
39 explicit TDecompBK(Int_t nrows);
40 TDecompBK(Int_t row_lwb,Int_t row_upb);
41 TDecompBK(const TMatrixDSym &m,Double_t tol = 0.0);
42 TDecompBK(const TDecompBK &another);
43 virtual ~TDecompBK() {if (fIpiv) delete [] fIpiv; fIpiv = 0; }
44
45 virtual Int_t GetNrows () const { return fU.GetNrows(); }
46 virtual Int_t GetNcols () const { return fU.GetNcols(); }
47 const TMatrixD &GetU () { if ( !TestBit(kDecomposed) ) Decompose();
48 return fU; }
49
50 virtual void SetMatrix (const TMatrixDSym &a);
51
52 virtual Bool_t Decompose ();
53 virtual Bool_t Solve ( TVectorD &b);
54 virtual TVectorD Solve (const TVectorD& b,Bool_t &ok) { TVectorD x = b; ok = Solve(x); return x; }
55 virtual Bool_t Solve ( TMatrixDColumn &b);
56 virtual Bool_t TransSolve ( TVectorD &b) { return Solve(b); }
57 virtual TVectorD TransSolve (const TVectorD& b,Bool_t &ok) { TVectorD x = b; ok = Solve(x); return x; }
58 virtual Bool_t TransSolve ( TMatrixDColumn &b) { return Solve(b); }
59 virtual void Det (Double_t &/*d1*/,Double_t &/*d2*/)
60 { MayNotUse("Det(Double_t&,Double_t&)"); }
61
63 TMatrixDSym Invert (Bool_t &status);
64 TMatrixDSym Invert () { Bool_t status; return Invert(status); }
65
66 void Print(Option_t *opt ="") const; // *MENU*
67
68 TDecompBK &operator= (const TDecompBK &source);
69
70 ClassDef(TDecompBK,1) // Matrix Decomposition Bunch-Kaufman
71};
72
73#endif
#define b(i)
Definition: RSha256.hxx:100
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
double Double_t
Definition: RtypesCore.h:55
const char Option_t
Definition: RtypesCore.h:62
#define ClassDef(name, id)
Definition: Rtypes.h:326
The Bunch-Kaufman diagonal pivoting method decomposes a real symmetric matrix A using.
Definition: TDecompBK.h:27
TDecompBK & operator=(const TDecompBK &source)
Assignment operator.
Definition: TDecompBK.cxx:660
TDecompBK()
Default constructor.
Definition: TDecompBK.cxx:64
const TMatrixD & GetU()
Definition: TDecompBK.h:47
TMatrixDSym Invert()
Definition: TDecompBK.h:64
virtual TVectorD Solve(const TVectorD &b, Bool_t &ok)
Definition: TDecompBK.h:54
virtual Bool_t Solve(TVectorD &b)
Solve Ax=b assuming the BK form of A is stored in fU . Solution returned in b.
Definition: TDecompBK.cxx:339
virtual TVectorD TransSolve(const TVectorD &b, Bool_t &ok)
Definition: TDecompBK.h:57
Int_t fNIpiv
Definition: TDecompBK.h:30
virtual Bool_t TransSolve(TMatrixDColumn &b)
Definition: TDecompBK.h:58
void Print(Option_t *opt="") const
Print the class members.
Definition: TDecompBK.cxx:648
virtual const TMatrixDBase & GetDecompMatrix() const
Definition: TDecompBK.h:34
virtual ~TDecompBK()
Definition: TDecompBK.h:43
TMatrixD fU
Definition: TDecompBK.h:32
virtual void SetMatrix(const TMatrixDSym &a)
Set the matrix to be decomposed, decomposition status is reset.
Definition: TDecompBK.cxx:314
virtual Int_t GetNrows() const
Definition: TDecompBK.h:45
virtual Bool_t TransSolve(TVectorD &b)
Definition: TDecompBK.h:56
virtual Bool_t Decompose()
Matrix A is decomposed in components U and D so that A = U*D*U^T If the decomposition succeeds,...
Definition: TDecompBK.cxx:131
virtual Int_t GetNcols() const
Definition: TDecompBK.h:46
virtual void Det(Double_t &, Double_t &)
Matrix determinant det = d1*TMath::Power(2.,d2)
Definition: TDecompBK.h:59
Int_t * fIpiv
Definition: TDecompBK.h:31
Decomposition Base class.
Definition: TDecompBase.h:34
Int_t GetNrows() const
Definition: TMatrixTBase.h:124
Int_t GetNcols() const
Definition: TMatrixTBase.h:127
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:172
void MayNotUse(const char *method) const
Use this method to signal that a method (defined in a base class) may not be called in a derived clas...
Definition: TObject.cxx:933
Double_t x[n]
Definition: legend1.C:17
void inv(rsa_NUMBER *, rsa_NUMBER *, rsa_NUMBER *)
Definition: rsaaux.cxx:949
auto * m
Definition: textangle.C:8
auto * a
Definition: textangle.C:12