# A Modification of Minimal Residual Iterative Method to Solve Linear Systems

Authors: Xingping Sheng; Youfeng Su; Guoliang Chen
Journal: Mathematical Problems in Engineering
 (2009)

---

## Abstract

We give a modification of minimal residual iteration (MR), which is 1V-DSMR to solve the linear systemAx=b. By analyzing, we find the modifiable iteration to be a projection
technique; moreover, the modification of which gives a better (at least the same) reduction of the residual
error than MR. In the end, a numerical example is given to demonstrate the reduction of the residual
error between the 1V-DSMR and MR.


---

## Body

## 1. Introduction

One of the important computational problems in the
applied science and engineering is the solution ofn×n nonsingular
linear systems of equations:(1.1)Ax=b,where A∈Rn×n is a symmetric
positive definite matrix (referred to as an SPD matrix), b∈Rn is given, and x∈Rn is unknown. To
solve this problem, usually an iterative method is spurred by demands,  which
can be found in excellent papers [1, 2]. Most of the existing practical
iterative techniques for solving larger linear systems of  (1.1)
utilize a projection process in one way or another; see, for example, [3–9].Projection techniques are presented in different forms
in many other areas of scientific computing, which can be formulated in
abstract Hilbert functional spaces and finite element spaces. Furthermore,
projection techniques are the process in which one attempts to solve a set of
equations by solving each separate equation by correcting so that it is small
in some norm. The idea of projection process is to extract an approximate
solution to (1.1) from a subspace of Rn. Denote 𝒦 and ℒ the search
subspace and the constraints subspace, respectively. Let m be their
dimension and x0∈Rn be an initial
guess to the solution of (1.1). A projection method onto the subspace 𝒦 and orthogonal
to ℒ is a process
which finds an approximation solution x∈Rn to (1.1) by
imposing the Petrov-Galerkin conditions that x belongs to the
affine space x0+𝒦 and the new
residual vector is orthogonal to ℒ, that is,(1.2)findx∈x0+𝒦,suchthatb−Ax⊥ℒ.From this point of view, the basic iterative methods
for solving (1.1), such as Gauss-Seidel Iteration (GS), Steepest Descent
Iteration (SD), Minimal Residual Iteration (MR), and Residual Norm Steepest Descent
Iteration (RNSD), all can be viewed as a special case of the
projection techniques.In [2], Ujević obtained a new
iterative method for solving (1.1), which is considered as a modification of
Gauss-Seidel method. In [10], Jing and Huang pointed that this iterative method
is also a projection process and named this method as “one-dimensional double
successive projection method” (referred to as 1D-DSPM). In the same paper
[10], the authors obtained another iterative method, which is named as
“two-dimensional double successive projection method” (referred to as
2D-DSPM). The theory indicates that 2D-DSPM gives a better reduction of error
than 1D-DSPM.


## 2. Notations and Preliminaries

In this paper, we will consider the following linear
system of equations:(2.1)Ax=b,where A∈Rn×n is not a
symmetric but a positive definite matrix of order n,b∈Rn is a given
element, and x unknown. For
linear systems (2.1), we can use the classical minimal residual iteration to
solve, which can be found in [11]. Here, we will give a modification of minimal
residual iteration to solve linear system (2.1). We call the modification
method as one vector double successive MR (abbreviated 1V-DSMR) and compare
reduction of the residual error at step k+1 between the
modification iteration and original MR. Hence, we find that the modification
iteration gives a better reduction of the residual error than the original MR.
〈x,y〉=yTx denotes a
vector inner product between the vector x,y∈Rn.We define the inner products as(2.2)a=〈Av1,Av1〉,c=〈Av1,Av2〉=〈Av2,Av1〉,d=〈Av2,Av2〉,(2.3)p=〈b−Axk,Av1〉=〈rk,Av1〉,q=〈b−Axk,Av2〉=〈rk,Av2〉.In this subsection, we will recall minimal residual
iterative method and give some properties of this iteration.For the linear system (2.1), we can use the following
algorithm which is called minimal residual iteration, viewed in [11].Algorithm 2.1.
(1) Choose an initial guess solutionx0∈Rn to (2.1), k:=0.
(2) Calculate(2.4)r0=b−Ax0,α0=〈Ar0,r0〉〈Ar0,Ar0〉.
(3) Ifrk=0, then stop; else,
(4) calculate(2.5)xk+1=xk+αkrk,rk+1=b−Axk+1,αk+1=〈Ark+1,rk+1〉〈Ark+1,Ark+1〉,k:=k+1.
(5) Go to step (3).The minimal residual iteration can be interpreted with
projection techniques. Here we represent the principles of this method in our
uniform notation as follows:(2.6)xk+1=xk+αv1,where α=p/a.If we choose𝒦=span{v1} and ℒ=span{Av1}, then (1.2) turns to find(2.7)xk+1∈xk+𝒦,suchthatb−Axk+1⊥ℒ,where xk+1=xk+αv1.Equation (2.7) can be represented in terms of inner product
as(2.8)〈b−Axk+1,Av1〉=0,which is(2.9)〈b−Axk−Aαv1,Av1〉=〈b−Axk,Av1〉−α〈Av1,Av1〉=p−αa=0,giving rise to α=p/a, which is the same as in (2.6).If we choose a specialv1=rk, then (2.6) is the minimal residual iteration (MR); up
to now, it is clear that MR is a special case of projection methods.For the MR, we have the following property.Lemma 2.2.
Let{xk} and {rk} be generated by Algorithm
2.1, then we have(2.10)∥rk+1∥2=∥rk∥2−p2a,where a,p are defined in
(2.2) and (2.3), respectively.Proof.
Using (2.6), we obtain(2.11)rk+1=b−Axk+1=b−Axk^−αAv1=rk−αAv1,where v1=rk. Hence we have(2.12)∥rk+1∥2=〈rk+1,rk+1〉=〈b−Axk−αAv1,b−Axk−αAv1〉=∥rk∥2−α〈rk,Av1〉−α〈Av1,rk〉+α2〈Av1,Av1〉=∥rk∥2−2αp+α2a=∥rk∥2−p2a.
From Lemma2.2, easily, we can get the reduction of
the residual error of MR as follows:(2.13)∥rk∥2−∥rk+1∥2=p2a.


