Tài liệu Giải số mô hình động học rừng bằng phương pháp adi

http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN MỤC LỤC M Lời mở đầu Chương 1. Giới thiệu mô hỡnh và một số kết quả lý thuyết đạt được 1.1 Giới thiệu .CO 1.2 Nghiệm địa phương 1.3 Nghiệm không âm 1.4 Nghiệm toàn cục 1.5 Hệ động lực OKS 1.6 Hàm Lyapunov 1.7 Sự ổn định của nghiệm dừng thuần nhất Chương 2. Phương pháp ADI và một số thuật toán liên quan 2.1 Giới thiệu về sai phân và phương pháp ADI 2.1.1 Giới thiệu 2.1.2 Phân loại phương trình tuyến tính cấp hai OBO 2.1.3 Các bước chính của phương pháp sai phân 2.1.4 Phương pháp ADI 2.2 Một số thuật toán liên quan 2.2.1 Thuật toán Thomas 2.2.2 Thuật toán Newton Chương 3. Giải số mô hình KIL 3.1 Giải số mô hình 3.2 Một số kết quả tính toán được Kết Luận Tài liệu tham khảo Phụ Lục http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN Lời mở đầu Vấn đề nghiên cứu sự bảo tồn và phát triển của rừng có ý nghĩa rất quan M trọng trong việc bảo vệ môi trường cũng như phát triển kinh tế, xã hội. Để quan sát sự phát triển của rừng cần thời gian rất dài và chi phí lớn. Hiện nay, việc sử máy tính là phương pháp rất quan trọng. .CO dụng các mô hỡnh toỏn học và nghiờn cứu cỏc mụ hỡnh đó với sự trợ giúp của Vào năm 1994, Kuznetsov cùng các tác giả khác (xem [6]) đó đưa ra mô hỡnh Cấu Trỳc Tuổi để mô tả sự phát triển của hệ động học rừng. Họ nghiên cứu mô hỡnh rừng đơn loài và chỉ cú hai lớp, lớp cõy non và lớp cõy già. Quỏ OKS trỡnh tỏi tạo của rừng được mô tả thông qua sự tạo hạt của cây già và sự nảy mầm phát triển thành cây non của hạt. Mô hỡnh Cấu Trỳc Tuổi được biểu diễn bởi hệ phương trỡnh sau: OBO  ∂u  ∂t = βδw − γ ( v )u − fu,   ∂v  ∂t = fu − hv,   ∂w = αv − βw + d ∂ 2 w .  ∂t ∂x 2 (0.1) Trong đó ẩn hàm u (t , x ) và v (t , x ) lần lượt là mật độ cây non và cây già ở tại vị KIL trí x ∈Ω và tại thời gian t ∈ (0, ∞) . Ẩn hàm w(t , x ) là mật độ hạt trong không khí tại x ∈Ω và t ∈ (0, ∞) . Phương trỡnh thứ ba của (0.1) mụ tả sự biến đổi hạt, d > 0 là hằng số khuếch tỏn của hạt; α > 0 và β > 0 lần lượt là tốc độ tạo hạt của cây già và tốc độ lắng đọng của hạt. Các phương trỡnh thứ nhất và thứ hai của (0.1) mụ tả sự phỏt triển của cõy non và cõy già; 0 < δ ≤ 1 là tỷ lệ nảy mầm của hạt, γ (v) > 0 là tỷ lệ chết của cây non, phụ thuộc vào mật độ cây già v ; f > 0 là tốc độ phát triển của cây non và h > 0 là tỷ lệ chết của cõy già. 1 http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN Mụ hỡnh này đó được nhiều tác giả quan tâm nghiên cứu. Trong [3], và [4], L. H. Chuẩn, A. Yagi và T. Tsujikawa đó nghiờn cứu mụ hỡnh trong trường hợp hai chiều. Bằng cách biến đổi hệ phương trỡnh (0.1) về phương trỡnh tiến M húa dạng parabolic và sử dụng cụng cụ nửa nhúm giải tớch, cỏc tỏc giả đó xõy dựng được nghiệm tổng quát cho mô hỡnh, nghiờn cứu sự ổn định và dáng điệu .CO tiệm cận của nghiệm khi thời gian đủ lớn. Dựa trờn cỏc kết quả lý thuyết đó đạt được, luận văn nghiên cứu việc xây thuận toán và viết chương trỡnh để giải số mô hỡnh Cấu Trỳc Tuổi trong trường hợp hai chiều. Hệ (0.1) là một hệ phương trỡnh vi phõn đạo hàm riêng khá phức tạp không có phương pháp giải đúng. Phương pháp giải gần đúng được sử dụng OKS trong luận văn là phương pháp sai phân, kết hợp với phương pháp ADI (Alternating Direction Implicit) để giảm thời gian tớnh toỏn và tăng độ chớnh xỏc. Luận văn gồm cú ba chương: Chương 1. Giới thiệu mụ hỡnh Cấu Trỳc Tuổi mụ tả sự phỏt triển của OBO rừng do Kuznetsov cựng cỏc tỏc giả khỏc đưa ra và trỡnh bày một số kết quả lý thuyết đó được chứng minh. Chương 2. Giới thiệu cỏc khỏi niệm cơ bản của phương phỏp sai phõn, phương phỏp ADI (Alternating Direction Implicit) và phương pháp Thomas ứng dụng cho việc giải gần đúng phương trỡnh đạo hàm riêng. Chương 3. Trỡnh bày quỏ trỡnh sai phõn mụ hỡnh Cấu Trỳc Tuổi, sơ đồ KIL khối mụ tả thuật toỏn để giải nghiệm gần đỳng. Đồng thời trỡnh bày một số kết quả tớnh toỏn thu được. 2 http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN Chương 1 GIỚI THIỆU Mễ HèNH CẤU TRÚC TUỔI VÀ MỘT SỐ KẾT QUẢ Lí 1.1. M THUYẾT Giới thiệu .CO Vào năm 1994, Kuznetsov, Antonovsky, Biktashev và Aponina (xem [6]) đó đưa ra mô hỡnh Cấu Trỳc Tuổi để mô tả sự phát triển của hệ động học rừng (xem (0.1)). Trong [3] và [4], L. H. Chuẩn cùng các tác giả khác đó nghiờn cứu mụ hỡnh đó trong trường hợp hai chiều như sau:  ∂u = βδw − γ (v )u − fu  ∂t  ∂v  = fu − hv  ∂t  ∂w = d ∆w − βw + αv  ∂t  ∂w =0   ∂n u ( x,0) = u0 ( x),  v ( x,0) = v0 ( x ), w( x,0) = w0 ( x ) OKS trên Ω × (0,∞), trên Ω × (0,∞), trên Ω × (0,∞), (1.1) OBO trên ∂Ω × (0,∞), trên Ω. Trong đó Ω là miền lồi hoặc cú biờn C 2 và bị chặn trong R 2 . Cỏc ẩn hàm u(t , x ) và v (t , x ) lần lượt là mật độ cây non và cây già ở tại vị trí x ∈Ω và tại thời gian t ∈ (0, ∞) . Ẩn hàm w(t , x ) là mật độ hạt trong không khí tại x ∈Ω và t ∈ (0, ∞) . Phương trỡnh thứ ba của (1.1) mụ tả sự biến đổi hạt, d > 0 là hằng số KIL khuếch tỏn của hạt; α > 0 và β > 0 lần lượt là tốc độ tạo hạt của cây già và tốc độ lắng đọng của hạt. Các phương trỡnh thứ nhất và thứ hai của (1.1) mụ tả sự phỏt triển của cõy non và cõy già; 0 < δ ≤ 1 là tỷ lệ nảy mầm của hạt, γ ( v ) > 0 là tỷ lệ chết của cây non, phụ thuộc vào mật độ cây già v ; f > 0 là tốc độ phát triển của cây non và h > 0 là tỷ lệ chết của cõy già. Hàm w(t , x ) thỏa món điều kiện Neumann trên biên ∂Ω . Các giá trị ban đầu u0 ( x ) ≥ 0 , v0 ( x ) ≥ 0 và w0 ( x ) ≥ 0 3 http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN được cho trên Ω . Giả sử rằng tỷ lệ chết của cây non được cho bởi một phương trỡnh bậc hai cú dạng γ (v ) = a (v − b ) 2 + c, M trong đó a, b và c là các hằng số dương, tức là ta có thể thấy tỷ lệ chết của cây non phụ thuộc vào cây già và đạt giá trị nhỏ nhất khi mật độ cây già là b . .CO Sau đây ta sẽ trỡnh bày một số kết quả lý thuyết đó đạt được (xem [3] và [4]). 1.2. Nghiệm địa phương Trong mục này, ta đi xây dựng nghiệm địa phương cho (1.1) trên không gian nền và không gian các giá trị ban đầu OKS  u   X =  v  ; u ∈ L∞ (Ω), v ∈ L∞ (Ω), w ∈ L2 (Ω)  w    OBO  u0    K =  v0  ∈ X ; u0 ≥ 0, v0 ≥ 0, w0 ≥ 0    w0    ,     .   Toỏn tử tuyến tớnh A được định nghĩa bởi:  f 0 0   A = 0 h 0  0 0 Λ   với  u   D ( A) =  v  ; u, v ∈ L∞ (Ω), w ∈ H N2 (Ω )  w      ,   KIL trong đó Λ là toán tử đóng liên kết với toán tử Laplace − d ∆ + β trong L2 (Ω) với điều kiện biên Neumann. Miền xác định của Λ là khụng gian D (Λ) = H N2 (Ω) = {u ∈ H 2 (Ω); ∂u = 0 trờn ∂Ω }. ∂n Toỏn tử phi tuyến F được định nghĩa bởi  βδ w − γ (v)u    F (U ) =  fu ,   αv   u    U = v ∈ X .  w   4 http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN Khi đó, hệ phương trỡnh (1.1) được viết lại dưới dạng sau:  dU + AU = F (U ),   dt U (0) = U 0 , (1.2) M 0 < t < ∞. Bằng cỏch kiểm tra một số tớnh chất của toỏn tử A và toỏn tử F , ta chứng .CO minh được kết quả sau: Định lý 1.1. Với mỗi giá trị ban đầu U 0 = (u0 , v0 , w0 ) ∈ K , phương trỡnh (1.2) cú nghiệm duy nhất (u, v, w) trong khụng gian hàm (1.3) OKS u, v ∈ C ([0, T1 ]; L∞ (Ω)) ∩ C 1 ((0, T1 ]; L∞ (Ω)),  2 1 2 2  w ∈ C ([0, T1 ]; L (Ω)) ∩ C ((0, T1 ]; L (Ω)) ∩ C ((0, T1 ]; H N (Ω )). Trong đó T1 > 0 được xác định chỉ bởi U 0 . 1.3. Nghiệm khụng õm Định lý 1.2. Cho U 0 ∈ K và giả sử U = (u, v, w) là nghiệm địa phương của (1.2) thu được trong Định lý 1.1. Khi đó u(t ) ≥ 0, v (t ) ≥ 0, w(t ) ≥ 0 trờn Ω với mọi Chứng minh:  ∂u% = βδw% − γ (v% )u% − fu% trên Ω × (0,∞),  ∂t  ∂v% trên Ω × (0,∞),  = fu% − hv% ∂ t   ∂w% trên Ω × (0,∞),  ∂t = d ∆w% − βw% + αχ (Re v% )   ∂w% = 0 trên ∂Ω × (0,∞),  ∂n u% ( x,0) = u ( x), v% ( x,0) = v ( x ), w% ( x,0) = w ( x). 0 0 0   KIL Xột bài toỏn phụ: OBO 0 ≤ t ≤ T0 , trong đó 0 < T0 ≤ T1 được xác định bởi U 0 . (1.4) 5 http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN Trong đó χ (v% ) là hàm cho bởi v% 0 χ (v% ) =  khi v% ≥ 0, khi v% < 0. M Dễ thấy | χ ( v%2 ) − χ ( v%1 ) | ≤ | v%2 − v%1 | ∀v%1 , v%2 ∈ R . Sử dụng Định lý 1.1, ta suy ra bởi U 0 . .CO % trong [0, T2 ] với T2 được xác định (1.4) có nghiệm địa phương duy nhất (u% , v%, w) Ta sẽ chứng minh rằng u%, v%, w% khụng õm. Xột hàm H ( τ ) định nghĩa bởi τ ≥0 khi khi τ <0 OKS 0  H (τ ) = τ 2  2 và đặt Ψ1 (t ) = ∫ H ( w% ( x, t ))dx, 0 ≤ t ≤ T2 . Dễ thấy Ψ1 (t ) là hàm liên tục không âm có Ω đạo hàm % + α ∫ H ' ( w% ) χ (v% )dx ≤ 0. Ψ1' (t ) = − d ∫ H '' ( w% ) | ∇w% |2 dx − β ∫ H ' ( w% ) wdx Ω Ω Ω Tương tự, đặt OBO Vỡ Ψ1 (0) = 0 , ta suy ra Ψ1 (t ) = 0 với mọi t ∈ [0, T2 ] . Do đó w% (t ) ≥ 0 trờn [0, T2 ] . Ψ 2 (t ) = ∫ H (u% ( x, t ))dx, 0 ≤ t ≤ T2 . Ω Ta cũng được Ψ 2 (t ) là hàm liên tục, không âm, có đạo hàm % − ∫ (γ (v% ) + f ) H ' (u% )udx % Ψ '2 (t ) = βδ ∫ H ' (u% )vdx Ω KIL Ω Vỡ w% ≥ 0 , Ψ 2 (0) = 0 nờn u% (t ) ≥ 0 trong [0, T2 ] . Tương tự ta có v% ≥ 0 . Vỡ v% ≥ 0 nờn χ (v% ) = v% , chứng tỏ ( u%, v%, w% ) cũng là nghiệm địa phương của (1.2) trên [0, T0 ] . Theo tớnh duy nhất nghiệm ta suy ra rằng (u , v, w) = (u% , v% , w% ) trờn [0, T0 ] , trong đó T0 = min {T1 , T2 } . Như vậy (1.2) có nghiệm địa phương không âm trong không gian hàm (1.3), trong đó T0 được xác định bởi U 0 . 6 http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN 1.4. Nghiệm toàn cục M Để xây dựng nghiệm toàn cục của bài toán, ta cần đánh giá tiên nghiệm sau: Định lý 1.3. Cho U 0 ∈ K và giả sử rằng U = (u, v, w) là nghiệm địa phương .CO của phương trỡnh (1.2) trờn khoảng [0, TU ] như sau: ∞ 1 ∞ 0 ≤ u, v ∈ C ([0, TU ); L (Ω)) ∩ C ((0, TU ); L (Ω),  2 1 2 2 0 ≤ w ∈ C ([0, TU ); L (Ω)) ∩ C ((0, TU ); L (Ω )) ∩ C ((0, TU ); H N (Ω )). Khi đó ta có đánh giá X ≤ C{e − ρt U 0 X + 1}, 0 ≤ t < TU , OKS U trong đó các hằng số C > 0 và ρ > 0 khụng phụ thuộc U . Sử dụng đánh giá tiên nghiệm trên, ta chứng minh được sự tồn tại và duy nhất nghiệm của hệ phương trỡnh (1.1). Định lý 1.4. Với mọi giá trị ban đầu U 0 ∈ K , hệ phương trỡnh (1.1) cú duy OBO nhất nghiệm toàn cục trong khụng gian hàm 0 ≤ u , v ∈ C ([0, ∞); L∞ (Ω)) ∩ C1 ((0, ∞); L∞ (Ω))  2 1 2 2 0 ≤ w ∈ C ([0, ∞); L (Ω)) ∩ C ((0, ∞); L (Ω)) ∩ C ((0, ∞); H N (Ω)) Chứng minh: Theo Định lý 1.2, bài toán luôn có duy nhất nghiệm địa phương U trong KIL [0, T0 ] . Theo Định lý 1.3, U (T0 ) được xác định bởi U 0 . Do đó nghiệm U cú thể thỏc triển thành nghiệm địa phương trên [0, T0 + τ ] với τ > 0 được xác định bởi U (T0 ) , tức là chỉ phụ thuộc vào U 0 . Tiếp tục quỏ trỡnh thỏc triển đó, ta sẽ thu được nghiệm toàn cục của bài toỏn. 1.5. Hệ động lực 7 http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN Theo Định lý 1.4, với mỗi U 0 ∈ K, tồn tại duy nhất nghiệm toàn cục U (t ,U 0 ) = (u (t ), v (t ), w(t )) của (1.1) và nghiệm này liên tục theo giá trị ban đầu S (t )U 0 = U (t , U 0 ) . M U 0 . Do đó, ta có thể định nghĩa nửa nhóm { S(t) }t≥ 0 trờn K bởi Do đó ta xây dựng được một hệ động lực ( S (t ), K , X ) sinh ra bởi hệ phương 1.6. .CO trỡnh (1.1). Hàm Lyapunov Trong phần này ta sẽ đi xây dựng hàm Lyapunov cho (1.1). Giả sử OKS (u, v, w) là nghiệm toàn cục của (1.1) với giá trị ban đầu U 0 = (u0 , v0 , w0 ) ∈ K . Đặt ϕ (t ) = fu (t ) − hv (t ), 0 ≤ t < ∞ . Từ hai phương trỡnh đầu của (1.1) ta có Nhõn (1.5) với ϕ (t ) = OBO ∂ϕ = f βδ w - {γ (v ) + f + h}ϕ − h{γ (v )v + fv}, 0 < t < ∞. ∂t (1.5) ∂v và tớnh tớch phõn trờn miền Ω ∂t 1 d d ∂v ∂v ϕ 2 dx + h ∫ Γ(v)dx − f βδ ∫ wdx = − ∫ {γ (v) + f + h}( )2 dx , ∫ 2 dt Ω dt Ω ∂t ∂t Ω Ω v trong đó Γ(v) = ∫ (γ (v)v + fv)dv . KIL 0 Tương tự nhân phương trỡnh thứ ba của (1.1) với ∂w và lấy tớch phõn ∂t trờn miền Ω ta được d d β d ∂w ∂w |∇w |2 dx + w2dx − α ∫ v dx = − ∫ ( ) 2 dx . ∫ ∫ 2 dt Ω 2 dt Ω ∂t ∂x Ω Ω Từ hai đẳng thức trên cho ta 8 http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN d α 2 df βδ f β 2δ 2 2 ϕ α [ + | ∇ w | + h Γ ( v ) + w − ( f αβδ )vw] dx dt Ω∫ 2 2 2 = − ∫ [α {γ (v ) + f + h}( Chứng tỏ, hàm số Ω α 2 ( fu − hv) 2 + df βδ f β 2δ 2 | ∇w |2 + hαΓ(v ) + w − ( f αβδ )uv]dx 2 2 .CO Ψ (U ) = ∫ [ là hàm Lyapunov của (1.1). 1.7. 0 < t < ∞. M Ω ∂v 2 ∂w ) + f βδ ( )2 ] dx ≤ 0, ∂t ∂t Sự ổn định của nghiệm dừng thuần nhất OKS Giải sử (u , v , w) là nghiệm dừng, thuần nhất và không âm của (1.1). Khi đó u ≥ 0 , v ≥ 0 và w ≥ 0 thỏa món hệ phương trỡnh  βδw − γ ( v )u − fu = 0,   fu − hv = 0,  − βw + αv = 0.  OBO (1.6) Ta chứng minh được rằng: • Nếu 0 < h ≤ fαδ ab + c + f 2 thỡ (1.6) cú hai nghiệm O = (0,0,0) và fαδ − ( c + f )h P+ =  h (b + D ), b + D , α (b + D )  , trong đó D = , và β ah  f  • Nếu KIL cả hai nghiệm đó đều ổn định. fαδ fαδ 0, cú dạng Parabolic nếu D = 0 và Cỏc dạng phổ biến OBO Hyperbolic nếu D < 0. ● Phương trỡnh Laplace: ∆u = δ 2u δ 2u + = 0. δ x2 δ y 2 ● Phương trỡnh truyền nhiệt: δu δ 2u = a2 2 . δt δx KIL ● Phương trỡnh dõy rung: 2 δu 2 δ u = a . δt2 δ x2 2.1.3. Các bước chính của phương pháp sai phân ● Rời rạc hoỏ miền Ω. ● Thay toỏn tử vi phõn bằng toỏn tử sai phõn. ● Giải hệ phương trỡnh đại số tuyến tính thu được. 11 http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN ● Khảo sát sự hội tụ và ổn định của lược đồ sai phân. .CO M a) Rời rạc hoỏ miền Ω OKS Hỡnh 2.1. Lưới đồng dạng sử dụng cho sai phân hữu hạn. Trong phạm vi luận văn, ta chỉ xét miền Ω là miền hỡnh chữ nhật [0,a] × [0, b] , như vậy ta có thể rời rạc hoá miền Ω bằng các đường thẳng, ta chia đoạn [0,a] và đoạn [0,b] thành J và K đoạn bởi các điểm chia x j = j∆x, j = 0...J , và yk = k ∆y, k = 0...K trong đó ∆x = a b , ∆y = . J K OBO Giao điểm của những đường thẳng đó gọi là những điểm lưới, điểm lưới (j,k) có toạ độ là ( x j , xk ) . Điểm kề của (j, k) là các điểm ( j ± 1, k ) và ( j , k ± 1) . Ký hiệu Ω h = { ( j , k ) : 0 ≤ j ≤ J , 0 ≤ k ≤ K } . Tập hợp các điểm trong Ω h = { ( j, k ) : 0 < j < J , 0 < k < K } , KIL và Γ h = Ω h \ Ω h được gọi là các điểm biên. b) Thay toỏn tử vi phõn bằng toỏn tử sai phõn Ta xột biến thời gian t trờn miền T = (0, +∞). Rời rạc hoỏ miền T bằng ∆t, như vậy ta có: t0 = 0, t1 = ∆t , t2 = 2∆t …và ký hiệu u nj,k = u ( j∆x, k ∆y, n∆t ) . Trong đó u 0j ,k = u0 ( j ∆x, k ∆y ) là giá trị ban đầu tại t0 . Ta có thể biểu diễn các toán tử vi phân bằng sai phân như sau: Sai phân tiến theo hướng x 12 http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN ∂u nj ,k = ∂x u nj +1,k − u nj ,k , ∆x (2.2) ∂u nj ,k u nj ,k − u nj −1,k = ∂x M sai phân lùi theo hướng x . ∆x (2.3) ∂u nj ,k = ∂x u nj +1,k − u nj −1,k 2 ∆x Tương tự sai phân theo hướng y ta cú ∂y = u nj ,k +1 − u nj ,k ∆y và sai phõn trung tõm theo y là ∂u nj ,k ∂y = = . u nj ,k − u nj ,k −1 OKS ∂u nj ,k .CO Từ (2.2) và (2.3) ta có sai phân trung tâm theo hướng x u nj ,k +1 − u nj ,k −1 2∆y ∆y , . OBO Ta sử dụng sai phân trung tâm đối với đạo hàm cấp 2 ∂ 2u nj ,k ∂x 2 ∂ 2u nj ,k ∂y 2 = = u nj +1,k − 2u nj ,k + u nj −1,k , ∆x 2 u nj ,k +1 − 2u nj ,k + u nj ,k −1 ∆y 2 . c) Sau khi được biểu diễn dưới dạng sai phân, các phương trỡnh đạo hàm riêng có dạng một hệ phương trỡnh đại số tuyến tính Ax = b, trong đó A thường Thomas. KIL có dạng ma trận chéo, thưa. Để giải hệ này ta sử dụng thuật toán thuật toán Xét phương trỡnh sau L(u ) := a ∂ 2u ∂ 2u ∂u ∂u + b +c +d + gu = f , 2 2 ∂x ∂y ∂x ∂y trong đó a, b, c, d, g, f là cỏc hàm của x, y và ∀( x, y ) ∈ Ω, D := ab > 0. 13 (2.4) http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN Thay đạo hàm trong (2.4) bằng sai phân tương ứng, ta được phương trỡnh sai phõn sau: u j +1,k − 2u j ,k + u j −1,k ∆x 2 + b j ,k + d j ,k u j ,k +1 − 2u j ,k + u j ,k −1 2 ∆y u j ,k +1 − u j ,k −1 2 ∆y + c j ,k u j +1,k − u j −1,k M Lh (u j ,k ) := a j ,k 2 ∆x (2.5) + g j ,k u j ,k = f j , k , ( j , k ) ∈ Ω h . .CO Ta gọi (2.5) là sơ đồ 5 điểm ( j , k ), ( j ± 1, k ), ( j , k ± 1) . Gộp các số hạng đồng dạng trong (2.5) ta được Lh (u j ,k ) := Aj , k u j +1,k + B j ,k u j −1, k + C j , k u j ,k +1 + D j ,k u j ,k −1 − E j ,k u j ,k = f j ,k , trong đó a j .k ∆x 2 + D j ,k c j ,k , a j ,k b d c , C j , k = j ,k2 + j ,k , ∆x 2∆x ∆y 2∆y 2a j ,k 2b j ,k E j ,k = + − g j ,k . ∆x 2 ∆y 2 OKS Aj ,k = B j ,k = 2∆x b j ,k d j ,k = 2− , ∆y 2 ∆y 2 − Giả sử các điều kiện sau được thoả món ●g ≤ 0 OBO ● a,b > 0 ∀( x, y ) ∈ Ω , (điều kiện để L(u ) là toỏn tử eliptic). ∀( x, y ) ∈ Ω , (điều kiện để (2.4) có nghiệm). Để chứng minh hệ phương trỡnh sai phõn  L h ( u j , k ) = f j , k   u j , k | ∂ Ω h = ϕ j , k tương ứng KIL giải được duy nhất, ta cần chứng tỏ hệ phương trỡnh (tuyến tính) thuần nhất  L h ( u j , k ) = 0   u j , k | ∂ Ω h = 0 chỉ có nghiệm tầm thường. Nguyờn lý maximum: Giả sử dóy {v j ,k } ≠ const xác định trên Ω h và thoả món điều kiện 14 http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN Lh (v j ,k ) ≥ 0, ∀( j , k ) ∈ Ω h , hoặc M Lh (v j ,k ) ≤ 0, ∀( j , k ) ∈ Ω h . Khi đó {v j ,k } không đạt max dương (hoặc min âm) trong Ωh . .CO Chứng minh: Giả sử Lh (v j ,k ) ≥ 0 và M := max( v j ,k ) > 0 . Khi đó tỡm được (i0 , j0 ) ∈ Ω h sao cho Ω h vi0 , j0 = M > 0 và một điểm kề ( j, k ) sao cho v j ,k < M (hoặc ta cú thể chứng tỏ rằng {v j ,k } = const ). Ta cú OKS 0 ≤ Lh ( v j0 , k0 ) < M { A j0 , k0 + B j0 , k0 + C j0 , k0 + D j0 , k0 − E j0 , k0 } = Mg j0 , k0 ≤ 0, vỡ g ≤ 0 . Điều mâu thuẫn này chứng tỏ max v j ,k ≤ 0 . □ Ω h Áp dụng nguyờn lý maximum, dễ dàng chứng minh được hệ phương trỡnh thuần nhất chỉ cú nghiệm tầm thường u j ,k ≡ 0, ( j, k ) ∈ Ω h . Từ đó suy ra hệ OBO phương trỡnh sai phõn  L h ( u j , k ) = f j , k   u j , k | ∂ Ω h = ϕ j , k là giải được duy nhất. d) Khảo sát sự ổn định và hội tụ của lược đồ sai phân trỡnh KIL Để đơn giản ta gộp cả điều kiện biên và điều kiện ban đầu đưa vào phương Lu = f . (2.6) Trong đó L là toán tử tuyến tính đưa không gian tuyến tính định chuẩn (∪,|| . ||∪ ) vào không gian tuyến tính định chuẩn ( Fh ,|| . ||U ) . h Song song với bài toỏn liờn tục ta xột bài toỏn “rời rạc “ sau 15 http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN Lh (uh ) = f h . (2.7) trong đó Lh là toán tử tuyến tính đưa không gian tuyến tính định chuẩn M (∪h ,|| . ||∪h ) vào không gian tuyến tính định chuẩn ( Fh ,|| . ||Fh ) . .CO Ta có sơ đồ sau L F U τh πh OKS Lh Fh Uh trong đó π h và τ h là cỏc toỏn tử rời rạc hoỏ thoả món điều kiện tương thích chuẩn sau OBO ∀u ∈ U , || π h u ||→ || u ||U (| h |→ 0) , ∀f ∈ F , || τ h f ||→|| f ||F , (| h |→ 0). Nếu chỉ số h ∈ R K thỡ |h| là một chuẩn nào đó của h. Định nghĩa 1. i) Lược đồ (2.7) xấp xỉ bài toán (2.6) tại nghiệm u* của (2.6) nếu KIL || Lhπ h u * −τ h f ||Fh → 0, (| h |→ 0). ii) Lược đồ (2.7) xấp xỉ bậc k bài toỏn (2.6) tại nghiệm u* nếu || Lhπ hu * −τ h f ||Fh ≤ c1 | h |k . Định nghĩa 2. Lược đồ sai phân (2.7) là ổn định nếu i) (2.7) cú nghiệm duy nhất với mọi vế phải f h . ii) Nghiệm phụ thuộc liên tục vào điều kiện ban đầu và vế phải. 16 http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN Định nghĩa 3. Nghiệm của bài toỏn sai phõn hội tụ tới nghiệm u* của bài toỏn vi phõn nếu M || π hu * −uh ||∪h → 0, (| h |→ 0). hội tụ bậc k, nếu || π hu * −uh ||∪ ≤ C | h |h . h Định lý Lax. Nếu lược đồ (2.7) ổn định và xấp xỉ (bậc k) bài toỏn (2.6) thỡ .CO nghiệm của (2.7) hội tụ (bậc k) tới nghiệm của (2.6). Ta đi xét một phương pháp khảo sát sự ổn định của lược đồ sai phân sau: Phương pháp phổ Newmann Phương pháp này cho điều kiện cần để lược đồ sai phân ổn định theo giá OKS trị ban đầu. Ta gọi phổ của bài toán sai phân là tập hợp cỏc giỏ trị λ = λ (α , h) sao cho nghiệm của phương trỡnh sai phõn thuần nhất tương ứng có dạng u j ,k [λ (α , h)]k u j ,0eiα j Để lược đồ sai phân ổn định cần (j = 0, ±1…). sup | u j ,k |≤ C sup | u j ,0 | (k = 0…K). j OBO j Từ đây suy ra | λ (α , h) |k ≤ c1 ∀k = 0...K . Nếu λ không phụ thuộc vào các đường j, k chia lưới thỡ điều này thoả món khi | λ ( α )| ≤ 1 với ∀ α . Vớ dụ: Xét phương trỡnh KIL ∂u ∂ 2u ∂ 2u  = 2 + 2 (t > 0), ∂y  ∂t ∂x u ( x, y ,0) = ϕ ( x, y ) (| x |,| y |< ∞ ).  Lược đồ sai phân có dạng  u nj ,+k1 − u nj ,k u nj +1,k − 2u nj ,k + u nj −1,k u nj ,k +1 − 2u nj ,k + u nj ,k −1 = − ,  ∆t ∆x 2 ∆y 2  u 0 = ϕ . j ,k  j ,k 17 http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN Từ u nj ,k = λ nu0 ei (α j + β k ) nghiệm của phương trỡnh thuần nhất, ta cú λ = 1 − 4r sin 2 α 2 − 4r sin 2 β 2 , (r = τ ∆x 2 = τ ∆y 2 = const > 0) . 1 1 . Vậy r > thỡ lược đồ sai phân không ổn định. 4 4 .CO r≤ M Từ đây suy ra điều kiện cần để lược đồ sai phân ổn định là −1 ≤ 1 − 8r ≤ λ ≤ 1 hay Nguyờn tắc maximum Xột bài toỏn ∂u ∂ 2u  − 2 = f ( x, t )  ∂t ∂x u ( x, 0) = ϕ ( x)  OKS ( x < ∞;0 ≤ t ≤ T ) . Lược đồ sai phân có dạng OBO  u j , p +1 − u j , p u j +1. p − 2u j , p + u j −1, p − = f j, p ,  τ h2   T  u j ,0 = ϕ j ( j = 0, ±1,...; p = 0,...,   − 1), τ    τ  r := 2 = const. h  Biểu diễn u j , p +1 qua u j , p tớnh ở lớp p (2.8) u j , p +1 = ru j +1, p + (1 − 2 r )u j , p + ru j −1, p + τ f j , p . Thay u j , p = λ p u0 eiα j vào phương trỡnh thuần nhất, ta tỡm được KIL λ = 1 − 4 r sin 2 α 2 ≤1. 1 2 Để λ ≥ −1 với mọi α ta phải cú r ≤ . Như vậy với r > 1 theo Neumann thỡ lược 2 đồ sai phân không ổn định. Nếu r ≤ 1 thỡ từ (2.8) ta suy ra 2 u j , p +1 ≤ [ r + (1 − 2r ) + r ] sup u j , p + τ f h j h = sup u j , p + τ f h j 18 h ≤ ... ≤ sup u j ,0 + T f h h . j http://kilobooks.com THÖ VIEÄN ÑIEÄN TÖÛ TRÖÏC TUYEÁN Như vậy u h h ϕ ≤ h + T fh h , M hay lược đồ sai phân là ổn định. 2.1.4. Phương pháp ADI (Alternating Direction Implicit) .CO Ta xột bài toỏn sau  u = σ ( u xx + u yy )   u ( x , y , 0) = φ ( x , y ) u( x, y , t ) = γ ( x, y , t )  trên Ω × (0,∞ ), trên Ω , trên ∂Ω × (0,∞ ). (2.9) Chia miền thời gian t ∈ (0, +∞) thành cỏc khoảng ∆t , 2 ∆t , 3 ∆t …, miền lưới các bước chia ∆x và ∆y như sau ∆x = OKS hỡnh chữ nhật Ω có kích thước a × b được trỡnh bày trờn. Ta chia miền Ω với a b , ∆y = , u nj,k = u ( j ∆x, k ∆y, n∆t ). J K Như vậy phương trỡnh (2.9) được viết dưới dạng sai phân  δ x2u nj ,k δ y2u nj ,k  =σ  + , 2 ∆y 2   ∆x OBO u nj ,+k1 − u nj ,k ∆t trong đó δ x , δ y là cỏc toỏn tử sai phõn trung tõm δ x2 u nj ,k = u nj −1,k − 2 u nj ,k + u nj +1,k , Như vậy KIL δ y2 u nj ,k = u nj ,k −1 − 2 u nj , k + u nj ,k +1 . u nj ,+k1 = rx (u nj −1, k + u nj +1, k ) + ry (u nj , k −1 + u nj , k +1 ) + (1 − 2 rx − 2 ry )u nj , k , trong đó rx = σ∆t ∆x 2 , ry = σ∆t ∆y 2 , u 0j ,k = φ ( j ∆x, k ∆y ). Ta có thể biểu diễn lược đồ của công thức (2.10) như sau 19 (2.10)