# What does LinearSolve Preconditioner do?

I have looked at a few places in order to find an answer to my question - specifically, I have looked at this old Mathematica documentation page, at the Preconditioner option for NDSolve, as explained here, also at this older question and at this unanswered post here, but I still don't understand what the Preconditioner option does. Does it multiply the matrix and the vector in the linear solve by something (or by its inverse)? How and from where does it take arguments? How can I specify a user-defined function and what should be the form of that?

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 SparseArrayKrylovLinearSolve.

More or less,

SparseArrayKrylovLinearSolve[A,b, "Method" -> "BiCGSTAB", "Preconditioner" -> f]


should be equivalent to

LinearSolve[A,b, Method -> {
"Krylov",
"Method" -> "BiCGSTAB", "Preconditioner" -> "ILU0"
}]


However, calling SparseArrayKrylovLinearSolve 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 = SparseArraySparseMatrixILU[A, "Method" -> "ILU0"]
f = x \[Function] SparseArraySparseMatrixApplyILU[precdata, x]


Other accepted values for the option "Method" in SparseArraySparseMatrixILU are "ILUT" and "ILUTP".

Another example is this:

f = x \[Function] Evaluate[1/Diagonal[A] 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.

• I had a look at the page you have linked adn the LinearSolve documentation, before posting my question. What is x in your example? f = x [Function] Evaluate[1/Diagonal[A] x] Also you said "Using SparseArrayKrylovLinearSolve", but this is just a different linear solver, isn't it? How can I build preconditioners with it? I can't find any documentation on it. And finally, if I understand correctly "In constrast to the mathematical nomenclature, a preconditioner itself gets applied, not its inverse, in Mathematica", means I have to pass Mathematica the inverse of P (Wiki notation) – ThunderBiggi May 1 '18 at 13:44
• x is a certain vector that occurs in the iterative algorithms (some sort of residual). SparseArrayKrylovLinearSolve is the backend of LinearSolve[A,b,Method -> "Krylov"]. The only built-in precondtioners I know of are those generated by SparseArraySparseMatrixILU. And yes, there is no documentation about it. "And finally...": You are correct. – Henrik Schumacher May 1 '18 at 13:59
• Thank you! Can you give me an example how to use SparseArray""KrylovLinearSolve in order to build a preconditioner? – ThunderBiggi May 1 '18 at 14:31
• I've already given an example (it is SparseArraySparseMatrixApplyILU that creates data for a preconditioner). – Henrik Schumacher May 1 '18 at 14:32
• The preconditioners in Mathematica are meant to be actions on a vector. This particular function is a fast version of f = x \[Function] Evaluate[DiagonalMatrix[1/Diagonal[A]].x]`... ;) – Henrik Schumacher May 1 '18 at 14:57