# A Polynomial Preconditioner for the CMRH Algorithm

Authors: Jiangzhou Lai; Linzhang Lu; Shiji Xu
Journal: Mathematical Problems in Engineering
 (2011)

---

## Abstract

Many large and sparse linear systems can be solved efficiently by restarted GMRES and
CMRH methods Sadok 1999. The CMRH(m) method is less expensive and requires slightly less storage
than GMRES(m). But like GMRES, the restarted CMRH method may not converge. In
order to remedy this defect, this paper presents a polynomial preconditioner for CMRH-based
algorithm. Numerical experiments are given to show that the polynomial preconditioner is
quite simple and easily constructed and the preconditioned CMRH(m) with the polynomial
preconditioner has better performance than CMRH(m).


---

## Body

## 1. Introduction

In many practical applications, we have to solve a large and sparse, nonsymmetric linear system of equations(1.1)Ax=b,
where A∈ℛn×n is nonsingular, and x,b∈ℛn. The linear system (1.1) can be solved efficiently by iterative methods, especially those based on the Krylov subspace(1.2)Kk(r0,A)=span⁡{r0,Ar0,…,Ak-1r0}.
Although the GMRES method [1] and CMRH method [2] are Krylov-type algorithms, they construct different Krylov subspace bases in different ways; the GMRES method uses the Arnoldi process, while CMRH uses the Hessenberg process. As the number of iterations increases, the methods become impractical because of growth of memory and computational requirement as k increases, so they must be remedied by restarting; thus, the GMRES(m) and CMRH(m) algorithms were developed. It is concluded in [2] that the CMRH(m) algorithm only requires half the arithmetical operations per iteration and slightly less storage than GMRES(m). However, like GMRES(m), the CMRH(m) method may not converge in some cases. In order to remedy this defect, van Gijzen [3] presented a polynomial preconditioner to improve the GMRES(m) method. The main idea is finding the low-degree polynomials pm(x) satisfying pm(A)≈A-1 and transferring the original linear system of (1.1) into(1.3)pm(A)Ax=pm(A)b,
then applying the iterative method based on the following Krylov subspace:(1.4)Kk(pm(A)A;r0)=span⁡{r0,pm(A)Ar0,…,(pm(A)A)k-1r0}
to the new linear system (1.3), obviously, Kk(pm(A)A;r0)∈K(k+1)(m-1)+1(A;r0). This method becomes very effective when large computation of inner product is needed. In this paper, we apply the same technique to the CMRH(m) algorithm. Our method is efficient since the CMRH process itself is utilized to construct the polynomial preconditioner. Numerical experiments are given to show that the CMRH(m) method with a polynomial preconditioner has better performance, in the sense of fewer iterations, than the original CMRH(m) method.The rest of the paper is arranged as follows. The first part of Section2 reviews the CMRH(m) method. The construction of the polynomial preconditioner is described in the second part of Section 2, and in the third part, the PCMRH(m) method is developed. We present some numerical experiments to show that the polynomial preconditioner is effective for the CMRH(m) method in Section 3. Section 4 is a simple conclusion.Throughout the paper, all vectors and matrices are assumed to be real. For a vectorv, ∥v∥ denotes the Euclidean norm ∥v∥=vTv, and ∥v∥∞ denotes the maximum norm ∥v∥∞=maxi=1,…n|vi|, where vi is the ith component of a vector v, and |vi| denotes the absolute value(modulus) of a real(complex) number vi. For a matrix A, ∥A∥ denotes the 2-norm. Furthermore, we denote by ek(n) the kth canonical vector of Rn,(1.5)ek(n)=(0,…,1,0,…,0)T,
and Ik(n) the first kth columns of an n-dimensional identity matrix, and we denote by H† the pseudoinverse of a matrix H.


## 2. PCMRH(m) Algorithm

### 2.1. CMRH(m) Algorithm

