›
›

Bài giảng Phương pháp số - Bài 3: Ma trận và hệ phương trình tuyến tính

pdf 44 trang ngocly 80 Free
Download
Bạn đang xem 20 trang mẫu của tài liệu "Bài giảng Phương pháp số - Bài 3: Ma trận và hệ phương trình tuyến tính", để tải tài liệu gốc về máy bạn click vào nút DOWNLOAD ở trên

Tài liệu đính kèm:

  • pdfbai_giang_phuong_phap_so_bai_3_ma_tran_va_he_phuong_trinh_tu.pdf

Nội dung text: Bài giảng Phương pháp số - Bài 3: Ma trận và hệ phương trình tuyến tính

  1. BÀI 3 MA TRẬN VÀ HỆ PHƯƠNG TRÌNH TUYẾN TÍNH
  2. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (1) 1. HỆ PHƯƠNG a11x1 a12x2 a1n xn b1 TRÌNH TUYẾN TÍNH a 21x1 a 22x2 a 2n xn b2 gồm m phương trình n ẩn là một hệ cĩ dạng a m1x1 a m2x2 a mmxn bm Nếu đặt a1j b1 x1 a11 a12 a1n a 2j b x a a a v ,(j 1,2, ,n),b 2 , x 2 và A 21 22 2n , j .   a mj bm x n a m1 a m2 a mn thì hệ trên cịn cĩ thể viết ở dạng vectơ cột x1v1 +x2v2 + + xnvn = b hay dạng phương trình ma trận Ax = b PHƯƠNG PHÁP SỐ-Bài 3 2
  3. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (2) 2. HỆ DẠNG TAM GIÁC TRÊN VÀ CÁCH GIẢI Hệ dạng tam giác trên là hệ cĩ dạng a11x1 + a12x2 + + a1nxn = b1 trong đĩ a11, a22, ,ann ≠ 0 a22x2 + + a2nxn = b2 Cách giải: giải ngược từ dưới lên annxn = bn void heTGiac(vector > a, vector &x) { unsigned n = a.size(); vector y(n, 0); // y co n phan tu 0 y[n-1] = a[n-1][n]/a[n-1][n-1]; for (int i = n-2; i >= 0; i ) { double tong = 0.; for(unsigned j = i+1; j <= n-1; j++) tong = tong + a[i][j] * y[j]; y[i] = (a[i][n] - tong) / a[i][i]; } x = y; } PHƯƠNG PHÁP SỐ-Bài 3 3
  4. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (3) 3. MỘT SỐ ĐỊNH LÍ VỀ NGHIỆM Định lí 4.1: Hệ Ax = b cĩ nhiều nhất một nghiệm (tức là, nghiệm là duy nhất nếu tồn tại ) nếu và chỉ nếu hệ thuần nhất tương ứng Ax = 0 chỉ cĩ nghiệm “tầm thường” x= 0. Định lí 4.2: Bất kì hệ PTTT thuần nhất nào với số phương trình ít hơn số ẩn đều cĩ nghiệm khơng tầm thường (khác 0) Định lí 4.3 Nếu A là một ma trận cấp m × n và hệ Ax = b cĩ nghiệm với mọi vectơ m chiều b, thì m ≤ n Định lí 4.4 Cho A là ma trận cấp n × n. Các khẳng định sau đây là tương đương: (i) Hệ thuần nhất Ax = 0 chỉ cĩ nghiệm tầm thường x = 0. (ii) Với mọi vế phải b, hệ Ax = b luơn cĩ nghiệm. (iii) A là ma trận khả nghịch PHƯƠNG PHÁP SỐ-Bài 3 4
  5. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (4) 4. LỜI GIẢI BẰNG SỐ CỦA HỆ PTTT Ax = b Xét các hệ PTTT cĩ số phương trình đúng bằng số ẩn. Nếu A khả nghịch thì hệ cĩ duy nhất nghiệm đối với mọi b Cĩ hai loại phương pháp giải: Trực tiếp: là các phương pháp đưa ra nghiệm chính xác bởi một số hữu hạn các phép tốn số học sơ cấp (khử Gauss). Các phương pháp trực tiếp thực hiện trên máy tính thường khơng đưa đến các nghiệm chính xác do các sai số làm trịn. Lặp: xuất phát với một xấp xỉ ban đầu và dùng một thuật tốn đã được lựa chọn phù hợp, sẽ đưa ra các xấp xỉ liên tiếp ngày càng tốt hơn. Chỉ nhận được nghiệm gần đúng. PHƯƠNG PHÁP SỐ-Bài 3 5
  6. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (5) 5. CÁC PHÉP BIẾN ĐỔI TƯƠNG ĐƯƠNG HỆ PTTT Định lí 4.7 Cho Ax = b là hệ PTTT, các xử lí hệ này bằng các phép tốn sau đây dẫn tới các hệ tương đương: (i) Nhân một phương trình với một hằng số khác 0 (ii) Cộng một phương trình đã nhân với một hằng số vào một phương trình khác (iii) Đổi chỗ hai phương trình Nhận xét: Nếu với k và i xác định mà akk ≠ 0 thì chúng ta cĩ thể khử ẩn xk từ phương trình thứ i bất kì bằng cách cộng phương trình thứ k đã nhân với –(aik /akk) vào phương trình thứ i. Khi đĩ hệ phương trình hệ quả sau là tương đương với hệ ban đầu ~ Ãx b PHƯƠNG. PHÁP SỐ-Bài 3 6
  7. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (6) 6. PHÉP KHỬ GAUSS: đưa hệ Ax = b, về hệ tam giác trên: Gọi W là ma trận cấp n x (n + 1) chứa ma trận vuơng A ở n cột đầu và vectơ b ở cột cuối cùng. Với k = 1, , n–1 lặp các cơng việc: Tìm hàng trụ i ≥ k gần hàng k nhất sao cho wik ≠ 0 Nếu khơng cĩ hàng i như vậy thì dừng (A khơng khả nghịch) Nếu tìm được hàng i, thì đổi chỗ hàng i với hàng k trụ wkk Với i = k+1, , n, lặp các cơng việc sau (khử các pt dưới trụ) mi = wik / wkk là hệ số nhân cho hàng i Với j = k+1, , n+1, lặp các cơng việc sau (biến đổi hàng i) wij = w ij – mi x wkj Nếu wnn = 0 thì dừng (A khơng khả nghịch) w ij i j Nếu trái lại hệ tam giác trên Ux = y với uij  0 i j và yj = wi,n+1 i = 1, , n PHƯƠNG PHÁP SỐ-Bài 3 7
  8. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (7) 6. PHÉP KHỬ GAUSS (tiếp) Ví dụ : Giải hệ PTTT sau: 2x – 3y = 3 4x – 5y + z = 7 2x – y – 3z = 5 Giải: Bước 1: Dùng trụ đầu tiên của hệ để khử những hệ số bên dưới trụ đĩ: m2 = 4/2 = 2, m3 = 2/2 = 1 lấy pt 2 trừ đi 2 lần pt 1, và pt 3 trừ đi 1 lần pt 1 ta được hệ: 2x – 3y = 3 1 y + z = 1 2 y – 3z = 2 Bước 2: Dùng trụ thứ hai để khử hệ số bên dưới nĩ: m3 = 2/1 = 2 lấy pt 3 trừ đi 2 lần pt 2 ta được hệ 2x – 3y = 3 1y + z = 1 –5z = 0 giải ngược từ dưới lên ta được: z = 0, y = 1, x = 3 PHƯƠNG PHÁP SỐ-Bài 3 8
  9. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (8) void giaiGauss(vector > a) { unsigned i, j, k, n = a.size(); double t, m; for (k = 0; k < n–1; k++) { // BUOC KHU XUOI HE TAM GIAC TREN while (a[i] [k] == 0 && i < n) i++; // Tim hang i gan k nhat: a[i] [k]≠0 if (a[i] [k] != 0 && i != k) { for (j = i; j <= n; j++) { // Doi cho hang i va hang k t = a[i] [j]; a[i] [j] = a[k] [j]; a[k] [j] = t; } } else if(a[i] [k] == 0) { // Ma tran A suy bien cout<<"He khong co duy nhat nghiem"<< endl; exit(1); } for (i = k+1; i <= a.size()–1; i++) { // Tao ma tran dang tam giac tren m = a[i] [k] / a[k] [k]; for (j = k+1; j <= a.size(); j++) a[i] [j] = a[i] [j] – m * a[k] [j]; } } // BUOC QUET NGUOC - GIAI HE TAM GIAC TREN } Độ phức tạp thời gian của thuật tốn là O(n3) PHƯƠNG PHÁP SỐ-Bài 3 9
  10. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (9) 7. GIẢI HỆ PTTT DẠNG 3 ĐƯỜNG CHÉO Ma trận A = (aij) cấp n cĩ dạng ba đường chéo nếu aij = 0 bất cứ khi nào |i – j| > 1. Hệ PTTT 3 đường chéo cĩ dạng d1x1 u1x2 b1 l2x1 d2x2 u2x3 b2 l3 x2 d3x3 u3x4 b3     ln 1xn 2 dn 1xn 1 un 1xn bn 1 lnxn 1 dnxn bn 2 1 0 0 1 1 2 1 0 0 Ví dụ: hệ PTTT bên phải là hệ ba đường chéo 0 1 2 1 0 0 0 1 2 0 PHƯƠNG PHÁP SỐ-Bài 3 10
  11. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (10) 7. GIẢI HỆ PTTT DẠNG 3 ĐƯỜNG CHÉO Cho trước các hệ số li, di, ui và vế phải bi của hệ ba đường chéo lixi–1 + dixi + uixi+1 = bi, i = 1 , . . . , n (với l1 = un = 0) Bước 1: Với k = 2, , n lặp các cơng việc sau đây Nếu dk−1 = 0, đĩ là dấu hiệu sai và dừng Trái lại, m = lk / dk–1 và tiếp tục tính dk = dk – m*uk–1 bk = bk – m*bk–1 Nếu dn = 0, đĩ là dấu hiệu sai và dừng Trái lại: xn = bn / dn Bước 2: Với k = n–1, , 1 lặp các cơng việc sau b u *x x k k k 1 k dk Nhận xét: Độ phức tạp thời gian của thuật tốn là O(n) PHƯƠNG PHÁP SỐ-Bài 3 11
  12. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (11) 7. GIẢI HỆ PTTT DẠNG 3 ĐƯỜNG CHÉO (tiếp) for (unsigned k=1; k = 0; i ) x[i] = (b[i] - u[i] * x[i+1]) / d[i]; } for (unsigned i = 0; i < n; i++) cout << "x[" << i + 1 << "]=" << x[i] << endl; PHƯƠNG PHÁP SỐ-Bài 3 12
  13. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (12) 8. CHIẾN LƯỢC XOAY Nhận xét: Thuật tốn khử Gauss đưa ra lời giải sai cho ví dụ sau 0.0003 x1 + 1.566 x2 = 1.569 0.0003 x1 + 1.566 x2 = 1.569 0.3454 x1 – 2.436 x2 = 1.018 – 1805.42 x2 = – 1807.46 Hệ này cĩ nghiệm đúng là x1= 10, x2 = 1, nhưng phép khử Gauss đưa ra x2 =1.001; x1 = 4.780 do trụ a11 = 0.0003 “rất nhỏ” so với a12 Nếu lấy a12 làm trụ, ta cĩ nghiệm đúng! Phân tích Trong bước khử Gauss đưa hệ Ax = b về hệ tam giác trên Ux = y sử dụng ma trận W cấp n × (n + 1) ta đã tìm các hàng trụ ở các bước lặp k = 1, , n–1 theo tiêu chí: hàng trụ i ≥ k gần hàng k nhất sao cho wik ≠ 0 PHƯƠNG PHÁP SỐ-Bài 3 13
  14. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (13) 8. CHIẾN LƯỢC XOAY (tiếp) Xoay theo tỉ lệ riêng: Ở lần lặp thứ k của phép khử Gauss • Tính “cỡ” di của hàng i của ma trận A, với i = k, . . . , n. di max | aij | k j n • Chọn ra phương trình xoay là một trong các phương trình trong n – k ứng cử viên cịn lại mà hệ số của xk cĩ giá trị tuyệt đối lớn nhất so với cỡ của phương trình đĩ Trong thuật tốn Gauss, điều này cĩ nghĩa là số nguyên j được lựa chọn là số nguyên nằm giữa k và n (thường là số nhỏ nhất) sao cho tỉ lệ riêng | w jk | | w | ik với mọi i = k, , n d j di PHƯƠNG PHÁP SỐ-Bài 3 14
  15. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (14) 9. PHÂN TÍCH MA TRẬN THÀNH TÍCH − GIẢI HỆ PTTT Sử dụng phép khử Gauss, khơng đổi chỗ các hàng, giải hệ Ax = b Bước 1: Phân tích A = LU, với U là ma trận tam giác trên, L là ma trận tam giác dưới, chứa các hệ số nhân của lần lặp thứ k w m , i k ik ik (k 1) (k 1) mik aik / akk , lik 1, i k k 1,, n 1; i k 0, i k Bước 2: - Giải hệ tam giác dưới tìm được y Ly b * (b ), i 1,,n pi Vector p = (pi), i=1, ,n dùng để lưu thay đổi thứ tự của các hàng - Giải hệ tam giác trên Ux = y tìm được x PHƯƠNG PHÁP SỐ-Bài 3 15
  16. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (15) 9. PHÂN TÍCH MA TRẬN THÀNH TÍCH − GIẢI HỆ PTTT THUẬT TỐN THÉ XUƠI VÀ KHỬ NGƯỢC PHA 1: Biến đổi ma trận A ban đầu về tích các ma trận tam giác Với k = 1, , n–1 lặp các cơng việc: Tìm hàng trụ i ≥ k gần hàng k nhất sao cho nĩ cĩ tỉ lệ riêng lớn nhất Nếu khơng cĩ hàng i như vậy thì dừng (A khơng khả nghịch) Trái lại, thì đổi chỗ 2 hàng i và k, nếu i ≠ k, lưu thứ tự p mới của b Với i = k+1, , n, lặp các cơng việc (khử các phần tử dưới trụ): m = aik / akk (= Giá trị cần khử ở hàng i / Điểm trụ ở hàng trụ k) Với j = k+1, , n-1, lặp các cơng việc (giữ nguyên vế phải): aij = aij – m x akj Gán aik = m PHƯƠNG PHÁP SỐ-Bài 3 16
  17. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (16) 9. PHÂN TÍCH MA TRẬN THÀNH TÍCH - GiẢI HỆ PTTT THUẬT TỐN THÉ XUƠI VÀ KHỬ NGƯỢC PHA 2: Giải hệ tam giác dưới Ly = b* tìm y (L cĩ đường chéo =1) Với k = 1, , n lặp lại việc gán k 1 y b a y k pk  kj j j 1 Giải hệ tam giác trên Ux = y tìm được x Với k = n, n–1, , 1 lặp lại việc gán n yk akjx j j k 1 xk akk PHƯƠNG PHÁP SỐ-Bài 3 17
  18. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (17) 9. PHÂN TÍCH MA TRẬN THÀNH TÍCH − GIẢI HỆ PTTT Ví dụ 1: Giải hệ u + 2v = 5 1 2 A 4u + 9v = 21 4 9 với hệ số nhân là 4, lưu trong L. Vế phải được dùng để tìm y 1 0 5 5 Giải hệ tam giác dưới Ly = b: y y 4 1 21 1 1 2 5 3 Giải hệ tam giác trên Ux = y: x x 0 1 1 1 PHƯƠNG PHÁP SỐ-Bài 3 18
  19. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (18) 9. PHÂN TÍCH MA TRẬN THÀNH TÍCH − GIẢI HỆ PTTT Ví dụ 2: Giải hệ Ax = b PHA 1 2 3 0 3 3 | ai1 | | a21 | 4 k 1 xét A 4 5 1 7 d 5 max i 2 1 i 3 di d2 5 2 1 3 5 3 4 5 1 7 4 5 1 7 1 l21 2 A* 2 3 0 3 1 - 1 1 3 l 1 2 2 2 31 2 1 3 -7 2 1 3 5 2 2 2 5 4 5 1 7 1 0 0 k 2 hàng tru i 2 l - 3 1 - 1 1 3 L 1 1 0 32 2 2 2 2 1 1 2 - 3 -5 5 2 3 1 PHƯƠNG PHÁP SỐ-Bài 3 19
  20. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (19) 9. PHÂN TÍCH MA TRẬN THÀNH TÍCH − GIẢI HỆ PTTT Ví dụ 2: PHA 2 Giải hệ Ly = b* với 1 0 0 7 7 1 1 L | b * 2 1 0 3 y 2 1 - 3 1 5 0 2 Giải hệ Ux = y với 4 5 1 7 3 U | y 0 - 1 1 1 x 1   2 2 2 0 0 -5 0 0 PHƯƠNG PHÁP SỐ-Bài 3 20
  21. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (20) 9. PHÂN TÍCH MA TRẬN THÀNH TÍCH − GiẢI HỆ PTTT Nhận xét 1. Để nhận được nội dung cuối cùng của ma trận làm việc W, phép khử Gauss cần O(n3) phép cộng và bằng ấy phép nhân / chia. Thuật tốn phân tích A = LU giải hệ PTTT chỉ cần n2 phép nhân / chia và n(n–1) phép cộng 2. Thuật tốn phân tích A = LU giải hệ PTTT cĩ thể dùng giải hệ phương trình Ax = b với các vế phải khác nhau 3. Để đổi chỗ hàng i cho hàng j của ma trận W ta sử dụng phép nhân trái ma trận hốn vị Pij (nhận được từ ma trận đơn vị I cùng cấp, bằng cách đổi chỗ các hàng i và j của I) với W. Ví dụ P 1 0 0 1 1 1 0 0 2 4 1 2 4 1 32 3 5 0 0 1 0 0 3 0 6 5 0 0 1 và 0 1 0 5 3 0 1 0 0 6 5 0 0 3 PHƯƠNG PHÁP SỐ-Bài 3 21
  22. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (21) 10. PHÉP LẶP ĐIỂM CỐ ĐỊNH GIẢI HỆ PTTT Nếu đặt a11 a12 a1n b1 x1 a a a b 2 x2 A 21 22 2n , và b , x   bn xn a m1 a m2 a mn thì hệ Ax = b cịn cĩ thể viết ở dạng phương trình ma trận b – Ax = 0 Giả sử tồn tại ma trận C xấp xỉ ma trận nghich đảo A-1 theo nghĩa ||I – CA|| < 1 Khi đĩ b – Ax = 0  C(b – Ax) = 0  x = x + C(b – Ax) Chọn g(x) = x + C(b – Ax) = (I – CA)x + Cb PHƯƠNG PHÁP SỐ-Bài 3 22
  23. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (22) 10. PHÉP LẶP ĐIỂM CỐ ĐỊNH GIẢI HỆ PTTT Nhận xét: g(x) – g(y) = Cb + (I– CA)x – [Cb + (I – CA)y] = (I – CA)(x – y) Vậy ||g(x) – g(y)|| < ||I – CA|| ||x – y|| Chứng tỏ x-y được g co lại với tỉ lệ co K=||I – CA||< 1, phép lặp x(m +1) = x(m) + C(b – Ax(m)) = Cb + (I – CA)x(m) , m = 0, 1, 2, . khởi đầu từ bất kì x(0) nào, cũng sẽ hội tụ về nghiệm duy nhất ξ của PT Ax = b với sai số ở mỗi bước lặp giảm xuống ít nhất theo tỉ lệ K = || I – CA||. Tìm các ma trận C xấp xỉ ma trận nghịch đảo A-1? PHƯƠNG PHÁP SỐ-Bài 3 23
  24. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (23) 10. PHÉP LẶP ĐIỂM CỐ ĐỊNH GIẢI HỆ PTTT Thuật tốn Phép Lặp Điểm Cố Định đối với các hệ PTTT Cho hệ phương trình tuyến tính Ax = b bậc n. 1. Lấy một ma trận C bậc n thỏa mãn các điều kện (i) Đối với vectơ r cho trước, vectơ Cr tính được “một cách dễ dàng” (ii) Với một chuẩn ma trận nào đĩ, ||I – CA|| < 1 2. Lấy một vectơ n chiều x(0), (chẳng hạn, x(0) = 0) Với m = 0, 1, 2, , cho đến khi thỏa mãn điều kiện lặp, tính x(m+1) = x(m) + C(b – Ax(m)) PHƯƠNG PHÁP SỐ-Bài 3 24
  25. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (23) 10. PHÉP LẶP ĐIỂM CỐ ĐỊNH GIẢI HỆ PTTT Sự hội tụ và sai số: ||x(m+1) - x(m|| = ||I – CA|| ||x(m) - x(m-1) || = K ||x(m) - x(m-1) || = = Kn ||x(1) - x(0) | → 0 khi m → ∞. Bỏ qua sai số làm trịn, dãy x(0), x(1), x(2), hội tụ về nghiệm đúng ξ của hệ PTTT đã cho. Với K = ||I – CA|| < 1 ta suy ra được sai số ||x(m) – ξ|| = ||g(x(m-1))– g(ξ)|| ≤ K ||x(m-1) – ξ|| (*) Do ||x(m-1) – ξ|| ≤ ||x(m) – ξ|| + ||x(m-1) - x(m) || ( ) → ||x(m) – ξ|| ≤ ||x(m) – x(m –1)|| * K / (1-K) PHƯƠNG PHÁP SỐ-Bài 3 25
  26. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (24) 10. PHÉP LẶP ĐIỂM CỐ ĐỊNH GIẢI HỆ PTTT Tiêu chuẩn kết thúc lặp (a): ||x(m) – x(m –1)|| M với M cho trước Vì với K = ||I – CA|| < 1 yếu tố ||x(m) – x(m –1)|| < ε bao hàm điều ||x(m) – ξ|| < ε Tiêu chuẩn cuối luơn cĩ mặt trong bất cứ chương trình lặp nào để tránh việc lặp khơng kết thúc PHƯƠNG PHÁP SỐ-Bài 3 26
  27. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (25) 10. PHÉP LẶP ĐIỂM CỐ ĐỊNH GIẢI HỆ PTTT Phép lặp Jacobi: Giả sử A là ma trận chéo trội thực sự, nghĩa là với i = 1, 2, , n |a ii| |a ij| j i Gọi D = diag(a11, a22, . . . ,ann) là ma trận đường chéo của A. Xét  1 || I D A || max 1 |aij|/|aii| max|aij|/|aii| 1 i j  i j i Chứng tỏ D-1 là một nghịch đảo xấp xỉ của A. Sơ đồ lặp x(m +1) = x(m) + D–l (b – Ax(m)) ↔ x(m +1) = A’x(m) + b’ và cơng thức lặp tính các thành phần của nghiệm, j= 1, , n x(m 1) b a x(m) / a (*) i i  ij j ii j i PHƯƠNG PHÁP SỐ-Bài 3 27
  28. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (26) 10. PHÉP LẶP ĐIỂM CỐ ĐINH GIẢI HỆ PTTT Thuật tốn lặp Jacobi: từ cơng thức lặp (*) ta suy ra thuật tốn sau 1. Cấu tạo ma trận lặp A’ và b’ với a a b 0 12  1n 1 a a a 11 11 11 a a b 21 0  2n 2 A' b' a22 a22 a22      a a b n1 n2  0 n ann ann ann 2. Lặp x(m+1) = A’ x(m) + b’ với nghiệm xấp xỉ ban đầu x(0) tùy ý 3. Điều kiện dừng: khi một trong 3 tiêu chuẩn (a), (b), (c) thỏa mãn . PHƯƠNG PHÁP SỐ-Bài 3 28
  29. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (27) 10. PHÉP LẶP ĐIỂM CỐ ĐỊNH GIẢI HỆ PTTT Ví dụ, giải hệ sau bằng phép lặp Jacobi: (0) (0) (0) 10x1 + x2 + x3 = 12 x1 x2 x3 x1 + 10x2 + x3 = 12 0 0 0 x1 + x2 + 10x3 = 12 1.2 1.2 1.2 (0) với x = 0 và ma trận lặp 0.96 0.96 0.96 0 0.1 0.1 1.2 1.008 1.008 1.008 A' 0.1 0 0.1 và b' 1.2 0.9984 0.9984 0.9984 0.1 0.1 0 1.2 1.00032 1.00032 1.00032 Dãy này hội tụ về nghiệm đúng 0.999936 0.999936 0.999936 [1 1 l]T. Do 1 1 1 (6) || I D A || max  0.2 || x ξ || 0.000096 i 10 10 NHẬN XÉT: Sai số này quá lớn vì thực ra ||ξ – x(6)|| = 0.000064 PHƯƠNG PHÁP SỐ-Bài 3 29
  30. void lapJacobi (double eps, vector > a) { unsigned i, j, k , n = a.size(); vector xt(n,0), xs(n,0); // Khởi tạo vector nghiệm bool t; // Biến logic kiểm tra sai số k = 1; // Biến đếm số lần lặp do { for(i = 0; i <= n –1; i++) { xs[i] = 0; for (j = 0; j <= n – 1; j++) if (j != i) xs[i] = xs[i] – a[i][j] * xt[j]; xs[i] = (xs[i] + a[i][n]) / a[i][i]; } // Kiểm tra đ/k dừng t t = max {fabs(xs[i] – xt[i])} <= eps } while ( ( !t ) && (k <=100) ); }PHƯƠNG PHÁP SỐ-Bài 3 30
  31. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (28) 10. PHÉP LẶP ĐIỂM CỐ ĐỊNH GIẢI HỆ PTTT Phép lặp Gauss–Seidel, hay phương pháp thay thế liên tiếp: Phân tích A = L + D + U với L là ma trận tam giác dưới thực sự, U là ma trận tam giác trên thực sự. Khơng mất TQ, giả sử ma trận đường -1 chéo D khả nghịch, tức là aii≠ 0 với mọi i. Chọn C = (L + D) làm nghịch đảo xấp xỉ của A (nới lỏng đ/k của C), nếu || I (L D) 1A || || I D 1A || 1 Sơ đồ lặp x(m+1) = x(m) + ( L+ D)–1(b – Ax(m)) . Dx(m 1) Lx(m 1) Ux(m) b và cơng thức lặp tính các thành phần của nghiệm (m 1) (m) ai jx j ai jx j b i (m 1) j i j i xi , i 1, 2, , n a ii PHƯƠNG PHÁP SỐ-Bài 3 31
  32. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (29) 10. PHÉP LẶP ĐIỂM CỐ ĐỊNH GIẢI HỆ PTTT Thuật tốn lặp Gauss–Seidel : 1. Cấu tạo ma trận lặp A’ và vectơ b’ như trong phép lặp Jacobi, với giả thiết aii ≠ 0 với mọi i 2. Chọn một nghiệm xấp xỉ ban đầu x(0) 3. Với m = 1, 2, . . . , cho đến khi thỏa mãn, lặp các cơng việc sau: (0) (0) (0) Với i = 1, , n, lặp các cơng việc sau: x1 x2 x3 i 1 n 0 0 0 ( m ) ( m ) ( m 1 ) x i  a ' ij x j  a ' ij x j b ' i 1.2 1.08 0.972 j 1 j i 1 0.9948 1.0033 1.00019 Xét ví dụ giải hệ trên, ta cĩ bảng nghiệm 0.99965 1.000016 1.0080033 NHẬN XÉT: Phép lặp Gauss-Seidel hội tụ nếu ma trận hệ số A là chéo trội / tam giác dưới / xác định dương. Khi đĩ nĩ cĩ tốc độ hội tụ cao hơn phép lặp Jacobi PHƯƠNG PHÁP SỐ-Bài 3 32
  33. void lapSeidel (double eps, vector > a) { unsigned i, j, k, n = a.size(); vector xt(n,0), xs(n,0); // Khởi tạo vector nghiệm bool t; k = 1; // Đếm số lần lặp do { for (i = 0; i I ) xs[i] = xs[i] – a[i][j] * xt[j]; } xs[i] = (xs[i] + a[i][n]) / b[i][i]; } // Kiểm tra đ/k dừng t = max {fabs(xs[i] – xt[i])} <= eps } while ( ( !t ) && (k <=100) ); }PHƯƠNG PHÁP SỐ-Bài 3 33
  34. HỆ PHƯƠNG TRÌNH TUYẾN TÍNH (30) 11. TỔNG KẾT CÁC PHƯƠNG PHÁP SỐ GiẢI HỆ PTTT Trực tiếp (khử Gauss và các cải thiện) : 1. Là các phương pháp giải đúng khơng cĩ sai số thường khơng đưa đến các nghiệm chính xác do mắc phải sai số làm trịn số khi tính tốn. 2. Áp dụng được cho mọi hệ vuơng Lặp: 1. Là phương pháp lặp điểm cố định, xuất phát từ một xấp xỉ ban đầu và dùng một thuật tốn lặp phù hợp, đưa ra các nghiệm xấp xỉ liên tiếp ngày càng tốt hơn. 2. (i) Tốc độ hội tụ phụ thuộc vào đốn nhận ban đầu x(0) (ii) Phải xây dựng ma trận lặp thỏa mãn điều kiện hội tụ Thường áp dụng cho các hệ PTTT cĩ ma trận hệ số thưa hay cĩ các tính chất đặc biệt như chéo trội, xác định dương, PHƯƠNG PHÁP SỐ-Bài 3 34
  35. SỬ DỤNG SOLVER GIẢI HỆ PTTT (1) PHƯƠNG PHÁP SỐ-Bài 3 35
  36. SỬ DỤNG SOLVER GIẢI HỆ PTTT (2) PHƯƠNG PHÁP SỐ-Bài 3 36
  37. SỬ DỤNG SOLVER GIẢI HỆ PTTT (3) 1. Nhập các giá trị đốn nhận ban đầu cho x1, x2, , xn 2. Nhập các vế trái của các phương trình f1, f2, , fn vào các ơ riêng được liên kết với giá trị của các ẩn 3. Chọn cơng cụ Solver (Data Tab Solver) a. Đặt “Target Cell” cho ơ cĩ cơng thức = 0 b. Chọn “Equal to Value of 0” c. Đối với “By Changing Cells,” trỏ tới các ơ chứa giá trị ban đầu của x1, x2, , xn d. Thêm các ràng buộc (các phương trình) e. Nhấn nút “Options” , chọn “Assume linear model” f. Nhấn nút “Solve” PHƯƠNG PHÁP SỐ-Bài 3 37
  38. TÍNH MA TRẬN NGHỊCH ĐẢO (1) Phương pháp Gauss-Jordan: để giải AA-1= I tìm A-1, ta tìm -1 1 từng cột của A : AA Ax1 x2 x3  i1 i2 i3  I Chi phí tính tốn chỉ tăng gấp 3 do tất cả các hệ PTTT đều cĩ chung ma trận vế phải A Ví dụ: 2 1 0 1 0 0 [A i1 i2 i3 ] = 1 2 1 0 1 0 0 1 2 0 0 1 Đem 1/2 hàng 1 + hàng 2 Đem 2/3 hàng 2 + hàng 3 2 1 0 1 0 0 2 1 0 1 0 0 3 1 0 1 1 0 0 2 1 1 1 0 2 2 3 2 0 1 2 0 0 1 0 0 4 1 2 1 3 3 3 PHƯƠNG PHÁP SỐ-Bài 3 38
  39. TÍNH MA TRẬN NGHỊCH ĐẢO (2) Phương pháp Gauss-Jordan: tiếp tục phép khử - các hàng được cộng với các hàng trên chúng, để tạo ra các giá trị 0 trên các giá trị trụ 2 1 0 1 0 0 3 3 3 3 (4/3 x hàng 3 + hàng2) 0 0 2 4 2 4 4 1 2 0 0 1 3 3 3 2/3 x hàng 2 + hàng 1 chia mỗi hàng cho giá trị tại điểm trụ của nĩ 3 1 3 1 1 2 0 0 1 1 0 0 2 2 4 2 4 3 3 3 3 1 1 0 0 0 1 0 1  I x x x  2 2 1 2 3 2 4 2 4 4 1 2 1 1 3 0 0 1 0 0 1 3 3 3 4 2 4 PHƯƠNG PHÁP SỐ-Bài 3 39
  40. SỬ DỤNG HÀM MINVERSE TÌM A-1 CHÚ Ý: Cơng thức phải được nhập như cơng thức mảng: 1. Chọn mảng các ơ chứa A-1. 2. Nhập hàm MINVERSE tính A-1 với các tham số là các ơ chứa ma trận A 3. Nhấn phím F2, rồi nhấn tiếp CTRL+SHIFT+ENTER PHƯƠNG PHÁP SỐ-Bài 3 40
  41. SỬ DỤNG HÀM MINVERSE TÌM A-1 PHƯƠNG PHÁP SỐ-Bài 3 41
  42. TÍNH ĐỊNH THỨC Phương pháp khử Gauss: tiến hành phép khử xuơi ma trận A để n đưa nĩ về dạng tam giác trên U Khi đĩ det(A) = (−1) u11u22 unn Ví dụ 2 3 0 n là số lần A 4 -5 1 § ỉi_hµng 1 đổi hàng 2 -1 -3 4 5 1 4 5 1 2 -3 0 0 - 1 1 2 2 2 -1 -3 0 3 -7 2 2 4 5 1 § ỉi_hµng 1 0 - 1 1 det(A) ( 1) * 4 * * ( 5) 10 2 2 2 0 0 - 5 PHƯƠNG PHÁP SỐ-Bài 3 42
  43. SỬ DỤNG HÀM MDETERM TÌM det(A) PHƯƠNG PHÁP SỐ-Bài 3 43
  44. BÀI TẬP 1. Bài tập tính tốn: 4.2-5, 4.2-6, 4.2-8, 4.3-4, 4.4-2, 4.4-4, 5.3-7 2. Bài tập lập trình: Viết các hàm trong C++ giải hệ Ax = b bằng phép khử Gauss cĩ sử dụng chiến lược xoay theo tỉ lệ riêng và phân tích ma trận vuơng thành tích các ma trận tam giác; giải hệ cĩ ma trận chéo trội bằng các phép lặp Jacobi và Gauss-Seidel. PHƯƠNG PHÁP SỐ-Bài 3 44