## 3. An Interpretation of 1V-DSMR with Projection Technique

In this
section, we will give the modification of the minimal residual iterative
method, which is abbreviated to 1V-DSMR; we can present this method in our
uniform notation as follows:(3.1)x^k+1=xk^+α^v1+β^v2,where α^=p/a and β^=(aq−cp)/ad.We will have a two-step investigation of 1V-DSMR.The first step is to choose𝒦1=span{v1} and ℒ1=A𝒦1=span{Av1}, then it turns into the proceeding of MR, so we have α^=p/a.The next step is a similar way to choose𝒦2=span{v2} and ℒ2=A𝒦2=span{Av2}; denote xk+1=xk^+α^v1 and (1.2)
turns to find(3.2)x^k+1∈xk+1+𝒦2,suchthatb−Ax^k+1⊥ℒ2,where x^k+1=xk+1+β^v2.Equation (3.2) can be represented as in terms of inner
product as(3.3)〈b−Ax^k+1,Av2〉=0,which is(3.4)〈b−Axk+1−β^Av2,Av2〉=〈b−Axk−α^Av1−β^Av2,Av2〉=〈b−Axk,Av2〉−α^〈Av1,Av2〉−β^〈Av2,Av2〉=q−α^c−β^d=0.This gives rise to β^=(aq−cp)/ad, which is the same as in (3.1).If we choose a specialv1=rk, then (3.1) is a modification of the minimal residual
iteration, which is named as 1V-DSMR; up to now, it is clear that 1V-DSMR is
also a special case of projection methods.As 1V-DSMR, we have the following relation of residual
errors.Theorem 3.1.
Let{x^k+1} generated by (3.1) and r^k+1=b−Ax^k+1, then we have(3.5)∥r^k+1∥2=∥r^k∥2−1a2d[(ad+c2)p2+a2q2−2acpq],where
r^k
is the same  as
in Lemma 2.2, and a,c,d,p,q are defined in
(2.2) and (2.3).Proof.
Using (3.1), we have(3.6)r^k+1=b−Ax^k+1=b−Axk−α^Av1−β^Av2=r^k−α^Av1−β^Av2.By deduction, we get(3.7)∥r^k+1∥2=〈r^k−α^Av1−β^Av2,r^k−α^Av1−β^Av2〉=〈r^k,r^k〉−2α^〈r^k,Av1〉−2β^〈r^k,Av2〉+α^2〈Av1,Av1〉+2α^β^〈Av1,Av2〉+β^2〈Av2,Av2〉=∥r^k∥2−2α^p−2β^q+α^2a+2α^β^c+β^2d.If we substitute α^=p/a and β^=(aq−cp)/ad into (3.7),
then we obtain(3.8)∥r^k+1∥2=∥r^k∥2−1a2d[(ad+c2)p2+a2q2−2acpq].From Theorem 3.1, we also get a
reduction of residual error of 1V-DSMR as follows:(3.9)∥r^k∥2−∥r^k+1∥2=1a2d[(ad+c2)p2+a2q2−2acpq].Next we  will depict the comparison results with
respect to residual error reduction between 1V-DSMR and MR.Theorem 3.2.
1V-DSMR gives a better (at least the same) reduction
of the residual error than MR.Proof.
From the equalities (2.13) and (3.9), we
have(3.10)(∥r^k∥2−∥r^k+1∥2)−(∥rk∥2−∥rk+1∥2)=1a2d[(ad+c2)p2+a2q2−2acpq]−p2a=1a2d(cp−aq)2≥0,which proves the assertion of Theorem
3.2.Theorem3.2 implies that the residual error of MR is bigger
than that of 1V-DSMR at k+1 step if the
residual vectors rk and r^k at the kth iteration
are equal to each other.


