二維波動方程的有限差分法

1、生實驗報實驗課程名稱偏微分方程數(shù)值解數(shù)統(tǒng)學院開課實驗室數(shù)統(tǒng)年級 2013 專業(yè)班信計02班學生姓名開課時間 2015至2016學年第2學期總成績教師簽名數(shù)學與統(tǒng)計學院制開課學院、實驗室： 數(shù)統(tǒng)學院實驗時間：2016年6月20日實驗項目名稱二維波動方程的有限差分法實驗項目類型驗證演示綜合設計其他指導教師曾芳成績是一.實驗目的通過該實驗，要求學生掌握求解二維波動方程的有限差分法， 并能通過計算機語言編程 實現(xiàn)。二.實驗內(nèi)容考慮如下的初值問題:u x, y,0u x, y,t2u,x,y y0,1 2,tx,y,00,1.4sin xsin y,-u0, x,y ,t 0,1.40, x, y20,

2、1(1)在第三部分寫出問題（1）三層顯格式。根據(jù)你寫出的差分格式，編寫有限差分法程序。取 h 0.1,該問題的解析解為u x,y,t示，對數(shù)值結(jié)果進行簡單的討論。4.將所寫程序放到第四部分。0.1h，分別將t 0.5,1.0,1.4時刻的數(shù)值解畫圖顯示。cosJ2 tsin xs in y，將四個時刻的數(shù)值解的誤差畫圖顯三. 實驗原理、方法（算法）、步驟11 4網(wǎng)格劃分 h 0.1,0.1h，故 N - 10,M 一 140，Xiih, yjhjh, i,j 0,1,L ,10 ,tk k ，k0,1丄,140。在內(nèi)網(wǎng)點Xi,yj,tk ,利用二階中心差商，對（1 ）建立差分格式:k 1Ui,

3、jkk 12Ui,jUi,jkUi 1,jc k2Ui,jkUi 1,jkkkUi,j 12Ui,jUi,j 12h2h2整理得到：k 1ui,j2kr ui 1,jkui 1,jkUi,j 1kui,j 1CX 2kk 12 4rui,jU,j其中，i,j 1,2,L,9, k 1,2, L,139，網(wǎng)比r h0.1，局部截斷誤差為考慮邊界條件u x,y,t0, x,y,t 0,1.4，差分格式為：kkU0,0U0,NUN,0kUN ,N0,k0,1,L ,140考慮初始條件U x,y,0sin xsin y，差分格式為：u0. Ui, jsinxi sinyjsinih sinjh ,i,

4、j 0,1,L ,10考慮初始條件Ut X, y,020, x, y0,1 ，利用二階差商近似：1ui,j1ui,j0,i, j0,1,L,1020h2。設k 0時刻的點為內(nèi)點，貝W足差分格式（2），代入上式得到:120000C/201ui,j rui 1,j ui 1,j u,j 1 ui,j 12 4r ui,j ui,j(7)將（6）得到的結(jié)果Sj Ui,1代入（7、中，整理得到:1ui,j1 2 0 0 0 0 .20 2r Ui 1,J Ui 1,J Ui,j 1 Ui,j 11 2r Ui,j綜上（2）、(4 )、（5）、（8）得到三層顯格式的差分格式為:0Ui,jUi,jk 1u

5、i,jsin2r20.1 ,2 kkkkc/2kr Ui 1,j Ui 1,j Ui,j 1Ui,j 12 4r 嘔i,j 1,2,L ,9,k 1,2,L ,139kkkU0,0U0,NUN,0Xi sin yj sinuN,N 0,k0,1,L ,140ih sin jh ,i, j0,1,L0 0 0ui 1,jUi 1,j Ui,j 1局部截斷誤差為00“Ic 20 Ui,j 112rUi,j,i, J2 h2。k 1Ui,j,100,1L ,10四. 實驗環(huán)境(所用軟件、硬件等)及實驗數(shù)據(jù)文件Matlab%二維波動方程數(shù)值計算(關(guān)鍵：怎么運用i,j,k三個指標建立循環(huán))cic;%可以

