Log of /trunk/math/matrix/src/TDecompLU.cxx
Parent Directory
Revision
22885 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Fri Mar 28 13:57:25 2008 UTC (6 years, 9 months ago) by
rdm
File length: 26237 byte(s)
Diff to
previous 22039
move the following directories under the new "math" meta directory:
mathcore
mathmore
fftw
foam
fumili
genvector
matrix
minuit
minuit2
mlp
physics
smatrix
splot
unuran
quadp
Revision
22039 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Thu Feb 7 05:48:31 2008 UTC (6 years, 11 months ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 26237 byte(s)
Diff to
previous 20882
From Eddy:
- remove bug in Decompose of all TDecompXXX classes : It is checked
whether
the matrix was already decomposed
- Added to TVector and the different TMatrix flavors the option to use
"Use" on constant arrys and vectors/matrices.
Revision
16458 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Fri Oct 6 06:52:34 2006 UTC (8 years, 3 months ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 26236 byte(s)
Diff to
previous 15325
From Eddy Offermann:
Previously, the matrix package contained a large amount of ASSERT
statements which are a nuisance for programs analyzing a series of
independent events like in high-energy physics .
Assert's were issued when for instance a division by zero was requested
or a matrix was invalid . Most algorithms made matrices/vectors invalid
after an error occured in an operation , like inversion of a singular matrix .
Unfortunately, not all assert's were accompanied by error messages .
This situation has been completely overhauled :
- All error conditions in the algorithms are now accompanied by error
messages .
- In all algorithms it is still asserted that vectors/matrices are
valid BUT only in very few cases is a matrix/vector made invalid :
for instance if memory is allocated with incorrect parameters .
- In case of division by zero, the division is skipped . In case of a
singular matrix, the inversion routine returns the original matrix .
In the past the result of an inversion could be checked through the
value of the returned determinant or checking whether the inverted matrix
was valid .
Since from now on, we never make the matrix invalid in this operation, the
latter check will not indicate a singular matrix anymore .
The decompostion classes TDecomp... have a backward-compatible
change in the interface which makes detection of singularity easier :
old interface :
void Invert (TMatrixD &inv);
TMatrixD Invert ();
new interface :
Bool_t Invert (TMatrixD &inv);
TMatrixD Invert (Bool_t &status);
TMatrixD Invert () { Bool_t status; return Invert(status);
}
The returned status is kFALSE in case of singularity .
The old situation is easily reproduced by setting the ROOT variable
gErrorAbortLevel to kError . This cause an exception when there is an error
message (Error...) ., In the past the matrix would be made invalid which
would cause the next operation to throw an exception .
Revision
15300 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Fri Jun 2 05:11:20 2006 UTC (8 years, 7 months ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 26056 byte(s)
Diff to
previous 15176
From Eddy:
Worked more on the comments in the matrix package . Actually the
important
ones are all there .. but have to add thing slike "cop constructor" to
make the checker happy .
Added also namespaces TMatrixTCramer and TMatrixTSymCramerInv to the
dictionary
Revision
14745 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Wed Apr 19 08:22:26 2006 UTC (8 years, 9 months ago) by
rdm
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 24451 byte(s)
Diff to
previous 12591
Change the TError.h macros:
Assert -> R__ASSERT
Check -> R__CHECK
Change the TCollection.h macro:
ForEach -> R__FOR_EACH
This to avoid potential problems due too trivial macro names.
The old macros will be removed in the next release. Currently
they will print out warning messages with the advice to move
to the new macro names.
Revision
12111 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Wed Jun 22 12:37:09 2005 UTC (9 years, 7 months ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 24383 byte(s)
Diff to
previous 11489
From Eddy Offermann:
remove a warning in the calculation of the determinant in TDecompLU :
The idea was that in case the d2 was small/large, somewhere in the
calculation situations would appear of the type 1+epsilon where
epsilon < 2^-52 which will cause an underflow . Howevere, the current
test is not sufficient and will surely lead to erroneous warnings .
Revision
11119 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Tue Feb 15 16:17:10 2005 UTC (9 years, 11 months ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 24259 byte(s)
Diff to
previous 10711
From Eddy Offermann
Here a patch to remove some inconsistencies/bugs in the indexing of the
decomposed matrices and eigen-vectors/-values in case the row/column
index
did not start 0 .
The tests in stressLinear were modified/extended to check these cases .
Revision
10711 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Tue Dec 7 19:11:26 2004 UTC (10 years, 1 month ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 24070 byte(s)
Diff to
previous 10667
From Eddy Offermann:
Instead of warning of a under/over flow in TDecompLU::InvertLU and then
returning a determinant of zero , we now return the calculated value .
Now the user will not be confused anymore and think that the matrix
was singular .
Revision
10329 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Sat Oct 16 18:09:17 2004 UTC (10 years, 3 months ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 24133 byte(s)
Diff to
previous 9893
From Eddy Offermann:
- added new decomposition class TDecompBK: Bunch-Kaufman algorithm
this class is designed to decompose real symmetric matrices . Now
we can guarantee that inverse of a symmetric matrix is again
symmetric .
- As a consequence of the introduction of TDecompBK, we added
Invert(..) functionality to TMatrixD/FSym
Also added InvertFast to TMatrixD/FSym which uses Cramer
inversion for matrices <= (6 x 6) , resulting in large speedup
but less accuracy . The necessary routines are encapsulated
in the new TMatrixDSymCramerInv class .
Note of caution :
Unlike the other decompositions, Bunch-Kaufmann does not result
in a triangular matrix of which the determinant is trivially calculated
by multiplying the diagonal elements . As a consequence,
we have to invoke yet another decomposition when the user
requests the determinant .
So Invert() will be roughly twice faster than Double_t *det; Invert(det)
- Changed the algoritm for the Cholesky decomposition so that we
do not need to have a copy of the original matrix around !
After removal of class member fA, increased the version to 2
- Added sorting of eigen -values/-vectors to TMatrixDEigen .
They are now sorted according to Re^2+Im^2 of the (possible)
complex eigenvalues in descending order . TMatrixDSymEigen
already sorted . Since most users use TMatrixD instead of
TMatrixDSym, this will ensure backward compatiblity since
old root contained the symmetric algorithms in TMatrixD .
- Added rank 1 update algorihms : A += alpha * x y^T
(x and y are vectors)
to the non-sparse matrix classes (adding this to sparse
would make this matrix immediately dense :)
- Moved Invert fom TDecompBase to the different decomposition
classes . Reason is that some classes can decompose A(m,n)
where m > n . The result of Invert() is then a "pseudo-invert"
: A_pseudo = (A^T A)^-1 A^T .
- Removed an unfortunate bug in the TMatrixD/FSub class :
both operators
void operator+=(Double_t val);
void operator*=(Double_t val);
performed actually
void operator=(Double_t val);
~
Revision
9893 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Fri Sep 3 13:41:35 2004 UTC (10 years, 4 months ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 24346 byte(s)
Diff to
previous 9503
From Eddy Offermann
The matrix includes contained circular dependencies (which were of
course
stopped by the #ifdef's) . Situation has been cleaned up by inserting
several forward declarations in the header files and moving some
includes
to the source files . A change in a header file now causes less
compilation
in other directories like hist.
This cleanup made adding an include file in quadp/src/TQpVar necessary
and some additional include files in the tutorials and test files
Revision
9503 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Mon Jul 12 20:00:41 2004 UTC (10 years, 6 months ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 24346 byte(s)
Diff to
previous 9181
From Eddy Offermann:
Some code re-organization:
- made TDecompLU::DecomposeLUCrout and TDecompLU::DecomposeLUGauss
protected.
As a consequence changed TDecompLU::InvertLU so that it now calls
the decomposition step .
This change prohibits the user of calling InvertLU without applying
first
the decomposition step
- Routines in TMatrixD/FCramerInv return now a Bool_t instead of Int_t
.
More important, these routines check that the matrix is indeed of the
appropriate shape, since the user could call them directly instead
of going through TMatrixD/F::InvertFast
- Removed the redundant TDecompBase::MultiSolve(TMatrixDSym &B), it did
exactly what TDecompBase::MultiSolve(TMatrixD &B) did .
- void (Trans)Solve(TVectorD &b,Int_t ok) and
void (Trans)Solve(TMatrixDColumn &cb, Int_t ok) routines in TDecomp*
now make vector b or matrix b in case of an error (when ok was set to
false)
Revision
8892 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Wed May 12 10:39:29 2004 UTC (10 years, 8 months ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 22457 byte(s)
Diff to
previous 8456
From Eddy Offermann:
1) patch for the matrix directory (cvs diff -a matrix)
2) tar file containing new modules:
matrix/inc/TDecompSparse.h
matrix/inc/TMatrixDSparse.h
matrix/src/TDecompSparse.cxx
matrix/src/TMatrixDSparse.cxx
3) new stressLinear.cxx
Description of changes:
1) many updates of comments in code
2) possibility to construct random matrices, general
symmetric and positive definite
3) Introduction of new matrix class: TMatrixDSparse .
It defines a general sparse matrix in
Harwell-Boeing
format (sparse structure definition optimized for
matrix-vector multiplications).
Like all the other matrix classes, it derives from
TMatrixDBase
4) Introduction of a new decomposition class:
TDecompSparse . It allows decomposition (and
therefore equation solving) of sparse symmetric
matrices
Revision
8456 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Mon Mar 22 08:34:36 2004 UTC (10 years, 10 months ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 22450 byte(s)
Diff to
previous 8442
From Eddy:
Added to the decomposition class constructors in
which just the matrix size is given as argument .
Also a new state (kMatrixSet) was added .
Instead of
TMatrixD a;
TDecompLU lu(a);
one can do
TMatrixD a;
TDecompLU lu(n);
lu.SetMatrix(a);
The SetMatrix command resizes, if necessary, the
data members in the decomposition class .
This set up willl be usefull when one wants to
decompose
matrices in a loop
Revision
8442 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Fri Mar 19 14:20:40 2004 UTC (10 years, 10 months ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 21860 byte(s)
Diff to
previous 8189
From Eddy Offermann:
1. Removed streaming of matrix view classes (TMatrixDRow,TMatrixDColumn,TMatrixDDiag,TMatrixDFlat)
because a. we want these classes to be as light as possible
b. old implementation was wrong, some member elements were "(!)" out
2. Added a series of functions to support upcoming quadratic programming functionality:
TVectorD/F:
- AddElemMult(TVectorF &t,Float_t alpha,const TVectorF &s1,const TVectorF &s2);
t += alpha * s1*s2 , where "s1*s2 is an element-wise multiplication
- AddElemMult(TVectorF &t,Float_t alpha,const TVectorF &s1,const TVectorF &s2,const TVectorF &map);
t += alpha * s1*s2 , where "s1*s2 is an element-wise multiplication. It is only performed
for those elements i where map[i] != 0
- AddElemDiv (TVectorF &t,Float_t alpha,const TVectorF &s1,const TVectorF &s2);
t += alpha * s1\s2 , where "s1\s2 is an element-wise division
- AddElemDiv (TVectorF &t,Float_t alpha,const TVectorF &s1,const TVectorF &s2,const TVectorF &map);
t += alpha * s1\s2 , where "s1\s2 is an element-wise division. It is only performed
for those elements i where map[i] != 0
- ElementMult(TVectorF &t,const TVectorF &s,const TVectorF &map);
t += t*s , where "t*s" is an element-wise multiplication. It is only performed
for those elements i where map[i] != 0
- ElementDiv (TVectorF &t,const TVectorF &s,const TVectorF &map);
t += t/s , where "t/s" is an element-wise division. It is only performed
for those elements i where map[i] != 0
- TVectorD &SelectNonZeros (const TVectorD &select);
put elements for which select[i] == 0 to zero
- Bool_t MatchesNonZeroPattern(const TVectorD &select);
verify that elements are != 0 for which select[i] != 0
- Bool_t SomePositive (const TVectorD &select);
verify that all elements > 0 for which select[i] != 0
- void AddSomeConstant (Double_t val,const TVectorD &select);
Add val to those elements for which select[i] != 0
- TVectorF &TVectorF::Invert()
v[i] = 1/v[i]
TVectorD/F and TMatrixD/FBase:
- virtual Double_t Sum () const;
- virtual Double_t Min () const;
- virtual Double_t Max () const;
3. Added another data member fA to TDecompChol . Unlike the other decomposition schemes, it is not
possible to store the original input matrix in one of the decomposition members
4. Added method "SetMatrix(const TMatrix.. &)" to all decomposition classes . It allows to
change the matrix that should be decomposed. Of course it resets the class status to kInit
and triggers a new factorization.
5. Added private methods to TMatrixD and TMatrixDSym
const TMatrixD EigenVectors(TVectorD &eigenValues) const; // This function is now obsolete (and not implemented), you
// should use TMatrixDSymEigen or TMatrixDEigen .
6. Changed the name "Adopt" to "Use" to be in agreement with the convention that the memory management
of adopted members is now the reponsibility of the class
Revision
8118 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Wed Feb 4 17:12:44 2004 UTC (10 years, 11 months ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 20339 byte(s)
Diff to
previous 8112
From Eddy Offermann:
1. LinkDef.h : Added NormalEqn friend functions from
TDecompChol class
2. TMatrixDBase: GetDecompMatrix() returns a TMatrixD
instead of TMatrixDBase
Added Solve and TransSolve functions
that return the solution (others
replace
the input argument which is
inconvinient
when output vector size is smaller
than
input (in ncase of fit).
Also added Invert that returns matrix
instead of replacing the input
argument.
This is convinient for a generalized
(pseudo) inverse: in (m x n) out (n x
m)
3. TDecompChol: Solve/TrnasSolve had incorrect
forward/
backward substitution function ! Tests
have been added to stressLinear to
check
Added NormlEqn solver which can handle
weights in data
4. TDecompSVD: Bug in Solve function ! index was
stop-
ping one too early . Was not noticed
by
stressLinear because the linear
equation
test wasnumerically to difficult
causing accuracy tolerances to be not
strict enough (was 0.005 , is now
DBL_EPSILON on a diagonal dominant
matrix)
Revision
8062 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Modified
Wed Jan 28 07:39:18 2004 UTC (10 years, 11 months ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 19695 byte(s)
Diff to
previous 8046
From Eddy Offermann:
Reacting to Philippe's worry , I have removed
one TDecompLU constructor .
Answering George's Irwin question, I noticed
that I forgot to forward a change for the
inversion of symmetric matrices: I was
using the Cholesky decomp which wants the matrix
to be also positive definite
The patch rectifies these issues .
Revision
8013 -
(
view)
(
download)
(
as text)
(
annotate)
-
[select for diffs]
Added
Sun Jan 25 20:33:32 2004 UTC (11 years ago) by
brun
Original Path:
trunk/matrix/src/TDecompLU.cxx
File length: 20079 byte(s)
New Linear Algebra package from Eddy Offermann.
This new package reimplements the previous classes TMatrix and TMatrixD.
The new classes should be back compatible with the previous version except
the function GetElements.
New classes have been introduced for symmetric matrices,
lazy matrices.
New algorithms (LU, SVD) have been introduced.
A new test suite test/stressLinear.cxx is introduced.
A complete description of this package will be posted in the coming days.
The classes are well documented in the implementation headers.
This form allows you to request diffs between any two revisions of this file.
For each of the two "sides" of the diff,
enter a numeric revision.
