88 if (
a.GetNrows() !=
a.GetNcols() ||
a.GetRowLwb() !=
a.GetColLwb()) {
89 Error(
"TDecompLU(const TMatrixD &",
"matrix should be square");
130 Error(
"Decompose()",
"Matrix has not been set");
153 Error(
"GetMatrix()",
"Matrix is singular");
158 Error(
"GetMatrix()",
"Decomposition failed");
182 for (
Int_t i =
n-1; i >= 0; i--) {
187 for (
Int_t k = 0; k <
n; k++) {
206 if (
a.GetNrows() !=
a.GetNcols() ||
a.GetRowLwb() !=
a.GetColLwb()) {
207 Error(
"TDecompLU(const TMatrixD &",
"matrix should be square");
236 Error(
"Solve()",
"Matrix is singular");
241 Error(
"Solve()",
"Decomposition failed");
247 Error(
"Solve(TVectorD &",
"vector and matrix incompatible");
259 for (i = 0; i <
n ; i++) {
269 for (i = 0; i <
n; i++) {
283 for (i =
n-1; i >= 0; i--) {
303 Error(
"Solve()",
"Matrix is singular");
308 Error(
"Solve()",
"Decomposition failed");
314 Error(
"Solve(TMatrixDColumn &",
"vector and matrix incompatible");
324 for (i = 0; i <
n ; i++) {
327 Error(
"Solve(TMatrixDColumn &cb)",
"LU[%d,%d]=%.4e < %.4e",i,i,
pLU[
off_i+i],
fTol);
337 for (i = 0; i <
n; i++) {
353 for (i =
n-1; i >= 0; i--) {
374 Error(
"TransSolve()",
"Matrix is singular");
379 Error(
"TransSolve()",
"Decomposition failed");
385 Error(
"TransSolve(TVectorD &",
"vector and matrix incompatible");
397 for (i = 0; i <
n ; i++) {
406 for (i = 0; i <
n; i++) {
418 for (i =
n-1 ; i >= 0; i--) {
444 Error(
"TransSolve()",
"Matrix is singular");
449 Error(
"TransSolve()",
"Decomposition failed");
455 Error(
"TransSolve(TMatrixDColumn &",
"vector and matrix incompatible");
467 for (i = 0; i <
n ; i++) {
470 Error(
"TransSolve(TMatrixDColumn &cb)",
"LU[%d,%d]=%.4e < %.4e",i,i,
pLU[
off_i+i],
fTol);
476 for (i = 0; i <
n; i++) {
488 for (i =
n-1 ; i >= 0; i--) {
529 Error(
"Invert(TMatrixD &",
"Input matrix has wrong shape");
615 for (
Int_t i = 0; i <
n ; i++) {
623 scale[i] = (max == 0.0 ? 0.0 : 1.0/max);
629 for (
Int_t i = 0; i <
j; i++) {
632 for (
Int_t k = 0; k < i; k++) {
646 for (
Int_t i =
j; i <
n; i++) {
649 for (
Int_t k = 0; k <
j; k++) {
664 for (
Int_t k = 0; k <
n; k++ ) {
680 for (
Int_t i =
j+1; i <
n; i++) {
686 ::Error(
"TDecompLU::DecomposeLUCrout",
"matrix is singular");
725 for (
Int_t i =
j+1; i <
n; i++) {
737 for (
Int_t k = 0; k <
n; k++ ) {
751 for (
Int_t i =
j+1; i <
n; i++) {
756 for (
Int_t k =
j+1; k <
n; k++) {
763 ::Error(
"TDecompLU::DecomposeLUGauss",
"matrix is singular");
779 if (
lu.GetNrows() !=
lu.GetNcols() ||
lu.GetRowLwb() !=
lu.GetColLwb()) {
780 ::Error(
"TDecompLU::InvertLU",
"matrix should be square");
800 ::Error(
"TDecompLU::InvertLU",
"matrix is singular, %d diag elements < tolerance of %.4e",
nrZeros,
tol);
817 for (
j = 0;
j <
n;
j++) {
827 for (k = 0; k <=
j-1; k++) {
831 for (
Int_t i = 0; i <= k-1; i++) {
838 for (k = 0; k <=
j-1; k++) {
854 for (
j =
n-1;
j >= 0;
j--) {
857 for (
Int_t i =
j+1; i <
n; i++) {
885 for (
j =
n-1;
j >= 0;
j--) {
888 for (
Int_t i = 0; i <
n; i++) {
const char Option_t
Option string (const char)
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
#define R__ASSERT(e)
Checks condition e and reports a fatal error if it's false.
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t r
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t index
TMatrixTDiag_const< Double_t > TMatrixDDiag_const
TMatrixT< Double_t > TMatrixD
Decomposition Base class.
static void DiagProd(const TVectorD &diag, Double_t tol, Double_t &d1, Double_t &d2)
virtual Bool_t MultiSolve(TMatrixD &B)
Solve set of equations with RHS in columns of B.
TDecompBase & operator=(const TDecompBase &source)
Assignment operator.
void Print(Option_t *opt="") const override
Print class members.
virtual void Det(Double_t &d1, Double_t &d2)
Matrix determinant det = d1*TMath::Power(2.,d2)
Int_t GetNrows() const override
Bool_t TransSolve(TVectorD &b) override
Solve A^T x=b assuming the LU form of A^T is stored in fLU, but assume b has not been transformed.
Int_t GetNcols() const override
virtual void SetMatrix(const TMatrixD &a)
Set matrix to be decomposed.
static Bool_t DecomposeLUGauss(TMatrixD &lu, Int_t *index, Double_t &sign, Double_t tol, Int_t &nrZeros)
LU decomposition using Gaussian Elimination with partial pivoting (See Golub & Van Loan,...
TDecompLU()
Default constructor.
const TMatrixD GetMatrix()
Reconstruct the original matrix using the decomposition parts.
static 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.
void Det(Double_t &d1, Double_t &d2) override
Calculate determinant det = d1*TMath::Power(2.,d2)
void Print(Option_t *opt="") const override
Print internals of this object.
static Bool_t InvertLU(TMatrixD &a, Double_t tol, Double_t *det=nullptr)
Calculate matrix inversion through in place forward/backward substitution.
Bool_t Solve(TVectorD &b) override
Solve Ax=b assuming the LU form of A is stored in fLU, but assume b has not been transformed.
Bool_t Decompose() override
Matrix A is decomposed in components U and L so that P * A = U * L If the decomposition succeeds,...
TDecompLU & operator=(const TDecompLU &source)
assignment operator
void Print(Option_t *name="") const override
Print the matrix as a table of elements.
const TMatrixTBase< Element > * GetMatrix() const
const Element * GetMatrixArray() const override
TMatrixTBase< Element > & ResizeTo(Int_t nrows, Int_t ncols, Int_t=-1) override
Set size of the matrix to nrows x ncols New dynamic elements are created, the overlapping part of the...
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Returns x raised to the power y.
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
void inv(rsa_NUMBER *, rsa_NUMBER *, rsa_NUMBER *)
static uint64_t sum(uint64_t i)