mpkrir: Krylov-IR solver
MultiPrecisionArrays.mpkrir
— Methodmpkrir(AF::MPKFact, b; reporting=false, verbose=false, mpdebug=false, termparms=termparms_default)
Krylov-IR solver
I do not export this function. The idea is that you use mpglu
and do not touch either the constructor or the solver directly.
Use a multi-precision factorization to solve a linear system with IR-Krylov metods
This is the generic solver used by GMRES-IR and BiCGSTAB-IR. You use the correct MPKFact = Union{MPGFact,MPBFact} structure and mpkrir will do the right thing.
You should not be calling this directly. Use \
to solve linear systems with the multiprecision factorization and to use the optional kwargs.
We overload the backslash operator to call mpkrir for a multiprecision MPGFact factorization. So if MPA
is an MPGArray and
AF = mpglu!(MPA)
Then AF\b
maps to
mpkrir(AF, b)
which does the GMRES-IR solve. You can also get the multiprecision factoriztion directly with
AF = mpglu(A)
which builds the mutliprcision MPGArray and then factors the low preicsion copy.
Similarly if MPA
is an MPBArray. Then
AF = mpblu!(MPA)
Then AF\b
maps to
mpkrir(AF, b)
which does the BiCGSTAB-IR solve.
You can also use the \
operator to harvest iteration statistics.
You may not get exactly the same results for this example on different hardware, BLAS, versions of Julia or this package. I am still playing with the termination criteria and the iteration count could grow or shrink as I do that.
Example
```jldoctest julia> using MultiPrecisionArrays.Examples
julia> N=4096; A = I - 800.0 * Gmat(N); b=ones(N);
julia> AF=mpglu(A);
julia> solout=(AF, b; reporting=true);
Correct result?
julia> x=solout.sol; norm(b-A*x,Inf) 8.20877e-12
Look at the residual history
julia> solout.rhist
4-element Vector{Float64}: 1.00000e+00 1.25881e-10 8.58902e-12 8.20877e-12
and the correction norm history. No correction for the final iteration.
julia> solout.dhist 3-element Vector{Float64}: 2.00129e+02 1.05241e-10 2.28629e-11
Stagnation after the 2nd iteration. Now the Krylovs/iteration
julia> solout.khist 3-element Vector{Int64}: 4 5 5
4-5 Krylovs per iteration.
BiCGSTAB works the same way.