Fork me on GitHub

6. The Bi-Conjugate Gradient Stabilized Method

6.1. The bicgstab Module

class bicgstab.BiCGSTAB(op, **kwargs)

Bases: pykrylov.generic.generic.KrylovMethod

A pure Python implementation of the bi-conjugate gradient stabilized (Bi-CGSTAB) algorithm. Bi-CGSTAB may be used to solve unsymmetric systems of linear equations, i.e., systems of the form

A x = b

where the operator A is unsymmetric and nonsingular.

Bi-CGSTAB requires 2 operator-vector products, 6 dot products and 6 daxpys per iteration.

In addition, if a preconditioner is supplied, it needs to solve 2 preconditioning systems per iteration.

The original description appears in [VdVorst92]. Our implementation is a preconditioned version of that given in [Kelley].

Reference:

[VdVorst92]H. Van der Vorst, Bi-CGSTAB: A Fast and Smoothly Convergent Variant of Bi-CG for the Solution of Nonsymmetric Linear Systems, SIAM Journal on Scientific and Statistical Computing 13 (2), pp. 631–644, 1992.
solve(rhs, **kwargs)

Solve a linear system with rhs as right-hand side by the Bi-CGSTAB method. The vector rhs should be a Numpy array.

Keywords:
guess:Initial guess (Numpy array, default: 0)
matvec_max:Max. number of matrix-vector produts (2n)

6.2. Example

Here is an example using Bi-CGSTAB on a linear system. The coefficient matrix is read from file in Matrix Market format:

import numpy as np
from pykrylov.cgs import BiCGSTAB as KSolver
from pykrylov.linop import PysparseLinearOperator
from pysparse import spmatrix
from pysparse.pysparseMatrix import PysparseMatrix as sp

A = sp(matrix=spmatrix.ll_mat_from_mtx('jpwh_991.mtx'))
n = A.shape[0]
e = np.ones(n)
rhs = A*e

ks = KSolver( PysparseLinearOperator(A),
              matvec_max=2*n,
              verbose=False,
              reltol = 1.0e-5 )
ks.solve(rhs, guess = 1+np.arange(n, dtype=np.float))

print 'Number of matvecs: ', ks.nMatvec
print 'Initial/final res: %8.2e/%8.2e' % (ks.residNorm0, ks.residNorm)
print 'Direct error: %8.2e' % (np.linalg.norm(ks.bestSolution-e)/sqrt(n))

Running this script produces the following output:

Number of matvecs:  57
Initial/final res: 8.64e+03/5.18e-02
Direct error: 3.35e-03