Class SymmLQ
 java.lang.Object

 org.apache.commons.math4.legacy.linear.IterativeLinearSolver

 org.apache.commons.math4.legacy.linear.PreconditionedIterativeLinearSolver

 org.apache.commons.math4.legacy.linear.SymmLQ

public class SymmLQ extends PreconditionedIterativeLinearSolver
Implementation of the SYMMLQ iterative linear solver proposed by Paige and Saunders (1975). This implementation is largely based on the FORTRAN code by Pr. Michael A. Saunders, available here.
SYMMLQ is designed to solve the system of linear equations A · x = b where A is an n × n selfadjoint linear operator (defined as a
RealLinearOperator
), and b is a given vector. The operator A is not required to be positive definite. If A is known to be definite, the method of conjugate gradients might be preferred, since it will require about the same number of iterations as SYMMLQ but slightly less work per iteration.SYMMLQ is designed to solve the system (A  shift · I) · x = b, where shift is a specified scalar value. If shift and b are suitably chosen, the computed vector x may approximate an (unnormalized) eigenvector of A, as in the methods of inverse iteration and/or Rayleighquotient iteration. Again, the linear operator (A  shift · I) need not be positive definite (but must be selfadjoint). The work per iteration is very slightly less if shift = 0.
Preconditioning
Preconditioning may reduce the number of iterations required. The solver may be provided with a positive definite preconditioner M = P^{T} · P that is known to approximate (A  shift · I)^{1} in some sense, where matrixvector products of the form M · y = x can be computed efficiently. Then SYMMLQ will implicitly solve the system of equations P · (A  shift · I) · P^{T} · x_{hat} = P · b, i.e. A_{hat} · x_{hat} = b_{hat}, where A_{hat} = P · (A  shift · I) · P^{T}, b_{hat} = P · b, and return the solution x = P^{T} · x_{hat}. The associated residual is r_{hat} = b_{hat}  A_{hat} · x_{hat} = P · [b  (A  shift · I) · x] = P · r.
In the case of preconditioning, the
IterativeLinearSolverEvent
s that this solver fires are such thatIterativeLinearSolverEvent.getNormOfResidual()
returns the norm of the preconditioned, updated residual, P · r, not the norm of the true residual r.A default stopping criterion is implemented. The iterations stop when  rhat  ≤ δ  Ahat   xhat , where xhat is the current estimate of the solution of the transformed system, rhat the current estimate of the corresponding residual, and δ a userspecified tolerance.
Iteration count
In the present context, an iteration should be understood as one evaluation of the matrixvector product A · x. The initialization phase therefore counts as one iteration. If the user requires checks on the symmetry of A, this entails one further matrixvector product in the initial phase. This further product is not accounted for in the iteration count. In other words, the number of iterations required to reach convergence will be identical, whether checks have been required or not.
The present definition of the iteration count differs from that adopted in the original FOTRAN code, where the initialization phase was not taken into account.
The
x
parameter insolve(RealLinearOperator, RealVector, RealVector)
,solve(RealLinearOperator, RealLinearOperator, RealVector, RealVector)
},solveInPlace(RealLinearOperator, RealVector, RealVector)
,solveInPlace(RealLinearOperator, RealLinearOperator, RealVector, RealVector)
,solveInPlace(RealLinearOperator, RealLinearOperator, RealVector, RealVector, boolean, double)
,
Besides standard
DimensionMismatchException
, this class might throwNonSelfAdjointOperatorException
if the linear operator or the preconditioner are not symmetric. In this case, theExceptionContext
provides more information key
"operator"
points to the offending linear operator, say L,  key
"vector1"
points to the first offending vector, say x,  key
"vector2"
points to the second offending vector, say y, such that x^{T} · L · y ≠ y^{T} · L · x (within a certain accuracy).
NonPositiveDefiniteOperatorException
might also be thrown in case the preconditioner is not positive definite. The relevant keys to theExceptionContext
are key
"operator"
, which points to the offending linear operator, say L,  key
"vector"
, which points to the offending vector, say x, such that x^{T} · L · x < 0.
References
 Paige and Saunders (1975)
 C. C. Paige and M. A. Saunders, Solution of Sparse Indefinite Systems of Linear Equations, SIAM Journal on Numerical Analysis 12(4): 617629, 1975
 Since:
 3.0


