Class SymmLQ.State
- java.lang.Object
-
- org.apache.commons.math3.linear.SymmLQ.State
-
- Enclosing class:
- SymmLQ
private static class SymmLQ.State extends java.lang.Object
A simple container holding the non-final variables used in the iterations. Making the current state of the solver visible from the outside is necessary, because during the iterations,
x
does not exactly hold the current estimate of the solution. Indeed,x
needs in general to be moved from the LQ point to the CG point. Besides, additional upudates must be carried out in casegoodb
is set totrue
.In all subsequent comments, the description of the state variables refer to their value after a call to
update()
. In these comments, k is the current number of evaluations of matrix-vector products.
-
-
Field Summary
Fields Modifier and Type Field Description private RealLinearOperator
a
Reference to the linear operator.private RealVector
b
Reference to the right-hand side vector.private double
beta
The value of beta[k+1].private double
beta1
The value of beta[1].private boolean
bIsNull
The value ofb == 0
(exact floating-point equality).private double
bstep
The value of bstep[k-1].(package private) static double
CBRT_MACH_PREC
The cubic root ofMACH_PREC
.private double
cgnorm
The estimate of the norm of P * rC[k].private boolean
check
true
if symmetry of matrix and conditioner must be checked.private double
dbar
The value of dbar[k+1] = -beta[k+1] * c[k-1].private double
delta
The value of the custom tolerance δ for the default stopping criterion.private double
gammaZeta
The value of gamma[k] * zeta[k].private double
gbar
The value of gbar[k].private double
gmax
The value of max(|alpha[1]|, gamma[1], ..., gamma[k-1]).private double
gmin
The value of min(|alpha[1]|, gamma[1], ..., gamma[k-1]).private boolean
goodb
Copy of thegoodb
parameter.private boolean
hasConverged
true
if the default convergence criterion is verified.private double
lqnorm
The estimate of the norm of P * rL[k-1].private RealLinearOperator
m
Reference to the preconditioner, M.(package private) static double
MACH_PREC
The machine precision.private RealVector
mb
The value of M * b.private double
minusEpsZeta
The value of (-eps[k+1] * zeta[k-1]).private double
oldb
The value of beta[k].private RealVector
r1
The value of beta[k] * M^(-1) * P' * v[k].private RealVector
r2
The value of beta[k+1] * M^(-1) * P' * v[k+1].private double
rnorm
The value of the updated, preconditioned residual P * r.private double
shift
Copy of theshift
parameter.private double
snprod
The value of s[1] * ...private double
tnorm
An estimate of the square of the norm of A * V[k], based on Paige and Saunders (1975), equation (3.3).private RealVector
wbar
The value of P' * wbar[k] or P' * (wbar[k] - s[1] * ...private RealVector
xL
A reference to the vector to be updated with the solution.private RealVector
y
The value of beta[k+1] * P' * v[k+1].private double
ynorm2
The value of zeta[1]^2 + ...
-
Constructor Summary
Constructors Constructor Description State(RealLinearOperator a, RealLinearOperator m, RealVector b, boolean goodb, double shift, double delta, boolean check)
Creates and inits to k = 1 a new instance of this class.
-
Method Summary
All Methods Static Methods Instance Methods Concrete Methods Modifier and Type Method Description (package private) boolean
bEqualsNullVector()
Returnstrue
if the right-hand side vector is zero exactly.(package private) boolean
betaEqualsZero()
Returnstrue
ifbeta
is essentially zero.private static void
checkSymmetry(RealLinearOperator l, RealVector x, RealVector y, RealVector z)
Performs a symmetry check on the specified linear operator, and throws an exception in case this check fails.private static void
daxpbypz(double a, RealVector x, double b, RealVector y, RealVector z)
A BLAS-like function, for the operation z ← a · x + b · y + z.private static void
daxpy(double a, RealVector x, RealVector y)
A clone of the BLASDAXPY
function, which carries out the operation y ← a · x + y.(package private) double
getNormOfResidual()
Returns the norm of the updated, preconditioned residual.(package private) boolean
hasConverged()
Returnstrue
if the default stopping criterion is fulfilled.(package private) void
init()
Performs the initial phase of the SYMMLQ algorithm.(package private) void
refineSolution(RealVector x)
Move to the CG point if it seems better.private static void
throwNPDLOException(RealLinearOperator l, RealVector v)
Throws a newNonPositiveDefiniteOperatorException
with appropriate context.(package private) void
update()
Performs the next iteration of the algorithm.private void
updateNorms()
Computes the norms of the residuals, and checks for convergence.
-
-
-
Field Detail
-
CBRT_MACH_PREC
static final double CBRT_MACH_PREC
The cubic root ofMACH_PREC
.
-
MACH_PREC
static final double MACH_PREC
The machine precision.
-
a
private final RealLinearOperator a
Reference to the linear operator.
-
b
private final RealVector b
Reference to the right-hand side vector.
-
check
private final boolean check
true
if symmetry of matrix and conditioner must be checked.
-
delta
private final double delta
The value of the custom tolerance δ for the default stopping criterion.
-
beta
private double beta
The value of beta[k+1].
-
beta1
private double beta1
The value of beta[1].
-
bstep
private double bstep
The value of bstep[k-1].
-
cgnorm
private double cgnorm
The estimate of the norm of P * rC[k].
-
dbar
private double dbar
The value of dbar[k+1] = -beta[k+1] * c[k-1].
-
gammaZeta
private double gammaZeta
The value of gamma[k] * zeta[k]. Was calledrhs1
in the initial code.
-
gbar
private double gbar
The value of gbar[k].
-
gmax
private double gmax
The value of max(|alpha[1]|, gamma[1], ..., gamma[k-1]).
-
gmin
private double gmin
The value of min(|alpha[1]|, gamma[1], ..., gamma[k-1]).
-
goodb
private final boolean goodb
Copy of thegoodb
parameter.
-
hasConverged
private boolean hasConverged
true
if the default convergence criterion is verified.
-
lqnorm
private double lqnorm
The estimate of the norm of P * rL[k-1].
-
m
private final RealLinearOperator m
Reference to the preconditioner, M.
-
minusEpsZeta
private double minusEpsZeta
The value of (-eps[k+1] * zeta[k-1]). Was calledrhs2
in the initial code.
-
mb
private final RealVector mb
The value of M * b.
-
oldb
private double oldb
The value of beta[k].
-
r1
private RealVector r1
The value of beta[k] * M^(-1) * P' * v[k].
-
r2
private RealVector r2
The value of beta[k+1] * M^(-1) * P' * v[k+1].
-
rnorm
private double rnorm
-
shift
private final double shift
Copy of theshift
parameter.
-
snprod
private double snprod
The value of s[1] * ... * s[k-1].
-
tnorm
private double tnorm
An estimate of the square of the norm of A * V[k], based on Paige and Saunders (1975), equation (3.3).
-
wbar
private RealVector wbar
The value of P' * wbar[k] or P' * (wbar[k] - s[1] * ... * s[k-1] * v[1]) ifgoodb
istrue
. Was calledw
in the initial code.
-
xL
private final RealVector xL
A reference to the vector to be updated with the solution. Contains the value of xL[k-1] ifgoodb
isfalse
, (xL[k-1] - bstep[k-1] * v[1]) otherwise.
-
y
private RealVector y
The value of beta[k+1] * P' * v[k+1].
-
ynorm2
private double ynorm2
The value of zeta[1]^2 + ... + zeta[k-1]^2.
-
bIsNull
private boolean bIsNull
The value ofb == 0
(exact floating-point equality).
-
-
Constructor Detail
-
State
State(RealLinearOperator a, RealLinearOperator m, RealVector b, boolean goodb, double shift, double delta, boolean check)
Creates and inits to k = 1 a new instance of this class.- Parameters:
a
- the linear operator A of the systemm
- the preconditioner, M (can benull
)b
- the right-hand side vectorgoodb
- usuallyfalse
, except ifx
is expected to contain a large multiple ofb
shift
- the amount to be subtracted to all diagonal elements of Adelta
- the δ parameter for the default stopping criterioncheck
-true
if self-adjointedness of both matrix and preconditioner should be checked
-
-
Method Detail
-
checkSymmetry
private static void checkSymmetry(RealLinearOperator l, RealVector x, RealVector y, RealVector z) throws NonSelfAdjointOperatorException
Performs a symmetry check on the specified linear operator, and throws an exception in case this check fails. Given a linear operator L, and a vector x, this method checks that x' · L · y = y' · L · x (within a given accuracy), where y = L · x.- Parameters:
l
- the linear operator Lx
- the candidate vector xy
- the candidate vector y = L · xz
- the vector z = L · y- Throws:
NonSelfAdjointOperatorException
- when the test fails
-
throwNPDLOException
private static void throwNPDLOException(RealLinearOperator l, RealVector v) throws NonPositiveDefiniteOperatorException
Throws a newNonPositiveDefiniteOperatorException
with appropriate context.- Parameters:
l
- the offending linear operatorv
- the offending vector- Throws:
NonPositiveDefiniteOperatorException
- in any circumstances
-
daxpy
private static void daxpy(double a, RealVector x, RealVector y)
A clone of the BLASDAXPY
function, which carries out the operation y ← a · x + y. This is for internal use only: no dimension checks are provided.- Parameters:
a
- the scalar by whichx
is to be multipliedx
- the vector to be added toy
y
- the vector to be incremented
-
daxpbypz
private static void daxpbypz(double a, RealVector x, double b, RealVector y, RealVector z)
A BLAS-like function, for the operation z ← a · x + b · y + z. This is for internal use only: no dimension checks are provided.- Parameters:
a
- the scalar by whichx
is to be multipliedx
- the first vector to be added toz
b
- the scalar by whichy
is to be multipliedy
- the second vector to be added toz
z
- the vector to be incremented
-
refineSolution
void refineSolution(RealVector x)
Move to the CG point if it seems better. In this version of SYMMLQ, the convergence tests involve only cgnorm, so we're unlikely to stop at an LQ point, except if the iteration limit interferes.
Additional upudates are also carried out in case
goodb
is set totrue
.- Parameters:
x
- the vector to be updated with the refined value of xL
-
init
void init()
Performs the initial phase of the SYMMLQ algorithm. On return, the value of the state variables ofthis
object correspond to k = 1.
-
update
void update()
Performs the next iteration of the algorithm. The iteration count should be incremented prior to calling this method. On return, the value of the state variables ofthis
object correspond to the current iteration countk
.
-
updateNorms
private void updateNorms()
-
hasConverged
boolean hasConverged()
Returnstrue
if the default stopping criterion is fulfilled.- Returns:
true
if convergence of the iterations has occurred
-
bEqualsNullVector
boolean bEqualsNullVector()
Returnstrue
if the right-hand side vector is zero exactly.- Returns:
- the boolean value of
b == 0
-
betaEqualsZero
boolean betaEqualsZero()
Returnstrue
ifbeta
is essentially zero. This method is used to check for early stop of the iterations.- Returns:
true
ifbeta <
MACH_PREC
-
getNormOfResidual
double getNormOfResidual()
Returns the norm of the updated, preconditioned residual.- Returns:
- the norm of the residual, ||P * r||
-
-