• Source: Biconjugate gradient stabilized method

In numerical linear algebra, the biconjugate gradient stabilized method, often abbreviated as BiCGSTAB, is an iterative method developed by H. A. van der Vorst for the numerical solution of nonsymmetric linear systems. It is a variant of the biconjugate gradient method (BiCG) and has faster and smoother convergence than the original BiCG as well as other variants such as the conjugate gradient squared method (CGS). It is a Krylov subspace method. Unlike the original BiCG method, it doesn't require multiplication by the transpose of the system matrix.




Algorithmic steps






= Unpreconditioned BiCGSTAB

=
In the following sections, (x,y) = xT y denotes the dot product of vectors. To solve a linear system Ax = b, BiCGSTAB starts with an initial guess x0 and proceeds as follows:

r0 = b − Ax0
Choose an arbitrary vector r̂0 such that (r̂0, r0) ≠ 0, e.g., r̂0 = r0
ρ0 = (r̂0, r0)
p0 = r0
For i = 1, 2, 3, …
v = Api−1
α = ρi−1/(r̂0, v)
h = xi−1 + αpi−1
s = ri−1 − αv
If h is accurate enough, i.e., if s is small enough, then set xi = h and quit
t = As
ω = (t, s)/(t, t)
xi = h + ωs
ri = s − ωt
If xi is accurate enough, i.e., if ri is small enough, then quit
ρi = (r̂0, ri)
β = (ρi/ρi−1)(α/ω)
pi = ri + β(pi−1 − ωv)
In some cases, choosing the vector r̂0 randomly improves numerical stability.




= Preconditioned BiCGSTAB

=
Preconditioners are usually used to accelerate convergence of iterative methods. To solve a linear system Ax = b with a preconditioner K = K1K2 ≈ A, preconditioned BiCGSTAB starts with an initial guess x0 and proceeds as follows:

r0 = b − Ax0
Choose an arbitrary vector r̂0 such that (r̂0, r0) ≠ 0, e.g., r̂0 = r0
ρ0 = (r̂0, r0)
p0 = r0
For i = 1, 2, 3, …
y = K −12 K −11 pi−1
v = Ay
α = ρi−1/(r̂0, v)
h = xi−1 + αy
s = ri−1 − αv
If h is accurate enough then xi = h and quit
z = K −12 K −11 s
t = Az
ω = (K −11 t, K −11 s)/(K −11 t, K −11 t)
xi = h + ωz
ri = s − ωt
If xi is accurate enough then quit
ρi = (r̂0, ri)
β = (ρi/ρi−1)(α/ω)
pi = ri + β(pi−1 − ωv)
This formulation is equivalent to applying unpreconditioned BiCGSTAB to the explicitly preconditioned system

Ãx̃ = b̃
with à = K −11 AK −12 , x̃ = K2x and b̃ = K −11 b. In other words, both left- and right-preconditioning are possible with this formulation.




Derivation






= BiCG in polynomial form

=
In BiCG, the search directions pi and p̂i and the residuals ri and r̂i are updated using the following recurrence relations:

pi = ri−1 + βipi−1,
p̂i = r̂i−1 + βip̂i−1,
ri = ri−1 − αiApi,
r̂i = r̂i−1 − αiATp̂i.
The constants αi and βi are chosen to be

αi = ρi/(p̂i, Api),
βi = ρi/ρi−1
where ρi = (r̂i−1, ri−1) so that the residuals and the search directions satisfy biorthogonality and biconjugacy, respectively, i.e., for i ≠ j,

(r̂i, rj) = 0,
(p̂i, Apj) = 0.
It is straightforward to show that

ri = Pi(A)r0,
r̂i = Pi(AT)r̂0,
pi+1 = Ti(A)r0,
p̂i+1 = Ti(AT)r̂0
where Pi(A) and Ti(A) are ith-degree polynomials in A. These polynomials satisfy the following recurrence relations:

Pi(A) = Pi−1(A) − αiATi−1(A),
Ti(A) = Pi(A) + βi+1Ti−1(A).




= Derivation of BiCGSTAB from BiCG

=
It is unnecessary to explicitly keep track of the residuals and search directions of BiCG. In other words, the BiCG iterations can be performed implicitly. In BiCGSTAB, one wishes to have recurrence relations for

r̃i = Qi(A)Pi(A)r0
where Qi(A) = (I − ω1A)(I − ω2A)⋯(I − ωiA) with suitable constants ωj instead of ri = Pi(A)r0 in the hope that Qi(A) will enable faster and smoother convergence in r̃i than ri.
It follows from the recurrence relations for Pi(A) and Ti(A) and the definition of Qi(A) that

Qi(A)Pi(A)r0 = (I − ωiA)(Qi−1(A)Pi−1(A)r0 − αiAQi−1(A)Ti−1(A)r0),
which entails the necessity of a recurrence relation for Qi(A)Ti(A)r0. This can also be derived from the BiCG relations:

Qi(A)Ti(A)r0 = Qi(A)Pi(A)r0 + βi+1(I − ωiA)Qi−1(A)Pi−1(A)r0.
Similarly to defining r̃i, BiCGSTAB defines

p̃i+1 = Qi(A)Ti(A)r0.
Written in vector form, the recurrence relations for p̃i and r̃i are

