I'm developing an FEM solver for a coupled system. I have diffusion and potential equations which result in positive definite matrices for each equation, but the coupling makes the overall system matrix non-symmetric, non-positive-definite.

The solver is planned to run entirely on a GPU. Mesh and data representation is highly optimized for that and standard libraries are therefore not very useful, also because one would need to transfer between GPU and CPU, or had to use CPU-optimized data formats that would not perform well on the GPU.

I plan to use multigrid later on, so I currently use jacobi iterations to solve the equation. That works, let's say, semi-well. It converges not very fast and when errors get small, there are some signs of other problems (which I do not fully understand currently).

I'm now thinking about my next step. Should I try some CG-based algorithms? Or what would you try or investigate next?

  • $\begingroup$ You should use PETSc. Good preconditioners typically don't perform well on GPUs (compared to CPUs) but you can test this yourself. $\endgroup$
    – stali
    Oct 20, 2015 at 13:59
  • $\begingroup$ @stali: I clarified on that, see the bold printed section. $\endgroup$
    – Michael
    Oct 20, 2015 at 14:53
  • $\begingroup$ They why not just use NVIDIAs cuSPARSE (CG+ILU0)? $\endgroup$
    – stali
    Oct 20, 2015 at 15:58
  • $\begingroup$ @Michael there are a couple options if multigrid isn't doing great - you can wrap it into a very effective preconditioner for CG, or try a more problem-flexible AMG for the GPU (see for example Luke Olsen or sciencedirect.com/science/article/pii/S0898122114004143), explore GPU-friendly preconditioners (for example cc.gatech.edu/~echow/pubs/parilu-sisc.pdf). Perhaps a little more clarification on what your problems are might help? $\endgroup$
    – Jesse Chan
    Oct 20, 2015 at 17:43

1 Answer 1


The most successful solvers for coupled systems typically apply something like GMRES to the whole system, but then use the block decomposition to derive preconditioners that make use of the properties of individual blocks. If you want to use a multigrid solver, then use it as a solver in the preconditioner.

Since you don't state how the system you consider looks like, it's hard to give concrete advice. That said, the general idea of building preconditioners for block systems is discussed in lecture 38 of my video lectures here: http://www.math.tamu.edu/~bangerth/videos.html


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.