Constructor Summary
Constructors Constructor Description SymmLQ(int maxIterations, double delta, boolean check)
Creates a new instance of this class, with default stopping criterion.SymmLQ(IterationManager manager, double delta, boolean check)
Creates a new instance of this class, with default stopping criterion and custom iteration manager.

Method Summary
All Methods Instance Methods Concrete Methods Modifier and Type Method Description boolean
getCheck()
Returnstrue
if symmetry of the matrix, and symmetry as well as positive definiteness of the preconditioner should be checked.RealVector
solve(RealLinearOperator a, RealLinearOperator m, RealVector b)
Returns an estimate of the solution to the linear system A · x = b.RealVector
solve(RealLinearOperator a, RealLinearOperator m, RealVector b, boolean goodb, double shift)
Returns an estimate of the solution to the linear system (A  shift · I) · x = b.RealVector
solve(RealLinearOperator a, RealLinearOperator m, RealVector b, RealVector x)
Returns an estimate of the solution to the linear system A · x = b.RealVector
solve(RealLinearOperator a, RealVector b)
Returns an estimate of the solution to the linear system A · x = b.RealVector
solve(RealLinearOperator a, RealVector b, boolean goodb, double shift)
Returns the solution to the system (A  shift · I) · x = b.RealVector
solve(RealLinearOperator a, RealVector b, RealVector x)
Returns an estimate of the solution to the linear system A · x = b.RealVector
solveInPlace(RealLinearOperator a, RealLinearOperator m, RealVector b, RealVector x)
Returns an estimate of the solution to the linear system A · x = b.RealVector
solveInPlace(RealLinearOperator a, RealLinearOperator m, RealVector b, RealVector x, boolean goodb, double shift)
Returns an estimate of the solution to the linear system (A  shift · I) · x = b.RealVector
solveInPlace(RealLinearOperator a, RealVector b, RealVector x)
Returns an estimate of the solution to the linear system A · x = b.
Methods inherited from class org.apache.commons.math4.legacy.linear.PreconditionedIterativeLinearSolver
checkParameters

Methods inherited from class org.apache.commons.math4.legacy.linear.IterativeLinearSolver
checkParameters, getIterationManager




Constructor Detail

SymmLQ
public SymmLQ(int maxIterations, double delta, boolean check)
Creates a new instance of this class, with default stopping criterion. Note that settingcheck
totrue
entails an extra matrixvector product in the initial phase. Parameters:
maxIterations
 the maximum number of iterationsdelta
 the δ parameter for the default stopping criterioncheck
true
if selfadjointedness of both matrix and preconditioner should be checked

SymmLQ
public SymmLQ(IterationManager manager, double delta, boolean check)
Creates a new instance of this class, with default stopping criterion and custom iteration manager. Note that settingcheck
totrue
entails an extra matrixvector product in the initial phase. Parameters:
manager
 the custom iteration managerdelta
 the δ parameter for the default stopping criterioncheck
true
if selfadjointedness of both matrix and preconditioner should be checked


Method Detail

getCheck
public final boolean getCheck()
Returnstrue
if symmetry of the matrix, and symmetry as well as positive definiteness of the preconditioner should be checked. Returns:
true
if the tests are to be performed

solve
public RealVector solve(RealLinearOperator a, RealLinearOperator m, RealVector b) throws NullArgumentException, NonSquareOperatorException, DimensionMismatchException, MaxCountExceededException, NonSelfAdjointOperatorException, NonPositiveDefiniteOperatorException, IllConditionedOperatorException
Returns an estimate of the solution to the linear system A · x = b. Overrides:
solve
in classPreconditionedIterativeLinearSolver
 Parameters:
