Phương pháp sai phân

Mục Lục Tổng kết phương pháp sai phân Tổng kết các vấn đề lý thuyết Mục tiêu của phương pháp: Phương pháp sai phân là 1 trong 3 phương pháp để giải phương trình vi phân Bao gồm : Giải đúng nhờ phương pháp biến đổi giải tích sử dụng các công cụ vi tích phân và biến đổi tương ứng. Giải gần đúng bằng phương pháp sai phân Giải gần đúng bằng phương pháp phần tử hữu hạn Thực chất của phương phấp sai phân là việc xấp xỉ các giá trị của hàm tại các điểm đặc biệt ( thường chọn là lưới sai phân đều cho đơn giản) : Như vậy để xây dựng được 1 nghiệm sai phân đúng theo nghĩa sai phân cần tính được các giá trị của nghiệm phương trình vi phân tại các đầu mút quan trọng, có tất cả là m*n ( m là kích thước lưới sai phân, t là thời gian trong trường hợp 1 chiều) các giá trị này với 1 số giá trị biên đã biết còn lại ta cần tìm các giá trị còn lại điều đó đưa đến nhận xét là : Thực chất của phương pháp sai phân chính là giải hệ phương trình đại số tuyến tính khi đã biết 1 số giá trị đặc biệt. Chỉ cần lập được hệ đại số tuyến tính để tính các nghiệm chính là đã giải quyết được 70% công việc của bài toán. Việc còn lại chỉ là tìm phương pháp giải phù hợp ( Euler,Luân phương ẩn, Seidel co dãn, truy hồi 3 đường chéo…..) Giải thích thêm : Ở đây sở dĩ ta phải giải hệ vì vịj thì các chỉ số i,j chạy tương ứng vì vậy ta thu được 1 hệ phương trình chứ ko phải là 1 phương trình(biểu thức) Còn phương pháp phần tử hữu hạn là xấp xỉ các hàm trên các miền hữu hạn thành các toán tử tuyến tính để xấp xỉ được như vậy thực chất là giải bài toán xấp xỉ hàm với sai số đã chọn Cách làm: Để xấp xỉ được giá trị của hàm thì ta dung cách thay thế các đạo hàm bằng các sai phân tương ứng của nó và đưa về giải hệ phương trình. Cụ thể uxi=u’(xi)=uẋi= uẋi+1 (Sai phân cấp 1) =(Sai phân cấp 2) Còn rất nhiều cách thay khác có thể tham khảo trong sách đánh này ko quen lâu quá J. Từ các công thức đó ta thay vào các đạo hàm tương ứng thì ta có một hệ phương trình sai phân với các giá trị nút lưới tương ứng. Giải hệ này ta có đáp số.Cụ thể từng bài toán xem xét dưới đây Các dạng phương trình vi phân thường sử dụng phương pháp sai phân để giải Bài toán truyền nhiệt Truyền nhiệt dừng 1 chiều ( tức là biểu thức =0) Lu = (k(x)u’)’-q(x)u(x)=-f(x) Bài toán biên loại 1: u(a)=A, u(b)=B với a,b là 2 biên Như đã nói ở trên ta đưa về hệ phương trình tuyến tính: aiyi-1-(ai+1+ai+h2qi)yi+ai+1yi+1= -h2fi ( Để có được hệ thức này như đã nói ở trên ta thay các đạo hàm bằng sai phân tương ứng và rút gọn). Ở đây các chỉ số i chạy nên ta thu được hệ phương trình với mỗi phương trình với 3 ẩn yi-1, yi, yi+1 với i thay lần lượt từ 1 đến N-1 -C1y1+B1y2=-(h2f1+ A1y0)=-d1 A2y1-C2y2+B2y3=-h2f2=-d2 ………………………. Ai-1yi-1-Ciyi+Biyi+1=-h2fi=-di-1 AN-1yN-2-CN-1yN-1=-(h2fN-1+BN-1yN) =-dN-1 Với : Ai=ai=k(xi-h/2) Bi=k(xi+h/2) Ci= Ai+ Bi+ h2qi Cần chú ý các công thức tính A,B,C để tính cho nhanh J Đưa về dạng ma trận : AY=-d như trong vở ( ở đây dài quá nên đưa công thức nghiệm luôn để thay như sau : ) Ai=ai=k(xi-h/2) Bi=k(xi+h/2) Ci= Ai+ Bi+ h2qi α1=0 β1=A( Điều kiện đầu) α2= β2= αi+1= βi+1= y0,yN đã cho từ điều kiện biên yN-1= yi= αi+1yi+1+ βi+1 Chú ý quan trọng : di=h2fi phải chú ý dấu của fi cho trong đề bài ( bài giữa kì cả lớp sai vì lỗi này nên cần chú ý đặc biệt) Để đỡ nhầm lẫn nên thay ngay các giá trị trung gian có thể tính từ ban đầu luôn như Ai, Bi,Ci, di=h2fi rồi mới đi giải hệ bằng công thức trên Với lưới sai phân nhỏ ( từ khoảng 4 đến 5 điểm cần tính) có thể dung phương pháp khử Gauss-Jordan như vậy sẽ tiết kiệm khối lượng tính ít nhất 2.5 lần . Ví dụ : Lấy lại ví dụ trong vở thầy cho:Dùng công thức trên ta tính được xi 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 Ai 1.05 1.15 1.25 1.35 1.45 1.55 1.65 1.75 1.85 Bi 1.15 1.25 1.35 1.45 1.55 1.65 1.75 1.85 1.95 h2qi 1.10-3 2.10-3 3.10-3 4.10-3 5.10-3 6.10-3 7.10-3 8.10-3 9.10-3 Ci 2.201 2.402 2.603 2.804 3.005 3.206 3.407 3.608 3.809 di 1.004 -0.0523 -0.0587 -0.0651 -0.0717 -0.0782 -0.0848 -0.0913 7.153 αi 0 0.522 0.694 0.778 0.827 0.858 0.879 0.894 0.905 βi 1 0.456 0.262 0.155 0.082 0.026 -0.02 -0.06 -0.096 y9===3.268 y9=3.268 y8= = 2.862 y7= 2.5034 y6= 2.1820 y5= 1.8988 y4= 1.6520 y3= 1.4401 y2= 1.2617 y1= 1.1154 Truyền nhiệt trong không gian 1 chiều ko dừng Với điều kiện ban đầu : u(x,0)=g(x) Và điều kiện biên :u(a,t)=α(t) :u(b,t)=β(t) Thay các đạo hàm bằng sai phân tương ứng ta xây dựng được 3 phương pháp giải tương ứng là : sai phân hiện, sai phân ẩn và 6 điểm đối xứng Crank-Nicolson: Phương pháp sai phân hiện. Thay Chú ý rằng đạo hàm theo x thì chỉ số j ( tức là của t giữ nguyên ) Và đạo hàm theo t thì chỉ số i ( tức là của x giữ nguyên) Cho dễ nhớ Bài toán đưa về: =+ Hay : =[]+ Chú ý: Ở đây do đã biết lớp thứ j-1 ( lớp 0 đầu tiên suy ra từ điều kiện ban đầu u(x,0)=g(x), ngoại trừ u(a,0) và u(b,0 ) xác định theo điều kiện đầu) nên ta có thể tính được lớp thứ j từ lớp j-1 mà không cần giải hệ . Nếu bài toán có nhân them hằng số hoặc 1 hàm số k(x,t) thì tương tự ta thay các hằng, hoặc hàm số theo sai phân tương ứng là được. Ví dụ: Lấy lại ví dụ trong như trong vở f(x,t)=2; u(0,t)=et; u(1,t)=e1+t+1; u(x,0)= ex+x2 Chọn =0.005, h=0.01 =et (Điều kiện biên) =e1+t+1(Điều kiện biên) vi0=ex+x2 (Điều kiện đầu) Sử dụng lược đồ sai phân ta có công thức tổng quát =[]+ Với =0.5( Ở đây có thể làm tổng quát cho bất kì đặt =r rồi thế vào các biểu thức ta có đáp số, tuy nhiên ở đây chọn luôn cho đơn giản và tiết kiệm công gõ J ). =et =e1+t+1 vi0=ex+x2 Lớp thứ nhất j=1,t=0.005 =0.5[-2+]++0.005. =0.5[e0.005 -2(ex+x2)+ ex+x2]+ ex+x2+0.005*2= =0.5[e0.005-2(e0.1+0.12)+ e0.2+0.22] + e0.1+0.12+0.005*2 =0.5(1-2*1.1152+1.261)+1.1152+0.01=1.141 =0.5[-2+]++0.005 =0.5[(ex+x2)-2(ex+x2)+ (ex+x2)]+ (ex+x2)+0.005*2 =0.5[(e0.1+0.12)-2(e0.2+0.22)+ (e0.3+0.32)]+ (e0.2+0.22)+0.005*2 =0.5*0+1.271= 1.2875 ………………………. =0.5[-2+]++0.005*2 =0.5[(ex+x2)-2(ex+x2)+ e1+t+1]+ ex+x2+0.01 =0.5[(e0.8+0.82)-2(e0.9+0.92)+ e1+0.005+1]+ e0.9+0.92+0.01 =3.309 Các bước 2,3,4… tương tự nhau chỉ có bước 1, và 9 là cần chú ý vì điều kiện biên. Ở đây kết quả xấp xỉ “ không được tốt lắm với nghiệm đúng do ta chọn tỉ số tương đối lớn, và các bước làm tròn nhiều. Ta có bảng sau: v11 v21 v31 v41 v51 v61 v71 v81 v91 Lớp thứ 2 j=2,t=0.01 =0.5[-2+]++0.005*2 =0.5[e0.01-2*1.141+1.2875]+ 1.141+0.01 =0.5*0.0155+1.141+0.01 =1.159 Các nút khác tính tương tự Ta có bảng sau v12 v22 v32 v42 v52 v62 v72 v82 v92 Lớp thứ 3 j=3,t=0.015 Tính tương tự như 2 phần trên chỉ thay t=0.015 là được Ta có bảng sau v13 v23 v33 v43 v53 v63 v73 v83 v93 Phương pháp sai phân ẩn Tương tự như phương pháp sai phân hiện phương pháp sai phân ẩn cũng sử dụng việc thay các đạo hàm bằng các sai phân tương ứng nhưng thay vì sử dụng các sai phân ở lớp dưới (j) ta sử dụng luôn sai phân ở lớp trên (j+1) để từ 1 phương trình ta thu được 3 giá trị xấp xỉ ngay lập tức. Tuy nhiên để từ lớp j+1 xấp xỉ được chính nó ta cần đưa đến việc giải hệ phương trình với mỗi phương trình trong hệ gồm 3 ẩn ở lớp j+1 và các hằng số biết trước do đã tính được lớp j Cụ thế như sau: Thay vào phương trình vi phân ban đầu ta được -= =g(xi) Rút gọn phương trình trên ta được hệ phương trình sau : -(1+2)+=--(1) Như đã nói ở trên biểu thức trên xác định 1 hệ phương trình do các chỉ số i,j chạy . Mỗi phương trình của hệ chứa 3 ẩn ở lớp trên ,, và 1 ẩn ở lớp dưới ( ẩn này có thể tính trực tiếp truy hồi thong qua điều kiện biên và biểu thức trên) nên thực tế mỗi phương trình chỉ có 3 ẩn. Ví dụ: Lấy lại ví dụ trên f(x,t)=2; u(0,t)=et; u(1,t)=e1+t+1; u(x,0)= ex+x2 Ở đây chọn h=0.2 cho ít phương trình Chọn = 0.001 thế thì =0.025 Từ hệ 1 viết lại dạng tường minh ta được: Với j=0,i=1 0.025-(1+2*0.025)+ 0.025=-- Với =e0.001= 1.001, =e0.2 +0.22=1.261 Rút gọn ta đưa về được -1.05+0.025=-1.288 Tương tự với j=0,i=2,3,4 Ta thuc được hệ như sau: -1.05 0.025 0 0 -1.288 0,025 -1.05 0.025 0 -1.654 0 0.025 -1.05 0.025 = -2.184 0 0 0.025 -1.05 -2.868 . Để giải hệ này về mặt tổng quát thì dung truy đuổi 3 đường chéo Ai= Bi= Ci= (1+2) di=+(i>1) d1=++ α1= β1= αi+1= βi+1= y0,yN đã cho từ điều kiện biên yN-1= yi= αi+1yi+1+ βi+1 Nếu hệ ít phương trình có thể dùng phương pháp khử Gauss-Jordan thì nhanh hơn. Ở bài toán trên ta có nghiệm là (1.2661 1.6574 2.1857 2.7835) Lớp thứ 2,3 4…. tính tương tự Phương pháp Crank-Nicolson Lý luận tương tự như phương pháp sai phân ẩn cho cách tìm ra nghiệm Ở đây để tiết kiệm thời gian ta đi luôn và hệ tổng quát tức là sau quá trình biến đổi và thay thế sai phân =g(xi) Lý luận và khai triển như trên ta thấy đây cũng là 1 hệ phương trình 3 đường chéo chính nên có thể giải bằng công thức truy đuổi Ai=0.5 Bi=0.5 Ci= (1+) di=+(i>1) d1= α1= β1= αi+1= βi+1= y0,yN đã cho từ điều kiện biên yN-1= yi= αi+1yi+1+ βi+1 Hoặc hệ nhỏ có thể giải Gauss Ví dụ: Lấy lại ví dụ trên f(x,t)=2; u(0,t)=et; u(1,t)=e1+t+1; u(x,0)= ex+x2 Ở đây chọn h=0.2 cho ít phương trình Chọn = 0.001 thế thì =0.025 Ta thiết lập được hệ dựa vào công thức trên và phân tích đã có từ trước: -2.1 0.025 0 0 2.555 0.025 -2.1 0.025 0 = = 3.311 0 0.025 -2.1 0.025 4.372 0 0 0.025 -2.1 5.832 Có nghiệm là : 1.2359 1.6168 2.1345 2.8026 Truyền nhiệt dừng trong không gian 2 chiều Với điều kiện biên u(x,y)=g(x,y) trên biên của miền đang xét D Thay các đạo hàm bằng sai phân tương ứng ta có: =f(xi,yj) (2) vij=g(xi,yj) với (xi,yj) thuộc biên Đặt phương trình (2) được viết lại thành : ++(+)-2(1+)=h2fij Hay = Ở đây phương trình có tới 4 ẩn ,,,Vì vậy việc giải nghiệm đúng là vô cùng phức tạp, do đó chỉ có thể chọn phương pháp giải gần đúng. ( Đây là lí do vì sao ở phương pháp này ta ko dung 3 đường chéo, hay Gauss để giải). Có rất nhiều phương pháp giải gần đúng hệ này xong ở đây ta quan tâm đến phương pháp giải bằng lặp Seidel đã học trong giải tích số. Chọn xấp xỉ đầu bất kì ( thường chọn là bằng 0). Sử dụng phương pháp lặp Seidel ta có công thức tổng quát j=1,2…M-1; i=1,2…N-1 Sau mỗi bước lặp ta tiến hành phép co giãn để tăng tốc độ hội tụ =λ+(1- λ)( Chọn 1<λ<2) Ví dụ: Lấy ví dụ trong sách giáo khoa ( trang 118) Truyền nhiệt không dừng trong không gian 2 chiều (1) (x,y)D, 0<tT u(x,y,0)=g(x,y) (2) (x,y) = ={x=a,x=b,y=0,y=d} u(a,y,t)=ga(y,t) u(b,y,t)=gb(y,t) u(x,c,t)=gc(y,t) u(x,d,t)=gd(y,t) Lược đồ sai phân hiện: Thay các đạo hàm bằng sai phân tương ứng: (1) (2) (3) Ta được : =++ (4) Hay (*) Từ ta có được Và từ đó ta tính được ,……..Theo công thức (*) Và các giá trị sai phân xấp sỉ trên biên được tính theo =ga(yj,tn) =gb(yj,tn) =gc(xi,tn) =gd(xi,tn) Ví Dụ Lược đồ sai phân ẩn . Cách làm tương tự như trên nhưng ở đây ta thay các xấp sỉ sai phân bằng các sai phân của lớp trên n+1 chứ ko phải lớp dưới n-1 nữa và từ đó ta có : (1) =ga(yj,tn+!) =gb(yj,tn+1) (2) =gc(xi,tn+1) =gd(xi,tn+1) =g(xi,yj) (3) Với t=0 ta tính theo (3) = g(xi,yj) t>0 ta thiết lập được các phương trình, và mỗi phương trình gồm 5 ẩn ,,,, (5) Như phân tích với bài toán truyền nhiệt dừng trong không gian 2 chiều ở trên,hệ này mỗi phương trình có tới 5 ẩn, vì vậy để giải hệ một cách có hiệu quả thì ta dung một phương pháp giải gần đúng nghiệm mà thôi. Cũng như trên ở đây ta chọn phương pháp lặp Seidel. Công thức lặp cụ thể như sau: Chọn xấp xỉ đầu = Với n=0 thì đã biết và tính được từ các điều kiện ban đầu. Với n=1,2…. thì được tính truy hồi theo các bước trước. Sau khi chọn xấp sỉ đầu ta có công thức lặp như sau: Việc tính toán lặp và giải hệ đại số này là tương đối phức tạp nên thường phương pháp này ít dùng. Và ta nghiên cứu một số phương pháp cải tiến như: luân phương ẩn, một chiều địa phương, nhảy ô Lược đồ luân phương ẩn. Để đi tới lược đồ luân phương ẩn ta thay các đạo hảm bằng các sai phân tương ứng ở lớp n+1/2 thay vì lớp n+1 và lớp n( n là lớp theo thời gian, tức là nếu từ t=0 ta xét đến lớp t=/2 trước khi đến lớp t=). Cụ thể sau khi thay thế ta có: Bài toán thứ nhất : Cho phép tính từ Và bài toán thứ 2: Cho phép tính từ Các bước giải tuần tự như sau: Bước 1: Từ = gij=g(xi,yj) ta tính được theo hệ sau: Bước 2: Từ ta tính ra theo hệ sau: Bước 3: Từ ta tính ra theo hệ sau: Bước 4: Từ ta tính ra theo hệ sau: Phương pháp một chiều địa phương. Đúng như tên gọi của phương pháp, phương pháp 1 chiều địa phương dựa trên ý tưởng là tách bài toán truyền nhiệt ko dừng 2 chiều về thành 2 bài toán truyền nhiệt không dừng 1 chiều như sau: (1) Tách thành 2 bài toán: ( Theo chiều của x) (2) (Theo chiều của y) (3) Như đã xét ở tiết trước ta lập lược đồ sai phân cho lần lượt 2 bài toán này: Bài toán 2 đối với x (*) Bài toán 3 đối với y: (**) Các bước giải chi tiết: Bước 1: Bắt đầu từ n=0, Ta tính được Từ hệ (*) Bước 2: Từ ta tính được từ hệ (**) Bước 3: Từ ta tính được từ hệ (*) Bước 4: Từ ta tính được từ hệ (**) Giải các hệ trên bằng công thức truy đuổi 3 đường chéo hoặc khủ Gauss ta có kết quả cần tìm Phương pháp nhảy ô Phương pháp này trong chương trình thầy ko dậy chắc ko phải thi II. Phương trình truyền sóng ( dây rung) Ở đây chỉ xét bài toán dây rung 1 chiều (1) (2) (3) (4) Bài toán sai phân: Chú ý quan trọng : Để giải hệ phương trình đại số tuyến tính có rất nhiều các không nhất thiết phải dùng truy đuổi hay Gauss Khi dùng phương pháp truy đuổi cần chú ý công thức truy đuổi đầu tiên với , 2 số hạng này dùng để tính y0 chứ ko phải y1 trong công thức truy đuổi từ phải vì vậy phải hết sức cẩn thận. Có những bài ( ví dụ như ví dụ đầu tiên ở trên ta lại không dùng , mà tính ngay , Về mặt kiểm tra kết quả tính toán thì có thể nhẩm nghiệm đúng của phương trình vi phân ban đầu trước rồi thay số vào, việc nhẩm nghiệm và kiểm tra tính đúng của các xấp sỉ thì ta sử dụng tính chất đơn điệu của hàm số. Nếu nghiệm đúng là đa thức bậc p thì ta có thể nhẩm ra nghiệm chính xác là 1 hàm số, tuy nhiên với các nghiệm vô tỉ thì ko thể L. Có lẽ các phần điều kiện biên loại 3, miền bất kì là sẽ không thi vì công thức tương đối phức tạp. Về lý thuyết nên học các xây dựng lược đồ sai phân, các chứng minh định lý, bổ đề về sự hội tụ.