With SPD matrix, call it $A$, PCG is a good method.
The kind of preconditioner depend a lot by the matrix structure (read kind of matrix) that you have got and there is not a unique answer.
The basic criteria to choose the preconditioner, $P$, are more than one.
Considering the PCG the three most important are:
- $AP^{-1} \approx I$
- $P$ is cheap to built
- calculated $P^{-1}r_{k+1}$ is cheap
In general for problem like laplacian you can try the two follow.
Diagonal $D^{-1}$ The diagonal preconditioner is very easy to implement and cheap to apply, but is not very efficient. The matrix contains the verse of $A$ diagonal terms. Better if it is used with matrix $A$ diagonally dominant, but when the extra diagonal terms are big you can find problem.
Incomplete Cholesky $\tilde{L}\tilde{L^T}$ Whit $A$ SPD is defined the Cholesky decomposition. The idea of this preconditioner is to use the decomposition but to conserve the sparsity pattern of $A$, i.e you use the term calculated by the decomposition only for the terms of $A$ different from zero. The existence of this incomplete decomposition depend by the structure of $A$, for example if $A$ is a $M$ matrix it is ok. This kind of preconditioner can improve if you reorder the nodes of $A$ so that the mustache of $A$ are shorter (here you can see an example that the structure of $A$ is important).
This preconditioned is apply with the method forward and back (as for the $LU$ decomposition)
When you talk about parallel execution here is important to know the kind of hardware. Anyway you must parallelize the preconditioner application, for example with $D^{-1}$ you parallelize the $\frac{v_i}{d_i}$ operation for the whole vector.
Since you mentioned the Stokes's problem I add some rows about parallelition. If your goal is to use heavy parallelization as gpu is important to use preconditioner with an easy application for gpu, as the product matrix vector sparse (SpMV). In the article of Li and Saad [1] there is a comparison of various preconditioners for gpu.
In the case of gpu polynomial preconditioner is good.
I personally use this kind with good performance respect other traditional gpu preconditioner, see this link to CUSP group where I give some detail, with reference, and also the code with some result test.
This is a least square preconditioner, it uses Chebyshev polynomials so explicit quadrature formula are not involved, detail are in the article by J. Erhel, F. Guyomarc, and Y. Saad [2].
[1] Ruipeng Li and Yousef Saad. Gpu-accelerated preconditioned iterative linear solvers. The Journal of Supercomputing, 63(2):443–466, 2013.
[2] J. Erhel, F. Guyomarc, and Y. Saad, Least-squares polynomial filters for ill-conditioned lin- ear systems, Tech. Report umsi-2001-32, Minnesota Supercomputer Institute, University of Minnesota, Minneapolis, MN, 2001.