The principle of CMRH [2] method is to use the Hessenberg process with pivoting to construct a basis of the Krylov subspace Kk(A,r0): l1,l2,…,lk. Precisely, the Hessenberg process produces an n×k unit lower trapezoidal matrix Lk=(l1,…,lk) and a (k+1)×k upper Hessenberg matrix H¯k∈ℛ(k+1)×k, such that(2.1)ALk=Lk+1H¯k.
Then, we find a correction vector z∈Kk(A,r0), such that(2.2)min⁡z∈Kk(A,r0)‖b-A(x0+z)‖=min⁡z∈Kk(A,r0)‖r0-Az‖.
Since Kk(A,r0)=span{l1,l2,…,lk}, we let z=Lkyk and have(2.3)min⁡z∈Kk(A,r0)‖b-A(x0+zk)‖=min⁡yk∈Rk‖r0-ALkyk‖=min⁡yk∈Rk‖r0-Lk+1H¯kyk‖=min⁡yk∈Rk‖e1(k+1)‖r0‖e1(k+1)-H¯kyk‖.
Therefore, a least-squares solution could be obtained from the Krylov subspace by minimizing ∥∥r0∥e1(k+1)-H¯kyk∥ over y, just as in the GMRES method, except that the Hessenberg process is replaced by the Arnoldi process; this is the basic idea of the CMRH method. Since the Hessenberg process will be breakdown if (lj)j is zero, Sadok [2] uses a pivoting strategy such as in the Gaussian elimination method to avoid such a breakdown, and then the approximating solution can be computing by (2.4)xk=x0+LkH†(r0)i0e1k+1,
where i0 satisfies |(r0)i0|=∥r0∥∞, Lk=(l1,…,lk) is an n×k unit lower trapezoidal matrix, and H† denotes the pseudoinverse of a matrix H for more details, see (Algorithm 1) [2].Algorithm 1:Hessenberg process with pivoting [2].
setp=(1,…,n)T.Determinei0 so |(v)i0|=∥v∥∞; l1=v0(v)i0; p(1)↔p(i0) (interchange p(1)  and p(i0));iterate:fork=1,2,…,nu=Alkforj=1,2,…,kc=up(j); hj,k=c; up(j)=0forl=j+1,…,nup(l)=up(l)-c*(lj)p(l)endendifk<ndeterminei0  so  |(u)i0|=∥u∥∞; p(k+1)↔p(i0)hk+1,k=ui0; lk+1=u/(u)i0endendThe notation “↔” means as “swap contents”: α↔β⇔γ=α;α=β;β=γ. As k increases, the CMRH method(Hessenberg process with pivot) becomes impractical because of the growth of memory and computational requirement, so it must be remedied by restarting, thus the CMRH(m) algorithm is developed, which is described in (Algorithm 2) [2].Algorithm 2:CMRH(m) [2].
Input: m: the dimension of a Krylov subspace and the approximation precisionεstart:x0=0,r0=b-Ax0, and set p(i)=i, i=1,…,nDeterminei0  so|(r0)i0|=∥r0∥∞; l1=r0(r0)i0; p(1)↔p(i0) (interchange p(1) and p(i0));iterate:fork=1,2,…,mu=Alkforj=1,2,…,kc=up(j); hj,k=c; up(j)=0forl=j+1,…,nup(l)=up(l)-c*(lj)p(l)endendifk<ndeterminei0so|(u0)i0|=∥u0∥∞; p(k+1)↔p(i0)hk+1,k=ui0; lk+1=u/(u)i0endif (an estimate of)∥b-Axk∥ is small enough ork=n thenxk=x0+Lkyk, where yk minimize ∥H¯ky-(r0)i0e1(k+1)∥stop iterationendendxm=x0+Lmym, where ym∈ℛm minimizes ∥H¯my-(r0)i0e1(m+1)∥.x0:=xm, go to startThe main body of Algorithm2 is the Hessenberg process with pivoting [2], in the last, we get xm=x0+Lmym, where ym∈ℛm minimizes ∥H¯my-(r0)i0e1(m+1)∥; this is a little different from the GMRES method because of the pivot strategy used in the Hessenberg process.


