›
›

Bài giảng Phương pháp số trong công nghệ hóa học - TS. Nguyễn Đặng Bình Thành

pdf 92 trang phuongnguyen 2140
Download
Bạn đang xem 20 trang mẫu của tài liệu "Bài giảng Phương pháp số trong công nghệ hóa học - TS. Nguyễn Đặng Bình Thà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_trong_cong_nghe_hoa_hoc_ts_nguyen_d.pdf

Nội dung text: Bài giảng Phương pháp số trong công nghệ hóa học - TS. Nguyễn Đặng Bình Thành

  1. MỞ ĐẦU PHƯƠNG PHÁP SỐ Phương pháp số được dùng để phân tích và TRONG CÔNG NGHỆ HÓA HỌC giải gần đúng các bài toán với sai số nằm trong giới hạn cho phép. Mã học phần: CH3454 bởi vì hầu hết các bài toán khoa học kỹ thuật đều không có các lời giải chính xác. TS. Nguyễn Đặng Bình Thành BM:Máy & TBCN Hóa chất Phương pháp số thường được bắt đầu từ việc xây dựng mô hình, lựa chọn thuật toán, và đưa ra các đáp số gần đúng. Numerical Methods in Chemical Engineering 1 2 MỞ ĐẦU MỞ ĐẦU Phương pháp số có vai trò quan trọng Phương pháp số trong Kỹ thuật hóa học: trong nhiều lĩnh vực như: Thiên văn Mô tả bằng toán học các quá trình và và thiết học, nông nghiệp, kiến trúc, bị trong công nghệ hóa học. Tính toán thiết kế các quá trình và thiết bị hoạt động trong lĩnh vực kỹ thuật hóa học. Và tất nhiên rất quan trọng trong kỹ Tính toán tối ưu hóa các điều kiện làm việc thuật. và kết cấu các thiết bị hóa chất. Xác định các hằng số thực nghiệm bằng phương pháp hồi quy. 3 4 1
  2. NỘI DUNG NỘI DUNG Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình. trình và hệ phương trình Chương 2. Phương pháp tính tích phân 1.1 Phương pháp giải hệ phương trình tuyến tính Chương 3. Phương trình và hệ phương trình và ứng dụng vi phân Chương 4. Tối ưu hóa 1.1.1 Giải hệ phương trình tuyến tính bằng phương pháp Gauss và phương pháp nghịch đảo ma trận 1.1.2 Ứng dụng để tính toán cân bằng vật chất của hệ thống CNHH 5 6 NỘI DUNG NỘI DUNG Chương 1. Các phương pháp giải phương Chương 2. Phương pháp tính tích phân trình và hệ phương trình 2.1 Tính tích phân xác định bằng phương pháp 1.2 Phương pháp giải phương trình và hệ phương hình thang trình phi tuyến 2.2 Tính tích phân xác định bằng phương pháp Simpson 1.2.1 Giải phương trình phi tuyến bằng phương pháp lặp đơn giản và phương pháp Newton-Raphson 2.3 Ứng dụng 1.2.2 Giải hệ phương trình phi tuyến bằng phương pháp 2.3.1 Tính toán tháp chưng luyện lặp đơn giản và phương pháp Newton-Raphson 2.3.2 Tính toán tháp hấp thụ 1.3 Ứng dụng 7 8 2
  3. NỘI DUNG NỘI DUNG Chương 3. Phương trình và hệ phương Chương 4. Tối ưu hóa trình vi phân 4.1 Tìm cực trị hàm một biến: phương pháp điểm 3.1 Giải phương trình vi phân bằng phương pháp vàng, phương pháp gradien Euler 4.2 Tìm cực trị hàm nhiều biến: phương pháp 3.2 Giải phương trình vi phân bằng phương pháp gradien, phương pháp đơn hình Runge-Kutta 3.3 Giải hệ phương trình vi phân bằng phương 4.3 Cực trị có ràng buộc: phương pháp hàm phạt pháp Euler 3.4 Giải hệ phương trình phi phân bằng phương pháp Runge-Kutta 3.5 Ứng dụng tính toán hệ phản ứng hóa học 9 10 TÀI LIỆU THAM KHẢO Ngôn ngữ lập trình [1] Sổ tay quá trình và công nghệ hóa chất T1, 2, NXB Có nhiều ngôn ngữ lập trình có thể KHKT, 2004. ứng dụng để tính toán các quá [2] Nguyễn Bin. Các quá trình và thiết bị công nghệ hóa trình công nghệ hóa học: chất T1, 2, 3, NXB KHKT, 2001. Matlab; C; C++; Visual Basic; Delphi; [3] R. Perry. Chemical Engineers’ Handbook, 7th Ed., Pascal; Mc. Graw Hill, 2007. Các phầm mềm ứng dụng khác trong [4] K. Johnson. Numerical Methods in Chemistry, Mc. công nghệ hóa học: Aspen Plus; Graw Hill, 1978. gProms; [5] Nguyễn Minh Tuyển, Phạm Văn Thiêm. Kỹ thuật hệ thống trong CN Hóa học, T2, NXB KHKT, 2001. 11 12 3
  4. Ngôn ngữ lập trình Nhắc lại các kiến thức về lập trình Pascal Tìm hiểu bản chất của quá trình và Cấu trúc chương trình các ứng dụng các thuật toán Do đó: Pascal được sử dụng chính trong môn học này! 13 14 Nhắc lại các kiến thức về lập trình Pascal Nhắc lại kiến thức về lập trình Pascal Tên chương trình Sử dụng các thư viện: CRT, GRAPH, Khai báo theo kiểu Khai báo nhãn (khi dùng lệnh goto) Khai báo các hằng số (một giá trị cụ thể) Số thực: “real” Khai báo các kiểu dữ liệu đặc biệt như ma trận, Var Khai báo các biến số cùng với các kiểu tương ứng a,b,c: real; Số nguyên: “integer” Chương trình chính (Begin End.) Var i,j,k,n: integer; Chú ý: Cần phải tuân thủ nghiêm ngặt trình15 tự! 16 4
  5. Nhắc lại kiến thức về lập trình Pascal Nhắc lại kiến thức về lập trình Pascal Khai báo theo kiểu Khai báo theo kiểu Mảng hay ma trận: “array” Mảng hay ma trận: “array” Đối với kiểu này, trước hết phải khai báo kiểu trước! Tuy nhiên cũng có thể khai báo trực tiếp Type Var mx = array [1 50] of real; x: array [1 50] of real; ma = array [1 50, 1 100] of real; a: array [1 50, 1 100] of real; Var x: mx; a: ma; 17 18 Nhắc lại kiến thức về lập trình Pascal Nhắc lại kiến thức về lập trình Pascal Các loại chương trình con Các loại chương trình con Dùng chương trình con khi cần thực hiện một đoạn Trong Pascal có hai loại chương trình con: chương trình lặp đi lặp lại nhiều lần. -Hàm (function) Do đó: khi cần đến những đoạn chương trình như Hàm chỉ trả lại một kiểu dữ liệu và một giá trị duy vậy thì chỉ cần gọi tên chương trình con đó. nhất Thuận lợi: -Thủ tục (procedure) -Chương trình chính đơn giản Thủ tục có thể trả lại nhiều kiểu dữ liệu khác nhau và -Mức độ khái quát hóa chương trình cao có thể trả lại nhiều giá trị -Dễ kiểm tra lỗi cho toàn bộ chương trình -Thuận lợi cho người sử dụng 19 20 5
  6. Nhắc lại kiến thức về lập trình Pascal Nhắc lại kiến thức về lập trình Pascal Các loại chương trình con Các loại chương trình con -Hàm (function) -Thủ tục (procedure) 2 2 Ví dụ: xác định giá Nhượctrị của hàm điểm: số y = 2x + 3x – 5 Ví dụ: xác định giá trị của 2 hàm số y1 = 2x + 3x – 5; Khi tính toán cho nhiều biểu thức, ví dụ: y = 5x3 – 2x – 6. Function F(x: real): real; 2 y = 2x2 + 3x – 5 Khi gọi thủ tục: 1 Procedure HAM; Begin y = 5x3 – 2x – 6 2 Chú ý: {Tính giá trị tại x:=x1} Cần phải dùng: Begin F:=2*x*xĐằng + sau 3*x “End” – 5; của chương trình con dùng dấu 2 chương trình con kiểu hàm x:=x1; “;” thay vì dấu “.” trong chương trình chính. y1:=2*sqr(x) + 3*x – 5; End; HAM; y1:=5*x*sqr(x) – 2*x – 6; Khi gọi hàm trong chương trình chính: {Có thể gán giá trị cho y , y } End; 3 4 y:= F(x1); 21 y3:=y1; y4:=y2; 22 Có cần chương trình con??? Nhắc lại kiến thức về lập trình Pascal Ví dụ áp dụng Các loại chương trình con Thuật toán? Ví dụ 1 a) Tính tổng -Thủ tục (procedure) có tham trị và tham biến hình thức Cho dãy số thực i:=0;S:=0;x1, x2, ,xn Ví dụ: a) Hãy tính tổngi:=i S +của 1; dãy S:=S số + trên x[i]; Procedure HAM(x: real; Var y1,y2:real) b) Tìm giá trị nhỏb) Tìm nhất giá của trị dãy nhỏ nhất i:=1; Begin Khi gọi thủ tục: xmin:= một số rất lớn??? y1:=2*sqr(x) + 3*x – 5; {Tính giá trị tại x:=x1} Nếu x[i] <= xmin thì xmin:=x[i]; y1:=5*x*sqr(x) –HAM(x1,y1,y2);2*x – 6; i:=i+1; End; Hoặc cách khác??? {Có thể gán giá trị cho y3, y4} HAM(x1,y3,y4) 23 24 6
  7. 25 26 27 28 7
  8. Có cần chương trình con??? Ví dụ áp dụng Thuật toán? Ví dụ 2 a) Tính tổng Cho dãy số thực x , x , ,x và hàm số 1 i:=0;S:=0;2 n y = 2x 2 + 3x – 7 i i i i:=i + 1; y[i]:= ; a) Hãy tính tổng S của các giá trị yi S:=S + y[i]; b) Tìm giá trị nhỏ nhất của dãy yi b) Tìm giá trị nhỏ nhất i:=1; ymin:= một số rất lớn??? Nếu y[i] <= ymin thì ymin:=x[i]; i:=i+1; Hoặc cách khác??? 29 30 31 32 8
  9. Bài tập!!! 1) Thực hiện lại ví dụ 1 và ví dụ 2 bằng chương trình Pascal 2) Cho dãy số thực tăng dần x1, x2, ,xn và một số thực xs. Hãy xác định vị trí của xs trong dãy số trên, biết x1 < xs <xn 33 34 Ứng dụng đơn giản: Vấn đề nội suy trong Ứng dụng đơn giản: Vấn đề nội suy trong kỹ thuật hóa học kỹ thuật hóa học Tại sao phải sử Tính toán? dụng nội suy trong tính toán các quá trình CN Hóa học??? Các đường cong này được xây dựng từ??? 35 36 9
  10. Ứng dụng đơn giản: Vấn đề nội suy trong Ứng dụng đơn giản: Vấn đề nội suy trong kỹ thuật hóa học kỹ thuật hóa học Các thuật toán nội suy: Tuyến tính, Lagrance, Nội suy tuyến tính Newton, y y Phương trình đường y Nhưng thẳng đi qua hai n yn điểm (x ,y ) và Không có số liệu k-1 k-1 yk y thực nghiệm!!! (xk,yk): k ys y s yk-1 Đó là??? yk-1 Giả thiết đường cong y Nội suy nối giữa hai điểm là 1 tuyến tính!!! y1 đường thẳng 0 x1 xk-1 xs xk xn x 0 x1 xk-1 xs xk xn x 37 38 Ứng dụng đơn giản: Vấn đề nội suy trong Ứng dụng đơn giản: Vấn đề nội suy trong kỹ thuật hóa học kỹ thuật hóa học Nội suy tuyến tính Nội suy tuyến tính Procedure NOSUY(xs:real;VAR ys:real;Y,X:mX); Thuật toán: Begin 1. Chỉ ra khoảng (xk-1,xk) chứa giá trị xs giá trị của k k:=0; {so sánh xs với các giá trị x1, , xn} 2. Đưa giá trị của k tìm được vào biểu thức nội suy tuyến tính {Số vòng lặp sẽ là không xác định!!!?} {Sử dụng cấu trúc:} {Repeat Until hoặc While End} Repeat k:=k+1; 39 Until xs < x[k]; 40 10
  11. Ứng dụng đơn giản: Vấn đề nội suy trong Ứng dụng đơn giản: Vấn đề nội suy trong kỹ thuật hóa học kỹ thuật hóa học Nội suy tuyến tính Ví dụ 1: Procedure NOSUY(xs:real;VAR ys:real;Y,X:mX); Cho hỗn hợp lỏng Bezne – Toluen, biết hàm lượng Begin Benzen trong pha lỏng x = 0,4 (phần mol). Hãy xác định k:=0; hàm lượng Benze trong pha hơi ở trạng thái cân bằng. Repeat k:=k+1; yCB = ? Until xs < x[k]; {Ra khỏi vòng lặp trên đã tìm được giá trị k} ys:=y[k-1]+(y[k]-y[k-1])*(xs-x[k-1]) x = 0,4 /(x[k]-x[k-1]); End; 41 42 Ứng dụng đơn giản: Vấn đề nội suy trong Ứng dụng đơn giản: Vấn đề nội suy trong kỹ thuật hóa học kỹ thuật hóa học Ví dụ 1: Ví dụ 1: Dữ liệu thực nghiệm về cân bằng pha: Program CB1; x y T x y T uses crt; 0 0 110,6 50 71,2 92,1 type 5 11,8 108,3 60 79 89,4 mX=array [1 50] of real; 10 21,4 106,1 70 85,4 86,8 var 20 38 102,2 80 91 84,4 X,Y:mX; 30 51,1 98,6 90 95,9 82,3 xs,ys:real; 40 61,9 95,2 100 100 80,2 n,i,j,k:integer; Procedure NOSUY(xs:real;VAR ys:real;Y,X:mX); 43 Begin 44 11
  12. Ứng dụng đơn giản: Vấn đề nội suy trong Ứng dụng đơn giản: Vấn đề nội suy trong kỹ thuật hóa học kỹ thuật hóa học Ví dụ 1: Ví dụ 1: Procedure NOSUY(xs:real;VAR ys:real;Y,X:mX); {Chương trình chính} Begin BEGIN k:=0; clrscr; Repeat writeln (‘Nhập số điểm thực nghiệm n = ’); k:=k+1; readln (n); Until xs<x[k]; {Nhập các số liệu của pha lỏng x[i]} ys:=y[k-1]+ (y[k]-y[k-1])*(xs-x[k-1]) For i:=1 to n do /(x[k]-x[k-1]); Begin End; writeln (‘x[‘,i,’] =‘);readln (x[i]); {Chương trình chính} 45 End; 46 Ứng dụng đơn giản: Vấn đề nội suy trong Ứng dụng đơn giản: Vấn đề nội suy trong kỹ thuật hóa học kỹ thuật hóa học Ví dụ 1: Ví dụ 1: {Chương trình chính} {Chương trình chính} BEGIN BEGIN {Nhập các số liệu của pha hơi ở TTCB y[i]} {Tìm hàm lượng pha hơi cân bằng với xs=0.4} For i:=1 to n do writeln (‘Nhập giá trị xs =’);readln(xs); Begin NOISUY (xs,ys,Y,X); writeln (‘y[‘,i,’] =‘);readln (y[i]); {Hiển thị kết quả} End; writeln (‘Ham lượng y cân bằng, yCB =’,ys); {Tìm hàm lượng pha hơi cân bằng với xs=0.4} readln; 47 END. 48 12
  13. Ứng dụng đơn giản: Vấn đề nội suy trong Ứng dụng đơn giản: Vấn đề nội suy trong kỹ thuật hóa học kỹ thuật hóa học Ví dụ 2: Ví dụ 2: Cho hỗn hợp lỏng Bezne – Toluen, biết hàm lượng Dữ liệu thực nghiệm về cân bằng pha: Benzen trong pha khí (hơi) y = 0,6 (phần mol). Hãy xác x y T x y T định hàm lượng Benze trong pha lỏng ở trạng thái cân 0 0 110,6 50 71,2 92,1 bằng. 5 11,8 108,3 60 79 89,4 10 21,4 106,1 70 85,4 86,8 y = 0,6 20 38 102,2 80 91 84,4 30 51,1 98,6 90 95,9 82,3 40 61,9 95,2 100 100 80,2 xCB = ? 49 50 Ứng dụng đơn giản: Vấn đề nội suy trong Ứng dụng đơn giản: Vấn đề nội suy trong kỹ thuật hóa học kỹ thuật hóa học Ví dụ 2: Ví dụ 2: Program CB2; Procedure NOSUY(xs:real;VAR ys:real;Y,X:mX); uses crt; Begin type k:=0; mX=array [1 50] of real; Repeat var k:=k+1; X,Y:mX; Until xs<x[k]; xs,ys:real; ys:=y[k-1]+ (y[k]-y[k-1])*(xs-x[k-1]) n,i,j,k:integer; /(x[k]-x[k-1]); Procedure NOSUY(xs:real;VAR ys:real;Y,X:mX); End; Begin 51 {Chương trình chính} 52 13
  14. Ứng dụng đơn giản: Vấn đề nội suy trong Ứng dụng đơn giản: Vấn đề nội suy trong kỹ thuật hóa học kỹ thuật hóa học Ví dụ 2: Ví dụ 2: {Chương trình chính} {Chương trình chính} BEGIN BEGIN clrscr; writeln (‘Nhập số điểm thực nghiệm n = ’); {Nhập các số liệu của pha hơi ở TTCB y[i]} readln (n); For i:=1 to n do {Nhập các số liệu của pha lỏng x[i]} Begin For i:=1 to n do writeln (‘y[‘,i,’] =‘);readln (y[i]); Begin End; writeln (‘x[‘,i,’] =‘);readln (x[i]); {Tìm hàm lượng pha hơi cân bằng với xs=0.4} End; 53 54 Ứng dụng đơn giản: Vấn đề nội suy trong Ứng dụng đơn giản: Vấn đề nội suy trong kỹ thuật hóa học kỹ thuật hóa học Ví dụ 2: Ví dụ 3: {Chương trình chính} Cho hệ hơi nước bão hòa, biết nhiệt độ của hệ là BEGIN o T = 119,6 C. Hãy xác định áp suất hơi bão hòa pS của hệ? {Tìm hàm lượng pha lỏng cân bằng với ys=0.6} writeln (‘Nhập giá trị ys =’);readln(ys); NOISUY (ys,xs,X,Y); pS = ? {Hiển thị kết quả} T = 119,6oC writeln (‘Ham lượng x cân bằng, yCB =’,xs); readln; END. 55 56 14
  15. Ứng dụng đơn giản: Vấn đề nội suy trong Ứng dụng đơn giản: Vấn đề nội suy trong kỹ thuật hóa học kỹ thuật hóa học Ví dụ 3: Ví dụ 3: {Chương trình chính} BEGIN {Tìm áp suất hơi bão hòa tại ts=119.6} writeln (‘Nhập giá trị ts =’);readln(ts); NOISUY (ts,ps,P,T); {Hiển thị kết quả} writeln (‘Áp suất hơi bão hòa, pS =’,ps); readln; 57 END. 58 Ứng dụng đơn giản: Vấn đề nội suy trong Ứng dụng đơn giản: Vấn đề nội suy trong kỹ thuật hóa học kỹ thuật hóa học Ví dụ 4: Ví dụ 4: Cho hệ hơi nước bão hòa, biết áp suất hơi bão hòa của {Chương trình chính} hệ là pS = 2 atm. Hãy xác định nhiệt độ của hệ? BEGIN {Tìm nhiệt độ của hệ tại ps=2.0} writeln (‘Nhập giá trị ps =’);readln(ps); pS = 2 atm NOISUY (ps,ts,T,P); T = ? {Hiển thị kết quả} writeln (‘Nhiệt độ của hệ, tS =’,ts); readln; 59 END. 60 15
  16. Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Sơ đồ chưng đơn giản 61 Sơ đồ chưng luyện liên tục 62 Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Lượng lỏng chảy xuống Lượng lỏng đi xuống nhận từ đĩa trên trao đổi nhiệt nhiệt từ dòng hơi và bay và chất với lượng hơi từ hơi một phần cấu tử dễ đĩa dưới đi lên bay hơi Do đó: Lượng lỏng đi xuống và Trong pha hơi ở đĩa trên lượng hơi đi lên từ một có hàm lượng cấu tử dễ đĩa ở trạng thái cân bằng bay hơi cao hơn đĩa dưới Truyền nhiệt và chuyển khối 63 Truyền nhiệt và chuyển khối trên từng đĩa 64 16
  17. Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Đi từ đáy tháp đến đỉnh Dòng hơi mất nhiệt ngưng tháp làm lượng cấu tử dễ tụ một phần cấu tử khó bay hơi tăng dần và bay hơi ngưng tụ cho sản phẩm đỉnh Do đó: Đi từ đỉnh tháp xuống đáy Trong pha lỏng ở đĩa dưới tháp hàm lượng cấu tử khó có hàm lượng cấu tử khó bay hơi tăng dần và cho bay hơi cao hơn đĩa trên sản phẩm đáy Truyền nhiệt và chuyển khối trên từng đĩa 65 Chưng luyện liên tục 66 Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Xác định số đĩa lý thuyết bằng phương pháp đồ thị Một đĩa lý (Phương pháp McCabe-Thiele) thuyết Dữ liệu cân bằng pha lấy từ các thực nghiệm Đường cân bằng được xác định từ cân bằng vật chất 67 68 17
  18. Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Các giả thiết: Đường làm 1- Ẩn nhiệt hóa hơi của hai cấu tử bằng nhau và là hằng việc đoạn số. luyện 2- Nhiệt lượng làm thay đổi nhiệt độ các cấu tử trong tháp là rất nhỏ so với ẩn nhiệt hóa hơi Với các giả thiết đã nêu thì 3- Hỗn hợp lỏng hai cấu tử là hỗn hợp lý tưởng đường làm việc của đoạn Đường làm 4- Không có sự trao đổi nhiệt với môi trường xung quanh chưng và đoạn luyện là đường thẳng việc đoạn 5- Áp suất làm việc là hằng số tại mọi điểm trong tháp chưng 6- Hỗn hợp dòng vào ở trạng thái bão hòa 69 70 Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Phương trình đường làm việc Chỉ số hồi lưu: Đoạn luyện: Hồi lưu hoàn toàn (Rmax) Đoạn chưng: F, P,W: mol/h hoặc kmol/h R: Chỉ số hồi lưu (-) 71 72 18
  19. Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Chỉ số hồi lưu: Xác định số đĩa lý thuyết bằng phương pháp số Hồi lưu nhỏ nhất (Rmin) Một đĩa lý thuyết 73 74 Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Xác định số đĩa lý thuyết bằng phương pháp số Xác định số đĩa lý thuyết bằng phương pháp số Thuật toán? Một đĩa lý thuyết Thuật toán? Một đĩa lý thuyết ys:=xp Tìm yf qua đường làm việc Nội suy tìm xs ys:=yf Tìm ys qua phương Nội suy tìm xs trình đường làm việc Tìm ys qua đường làm việc Lặp lại cho tới khi xF xs xf Lặp lại cho tới khi xs xw Chú ý: Đối với đoạn luyện tính từ xp 75 Chú ý: Đối với đoạn luyện tính từ xf 76 19
  20. Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Chương trình (Đoạn luyện) Chương trình (Đoạn chưng) ys:=xp; yf:=DLVC(xf); NLTL:=0; ys:=yf; Repeat NLTC:=0; NOISUY(ys,xs,X,Y); Repeat ys:=DLVL(xs); NOISUY(ys,xs,X,Y); NLTL:=NLTL+1; ys:=DLVC(xs); Until xs <= xf; NLTC:=NLTC+1; Until xs <= xw; 77 78 Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Chương trình Chương trình Program chungluyen_1; Program chungluyen_1; uses crt; type {Các chương trình con} Procedure NOSUY(xs:real;VAR ys:real;Y,X:mX); mX = array [1 50] of real; Begin var k:=0; Repeat X,Y:mX; k:=k+1; F,P,W,xF,xP,xW,R:real; Until xs<x[k]; ys:=y[k-1]+ (y[k]-y[k-1])*(xs-x[k-1]) NLT,NLTC,NLTL,n,i,k:integer; /(x[k]-x[k-1]); {Các chương trình con} End; 79 80 20
  21. Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Chương trình Chương trình Program chungluyen_1; Program chungluyen_1; {Các chương trình con} Function DLVL(xs:real):real; {Chương trình chính} Begin BEGIN DLVL:=R/(R+1)*xs+xP/(R+1); clrscr; End; {số liệu đầu} Function DLVC(xs:real):real; write (‘Nhập số điểm thực nghiệm n = ’); Begin readln(n); DLVC:=(R+F/P)*xs/(R+1)+(1-F/P)/(R+1)*xW; {Nhập các giá trị thực nghiệm} End; 81 82 Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Chương trình Chương trình Program chungluyen_1; Program chungluyen_1; BEGIN BEGIN {Nhập các giá trị thực nghiệm} {Nhập các giá trị thực nghiệm} For i:=1 to n do For i:=1 to n do Begin Begin write (‘X[’,i,‘] = ’);readln(x[i]); write (‘Y[’,i,‘] = ’);readln(Y[i]); End; End; 83 84 21
  22. Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Cần có cân bằng chất cho toàn Chương trình Chương trình tháp: Program chungluyen_1; Program chungluyen_1; FP:=F*(xP = P + W -xF)/(xP-xW); F.xF = P.xP + W.xW BEGIN W:=F-P; BEGIN {Tính toán số đĩa lý thuyết đoạn luyện} {Nhập các số liệu yêu cầu} ys:=xP; write (‘Nhập F = ’);readln(F); NLTL:=0; Repeat write (‘Nhập xF = ’);readln(xF); NOISUY(ys,xs,X,Y); write (‘Nhập xP = ’);readln(xP); ys:=DLVL(xs); write (‘Nhập xW = ’);readln(xF); NLTL:=NLTL+1; Until xs <= xF; write (‘Nhập R = ’);readln(R); 85 86 Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Chương trình Chương trình Program chungluyen_1; Program chungluyen_1; BEGIN BEGIN {Tính toán số đĩa lý thuyết đoạn chưng} {Hiển thị kết quả} yf:=DLVC(xf); writeln (‘NLTL = ’,NLTL); ys:=yf; writeln (‘NLTC = ’,NLTC); NLTC:=0; {Xác định số đĩa lý thuyết cho toàn tháp} Repeat NLT = NLTL + NLTC; NOISUY(ys,xs,X,Y); writeln (‘NLT = ’,NLT); ys:=DLVC(xs); readln; NLTC:=NLTC+1; END. Until xs <= xw; 87 88 22
  23. Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Chương trình Chương trình Procedure DiaLT(R:real; var NLT:integer); Trong một số chương trình việc xác định số đĩa lý thuyết Begin cần phải lặp đi lặp lại nhiều lần. {Xác định số đĩa lý thuyết đoạn luyện} Do đó: {Xác định số đĩa lý thuyết đoạn chưng} Nên xây dựng chương trình con xác định số đĩa lý thuyết. NLT:= NLTL+NLTC; End; 89 90 Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Phương pháp xác định chỉ số hồi lưu thích hợp Phương pháp xác định chỉ số hồi lưu thích hợp Phương trình đường làm việc đoạn luyện Và khi nào thì R min? Điều gì xảy ra khi R ? 91 92 23
  24. Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Phương pháp xác định chỉ số hồi lưu thích hợp Phương pháp xác định chỉ số hồi lưu thích hợp Source: Richardson & Coulson. Chemical Engineering, vol.2 93 Source: Richardson & Coulson. Chemical Engineering, vol.2 94 Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Phương pháp xác định chỉ số hồi lưu thích hợp Phương pháp xác định chỉ số hồi lưu thích hợp Thuật toán? Vậy: Xác định R Xác định Rmin: NOISUY(xF,yFCB,Y,X) thích hơp thế nào? Cho R tăng dần: R = Rmin + R Mối quan hệ thể tích tháp, Xác định NLT với R đã biết DiaLT(R,NLT) Lặp NLt, và R: lại Xác định V: nhiều Tuy nhiên có thể lấy: Tìm giá trị nhỏ nhất của V lần NLT và Reff 95 96 Source: Richardson & Coulson. Chemical Engineering, vol.2 Trong khoảng (Rmin – 5Rmin) 24
  25. Ứng dụng nội suy trong tính toán số đĩa Ứng dụng nội suy trong tính toán số đĩa lý thuyết và chiều cao tháp chưng luyện lý thuyết và chiều cao tháp chưng luyện Phương pháp xác định chỉ số hồi lưu thích hợp Phương pháp xác định chỉ số hồi lưu thích hợp Chương trình Chương trình Repeat NOISUY(xF,yFCB,Y,X); R:=R+deltaR; Rmin:=(xP-yFCB)/(yFCB-xF); DiaLT(R,NLT); deltaR:=0.05; V:=NLT*(R+1); R:=Rmin; if V =5*Rmin; 98 Ứng dụng nội suy trong tính toán số đĩa Chương 1. Các phương pháp giải phương lý thuyết và chiều cao tháp chưng luyện trình và hệ phương trình Bài tập 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng 1. Cho hỗn hợp Benzen – Toluen biết Hệ phương trình tuyến tính F = 300 kmol/h; xF = 0.79 Xác định số đĩa lý thuyết của tháp để có được: xP = 0.99; xW = 0.01 tại các giá trị của R: 2. Xác định chỉ số hồi lưu thích hợp Reff và số đĩa lý thuyết tương ứng cho quá trình chưng luyện trên. 99 100 25
  26. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Dạng rút gọn (dạng vector): 101 102 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss Phương pháp khử Gauss Là phương pháp khử dần các ẩn để đưa hệ phương trình đã cho về dạng tam giác trên rồi giải hệ này từ dưới lên không phải tính định thức 103 104 26
  27. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss Phương pháp khử Gauss Các bước thực hiện: 1. Quá trình xuôi 1. Quá trình xuôi Bước 0: Dùng pt đầu tiên để khử x1 trong n-1 pt còn lại. 2. Quá trình ngược Để khử x1 ở hàng thứ k (k = 2,3, ,n) tính lại các hệ số ak,j ở hàng thứ k (j = 1,2, ,n): ak,j = ak,j – a1,j*ak,1/a1,1 và tính lại hệ số bk ở hàng thứ k: bk = bk – b1*ak,1/a1,1 105 106 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss Phương pháp khử Gauss 1. Quá trình xuôi 1. Quá trình xuôi Bước 1: Dùng pt thứ 2 để khử x2 trong n-2 pt còn lại phía sau. Bước i: Dùng pt thứ i để khử xi trong (n-i) pt còn lại phía sau. Để khử x2 ở hàng thứ k (k = 3,4, ,n) Để khử xi ở hàng thứ k (k = i+1,i+2, ,n) tính lại các hệ số ak,j ở hàng thứ k (j = 2,3, ,n): tính lại các hệ số ak,j ở hàng thứ k (j = i,i+1, ,n): ak,j = ak,j – a2,j*ak,2/a2,2 ak,j = ak,j – ai,j*ak,i/ai,i và tính lại hệ số bk ở hàng thứ k: và tính lại hệ số bk ở hàng thứ k: bk = bk – b2*ak,2/a2,2 bk = bk – bi*ak,i/ai,i 107 108 27
  28. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss Phương pháp khử Gauss 1. Quá trình xuôi 1. Quá trình xuôi: Sau khi khử hệ phương trình có dạng Bước n-1: Dùng pt thứ i để khử xn-1 trong pt thứ n. Dạng 1: Nếu tại các bước (bước i) không chia cho hệ số ai,i trước khi thực hiện quá trình khử. Để khử xn-1 ở hàng thứ n tính lại các hệ số an,j ở hàng thứ n (j = n-1,n): an,j = an,j – an-1,j*an,n-1/an-1,n-1 và tính lại hệ số bn ở hàng thứ n: bn = bn – bn-1*an,n-1/an-1,n-1 109 110 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss Phương pháp khử Gauss 1. Quá trình xuôi: Sau khi khử hệ phương trình có dạng 2. Quá trình ngược Dạng 2: Nếu tại các bước (bước i) chia cho hệ số ai,i trước khi thực Xuất phát từ pt thứ n ở các hệ pt dạng 1 hoặc dạng 2 lần lượt xác hiện quá trình khử. định được các giá trị xi thông qua các biểu thức: 111 112 28
  29. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss Phương pháp khử Gauss Ví dụ Ví dụ Bước 0: Dùng pt đầu tiên để khử x1 trong 2 pt còn lại. Để khử x1 ở hàng thứ k (k = 2,3) các hệ số ak,j ở hàng thứ k (j = 1,2,3): ak,j = ak,j – a1,j*ak,1/a1,1 hệ số bk ở hàng thứ k: bk = bk – b1*ak,1/a1,1 113 114 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss Phương pháp khử Gauss Ví dụ Ví dụ Bước 1: Dùng pt thứ 2 để khử x2 trong pt thứ 3. Để khử x2 ở hàng thứ 3 các hệ số a3,j ở hàng thứ 3 (j = 2,3): a3,j = a3,j – a2,j*a3,2/a2,2 hệ số b3 ở hàng thứ n: b3 = b3 – b2*a3,2/a2,2 115 116 29
  30. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss Phương pháp khử Gauss Ví dụ Chương trình Procedure GAUSS(A:ma;B:mX;Var X:mX;nF:integer); Begin End; Để giải hệ phương trình trước hết cần biết: -Số phương trình và ẩn số nF -Giá trị các phần tử của ma trận hệ số A -Giá trị các phần tử của ma trân hệ số tự do B 117 118 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss Phương pháp khử Gauss Ví dụ Ví dụ Program HTT1; Program HTT1; uses crt; Type Procedure GAUSS(A:ma;B:mX;Var X:mX;nF:integer); mX= Begin ma= Var End; X,B:mX; {Chương trình chính} A:ma; BEGIN nF,i,j,k:integer; clrscr; 119 {Nhập số ẩn số và phương trình} 120 30
  31. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss Phương pháp khử Gauss Ví dụ Ví dụ Program HTT1; Program HTT1; {Chương trình chính} {Chương trình chính} BEGIN BEGIN clrscr; {Nhập số ẩn số và phương trình} {Nhập ma trận hệ số A} write (‘Số ẩn số nF = ’);readln(nF); For i:=1 to nF do {Nhập ma trận hệ số tự do} For j:=1 to nF do For i:=1 to nF do readln(A[i,j]); readln(B[i]); 121 122 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss Phương pháp khử Gauss-Jordan Ví dụ Program HTT1; Là phương pháp khử dần các ẩn để đưa hệ phương trình đã cho về dạng ma trận đường chéo rồi giải. {Chương trình chính} không phải tính định thức BEGIN GAUSS(A,B,X,nF); {In kết quả} For i:=1 to nF do writeln (‘X[’,i,‘] = ’,X[i]:8:4); readln; END. 123 124 31
  32. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss-Jordan Phương pháp khử Gauss-Jordan Các bước thực hiện Bước 1: Dùng pt đầu tiên để khử x1 trong n-1 pt còn lại. Để khử x1 ở hàng thứ k (k = 2,3, ,n) tính lại các hệ số ak,j ở hàng thứ k (j = 1,2, ,n): ak,j = ak,j – a1,j*ak,1/a1,1 và tính lại hệ số bk ở hàng thứ k: bk = bk – b1*ak,1/a1,1 125 126 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss-Jordan Phương pháp khử Gauss-Jordan Các bước thực hiện Các bước thực hiện Bước i: Dùng pt thứ i để khử xi trong (n-1) pt còn lại. Bước n: Dùng pt thứ n để khử xn trong (n-1) pt còn lại. Để khử xi ở hàng thứ k (k = 1,2, ,i-1,i+1,i+2, ,n) Để khử xn ở hàng thứ k (k = 1,2, ,n-1) tính lại các hệ số ak,j ở hàng thứ k (j = i,i+1, ,n): tính lại các hệ số ak,j ở hàng thứ k (j = n): ak,j = ak,j – ai,j*ak,i/ai,i ak,j = ak,j – an,j*ak,n/an,n và tính lại hệ số bk ở hàng thứ k: và tính lại hệ số bk ở hàng thứ k: bk = bk – bi*ak,i/ai,i bk = bk – bn*ak,n/an,n 127 128 32
  33. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss-Jordan Phương pháp khử Gauss-Jordan Sau khi khử hệ phương trình có dạng Sau khi khử hệ phương trình có dạng Dạng 1: Nếu tại các bước (bước i) không chia cho hệ số ai,i trước Dạng 2: Nếu tại các bước (bước i) chia cho hệ số ai,i trước khi thực khi thực hiện quá trình khử. hiện quá trình khử. 129 130 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss-Jordan Phương pháp khử Gauss-Jordan Ví dụ Ví dụ Bước 1: Dùng pt đầu tiên để khử x1 trong 2 pt còn lại. Để khử x1 ở hàng thứ k (k = 2,3) các hệ số ak,j ở hàng thứ k (j = 1,2,3): ak,j = ak,j – a1,j*ak,1/a1,1 hệ số bk ở hàng thứ k: bk = bk – b1*ak,1/a1,1 131 132 33
  34. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss-Jordan Phương pháp khử Gauss-Jordan Ví dụ Ví dụ Bước 2: Dùng pt thứ 2 để khử x2 trong 2 pt còn lại. Để khử x2 ở hàng thứ k (k = 1,3) các hệ số ak,j ở hàng thứ k (j = 2,3): ak,j = ak,j – a2,j*ak,2/a2,2 hệ số bk ở hàng thứ k: bk = bk – b2*ak,2/a2,2 133 134 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss-Jordan Phương pháp khử Gauss-Jordan Ví dụ Chương trình Bước 3: Dùng pt thứ 3 để khử x trong 2 pt còn lại. Procedure GauJor(A:ma;Var X:mX;nF:integer); 3 Begin Để khử xn ở hàng thứ k (k = 1,2) End; các hệ số ak,j ở hàng thứ k (j = n): ak,j = ak,j – a3,j*ak,3/a3,3 hệ số b ở hàng thứ k: b = b – b *a /a k k k 3 k,3 3,3 Để giải hệ phương trình trước hết cần biết: -Số phương trình và ẩn số nF -Giá trị các phần tử của ma trận hệ số mở rộng A[i,nF+1] 135 136 34
  35. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp khử Gauss-Jordan Phương pháp khử Gauss-Jordan Ví dụ: Ví dụ: Program HTT2; Program HTT2; uses crt; {Chương trình chính} procedure GauJor(A:ma;Var X:mX;n:integer); BEGIN Begin For i:=1 to nF do End; For j:=1 to (nF+1) do {Chương trình chính} readln(A[i,j]); BEGIN GauJor(A,X,nF); 137 138 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp nghịch đảo ma trận Phương pháp nghịch đảo ma trận Tìm ma trận nghịch đảo bằng phương pháp Gauss-Jordan: Bước 1: Viết thêm ma trận đơn vị vào bên phải ma trận hệ số: Ma trận nghịch đảo có thể được xác định bằng phương pháp khử Gauss-Jordan 139 140 35
  36. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp nghịch đảo ma trận Phương pháp nghịch đảo ma trận Tìm ma trận nghịch đảo bằng phương pháp Gauss-Jordan: Tìm ma trận nghịch đảo bằng phương pháp Gauss-Jordan: Bước 2: Sử dụng phương pháp Gauss-Jordan thực hiện phép biến đổi Bước 3: Ma trận bên phải sẽ là ma trận nghịch đảo của ma trận hệ cơ sở trên các hàng của ma trận đến khi ma trận có dạng: số: 141 142 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp nghịch đảo ma trận Phương pháp nghịch đảo ma trận Chương trình Chương trình Procedure DMT(Var A:ma;nF:integer); Do đó: để tìm nghiệm phải thực hiện phép nhân hai ma trận: Begin End; Để giải hệ phương trình trước hết cần biết: -Số phương trình và ẩn số nF -Giá trị các phần tử của ma trận hệ số A -1 Kết quả trả lại là ma trận nghịch đảo A 143 144 36
  37. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp nghịch đảo ma trận Phương pháp nghịch đảo ma trận Chương trình nhân hai ma trận Ví dụ: Program HTT3; For i:=1 to nF do use crt; Begin x[i]:=0.0; Procedure DMT(Var A:ma;nF:integer); For j:=1 to nF do Begin x[i]:=x[i]+A[i,j]*B[j]; End; End; {Chương trình chính} BEGIN 145 146 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Phương pháp nghịch đảo ma trận Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Ví dụ: Program HTT3; {Chương trình chính} BEGIN {Nhập nF, A, B} DMT(A,nF); {Nhân hai ma trận} Hồi quy thực nghiệm: {In kết quả} Xây dựng hàm toán học tường minh mô tả END. 147 chính xác nhất bộ số liệu thực nghiệm 148 37
  38. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Các dạng hàm số thường xuất hiện trong kỹ thuật hóa học Hồi quy tuyến tính: Xây dựng hàm tuyến tính mô tả chính xác nhất bộ số liệu thực nghiệm. Phương pháp xây dựng: Phương pháp bình phương tối thiểu (least square method) 149 150 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Hồi quy tuyến tính Tổng bình phương sai số giữa Hồi quy tuyến tính dự đoán và thực nghiệm Các hệ số a và b thích hợp nhất khi tổng bình phương sai số là nhỏ nhất Phương pháp tìm cực tiểu: 151 152 38
  39. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Hồi quy tuyến tính Hồi quy tuyến tính 153 154 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Hồi quy đa thức Hồi quy đa thức Tổng bình phương sai số giữa dự đoán và thực nghiệm 155 156 39
  40. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Hồi quy đa thức Hồi quy đa thức 157 158 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Hồi quy đa thức Hồi quy đa thức Hệ số tương quan: Là hệ số đánh giá tính tương hợp của hàm toán được xây dựng -Không tương hợp: r2 0,8 -Thích hợp: r2 1 159 160 40
  41. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Ví dụ: Xây dựng hàm toán học mô tả bộ số liệu thực nghiệm Ví dụ: Xây dựng hàm toán học mô tả bộ số liệu thực nghiệm 161 162 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng ứng dụng Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Ví dụ: Xây dựng hàm toán học mô tả bộ số liệu thực nghiệm Ví dụ: Xây dựng hàm toán học mô tả bộ số liệu thực nghiệm 163 164 41
  42. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.1 Phương pháp giải hệ phương trình tuyến tính và ứng dụng Ứng dụng: Xây dựng hàm hồi quy thực nghiệm Ví dụ: Xây dựng hàm toán học mô tả bộ số liệu thực nghiệm 165 166 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Nghiệm thực của phương trình – Ý nghĩa hình học. Sự tồn tại của nghiệm thực f(x) = 0; ( 1 ) y Định lý. Nếu có hai số thực a, b y f – hàm cho trước của đối số x f(x) (a < b) sao cho f(a) và f(b) trái B dấu, tức là α - nghiệm thực của ( 1 ) O M f(a).f(b) < 0 ( 3 ) O a f(α) = 0; ( 2 ) α x đồng thời f(x) liên tục trên [a, b] b x - Vẽ đồ thị y = f(x) y g(x) thì trong khoảng [a, b] ít nhất có A - hoặc (1)~ g(x) = h(x) M một nghiệm thực của phương h(x) đồ thị y1 = g(x) và y2 = h(x) trình f(x) = 0. O Hoành độ điểm M nghiệm α. α 167 x 168 42
  43. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Khoảng phân ly nghiệm (tách nghiệm) Các phương pháp xác định gần đúng nghiệm thực của phương trình phi tuyến Định nghĩa: Khoảng [a, b] nào y đó gọi là khoảng phân ly nghiệm B của phương trình f(x) = 0 nếu nó 1. Phương pháp đồ thị. 2. Phương pháp thử. chứa một và chỉ một nghiệm O a của phương trình đó. b x 3. Phương pháp chia đôi. - hàm f(x) đơn điệu 4. Phương pháp lặp. trong [a, b] : A f’(x) không đổi dấu 5. Phương pháp tiếp tuyến (phương pháp Newton-Raphson). Định lý: Nếu hàm số f(x) liên tục và đơn điệu trên khoảng [a, b], 6. Phương pháp dây cung. đồng thời f(a) và f(b) trái dấu thì [a, b] là khoảng phân ly nghiệm của phương trình f(x) = 0. 169 170 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Phương pháp Newton-Raphson (phương pháp tiếp tuyến) Phương pháp Newton-Raphson (phương pháp tiếp tuyến) - Giả sử f(x) =0 : - Có nghiệm thực α phân ly trong [a, b]; Cơ sở : khai triển Taylor: - Có đạo hàm f’(x) ≠ 0 tại x [a, b]; - Hàm F(x) xác định và có đạo hàm đến cấp n+1 tại xo và lân - Có đạo hàm cấp hai f’’(x) tại x [a, b]; cận xo. - Chọn xo [a,b] khai triển Taylo bậc nhất của f(x) tại xo: - Khai triển Taylo bậc n của F(x) tại xo: Bỏ qua số hạng cuối . . . 171 172 43
  44. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Phương pháp Newton-Raphson (phương pháp tiếp tuyến) Phương pháp Newton-Raphson (phương pháp tiếp tuyến) Ý nghĩa hình học: thay đường cong y = f(x) bằng tiếp tuyến kẻ từ A(a,f(a)) hay B(b,f(b)), hoành độ giao điểm của tiếp tuyến với trục hoành là nghiệm gần đúng của phương trình. y A Đặt: - xo = a, nếu tiếp tuyến kẻ từ A; Tiếp tục vẽ tiếp tuyến với điểm [ (x1, f(x1) ] - xo = b, nếu tiếp tuyến kẻ từ B; . . . Phương trình tiếp tuyến của y = f(x) tại [xo, f(xo)] : O α b ( a ) xo=a x1 x2 x B Giao điểm với trục hoành (x1, y1 = 0) ( b173 ) 174 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Phương pháp Newton-Raphson (phương pháp tiếp tuyến) Phương pháp Newton-Raphson (phương pháp tiếp tuyến) * Định lý về sự hội tụ. Giả sử: - [a, b] là khoảng phân ly nghiệm của f(x) = 0; y - f có đạo hàm f’, f” với f’ và f” : + liên tục trên [a, b]; + không đổi dấu trên [a, b]; Tiếp tục vẽ tiếp tuyến với điểm [ (x1, f(x1) ] - xấp xỉ xo chọn f(xo).f”(xo) > 0; . . . - xn α khi A x O 2 x1 xo=b * Sai số. Lấy xn nghiệm gần đúng sai số: a α x với Trong thực tế, thường dừng quá trình tính khi: B 175 xn – xn-1 < ε (sai số cho phép) 176 44
  45. Ví dụ. Tìm nghiệm gần đúng của phương trình Chương 1. Các phương pháp giải phương f(x) = x3 – x – 1 = 0; Với sai số cho phép ε =10-3 trình và hệ phương trình 1. Điều kiện phương trình có nghiệm thực và tìm khoảng phân ly nghiệm. 1.2 Phương pháp giải phương trình và hệ phương trình -Hàm số xác định và liên tục phi tuyến x tại mọi x Phương pháp Newton-Raphson (phương pháp tiếp tuyến) f’(x) + 0 _ 0 + - f’(x) = 3x2 – 1 = 0 tại Sơ đồ tóm tắt các bước giải: M f(x) m 1/ Cho phương trình f(x) = 0; 2/ Ấn định sai số cho phép ε; - Bảng biến thiên hàm số: 3/ Tìm khoảng phân ly nghiệm; 4/ Chọn giá trị đầu xo: f(xo).f”(xo) > 0; đồ thị hàm số chỉ cắt trục hoành tại một điểm, phương trình có một nghiệm thực trong khoảng 5/ Tính toán sai số - Chọn khoảng chứa nghiệm [1, 2] e 0 chứa nghiệm. x = α 177 178 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Phương pháp Newton-Raphson (phương pháp tiếp tuyến) Phương pháp Newton-Raphson (phương pháp tiếp tuyến) Lập bảng tính: Phương trình: x f(x)=x3-x-1 f (x)=3x2-1 f’(x) = 3x2 – 1 > 0 trong khoảng [1, 2] 2,0 5,0 11,0 1,5454545 f”(x) = 6x > 6 trong khoảng [1, 2] 1,5454545 1,145755 6,165288 1,3596148 f(1) = -1; f(2) = 5; 1,3596148 0,153704 4,545657 1,3258015 1,3258015 0,004625 4,273245 1,3247190 1,3247190 0,0000034 4,264641 1,3247182 f(2).f”(2) > 0 Chọn đầu tính x = 2. 179 1,3247182 0,0000010 180 45
  46. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Phương pháp Newton-Raphson (phương pháp tiếp tuyến) Phương pháp Newton-Raphson (phương pháp tiếp tuyến) Chương trình: Ví dụ: Function F(x:real):real; f(x) = x3 – x – 1 = 0 Program PT1; Begin Uses crt; f’(x) = 3x2 – 1 = 0 F:=x*sqr(x)-x-1; Var End; n,i,j,k:integer; Function dF(x:real):real; x,x0,eps,ss:real; Begin Function F(x:real):real; dF:=3*sqr(x)-1; Begin End; F:=x*sqr(x)-x-1; End; 181 182 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Phương pháp Newton-Raphson (phương pháp tiếp tuyến) Phương pháp Newton-Raphson (phương pháp tiếp tuyến) Ví dụ: Ví dụ: Program PT1; Program PT1; Function dF(x:real):real; {Chương trình chính} Begin BEGIN dF:=3*sqr(x)-1; clrscr; End; write (‘Nhập x0 = ’);readln(x0); {Chương trình chính} write (‘Nhập eps = ’);readln(eps); BEGIN x:=x0;k:=0; clrscr; Repeat 183 184 46
  47. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Phương pháp Newton-Raphson (phương pháp tiếp tuyến) Phương pháp Newton-Raphson (phương pháp tiếp tuyến) Ví dụ: Ví dụ: Program PT1; Program PT1; BEGIN BEGIN Repeat x:=x0;k:=0; Repeat Until ss<=eps; x:=x-F(x)/dF(x); {In kết quả} ss:=abs(x-x0); write (‘Nghiệm x = ’,x:8:4); x0:=x;k:=k+1; readln; Until ss<=eps; 185 END. 186 187 188 47
  48. 189 190 Chương 1. Các phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến Phương pháp lặp Cho phương trình f(x) = 0 có nghiệm thực trong khoảng [a,b]; - Viết lại x + f(x) – x = 0; Đặt φ(x) = x + f(x); x = φ(x); - Chọn sơ bộ giá trị gần đúng đầu tiên của nghiệm: x1 = φ(xo); x = φ(x ); . 2. . . . . .1 . . xn = φ(xn-1); n = 1, 2, . . . (*) - Hàm φ(x) gọi là hàm lặp. 191 192 48
  49. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Phương pháp lặp Phương pháp lặp Sự hội tụ của quá trình tính toán Ý nghĩa hình học x1 = φ(xo); - Giả sử khi n ; xn nghiệm α của phương trình f(x) = 0; x + f(x) – x = 0; x = φ(x); x2 = φ(x1); phương pháp lặp hội tụ, có thể coi xn là nghiệm gần đúng . . . . . . . . . . Đặt φ(x) = x + f(x); xn = φ(xn-1); n = 1, 2, . . . y -Quá trình tính cũng có thể phân kỳ, xn ngày càng đi xa khỏi y nghiệm. y=φ(x) O O 193 α x3x2 x1 xo x α x0 x1 x2 x3 x194 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Phương pháp lặp Phương pháp lặp Sai số của phép tính: Định lý về sự hội tụ. Giả sử: - [a, b] là khoảng phân ly nghiệm α của phương trình f(x) = 0; Ước lượng sai số: - Mọi xn tính theo (*) đều - Hàm φ(x) có đạo hàm hạng nhất thoả mãn điều kiện: Chú ý: - Nếu φ’(x) > 0; Có thể chọn xo bất kỳ - Nếu φ’(x) < 0: q - hằng số; khi - Phương pháp lặp hội tụ với mọi x xét dấu x α khi n n khi 195 196 49
  50. Ví dụ. Tìm nghiệm gần đúng của phương trình Chương 1. Các phương pháp giải phương f(x) = x3 – x – 1 = 0; Với sai số cho phép ε =10-3 trình và hệ phương trình 1. Điều kiện phương trình có nghiệm thực và tìm khoảng phân ly nghiệm. 1.2 Phương pháp giải phương trình và hệ phương trình -Hàm số xác định và liên tục phi tuyến x tại mọi x Phương pháp lặp f’(x) + 0 _ 0 + - f’(x) = 3x2 – 1 = 0 tại f(x) M x=xo m - Bảng biến thiên hàm số: y = φ(x) đồ thị hàm số chỉ cắt trục hoành tại một điểm, phương trình s có một nghiệm thực trong khoảng x = y y – x 0 α = y 197 chứa nghiệm.198 3 f(x) = x – x – 1 = 0; Lập bảng tính: - Đặt x = φ(x) = x3 – 1 ; 2 φ’(x) = 3x ≥ 3 tại mọi x trong khoảng [1, 2] không No x φ(x) x – x đảm bảo điều kiện hội tụ. i+1 i x 1,0 1,259921 - Đặt x = φ(x) = x – λf(x) với o -0,259921 x 1,259921 1,3122938 1 -0,052373 2 x 1,3122938 1,3223538 f’(x) = 3x -1 f’(2) = 11; 2 -0,010060 x 1,3223538 1,3242687 3 -0,001915 x 1,3242687 1,3242826 4 -0,00185 x = 1 φ’(x) = 1-(2/11) <1; x = 2 φ’(x) = 0 <1 x 1,3242826 1,3246326 5 -0,000364 3 - Hoặc đặt x = x + 1; x6 1,3246326 Lấy α = 1,3246326 sai số sẽ < ε = 10-3. tại mọi đảm bảo điều kiện hội tụ. 199 200 50
  51. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Phương pháp lặp Phương pháp lặp Chương trình: Ví dụ: Function Fi(x:real):real; f(x) = x3 – x – 1 = 0 Program PT2; Begin Uses crt; 3 Fi:=exp(1/3*ln(x+1)); x = φ(x) = x – 1 Var End; n,i,j,k:integer; x,x0,eps,ss:real; Function Fi(x:real):real; Begin Fi:=exp(1/3*ln(x+1)); End; 201 202 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Phương pháp lặp Phương pháp lặp Ví dụ: Ví dụ: Program PT2; Program PT2; {Chương trình chính} {Chương trình chính} BEGIN BEGIN clrscr; clrscr; write (‘Nhập x0 = ’);readln(x0); write (‘Nhập eps = ’);readln(eps); x:=x0;k:=0; Repeat 203 204 51
  52. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Phương pháp lặp Phương pháp lặp Ví dụ: Ví dụ: Program PT2; Program PT2; BEGIN BEGIN Repeat x:=x0;k:=0; Repeat Until ss<=eps; x:=F(x); {In kết quả} ss:=abs(x-x0); write (‘Nghiệm x = ’,x:8:4); x0:=x;k:=k+1; readln; Until ss<=eps; 205 END. 206 207 208 52
  53. 209 210 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Giải hệ phương trình phi tuyến bằng phương pháp Newton Giải hệ phương trình phi tuyến bằng phương pháp Newton Phương pháp Newton có thể tổng quát hóa để giải hệ phương trình phi Công thức Newton với phương trình 1 biến: tuyến có dạng: Hay: Với: Dạng ma trận: Trong đó: 211 212 53
  54. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Giải hệ phương trình phi tuyến bằng phương pháp Newton Giải hệ phương trình phi tuyến bằng phương pháp Newton Đối với hệ phương trình phi tuyến, công thức Newton tổng quát: Phương pháp Newton giải hệ phương trình phi tuyến là phương pháp tuyến tính hóa hệ phương trình đã cho thành một hệ phương trình tuyến tính mà biến số của hệ là X. Như vậy ở mỗi bước lặp (bước thứ i), cần phải giải một hệ phương trình Trong đó J(Xi) là ma trận (toán tử) Jacobi. Nó là ma tuyến tính với biến số là Xi cho đến khi được nghiệm gần đúng. trận cấp n có dạng: Và: Vì vậy: việc giải hệ phương trình phi tuyến bằng phương pháp Newton là lặp lại việc giải hệ phương trình tuyến tính: 213 214 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Giải hệ phương trình phi tuyến bằng phương pháp Newton Giải hệ phương trình phi tuyến bằng phương pháp Newton - Do đó việc giải một hệ phi tuyến bằng phương pháp Newton, chính là Thuật toán: việc giải hệ phương trình tuyến tình với: 1 Chọn giá trị đầu X0: 215 216 54
  55. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Giải hệ phương trình phi tuyến bằng phương pháp Newton Giải hệ phương trình phi tuyến bằng phương pháp Newton Thuật toán: Thuật toán: 2 Giải hệ phương trình tuyến tính (Gauss hoặc Gauss-Jordan): Procedure HAM(X:mX; nF:integer; Var F:mX); Begin F[1]:= ; F[2]:= ; F[nF]:= ; End; 3 Kiểm tra sai số: 217 218 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Giải hệ phương trình phi tuyến bằng phương pháp Newton Giải hệ phương trình phi tuyến bằng phương pháp Newton Thuật toán: Ví dụ: Procedure DHAM(X:mX; Var A:mA); Giải hệ phương trình phi tuyến Begin A[1,1]:= ; A[1,2]:= ; A[1,nF]:= ; A[nF,1]:= ; A[nF,2]:= ; A[nF,nF]:= ; End; 219 220 55
  56. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Giải hệ phương trình phi tuyến bằng phương pháp Newton Giải hệ phương trình phi tuyến bằng phương pháp Newton Ví dụ: Ví dụ: Program HFT1; Program HFT1; uses crt; Type Procedure HAM(X:mX; nF:integer; Var F:mX); mX = ; mA = ; Begin Var F[1]:=-2*x[1]+exp(x[2]); X0,X,dX,B,F:mX; F[2]:=-exp(-x[1])-5*x[2]; A:mA; End; nF,i,j,k:integer; dXmax,eps:real; 221 222 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Giải hệ phương trình phi tuyến bằng phương pháp Newton Giải hệ phương trình phi tuyến bằng phương pháp Newton Ví dụ: Ví dụ: Program HFT1; Program HFT1; Procedure DHAM(X:mX; Var A:mA); BEGIN Begin clrscr; A[1,1]:=2; A[1,2]:=-exp(x[2]); writeln (‘Nhập các giá trị đầu X0:’); A[2,1]:=-exp(-x[1]); A[2,2]:=5; For i:=1 to nF do End; readln(x0[i]); {Chương trình chính} For j:=1 to nF do x[j]:=x0[j]; 223 224 56
  57. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình phi tuyến phi tuyến Giải hệ phương trình phi tuyến bằng phương pháp Newton Giải hệ phương trình phi tuyến bằng phương pháp Newton Ví dụ: Ví dụ: Program HFT1; Program HFT1; k:=0; k:=0; Repeat Repeat HAM(X,nF,F); For j:=1 to nF do For j:=1 to nF do B[j]:=F[j]; X[j]:=X[j]+dX[j]; DHAM(X,A); For j:=1 to nF do GAUSS(A,B,dX,nF); if dX[j]>=dXmax then 225 dXmax:=dX[j]; 226 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.2 Phương pháp giải phương trình và hệ phương trình 1.3 Ứng dụng phi tuyến Giải hệ phương trình phi tuyến bằng phương pháp Newton Tính toán thiết bị trao đổi nhiệt Ví dụ: Program HFT1; Repeat k:=k+1; Until dXmax<=eps; {In kết quả} For i:=1 to nF do writeln (‘X[’,i,‘] = ’,X[i]); readln; END. 227 228 57
  58. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Dung dịch Tính toán thiết bị trao đổi nhiệt Tính toán thiết bị trao đổi nhiệt sau gia nhiệt Hơi nước bão hòa Nước ngưng Dung dịch cần gia 229 nhiệt230 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán thiết bị trao đổi nhiệt Tính toán thiết bị trao đổi nhiệt Phương trình trao đổi nhiệt cơ bản: Cơ chế trao đổi nhiệt Biến thiên nhiệt độ trong thiết bị TĐN 231 232 58
  59. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán thiết bị trao đổi nhiệt Tính toán thiết bị trao đổi nhiệt Xác định diện tích bề mặt trao đổi Để thuận tiện trong tính toán thiết nhiệt và số ống của chùm ống: bị trao đổi nhiệt, sử dụng đại lượng mật độ dòng nhiệt: Giả thiết trao đổi nhiệt ổn định: Biến thiên nhiệt độ trong Cơ chế trao đổi nhiệt thiết bị TĐN 233 234 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán thiết bị trao đổi nhiệt Tính toán thiết bị trao đổi nhiệt Giả thiết trao đổi nhiệt ổn định: Phương pháp tính: - Chọn sơ bộ kích thước ống chùm: Chiều dài , đường kính -Chọn sơ bộ chế độ chuyển động: Thường chọn trước Re ở chể độ chảy rối (xoáy) -Xác định các hệ số cấp nhiệt trên Cơ chế trao đổi nhiệt hai bề mặt của ống. Cơ chế trao đổi nhiệt 235 236 59
  60. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán thiết bị trao đổi nhiệt Tính toán thiết bị trao đổi nhiệt Phương pháp tính: Phương pháp tính: - Chọn sơ bộ kích thước ống chùm: Chọn giá trị đầu: Chiều dài , đường kính Xác định: -Chọn sơ bộ chế độ chuyển động: Xác định: theo: Thường chọn trước Re ở chể độ chảy rối (xoáy) -Xác định các hệ số cấp nhiệt trên hai bề mặt của ống. Cơ chế trao đổi nhiệt Xác định: Cơ chế trao đổi nhiệt 237 Kiểm tra sai số: 238 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán thiết bị trao đổi nhiệt Tính toán thiết bị trao đổi nhiệt Phương pháp tính: Phương pháp tính: Hệ số cấp nhiệt 1 phía hơi nước bão hòa Hệ số cấp nhiệt 1 phía hơi nước bão hòa Ngưng tụ dọc Ngưng tụ dọc theo vách đứng theo vách đứng Chú ý: Tra sổ tay (Dùng thuật toán nội suy) Ẩn nhiệt (J/kg), khối lượng riêng (kg/m3), hệ số dẫn nhiệt (W/m), độ nhớt (Ns/m2) của nước ngưng Nhiệt độ màng nước ngưng 239 240 60
  61. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán thiết bị trao đổi nhiệt Tính toán thiết bị trao đổi nhiệt Phương pháp tính: Thuật toán: Hệ số cấp nhiệt 2 phía dung dịch 1) Xác định q1: Giả thiết: 241 242 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán thiết bị trao đổi nhiệt Tính toán thiết bị trao đổi nhiệt Thuật toán: Thuật toán: 2) Xác định q2: 3) So sánh q1 và q2: Giả thiết: 4) Xác định mật độ truyền nhiệt q: 243 244 61
  62. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán thiết bị trao đổi nhiệt Tính toán thiết bị trao đổi nhiệt Thuật toán: Chương trình: NOISUY(tm,lamNN,lamN,tlamN); 5) Xác định diện tích bề mặt NOISUY(tm,RoNN,RoN,tRoN); trao đổi nhiệt F: NOISUY(tm,MuyNN,MuyN,tMuyN); NOISUY(t1,rNN,rN,tRN); 6) Xác định số ống chùm n: 245 246 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán thiết bị trao đổi nhiệt Tính toán thiết bị trao đổi nhiệt Chương trình: Chương trình: Program TBN1; Program TBN1; uses crt; type Procedure ALPHA1(var Alfa1:real; deltat1:real); mX= var Procedure ALPHA2(var Alfa2:real; deltat2:real); lamN,tlamN, :mX; {Chương trình chính} lamNN, ,tm,G,Cpdd, ,q1,q2,qtb:real; BEGIN nCB,i,j,k:integer; clrscr; Procedure NOISUY(xs:real; var ys: real; Y,X:mX); {Nhập G, t2d,t2c,t1,deltaT,lamdaT,Cpdd,r1,r2} 247 248 62
  63. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán thiết bị trao đổi nhiệt Tính toán thiết bị trao đổi nhiệt Chương trình: Chương trình: Program TBN1; Program TBN1; rT:=r1+r2+deltaT/lamdaT; Repeat deltaTtb:=((t1-t2d)-(t1-t2c))/ln((t1-t2d)/(t1- ALPHA1(Alfa1,deltat1); t2c)); q1:=Alfa1*deltat1; Q:=G*Cpdd*(t2c-t2d); tT1:=t1-deltat1; t2:=t1-deltaTtb; deltaT:=q1*rt; deltat1:=0.5; tT2:=tT1-deltaT; Repeat deltat2:=tT2-t2; 249 250 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán thiết bị trao đổi nhiệt Tính toán thiết bị trao đổi nhiệt Chương trình: Chương trình: Program TBN1; Program TBN1; Repeat F:=Q/q; no:=F/(3.14*dotb*H); ALPHA2(Alfa2,deltat2); {In kết quả: Diện tích F, số ống no} q2:=Alfa2*deltat2; readln; if q1<q2 then deltat1:=deltat1+0.05 END. else deltat1:=deltat1-0.05; Until abs(q1-q2)/q1<=0.03; q:=(q1+q2)/2; 251 252 63
  64. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán thiết bị trao đổi nhiệt Tính toán thiết bị trao đổi nhiệt Chương trình: Bài tập: Procedure ALPHA1(Var Alfa1:real;deltat1:real); Viết chương trình xác định diện tích bề mặt trao đổi nhiệt của thiết bị Begin trao đổi nhiệt ống chùm kiểu đứng để đun nóng dung dịch Benzene o NOISUY(tm,lamNN,lamN,tlamN); với năng suất 5,5 tấn/h từ 25 đến 75 C bằng hơi nước bão hòa ở nhiệt độ 120oC. NOISUY(tm,RoNN,RoN,tRoN); Cho trước: NOISUY(tm,MuyNN,MuyN,tMuyN); 1- Chiều ống chùm: H=2m; Chiều dày ống: T = 2mm; Vật liệu: thép CT3. NOISUY(t1,rNN,rN,trN); 2 2- Hệ số cấp nhiệt phía dung dịch: 2 = 397 W/m TU:=rNN*sqr(RoNN)*lamNN*sqr(lamNN); Các số liệu khác tra sổ tay Quá trình và thiết bị công nghệ hóa chất T1,2. MAU:=MuyNN*deltat1*H; Alfa1:=2.04*sqrt(sqrt(TU/MAU)); End; 253 254 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán hệ thống thiết bị cô đặc Tính toán hệ thống thiết bị cô đặc 255 256 64
  65. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán hệ thống thiết bị cô đặc Tính toán hệ thống thiết bị cô đặc Mục đích của việc tính toán hệ thống cô đặc (nhiều nồi liên tiếp – Cơ sở tính toán: multi-effect evaporation): - Xây dựng phương trình cân bằng chất cho từng nồi và cho - Xác định các đại lượng D, W1, W2, W3 để đảm bảo hệ thống 1) Nâng cao nồng độ dung dịch cần cô đặc từ ađ đến ac - Xây dựng phương trình cân bằng nhiệt (năng lượng) cho từng nồi và cho hệ thống 2) Đảm bảo đủ khả năng trao đổi nhiệt từ hơi đốt D và hơi thứ Wi trong từng thiết bị cô đặc. - Kết hợp với một số giả thiết nhằm đơn giản hóa mô hình - Dựa vào hai lựa chọn chính: 1) Diện tích bề mặt trao đổi nhiệt trong các thiết bị là bằng nhau Hệ phương trình tuyến tính với các ẩn số: D, Wi 2) Tổng diện tích bề mặt trao đổi nhiệt là nhỏ nhất 257 258 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán hệ thống thiết bị cô đặc Tính toán hệ thống thiết bị cô đặc Cân bằng chất: Cân bằng chất: TB1 (n=1) TB2 (n=2) 259 260 65
  66. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán hệ thống thiết bị cô đặc Tính toán hệ thống thiết bị cô đặc Cân bằng chất: Cân bằng nhiệt: TB3 (n=3) TB1 (n=1) Thông thường dung dịch được gia nhiệt đến nhiệt độ sôi trước khi đưa vào cô đặc: 261 262 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán hệ thống thiết bị cô đặc Tính toán hệ thống thiết bị cô đặc Cân bằng nhiệt: Kết hợp cân bằng chất và cân bằng nhiệt: TB2 (n=2) TB3 (n=3) 263 264 66
  67. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán hệ thống thiết bị cô đặc Tính toán hệ thống thiết bị cô đặc Để giải được hệ phương trình cân bằng vật chất và năng lượng: 1. Xác định áp suất làm 2. Xác định lượng hơi 1. Xác định áp suất làm việc trong từng thiết bị việc trong từng thiết thứ trong từng thiết bị P1, P2, P3 bị W1, W2, W3 2. Xác định lượng hơi thứ trong từng thiết bị Xác định nhiệt độ hơi Xác định nồng độ chất thứ trong từng thiết bị tantrong từng thiết bị tHT1, tHT2, tHT3 a1, a2, a3 265 266 Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình 1.3 Ứng dụng 1.3 Ứng dụng Tính toán hệ thống thiết bị cô đặc Tính toán hệ thống thiết bị cô đặc 1. Xác định áp suất làm Giả thiết ban đầu việc trong từng thiết bị P1, P2, P3 Ảnh hưởng Giả thiết lại của nồng độ Hiệu nhiệt độ hữu ích Kiểm tra giả thiết Ảnh hưởng của áp suất Ảnh hưởng do thủy tĩnh ma sát trên đường ống 267 268 67
  68. Chương 1. Các phương pháp giải phương Chương 1. Các phương pháp giải phương trình và hệ phương trình trình và hệ phương trình Nhập số liệu đầu 1.3 Ứng dụng 1.3 Ứng dụng Giả thiết phân bố áp suất Tính toán hệ thống thiết bị cô đặc Tính toán hệ thống thiết bị cô đặc P1, P2, P3 Giả thiết ban đầu Giả thiết phân bố hơi 2. Xác định lượng hơi thứ W , W , W thứ trong từng thiết 1 2 3 bị W1, W2, W3 Xác định các nhiệt độ và thông số hóa lý tương Giả thiết lại ứng Các thông số Kiểm tra giả thiết Tính toán lại W1, W2, W3 Từ việc giải hệ Tính toán trao đổi nhiệt phương trình Kiểm tra hiệu CBC và NL nhiệt độ hữu ích 269 270 Chương 1. Các phương pháp giải phương Chương 2 trình và hệ phương trình Các phương pháp tính tích phân xác định 1.3 Ứng dụng Tính gần đúng các tích phân xác định Tính toán hệ thống thiết bị cô đặc - Xét tích phân xác định: Tính toán hệ thông cô đặc hai nồi xuôi chiều để cô đặc dung dịch đường sucrose: - Nếu f(x) liên tục trên [a, b] và có nguyên hàm là F(x) -Năng suất Gđ = 3000 kg/h -Nồng độ đầu ađ = 12% + thường khó khăn khi tìm nguyên hàm -Nồng độ cuối a = 60% - Thực tế: c + Hàm f(x) được cho dưới dạng bảng số. Sử dụng hơi bão hòa ở áp suất PH = 2atm Áp suất tại baromet PB = 0,2 atm 271 272 68
  69. Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định Tính gần đúng các tích phân xác định Tính gần đúng các tích phân xác định Đa thức xấp xỉ trực tiếp: -Tính gần đúng giá trị của tích phân thay hàm dưới dấu b tích phân bằng một đa thức xấp xỉ. a 273 274 Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định Tính gần đúng các tích phân xác định Tính gần đúng các tích phân xác định Đa thức Newton thứ nhất (Newton tiến): (với dx = hdt) - Chọn điểm cơ sở là điểm a (x0 = a) thì tại đó t(a) = 0 và x = b ứng với t = k; x = x0 + ht x = x0 + ht Bậc của đa thức được chọn công thức tính tương ứng. - Chia [a, b] thành n đoạn con bằng nhau bởi các nút xi: n = 0 công thức hình chữ nhật; n = 1 công thức hình thang; n = 2 công thức Simpson 1/3; xi = a + ih ; 275 n = 3 công thức Simpson 3/8; 276 69
  70. Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.1 Tính tích phân xác định bằng phương pháp hình 2.1 Tính tích phân xác định bằng phương pháp hình thang thang - Ý nghĩa hình học của công thức: Thay diện tích hình thang cong bằng - Thay f(x) bằng đa thức nội suy P (x). n diện tích của hình thang thường. - Công thức hình thang n = 1 - Đổi biến: x = x + ht dx = hdt 0 M x = x t = 0; x = x t = 1 Tích phân thứ 1: 1 0 1 M t=1 0 t=0 x0 x1 277 278 Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.1 Tính tích phân xác định bằng phương pháp hình 2.1 Tính tích phân xác định bằng phương pháp hình thang thang Tích phân thứ i+1: Ví dụ: Tính tích phân - Đã chứng minh được sai số của công thức là 279 280 70
  71. Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.1 Tính tích phân xác định bằng phương pháp hình 2.1 Tính tích phân xác định bằng phương pháp hình thang thang Ví dụ: Ví dụ: Program HT1; Program HT1; Uses crt; Var {Chương trình chính} TP,a,b,x,x0,S1,S0,hx:real; BEGIN n,i,j,k:integer; clrscr; Function F(x:real):real; {nhập a,b,n} Begin hx:=(b-a)/n; F:=1/(1+sqr(x)); S0:=(F(a)+F(b))/2; End; S1:=0; 281 282 Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.1 Tính tích phân xác định bằng phương pháp hình 2.2 Tính tích phân xác định bằng phương pháp Simpson thang Ví dụ: - Chia [a, b] thành 2n đoạn con bởi các nút xi. Program HT1; For i:=1 to (n-1) do Begin - Cho hàm f(x): x:=a+i*hx; S1:=S1+F(x); End; TP:=hx*(S0+S1); {In kết quả} END. 283 284 71
  72. Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.2 Tính tích phân xác định bằng phương pháp Simpson 2.2 Tính tích phân xác định bằng phương pháp Simpson - f(x) đa thức nội suy Newton bậc 2: 2 0 - Đổi biến: x = x0 + ht; dx = hdt; x = x0 t = 0; x = x2 t = 2; 285 286 Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.2 Tính tích phân xác định bằng phương pháp Simpson 2.2 Tính tích phân xác định bằng phương pháp Simpson - Các tích phân sau cũng tính tương tự - Sai số: Cộng tất cả: với với 287 288 72
  73. Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.2 Tính tích phân xác định bằng phương pháp Simpson 2.2 Tính tích phân xác định bằng phương pháp Simpson Ví dụ: Ví dụ: Tính tích phân Program Simpson1; Uses crt; Var TP,a,b,x,x0,S0,S1,S2,hx:real; n,i,j,k:integer; Function F(x:real):real; Begin F:=1/(1+sqr(x)); End; 289 290 Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.2 Tính tích phân xác định bằng phương pháp Simpson 2.2 Tính tích phân xác định bằng phương pháp Simpson Ví dụ: Ví dụ: Program Simpson1; Program Simpson1; {Chương trình chính} k:=-1; S1:=0; BEGIN For i:=1 to n do clrscr; Begin {nhập a,b,n} k:=k+2;{Tạo dãy lẻ 1,3, ,2n-1} hx:=(b-a)/(2*n); x:=a+k*hx; S0:=F(a)+F(b); S1:=S1+F(x); End; 291 292 73
  74. Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.2 Tính tích phân xác định bằng phương pháp Simpson 2.2 Tính tích phân xác định bằng phương pháp Simpson Ví dụ: Ví dụ: Program Simpson1; Program Simpson1; k:=0; S2:=0; TP:=hx*(S0+4*S1+2*S2)/3; For i:=1 to n do {In kết quả} Begin writeln (‘Tích phân I = ’,TP:8:4); k:=k+2; readln; x:=a+k*hx; END. S2:=S2+F(x); End; 293 294 Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.2 Tính tích phân xác định bằng phương pháp Simpson 2.2 Tính tích phân xác định bằng phương pháp Simpson Xác định số khoảng chia thích hợp Program SS2; Xác định số khoảng chia thích hợp Uses crt; Procedure SS(Var TP:real;n:integer); Var Begin hx:=(b-a)/(2*n); Function F(x:real):real; S0:=F(a)+F(b); Begin k:=-1; S1:=0; End; k:=0; S2:=0; Procedure SS(Var TP:real;n:integer); TP:=hx*(S0+4*S1+2*S2)/3; Begin End; End; 295 296 74
  75. Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.2 Tính tích phân xác định bằng phương pháp Simpson 2.2 Tính tích phân xác định bằng phương pháp Simpson Xác định số khoảng chia thích hợp Xác định số khoảng chia thích hợp Program SS2; Program SS2; {Chương trình chính} {In kết quả} BEGIN writeln (‘Tich phan I = ’,TP1:8:4); clrscr; writeln (‘số khoảng chia thích hợp n = ’,n); {nhập a,b,n} readln; Repeat END. SS(TP1,n); n:=2*n; SS(TP2,n); Until abs(TP1-TP2)<=eps;{ví dụ eps = 0,01} 297 298 Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.2 Ứng dụng 2.2 Ứng dụng Tính toán quá trình chưng đơn giản (batch distillation) Tính toán quá trình chưng đơn giản (batch distillation) Lượng sản phẩm đỉnh và sản phẩm đáy thay đổi theo thời gian Lượng sản phẩm đỉnh và sản phẩm đáy thay đổi theo thời gian Cân bằng vật chất viết cho Khi quá trình chưng diễn ra hệ thống trong một khoảng trong khoảng thời gian t để thời gian vi phân dt: có được xD và xW: 299 300 75
  76. Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.2 Ứng dụng 2.2 Ứng dụng Tính toán quá trình chưng đơn giản (batch distillation) Tính toán quá trình chưng đơn giản (batch distillation) Thuật toán 1. Chia đoạn [xW,xF] thành n khoảng đều nhau 2. Xác định giá trị FFi tại xi tương ứng x y* T x y* T 3. Tính tích phân bằng phương pháp hình thang hoặc Simpson 0 0 110,6 50 71,2 92,1 5 11,8 108,3 60 79 89,4 n 0 1 2 n 10 21,4 106,1 70 85,4 86,8 x xW x1 x2 xF * * * * * 20 38 102,2 80 91 84,4 y y 0 y 1 y 2 y n * 30 51,1 98,6 90 95,9 82,3 = y x 0 1 2 n Phải tính I bằng phương pháp số 301 302 40 61,9 95,2 100 100 80,2 FF = 1/ FF0 FF1 FF2 FFn Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.2 Ứng dụng 2.2 Ứng dụng Tính toán quá trình chưng đơn giản (batch distillation) Tính toán quá trình chưng đơn giản (batch distillation) Chương trình Chương trình Procedure HAM(Var FF:mX;n:integer); Program CDG; Begin Uses crt; hx:=(xF-xW)/n; Type For i:=0 to n do mX=array[1 50] of real; Begin Var xs:=xW+i*hx; FF,X,Y:mX; NOISUY(xs,yCB,Y,X); P,W,F,xF,xP,xW,hx,S0,S1,S2:real; FF[i]:=1/(yCB-xs); n,nCB,i,j,k:integer; End; Procedure NOISUY( ); End; Procedure HAM( ); 303 304 76
  77. Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.2 Ứng dụng 2.2 Ứng dụng Tính toán quá trình chưng đơn giản (batch distillation) Tính toán quá trình chưng đơn giản (batch distillation) Chương trình Chương trình Program CDG; Program CDG; {Chương trình chính} For i:=1 to (n-1) do BEGIN S1:=S1+FF[i]; clrscr; TP:=hx*(S0+S1); {Nhập F,xF,xW,nCB,X,Y}; W:=F/exp(TP); HAM(FF,n); P:=F-W; {Tính tích phân:phương pháp hình thang} xP:=(F*xF-W*xW)/P; S0:=(FF[0]+FF[n])/2; {In kết quả: P,xP} S1:=0; readln; 305 END. 306 Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.2 Ứng dụng 2.2 Ứng dụng Tính toán số đơn vị chuyển khối và chiều cao Tính toán số đơn vị chuyển khối và chiều cao tháp chưng luyện và hấp thụ loại đệm tháp chưng luyện và hấp thụ loại đệm Chương trình Số đơn vị chuyển khối: Procedure HAML(Var FF:mX;nL:integer); 1. Chưng luyện Đoạn luyện Đoạn chưng Begin hxL:=(xP-xF)/nL; yF:=DLVL(xF); hyL:=(xP-yF)/nL;{Chia khoảng cho y} For i:=1 to nL do Begin 2. Hấp thụ xs:=xF+i*hxL; NOISUY(xs,yCB,Y,X); ys:=DLVL(xs); FF[i]:=1/(yCB-ys); End; End; 307 308 77
  78. Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.2 Ứng dụng 2.2 Ứng dụng Tính toán số đơn vị chuyển khối và chiều cao Tính toán số đơn vị chuyển khối và chiều cao tháp chưng luyện và hấp thụ loại đệm tháp chưng luyện và hấp thụ loại đệm Chương trình Chương trình Procedure HAMC(Var FF:mX;nL:integer); Program TDem1; Begin Uses crt; hxC:=(xF-xW)/nC; Type yF:=DLVC(xF); mX= hyC:=(yF-xW)/nC;{Chia khoảng cho y} Var For i:=1 to nC do FF,X,Y:mX; Begin No,NoL,NoC,xs,ys,F,P,W,xF,xP,xW, :real; xs:=xW+i*hxC; nCB,nL,nC,n,i,j,k:integer; NOISUY(xs,yCB,Y,X); Procedure NOISUY( ); ys:=DLVC(xs); Function DLVL( ):real; FF[i]:=1/(yCB-ys); Function DLVC( ):real; End; Procedure HAML( ); End; Procedure HAMC( ); 309 310 Chương 2 Chương 2 Các phương pháp tính tích phân xác định Các phương pháp tính tích phân xác định 2.2 Ứng dụng 2.2 Ứng dụng Tính toán số đơn vị chuyển khối và chiều cao Tính toán số đơn vị chuyển khối và chiều cao tháp chưng luyện và hấp thụ loại đệm tháp chưng luyện và hấp thụ loại đệm Chương trình Chương trình Program TDem1; Program TDem1; {Chương trình chính} BEGIN HAMC(FF,nC); {Nhập số liệu đầu: F,xF,xP,xW,R,nCB,X,Y,nL,nC, } S0:=(FF[0]+FF[nC])/2; {Cân bằng chất} S1:=0; P:=F*(xF-xW)/(xP-xW); For i:=1 to (nC-1) do W:=F-P; S1:=S1+FF[i]; HAML(FF,nL); NoC:=hyC*(S0+S1); S0:=(FF[0]+FF[nL])/2; No:=NoL+NoC; S1:=0; {In kết quả} For i:=1 to (nL-1) do writeln (‘Số đơn vị chuyển khối No = ’,No:8:4); S1:=S1+FF[i]; readln; NoL:=hyL*(S0+S1); END. 311 312 78
  79. Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân Mở đầu. Các bài toán thường gặp có thể 2 loại: Bài toán Côsi đối với phương trình vi phân cấp 1: * Bài toán Côsi : là bài toán dạng phương trình vi phân với - Cho khoảng [x0, X] điều kiện bổ sung (điều kiện ban đầu) đã cho tại không quá - Tìm hàm số y = y(x) xác định trên [x0, X] thoả mãn: một điểm. y’ = f(x,y); ( 1 ) Ví dụ: Cho phương trình vi phân cấp 1: y’ = 2x + 1; (a) y(x0) = η ; ( 2 ) - Nghiệm tổng quát : y = x2 + x + C; (b) Trong đó f(x, y) – hàm đã biết; η - số thực cho trước C - hằng số tích phân, phụ thuộc điều kiện ban đầu ( 2 ) - điều kiện Côsi hay điều kiện ban đầu. - Mỗi giá trị của C 1 nghiệm xác định. - Xác định C cần biết thêm 1 điều kiện ban đầu, ví dụ y(x=1) = 2; (c) (b) C = 0; Nghiệm của (a) là y = x2 + x thoả mãn (a) và (c). Bài toán tìm hàm số y(x) thoả mãn p/t vi phân (a) và điều kiện ban đầu (c) bài toán Côsi. 313 314 Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân * Bài toán biên. Giải bài toán Côsi. Bài toán giải phương trình vi phân với điều kiện bổ sung được Phương pháp chuỗi Taylo. y’ = f(x, y); cho tại nhiều hơn 1 điểm. y(x0) = η ; - Cho khoảng [a, b]; Khai triển nghiệm y(x) tại x = x0: - Tìm hàm y = y(x) trên [a, b] thoả mãn: ( 5 ) y’ + p(x)y’ +q(x,y) = f(x); ( 3 ) với điều kiện y(a) = α; y(b) = β ( 4 ) ( 6 ) Trong nhiều trường hợp giải gần đúng . ( 7 ) (3) 315 Tương tự y’” y (x0) chuỗi ( 5 ). 316 79
  80. Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân Đã CM được rằng: 3.1 Giải phương trình vi phân bằng phương pháp Euler Với đủ bé, chuỗi ( 5 ) nghiệm của ( 1 ), ( 2 ) - Là phương pháp số; tổng Sn(x) của n số hạng đầu của ( 5 ) nghiệm xấp xỉ của ( 1 ) , ( 2 ); n càng lớn độ chính xác càng cao. - Xác định từng giá trị của y(x) theo giá trị cụ thể của x bảng các giá trị x và y(x) tương ứng. Nội dung: - Chia [x0, X] n đoạn bằng các nút xi cách đều. ( 5 ) xi = x0 + ih; i = 0, 1, 2, . . ., n; xi Lưới sai phân trên [x0, X] x – nút của lưới; h - bước của lưới: h = const; y’ = f(x,y); i ( 1 ) - y(x) nghiệm đúng của (1), (2) y(x ) = η ; ( 2 ) 0 y’ = f(x,y); ( 1 ) y(x0) = η ; ( 2 ) 317 318 Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.1 Giải phương trình vi phân bằng phương pháp Euler 3.1 Giải phương trình vi phân bằng phương pháp Euler Thành lập công thức tính: y’ = f(x,y); ( 1 ) - y(xi) – giá trị đúng của y(x) tại xi; y(x0) = η ; ( 2 ) - yi – giá trị gần đúng tính được của y(xi); - Điều kiện ban đầu y = η - Giả sử đã biết yi, cần tính yi+1 tại xi+1. 0 - Khai triển Taylor tại xi; h đủ nhỏ bỏ qua các số hạng cuối. ( 8 ) . . . . . . . . . . . . . . . . . (8 a) Nhận xét: - Đơn giản, không phải giải p/trình nào, thuận tiện lập trình giải trên máy tính - biết yi ( 9 ) 319 - Độ chính xác không cao. 320 80
  81. Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.1 Giải phương trình vi phân bằng phương pháp Euler 3.1 Giải phương trình vi phân bằng phương pháp Euler Các bước tính: Cho y’ = f(x, y); y(x0) = η ; - Đánh giá sai số: Sau khi tính được u tại xi với bước h: u(xi,h) tính u(xi, h/2) nghiệm sai số : - Ấn định số khoảng chia n; - Tính h = (x – x0)/n ; (10) - Tính xi = x0 + ih; - Đặt y0 = η - Tính yi+1 = yi + h.f(xi,yi) với i = 0, 1, 2, . . ., n ; - Đặt y(xi, h) = yi; thay h = h/2 tính lại; - Đặt y(xi,h/2) = yi; y(xi, h/2) ~ y(xi) - Sai số: 321 322 Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.1 Giải phương trình vi phân bằng phương pháp Euler 3.1 Giải phương trình vi phân bằng phương pháp Euler Chương trình Chương trình Program Euler1; Program Euler1; Uses crt; Var h:=(xc-x0)/n; x0,x,y0,y,h:real; x:=x0; n,i,j,k:interger; y:=y0; Function F(x,y:real):real; For i:=1 to n do Begin Begin F:=y-(2*x/y); x:=x+i*h; End; y:=y+h*(F(x,y); {Chương trình chính} yf[i]:=y; BEGIN End; clrscr; {In kết quả} {Nhập x0,xc,y0,n} 323 324 81
  82. Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.1 Giải phương trình vi phân bằng phương pháp Euler 3.2 Giải phương trình vi phân bằng phương pháp hình thang Chương trình Program Euler1; Phương pháp Ơle có độ chính xác không cao. Thay (9) ( 9 ) For i:=1 to n do writeln (‘y[‘,i,’] = ‘,yf[i]:8:4); readln; (11) END. u*i+1 trong (11) được tính theo công thức Ơle: Đã chứng minh được phương pháp này có độ chính xác cấp 2: 325 326 Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.2 Giải phương trình vi phân bằng phương pháp hình 3.2 Giải phương trình vi phân bằng phương pháp thang Runge-Kutta (11) Bậc 2: (14) Để tính nghiệm (11) chính xác hơn khi đã biết ui, có thể dùng phép lặp đơn: trong đó (12) (15) Quá trình lặp dừng lại ở bước k khi ε – sai số cho phép. (13) 327 328 82
  83. Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.2 Giải phương trình vi phân bằng phương pháp 3.2 Giải phương trình vi phân bằng phương pháp Runge-Kutta Runge-Kutta Bậc 3: Bậc 4: (16) (18) trong đó trong đó (17) (19) 329 330 Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.2 Giải phương trình vi phân bằng phương pháp 3.2 Giải phương trình vi phân bằng phương pháp Runge-Kutta Runge-Kutta Chương trình Chương trình Program RunKut4; Program RunKut4; Uses crt; Var h:=(xc-x0)/n; x0,x,y0,y,h,k1,k2,k3,k4:real; x:=x0; n,i,j,k:interger; y:=y0; Function F(x,y:real):real; For i:=1 to n do Begin Begin F:=y-(2*x/y); k1:=h*F(x,y); End; x:=x0+h/2; {Chương trình chính} y:=y0+k1/2; BEGIN k2:=h*F(x,y); clrscr; y:=y0+k2/2; {Nhập x0,xc,y0,n} k3:=h*F(x,y); 331 332 83
  84. Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.2 Giải phương trình vi phân bằng phương pháp 3.3 Giải hệ phương trình vi phân bằng phương pháp Runge-Kutta Euler Chương trình Program RunKut4; Cho khoảng [x0,X], tìm các hàm một biến y1(x), y2(x), , yn(x) xác định trên khoảng [x0,X] và thoả mãn hệ phương trình: x:=x0+h; y:=y0+k3; k4:=h*F(x,y); y:=y0+(k1+2*k2+2*k3+k4)/6; ( 20 ) yf[i]:=y; . . . . . . . . . . . . . . . . . x0:=x; y0:=y; End; với ( 21 ) {In kết quả} For i:=1 to n do trong đó f1, f2, , fn – các hàm đã biết. writeln (‘y[‘,i,’] = ‘,yf:8:4); readln; 333 334 END. Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.3 Giải hệ phương trình vi phân bằng phương pháp 3.3 Giải hệ phương trình vi phân bằng phương pháp Euler Euler Công thức Euler: Chương trình Program Euler2; Chia khoảng [x0, X] thành n đoạn con đều nhau bởi các nút uses crt Type ( 22 ) Var Procedure HAM(x:real;y:mX;var F:real); Begin F[1]:= F[n]:= ( 23 ) End; {Chương trình chính} BEGIN clrscr; 335 336 84
  85. Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.3 Giải hệ phương trình vi phân bằng phương pháp 3.3 Giải hệ phương trình vi phân bằng phương pháp Euler Euler Chương trình Chương trình Program Euler2; Program Euler2; {Nhập nF,n,x0,xc} Begin writeln (‘Nhập các giá trị y0[i]’); HAM(x,y,F); For i:=1 to nF do y[j]:=y[j]+h*F[j]; readln(y0[i]); yf[i,j]:=y[j]; h:=(xc-x0)/h; End; For i:=1 to nF do {In kết quả} y[i]:=y0[i]; For i:=1 to n do x:=x0; For j:=1 to nF do For i:=1 to n do Begin Begin write (‘y[‘,i,’,’,j,’] = ‘,yf[i,j]:8:4);writeln; x:=x+i*h; End; For j:=1 to nF do readln; 337 338 END. Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.3 Giải hệ phương trình vi phân bằng phương pháp 3.4 Giải hệ phương trình vi phân bằng phương pháp Euler Runge-Kutta Bài tập: Program RunKut2; Giải hệ phương trình vi phân thường cấp 1 Uses crt; Type mX = array[1 50] of real; Var k1,k2,k3,k4,F,y,y0:mX; h,x0,xc:real; nF,n,i,j,k:integer; Procedure HAM(x:real;y:mX;var F:mX); Begin F[1]:= ; F[nF]:= ; End; 339 340 85
  86. Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.4 Giải hệ phương trình vi phân bằng phương pháp 3.4 Giải hệ phương trình vi phân bằng phương pháp Runge-Kutta Runge-Kutta Program RunKut2; Program RunKut2; {Chương trình chính} {Chương trình chính} BEGIN BEGIN clrscr; {Nhập nF,n,x0,xc} For i:=1 to n do {Nhập yo[i]} Begin For i:=1 to nF do HAM(x,y,F); readln(y0[i]); For j:=1 to nF do h:=(xc-x0)/n; k1[j]:=h*F[j]; x:=x0; x:=x0+h/2; For i:=1 to nF do For j:=1 to nF do y[i]:=y0[i]; y[j]:=y0[j]+k1[j]/2; HAM(x,y,F); 341 342 Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.4 Giải hệ phương trình vi phân bằng phương pháp 3.4 Giải hệ phương trình vi phân bằng phương pháp Runge-Kutta Runge-Kutta Program RunKut2; Program RunKut2; For j:=1 to nF do For j:=1 to nF do k2[j]:=h*F[j]; k4[j]:=h*F[j]; x:=x0+h/2; For j:=1 to nF do For j:=1 to nF do y[j]:= y0[j]+(k1[j]+2*k2[j]+2*k3[j]+k4[j])/6; y[j]:=y0[j]+k2[j]/2; For j:=1 to nF do HAM(x,y,F); yf[i,j]:=y[j]; For j:=1 to nF do x0:=x; k3[j]:=h*F[j]; For j:=1 to nF do x:=x0+h; y0[j]:=y[j]; For j:=1 to nF do End; y[j]:=y0[j]+k3[j]; End; HAM(x,y,F); {In kết quả} 343 344 86
  87. Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.4 Giải hệ phương trình vi phân bằng phương pháp 3.5 Ứng dụng Runge-Kutta Khảo sát động học của một hệ phản ứng phức tạp Program RunKut2; For i:=1 to n do For j:=1 to nF do writeln(yf[i,j]); readln; Biết: END. Tại t = 0: CA = CA,0, CB = CB,0, CC = CC,0, CD = CD,0, CE = CE,0. Các hằng số tốc độ phản ứng: 345 346 Chương 3 Chương 3 Phương trình và hệ phương trình vi phân Phương trình và hệ phương trình vi phân 3.5 Ứng dụng 3.5 Ứng dụng Khảo sát động học của một hệ phản ứng phức tạp Khảo sát động học của một hệ phản ứng phức tạp Biết: Tại t = 0: CA = CA,0, CB = CB,0, CC = CC,0, CD = CD,0, CE = CE,0, CF = CF,0. Các hằng số tốc độ phản ứng: 347 348 87
  88. Chương 4 Chương 4 Tối ưu hóa (Optimization) Tối ưu hóa (Optimization) Đặt vấn đề Đặt vấn đề Thanh Nguyen (2008), Energy & Fuel, ACS Chỉ số hồi lưu trong chưng luyện liên tục khi tổng Hiệu chi phí và chi phí suất xử đầu tư nhỏ nhất lý khí thải NO trong khí lò đốt Chemical Engineering Vol2 349 350 Chương 4 Chương 4 Tối ưu hóa (Optimization) Tối ưu hóa (Optimization) Đặt vấn đề Thanh Nguyen (2010), Fuel, Elsevier 4.1 Tìm cực trị hàm một biến số: phương pháp điểm vàng, phương pháp gradient Độ chuyển Cực trị toàn cục, và hóa các bon cực trị địa phương trong quá trình khí hóa than 351 352 88
  89. Chương 4 Chương 4 Tối ưu hóa (Optimization) Tối ưu hóa (Optimization) 4.1 Tìm cực trị hàm một biến số: phương pháp điểm 4.1 Tìm cực trị hàm một biến số: phương pháp điểm vàng, phương pháp gradient vàng, phương pháp gradient Phương pháp điểm vàng (Golden-section search or Golden ratio) Miền tìm cực trị (Khoảng phân ly nghiệm) 353 354 Chương 4 Chương 4 Tối ưu hóa (Optimization) Tối ưu hóa (Optimization) 4.1 Tìm cực trị hàm một biến số: phương pháp điểm 4.1 Tìm cực trị hàm một biến số: phương pháp điểm vàng, phương pháp gradient vàng, phương pháp gradient Phương pháp điểm vàng (Golden-section search or Golden ratio) Phương pháp điểm vàng (Golden-section search or Golden ratio) 355 356 89
  90. Chương 4 Chương 4 Tối ưu hóa (Optimization) Tối ưu hóa (Optimization) 4.1 Tìm cực trị hàm một biến số: phương pháp điểm 4.1 Tìm cực trị hàm một biến số: phương pháp điểm vàng, phương pháp gradient vàng, phương pháp gradient Phương pháp điểm vàng (Golden-section search or Golden ratio) Phương pháp điểm vàng (Golden-section search or Golden ratio) 357 358 Chương 4 Chương 4 Tối ưu hóa (Optimization) Tối ưu hóa (Optimization) 4.1 Tìm cực trị hàm một biến số: phương pháp điểm 4.1 Tìm cực trị hàm một biến số: phương pháp điểm vàng, phương pháp gradient vàng, phương pháp gradient Phương pháp điểm vàng (Golden-section search or Golden ratio) Phương pháp điểm vàng (Golden-section search or Golden ratio) 359 360 90
  91. Chương 4 Chương 4 Tối ưu hóa (Optimization) Tối ưu hóa (Optimization) 4.1 Tìm cực trị hàm một biến số: phương pháp điểm 4.1 Tìm cực trị hàm một biến số: phương pháp điểm vàng, phương pháp gradient vàng, phương pháp gradient Phương pháp điểm vàng (Golden-section search or Golden ratio) Phương pháp điểm vàng (Golden-section search or Golden ratio) G:= (sqrt(5)-1)/2; d:= G*(xC-x0); xL:=xC-d;yL:=f(xL); xR:=x0+d;yR:=f(xR); Repeat If yL >= yR then Begin x0:=xL; xL:=xR; yL:=yR; d:=G*(xC-x0); xR:=x0+d; yR:=F(xR); 361 End; 362 Chương 4 Chương 4 Tối ưu hóa (Optimization) Tối ưu hóa (Optimization) 4.1 Tìm cực trị hàm một biến số: phương pháp điểm 4.1 Tìm cực trị hàm một biến số: phương pháp điểm vàng, phương pháp gradient vàng, phương pháp gradient Phương pháp điểm vàng (Golden-section search or Golden ratio) Phương pháp điểm vàng (Golden-section search or Golden ratio) If yL <= yR then Program Gold1; Begin Uses crt; xC:=xR; Var xR:=xL; G,d,x0,xC,xL,xR,y0,yL,yR,eps:real; yR:=yL; i,j,k:integer; d:=G*(xC-x0); Function F(x:real):real; xL:=xC-d; Begin yL:= F(xL); F:=x*x; End; End; Until abs(xC-x0)<=eps; {Chuong trinh chinh} x0:=(x0+xC)/2; BEGIN clrscr; 363 364 91
  92. Chương 4 Chương 4 Tối ưu hóa (Optimization) Tối ưu hóa (Optimization) 4.1 Tìm cực trị hàm một biến số: phương pháp điểm 4.1 Tìm cực trị hàm một biến số: phương pháp điểm vàng, phương pháp gradient vàng, phương pháp gradient Phương pháp điểm vàng (Golden-section search or Golden ratio) Phương pháp điểm vàng (Golden-section search or Golden ratio) Program Gold1; Program Gold1; {Nhập x0,xC,eps, } Repeat G:=(sqrt(5)-1)/2; If yL >= yR then d:=G*(xC-x0); Begin xL:=xC-d; x0:=xL; yL:=F(xL); xL:=xR; xR:=x0+d; yL:=yR; yR:=F(xR); d:=G*(xC-x0); xR:=x0+d; yR:=F(xR); End; 365 366 Chương 4 Tối ưu hóa (Optimization) 4.1 Tìm cực trị hàm một biến số: phương pháp điểm vàng, phương pháp gradient Phương pháp điểm vàng (Golden-section search or Golden ratio) Program Gold1; If yL <= yR then Begin xC:=xR; xR:=xL; d:=G*(xC-x0); xL:=xC-d; yL:= F(xL); End; Until abs(xC-x0)<=eps; x0:=(x0+xC)/2; y0:=F(x0); END. 367 92