Giáo trình Phương pháp tính - Chương 3: Hệ phương trình đại số tuyến tính

n
−
−

0 0 · · · 0
−a21 0 · · · 0
· · · · · · · · · · · ·
−an1 −an2 · · · 0
 −

0 −a12 · · · −a1n
0 0 · · · −a2n
· · · · · · · · · · · ·
0 0 · · · −ann
 =
= D − L− U
Chú ý rằng do aii 6= 0, ∀i = 1, n nên detD 6= 0. Và như vậy tồn tại ma
trận nghịch đảo:
D−1 =

1
a11
0 · · · 0
0
1
a22
· · · 0
· · · · · · · · · · · ·
0 0 · · · 1
ann

Khi đó hệ
Ax = b⇐⇒ (D − L− U )x = b (3.16)
Bây giờ chúng ta sẽ xét một vài phương pháp để chuyển hệ phương
trình (3.1) về dạng x = Tx+ c.
Từ hệ (3.16) ta có Dx = (L + U )x + b. Do tồn tại D−1 nên
x = D−1(L+ U )x+D−1b. Ký hiệu Tj = D−1(L+ U ) và cj = D−1b. Khi
đó công thức lặp theo (3.12) sẽ có dạng
x(m) = Tjx(m−1) + cj, m = 1, 2, 3, . . . (3.17)
Phương pháp lặp dựa trên công thức lặp (3.17) được gọi là phương
pháp Jacobi. Dạng tường minh của công thức (3.17) như sau:
x
(m)
i =
1
aii
− i−1∑
j=1
aijx
(m−1)
j −
n∑
j=i+1
aijx
(m−1)
j + bi
 (3.18)
3.5 Phương pháp lặp 57
với i = 1, 2, . . ., n. Ta có
‖TJ‖∞ = ‖D−1(L + U )‖∞ = max
i=1,n
n∑
j=1,j 6=i
∣∣∣∣aijaii
∣∣∣∣ = max
i=1,n
n∑
j=1,j 6=i
|aij |
|aii| < 1
do A là ma trận đường chéo trội nghiêm ngặt. Vậy ‖TJ‖∞ < 1, nghĩa
là phương pháp Jacobi luôn hội tụ với mọi vectơ lặp ban đầu x(0).
Ví dụ 3.10. Xét hệ phương trình 10x1 + x2 − x3 = 7x1 + 10x2 + x3 = 8−x1 + x2 + 10x3 = 9
Với vectơ lặp ban đầu x(0) = (0, 0, 0)T , hãy tính vectơ x(3) và đánh
giá sai số của nó. Ta có
Tj =
 0.0 −0.1 0.1−0.1 0.0 −0.1
0.1 −0.1 0.0
 và cj =
 0.70.8
0.9

Do đó:
x(1) = Tjx(0) + cj =
 0.70.8
0.9
; x(2) = Tjx(1) + cj =
 0.710.64
0.89
;
x(3) = Tjx(2) + cj =
 0.7250.640
0.907
. Ta có ‖Tj‖∞ = 0.2. Vì vậy
‖x(3) − x‖∞ 6 0.21− 0.2‖x
(3) − x(2)‖∞ = 0.0043
Ví dụ 3.11. Hệ phương trình Ax = b cho bỡi
10x1 − x2 + 2x3 = 6
−x1 + 11x2 − x3 + 3x4 = 25
2x1 − x2 + 10x3 − x4 = −11
3x2 − x3 + 8x4 = 15
58 HỆ PHƯƠNG TRÌNH ĐẠI SỐ TUYẾN TÍNH
có nghiệm duy nhất x = [1, 2,−1, 1]T . Để chuyển từ hệ Ax = b về
dạng x = Tjx+ cj , ta biến đổi như sau
x1 = 0.100x2 − 0.200x3+ 0.600
x2 = 0.091x1+ 0.091x3 − 0.273x4+ 2.273
x3 = 0.200x1+ 0.100x2 + 0.100x4− 1.100
x4 = −0.375x2 + 0.125x3+ 1.875
Khi đó ma trận Tj và vectơ cj có dạng:
Tj =