### 2.2. Construction of the Polynomial Preconditioner

The main idea of the polynomial preconditioned method is to construct a polynomials(t) satisfying s(A)≈A-1 and then solve the linear system of equations s(A)Ax=s(A)b instead of solving the linear system (1.1). Generally, we do not need to compute s(A)A really, since what we need is only that the matrix s(A)A approaches the unit matrix. In practice, if a polynomial qn is constructed, such that xk+1=x0+qk(A)r0 is an approximation of the solution of the linear system (1.1), then when x0=0, the following relationship holds: (2.5)A-1b≈xk+1=qk(A)b.
That is, qk(A) is a polynomial preconditioner.In the following, we utilize the residual polynomial produced from the CMRH process to construct a polynomial preconditioner.LetKk be an n×k matrix formed by Krylov vectors, that is, (2.6)Kk=(r0,Ar0,…,Ak-1r0).
From (2.1), we get(2.7)lk+1=hk+1,k-1(Alk-Lkhk),
where hk=(h1k,…,hkk)T is the kth column of H¯k. Since the columns of Lk span the same space as the columns of Kk, the following relationship holds:(2.8)Lk=KkCk,
where Ck is an upper triangular matrix. Thus, (2.9)Lkhk=KkCkhk=(KkAkr0)(Ckhk0),Alk=AKk(c1k⋮ckk)=(Ar0A2r0⋯Akr0)(c1k⋮ckk)=Kk+1(0c1k⋮ckk).
Hence, from (2.7), we have(2.10)lk+1=hk+1,k-1Kk+1((0c1k⋮ckk)-(Ckhk0)).
On the other hand, from Lk+1=Kk+1Ck+1, we have lk+1=Kk+1(c1,k+1⋮ck+1,k+1). Since the columns of Kk+1 are linear independent, we obtain(2.11)(c1,k+1c2,k+1⋮ck+1,k+1)=hk+1,k-1(0c1k⋮ckk)-hk+1,k-1(Ckhk0).
Thus, from the CMRH algorithm, we can obtain the approximate solution (2.12)xk+1=x0+Lk+1y=x0+Kk+1Ck+1y.
Since Lk+1y=x0+Kk+1Ck+1y, we know that the entries of Ck+1y are the coefficients of qk(z), that is, if Ck+1y=(α0,α1,…,αk)T, then qk(z)=α0+α1z+⋯+αkzk. Note that pk+1(z)=1-zqk(z) is the residual polynomial, so we have qk(z)=(1-pk+1(z))/z and qk(A)≈A-1 [4]; therefore, qk(A) becomes a preconditioner.In the next subsection, we present the preconditioned CMRH(m) with the polynomial preconditioner described above.


### 2.3. PCMRH(m) Algorithm

The polynomial preconditioned CMRH(m) (PCMRH(m)) method is given in the following, which is composed of two parts, in the first part, the polynomial preconditioner is constructed, and in the second part, the CMRH(m) method is applied to the new linear systems (Algorithm3).Algorithm 3:PCMRH(m).
Input:  the dimension ofKrylov subspace m; the prescribed stopping criterion ɛ;the degree of the polynomial preconditionerkk; the initial approximate solution x0;the residual vectorr0=b-Ax0.1. The procedure of constructing polynomial preconditionerqkk(1)  Let p(i)=i,i=1,2,…,n;determinei0so|(r0)|i0=∥r0∥∞;l1=r0(r0)i0;p(1)↔p(i0),C(1,1)=1(r0)i0(2)  for k=1,2,…,kkdou=Alkforj=1,2,…,kdoc=up(j),hj,k=c,up(j)=0forl=j+1,…,ndoup(l)=up(l)-c×(lj)p(l)endifk<ncomputei0so|(u0)|i0=∥u0∥∞;p(k+1)↔p(i0)hk+1,k=ui0;lk+1=u(u)i0lethj=(h1,j⋮hj,j),Cj=(c11⋯c1j0⋱⋮00cjj)compute(c1,j+1c2,j+1⋮cj+1,j+1)=hj+1,j-1(0c1j⋮cjj)-hj+1,j-1(Cjhj0)endend(3) Form the approximate solution xk=x0+Lkyk, where yk∈CkminimizesF(y)=∥H¯ky-(r0)i0e1∥.(4) Form the polynomial preconditioner qkk(z)Ckyk=(α0⋯αk-1)T;qkk(x)=α0+α1x+⋯+αk-1xk-1.2. Using the CMRH(m) algorithm for qkk(A)Ax=qkk(A)b.