6、將代碼換成函數(shù)m文件h=0.1;tau=0.1*h;%定義步長空間網(wǎng)格剖分精確解計算r=tau/h;% 網(wǎng)比x,y,t=meshgrid(0:h:1,0:h:1,0:tau:1.4);% uu=cos(sqrt (2)*p i*t).*si n(pi *x).*si n(p );%第一層網(wǎng)點計算u=si n(pi *x).*si n(pi *y);%初始條件u1=u(:.:.1);%因為此時得到的u為11x11x141，故只取第一層%第二層網(wǎng)點計算 for i=2:10for j=2:10u(i.j.2)=0.5*r2*(u(i+1.j.1)+u(i-1.j.1)+u(i.j+1.1)+u(i

7、.j-1.1)+(1-2*r2)*u(i.j.1); u(11.:.2)=0;u(:.11.2)=0;end end u2=u(:,:,2);%第3-141層網(wǎng)點計算 for k=2:140for i=2:10for j=2:10u(i.j.k+1)=rA2*(u(i+1.j.k)+u(i-1.j.k)+u(i.j+1.k)+u(i.j-1.k)+(2-4*rA2)*u(i.j.k)-u(i.j.k-1):u(11,:,k+1)=0;u(:,11,k+1)=0;endendend%結(jié) 果 分圖%wucha=abs(u-uu);% 求絕對誤差矩陣 11x11x141wucha1=wucha(:,

8、:,11);%計算t=0.1時刻的絕對誤差矩陣11x11wucha2=wucha(:,:,51);%計算t=0.5時刻的絕對誤差矩陣11x11wucha3=wucha(:,:,101);%計算t=1.0時刻的絕對誤差矩陣11x11wucha4=wucha(:,:,141);%計算t=1.4時刻的絕對誤差矩陣11x11x0=0:h:1;y0=0:h:1;%析%作t=0.1時刻的絕對誤差圖sub plot(2,2,1);mesh(x0,y0,wucha1);title('t=0.1時刻的絕對誤差');xlabel('x 變量');ylabel('y 變量&#

9、39;);zlabel('%作t=0.5時刻的絕對誤差圖絕對誤差值'）;sub plot(2,2,2);mesh(x0,y0,wucha2);title('t=0.5時刻的絕對誤差');xlabel('x 變量');ylabel('y 變量');zlabel('%作t=1.0時刻的絕對誤差圖絕對誤差值'）;sub plot(2,2,3);mesh(x0,y0,wucha3);title（'t=1.0時刻的絕對誤差'）;xlabel（'x 變量'）;ylabel（'y 變量&#

10、39;）;zlabel（'%作t=1.4時刻的絕對誤差圖絕對誤差值'）;sub plot(2,2,4);mesh(x0,y0,wucha4);title（'t=1.4時刻的絕對誤差'）;xlabel（'x 變量'）;ylabel（'y 變量'）;zlabel（'絕對誤差值'）;% 四個時刻數(shù)值解、精確 解%作t=0.1、0.5時刻的數(shù)值解與精確解 sub plot（2,2,1）;mesh（x0,y0,u（:,:,11）;% 作 t=0.1 時刻的數(shù)值解 title（'t=0.1時刻的數(shù)值解'）;xl

11、abel（'x 變量'）;ylabel（'y 變量'）;zlabel（'u 值'）; sub plot（2,2,2）;mesh（x0,y0,uu（:,:,11）;% 作 t=0.1 時刻的精確解 title（'t=0.1時刻的精確解'）;xlabel（'x 變量'）;ylabel（'y 變量'）;zlabel（'u 值'）;%作t=0.5時刻的數(shù)值解與精確解 sub plot（2,2,3）;mesh（x0,y0,u（:,:,51）;% 作 t=0.5 時刻的數(shù)值解 title（

12、9;t=0.5 時刻的數(shù)值解'）;xlabel（'x 變量'）;ylabel（'y 變量'）;zlabel（'u 值'）;sub plot（2,2,4）;mesh（x0,y0,uu（:,:,51）;% 作 t=0.5 時刻的精確解 title（'t=0.5時刻的精確解'）;xlabel（'x 變量'）;ylabel（'y 變量'）;zlabel（'u 值'）;% 分別復希 y 粘貼運 行%作t=1.0、1.4時刻的數(shù)值解與精確解 sub plot（2,2,1）;mesh（x0

13、,y0,u（:,:,101）;% 作 t=1.0 時刻的數(shù)值解 title（'t=1.0 時刻的數(shù)值解'）;xlabel（'x 變量'）;ylabel（'y 變量'）;zlabel（'u 值'）; sub plot（2,2,2）;mesh（x0,y0,uu（:,:,101）;% 作 t=1.0 時刻的精確解 title（'t=1.0 時刻的精確解）；xlabel（'x 變量'）;ylabel（'y 變量'）;zlabel（'u 值'）;%作t=1.4時刻的數(shù)值解與精確解 su

14、b plot（2,2,3）;mesh（x0,y0,u（:,:,141）;% 作 t=1.4 時刻的數(shù)值解 title（'t=1.4時刻的數(shù)值解'）;xlabel（'x 變量'）;ylabel（'y 變量'）;zlabel（'u 值'）;sub plot（2,2,4）;mesh（x0,y0,uu（:,:,141）;% 作 t=1.4 時刻的精確解 title（'t=1.4時刻的精確解'）;xlabel（'x 變量'）;ylabel（'y 變量'）;zlabel（'u 值'

15、;）;五. 實驗結(jié)果及實例分析1、t 0.1、0.51.01.4時刻的數(shù)值解與精確解圖a丄一v*t+Hfl0107=D4圖 1 t=0.1asat04Q0H-ArTHOJ(Siff、0.5時刻的數(shù)值解、精確解.1iQi b-1 fi."ri 'bi' *區(qū)打穿號<12fl 1 站 i t '.ih HK0204.厶鯉4fl Op：*- 卩諾曲利6J41D4一404DBcert詁圖2 t=1.0、1.4時刻的數(shù)值解、0.17Pp代表維數(shù)，本文P 2注：上兩圖為四個時刻的數(shù)值解與精確解，r -h，三層顯格式達二階收斂，不難看出，收斂效果很好，符合理論。下圖是

16、四個時刻的絕對誤 差圖像，從圖中看出，絕對誤差較小，且經(jīng)過計算得到，收斂階近似于2，正好符合理論值。2、t 0.1、0.51.01.4時刻的絕對誤差圖f1Du_fltil* 出”口.斧時 ilu&i>ieC'DQ-T；«訥.mon鳥- A亠、JS=.40D*O'«Q»0i*¥圖3四個時刻的絕對誤差3、四個時刻（t=0.1、0.5、1.0、1.4 ）的絕對誤差表t=o.i時刻的絕對誤差0.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.0001

17、0.00010.00020.00020.00020.00020.00020.00010.00010.00000.00010.00030.00040.00040.00050.00040.00040.00030.00010.00000.00020.00040.00050.00060.00060.00060.00050.00040.00020.00000.00020.00040.00060.00070.00070.00070.00060.00040.00020.00000.00020.00050.00060.00070.00080.00070.00060.00050.00020.00000.0002

18、0.00040.00060.00070.00070.00070.00060.00040.00020.00000.00020.00040.00050.00060.00060.00060.00050.00040.00020.00000.00010.00030.00040.00040.00050.00040.00040.00030.00010.00000.00010.00010.00020.00020.00020.00020.00020.00010.00010.00000.00000.00000.00000.00000.00000.00000.00000.00000.00000.0000.0000.

19、0000.0000.0000.0000.0000.0000.0000.0000.00000000000000t=0.5時刻的絕對誤差0.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.000000000000000.0000.0000.0010.0010.0020.0020.0020.0010.0010.0000.000073812183700.0000.0010.0020.0030.0040.0040.0040.0030.0020.0010.000035402045300.0000.0010.0030.0040.0050.0050.0050

20、.0040.0030.0010.000084758574800.0000.0020.0040.0050.0060.0060.0060.0050.0040.0020.000010558550100.0000.0020.0040.0050.0060.0070.0060.0050.0040.0020.000022881882200.0000.0020.0040.0050.0060.0060.0060.0050.0040.0020.000010558550100.0000.0010.0030.0040.0050.0050.0050.0040.0030.0010.000084758574800.0000

21、.0010.0020.0030.0040.0040.0040.0030.0020.0010.000035402045300.0000.0000.0010.0010.0020.0020.0020.0010.0010.0000.000073812183700.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.00000000000000t=1.0時刻的絕對誤差0.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.000000000000000.0000.0010.0030.0040.0050.005

22、0.0050.0040.0030.0010.000061313131600.0000.0030.0050.0080.0090.0100.0090.0080.0050.0030.000019261629100.0000.0040.0080.0110.0130.0130.0130.0110.0080.0040.000032329232300.0000.0050.0090.0130.0150.0160.0150.0130.0090.0050.000016264626100.0000.0050.0100.0130.0160.0170.0160.0130.0100.0050.00003194249130

23、0.0000.0050.0090.0130.0150.0160.0150.0130.0090.0050.000016264626100.0000.0040.0080.0110.0130.0130.0130.0110.0080.0040.000032329232300.0000.0030.0050.0080.0090.0100.0090.0080.0050.0030.000019261629100.0000.0010.0030.0040.0050.0050.0050.0040.0030.0010.000061313131600.0000.0000.0000.0000.0000.0000.0000.0000