0 0.100 −0.200 0
0.091 0 0.091 −0.273
−0.200 0.100 0 0.100
0 −0.375 0.125 0
 , cj =

0.600
2.273
−1.100
1.875

Chọn chuẩn vô cùng và ta có ‖Tj‖∞ = 12 < 1. Do đó phương pháp
lặp hội tụ. Chọn x(0) = [0, 0, 0, 0]T. Bảng sau đây cho chúng ta kết
quả tính toán sau 10 lần lặp.
m 1 2 3 4 5
x
(m)
1 0.6000 1.0473 0.9326 1.0152 0.9890
x
(m)
2 2.2727 1.7159 2.0533 1.9537 2.0114
x
(m)
3 −1.1000 −0.8052 −1.0493 −0.9681 −1.0103
x
(m)
4 1.8750 0.8852 1.1309 0.9739 1.0214
m 6 7 8 9 10
x
(m)
1 1.0032 0.9981 1.0006 0.9997 1.0001
x
(m)
2 1.9922 2.0023 1.9987 2.0004 1.9998
x
(m)
3 −0.9945 −1.0020 −0.9990 −1.0004 −0.9998
x
(m)
4 0.9944 1.0036 0.9989 1.0006 0.9998
Quá trình lặp dừng lại dựa theo đánh giá:
‖x(10)− x‖∞ 6 1/21− 1/2‖x
(10)− x(9)‖∞ = 8.0× 10−4 < 10−3
Trong khi sai số thực sự là ‖x(10) − x‖∞ = 0.0002.
3.5 Phương pháp lặp 59
Phương pháp lặp Jacobi được thể hiện trong Chương trình 3.8.
Đối số của chương trình gồm: N là cấp của ma trận, a là ma trận
hệ số cấp N ×(N+1 ), x0 là vectơ lặp ban đầu, eps là sai số (giá trị
mặc định là 10−6) và maxit là số lần lặp tối đa cho phép. Kết quả
trả về của chương trình là vectơ lặp x , sai số ss và số lần lặp thực
tế n .
Chương trình 3.8. - c3jacobi : Phương pháp Jacobi.
function [x,ss,n]=c3jacobi(N,a,x0,eps,maxit)
if nargin < 5, maxit = 100; end;
if nargin < 4, eps = 1.0E-6; end;
if nargin < 3, error('Hàm có ít nhất 3 đối số.'); end;
n=0;
for l=1:maxit
n=n+1;
for k=1:N
sum=0;
for j=1:N
if j∼=k, sum=sum+a(k,j)*x0(j);end;
end;
x(k)=(a(k,N+1)-sum)/a(k,k);
end;
ss=0;
for k=1:N
if abs(x(k)-x0(k))>ss
ss=abs(x(k)-x0(k));
end;
end;
if ss<eps, break; end;
for k=1:N, x0(k)=x(k); end;
end;
Trong công thức (3.18), để tính các toạ độ của vectơ lặp x(m),
60 HỆ PHƯƠNG TRÌNH ĐẠI SỐ TUYẾN TÍNH
chúng ta chỉ sử dụng các toạ độ của x(m−1). Tuy nhiên với
i > 1, x(m)1 , x
(m)
2 , . . . , x
(m)
i−1 đã được tính và xấp xỉ nghiệm chính xác
x1, x2, . . . , xi−1 tốt hơn x(m−1)1 , x
(m−1)
2 , . . . , x
(m−1)
i−1 . Do đó khi tính x
(m)
i
chúng ta nên sử dụng các giá trị vừa tính xong x(m)1 , x
(m)
2 , . . . , x
(m)
i−1. Ta
thu được
x
(m)
i =
1
aii
− i−1∑
j=1
aijx
(m)
j −
n∑
j=i+1
aijx
(m−1)
j + bi
 (3.19)