## 2.1. CMRH(m) Algorithm

The principle of CMRH [2] method is to use the Hessenberg process with pivoting to construct a basis of the Krylov subspace Kk(A,r0): l1,l2,…,lk. Precisely, the Hessenberg process produces an n×k unit lower trapezoidal matrix Lk=(l1,…,lk) and a (k+1)×k upper Hessenberg matrix H¯k∈ℛ(k+1)×k, such that(2.1)ALk=Lk+1H¯k.
Then, we find a correction vector z∈Kk(A,r0), such that(2.2)min⁡z∈Kk(A,r0)‖b-A(x0+z)‖=min⁡z∈Kk(A,r0)‖r0-Az‖.
Since Kk(A,r0)=span{l1,l2,…,lk}, we let z=Lkyk and have(2.3)min⁡z∈Kk(A,r0)‖b-A(x0+zk)‖=min⁡yk∈Rk‖r0-ALkyk‖=min⁡yk∈Rk‖r0-Lk+1H¯kyk‖=min⁡yk∈Rk‖e1(k+1)‖r0‖e1(k+1)-H¯kyk‖.
Therefore, a least-squares solution could be obtained from the Krylov subspace by minimizing ∥∥r0∥e1(k+1)-H¯kyk∥ over y, just as in the GMRES method, except that the Hessenberg process is replaced by the Arnoldi process; this is the basic idea of the CMRH method. Since the Hessenberg process will be breakdown if (lj)j is zero, Sadok [2] uses a pivoting strategy such as in the Gaussian elimination method to avoid such a breakdown, and then the approximating solution can be computing by (2.4)xk=x0+LkH†(r0)i0e1k+1,
where i0 satisfies |(r0)i0|=∥r0∥∞, Lk=(l1,…,lk) is an n×k unit lower trapezoidal matrix, and H† denotes the pseudoinverse of a matrix H for more details, see (Algorithm 1) [2].Algorithm 1:Hessenberg process with pivoting [2].
setp=(1,…,n)T.Determinei0 so |(v)i0|=∥v∥∞; l1=v0(v)i0; p(1)↔p(i0) (interchange p(1)  and p(i0));iterate:fork=1,2,…,nu=Alkforj=1,2,…,kc=up(j); hj,k=c; up(j)=0forl=j+1,…,nup(l)=up(l)-c*(lj)p(l)endendifk<ndeterminei0  so  |(u)i0|=∥u∥∞; p(k+1)↔p(i0)hk+1,k=ui0; lk+1=u/(u)i0endendThe notation “↔” means as “swap contents”: α↔β⇔γ=α;α=β;β=γ. As k increases, the CMRH method(Hessenberg process with pivot) becomes impractical because of the growth of memory and computational requirement, so it must be remedied by restarting, thus the CMRH(m) algorithm is developed, which is described in (Algorithm 2) [2].Algorithm 2:CMRH(m) [2].
Input: m: the dimension of a Krylov subspace and the approximation precisionεstart:x0=0,r0=b-Ax0, and set p(i)=i, i=1,…,nDeterminei0  so|(r0)i0|=∥r0∥∞; l1=r0(r0)i0; p(1)↔p(i0) (interchange p(1) and p(i0));iterate:fork=1,2,…,mu=Alkforj=1,2,…,kc=up(j); hj,k=c; up(j)=0forl=j+1,…,nup(l)=up(l)-c*(lj)p(l)endendifk<ndeterminei0so|(u0)i0|=∥u0∥∞; p(k+1)↔p(i0)hk+1,k=ui0; lk+1=u/(u)i0endif (an estimate of)∥b-Axk∥ is small enough ork=n thenxk=x0+Lkyk, where yk minimize ∥H¯ky-(r0)i0e1(k+1)∥stop iterationendendxm=x0+Lmym, where ym∈ℛm minimizes ∥H¯my-(r0)i0e1(m+1)∥.x0:=xm, go to startThe main body of Algorithm2 is the Hessenberg process with pivoting [2], in the last, we get xm=x0+Lmym, where ym∈ℛm minimizes ∥H¯my-(r0)i0e1(m+1)∥; this is a little different from the GMRES method because of the pivot strategy used in the Hessenberg process.


