I tested my conjugate gradients implementation with float and double precision and contrary to my guess the double code was twice faster than the single precision code. The reason is that I need many more iterations to achieve the same residual. Mind you, that was already on a pretty small system of about just 2000 unknowns (though the matrix is ill-conditioned). Considering this, are there any remedies to this? For example I looked here and they simply implement the dot product in double precision. Are there any approaches that are accurate and efficient in practice (not just theoretically) that work with single precision (e.g. if I want to use this on the GPU where double comes at a premium)? Or maybe there are some reformulations of the conjugate gradients algorithm that make it more stable (reorthogonalization?).
$\begingroup$
$\endgroup$
3
-
2$\begingroup$ This is an interesting question in my current scope. I will do some more research and try to get back with more information. But current wisdom is to use GMRES with modified (or block modified) Gram-Schmidt instead of CG. It uses more space and computation for sure, but for symmetric matrices GMRES iterations are equivalent to MINRES iterations. In addition, CG and MINRES performance (in terms of iterations) tend to be similar (see Fong, Saunders - CG versus MINRES: an Empirical Comparison). I haven't seen much on stabilizing CG in lower precisions, maybe Erin Carson published something? $\endgroup$– Abdullah Ali SivasNov 17 at 16:20
-
1$\begingroup$ @AbdullahAliSivas I am aware of the paper about CG vs MINRES, but I have been sceptical about it since I ran some experiments several years ago. At the very least in those CG seemed to minimize the error more efficiently than MINRES. Instead MINRES was minimizing the residual better, so I guess it depends on what you want to minimize. I'll look into GMRES - any recommendations for reference discussing practical implementation details? The only thing I have read is the part from Saad's book but from my experience the pseudo-codes from his book needed some modifications in other cases. $\endgroup$– lightxbulbNov 17 at 21:01
-
2$\begingroup$ @lightxbulb I am sure you know about the PETSc implementation of GMRES, which runs on GPUs already (given you use GPU friendly datastructures). I know "Sparse systems solving on GPUs with GMRES" by Couturier and Domas which shows their code. Otherwise, probably there are multiple papers that need to be read at the same time to figure out a practical implementation (at least a PETSc developer told me that their implementation takes pieces from many papers some years ago) $\endgroup$– Abdullah Ali SivasNov 18 at 3:14
Add a comment
|