Giáo trình tính toán khoa học - Chương 3 pot - Pdf 17


62

Chương 3
GIẢI TÍCH MA TRẬN VÀ ĐẠI SỐ TUYẾN TÍNH

3.1 GIẢI TÍCH MA TRẬN
Để các bạn làm quen với các công cụ của Matlab trong đại số tuyến tính,
trước hết chúng tôi cần nhắc lại một số khái niệm về ma trận và các phép toán
trên ma trận.
3.1.1 Chuyển vị ma trận
Cho ma trận A=(a
ij
)
mxn
. Ma trận chuyển vị của A là ma trận A' =(a'
ij
)
nxm
sao
cho a'
ij
=a
ji
. Nếu A’=A thì A được gọi là một ma trận đối xứng.
Thí dụ 1.
Nếu
1 2 3 4
5 6 7 8
9 10 11 12
A

với các phần tử a
23
và a
12
của nó:
23 12
1 2 3
1 2 4 6
4 5 6 , ,
7 8 7 9
7 8 9
A M M
 
   
 
  
   
 
   
 
.
Định thức của ma trận A vuông cấp n, gọi là định thức cấp n, được định
nghĩa theo phương pháp qui nạp như sau:
- Nếu A là ma trận cấp 1 hay A=(a
11
), thì det(A)=a
11
;

63

det( )
A

a
11
det(M
11
)-a
12
det(M
12
)+ + (-1)
1+n
a
1n
det(M
1n
),
hay
1
1 1
1
det( ) ( 1) det( )
n
j
j j
j
A a M



c c c
A C
A A
c c c

 
 
 
 
 
 
 
. (3.2)
trong đó C=(c
ij
)
nxn
, với c
ij
=(-1)
i+j
det(M
ij
) gọi là phần phụ đại số của phần tử a
ij

của ma trận A.
Chú ý rằng với khái niệm về phần phụ đại số ta có:
1 1
det( ) .


Tính giá trị nhỏ nhất m trong các toạ độ của vector x và vị trí
k đạt min. Nếu x là ma trận thì kết quả là vector hàng gồm giá
trị min của các cột.
[M,k]=max(
x)
Tính giá trị lớn nhất M trong các toạ độ của vector x và vị trí
k đạt max. Nếu x là ma trận thì kết quả là vector hàng gồm giá
trị max của các cột.
mean(x)
Tính giá trị trung bình của các phần tử của vector x. Nếu x là
ma trận thì kết quả là vector hàng gồm giá trị trung bình của
các cột.
[y,k]=sort(x)

Sắp xếp lại các phần tử của x theo thứ tự tăng dần, kết quả trả
cho vector y. Vector k là vector số thứ tự cũ trong x của các
phần tử trong y. Nếu x là ma trận thì các cột của x được sắp
xếp tăng dần.
sum(x)
Tính tổng các phần tử của vector x. Nếu x là một ma trận thì
kết quả là 1 vector hàng, mà mỗi phần tử của vector là tổng
các phần tử của một cột tương ứng.
cumsum(x)
Cộng dồn các phần tử của vector x.
prod(x)
Tính tích các phần tử của vector x. Nếu x là một ma trận thì
kết quả là 1 vector hàng, mà mỗi phần tử của vector là tích các
phần tử của một cột tương ứng.
cumprod

>> x= [ 0 0 1 -2 3 4 5 -6 -7 -8];
>> [M,k]=max(x)
M =

66
5
k =
7
Chú ý k là vị trí đầu tiên của phần tử đạt max trong x.Tương tự ta có:
>> [m,p]=min(x)
m =
-8
p =
9
>> mean(x)
ans =
-10/9
>> sort(x)
ans =
-8 -7 -6 -2 0 0 1 3 4 5
>> [y,k]=sort(x)
y =
-8 -7 -6 -2 0 0 1 3 4 5
k =
10 9 8 4 1 2 3 5 6 7
>> D= sort(A)
D =
1 2 1
1 2 1
3 5 3

1
,
n
n n i i
i
x y x y x y x y x y

     

