Bài giảng Mô hình hóa, mô phỏng và tối ưu hóa các quá trình hóa học - Phần 2: Mô phỏng máy tính

Phần 2: mô phỏng máy tính Modeling, simulation and optimization for chemical process Instructor: Hoang Ngoc Ha Email: ha.hoang@hcmut.edu.vn Bộ môn QT&TB CuuDuongThanCong.com https://fb.com/tailieudientucntt Introduction Numerical Analysis Computer Programming SIMULATION CuuDuongThanCong.com https://fb.com/tailieudientucntt Computer simulation „ Some simulation techniques for solving some of the systems of equations ‰ Solution of (nonlinear) algebraic equations ‰ Ordinary differential equations (ODEs) ‰ Partial differential equations (PDEs) „ Numerical methods ‰ Iterative methods ‰ Discrete difference methods ‰ Femlab, Fortran, Ansys Matlab/Simulink CuuDuongThanCong.com https://fb.com/tailieudientucntt Computer simulation „ Computer programming ‰ Assume that you know some computer programming language ‰ We are not interested in generating the most efficient and elegant code but in solving problems (from point of view of engineers) „ Including extensive comment statements „ Use of symbols (the same ones in the equations describing the systems) „ Debugging (for mistakes in coding and/or in logic) „ CuuDuongThanCong.com https://fb.com/tailieudientucntt Computer simulation Example: CuuDuongThanCong.com https://fb.com/tailieudientucntt Computer simulation CuuDuongThanCong.com https://fb.com/tailieudientucntt Computer simulation CuuDuongThanCong.com https://fb.com/tailieudientucntt Computer simulation „ Interval halving (chia đôi khoảng) CuuDuongThanCong.com https://fb.com/tailieudientucntt Computer simulation „ This problem can be formulated under the following form: „ The goal is to find the solution of this nonlinear equations (in ONE VARIABLE) „ Tools (Iterative methods) ‰ Bisection method (phương pháp phân đoạn) ‰ Newton’s (or Newton-Raphson) method f(x) = 0, x ∈ R CuuDuongThanCong.com https://fb.com/tailieudientucntt Iterative method „ Intermediate value theorem ‰ If f is a real-valued continuous function on the interval [a, b], and u is a number between f(a) and f(b), then there is a such that f(c) = uc ∈ [a, b] If f(a) and f(b) are of opposite sign, there exist a number p in [a, b] with f(p)=0 CuuDuongThanCong.com https://fb.com/tailieudientucntt Iterative method „ Bisection method Computer programming: Matlab CuuDuongThanCong.com https://fb.com/tailieudientucntt Iterative method „ Newton’s method Numerical solutions of nonlinear systems of equations (of SEVERAL VARIABLES)Î (See Ref.) CuuDuongThanCong.com https://fb.com/tailieudientucntt Computer simulation „ Interpolation and polynomial approximation ‰ Interpolation and the Lagrange polynomial ‰ Cubic spline interpolation ‰ „ Numerical differentiation and intergration ‰ Numerical differentiation ‰ Richardson’s extrapolation ‰ (See Ref.) CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) y(t) ttN x y(tN )x t0 y(t0) x t1 y(t1) Interpolation CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) „ Tools: ‰ Euler’s method ‰ Higher-Order Taylor methods ‰ Runge-Kutta methods ‰ CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) „ Euler’s method CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) „ Euler’s method CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) „ Example Exact solution? y(t) = −0.5 exp(t) + (t+ 1)2 P/p Euler n=10? y0 = y − t2 + 1, t ∈ [0 2] y(0) = 0.5 Approximate solution? n = 10⇒ h = b−a n = 0.2 Computer programming: Matlab CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) „ Local truncation error The local truncation error in Euler’s method is O(h) Definition CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) „ Higher-Order Taylor methods CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) „ Higher-Order Taylor methods CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) „ Runge-Kutta methods CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) „ Runge-Kutta methods ‰ Xây dựng công thức tính theo mà không phải đạo hàm « tay », cần xấp xỉ mà không dùng đạo hàm với „ Minh họa qua k=2 wi+1 wi T (k) O(hk) f(t, y) T (2)(t, y) = f(t, y) + h 2 f 0(t, y) f 0(t, y) = f 0t(t, y) + f 0 y(t, y)y 0(t) CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) Như vậy Cần tìm với sai số đểa1,α1,β1 a1f(t+ α1, y + β1) ' T (2)(t, y) O(h2) T (2)(t, y) = f(t, y) + h 2 f 0t(t, y) + h 2 f 0y(t, y)f(t, y) f(t+ α1, y + β1) ' f(t, y) + f 0t(t, y)α1 +f 0y(t, y)β1 CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) Cần chọn Đồng nhất hai vế a1f(t, y) + a1α1f 0t(t, y) + a1β1f 0 y(t, y) = f(t, y) + h2 f 0 t(t, y) + h 2 f 0 y(t, y)f(t, y) ⎧ ⎨ ⎩ a1 = 1 a1α1 = h2 a1β1 = h2 f(t, y) ⎧ ⎨ ⎩ a1 = 1 α1 = h2 β1 = h2 f(t, y) CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) Sơ đồ trung điểm (R_K bậc 2) wi+1 = wi + hT (2)(ti, wi) wi+1 = wi + h h f(ti + h 2 , wi + h 2 f(ti, wi)) i ⎧ ⎪⎪⎨ ⎪⎪⎩ w0 = α k1 = h 2 f(ti, wi) k2 = hf(ti + h 2 , wi + k1) wi+1 = wi + k2 Computer programming: Matlab CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Ordinary Differential Equations (ODEs) Sơ đồ R_K bậc 4 ⎧ ⎪⎪⎪⎪⎪⎪⎨ ⎪⎪⎪⎪⎪⎪⎩ w0 = α k1 = hf(ti, wi) k2 = hf(ti + h 2 , wi + k1 2 ) k3 = hf(ti + h 2 , wi + k2 2 ) k4 = hf(ti + h,wi + k3) wi+1 = wi + 1 6 (k1 + 2k2 + 2k3 + k4) Về nhà tự đọc R_K cho hệ và viết chương trình R_K cho hệ CuuDuongThanCong.com https://fb.com/tailieudientucntt Numerical intergration of Partial Differential Equations (PDEs) „ Click here CuuDuongThanCong.com https://fb.com/tailieudientucntt Giới thiệu chung u = u(x . . . , t) một đại lượng vật lý của hệ khảo sát n chemical species Inlet material and/or energetic flux Outlet material and/or energetic flux P k νkSk = 0 dV {C, T . . . } HỆ PHÂN BỐ Profile Local observation CuuDuongThanCong.com https://fb.com/tailieudientucntt Ba dạng phương trình đạo hàm riêng cơ bản Động học biến hệ thống có thể thuộc về các dạng phương trình sau: u = u(x . . . , t) „ Phương trình elliptic (tĩnh-static) ∂2u ∂x2 + ∂2u ∂y2 = f(x, y) „ P/t parabolic (b/toán truyền nhiệt) ∂u ∂t = a 2 ∂2u ∂x2 „ Phương trình hyperbolic (b/toán truyền sóng) ∂2u ∂t2 = a 2 ∂2u ∂x2 CuuDuongThanCong.com https://fb.com/tailieudientucntt Ba dạng phương trình đạo hàm riêng cơ bản „ Phương pháp tìm nghiệm ‰ Phương pháp giải tích ‰ Phương pháp số Ý tưởng: xấp xỉ sai phân các đạo hàm riêng tại các điểm rời rạc (kg,tg) và tính giá trị của tại đóu = u(x . . . , t) ∂u(x,t) ∂t ≈ u(x,t+∆t)−u(x,t) ∆t x • ∆t (x, t) • (x,∆t+ t)t u(x, t) CuuDuongThanCong.com https://fb.com/tailieudientucntt Xấp xỉ sai phân ∂2u(x,t) ∂x2 ≈ u(x+∆x,t)−2u(x,t)+u(x−∆x,t) (∆x)2 (x, t) • •(x+∆x, t)(x−∆x, t) t • ∆x x CuuDuongThanCong.com https://fb.com/tailieudientucntt Xấp xỉ sai phân • •• x y P1 P = (x, y) ∆x P • • ∆y P3 P2 P4 ∂2u(P ) ∂x2 + ∂2u(P ) ∂y2 ≈ u(P1)+u(P2)+u(P3)+u(P4)−4u(P ) h2 h = ∆x = ∆y CuuDuongThanCong.com https://fb.com/tailieudientucntt BÀI TOÁN ELLIPTIC „ Bài toán elliptic với điều kiện biên Dirichlet( ∆u = ∂ 2u ∂x2 + ∂2u ∂y2 = f(x, y), (x, y) ∈ Ω ⊂ R2 u(x, y) = g(x, y), (x, y) ∈ Γ CuuDuongThanCong.com https://fb.com/tailieudientucntt BÀI TOÁN ELLIPTIC „ Phân hoạch và tạo lưới : chia nhỏ bởi các đường thẳng // với Ox và Oy Ω Ω h = ∆x = ∆y Tạo lưới bước chia cách đều Kí hiệu P1,P2, P3 và P4 là 4 điểm rời rạc x/q P ∆u(P ) ≈ u(P1)+u(P2)+u(P3)+u(P4)−4u(P ) h2 Lần lượt thay vào phương trình elliptic, sử dụng công thức xấp xỉ & điều kiên biên Pk = (xk, yk) Hệ PTTT ẩn uk = u(Pk) CuuDuongThanCong.com https://fb.com/tailieudientucntt BÀI TOÁN ELLIPTIC „ Example Lưới 4 nút ẩn đánh số 4 giá trị cần tìm. SD PTSP & Giá trị trên biên Hệ PTTT Mô phỏng Matlab CuuDuongThanCong.com https://fb.com/tailieudientucntt BÀI TOÁN PARABOLIC (mô hình truyền nhiệt) CuuDuongThanCong.com https://fb.com/tailieudientucntt BÀI TOÁN PARABOLIC (mô hình truyền nhiệt) CuuDuongThanCong.com https://fb.com/tailieudientucntt Mức thứ k: t = k∆t Mức thứ 0: đã biếtu(0) Mức thứ 1: chưa biếtu(1) BÀI TOÁN PARABOLIC (mô hình truyền nhiệt) „ Phân hoạch và xấp xỉ Chia cách đều bước và∆x ∆t u(k) = (. . . uki . . . ) Xấp xỉ đạo hàm riêng cấp 1: TIẾN - LÙI∂u(x,t) ∂t CuuDuongThanCong.com https://fb.com/tailieudientucntt BÀI TOÁN PARABOLIC (mô hình truyền nhiệt) „ Ví dụ Sơ đồ hiện (s/p tiến) Sơ đồ ẩn (s/p lùi) CuuDuongThanCong.com https://fb.com/tailieudientucntt BÀI TOÁN PARABOLIC (mô hình truyền nhiệt tổng quát) Miền Ω = {(x, t)|0 ≤ x ≤ 1, t ≥ 0} Sai phân tiến & Sai phân lùi Xấp xỉ ∂u ∂t , ∂u ∂x , ∂2u ∂x2 & ĐK biên, đầu Giá trị u tại các điểm lưới CuuDuongThanCong.com https://fb.com/tailieudientucntt BÀI TOÁN PARABOLIC (mô hình truyền nhiệt) CuuDuongThanCong.com https://fb.com/tailieudientucntt BÀI TOÁN PARABOLIC (mô hình truyền nhiệt) Kí hiệu u(k) = [uk1 . . . u k n], f (k) = [fk1 . . . f k n ] CuuDuongThanCong.com https://fb.com/tailieudientucntt BÀI TOÁN PARABOLIC (mô hình truyền nhiệt) „ Sơ đồ Crank-Nicholson CuuDuongThanCong.com https://fb.com/tailieudientucntt Outline „ General introduction ‰ Structure and operation of chemical engineering systems ‰ What is a chemical process? ‰ Motivation examples „ Part I: Process modeling „ Part II: Computer simulation „ Part III: Optimization of chemical processes CuuDuongThanCong.com https://fb.com/tailieudientucntt