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?


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)


1 Answer 1


You certainly can try using Gauss-Seidel based preconditioning that is relatively easy to construct and is cheap enough (by my assessment) to give it a try.

Other choices of algebraic preconditioners (say, incomplete LU) might be too heavy relative to their potential impact.

I would also suggest looking into a full sparse direct solution of your problem. There are quite efficient implementations of sparse direct solvers nowadays, so it might be a reasonable alternative. By quick googling I was able to find at least this publicstion; however, I bet there are more canonical papers of this flavor.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.