## 2.2. Construction of the Polynomial Preconditioner

The main idea of the polynomial preconditioned method is to construct a polynomials(t) satisfying s(A)≈A-1 and then solve the linear system of equations s(A)Ax=s(A)b instead of solving the linear system (1.1). Generally, we do not need to compute s(A)A really, since what we need is only that the matrix s(A)A approaches the unit matrix. In practice, if a polynomial qn is constructed, such that xk+1=x0+qk(A)r0 is an approximation of the solution of the linear system (1.1), then when x0=0, the following relationship holds: (2.5)A-1b≈xk+1=qk(A)b.
That is, qk(A) is a polynomial preconditioner.In the following, we utilize the residual polynomial produced from the CMRH process to construct a polynomial preconditioner.LetKk be an n×k matrix formed by Krylov vectors, that is, (2.6)Kk=(r0,Ar0,…,Ak-1r0).
From (2.1), we get(2.7)lk+1=hk+1,k-1(Alk-Lkhk),
where hk=(h1k,…,hkk)T is the kth column of H¯k. Since the columns of Lk span the same space as the columns of Kk, the following relationship holds:(2.8)Lk=KkCk,
where Ck is an upper triangular matrix. Thus, (2.9)Lkhk=KkCkhk=(KkAkr0)(Ckhk0),Alk=AKk(c1k⋮ckk)=(Ar0A2r0⋯Akr0)(c1k⋮ckk)=Kk+1(0c1k⋮ckk).
Hence, from (2.7), we have(2.10)lk+1=hk+1,k-1Kk+1((0c1k⋮ckk)-(Ckhk0)).
On the other hand, from Lk+1=Kk+1Ck+1, we have lk+1=Kk+1(c1,k+1⋮ck+1,k+1). Since the columns of Kk+1 are linear independent, we obtain(2.11)(c1,k+1c2,k+1⋮ck+1,k+1)=hk+1,k-1(0c1k⋮ckk)-hk+1,k-1(Ckhk0).
Thus, from the CMRH algorithm, we can obtain the approximate solution (2.12)xk+1=x0+Lk+1y=x0+Kk+1Ck+1y.
Since Lk+1y=x0+Kk+1Ck+1y, we know that the entries of Ck+1y are the coefficients of qk(z), that is, if Ck+1y=(α0,α1,…,αk)T, then qk(z)=α0+α1z+⋯+αkzk. Note that pk+1(z)=1-zqk(z) is the residual polynomial, so we have qk(z)=(1-pk+1(z))/z and qk(A)≈A-1 [4]; therefore, qk(A) becomes a preconditioner.In the next subsection, we present the preconditioned CMRH(m) with the polynomial preconditioner described above.


## 2.3. PCMRH(m) Algorithm

