搜索
传热大作业二维导热物体温度场的数值模拟姓名：刘璇班级：能动A02学号：100310961一．物理问题有一个用砖砌成的长方形截面的冷空气通道，其截面尺寸如下图所示，假设在垂直于纸面方向上用冷空气及砖墙的温度变化很小，可以近似地予以忽略。在下列两种情况下试计算：（1）砖墙横截面上的温度分布；（2）垂直于纸面方向的每米长度上通过砖墙的导热量。第一种情况：内外壁分别均与地维持在0℃及30℃；第二种情况：内外壁均为第三类边界条件，且已知：𝑡∞1=30℃,ℎ1=10𝑤𝑚2∙℃𝑡∞2=10℃,ℎ2=4𝑤𝑚2∙℃砖墙的导热系数λ=0.53𝑊𝑚∙℃2二．数学描写由对称的界面必是绝热面，可取左上方的四分之一墙角为研究对象，该问题为二维、稳态、无内热源的导热问题，其控制方程和边界条件如下：𝜕𝑡2𝜕𝑥2+𝜕𝑡2𝜕𝑦2=0边界条件（情况一）t(x,0)=300≤x≤1.5t(0,y)=300≤y≤1.1t(0.5,y)=00.5≤y≤1.1t(x,0.5)=00.5≤x≤1.5𝜕𝑡(1.5,𝑦)𝜕𝑦=00≤y≤0.5∂t(x,1.1)∂x=00≤x≤0.5边界条件（情况二）λ𝜕𝑡𝜕𝑥=ℎ1(𝑡−𝑡𝑓1)𝑥=0,0≤𝑦≤1.1−λ𝜕𝑡𝜕𝑥=ℎ2(𝑡−𝑡𝑓2)𝑥=0.5,0.5≤𝑦≤1.1λ𝜕𝑡𝜕𝑦=ℎ1(𝑡−𝑡𝑓1)𝑦=0,0≤𝑥≤1.5−λ𝜕𝑡𝜕𝑦=ℎ2(𝑡−𝑡𝑓2)𝑦=0,0.5≤𝑥≤1.5∂t(1.5,y)∂y=00≤y≤0.5∂t(x,1.1)∂x=00≤x≤0.5三．网格划分网格划分与传热学实验指导书中“二维导热物体温度场的电模拟实验”一致，如下图所示：3四．方程离散对于内节点，离散方程t[i][j]=0.25*(t[i+1][j]+t[i-1][j]+t[i][j+1]+t[i][j-1])对于边界节点，则应对一、二两种情况分开讨论：情况一:绝热平直边界点：t[15][j]=0.25*(2*t[14][j]+t[15][j-1]+t[15][j+1])1≤j≤4t[i][11]=0.25*(2*t[i][10]+t[i-1][11]+t[i+1][11])1≤i≤4外等温边界点：t[i][j]=30内等温边界点：t[i][j]=0情况二：（Bi1,Bi2为网格Bi数，Bi1=ℎ1∆𝑥λBi2=ℎ2∆𝑥λ）绝热平直边界点：t[15][j]=0.25*(2*t[14][j]+t[15][j-1]+t[15][j+1])1≤j≤4t[i][11]=0.25*(2*t[i][10]+t[i-1][11]+t[i+1][11])1≤i≤4外侧对流平直边界：4t[i][0]=(2*t[i][1]+t[i+1][0]+t[i-1][0]+2*Bi1*tf1)/(2*Bi1+4)1≤i≤14t[0][j]=(2*t[1][j]+t[0][j+1]+t[0][j-1]+2*Bi1*tf1)/(2*Bi1+4)1≤j≤10内侧对流平直边界：t[i][5]=(2*t[i][4]+t[i+1][5]+t[i-1][5]+2*Bi2*tf2)/(2*Bi2+4)6≤i≤14t[5][j]=(2*t[4][j]+t[5][j+1]+t[5][j-1]+2*Bi2*tf2)/(2*Bi2+4)6≤j≤10特殊点：a点t[15][0]=(t[14][0]+t[15][1]+tf1*Bi1)/(Bi1+2)b点t[15][5]=(t[14][5]+t[15][4]+tf2*Bi2)/(Bi2+2)c点t[5][5]=(2*t[4][5]+2*t[5][4]+t[5][6]+t[6][5]+3*Bi2*tf2)/(2*Bi2+6)d点t[5][11]=(t[5][10]+t[4][11]+tf2*Bi2)/(Bi2+2)e点t[0][11]=(t[0][10]+t[1][11]+tf1*Bi1)/(Bi1+2)f点t[0][0]=(t[0][1]+t[1][0]+tf1*Bi1*2)/(2*Bi1+2)五．编程思路及流程图编程思路为设定两个二维数组t[i][j]、ta[i][j]分别表示本次迭代和上次迭代各节点的温度值，iter表示迭代进行的次数,daore_in、daore_out分别表示内外边界的散热量。开始时，给t[i][j]、ta[i][j]赋相同的初始值，t[i][j]根据内节点和各边界节点的离散方程进行迭代，迭代后比较t[i][j]、ta[i][j]各个节点之间温度之差，若两个温度之差小给定的精度，则此时迭代完成，t[i][j]就是所求的温度场分布，若两温度之差不满足精度要求，则将t[i][j]的值赋给ta[i][j]，t[i][j]继续迭代，直到二者各个点的温度之差满足精度要求，记下此时的迭代次数，并根据所得到的温度场分布计算内外边界上散热量以及偏差。5输入已知参数说明边界条件Iter=1计算新的内节点温度及新的边界点温度t[i][j]取定初始试探值t[i][j]=0ta[i][j]=t[i][j]比较所有节点|t[i][j]—ta[i][j]|精度ta[i][j]=t[i][j]计算内外边界上散热量及其平均值、偏差输出t[i][j]、iter平均导热量及偏差结束开始6六．实验结果等温边界程序运行结果：对流边界程序运行结果：7等温边界节点温度分布图8对流边界节点温度分布图：9七．结果讨论1.由实验结果可知：等温边界下，内外边界散热量分别为241.52,242.12，平均值241.82，偏差0.25%；对流边界下，内外边界散热量分别为118.13,111.16，平均值114.14，偏差5.9%，这与“二维导热物体温度场的电模拟实验“结果相似，说明了数值解法分析问题的可行性。用数值解法仅用计算机模拟就能解决某些复杂的工程问题，为复杂工程问题的求解提供了极大的便利。2.在实验中，内外边界散热量存在偏差，这在很大程度上是由于用数值计算分析问题时，采用离散平均的思想，用节点中心的温度代替节点的平均温度从而产生误差。不断提高所划分的网格数目，实验偏差会得到不断改善。3.由所做的温度分布图可知，温度分布大致对称与对角线fc,这是由于对角线两侧几何对称，而且两侧的控制方程以及边界条件都一致。4.由所做的温度分布图可知，等温线垂直于绝热边界，在绝热边界处，没有热流量，温度梯度为零。附源程序：1.等温边界（情况一）#includestdio.h#includemath.hintmain(){10intiter=0,n=0;doublet[16][12]={0},ta[16][12]={0};doubleepsilon=5.0e-3;doublelambda=0.53,error=0;doubledaore_in=0,daore_out=0,daore=0;FILE*fp;fp=fopen(data3,w);for(inti=0;i=15;i++)for(intj=0;j=11;j++){if((i==0)||(j==0))ta[i][j]=30;if(i==5)if(j=5&&j=11)ta[i][j]=0;if(j==5)if(i=5&&i=15)ta[i][j]=0;}for(inti=0;i=15;i++)for(intj=0;j=11;j++)t[i][j]=ta[i][j];n=1;while(n0){n=0;for(intj=1;j=4;j++)t[15][j]=0.25*(2*t[14][j]+t[15][j-1]+t[15][j+1]);for(inti=1;i=4;i++)t[i][11]=0.25*(2*t[i][10]+t[i-1][11]+t[i+1][11]);for(inti=1;i=14;i++)for(intj=1;j=4;j++)t[i][j]=0.25*(t[i+1][j]+t[i-1][j]+t[i][j+1]+t[i][j-1]);for(inti=1;i=4;i++)for(intj=5;j=10;j++)t[i][j]=0.25*(t[i+1][j]+t[i-1][j]+t[i][j+1]+t[i][j-1]);for(inti=0;i=15;i++)for(intj=0;j=11;j++)if(fabs(t[i][j]-ta[i][j])epsilon)n++;for(inti=0;i=15;i++)for(intj=0;j=11;j++)ta[i][j]=t[i][j];11iter++;//printf(%d\n,iter);}for(intj=0;j=5;j++){for(inti=0;i=15;i++){printf(%4.1f,t[i][j]);fprintf(fp,%4.1f,t[i][j]);}printf(\n);fprintf(fp,\n);}for(intj=6;j=11;j++){for(inti=0;i=5;i++){printf(%4.1f,t[i][j]);fprintf(fp,%4.1f,t[i][j]);}fprintf(fp,\n);printf(\n);}for(inti=1;i=14;i++)daore_out+=(30-t[i][1]);for(intj=1;j=10;j++)daore_out+=(30-t[1][j]);daore_out=4*(lambda*(daore_out+0.5*(30-t[1][11])+0.5*(30-t[15][1])));for(inti=5;i=14;i++)daore_in+=t[i][4];for(intj=5;j=10;j++)daore_in+=t[4][j];daore_in=4*(lambda*(daore_in+0.5*t[4][11]+0.5*t[15][4]));error=abs(daore_out-daore_in)/(0.5*(daore_in+daore_out));daore=(daore_in+daore_out)*0.5;printf(iter=%d\ndaore_in=%f\ndaore_out=%f\ndaore=%f\nerror=%f\n,iter,daore_in,daore_out,daore,error);}2.对流边界（情况二）#includestdio.h#includemath.hintmain(){doubletf1=30,tf2=10,h1=10,h2=4,dy=0.1,dx=0.1,lambda=0.53,Bi1,Bi2;12doublet[16][12],ta[16][12];intiter=0,n=0;doubleepsilon=5e-3;doubleerror=0;doubledaore_in=0,daore_out=0,daore=0;Bi1=h1*dx/lambda;Bi2=h2*dx/lambda;for(inti=0;i=15;i++)for(intj=0;j=11;j++){t[i][j]=0;ta[i][j]=0;}FILE*fp;fp=fopen(data2,w);n=1;while(n0){n=0;t[15][0]=(t[14][0]+t[15][1]+tf1*Bi1)/(Bi1+2);t[15][5]=(t[14][5]+t[15][4]+tf2*Bi2)/(Bi2+2);t[0][11]=(t[0][10]+t[1][11]+tf1*Bi1)/(Bi1+2);t[5][11]=(t[5][10]+t[4][11]+tf2*Bi2)/(Bi2+2);t[0][0]=(t[0][1]+t[1][0]+tf1*Bi1*2)/(2*Bi1+2);t[5][5]=(2*t[4][5]+2*t[5][4]+t[5][6]+t[6][5]+3*Bi2*tf2)/(2*Bi2+6)