a
 the linear operator A of the systemm
 the preconditioner, M (can benull
)b
 the righthand side vector Returns:
 a new vector containing the solution
 Throws:
NonSelfAdjointOperatorException
 ifgetCheck()
istrue
, anda
orm
is not selfadjointNonPositiveDefiniteOperatorException
 ifm
is not positive definiteIllConditionedOperatorException
 ifa
is illconditionedNullArgumentException
 if one of the parameters isnull
NonSquareOperatorException
 ifa
orm
is not squareDimensionMismatchException
 ifm
orb
have dimensions inconsistent witha
MaxCountExceededException
 at exhaustion of the iteration count, unless a customcallback
has been set at construction of theIterationManager

solve
public RealVector solve(RealLinearOperator a, RealLinearOperator m, RealVector b, boolean goodb, double shift) throws NullArgumentException, NonSquareOperatorException, DimensionMismatchException, MaxCountExceededException, NonSelfAdjointOperatorException, NonPositiveDefiniteOperatorException, IllConditionedOperatorException
Returns an estimate of the solution to the linear system (A  shift · I) · x = b.If the solution x is expected to contain a large multiple of
b
(as in Rayleighquotient iteration), then better precision may be achieved withgoodb
set totrue
; this however requires an extra call to the preconditioner.shift
should be zero if the system A · x = b is to be solved. Otherwise, it could be an approximation to an eigenvalue of A, such as the Rayleigh quotient b^{T} · A · b / (b^{T} · b) corresponding to the vector b. If b is sufficiently like an eigenvector corresponding to an eigenvalue near shift, then the computed x may have very large components. When normalized, x may be closer to an eigenvector than b. Parameters:
a
 the linear operator A of the systemm
 the preconditioner, M (can benull
)b
 the righthand side vectorgoodb
 usuallyfalse
, except ifx
is expected to contain a large multiple ofb
shift
 the amount to be subtracted to all diagonal elements of A Returns:
 a reference to
x
(shallow copy)  Throws:
NullArgumentException
 if one of the parameters isnull
NonSquareOperatorException
 ifa
orm
is not squareDimensionMismatchException
 ifm
orb
have dimensions inconsistent witha
MaxCountExceededException
 at exhaustion of the iteration count, unless a customcallback
has been set at construction of theIterationManager
NonSelfAdjointOperatorException
 ifgetCheck()
istrue
, anda
orm
is not selfadjointNonPositiveDefiniteOperatorException
 ifm
is not positive definiteIllConditionedOperatorException
 ifa
is illconditioned

solve
public RealVector solve(RealLinearOperator a, RealLinearOperator m, RealVector b, RealVector x) throws NullArgumentException, NonSquareOperatorException, DimensionMismatchException, NonSelfAdjointOperatorException, NonPositiveDefiniteOperatorException, IllConditionedOperatorException, MaxCountExceededException
Returns an estimate of the solution to the linear system A · x = b. Overrides:
solve
in classPreconditionedIterativeLinearSolver
 Parameters:
x
 not meaningful in this implementation; should not be considered as an initial guess (more)a
 the linear operator A of the systemm
 the preconditioner, M (can benull
)b
 the righthand side vector Returns:
 a new vector containing the solution
 Throws:
NonSelfAdjointOperatorException
 ifgetCheck()
istrue
, anda
orm
is not selfadjointNonPositiveDefiniteOperatorException
 ifm
is not positive definiteIllConditionedOperatorException
 ifa
is illconditionedNullArgumentException
 if one of the parameters isnull
NonSquareOperatorException
 ifa
orm
is not squareDimensionMismatchException
 ifm
,b
orx0
have dimensions inconsistent witha
MaxCountExceededException
 at exhaustion of the iteration count, unless a customcallback
has been set at construction of theIterationManager

