资源描述:
《陶文铨 数值传热学 第二版 第五章 5-2》由会员上传分享,免费在线阅读,更多相关内容在行业资料-天天文库。
1、精确解:p=[1,5,10];x=0:1/19:1;fori=1:1:3forj=1:1:20y(i,j)=(exp(p(1,i)*19*x(1,j))-1)/(exp(p(1,i)*19)-1);endplot(x,y(i,:));holdon;end由题对中心差分、一阶迎风、混合格式进行模块编程:他们之间可以通用,只需更改ae关于p的函数即可:程序如下:(1)中心差分p=[1,5,10];fori=1:1:3ae=1-0.5*p(1,i);aw=p(1,i)+ae;ap=ae+aw;fori=1:1:18forj=1:1:20a(i,j)=
2、0;endendfori=1:1:18j=i;a(i,j)=aw;a(i,j+1)=-ap;a(i,j+2)=ae;endfori=1:1:17n=i+1;form=i:-1:1b(1,1)=a(m,n);a(m,n)=-a(i+1,n)/a(i+1,n)*b(1,1)+a(m,n);a(m,n+1)=-a(i+1,n+1)/a(i+1,n)*b(1,1)+a(m,n+1);a(m,n+2)=-a(i+1,n+2)/a(i+1,n)*b(1,1)+a(m,n+2);endendF(1)=0;F(20)=1;F(19)=(-a(1,20)*F(2
3、0)-a(1,1)*F(1))/a(1,19);fori=2:1:18F(i)=(-a(i,20)*F(20)-a(i,19)*F(19))/a(i,i);endx=0:1/19:1;y(1,:)=F;plot(x,y);holdonend(2)一阶迎风更改ae表达式:ifp(1,i)<0ae=1+p(1,i);elseae=1;end(3)混合格式更改ae表达式:ifp(1,i)>2ae=0;elseif(p(1,i)<=2)&&(p(1,i)>0)ae=1-p(1,i)/2;elseae=-p(1,i);end(4)quick格式用matl
4、ab编quick格式:程序如下:p1=[1,5,10];fori=1:1:3p=p1(1,i);F(1,1)=0;F(1,20)=1;fori=2:1:19F(1,i)=0;endF1(1,1)=0;F1(1,20)=1;fori=2:1:19F1(1,i)=0;endwhile(F1(1,2)==0)
5、
6、max>0.00000001fori=1:1:20F(1,i)=F1(1,i);endA=1;B=1/(2+p);C=(1+p)/(2+p);D(1,2)=p*((6*F(1,2)-3*F(1,1)-3*F(1,3))/8)/(2+p);fo
7、ri=3:1:19D(1,i)=p*((5*F(1,i)-F(1,i-1)-F(1,i-2)-3*F(1,i+1))/8)/(2+p);endP(1,2)=B;Q(1,2)=C*F(1,1)+D(1,2);fori=3:1:19P(1,i)=B/(A-C*P(1,i-1));Q(1,i)=(D(1,i)+C*Q(1,i-1))/(A-C*P(1,i-1));endfori=19:-1:2F1(1,i)=P(1,i)*F1(1,i+1)+Q(1,i);endmax=abs(F1(1,2)-F(1,2));fori=3:1:19ifabs(F1(1
8、,i)-F(1,i))>maxmax=abs(F1(1,i)-F(1,i));endendendx=0:1/19:1;plot(x,F1(1,:));holdon;end