. (3.3)
3.1.5 Chuẩn của vector
Cho một hàm vector hàm f(x) xác định trên không gian vector V, ký hiệu là
x
. Nếu nó thoả mãn các tính chất sau đây với

x

V:
1) 00,0  xxx ,
2) x .

x , R,
3) yxyx  (gọi là bất đẳng thức tam giác),
thì
x
được gọi là chuẩn của vector x.

68
Từ định nghĩa về chuẩn vector, với mỗi số p>0 có một chuẩn loại p tương
ứng trong R

;
-
1 2
1
= x
n
x x x
   ;

-
 
1
2 2 2
2
1 2
2
x
n
x x x
   
.
Từ đó có 3 loại chuẩn của ma trận tương thích với chuẩn vector. Cho A là
một ma trận vuông A =( a
ij
)
n

n
, khi đó:
-



2
= max
j
j
A

, với

j
là các trị riêng của ma trận đối xứng A
T
A.
j


còn được gọi là các trị kì dị (single value) của ma trận A.
Các chuẩn của ma trận cũng thỏa mãn các tính chất 1)-3) của chuẩn vector.
Ngoài ra, giữa chuẩn loại p=0,1,2 của vector và chuẩn loại tương ứng của ma trận
thoả mãn tính chất sau:
p p p
Ax A x

.
Bất đẳng thức trên gọi là tính tương thích của chuẩn ma trận đối với chuẩn
vector. Chuẩn loại 2 cúa vector trong không gian R
n
:
2 2 2

24
>> Chuan1 =max(sum(abs(A)))
Chuan1 =
18
 Hàm NORM
Cú pháp:
norm(V,p)
Giải thích. Hàm NORM tính chuẩn loại p của ma trận và vector V.
- Nếu V là một vector p > 0 hoặc  inf thì
norm(V,p) = sum(abs(V).^p)^(1/p) : Chuẩn loại p;
norm(V) = norm(V,2) : Chuẩn Euclide;

