krypy.recycling
- Recycling Linear Systems Solvers¶
The recycling
module provides functions for the solution of sequences of
linear algebraic systems. Once a linear system has been solved, the
generated data is examined and a deflation space is determined
automatically for the solution of the next linear system. Several selection
strategys are available.
-
class
krypy.recycling.
RecyclingCg
(*args, **kwargs)¶ Bases:
krypy.recycling.linsys._RecyclingSolver
Recycling preconditioned CG method.
See
_RecyclingSolver
for the documentation of the available parameters.
-
class
krypy.recycling.
RecyclingMinres
(*args, **kwargs)¶ Bases:
krypy.recycling.linsys._RecyclingSolver
Recycling preconditioned MINRES method.
See
_RecyclingSolver
for the documentation of the available parameters.
-
class
krypy.recycling.
RecyclingGmres
(*args, **kwargs)¶ Bases:
krypy.recycling.linsys._RecyclingSolver
Recycling preconditioned GMRES method.
See
_RecyclingSolver
for the documentation of the available parameters.
-
class
krypy.recycling.linsys.
_RecyclingSolver
(DeflatedSolver, vector_factory=None)¶ Base class for recycling solvers.
Initialize recycling solver base.
Parameters: - DeflatedSolver – a deflated solver from
deflation
. - vector_factory –
(optional) An instance of a subclass of
krypy.recycling.factories._DeflationVectorFactory
that constructs deflation vectors for recycling. Defaults to None which means that no recycling is used.Also the following strings are allowed as shortcuts:
'RitzApproxKrylov'
: uses the approximate Krylov subspace bound evaluatorkrypy.recycling.evaluators.RitzApproxKrylov
.'RitzAprioriCg'
: uses the CG \(\kappa\)-bound (krypy.utils.BoundCG
) as an a priori bound withkrypy.recycling.evaluators.RitzApriori
.'RitzAprioriMinres'
: uses the MINRES bound (krypy.utils.BoundMinres
) as an a priori bound withkrypy.recycling.evaluators.RitzApriori
.
After a run of the provided
DeflatedSolver
viasolve()
, the resulting instance of theDeflatedSolver
is available in the attributelast_solver
.-
last_solver
= None¶ DeflatedSolver
instance from last run ofsolve()
.Instance of
DeflatedSolver
that resulted from the last call tosolve()
. Initialized withNone
before the first run.
-
solve
(linear_system, vector_factory=None, *args, **kwargs)¶ Solve the given linear system with recycling.
The provided vector_factory determines which vectors are used for deflation.
Parameters: - linear_system – the
LinearSystem
that is about to be solved. - vector_factory – (optional) see description in constructor.
All remaining arguments are passed to the
DeflatedSolver
.Returns: instance of DeflatedSolver
which was used to obtain the approximate solution. The approximate solution is available under the attributexk
.- linear_system – the
- DeflatedSolver – a deflated solver from