solve
public RealVector solve(RealLinearOperator a, RealVector b) throws NullArgumentException, NonSquareOperatorException, DimensionMismatchException, NonSelfAdjointOperatorException, IllConditionedOperatorException, MaxCountExceededException
Returns an estimate of the solution to the linear system A · x = b. Overrides:
solve
in classPreconditionedIterativeLinearSolver
 Parameters:
a
 the linear operator A of the systemb
 the righthand side vector Returns:
 a new vector containing the solution
 Throws:
NonSelfAdjointOperatorException
 ifgetCheck()
istrue
, anda
is not selfadjointIllConditionedOperatorException
 ifa
is illconditionedNullArgumentException
 if one of the parameters isnull
NonSquareOperatorException
 ifa
is not squareDimensionMismatchException
 ifb
has dimensions inconsistent witha
MaxCountExceededException
 at exhaustion of the iteration count, unless a customcallback
has been set at construction of theIterationManager

solve
public RealVector solve(RealLinearOperator a, RealVector b, boolean goodb, double shift) throws NullArgumentException, NonSquareOperatorException, DimensionMismatchException, NonSelfAdjointOperatorException, IllConditionedOperatorException, MaxCountExceededException
Returns the solution to the system (A  shift · I) · x = b.If the solution x is expected to contain a large multiple of
b
(as in Rayleighquotient iteration), then better precision may be achieved withgoodb
set totrue
.shift
should be zero if the system A · x = b is to be solved. Otherwise, it could be an approximation to an eigenvalue of A, such as the Rayleigh quotient b^{T} · A · b / (b^{T} · b) corresponding to the vector b. If b is sufficiently like an eigenvector corresponding to an eigenvalue near shift, then the computed x may have very large components. When normalized, x may be closer to an eigenvector than b. Parameters:
a
 the linear operator A of the systemb
 the righthand side vectorgoodb
 usuallyfalse
, except ifx
is expected to contain a large multiple ofb
shift
 the amount to be subtracted to all diagonal elements of A Returns:
 a reference to
x
 Throws:
NullArgumentException
 if one of the parameters isnull
NonSquareOperatorException
 ifa
is not squareDimensionMismatchException
 ifb
has dimensions inconsistent witha
MaxCountExceededException
 at exhaustion of the iteration count, unless a customcallback
has been set at construction of theIterationManager
NonSelfAdjointOperatorException
 ifgetCheck()
istrue
, anda
is not selfadjointIllConditionedOperatorException
 ifa
is illconditioned

solve
public RealVector solve(RealLinearOperator a, RealVector b, RealVector x) throws NullArgumentException, NonSquareOperatorException, DimensionMismatchException, NonSelfAdjointOperatorException, IllConditionedOperatorException, MaxCountExceededException
Returns an estimate of the solution to the linear system A · x = b. Overrides:
solve
in classPreconditionedIterativeLinearSolver
 Parameters:
x
 not meaningful in this implementation; should not be considered as an initial guess (more)a
 the linear operator A of the systemb
 the righthand side vector Returns:
 a new vector containing the solution
 Throws:
NonSelfAdjointOperatorException
 ifgetCheck()
istrue
, anda
is not selfadjointIllConditionedOperatorException
 ifa
is illconditionedNullArgumentException
 if one of the parameters isnull
NonSquareOperatorException
 ifa
is not squareDimensionMismatchException
 ifb
orx0
have dimensions inconsistent witha
MaxCountExceededException
 at exhaustion of the iteration count, unless a customcallback
has been set at construction of theIterationManager

solveInPlace
public RealVector solveInPlace(RealLinearOperator a, RealLinearOperator m, RealVector b, RealVector x) throws NullArgumentException, NonSquareOperatorException, DimensionMismatchException, NonSelfAdjointOperatorException, NonPositiveDefiniteOperatorException, IllConditionedOperatorException, MaxCountExceededException
Returns an estimate of the solution to the linear system A · x = b. The solution is computed inplace (initial guess is modified). Specified by:
solveInPlace
in classPreconditionedIterativeLinearSolver
 Parameters:
