What's the difference between conjugate gradient method and biconjugate gradient method

Question

What's the difference between these two methods? Can a problem be solved by one method will be able to solved by the other? Can both/or one of them be parallelized with OpenMP and/or MPI?

Answer 1

The conjugate gradient method is the provably fastest iterative solver, but only for symmetric, positive-definite systems. What would be awfully convenient is if there was an iterative method with similar properties for indefinite or non-symmetric matrices.

The CG method seeks approximate solutions at each step $k$ within the Krylov subspace

$K_k(A,b) = \{b, Ab, A^2b,\ldots,A^kb\}$.

The essential idea of the biconjugate gradient method is to maintain a second Krylov subspace

$K_k^*(A,b) = \{b, A^*b, (A^*)^2b,\ldots,(A^*)^kb\}$

and seeking a recurrence with similar orthogonality properties to that of CG, but without the stability issues of solving $A^*Ax = A^*b$.

Unfortunately, that fails if you apply it naively. However, by performing one step of the generalized minimum residual (GMRES) algorithm after each BiCG step, the resulting iteration is stable; this is usually referred to as BiCG-Stab.

So, BiCG-Stab is (in principle) a more general solver than CG but suffers worse efficiency when applied to the problems for which CG was intended. BiCG or BiCG-stab require more matrix-vector multiplications and more dot products, so if you parallelize them via distributed-memory multiprocessing you'll incur more communication overhead, but nonetheless they can be scaled up as much as you like.

There are two things worth noting here which are more important than all that other junk I just said:

1. For every iterative method (BiCG, GMRES, QMR...), there is a matrix that will make it fail to converge in finite-precision arithmetic.

2. Therefore, coming up with a good preconditioner for your specific matrix is probably more important than using the optimal outer-level iterative solver.

EDIT: For open-source libraries, the two most popular are PETSc and Trilinos. I highly recommend you also get the Python bindings, respectively petsc4py and PyTrilinos. You can also try Eigen. On the one hand, it doesn't have many features, but on the other hand, it has just what you need and no more; if you intend to read the code rather than just use it, Eigen might be the easiest.

${}$

See also: Yousef Saad, Iterative Methods for Sparse Linear Systems; Nachtigal et al, How Fast Are Nonsymmetrix Matrix Iterations?

Answer 2

The conjugate gradient method only works to solve the system

$ A x = b $

if $A$ is symmetric and positive-definite (also works for negative definite). The reason it must be symmetric is that conjugate gradient works by minimizing (or maximizing) the function

$ f(x) = \frac{1}{2} x^T A x - b^T x $

Note that the derivative is

$ f'(x) = \frac{1}{2} A^T x + \frac{1}{2} A x - b$

and if $A$ is symmetric $A^T = A$ so the above reduces to

$ f'(x) = A x - b $

At the minimum, $f'(x) = 0$ and $x$ is the solution to your system. This last step should make it obvious why $A$ must be symmetric. The positive/negative definite property is more subtle, but it is required so that the extrema exists.

The biconjugate gradient method will work for any system. It does so by solving both

$ A x = b$

along with

$ A^T x = b$

simultaneously. The details of which I am not familiar with, so I won't pretend to know. It is sufficient to know that biconjugate gradient is the more general of the two. It does have stability issues, and if $A$ is symmetric, then conjugate gradient will perform less work to get to the solution.

As far a parallelism, if you systems are large, you can gain a lot of parallelism in the linear algebra involved in the solver iterations, so there should be no reason a parallel linear algebra library wouldn't lead to performance gains.