p̃i = r̃i−1 + βi(I − ωi−1A)p̃i−1,
r̃i = (I − ωiA)(r̃i−1 − αiAp̃i).
To derive a recurrence relation for xi, define

si = r̃i−1 − αiAp̃i.
The recurrence relation for r̃i can then be written as

r̃i = r̃i−1 − αiAp̃i − ωiAsi,
which corresponds to

xi = xi−1 + αip̃i + ωisi.




= Determination of BiCGSTAB constants

=
Now it remains to determine the BiCG constants αi and βi and choose a suitable ωi.
In BiCG, βi = ρi/ρi−1 with

ρi = (r̂i−1, ri−1) = (Pi−1(AT)r̂0, Pi−1(A)r0).
Since BiCGSTAB does not explicitly keep track of r̂i or ri, ρi is not immediately computable from this formula. However, it can be related to the scalar

ρ̃i = (Qi−1(AT)r̂0, Pi−1(A)r0) = (r̂0, Qi−1(A)Pi−1(A)r0) = (r̂0, ri−1).
Due to biorthogonality, ri−1 = Pi−1(A)r0 is orthogonal to Ui−2(AT)r̂0 where Ui−2(AT) is any polynomial of degree i − 2 in AT. Hence, only the highest-order terms of Pi−1(AT) and Qi−1(AT) matter in the dot products (Pi−1(AT)r̂0, Pi−1(A)r0) and (Qi−1(AT)r̂0, Pi−1(A)r0). The leading coefficients of Pi−1(AT) and Qi−1(AT) are (−1)i−1α1α2⋯αi−1 and (−1)i−1ω1ω2⋯ωi−1, respectively. It follows that

ρi = (α1/ω1)(α2/ω2)⋯(αi−1/ωi−1)ρ̃i,
and thus

βi = ρi/ρi−1 = (ρ̃i/ρ̃i−1)(αi−1/ωi−1).
A simple formula for αi can be similarly derived. In BiCG,

αi = ρi/(p̂i, Api) = (Pi−1(AT)r̂0, Pi−1(A)r0)/(Ti−1(AT)r̂0, ATi−1(A)r0).
Similarly to the case above, only the highest-order terms of Pi−1(AT) and Ti−1(AT) matter in the dot products thanks to biorthogonality and biconjugacy. It happens that Pi−1(AT) and Ti−1(AT) have the same leading coefficient. Thus, they can be replaced simultaneously with Qi−1(AT) in the formula, which leads to

αi = (Qi−1(AT)r̂0, Pi−1(A)r0)/(Qi−1(AT)r̂0, ATi−1(A)r0) = ρ̃i/(r̂0, AQi−1(A)Ti−1(A)r0) = ρ̃i/(r̂0, Ap̃i).
Finally, BiCGSTAB selects ωi to minimize r̃i = (I − ωiA)si in 2-norm as a function of ωi. This is achieved when

((I − ωiA)si, Asi) = 0,
giving the optimal value

ωi = (Asi, si)/(Asi, Asi).




Generalization


BiCGSTAB can be viewed as a combination of BiCG and GMRES where each BiCG step is followed by a GMRES(1) (i.e., GMRES restarted at each step) step to repair the irregular convergence behavior of CGS, as an improvement of which BiCGSTAB was developed. However, due to the use of degree-one minimum residual polynomials, such repair may not be effective if the matrix A has large complex eigenpairs. In such cases, BiCGSTAB is likely to stagnate, as confirmed by numerical experiments.
One may expect that higher-degree minimum residual polynomials may better handle this situation. This gives rise to algorithms including BiCGSTAB2[1] and the more general BiCGSTAB(l)[2]. In BiCGSTAB(l), a GMRES(l) step follows every l BiCG steps. BiCGSTAB2 is equivalent to BiCGSTAB(l) with l = 2.




See also


Biconjugate gradient method
Conjugate gradient squared method
Conjugate gradient method




References



Van der Vorst, H. A. (1992). "Bi-CGSTAB: A Fast and Smoothly Converging Variant of Bi-CG for the Solution of Nonsymmetric Linear Systems". SIAM J. Sci. Stat. Comput. 13 (2): 631–644. doi:10.1137/0913035. hdl:10338.dmlcz/104566.
Saad, Y. (2003). "§7.4.2 BICGSTAB". Iterative Methods for Sparse Linear Systems (2nd ed.). SIAM. pp. 231–234. ISBN 978-0-89871-534-7.
^ Gutknecht, M. H. (1993). "Variants of BICGSTAB for Matrices with Complex Spectrum". SIAM J. Sci. Comput. 14 (5): 1020–1033. doi:10.1137/0914062.
^ Sleijpen, G. L. G.; Fokkema, D. R. (November 1993). "BiCGstab(l) for linear equations involving unsymmetric matrices with complex spectrum" (PDF). Electronic Transactions on Numerical Analysis. 1. Kent, OH: Kent State University: 11–32. ISSN 1068-9613.

Kata Kunci Pencarian:

• Biconjugate gradient stabilized method
• Biconjugate gradient method
• Gradient method
• Conjugate gradient squared method
• List of numerical analysis topics
• Lis (linear algebra library)
• Krylov subspace