The polynomial preconditioned CMRH(m) (PCMRH(m)) method is given in the following, which is composed of two parts, in the first part, the polynomial preconditioner is constructed, and in the second part, the CMRH(m) method is applied to the new linear systems (Algorithm3).Algorithm 3:PCMRH(m).
Input:  the dimension ofKrylov subspace m; the prescribed stopping criterion ɛ;the degree of the polynomial preconditionerkk; the initial approximate solution x0;the residual vectorr0=b-Ax0.1. The procedure of constructing polynomial preconditionerqkk(1)  Let p(i)=i,i=1,2,…,n;determinei0so|(r0)|i0=∥r0∥∞;l1=r0(r0)i0;p(1)↔p(i0),C(1,1)=1(r0)i0(2)  for k=1,2,…,kkdou=Alkforj=1,2,…,kdoc=up(j),hj,k=c,up(j)=0forl=j+1,…,ndoup(l)=up(l)-c×(lj)p(l)endifk<ncomputei0so|(u0)|i0=∥u0∥∞;p(k+1)↔p(i0)hk+1,k=ui0;lk+1=u(u)i0lethj=(h1,j⋮hj,j),Cj=(c11⋯c1j0⋱⋮00cjj)compute(c1,j+1c2,j+1⋮cj+1,j+1)=hj+1,j-1(0c1j⋮cjj)-hj+1,j-1(Cjhj0)endend(3) Form the approximate solution xk=x0+Lkyk, where yk∈CkminimizesF(y)=∥H¯ky-(r0)i0e1∥.(4) Form the polynomial preconditioner qkk(z)Ckyk=(α0⋯αk-1)T;qkk(x)=α0+α1x+⋯+αk-1xk-1.2. Using the CMRH(m) algorithm for qkk(A)Ax=qkk(A)b.


## 3. Numerical Experiments

