Applied Mathematics and Computation 181 (2006) 1499–1504 www.elsevier.com/locate/amc
A modified simple iterative method for nonsymmetric algebraic Riccati equations arising in transport theory q Liang Bao
a,b
, Yiqin Lin
c,* ,
Yimin Wei
a,b
a
b
Institute of Mathematics, School of Mathematical Sciences, Fudan University, Shanghai 200433, PR China Key Laboratory of Mathematics for Nonlinear Sciences (Ministry of Education), Fudan University, Shanghai 200433, PR China c School of Mathematics and Computational Science, Sun Yat-Sen University, Guangzhou 510275, PR China
Dedicated to Prof. Zhi-hao Cao on the occasion of his 70th birthday
Abstract In this paper we propose a modified simple iterative method for nonsymmetric algebraic Riccati equations arising in transport theory. It is similar to the method proposed by Lu in [L.Z. Lu, Solution form and simple iteration of a nonsymmetric algebraic Riccati equation arising in transport theory, SIAM J. Matrix Anal. Appl. 26 (2005) 679–685]. The theoretical results and numerical experiments show that the modified version is more efficient than its counterpart. 2006 Elsevier Inc. All rights reserved. Keywords: Nonsymmetric algebraic Riccati equations; Simple iterative method; Minimal positive solution
1. Introduction In this paper, we are interested in iteratively solving the following algebraic Riccati equation arising in transport theory (see [2,3] and the references therein): XCX XE AX þ B ¼ 0;
ð1Þ
where A; B; C; E 2 Rnn are given by A ¼ D eqT ;
B ¼ eeT ;
C ¼ qqT ;
E ¼ D qeT .
ci , Here e = (1, 1, . . . , 1)T, q = (q1, q2, . . . , qn)T with qi ¼ 2x i
q Y. Lin is supported by Scientific Research Startup Foundation of Sun Yat-Sen University for Young Teacher under Grant 2005-340001131042. Y. Wei is supported by the National Natural Science Foundation of China under Grant 10471027 and Shanghai Education Committee. * Corresponding author. E-mail addresses:
[email protected] (L. Bao),
[email protected] (Y. Lin),
[email protected] (Y. Wei).
0096-3003/$ - see front matter 2006 Elsevier Inc. All rights reserved. doi:10.1016/j.amc.2006.02.042
1500
L. Bao et al. / Applied Mathematics and Computation 181 (2006) 1499–1504
(
1 D ¼ diagðd1 ; d2 ; . . . ; dn Þ with di ¼ cxi ð1þaÞ ; 1 D ¼ diagðd 1 ; d 2 ; . . . ; d n Þ with d i ¼ cxi ð1aÞ
and 0 < c 6 1, 0 6 a < 1, 0 < xn < < x2 < x1 < 1, n X ci ¼ 1; ci > 0; i ¼ 1; 2; . . . ; n. i¼1
It has been shown in [2] and [3] that (1) has positive solutions (in the componentwise sense). Since only the minimal positive solution is physically meaningful, some iterative methods have been developed for computing the minimal positive solution of (1), see [2] and [3]. Recently, Lu [4] has shown that the matrix equation (1) is in essence equivalent to a ‘‘vector equation’’ and developed a simple but efficient iterative procedure to compute the minimal positive solution of (1). In this paper we present a modified version of this simple iterative procedure and prove that the modified simple iterative procedure is more efficient than the original one. This paper is organized as follows. In Section 2, we review the simple iterative procedure proposed by Lu [4] and present the modified simple iterative procedure. Some theoretical results are given in this section. Sections 3 and 4 give some numerical examples and conclusions, respectively. Let A = (ai,j), B = (bi,j) be n · n matrices with real entries. The Hadamard product of A and B is defined by A B = (ai,j Æ bi,j). 2. Simple iteration and modified simple iteration We rewrite (1) as DX þ XD ¼ ðXq þ eÞðqT X þ eT Þ and let v T ¼ qT X þ e T .
u ¼ Xq þ e;
ð2Þ
It is shown in [2] and [3] that the positive solutions of (1) are of the following form: X ¼ T ðuvT Þ ¼ ðuvT Þ T ; where
T ¼ ðti;j Þ ¼ X ¼ ðxi;j Þ;
ð3Þ
1 ; di þ d j
u ¼ ðu1 ; u2 ; . . . ; un Þ
T
T
v ¼ ðv1 ; v2 ; . . . ; vn Þ .
and
To find the minimal positive solution of (1), we need to find proper positive vectors u and v in (3). Substituting (3) into (2), we obtain u ¼ u ðPvÞ þ e; ð4Þ v ¼ v ðP~ uÞ þ e; where P ¼ ðpi;j Þ ¼
qj ; di þ d j
P~ ¼ ð~ pi;j Þ ¼
qj . dj þ d i
To get the positive vectors u and v from (4), Lu [4] has defined a simple iterative scheme 8 ðkþ1Þ > ¼ uðkÞ ðPvðkÞ Þ þ e;
: ð0Þ u ¼ vð0Þ ¼ 0 and the following results are given.
ð5Þ
L. Bao et al. / Applied Mathematics and Computation 181 (2006) 1499–1504
1501
Lemma 1 [4]. The minimal positive solution X* of (1) can be computed by X ¼ T ðu ðv ÞT Þ; where (u*, v*) is the limit of {(u(k), v(k))} defined by the simple iteration (5). Remark 2. Lu has showed that the sequence {(u(k), v(k))} is strictly monotonically increasing and bounded above. Before computing v(k+1), we have obtained u(k+1). It should be a better approximation to u* than u(k). We replace u(k) with u(k+1) in the second equation of (5), and obtain a modified simple iterative scheme 8 uðkþ1Þ ¼ ^ uðkÞ ðP^vðkÞ Þ þ e; > <^ ð6Þ ^vðkþ1Þ ¼ ^vðkÞ ðP~ ^ uðkþ1Þ Þ þ e; k ¼ 0; 1; . . . ; > : ð0Þ ð0Þ ^ u ¼ ^v ¼ 0: Then, we have the following result. Lemma 3. For all 0 < c 6 1 and 0 6 a < 1, the sequence fð^uðkÞ ; ^vðkÞ Þg defined by (6) is strictly monotonically increasing, bounded above and thus converging. Proof. It is easy to show by the induction that fð^ uðkÞ ; ^vðkÞ Þg is strictly monotonically increasing. To show that ðkÞ ðkÞ fð^ u ; ^v Þg is bounded above, we observe that ^uð0Þ < uð1Þ ;
^vð0Þ < vð1Þ
and
^ uð1Þ < uð2Þ ;
^vð1Þ < vð2Þ ;
^ð2Þ < uð3Þ ; ^vð2Þ < vð4Þ . Therefore, the statement it follows that u ^ uðiÞ < uð2i1Þ ;
^vðiÞ < vð2iÞ
ð7Þ
is true for i = 2. We now assume that (7) is true for i = k P 2. By (5) and (6) we have ^ uðkÞ ðP^vðkÞ Þ þ e < uð2kÞ ðPvð2kÞ Þ þ e ¼ uð2kþ1Þ uðkþ1Þ ¼ ^ and ^vðkþ1Þ ¼ ^vðkÞ ðP~ ^ uðkþ1Þ Þ þ e < vð2kþ1Þ ðP~ uð2kþ1Þ Þ þ e ¼ vð2kþ2Þ . We have thus proved that (7) is true for i = k + 1. Hence, by the principle of mathematical induction, (7) is true for all k P 2. Since {(u(k), v(k))} is bounded above (see [4]), so is fð^uðkÞ ; ^vðkÞ Þg. The sequence fð^uðkÞ ; ^vðkÞ Þg is now well defined, monotonically increasing, and bounded above. h Now we present our main result. Theorem 4. The minimal positive solution X* of (1) can be computed by X ¼ T ð^ u ð^v ÞT Þ; uðkÞ ; ^vðkÞ Þg defined by the modified simple iteration (6). where ð^ u ; ^v Þ is the limit of fð^ Proof. It is easy to show by the induction that uðkÞ < ^uðkÞ and vðkÞ < ^vðkÞ hold for all k P 3. By (7), we obtain uðkÞ < ^ uðkÞ < uð2k1Þ ;
vðkÞ < ^vðkÞ < vð2kÞ
ðk P 3Þ.
So the sequence fð^ uðkÞ ; ^vðkÞ Þg has the same limit (u*, v*) as the sequence {(u(k), v(k))}. By Lemma 1, we complete the proof. h Remark 5. From the proof of Theorem 4, we can see that the modified simple iteration (6) is more efficient than the simple iteration (5) since uðkÞ < ^ uðkÞ and vðkÞ < ^vðkÞ hold for all k P 3.
1502
L. Bao et al. / Applied Mathematics and Computation 181 (2006) 1499–1504
3. Numerical experiments In this section, we present some numerical experiments to illustrate the performance of the modified simple iterative method (6). Let SI and MSI denote the simple iterative method (5) and the modified simple iterative method (6), respectively. In the following examples, we compare MSI with SI. All numerical experiments are performed on an AMD 1.4 GHz PC with main memory 512 MB. The stopping criterion for two methods is kRk k1 6 1:0e 12; kR0 k1 where Rk denotes the kth residual matrix given by XkCXk XkE AXk + B and R0 = B = eeT.
Table 3.1 Example 1, MSI vs. SI (n = 64) c
a 8
6
Method
Iterations
CPU
Residual
10
1 10
SI MSI
17251 13572
27.94 20.38
9.9960e13 9.9969e13
0.001
0.999
SI MSI
746 588
1.08 0.80
9.8010e13 9.9428e13
0.01
0.99
SI MSI
246 195
0.34 0.27
8.9998e13 9.4635e13
0.5
0.5
SI MSI
19 17
0.03 0.02
8.3363e13 6.2768e13
0.85
0.1
SI MSI
8 7
0.01 0.01
2.3443e14 3.4609e13
0 SI MSI _2
log(||Rk||/||R0||)
_4
_6
_8
_10
_12
_14
0
2000
4000
6000
8000
10000
12000
14000
Number of iterations
Fig. 1. Comparison for Example 1 (a = 108, c = 1 106).
16000
18000
L. Bao et al. / Applied Mathematics and Computation 181 (2006) 1499–1504
1503
0 SI MSI _2
log(||Rk||/||R0||)
_4
_6
_8
_10
_12
_14
0
2000
4000
6000
8000
10000
12000
14000
16000
18000
Number of iterations
Fig. 2. Comparison for Example 2 (a = 108, c = 1 106).
Table 3.2 Example 2, MSI vs. SI (n = 128) a
c
Method
Iterations
CPU
Residual
108
1 106
SI MSI
17252 13574
149.34 115.65
9.9711e13 9.9985e13
0.001
0.999
SI MSI
746 588
6.42 4.97
9.8318e13 9.9843e13
0.01
0.99
SI MSI
246 195
2.11 1.68
9.0427e13 9.4950e13
0.5
0.5
SI MSI
19 17
0.16 0.13
8.3379e13 6.2264e13
0.85
0.1
SI MSI
8 7
0.07 0.06
2.3471e14 3.4595e13
Example 1. We consider (1) for n = 64. As in Example 5.2 in [1], the constants ci and xi are given by a numerical quadrature formula on the interval [0, 1], which is obtained by dividing [0, 1] into n/4 subinterval of equal length and applying a Gauss–Legendre quadrature with four nodes to each subinterval. We test five values of (a, c) taken to be (108, 1 106), (0.001, 0.999), (0.01, 0.99), (0.5, 0.5) and (0.85, 0.1). In Table 3.1, we give the results obtained by the two methods. We list the number of iterations, the CPUtime(seconds) required for convergence and the relative residual kRkk1/kR0k1. We can see that the number of iterations for MSI is about 0.75 time of that for SI in case (a, c) = (108, 1 106). See Fig. 1 for a graph of the convergence. Example 2. We now do the same experiment as in the previous example except that n is 128 instead of 64. See Fig. 2 for the computational results. From Table 3.2, we can see that the number of iterations is almost the same for n = 128 as that for n = 64.
1504
L. Bao et al. / Applied Mathematics and Computation 181 (2006) 1499–1504
4. Conclusion In this paper we propose a modified simple iterative method for the minimal positive solution of nonsymmetric Riccati equations arising in transport theory. Numerical experiments show the effectiveness of the proposed method. References [1] C.H. Guo, A.J. Laub, On the iterative solution of a class of nonsymmetric algebraic Riccati equations, SIAM J. Matrix Anal. Appl. 22 (2000) 376–391. [2] J. Juang, Existence of algebraic matrix Riccati equations arising in transport theory, Linear Algebra Appl. 230 (1995) 89–100. [3] J. Juang, W.W. Lin, Nonsymmetric algebraic Riccati equations and Hamiltonian-like matrices, SIAM J. Matrix Anal. Appl. 20 (1999) 228–243. [4] L.Z. Lu, Solution form and simple iteration of a nonsymmetric algebraic Riccati equation arising in transport theory, SIAM J. Matrix Anal. Appl. 26 (2005) 679–685.