一維對流方程在A、B、C三種差分格式

1、.一、上機(jī)目的用數(shù)值方法計(jì)算一維對流方程在A、B、C三種差分格式下的解。取為0.05. 取值為0.5，1，2。并作相關(guān)討論。二、實(shí)驗(yàn)原理三、上機(jī)要求：1.學(xué)會對MS-FORTRAN的基本操作。 2.用Fortran編寫程序計(jì)算一維對流方程在A、B、C三種格式下的解。3.討論各種格式下不同的值的差分格式解的特點(diǎn)。四、實(shí)驗(yàn)程序以A格式為例，對微分方程進(jìn)行離散化， 得出 A 格式的差分解的表達(dá)式：B、C格式同理可以寫出。由此編寫如下的Fortran程序。注：除了循環(huán)時(shí)間層的計(jì)算公式略有不同外，三個(gè)程序沒有區(qū)別，因此這里只用一個(gè)主程序，并根據(jù)格式選擇的不同，采用條件語句判斷執(zhí)行哪一部分循環(huán)體。!空間節(jié)

2、點(diǎn)321，dx=0.05 輸出依次為（時(shí)間，空間，數(shù)值）program mainimplicit nonereal dx_dt !定義x/t的值integer abc,r_t,i,j,k !定義變量，abc為格式類型，r_t為時(shí)間網(wǎng)格數(shù)，其余為循環(huán)變量real,allocatable:s(:,:) !定義存儲矩陣swrite(*,*) 輸入dx_dt=0.5,1,2read(*,*) dx_dtwrite(*,*) 選擇格式，A,B,C分別輸入1,2或3 read(*,*) abc !根據(jù)格式選擇生成相應(yīng)的文件 if(abc=1) then open(unit=8,file=out_a.csv

3、) elseif(abc=2) then open(unit=8,file=out_b.csv) elseif(abc=3) then open(unit=8,file=out_c.csv) endifr_t=160/dx_dt !計(jì)算時(shí)間網(wǎng)格總數(shù)allocate(s(r_t+1,321) !分配存儲矩陣的空間!第一層賦初值do i=1,140,1s(1,i)=0write(8,*)1,i,s(1,i)end dodo i=141,161,1s(1,i)=1+0.05*(i-161)write(8,*)1,i,s(1,i)end do do i=162,181,1s(1,i)=1-0.05*(

4、i-161)write(8,*)1,i,s(1,i)end do do i=182,321,1s(1,i)=0write(8,*)1,i,s(1,i)end do !循環(huán)時(shí)間層,根據(jù)格式的選擇來判斷執(zhí)行哪一部分if(abc=1) thendo i=2,r_t,1do j=i,322-i,1s(i,j)=s(i-1,j)-(s(i-1,j+1)-s(i-1,j-1)/(dx_dt*2) write(8,*)i,j,s(i,j) end do do k=1,i-1,1 !余下部分賦值0,下同 s(i,k)=0 write(8,*)i,k,s(i,k) end do do k=322-i,321,1

5、 s(i,k)=0 write(8,*)i,k,s(i,k) end do end do elseif(abc=2) then do i=2,r_t+1,1 do j=1,322-i,1 s(i,j)=s(i-1,j)-(s(i-1,j+1)-s(i-1,j)/dx_dt write(8,*)i,j,s(i,j) end do do k=322-i,321,1 s(i,k)=0 write(8,*)i,k,s(i,k) end do end do elseif(abc=3) then do i=2,r_t+1,1 do j=i,321,1 s(i,j)=s(i-1,j)-(s(i-1,j)-s

6、(i-1,j-1)/dx_dt write(8,*)i,j,s(i,j) end do do k=1,i-1,1 s(i,k)=0 write(8,*)i,k,s(i,k) end do end doendif!完成提示write(*,*)數(shù)據(jù)已輸出至源目錄pausestopend program五、實(shí)驗(yàn)結(jié)果及分析程序運(yùn)行后在對應(yīng)目錄下生成csv表格文件，根據(jù)輸入的xt的值不同生成對應(yīng)的網(wǎng)格并計(jì)算各節(jié)點(diǎn)數(shù)值。這里以xt=1為例，就A、B、C三種格式的結(jié)果進(jìn)行分析。首先我們用MATLAB軟件畫出初始波形，以便于對比（這里三個(gè)格式?jīng)]有區(qū)別）。A格式，根據(jù)輸出結(jié)果，選取不同時(shí)間節(jié)點(diǎn)用MATLAB軟件繪制波形圖像，這里從初始時(shí)刻開始每隔10單位取一點(diǎn)（即取t=1,11,21， 畫在同一張圖里，可明顯看出隨時(shí)間變化的趨勢）。波形畫出如下：（下一頁第一張）從圖中看出，隨著時(shí)間推移A格式呈現(xiàn)明顯的發(fā)散狀態(tài)。B格式圖像如下（下一頁第二張）：從圖中看出，隨著時(shí)間推移B格式也呈現(xiàn)明顯的發(fā)散狀態(tài)。C格式圖像如下：對比之下明顯可以看出，隨著時(shí)間的推移，C格式能較好的保證收斂性。同樣的，附上xt=0.5以及xt=2時(shí)的作圖結(jié)果，從圖中可以看出xt=0.5時(shí)三種格式都呈發(fā)散狀態(tài)，而xt=2時(shí)C格式能保持收斂性，A、B格式都呈發(fā)散狀態(tài)。首先給出xt=0.5時(shí)三種格式的波形圖（依次為