Question:
Is there some pre-conditioner for Conjugate-Gradient (CG) cheap enough, that it is worth using even if my operator is very local (i.e. already has a low number of non-zero elements), as it is when fitting B-spline functions?
or any other way how to improve the speed of CG convergence?
Background:
I fit function described on 2D/3D grid using localized basis function (e.g. B-splines) using Conjugate Gradient (CG) algorithm.
The least-square fitting reduce to soltion of this matrix problem:
$(B^T B) {\vec x}= B {\vec b}$
where $B$ is overlap-matrix between the basis function. Iterative method like CG use matrix-vector multiplications:
$\vec{y} = (B^T B) \vec{x}$
which in my case reduces into 2D/3D convolution with a stencil formed by tensor-product of 1D Cubic Spline (i.e. 3x3 or 3x3x3 tensor in 2D/3D case)
e.g. for 2D case:
$$[ [1,4,1], [4,16,4], [1,4,1], ]$$
You can see that this matrix is very sparse and relatively cheap (although not so cheap if the grid has ~1 million points).
Is there any preconditioner other than Diagonal Jacobi worth using in such a situation?
Currently CG takes 50 iterations ot achieve $10^{-6}$ preccision:
CG[1] err 0.126261
CG[2] err 0.052267
CG[3] err 0.0294687
CG[4] err 0.0193294
...
CG[47] err 4.449e-05
CG[48] err 4.24922e-05
CG[49] err 3.60459e-05
CG[50] err 3.32389e-05
Just in case code can be found here (1,2,3) (but not sure if it is easy to understand)