với i = 1, 2, . . . , n. Công thức (3.19) thường được gọi là công thức lặp
Gauss-Seidel. Bây giờ ta sẽ viết dạng ma trận của phương pháp
Gauss-Seidel.
Từ hệ phương trình (3.16) ta được (D−L)x = Ux+b. Ma trận D−L
cũng có ma trận nghịch đảo và do đó x = (D − L)−1Ux + (D − L)−1b.
Đặt Tg = (D − L)−1U và cg = (D − L)−1b. Khi đó công thức lặp có
dạng
x(m) = Tgx(m−1) + cg, m = 1, 2, 3, . . .
Ví dụ 3.12. Xét hệ phương trình trong ví dụ 3.11
10x1 − x2 + 2x3 = 6
−x1 + 11x2 − x3 + 3x4 = 25
2x1 − x2 + 10x3 − x4 = −11
3x2 − x3 + 8x4 = 15
Công thức lặp theo phương pháp Gauss-Seidel có dạng

x
(m)
1 = 0.100x
(m−1)
2 − 0.200x(m−1)3 + 0.600
x
(m)
2 = 0.091x
(m)
1 + 0.091x
(m−1)
3 − 0.273x(m−1)4 + 2.273
x
(m)
3 = 0.200x
(m)
1 + 0.100x
(m)
2 + 0.100x
(m−1)
4 − 1.100
x
(m)
4 = −0.375x(m)2 + 0.125x(m)3 + 1.875
Chọn x0 = [0, 0, 0, 0]T. Bảng sau đây cho chúng ta kết quả tính
3.5 Phương pháp lặp 61
toán sau 5 lần lặp.
m 1 2 3 4 5
x
(m)
1 0.6000 1.0300 1.0065 1.0009 1.0001
x
(m)
2 2.3272 2.0370 2.0036 2.0003 2.0000
x
(m)
3 −0.9873 −1.0140 −1.0025 −1.0003 −1.0000
x
(m)
4 0.8789 0.9844 0.9983 0.9999 1.0000
Ta thấy đến lần lặp thứ năm, nghiệm thu được bằng phương pháp
Gauss-Seidel tốt hơn nhiều so với phương pháp Jacobi.
Phương pháp lặp Gauss-Seidel được thể hiện trong Chương trình
3.9. Đối số của chương trình gồm: N là cấp của ma trận, a là ma
trận hệ số cấp N ×(N+1 ), x0 là vectơ lặp ban đầu, eps là sai số
(giá trị mặc định là 10−6) và maxit là số lần lặp tối đa cho phép.
Kết quả trả về của chương trình là vectơ lặp x , sai số ss và số lần
lặp thực tế n .
Chương trình 3.9. - c3seidel : Phương pháp Gauss-Seidel.
function [x,ss,n]=c3seidel(N,a,x0,eps,maxit)
n=0;
for l=1:maxit
n=n+1;
for k=1:N
sum=0;
for j=1:N
if j<k, sum=sum+a(k,j)*x(j);
elseif j>k, sum=sum+a(k,j)*x0(j);
end;
end;
x(k)=(a(k,N+1)-sum)/a(k,k);
end;
ss=0;
for k=1:N
62 HỆ PHƯƠNG TRÌNH ĐẠI SỐ TUYẾN TÍNH
if abs(x(k)-x0(k))>ss
ss=abs(x(k)-x0(k));
end;
end;
if ss<eps, break; end;
for k=1:N, x0(k)=x(k); end;
end;
3.6 BÀI TẬP
1. Sử dụng phương pháp phần tử trội giải các hệ phương trình sau
đây:
(a)

2x1 − 1.5x2 + 3x3 = 1
−x1 + 2x3 = 3
4x1 − 4.5x2 + 5x3 = 1
(b)
 2.13x1 + 3.45x2 − 6.21x3 = 1.450.43x1 + 4.24x2 − 5.05x3 = 2.232.67x1 − 1.13x2 + 3.27x3 = 3.21
(c)

x1 + x2 + x4 = 2
2x1 + x2 − x3 + x4 = 1
4x1 − x2 − 2x3 + 2x4 = 0
3x1 − x2 − x3 + 4x4 = −3
2. Dùng phương pháp Doolittle phân tích các ma trận sau thành
tích LU :
(a)
 4 1 −24 5 1
8 12 9
 (b)
 2 2 −1−1 2 1
−2 1 4

(c)

1 1 −3 2
−1 2 1 4
2 1 2 −2
2 2 −1 1

