非线性对流扩散方程不同解法稳定性比较(共7页)

1、PAGE 20132014学年第二(d r)学期Matlab编程技术作业(zuy) 专业(zhuny)班级 石工博13-02 研究方向 油气田开发 姓 名 王壮壮 学 号 B13020075 评 分 项 目一、选题的创新性或与专业相结合的程度(总分20分)二、选题的意义及难度(总分20分)三、程序的可读性、逻辑性(总分20分)四、程序运行情况（成果、效率等）(总分30分)五、其它（格式等）（10分）综合得分(总分100分)结合自己研究(ynji)方向，运用Matlab编写科学计算及可视化或其它相关(xinggun)程序。要求(yoqi)：将要解决的问题交代清楚（数学模型、目标等）；编写的程序的

2、关键语句要有注释说明；程序能顺利运行，运行结果和编写的m文件一并提交；独立完成。非线性对流扩散方程(fngchng)不同解法稳定性比较流体力学基本方程组本身(bnshn)就是非线性的对流扩散方程，非线性Burgers方程(fngchng)就是N-S方程很好的模型方程，它的一维形式如下： (1)边界条件为 (2)初始条件是任意可以给出的。我们知道，遇到对流项，我们用迎风格式是绝对没有问题，无论是一阶迎风还是二阶迎风格式都是能够解决非线性对流方程的，如果网格Peclet数允许的话，中心差分也是可以考虑的。不过，对于非线性对流，我们来看看另外两个G-S格式，一个是G-S型迎风半隐格式，另一个是G-S

3、型Samarskii半隐格式，对于任何类型的对流扩散方程，可以收敛到定常解，并且是绝对稳定的，特别适合于解决定常问题。对于式(1)这两个格式分别为 (3) (4)其中式(3)就是G-S型迎风半隐格式，它具有一阶精度，是从一阶迎风格式发展而来的；式(4)是G-S型Samarskii半隐格式，具有二阶精度，它是从Samarskii格式发展而来的。上面说过，它们只适用于求解定常解，因此上标的时层n可以看作是迭代步，可以说它们没有时间精度，如果想用这两个格式求解非定常解，那可是徒劳的。对于上两式，我们可以采用迭代法求解，把它们写成迭代式，分别为 (5) (6)这样我们可以看出来(ch li)，其实时间

4、步长就充当了一个(y )松弛因子的角色。首先(shuxin)式(1)是有定常精确解的，解为 (7)式中 (8) (9)程序3.3 非线性Burgers方程计算程序%程序开始 clear Re=100; %式(8) L=5; %计算域，可以随意设置 u0=1; %内边界条件，u0可随意 miu=u0*L/Re; %扩散系数 N=41; %网格节点数 t=2; %时间步长(松弛因子) h=L/(N-1); %网格步长 for i=1:N %网格划分 x(i)=(i-1)*h; end%G-S型迎风半隐格式计算 u=zeros(N,1); %迭代初值 u(1)=u0; %内边界条件 u1=u;tol

5、=1;p=0; while tol=1e-5 %收敛精度，0.00001 for i=2:N-1 R(i)=abs(u(i)*h/2/miu; u1(i)=(-t*h*(u(i)*(u(i+1)-u1(i-1)+2*t*miu*(1+R(i)*(u(i+1)+u1(i-1)+2*h2*u(i)/(2*h2+4*t*miu*(1+R(i); end tol=max(abs(u-u1); %计算两个(lin )迭代层的误差 u=u1; p=p+1; if p=5000 disp(G-S型迎风半隐格式(g shi)不收敛！) break end end%G-S型Samarskii型半隐格式(g sh

6、i)计算 tol=1;p=0; v=zeros(N,1); %迭代初值 v(1)=u0; %内边界条件 v1=v; while tol=1e-5 %收敛精度，0.00001 for i=2:N-1 R(i)=abs(v(i)*h/2/miu; v1(i)=(-t*h*(v(i)*(v(i+1)-v1(i-1)+2*t*miu*(1/(1+R(i)+R(i)*(v(i+1)+v1(i-1)+2*h2*v(i)/(2*h2+4*t*miu*(1/(1+R(i)+R(i); end tol=max(abs(v1-v); %计算量迭代层误差 v=v1; p=p+1; if p=5000 disp(G-

7、S型Samarskii型半隐格式不收敛！) break end end%牛顿迭代法计算 tol=1;U=2; %误差；迭代初值 while tol=1e-10 %收敛精度0.0000000001 f=(U-1)/(U+1)-exp(-U*Re); %式(9)建立的函数 fd=1/(U+1)-(U-1)/(U+1)2+Re/exp(U*Re); %式(9)建立的函数的导数 U1=U-f/fd; tol=abs(U-U1); %计算误差 U=U1; end%计算精确解 for i=1:N ue(i)=u0*U*(1-exp(U*Re*(x(i)/L-1)/(1+exp(U*Re*(x(i)/L-1); end plot(x,ue,-*,x,u,-or,x,v,-g) %作图Re=5Re=20Re=100式(8)定义(dngy)的参数就相当于流体力学中的雷诺数，如果雷诺数较大的话，表明(biomng)对流占优。上三个图就表示了随增大的各个(gg)解的情况。我们看到在比较小的时候，G-S迎风格式和G-S型Samarskii格式都与精确解十分接近，Sanaskii格式的精度更高一些。在较大的情况下，G-S迎风格式脱离精确解的情况比较严重，可以明显看出它伴随着很强的假扩散效应，而Samarskii格式却不存在这种现象，或者