x
 the vector to be updated with the solution;x
should not be considered as an initial guess (more)a
 the linear operator A of the systemm
 the preconditioner, M (can benull
)b
 the righthand side vector Returns:
 a reference to
x0
(shallow copy) updated with the solution  Throws:
NonSelfAdjointOperatorException
 ifgetCheck()
istrue
, anda
orm
is not selfadjointNonPositiveDefiniteOperatorException
 ifm
is not positive definiteIllConditionedOperatorException
 ifa
is illconditionedNullArgumentException
 if one of the parameters isnull
NonSquareOperatorException
 ifa
orm
is not squareDimensionMismatchException
 ifm
,b
orx0
have dimensions inconsistent witha
MaxCountExceededException
 at exhaustion of the iteration count, unless a customcallback
has been set at construction of theIterationManager

solveInPlace
public RealVector solveInPlace(RealLinearOperator a, RealLinearOperator m, RealVector b, RealVector x, boolean goodb, double shift) throws NullArgumentException, NonSquareOperatorException, DimensionMismatchException, NonSelfAdjointOperatorException, NonPositiveDefiniteOperatorException, IllConditionedOperatorException, MaxCountExceededException
Returns an estimate of the solution to the linear system (A  shift · I) · x = b. The solution is computed inplace.If the solution x is expected to contain a large multiple of
b
(as in Rayleighquotient iteration), then better precision may be achieved withgoodb
set totrue
; this however requires an extra call to the preconditioner.shift
should be zero if the system A · x = b is to be solved. Otherwise, it could be an approximation to an eigenvalue of A, such as the Rayleigh quotient b^{T} · A · b / (b^{T} · b) corresponding to the vector b. If b is sufficiently like an eigenvector corresponding to an eigenvalue near shift, then the computed x may have very large components. When normalized, x may be closer to an eigenvector than b. Parameters:
a
 the linear operator A of the systemm
 the preconditioner, M (can benull
)b
 the righthand side vectorx
 the vector to be updated with the solution;x
should not be considered as an initial guess (more)goodb
 usuallyfalse
, except ifx
is expected to contain a large multiple ofb
shift
 the amount to be subtracted to all diagonal elements of A Returns:
 a reference to
x
(shallow copy).  Throws:
NullArgumentException
 if one of the parameters isnull
NonSquareOperatorException
 ifa
orm
is not squareDimensionMismatchException
 ifm
,b
orx
have dimensions inconsistent witha
.MaxCountExceededException
 at exhaustion of the iteration count, unless a customcallback
has been set at construction of theIterationManager
NonSelfAdjointOperatorException
 ifgetCheck()
istrue
, anda
orm
is not selfadjointNonPositiveDefiniteOperatorException
 ifm
is not positive definiteIllConditionedOperatorException
 ifa
is illconditioned

solveInPlace
public RealVector solveInPlace(RealLinearOperator a, RealVector b, RealVector x) throws NullArgumentException, NonSquareOperatorException, DimensionMismatchException, NonSelfAdjointOperatorException, IllConditionedOperatorException, MaxCountExceededException
Returns an estimate of the solution to the linear system A · x = b. The solution is computed inplace (initial guess is modified). Overrides:
solveInPlace
in classPreconditionedIterativeLinearSolver
 Parameters:
x
 the vector to be updated with the solution;x
should not be considered as an initial guess (more)a
 the linear operator A of the systemb
 the righthand side vector Returns:
 a reference to
x0
(shallow copy) updated with the solution  Throws:
NonSelfAdjointOperatorException
 ifgetCheck()
istrue
, anda
is not selfadjointIllConditionedOperatorException
 ifa
is illconditionedNullArgumentException
 if one of the parameters isnull
NonSquareOperatorException
 ifa
is not squareDimensionMismatchException
 ifb
orx0
have dimensions inconsistent witha
MaxCountExceededException
 at exhaustion of the iteration count, unless a customcallback
has been set at construction of theIterationManager