70
norm(V,inf) = max(abs(V)) : Chuẩn loại vô cùng;
norm(V,-inf) = min(abs(V)).
- Nếu V là một ma trận thì p chỉ có thể 1, 2, inf hoặc 'fro' và
norm (V) = max(svd(V)) : Trị kì dị lớn nhất của V;
norm (V,2) = norm (V);
norm (V,1) = max(sum(abs(V)));
norm (V,inf) = max(sum(abs(V')));
norm (V,'fro')= sqrt(sum(diag(V'*V))) : chuẩn Frobenius.
3.1.6 Phân loại ma trận
Khi nghiên cứu các bài toán trong đại số tuyến tính, người ta thường chia
ma trận thành 2 loại:
- Ma trận lưu trữ được: Các ma trận mà các phần tử của chúng có thể lưu
trữ và xử lí được trong bộ nhớ của MTĐT. Kích thước của các ma trận này
thường không lớn lắm n = 0(10
3
).
- Ma trận thưa: Ma trận có kích thước rất lớn, nhưng phần lớn các phần tử

0
và m =
x
Ax
x
inf
071
thì dễ dàng chứng minh rằng M = A và nếu A không suy biến thì m =
1
1-
A

.
Khi đó, tỷ số
-1
A . A
M
m

.được gọi là số điều kiện của ma trận A và được ký
hiệu là cond(A).
Giả sử x là nghiệm đúng của hệ phương trình tuyến tính Ax = b và x+

x là
nghiệm của hệ phương trình xấp xỉ ( )
A x x b b
    


R : ta đều có cond(cA) = cond(A);
iv) Nếu D = diag


n
i
d
1
thì cond(D) =
max
min
i
i
d
d
.
Thí dụ 6. Giả sử A=diag


100
1
1,0 , khi đó det(A)=10
-100
nên A là ma trận gần
suy biến. Tuy nhiên do A=0,1E nên A là ma trận có số điều kiện rất tốt và hệ
phương trình Ax=1 dễ dàng giải được một cách chính xác là: x
i
=10,
1,100

Ma trận tích tensor Kronecker của ma trận x và y.
vander(x)
Ma trận Vandermonde của vector x
Ma trận ma phương là ma trận vuông cấp n sinh bởi các số 1, 2, 3, , n
2

các tổng của cột, các tổng hàng và các tổng đường chéo đều bằng nhau.
Thí dụ 8. Ma trận magic cấp 3:
magic(3) =
8 1 6
3 5 7
4 9 2
 
 
 
 
.
Ma trận Pascal là ma trận vuông cấp n chứa số của tam giác Pascal.
Thí dụ 9. Ma trận Pascal cấp 4:

73
pascal (4)=
1 1 1 1
1 2 3 4
1 3 6 10
1 4 10 20
 
 
 
 

Ma trận Vandermonde của vector x=(x
1
,x
2
, , x
n
) là ma trận:

74
vander(x) =















1

1
1
21


n
n
m m mn
a B a B a B
a B a B a B
a B a B a B
 
 
 
 
 
 
.
Bảng 3-3
Các hàm trích phần tử từ một ma trận.
Hàm Ý nghĩa
diag
(A,p)

Nếu A là 1 ma trận thì hàm sẽ tạo ra 1 vector cột từ đường chéo thứ
p của ma trận A. Nếu A là 1vector thì hàm sẽ tạo ra 1 ma trận vuông
có đường chéo thứ p là vector A, còn các phần tử khác đều bằng 0.
Mặc định của p là 0 (đường chéo chính).
triu
(A,p) Hàm sẽ tạo ra một ma trận cùng cỡ với ma trận A, có các phần tử từ
đường chéo thứ p trở lên được lấy từ A, còn các phần tử khác (phía
dưới đường chéo thứ p) đều bằng 0. Mặc định của p là 0.
tril
(A,p) Hàm sẽ tạo ra một ma trận cùng cỡ với ma trận A, có các phần tử từ

>> G = triu(A)
G =
1 2 3 4
0 6 7 8
0 0 11 12
>> H = triu(A,2)
H =
0 0 3 4

76
0 0 0 8
0 0 0 0
>> C = tril(A,-2)
C =
0 0 0 0
0 0 0 0
9 0 0 0

3.3 CÁC PHƯƠNG PHÁP GIẢI HỆ PHƯƠNG TRÌNH
TUYẾN TÍNH
Xét hệ phương trình đại số tuyến tính có dạng Ax=b, trong đó A là một ma
trận vuông cấp n và b là một vector cột n chiều. Để giải hệ phương trình Ax=b ta
có các phương pháp giải như sau:
3.3.1 Phương pháp Cramer
Nếu A là ma trận không suy biến thì hệ phương trình có duy nhất nghiệm
xác định bởi biểu thức:
det( )
; 1,
det( )
j

trên lại rất kém hiệu quả khi tính toán. Vì vậy việc giải hệ phương trình tuyến
tính trên MTĐT rất hiếm khi người ta sử dụng ma trận nghịch đảo A
-1
. Tuy nhiên
việc tính ma trận nghịch đảo rõ ràng là vẫn rất cần thiết bởi vì biểu diễn trên là
công thức “chính xác” để giải bài toán đại số tuyến tính. Nhiều khi chúng ta cần
quan tâm vector nghiệm x và các tính chất của nó.
3.3.3 Phương pháp Gauss
Xét hệ phương trình Ax =b trong trường hợp A là ma trận tam giác trên và
các hệ số trên đường chéo chính đều khác 0:
11 12 1, 1 1
1 1
22 2, 1 2
2 2
1 1
1, 1 1,

0 0 0
0 0 0
n n
n n
n n
n n n n
n n
nn
a a a a
x b

pháp thế ngược từ dưới lên theo công thức sau đây mà không phải tính một định
thức nào:
n
n
nn
b
x
a

,
1
, 1 1 , 2 2
- -
n
k kj j
j k
k k k k k k k kn n
k
kk kk
b a x
b a x a x a x
x
a a
 
   

 

,
với

a a


 

 

, với
2,
k n

 Phương pháp Gauss
Nếu A không phải một ma trận tam giác, có thể dùng các phép biến đổi sơ
cấp về hàng của ma trận để biến đổi dần ma trận mở rộng
A
=[A,b] của ma trận
A về dạng ma trận tam giác trên, rồi dùng phương pháp thế ngược nói trên để tính
vector x. Công việc đó được gọi là quá trình xuôi của phương pháp Gauss:
- Trước hết đặt A
(0)
=[A,b].
- Quá trình xuôi gồm n-1 bước. Tại bước thứ k=
1,1 n
thực hiện:
+ Đặt A
(k)
=A
(k-1)
;
+ Tại hàng i >k tính

a

đều khác không thì khi
kết thúc quá trình xuôi ta được ma trận A
(n-1)
có dạng tam giác trên.
 Phương pháp Gauss -Jordan (Còn gọi là phương pháp Gauss cải tiến
hay phương pháp Chọn trụ tối đại).
Trong quá trình tính toán của phương pháp Gauss, nếu gặp phần tử trụ


0
1

k
kk
a thì không thể tiến hành bước tiếp theo được; Hoặc nếu phần tử trụ có
giá trị quá nhỏ thì kết quả của phép chia cũng sẽ kém chính xác vì vậy cần phải
cải tiến phương pháp Gauss như sau:
- Đặt A
(0)
= [A,b].
- Quá trình xuôi gồm n bước. Tại bước thứ k=
n,1
thực hiện:
+ Tìm
 
1

k

, , 1
k
k
kj
kj
k
kk
a
a j k n
a


  
;

+ Tại hàng i >k tính
   




 
1 1
1
1
.
, , 1
k k
k k
ik kj

(n)
có các phần tử
trên đường chéo chính bằng 1. Có thể tính được số phép tính số học cần sử dụng
cho phương pháp Gauss (kể cả 2 quá trình xuôi và ngược) là:
 
3 2
G
1
n (4 9 7 )
6
N n n n
   .
3.3.4 Phương pháp phân tích QR và phương pháp phân tích LU
Để giải hệ phương trình tuyến tính trên MTĐT người ta thường phân tích
ma trận A thành tích của 2 ma trận, mỗi ma trận đó dễ nghịch đảo hơn ma trận A.
Do đó sẽ giảm nhẹ khối lượng tính toán. Có 2 phương pháp phân tích ma trận:
 Phương pháp phân tích QR
Phân tích ma trận A thành tích của ma trận trực giao Q (ma trận unitar) và
ma trận tam giác trên R. Ta có thể viết lại hệ phương trình gốc như sau:
Ax = QRx = b
Việc tính nghịch đảo của ma trận Q là tầm thường; Và khi đó hệ phương
trình trở thành:
Rx = Q
T
b
Hệ phương trình mới có dạng tam giác trên nên có thể dễ dàng giải bằng
quá trình thế ngược. Như vậy việc giải hệ không cần phải tính nghịch đảo của ma
trận A. Để xác định các ma trận Q và R có thể sử dụng hàm QR của Matlab có cú
pháp như sau:
[Q,R,P] =qr(A)

n-1
E
1
A= MA=U hay A =M
-1
U=LU.
Phép phân tích này thực hiện được nếu như tất cả các phần tử trụ


1
0 1
k
kk
a , k ,n

 
tức là các định thức con chính của ma trận A đều phải khác
không. Tuy nhiên trong khi tính toán ta có thể gặp trường hợp


1
0
k
kk
a


, giống
như phương pháp Gauss cải tiến, ta cần đổi chỗ hàng k với một hàng r nào đó ở
dưới. Vì vậy cũng cần lưu trữ các phép đổi chỗ bằng một ma trận giao hoán để sử

 Thủ tục Crout
Bước 1 . Đặt D
1
= A.
Bước k=2,3 n. Tại bước k-1 ma trận A đã được biến đổi thành:
11 12 1 1 1
21 22 2 1 2
1,1 1,2 -1, -1 1,
1 2 , -1
,
u k n
k n
k k k k k n
n n n k
u u u u
l u u u
l l u
l l l


  
 
 
 



1kk
T
kkk
Dl
uu

theo công thức:

1
k k
kk
l s
u


1 1

T
k k k k
D H l u
 
 
.
Bằng phép phân tích trên ta đã biểu diễn A thành tích của hai ma trận tam
giác dưới L và tam giác trên U. Vì vậy chúng ta có:
Ax = LUx = b
Đặt Ux = y, rồi giải hệ phương trình Ly = b. Vì L là tam giác dưới nên dễ
dàng giải bằng quá trình thế xuôi để tính y. Sau đó ta lại giải hệ phương trình

như trên hình vẽ:
R1 R3 R5

+
i
1
i
2
i
3

+

V1 -

R2 R4 - V2 Giả sử các điện trở và điện áp đã biết. Để tính các dòng điện trong mạch
điện ta sẽ sử dụng định luật Kitchhoff : tổng điện áp rơi trên mỗi mạch vòng phải
bằng không để mạch điện cân bằng. Ta có hệ phương trình:
- V
1
+ R
1
i
1
+ R
2
(i

3
+ V
2
= 0
Có thể viết lại hệ phương trình dưới dạng ma trận như sau:
1 2 2 1 1
2 2 3 4 4 2
4 4 5 3 2
0
0
0
R R R i V
R R R R R i
R R R i V
 
  
 
  
 
    
  
 
 
  
  
 
  
.
Sau đây là một chương trình Matlab cho bài toán trên:


y=b. Cần chú ý là s
ij
=0 nếu i>j và ma trận S
T
là ma trận
tam giác dưới. Nên có thể tính vector y như sau:
y
1
=
1
11

b
s

1
1
1
1 1 2 2 -1, -1
-
- -

k
k k j
j
k k k k k k
k
kk kk
b s y
b s y s y s y

k k k k k k k kn n
kk kk
b s x
y s x s x s x
s s
, k=n-1,n-2, 1.

84
Vấn đề còn lại là tìm cách tính ma trận S như thế nào. Từ phương trình
S
T
S=A ta có:
,
ij
1 1
n i
kj ki kj
ik
k k
a s s s s
 
 
 
a
ij
, với
i j

.
Từ đó suy ra:

i
ki kj ii
k
s a s s s


 
 
 
 
 

với j >i ; s
ij
= 0 với j <i.
Chú ý rằng nếu các phần tử trên đường chéo chính (s
ii
) của ma trận S đều
khác không thì det(A)=


2
1
ii
n
i
s
0 khi đó hệ phương trình có duy nhất nghiệm.
Phương pháp Cholesky tính toán được cả trong trường hợp s
ij

đạt được sai số cho phép.

85
Nếu B = q<1 thì phương pháp lặp đơn hội tụ với mọi vector xấp xỉ đầu
x
(0)
và có sai số của lời giải tại bước k được ước lượng theo công thức như sau:

     
1 0
1
k
k
q
x x x
q

  


hoặc
     
1
1
k k k
q
x x x
q




0
n
n
n
nn
n n n
b b b b
b
b b b
b b
B B
b b
b
b b b
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 

j j i
x b x b x


 
 
 
(3.7)
Ý tưởng của phương pháp Seidel là thông tin càng được khai thác sớm càng
tốt. Vì thế khi


1
k
i
x

vừa được tính xong, nó được được sử dụng ngay để tính


k
i
x
.
Nếu B =q<1 thì phương pháp lặp Seidel sẽ hội tụ với mọi vector xấp xỉ đầu x
(0)
và có thể sử dụng công thức sai số như của phương pháp lặp đơn.
 Phương pháp Gauss-Seidel
Ma trận vuông A được gọi là ma trận chéo trội nếu thoả mãn một trong hai
điều kiện:

a x a x b

 

với
i 1,n

được biến đổi thành
, i 1,n
ij
i
i j
ii ii
j i
a
b
x x
a a

   

. (3.8)
Khi đó, hệ có dạng x=Bx+g trong đó:
ij
0 khi
b
- khi
ij
ii
j i

ii
, rồi thay vào
(3.8) ta được:
, i 1,n
ij
i j i
jj
j i
a
z z b
a

   

. (3.10)
Hệ (3.10) có dạng z=Bz+g trong đó:
ij
ij
jj
0 khi j i
a
b
- khi j i
a







Nhờ tải bản gốc

Tài liệu, ebook tham khảo khác

Music ♫

Copyright: Tài liệu đại học © DMCA.com Protection Status