In this section, we give some numerical examples to compare the performance of PCMRH(m) with CMRH(m). All codes were written in Matlab 6.5 and run on AMD Athlon(tm) 64 X2 Dual Core Processor 2.20 GHz equipped with 1.5 G of RAM. In these numerical experiments, all the matrices are from Sadok's paper [2] we always take the initial vector x0=(0,0,…,0)T, b=(1,1,…,1)T, and m=20. We use the inequality norm(r)/norm(b)≤10-10 as stopping criteria. In the figures, solid lines(-) denote the CMRH(m) method, and dots(-.) denote the PCMRH(m) method. The x-axis denotes the number of iterations, and y-axis denotes the norm of the residual r: norm(r).Example 3.1.
The first example is taken from Brown [5]. The n×n matrix(3.1)A=(ɛ1-1ɛ1⋱⋱⋱-1ɛ1-1ɛ)
is used in three experiments, whose results are listed in Figures 1 and 2 and Table 1, respectively. In Figure 1, ɛ=0.1, n=40, and the degree of the preconditioned polynomial: kk=20. In Figure 2, ɛ=0.01, n=40, and kk=20. In Table 1, ɛ=0.01, n=100 and kk ranges from 1 to 20.
From Figures1 and 2, we can see that the PCMRH(m) method converges much faster than CMRH(m), when ɛ=0.1, n=40, the CMRH(20) method needs 107 restarts to make sure that the residual norm is 9.2925e-011, while PCMRH(20) can get 1.1143e-011 only in 3 restarts. When ɛ=0.01, n=40, the CMRH(20) method needs 840 restarts to make sure that the residual norm is 9.0075e-011, while PCMRH(20) can get 4.4771e-011 only by 6 restarts. The polynomial preconditioner is well done.
Letkk range from 1 to 20 and list PCMRH(20) in Table 1. We see that only when kk=1, the PCMRH(20) method does not converge, and when kk=19, PCMRH(20) has the best performance.Table 1kkNo. of restartsTime (sec)Residual normConvergence110002.75000.0011No21710.46909.4299e− 11Yes32890.78109.4596e− 11Yes41770.48409.9338e− 11Yes52360.65609.6960e− 11Yes6500.12509.6988e− 11Yes7720.20306.8890e− 11Yes8600.17208.5579e− 11Yes9750.21907.7776e− 11Yes10340.10905.8616e− 11Yes11400.12505.5309e− 11Yes12320.10906.1764e− 11Yes13280.07808.4069e− 11Yes14230.07803.6570e− 11Yes15230.09407.6101e− 11Yes16170.07804.8057e− 11Yes17180.07808.7496e− 12Yes18110.06207.8930e− 12Yes1970.04703.5654e− 12Yes20140.06201.3220e− 11YesFigure 1
ɛ=0.1, n=40.Figure 2
ɛ=0.01, n=40.Example 3.2.
The matrix(3.2)A=(111…11a111…11a1a21…11⋮⋮⋱⋱⋱⋮a1a2a3…an-11)withai=1+iɛ
is taken from [6]. In this experiment, we let ɛ=10-2, n=100, and kk=2.
From Figure3, we can also see that the PCMRH(20) method converges much faster than CMRH(20). The CMRH(20) method needs 317 restarts to make sure the residual norm is 9.9221e-011, while PCMRH(20) reaches 9.0007e-011 only in 37 restarts. So, we can conclude that the polynomial preconditioner makes CMRH(m) algorithm more efficient.Figure 3
ɛ=10-2, n=100.Example 3.3.
LetA=SDS-1∈ℛ1000×1000, where S is a bidiagonal matrix of the form(3.3)S=(10.910.9⋱⋱0⋱0.91)∈R1000×1000,
and D=diag(-10,-9,…,-1,1,2,…,990)∈ℛ1000×1000.
In this experiment, we compare the PCMRH(m) method with CMRH(m) and GMRES(m). In Figure4, we use the dot line to denote the GMRES(m) method.
It can be seen from Figure4 that the GMRES(20) method is stagnate since some eigenvalues of A are negative; it uses 1000 restarts or about 368 seconds to make the residual norm reach 1.8110e-006. Meanwhile, the CMRH(20) method uses 883 restarts or about 100 seconds to get the residual norm be 8.4761e-011; however, the PCMRH(20) method with kk=3 can reach 8.2944e-011 only by 481 restarts or about 57 seconds. These results are also listed in Table 2 and show that the polynomial preconditioner is effective.Table 2AlgorithmNo. of restartsTime (sec)Residual normConvergenceGMRES(20)1000368.03101.8110e − 6NoCMRH(20)883100.25008.4761e − 11YesPCMRH(20) (kk=3)48157.26508.2944e − 11YesFigure 4
A=SDS-1.Table3 lists the performance of PCMRH(20). It is seen that the PCMRH(20) method is convergent when the value of kk ranges from 1 to 20, but it can not conclude that the number of restarts and the time needed by PCMRH(20) decrease(increase) as kk increases(decreases). In fact, when kk=3, the PCMRH(20) method performs the best.Table 3kkNo. of restartsTime (sec)Residual normConvergence151860.48409.4131e − 11Yes249959.34409.6438e − 11Yes348157.26508.2944e − 11Yes462575.17209.7851e − 11Yes563977.79709.6175e−11Yes6934113.67109.7319e − 11Yes768084.71909.9634e − 11Yes864982.48409.8637e − 11Yes974794.92209.3972e − 11Yes10824104.37506.7137e − 11Yes11940118.85908.5893e − 11Yes1258579.57909.8813e − 11Yes1369093.14008.4271e − 11Yes1464089.28209.6628e − 11Yes1559984.95309.8641e − 11Yes16755104.54709.8972e − 11Yes1755381.84409.0119e − 11Yes18730104.39109.3932e − 11Yes1961190.29708.6073e − 11Yes2066097.04609.8196e − 11Yes


## 4. Conclusion

To accelerate the convergence, we take advantage of the CMRH process to construct a polynomial preconditioner for CMRH(m) algorithm. Numerical experiments show that the preconditioned CMRH(m) algorithm is more efficient than the CMRH(m) algorithm without the polynomial preconditioner. However, the degree of the polynomial used as preconditioner should be small, otherwise it possibly reduces the preconditioner's effectiveness.



---

*Source: 545470-2011-03-22.xml*