## 4. A Particular Method of 1V-DSMR

In this
section, particularv1 and v2 will be   chosen,
and an algorithm to interpret 1V-DSMR is obtained.Since 1V-DSMR is a modification of minimal residual iteration, takev1=r^k. In general, v2 may be chosen
in different ways. Here, we choose a particular v2=xk−1^, then from (3.1), each step of 1V-DSMR is as
follows:(4.1)x^k+1=xk^+α^r^k+β^x^k−1,fork=0,1,…,where r−1^=0 and(4.2)α^=pa=〈Ar^k,r^k〉〈Ar^k,Ar^k〉,β^=aq−cpad=〈Ar^k,Ar^k〉〈r^k,Ax^k−1〉−〈Ar^k,Ax^k−1〉〈r^k,Ar^k〉〈Ar^k,Ar^k〉〈Ax^k−1,Ax^k−1〉.In this case, the first step of
1V-DSMR is the same as minimal residual iteration.The above results provide the following algorithm of
1V-DSMR.Algorithm 4.1 (A particular implementation of 1V-DSMR in a generalized way).
(1) Choose an initial guess solutionx^0 of (2.1), k:=0.
(2) Calculate(4.3)r^0=b−Ax^0,α^0=〈Ar^0,r^0〉〈Ar^0,Ar^0〉.
(3) Calculate(4.4)x^1=x^0+α0^r^0=x^0+〈Ar^0,r^0〉〈Ar^0,Ar^0〉r^0,r^1=b−Ax^1=r^0−〈Ar^0,r^0〉〈Ar^0,Ar^0〉Ar^0.
(4) Ifr^k=0, then stop; else,
(5) calculate(4.5)α^k=〈Ar^k,r^k〉〈Ar^k,Ar^k〉,β^k=〈Ar^k,Ar^k〉〈r^k,Ax^k−1〉−〈Ar^k,Ax^k−1〉〈r^k,Ar^k〉〈Ar^k,Ar^k〉〈Ax^k−1,Ax^k−1〉,x^k+1=x^k+αk^r^k+βk^x^k−1,r^k+1=b−Ax^k+1,k:=k+1.
(6) Goto step (4).About Algorithm4.1, we have the following basic
property.Theorem 4.2.
IfA is a positive
matrix, then the sequence of the iterations of Algorithm 4.1 converges to the
solution of the linear system Ax=b.Proof.
From the equality (3.9), we get(4.6)∥r^k∥2−∥r^k+1∥2≥p2a=〈r^k,Ar^k〉2〈Ar^k,Ar^k〉≥λmin2∥r^k∥4λmax2∥r^k∥2=λmin2λmax2∥r^k∥2.
This means that the sequence∥r^k∥2 is a decreasing
and bounded one. Thus, the sequence in question is convergent implying that the
left-hand side tends to zero. Obviously, ∥r^k∥2 tends to zero,
and the proof is complete.


## 5. Numerical Examples

In this
section, we use examples to further examine the effectiveness and show the
advantages of 1V-DSMR over MR.We compare the numerical behavior of Algorithm4.1
with Algorithm 2.1. All the tests are performed by MATLAB 7.0. Because of the
influence of the error of roundoff, we regard the matrix A as zero matrix
if ∥A∥<10−10.For convenience of comparison, consider the
two-dimensional partial differential equations on the unit square regionΩ=[0,1]×[0,1] of the
form(5.1)−(pux)x−(quy)y+rux+(ru)x+suy+(su)y+tu=f,where p,q,r,s, and t are all given
real valued function of x and y, which are as follows: p=e−xy,q=exy,r=20(x+y),s=10(x+y), and t=1/(1+x+y).Here, we use a five-point finite difference scheme to
discretize the above problem with a uninform grid of mesh spacingΔx=Δy=1/m in x and y directions,
respectively; we can obtain a matrix of order m×m as m varies, which
is called PDE matrix. Now, we take m=30, then we get a matrix, which is called
PDE900, and denoted by P. It is easy to check that P is real
unsymmetrical and nonsingular.If we take the coefficient matrixA of linear
system (2.1) as P, the right vector b=(1,1,…,1)T, and initial iterative vector x0=b, then we use Algorithms 2.1 and 4.1 to
compute the linear system (2.1), respectively. The comparison results between
MR and 1V-DSMR are shown in Figure 1.Figure 1
Comparison of convergence curve of residual norm (1.1) between Algorithm 2.1 and Algorithm 4.1.From Figure1,  we can see that the convergence
velocity of Algorithm 4.1 is always faster than that of Algorithm 2.1. In fact,
when we use Algorithm 2.1 to compute the linear system (2.1), we only need to
iterate 814 steps, and the residual norm is ∥r814∥≤9.8625×10−11. While using Algorithm 4.1, we only need to
iterate 647 steps, and the residual norm is ∥r647∥≤9.8465×10−11.



---

*Source: 794589-2009-04-26.xml*