|
|||||||||
PREV CLASS NEXT CLASS | FRAMES NO FRAMES | ||||||||
SUMMARY: NESTED | FIELD | CONSTR | METHOD | DETAIL: FIELD | CONSTR | METHOD |
java.lang.Objectorg.apache.commons.math.linear.EigenDecompositionImpl
public class EigenDecompositionImpl
Calculates the eigen decomposition of a symmetric matrix.
The eigen decomposition of matrix A is a set of two matrices: V and D such that A = V D VT. A, V and D are all m × m matrices.
As of 2.0, this class supports only symmetric matrices,
and hence computes only real realEigenvalues. This implies the D matrix returned by
getD()
is always diagonal and the imaginary values returned getImagEigenvalue(int)
and getImagEigenvalues()
are always null.
When called with a RealMatrix
argument, this implementation only uses
the upper part of the matrix, the part below the diagonal is not accessed at all.
Eigenvalues are computed as soon as the matrix is decomposed, but eigenvectors
are computed only when required, i.e. only when one of the getEigenvector(int)
,
getV()
, getVT()
, getSolver()
methods is called.
This implementation is based on Inderjit Singh Dhillon thesis A New O(n2) Algorithm for the Symmetric Tridiagonal Eigenvalue/Eigenvector Problem, on Beresford N. Parlett and Osni A. Marques paper An Implementation of the dqds Algorithm (Positive Case) and on the corresponding LAPACK routines (DLARRE, DLASQ2, DLAZQ3, DLAZQ4, DLASQ5 and DLASQ6).
Nested Class Summary | |
---|---|
private static class |
EigenDecompositionImpl.Solver
Specialized solver. |
Field Summary | |
---|---|
private RealMatrix |
cachedD
Cached value of D. |
private RealMatrix |
cachedV
Cached value of V. |
private RealMatrix |
cachedVt
Cached value of Vt. |
private double |
dMin
Minimal value on current state of the diagonal. |
private double |
dMin1
Minimal value on current state of the diagonal, excluding last element. |
private double |
dMin2
Minimal value on current state of the diagonal, excluding last two elements. |
private double |
dN
Last value on current state of the diagonal. |
private double |
dN1
Last but one value on current state of the diagonal. |
private double |
dN2
Last but two on current state of the diagonal. |
private ArrayRealVector[] |
eigenvectors
Eigenvectors. |
private double |
eMin
Min value of off-diagonal elements in current segment. |
private double |
g
Shift ratio with respect to dMin used when tType == 6. |
private double[] |
imagEigenvalues
Imaginary part of the realEigenvalues. |
private double |
lowerSpectra
Lower bound of spectra. |
private double[] |
main
Main diagonal of the tridiagonal matrix. |
private double |
minPivot
Minimum pivot in the Sturm sequence. |
private int |
pingPong
Shift within qd array for ping-pong implementation. |
private double |
qMax
Max value of diagonal elements in current segment. |
private double[] |
realEigenvalues
Real part of the realEigenvalues. |
private double[] |
secondary
Secondary diagonal of the tridiagonal matrix. |
private double |
sigma
Current shift. |
private double |
sigmaLow
Low part of the current shift. |
private double |
splitTolerance
Split tolerance. |
private double[] |
squaredSecondary
Squared secondary diagonal of the tridiagonal matrix. |
private double |
tau
Shift increment to apply. |
private static double |
TOLERANCE
Tolerance. |
private static double |
TOLERANCE_2
Squared tolerance. |
private TriDiagonalTransformer |
transformer
Transformer to tridiagonal (may be null if matrix is already tridiagonal). |
private int |
tType
Type of the last dqds shift. |
private double |
upperSpectra
Upper bound of spectra. |
private double[] |
work
Work array for all decomposition algorithms. |
Constructor Summary | |
---|---|
EigenDecompositionImpl(double[] main,
double[] secondary,
double splitTolerance)
Calculates the eigen decomposition of the given tridiagonal symmetric matrix. |
|
EigenDecompositionImpl(RealMatrix matrix,
double splitTolerance)
Calculates the eigen decomposition of the given symmetric matrix. |
Method Summary | |
---|---|
private void |
computeGershgorinCircles()
Compute the Gershgorin circles for all rows. |
private void |
computeShiftIncrement(int start,
int end,
int deflated)
Compute the shift increment as an estimate of the smallest eigenvalue. |
private List<Integer> |
computeSplits()
Compute splitting points. |
private int |
countEigenValues(double t,
int index,
int n)
Count the number of realEigenvalues below a point. |
private void |
decompose()
Decompose a tridiagonal symmetric matrix. |
private void |
dqd(int start,
int end)
Perform a dqd step. |
private void |
dqds(int start,
int end)
Perform a dqds step, using current shift increment. |
private double[] |
eigenvaluesRange(int index,
int n)
Compute an interval containing all realEigenvalues of a block. |
private void |
findEigenvalues()
Find the realEigenvalues. |
private ArrayRealVector |
findEigenvector(double eigenvalue,
double[] d,
double[] l)
Find an eigenvector corresponding to an eigenvalue, using bidiagonals. |
private void |
findEigenVectors()
Find eigenvectors. |
private boolean |
flipIfWarranted(int n,
int step)
Flip qd array if warranted. |
RealMatrix |
getD()
Returns the block diagonal matrix D of the decomposition. |
double |
getDeterminant()
Return the determinant of the matrix |
RealVector |
getEigenvector(int i)
Returns a copy of the ith eigenvector of the original matrix. |
double |
getImagEigenvalue(int i)
Returns the imaginary part of the ith eigenvalue of the original matrix. |
double[] |
getImagEigenvalues()
Returns a copy of the imaginary parts of the eigenvalues of the original matrix. |
double |
getRealEigenvalue(int i)
Returns the real part of the ith eigenvalue of the original matrix. |
double[] |
getRealEigenvalues()
Returns a copy of the real parts of the eigenvalues of the original matrix. |
DecompositionSolver |
getSolver()
Get a solver for finding the A × X = B solution in exact linear sense. |
RealMatrix |
getV()
Returns the matrix V of the decomposition. |
RealMatrix |
getVT()
Returns the transpose of the matrix V of the decomposition. |
private int |
goodStep(int start,
int end)
Perform one "good" dqd/dqds step. |
private void |
initialSplits(int n)
Perform two iterations with Li's tests for initial splits. |
private boolean |
isSymmetric(RealMatrix matrix)
Check if a matrix is symmetric. |
private void |
ldlTDecomposition(double lambda,
int index,
int n)
Decompose the shifted tridiagonal matrix T-λI as LDLT. |
private void |
process1RowBlock(int index)
Find eigenvalue in a block with 1 row. |
private void |
process2RowsBlock(int index)
Find realEigenvalues in a block with 2 rows. |
private void |
process3RowsBlock(int index)
Find realEigenvalues in a block with 3 rows. |
private void |
processGeneralBlock(int n)
Find realEigenvalues using dqd/dqds algorithms. |
private void |
progressiveQuotientDifferenceWithShift(double[] d,
double[] l,
double lambda)
Decompose matrix LDLT - λ I as U-D-U-T. |
private void |
stationaryQuotientDifferenceWithShift(double[] d,
double[] l,
double lambda)
Decompose matrix LDLT - λ I as L+D+L+T. |
private void |
transformToTridiagonal(RealMatrix matrix)
Transform matrix to tridiagonal. |
private void |
updateSigma(double tau)
Update sigma. |
Methods inherited from class java.lang.Object |
---|
clone, equals, finalize, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait |
Field Detail |
---|
private static final double TOLERANCE
private static final double TOLERANCE_2
private double splitTolerance
private double[] main
private double[] secondary
private double[] squaredSecondary
private TriDiagonalTransformer transformer
private double lowerSpectra
private double upperSpectra
private double minPivot
private double sigma
private double sigmaLow
private double tau
private double[] work
private int pingPong
private double qMax
private double eMin
private int tType
private double dMin
private double dMin1
private double dMin2
private double dN
private double dN1
private double dN2
private double g
private double[] realEigenvalues
private double[] imagEigenvalues
private ArrayRealVector[] eigenvectors
private RealMatrix cachedV
private RealMatrix cachedD
private RealMatrix cachedVt
Constructor Detail |
---|
public EigenDecompositionImpl(RealMatrix matrix, double splitTolerance) throws InvalidMatrixException
matrix
- The symmetric matrix to decompose.splitTolerance
- tolerance on the off-diagonal elements relative to the
geometric mean to split the tridiagonal matrix (a suggested value is
MathUtils.SAFE_MIN
)
InvalidMatrixException
- (wrapping a ConvergenceException
if algorithm fails to convergepublic EigenDecompositionImpl(double[] main, double[] secondary, double splitTolerance) throws InvalidMatrixException
main
- the main diagonal of the matrix (will be copied)secondary
- the secondary diagonal of the matrix (will be copied)splitTolerance
- tolerance on the off-diagonal elements relative to the
geometric mean to split the tridiagonal matrix (a suggested value is
MathUtils.SAFE_MIN
)
InvalidMatrixException
- (wrapping a ConvergenceException
if algorithm fails to convergeMethod Detail |
---|
private boolean isSymmetric(RealMatrix matrix)
matrix
- matrix to check
private void decompose()
InvalidMatrixException
- (wrapping a ConvergenceException
if algorithm fails to convergepublic RealMatrix getV() throws InvalidMatrixException
V is an orthogonal matrix, i.e. its transpose is also its inverse.
The columns of V are the eigenvectors of the original matrix.
getV
in interface EigenDecomposition
InvalidMatrixException
public RealMatrix getD() throws InvalidMatrixException
D is a block diagonal matrix.
Real eigenvalues are on the diagonal while complex values are on 2x2 blocks { {real +imaginary}, {-imaginary, real} }.
getD
in interface EigenDecomposition
InvalidMatrixException
EigenDecomposition.getRealEigenvalues()
,
EigenDecomposition.getImagEigenvalues()
public RealMatrix getVT() throws InvalidMatrixException
V is an orthogonal matrix, i.e. its transpose is also its inverse.
The columns of V are the eigenvectors of the original matrix.
getVT
in interface EigenDecomposition
InvalidMatrixException
public double[] getRealEigenvalues() throws InvalidMatrixException
getRealEigenvalues
in interface EigenDecomposition
InvalidMatrixException
EigenDecomposition.getD()
,
EigenDecomposition.getRealEigenvalue(int)
,
EigenDecomposition.getImagEigenvalues()
public double getRealEigenvalue(int i) throws InvalidMatrixException, ArrayIndexOutOfBoundsException
getRealEigenvalue
in interface EigenDecomposition
i
- index of the eigenvalue (counting from 0)
InvalidMatrixException
ArrayIndexOutOfBoundsException
EigenDecomposition.getD()
,
EigenDecomposition.getRealEigenvalues()
,
EigenDecomposition.getImagEigenvalue(int)
public double[] getImagEigenvalues() throws InvalidMatrixException
getImagEigenvalues
in interface EigenDecomposition
InvalidMatrixException
EigenDecomposition.getD()
,
EigenDecomposition.getImagEigenvalue(int)
,
EigenDecomposition.getRealEigenvalues()
public double getImagEigenvalue(int i) throws InvalidMatrixException, ArrayIndexOutOfBoundsException
getImagEigenvalue
in interface EigenDecomposition
i
- index of the eigenvalue (counting from 0)
InvalidMatrixException
ArrayIndexOutOfBoundsException
EigenDecomposition.getD()
,
EigenDecomposition.getImagEigenvalues()
,
EigenDecomposition.getRealEigenvalue(int)
public RealVector getEigenvector(int i) throws InvalidMatrixException, ArrayIndexOutOfBoundsException
getEigenvector
in interface EigenDecomposition
i
- index of the eigenvector (counting from 0)
InvalidMatrixException
ArrayIndexOutOfBoundsException
EigenDecomposition.getD()
public double getDeterminant()
getDeterminant
in interface EigenDecomposition
public DecompositionSolver getSolver()
getSolver
in interface EigenDecomposition
private void transformToTridiagonal(RealMatrix matrix)
matrix
- matrix to transformprivate void computeGershgorinCircles()
private void findEigenvalues() throws InvalidMatrixException
InvalidMatrixException
- if a block cannot be diagonalizedprivate List<Integer> computeSplits()
private void process1RowBlock(int index)
In low dimensions, we simply solve the characteristic polynomial.
index
- index of the first row of the blockprivate void process2RowsBlock(int index) throws InvalidMatrixException
In low dimensions, we simply solve the characteristic polynomial.
index
- index of the first row of the block
InvalidMatrixException
- if characteristic polynomial cannot be solvedprivate void process3RowsBlock(int index) throws InvalidMatrixException
In low dimensions, we simply solve the characteristic polynomial.
index
- index of the first row of the block
InvalidMatrixException
- if diagonal elements are not positiveprivate void processGeneralBlock(int n) throws InvalidMatrixException
This implementation is based on Beresford N. Parlett and Osni A. Marques paper An Implementation of the dqds Algorithm (Positive Case) and on the corresponding LAPACK routine DLASQ2.
n
- number of rows of the block
InvalidMatrixException
- if block cannot be diagonalized
after 30 * n iterationsprivate void initialSplits(int n)
n
- number of rows of the matrix to processprivate int goodStep(int start, int end)
This implementation is based on Beresford N. Parlett and Osni A. Marques paper An Implementation of the dqds Algorithm (Positive Case) and on the corresponding LAPACK routine DLAZQ3.
start
- start indexend
- end index
private boolean flipIfWarranted(int n, int step)
n
- number of rows in the blockstep
- within the array (1 for flipping all elements, 2 for flipping
only every other element)
private double[] eigenvaluesRange(int index, int n)
index
- index of the first row of the blockn
- number of rows of the block
private int countEigenValues(double t, int index, int n)
t
- value below which we must count the number of realEigenvaluesindex
- index of the first row of the blockn
- number of rows of the block
private void ldlTDecomposition(double lambda, int index, int n)
A shifted symmetric tridiagonal matrix T can be decomposed as LDLT where L is a lower bidiagonal matrix with unit diagonal and D is a diagonal matrix. This method is an implementation of algorithm 4.4.7 from Dhillon's thesis.
lambda
- shift to add to the matrix before decomposing it
to ensure it is positive definiteindex
- index of the first row of the blockn
- number of rows of the blockprivate void dqds(int start, int end)
This implementation is a translation of the LAPACK routine DLASQ5.
start
- start indexend
- end indexprivate void dqd(int start, int end)
This implementation is a translation of the LAPACK routine DLASQ6.
start
- start indexend
- end indexprivate void computeShiftIncrement(int start, int end, int deflated)
This implementation is a translation of the LAPACK routine DLAZQ4.
start
- start indexend
- end indexdeflated
- number of realEigenvalues just deflatedprivate void updateSigma(double tau)
tau
- shift to apply to sigmaprivate void findEigenVectors()
private ArrayRealVector findEigenvector(double eigenvalue, double[] d, double[] l)
This method corresponds to algorithm X from Dhillon's thesis.
eigenvalue
- eigenvalue for which eigenvector is desiredd
- diagonal elements of the initial non-shifted D matrixl
- off-diagonal elements of the initial non-shifted L matrix
private void stationaryQuotientDifferenceWithShift(double[] d, double[] l, double lambda)
This method corresponds to algorithm 4.4.3 (dstqds) from Dhillon's thesis.
d
- diagonal elements of D,l
- off-diagonal elements of Llambda
- shift to applyprivate void progressiveQuotientDifferenceWithShift(double[] d, double[] l, double lambda)
This method corresponds to algorithm 4.4.5 (dqds) from Dhillon's thesis.
d
- diagonal elements of Dl
- off-diagonal elements of Llambda
- shift to apply
|
|||||||||
PREV CLASS NEXT CLASS | FRAMES NO FRAMES | ||||||||
SUMMARY: NESTED | FIELD | CONSTR | METHOD | DETAIL: FIELD | CONSTR | METHOD |