清华大学传热大作业-无限大平壁简单数值计算报告(共6页)

1、精选优质文档-倾情为你奉上传热大作业报告* 热动* *一、 大作业题目一厚度为0.1m的无限大平壁，两侧均为对流换热边界条件，初始时 两侧流体温度与壁内温度一致，tf1=tf2=t0=5 ；已知两侧对流换热系数分别为h1=11 W/m2K、h2=23W/m2K, 平壁材料的导热系数l=0.43W/mK，导温系数a=0.3437×10-6 m2/s。如果一侧的环境温度tf1突然升高为50并维持不变，计算在其它参数不变的条件下，平壁内温度分布及两侧壁面热流密度随时间的变化规律（用图形表示）。要求：将全部计算内容（包括网格的划分、节点方程组、计算框图、程序及计算结果）用A4纸打印。二、 网

2、格划分如图，将无限大平板作为一维处理，本题为一维非稳态导热问题，对流换热边界条件。l 空间网格划分：平板总厚度为delta=0.1m，定义空间步长为 dx=0.005m，则距离份数为N=delta/dx=20份。定义xn为以0为首项，以dx为公差的等差数列，尾项为delta=0.1m，共有N+1项，则xn中的每一项即表示一个沿平板厚度方向中的划分点。l 时间网格划分：设总时间长度为T= s，定义时间步长为dtao=20s，则时间份数为M=5000份。定义 taom是以0为首项，以dtao为公差的等差数列，尾项为T=s，共有M+1项,则taom中每一项即表示一个时刻。三、 计算框图l 程序中的各

3、个变量的名称及意义：1. 题设中各个常数lambda=0.43 导热系数；a=0.3437e-6 热扩散率；h1=11 边界对流换热系数；h2=23 边界对流换热系数2；t0=5 初始温度；tf1=50 初始流体温度；tf2=5 初始流体温度2；delta=0.1 总距离长度（无限大平板厚度）；2. 网格划分所设的变量T= 总时间长度（在T时间内考虑本问题）；dtao=20 定义时间步长；dx=0.005 定义距离步长；M=floor(T/dtao) 时间份数=总时间/时间步长（向下取整）；N=floor(delta/dx 距离份数=总厚度/距离步长(向下取整）；tao=0:dtao:T 定义

4、时间划分单元（以0为首项，以dtao为公差的等差数列，尾项为T），共有M+1项；x=0:dx:delta 定义距离划分单元（以0为首项，以dx为公差的等差数列，尾项为delta），共有N+1项；3. 判定稳定性的准则数Bi1=h1*dx/lambda 边界节点网格毕渥数；Bi2=h2*dx/lambda 边界节点网格毕渥数2；Fo=a*dtao/dx2 傅里叶数；l 程序计算框图开始输入delta,T,dtao,dx,M,N,tao,xFo>1/(2*Bi1+2)&&Fo>1/(2*Bi2+2)NoYes建立t(M+1,N+1)温度矩阵，令t(1,:)=t0,令m=

5、1m=m+1打印“不稳定”t(m,1)=2*Fo*(t(m-1,2)+Bi1*tf1)+(1-2*Bi1*Fo-2*Fo)*t(m-1,1)t(m,N+1)=2*Fo*(t(m-1,N)+Bi2*tf2)+(1-2*Bi2*Fo-2*Fo)*t(m-1,N+1)t(m,n)=Fo*(t(m-1,n-1)+t(m-1,n+1)+(1-2*Fo)*t(m-1,n)m>M+1NoYes输出温度矩阵t(M+1,N+1)和相应图象停机四、 程序代码本程序在MATLAB R2008a中运行通过，以下是源代码（%后为注释）：lambda=0.43;%导热系数a=0.3437e-6;%热扩散率h1=11

6、;%边界对流换热系数h2=23;%边界对流换热系数2t0=5;%初始温度tf1=50;%初始流体温度tf2=5;%初始流体温度2delta=0.1;%总距离长度（无限大平板厚度）T=;%总时间长度（在T时间内考虑本问题）dtao=20;%定义时间步长dx=0.005;%定义距离步长M=floor(T/dtao);%时间份数=总时间/时间步长（向下取整）N=floor(delta/dx);%距离份数=总厚度/距离步长(向下取整）tao=0:dtao:T;%定义时间划分单元（以0为首项，以dtao为公差的等差数列，尾项为T），共有M+1项x=0:dx:delta;%定义距离划分单元（以0为首项，以

7、dx为公差的等差数列，尾项为delta），共有N+1项Bi1=h1*dx/lambda;%边界节点网格毕渥数Bi2=h2*dx/lambda;%边界节点网格毕渥数2Fo=a*dtao/dx2;%傅里叶数if Fo>1/(2*Bi1+2)&&Fo>1/(2*Bi2+2)%判断稳定性，不稳定则显示毕渥数、傅里叶数 disp('不稳定'); disp(Bi1); disp(Bi2); disp(Fo); disp(1/(2*Bi1+2); disp(1/(2*Bi2+2);else%若稳定，则进行迭代计算 t=zeros(M+1,N+1);%建立一个(M+

8、1)*(N+1)的温度矩阵，M+1为时间节点个数，N+1为空间节点个数，以便进行迭代计算 q1=zeros(M+1,1);%根据题目要求算两壁面处热流密度 q2=zeros(M+1,1); t(1,:)=t0;%初始温度均为t0=5 for m=2:M+1%m=1时是初值上一行已计算出，则从m=2一直计算到m=M+1，m对应的时刻是tao=（m-1)dtao t(m,1)=2*Fo*(t(m-1,2)+Bi1*tf1)+(1-2*Bi1*Fo-2*Fo)*t(m-1,1);%首先计算一边界这个时刻温度 t(m,N+1)=2*Fo*(t(m-1,N)+Bi2*tf2)+(1-2*Bi2*Fo-2

9、*Fo)*t(m-1,N+1);%再计算另一边界这个时刻的温度 q1(m)=h1*(tf1-t(m,1); q2(m)=h2*(t(m,N+1)-tf2); for n=2:N%然后计算内部，n=1和n=N+1时是边界节点温度，上面两行已经计算出，n对应的坐标是x=(n-1)*dx t(m,n)=Fo*(t(m-1,n-1)+t(m-1,n+1)+(1-2*Fo)*t(m-1,n); end end %以下是画图 figure plot(x,t(1,:),x,t(11,:),x,t(21,:),x,t(51,:),x,t(101,:),x,t(1001,:),x,t(5001,:); lege

10、nd('t=0s','t=200s','t=400s','t=1000s','t=2000s','t=20000s','t=',0); title('一定时间下温度随距离的分布','fontsize',12,'fontweight','bold','fontname','楷体'); axis(0,0.1,0,40); figure plot(tao,t(:,1),tao,t(:,6),ta

11、o,t(:,11),tao,t(:,16),tao,t(:,21); legend('x=0','x=0.025','x=0.05','x=0.075','0.1',0); title('一定位置处温度随时间的分布','fontsize',12,'fontweight','bold','fontname','楷体'); axis(0,0,40); figure mesh(x,tao,t); title('温度随时间和空间的分布','fontsize',12,'fontweight','bold','fontname','楷体'); figure plot(tao,q1,tao,q2); legend('q1','q2'); title('两壁面热流密度随时间变化曲线','fontsize',12,'fontweight','bold','fon