一维对流方程在A、B、C三种差分格式

1、一、上机目的用数值方法计算一维对流方程在A、B、C三种差分格式下的解。取为0.05. 取值为0.5，1，2。并作相关讨论。二、实验原理三、上机要求：1.学会对MS-FORTRAN的基本操作。 2.用Fortran编写程序计算一维对流方程在A、B、C三种格式下的解。3.讨论各种格式下不同的值的差分格式解的特点。四、实验程序以A格式为例，对微分方程进行离散化， 得出 A 格式的差分解的表达式：B、C格式同理可以写出。由此编写如下的Fortran程序。注：除了循环时间层的计算公式略有不同外，三个程序没有区别，因此这里只用一个主程序，并根据格式选择的不同，采用条件语句判断执行哪一部分循环体。!空间节点

2、321，dx=0.05 输出依次为时间，空间，数值program mainimplicit nonereal dx_dt !定义x/t的值integer abc,r_t,i,j,k !定义变量，abc为格式类型，r_t为时间网格数，其余为循环变量real,allocatable:s(:,:) !定义存储矩阵swrite(*,*) 输入dx_dt=0.5,1,2read(*,*) dx_dtwrite(*,*) 选择格式，A,B,C分别输入1,2或3 read(*,*) abc !根据格式选择生成相应的文件 if(abc=1) then open(unit=8,file=out_a.csv) e

3、lseif(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 !计算时间网格总数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*(i-1

4、61)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 !循环时间层,根据格式的选择来判断执行哪一部分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 s(

5、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(i-

6、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(*,*)数据已输出至源目录pausestopend program五、实验结果及分析程序运行后在对应目录下生成csv表格文件，根据输入的xt的值不同生成对应的网格并计算各节点数值。这里以xt=1为例，就A、B、C三种格式的结果进行分析。首先我们用MATLAB软件画出初始波形，以便于比照这里三个格式没有区别。A格式，根据输出结果，选取不同时间节点用MATLAB软件绘制波形图像，这里从初始时刻开始每隔10单位取一点即取t=1,11,21， 画在同一张图里，可明显看出随时间变化的趋势。波形画出如下：下一页第一张从图中看出，随着时间推移A格式呈现明显的发散状态。B格式图像如下下一页第二张：从图中看出，随着时间推移B格式也呈现明显的发散状态。C格式图像如下：比照之下明显可以看出，随着时间的推移，C格式能较好的保证收敛性。同样的，附上xt=0.5以及xt=2时的作图结果，从图中可以看出xt=0.5时三种格式都呈发散状态，而xt=2时C格式能保持收敛性，A、B格式都呈发散状态。首先给出xt=0.5时三种格式的波形图依次为A、B、C格式