TL;DR Basically, a preconditioner is meant to speed up the convergence of iterative linear solvers. A preconditioner in Mathematica is the inverse of what is called preconditioner in numerical Mathematics.
The computational backend for LinearSolve[A,b,Method->"Krylov" ]
is SparseArray`KrylovLinearSolve
.
More or less,
SparseArray`KrylovLinearSolve[A,b, "Method" -> "BiCGSTAB", "Preconditioner" -> f]
should be equivalent to
LinearSolve[A,b, Method -> {
"Krylov",
"Method" -> "BiCGSTAB", "Preconditioner" -> "ILU0"
}]
However, calling SparseArray`KrylovLinearSolve
directly is usually a bit faster because LinearSolve
seems to have some overhead.
Other supported Krylov methods are "ConjugateGradient"
(only for symmetric positive-definite matrices) and "GMRES"
. See the documentation of LinearSolve
, section Options, subsection Methods, subsubsection "Krylov".
You can use an arbitrary function f
as preconditioner. Here an example for a built-in preconditioner (incomplete LU-factorization without fill-in):
precdata = SparseArray`SparseMatrixILU[A, "Method" -> "ILU0"]
f = x \[Function] SparseArray`SparseMatrixApplyILU[precdata, x]
Other accepted values for the option "Method"
in SparseArray`SparseMatrixILU
are "ILUT"
and "ILUTP"
.
Another example is this:
f = xWith[{invdiag \[Function]= Evaluate[11/Diagonal[A]Normal[Diagonal[A]]},
x] x |-> invdiag x
];
It is the well-known Jacobi preconditioner -- or rather its inverse.
In constrast to the mathematical nomenclature, a preconditioner itself gets applied, not its inverse, in Mathematica. That's very plausible from an algorithmic point of view.
See this and the standard literature about numerical linear equations for more information about what a preconditioner is. Basically, a preconditioner $f$ is meant to speed up the convergence of iterative linear solvers. Sloppily formulated, $f$ is a good preconditioner for the matrix $A$ if 1. $f(x)$ is easy to compute for each vector $x$ and 2. $A \,f(x)$ is close to $x$. So, a good preconditioner $f$ should best be almost an inverse of $A$ - but orders of magnitude faster than computing $A^{-1} x$ with direct methods.