Class RealSchur
- All Implemented Interfaces:
Decomposition<Matrix>
This class computes the Schur decomposition of a real dense square matrix.
That is, decompose a square matrix A
into A=UTU
T where U
is an orthogonal
matrix and T
is a block-upper triangular matrix called the real-Schur form of A
. T
is upper triangular
except for possibly 2x2 blocks along the diagonal. T
is similar to A
.
This code was adapted from the code found in the EJML library and the description of the Francis implicit double shifted QR algorithm given in Fundamentals of Matrix Computations 3rd Edition by David S. Watkins.
-
Field Summary
Modifier and TypeFieldDescriptionprotected double
Stores the scalar factor invalid input: '&alpha' for use in computation of the Householder reflectorP = I -
invalid input: '&alpha'vv
T .protected double
For computing the norm of a column for use when computing Householder reflectors.Fields inherited from class org.flag4j.linalg.decompositions.schur.Schur
checkFinite, computeU, DEFAULT_EXCEPTIONAL_ITERS, DEFAULT_MAX_ITERS_FACTOR, exceptionalThreshold, hess, householderVector, maxIterations, maxIterationsFactor, numExceptional, numRows, rng, shiftCol, sinceLastExceptional, T, temp, U, workArray
-
Constructor Summary
ConstructorDescriptionCreates a decomposer to compute the Schur decomposition for a real dense matrix.RealSchur
(boolean computeU) Creates a decomposer to compute the Schur decomposition for a real dense matrix where theU
matrix may or may not be computed.RealSchur
(boolean computeU, long seed) Creates a decomposer to compute the Schur decomposition for a real dense matrix.RealSchur
(long seed) Creates a decomposer to compute the Schur decomposition for a real dense matrix. -
Method Summary
Modifier and TypeMethodDescriptionprotected void
applyDoubleShiftReflector
(int i, boolean set) Applies reflector for the double shift.protected void
applyReflector
(int i, int shiftSize) Applies the constructed Householder reflector which has been constructed for the given shift size.protected void
applySingleShiftReflector
(int i, boolean set) Applies reflector for the double shift.protected int
checkConvergence
(int workingSize) Checks for convergence of lower 2x2 sub-matrix within working matrix to upper triangular or block upper triangular form.protected double
computeExceptionalShift
(int k) Computes a random shift to help theQR
algorithm converge if it gets stuck.protected void
computeImplicitDoubleShift
(int workingSize) Computes the shifts for a Francis double shift iteration.protected void
computeImplicitSingleShift
(int k, double shift) Computes the non-zero entries of the first column for the single shiftedQR
algorithm.Computes the Schur decomposition of the input matrix.protected boolean
makeReflector
(int i, double p1, double p2) Constructs a householder reflector given specified values for a column to apply the reflector to.protected boolean
makeReflector
(int i, double p1, double p2, double p3) Constructs a householder reflector given specified values for a column to apply the reflector to.protected void
performDoubleShift
(int workingSize) Performs a full iteration of the Francis implicit double shiftedQR
algorithm (this includes the bulge chase).protected void
performExceptionalShift
(int workingSize) Performs a full iteration of the single shiftedQR
algorithm (this includes the bulge chase) where the shift is chosen to be a random value with the same magnitude as the lower right element of the working matrix.protected void
performSingleShift
(int workingSize, double shift) Performs a full iteration of the implicit single shiftedQR
algorithm (this includes the bulge chase).CMatrix[]
Converts the real schur form computed in the last decomposition to the complex Schur form.setExceptionalThreshold
(int exceptionalThreshold) Sets the number of iterations of theQR
algorithm to perform without deflation before performing a random shift.setMaxIterationFactor
(int maxIterationFactor) Specify maximum iteration factor for computing the total number of iterations to run theQR
algorithm for when computing the decomposition.protected void
Initializes temporary work arrays to be used in the decomposition.Methods inherited from class org.flag4j.linalg.decompositions.schur.Schur
decomposeBase, getT, getU, setUp
-
Field Details
-
norm
protected double normFor computing the norm of a column for use when computing Householder reflectors. -
currentFactor
protected double currentFactorStores the scalar factor invalid input: '&alpha' for use in computation of the Householder reflectorP = I -
invalid input: '&alpha'vv
T .
-
-
Constructor Details
-
RealSchur
public RealSchur()Creates a decomposer to compute the Schur decomposition for a real dense matrix.
Note: This decomposer may use random numbers during the decomposition. If reproducible results are needed, set the seed for the pseudo-random number generator using
RealSchur(long)
-
RealSchur
public RealSchur(boolean computeU) Creates a decomposer to compute the Schur decomposition for a real dense matrix where the
U
matrix may or may not be computed.If the
U
matrix is not needed, passingcomputeU = false
may provide a performance improvement.By default if a constructor with no
computeU
parameter is called,U
WILL be computed.Note: This decomposer may use random numbers during the decomposition. If reproducible results are needed, set the seed for the pseudo-random number generator using
RealSchur(boolean, long)
- Parameters:
computeU
- Flag indicating if the unitaryU
matrix should be computed for the Schur decomposition. If true,U
will be computed. If false,U
will not be computed.
-
RealSchur
public RealSchur(long seed) Creates a decomposer to compute the Schur decomposition for a real dense matrix.- Parameters:
seed
- Seed to use for pseudo-random number generator when computing exceptional shifts during theQR
algorithm.
-
RealSchur
public RealSchur(boolean computeU, long seed) Creates a decomposer to compute the Schur decomposition for a real dense matrix.- Parameters:
seed
- Seed to use for pseudo-random number generator when computing exceptional shifts during theQR
algorithm.
-
-
Method Details
-
setExceptionalThreshold
Sets the number of iterations of the
QR
algorithm to perform without deflation before performing a random shift.That is, if
exceptionalThreshold = 10
, then at most 10 iterationsQR
algorithm iterations will be performed. If, by the 10th iteration, no convergence has been detected which allows for deflation, then aQR
algorithm iteration will be performed with a random (i.e. exceptional) shift.By default, the threshold is set to
Schur.DEFAULT_EXCEPTIONAL_ITERS
- Overrides:
setExceptionalThreshold
in classSchur<Matrix,
double[]> - Parameters:
exceptionalThreshold
- The new exceptional shift threshold. i.e. the number of iterations to perform without deflation before performing an iteration with random shifts.- Returns:
- A reference to this decomposer.
- Throws:
IllegalArgumentException
- IfexceptionalThreshold
is not positive.
-
setMaxIterationFactor
Specify maximum iteration factor for computing the total number of iterations to run the
QR
algorithm for when computing the decomposition. The maximum number of iterations is computed asmaxIteration = maxIterationFactor * src.numRows;
By default, this is computed as
maxIterations =
Schur.DEFAULT_MAX_ITERS_FACTOR
* src.numRows;
src
is the matrix being decomposed.- Overrides:
setMaxIterationFactor
in classSchur<Matrix,
double[]> - Parameters:
maxIterationFactor
- maximum iteration factor for use in computing the total maximum number of iterations to run theQR
algorithm for.- Returns:
- A reference to this decomposer.
- Throws:
IllegalArgumentException
- IfmaxIterationFactor
is not positive.
-
decompose
-
setUpArrays
protected void setUpArrays()Initializes temporary work arrays to be used in the decomposition.- Specified by:
setUpArrays
in classSchur<Matrix,
double[]>
-
performExceptionalShift
protected void performExceptionalShift(int workingSize) Performs a full iteration of the single shiftedQR
algorithm (this includes the bulge chase) where the shift is chosen to be a random value with the same magnitude as the lower right element of the working matrix. This can help theQR
converge for certain pathological cases where the double shift algorithm oscillates or fails to converge for repeated eigenvalues.- Specified by:
performExceptionalShift
in classSchur<Matrix,
double[]> - Parameters:
workingSize
- The current working size for the decomposition. I.e. all entries below this row have converged to an upper or possible 2x2 block upper triangular form.
-
computeExceptionalShift
protected double computeExceptionalShift(int k) Computes a random shift to help theQR
algorithm converge if it gets stuck.- Parameters:
k
- The current size of the working matrix. Specifically, the index of the lower right value in the working matrix is(k, k)
.- Returns:
- A shift in a random direction which has the same magnitude as the elements in the matrix.
-
computeImplicitSingleShift
protected void computeImplicitSingleShift(int k, double shift) Computes the non-zero entries of the first column for the single shiftedQR
algorithm.- Parameters:
k
- Size of current working matrix.shift
- The shift to use.
-
performSingleShift
protected void performSingleShift(int workingSize, double shift) Performs a full iteration of the implicit single shiftedQR
algorithm (this includes the bulge chase).- Parameters:
workingSize
- The current working size for the decomposition. I.e. all entries below this row have converged to an upper or possible 2x2 block upper triangular form.shift
- The shift to use in the implicit single shiftedQR
algorithm.
-
applySingleShiftReflector
protected void applySingleShiftReflector(int i, boolean set) Applies reflector for the double shift. This method can be used to apply either be the reflector constructed for the first column of the shifted matrix, or a reflector being used in the bulge chase of size 2 which arises from the first case.- Parameters:
i
- The starting row the reflector is being applied to.
-
performDoubleShift
protected void performDoubleShift(int workingSize) Performs a full iteration of the Francis implicit double shiftedQR
algorithm (this includes the bulge chase).- Specified by:
performDoubleShift
in classSchur<Matrix,
double[]> - Parameters:
workingSize
- The current working size for the decomposition. I.e. all entries below this row have converged to an upper or possible 2x2 block upper triangular form.
-
computeImplicitDoubleShift
protected void computeImplicitDoubleShift(int workingSize) Computes the shifts for a Francis double shift iteration. Specifically, the shifts are the generalized Rayleigh quotients of degree two.- Parameters:
workingSize
- Size of current working matrix.
-
applyDoubleShiftReflector
protected void applyDoubleShiftReflector(int i, boolean set) Applies reflector for the double shift. This method can be used to apply either be the reflector constructed for the first column of the shifted matrix, or a reflector being used in the bulge chase of size 2 which arises from the first case.- Parameters:
i
- The starting row the reflector is being applied to.
-
applyReflector
protected void applyReflector(int i, int shiftSize) Applies the constructed Householder reflector which has been constructed for the given shift size.- Parameters:
i
- The stating row the reflector is being applied to.shiftSize
- The size of the shift the reflector was constructed for.
-
makeReflector
protected boolean makeReflector(int i, double p1, double p2, double p3) Constructs a householder reflector given specified values for a column to apply the reflector to. This reflector is stored in indicesi
,i+1
, andi+2
ofSchur.householderVector
.- Parameters:
i
- Row of working matrix to construct reflector for.p1
- First entry to in column to apply reflector to.p2
- Second entry in column to apply reflector to.p3
- Third entry in column to apply reflector to.- Returns:
- True if a reflector needs to be constructed to return matrix to upper Hessenburg form. False if column is already in the correct form.
-
makeReflector
protected boolean makeReflector(int i, double p1, double p2) Constructs a householder reflector given specified values for a column to apply the reflector to. This reflector is stored in indicesi
andi+1
ofSchur.householderVector
.- Parameters:
i
- Row of working matrix to construct reflector for.p1
- First entry to in column to apply reflector to.p2
- Second entry in column to apply reflector to.- Returns:
- True if a reflector needs to be constructed to return matrix to upper Hessenburg form. False if column is already in the correct form.
-
checkConvergence
protected int checkConvergence(int workingSize) Checks for convergence of lower 2x2 sub-matrix within working matrix to upper triangular or block upper triangular form. If convergence is found, this will also zero out the values which have converged to near zero.- Specified by:
checkConvergence
in classSchur<Matrix,
double[]> - Parameters:
workingSize
- Size of current working matrix.- Returns:
- Returns the amount the working matrix size should be deflated. Will be zero if no convergence is detected, one if convergence to upper triangular form is detected and two if convergence to block upper triangular form is detected.
-
real2ComplexSchur
Converts the real schur form computed in the last decomposition to the complex Schur form.
That is, converts the real block upper triangular Schur matrix to a complex valued properly upper triangular matrix. If the unitary transformation matrix
U
was computed, the transformations will also be updated accordingly.This method was adapted from the code given by scipy.linalg.rsf2csf (v1.12.0).
- Returns:
- An array of length 2 containing the complex Schur matrix
T
from the last decomposition, and if computed, the complex unitary transformation matrixU
from the decomposition. IfU
was not computed, then the arrays second value will be null.
-