3. Sử dụng các phương pháp nhân tử LU (Doolittle) giải các hệ
phương trình sau:
3.6 Bài tập 63
(a)
 2x1 − 5x2 + 4x3 = 13x1 + 3x2 + 9x3 = 03x1 + 6x2 + 5x3 = 4
(b)
 2.2x1 + 0.3x2 + 0.2x3 = 1.50.3x1 + 3.4x2 + 0.2x3 = 2.40.2x1 + 0.2x2 + 4.1x3 = 3.2
(c)

x2 + x3 = 1
x1 − 2x2 − x3 = 0
x1 − x2 + x3 = −1
(d)

x1 + x2 − x3 + x4 = 1
x1 − x2 + 4x3 + 3x4 = 2
2x1 − x2 + 2x3 + 4x4 = 3
2x1 + x2 + 2x3 + 3x4 = 2
4. Cho A là ma trận vuông, ba đường chéo, cấp 10 với các hệ số
được xác định bỡi aii = 2, ai,i+1 = ai,i−1 = −1 với mọi i = 2, . . . , 9
và a11 = a10,10 = 2, a12 = a10,9 = −1. Cho b là vectơ cột cấp 10
được cho bỡi b1 = b10 = 2 và bi = 1 với mọi i = 2, . . . , 9. Hãy giải
hệ Ax = b sử dụng phương pháp nhân tử LU.
5. Tìm các giá trị của α để cho các ma trận sau đây là xác định
dương:
(a)
 α 1 −11 2 1
−1 1 4
 (b)
 1 α 2α 4 −1
2 −1 8
 (c)
 1 −1 α−1 3 1
α 1 5

6. Sử dụng phương pháp Choleski giải các hệ phương trình sau:
(a)
 2x1 − x2 = 2−x1 + 2x2 − x3 = 1− x2 + 2x3 = 2
(b)
 x1 + 3x2 − 2x3 = 13x1 + 4x2 − 2x3 = 4−2x1 − 2x2 + x3 = 3
(c)

4x1 + x2 + x3 + x4 = −1
x1 + 3x2 − x3 + x4 = 0
x1 − x2 + 2x3 = 1
x1 + x2 + 2x4 = 2
64 HỆ PHƯƠNG TRÌNH ĐẠI SỐ TUYẾN TÍNH
(d)

5.5x1 + 1.2x2 + 1.3x3 + 1.4x4 = 1.5
1.2x1 + 5.5x2 + 1.4x3 + 1.5x4 = 1.6
1.3x1 + 1.4x2 + 5.5x3 + 1.6x4 = 1.7
1.4x1 + 1.5x2 + 1.6x3 + 5.5x4 = 1.8
7. Tính các chuẩn ‖.‖1, ‖.‖∞ và số điều kiện theo các chuẩn một và
vô cùng của các ma trận sau:
(a)
(
5 −3
2 8
)
(b)
 3 1 −11 2 1
−1 1 4
 (c)
 2 −4 −12 2 1
−1 2 4

8. Tìm các giá trị của α > 0 và β > 0 để cho các ma trận sau đây
là đường chéo trội nghiêm ngặt:
(a)
 4 α 12β 5 4
β 2 α
 (b)
 3 2 βα 5 β
2 1 α
 (c)
 β 2 42 α β
4 −1 α

9. Sử dụng phương pháp Jacobi tìm nghiệm gần đúng của các hệ
phương trình sau với sai số 10−3, chọn chuẩn vô cùng:
(a)
 4x1 − x2 + x3 = 13x1 + 8x2 + 2x3 = 03x1 + 3x2 + 10x3 = 4
(b)

10x1 − x2 = 9
−x1 + 10x2 − 2x3 = 7
− 2x2 + 10x3 = 6
(c)

10x1 + 5x2 = 6
5x1 + 10x2 − 4x3 = 25
− 4x2 + 8x3 − x4 = −11
− x3 + 5x4 = −11
(d)

−4x1 + x2 + x3 = −2
x1 − 4x2 + x4 = −1
x1 + − 4x3 + x4 = 0
x2 + x3 − 4x4 = 1
10. Lặp lại bài tập 9 sử dụng phương pháp Gauss-Seidel.