1、偏微分方程组解法某厚度为 10cm 平壁原温度为 20 ，现其两侧面分别维持在 20 和 120 ，试求C C经过 8 秒后平壁内温度分布，并分析温度分布随时间的变化直至温度分布稳定为止。 2xtat式中 为导温系数， ； 。a/smc2a解：模型转化为标准形式： 21xta初始条件为： 0,t边界条件为：，12,0t 2,.函数：pdefun.m%偏微分方程(一维动态传热)function c,f,s=pdefun(x,t,u,dudx)c=1/2e-4;f=dudx;s=0;icbun.m%偏微分方程初始条件(一维动态传热)function u0=icbun(x)u0=20;bcfun.m
2、%偏微分方程边界条件(一维动态传热)function pl,ql,pr,qr=bcfun(xl,ul,xr,ur,t)pl=ul-120;ql=0;pr=ur-20;qr=0;命令：x=linspace(0,10,20)*1e-2;t=linspace(0,15,16);sol=pdepe(0,pdefun,icfun,bcfun,x,t);mesh(x,t,sol(:,:,1) %温度与时间和空间位置的关系图%画 1、2、4、6、8、15s 时刻温度分布图plot(x,sol(2,:,1) 1s 时刻， （因为本题 sol 第一行为 0 时刻）hold onplot(x,sol(3,:,1)
3、plot(x,sol(5,:,1)plot(x,sol(7,:,1)plot(x,sol(9,:,1)plot(x,sol(16,:,1)计算结果：%第 8 秒时温度分布xsol(9,:,1)经过 8 秒时的温度分布为：/cm 0 0.5263 1.0526 1.5789 2.1053 2.6316 3.1579/tC120.0000 112.5520 105.1653 97.8994 90.8100 83.9477 77.3562/cmx3.6842 4.2105 4.7368 5.2632 5.7895 6.3158 6.8421/ 71.0714 65.1202 59.5200 54.2
4、784 49.3930 44.8518 40.6338/cm 7.3684 7.8947 8.4211 8.9474 9.4737 10.0000/t36.7095 33.0419 29.5877 26.2982 23.1207 20.0000或者求第 8 秒时，x=0,2,4,6,8,10cm 处的温度uout,duoutdx=pdeval(0,x,sol(9,:,:),0,2,4,6,8,10*1e-2)120.0000 92.2279 67.5007 47.5765 32.3511 20.0000不同时刻温度分布图将上图的视角转至 xt 平面也得到本图，从本图可知当时间达到 15s 时平
5、壁内的温度分布已近稳定。某厚度为 20cm 钢板原温度为 20 ，现将其置于 1000 的炉中加热，平壁导热系CC数为 ，对流传热系数 ，导温系数为 ；CW/m8.34 W/m1742h /sm105.25a试分析温度分布随时间的变化及钢板表面温度达到 500 时所需的时间。2xtat解：模型转化为标准形式： 21xta初始条件为： 0,t边界条件为：（平壁中心坐标为 0，绝热） ，0,xt xttth,1.0,1.函数：pdefun1.m%偏微分方程(一维动态平壁两侧对流)function c,f,s=pdefun1(x,t,u,dudx)c=1/0.555e-5;f=dudx;s=0;ic
6、bun1.m%偏微分方程初始条件(一维动态平壁两侧对流)function u0=icbun1(x)u0=20;bcfun1.m%偏微分方程边界条件(一维动态平壁两侧对流)function pl,ql,pr,qr=bcfun1(xl,ul,xr,ur,t)pl=0;ql=1;pr=174*(ur-1000);qr=34.8; %平壁两侧置于同一流体中具有对流传热,平壁中心为绝热命令：%600s 内的温度分布变化x=linspace(0,10,20)*1e-2;t=0:60:600;sol=pdepe(0,pdefun1,icfun1,bcfun1,x,t);mesh(x,t,sol)%2160s
7、 内的温度分布变化t=0:60:2160;sol=pdepe(0,pdefun1,icfun1,bcfun1,x,t);mesh(x,t,sol)%60、120、180、240、300、 360、420s 时刻温度分布图plot(x,sol(2,:,1) 60s(t 网格为 0:60:2160,其时间间隔为 0,60,120,180,2160,第二点为 60s)hold onplot(x,sol(3,:,1)plot(x,sol(4,:,1)plot(x,sol(5,:,1)plot(x,sol(6,:,1)plot(x,sol(7,:,1)plot(x,sol(8,:,1)%1080、144
8、0、1800、2160s 时刻温度分布图plot(x,sol(19,:,1)hold onplot(x,sol(25,:,1)plot(x,sol(31,:,1)plot(x,sol(37,:,1)600s 内的温度分布变化2160s 内的温度分布变化不同时刻温度分布图，5.08.3417hBi 2.01.0365.252aFo根据以上两准数可知：该传热过程内外对流及导热阻力相当；当加热时间小于 360s 时为非正规阶段，加热时间大于 360s 后进入正规阶段，从上图也可得到该结论。由该图可知，当加热时间达到 2160s 时，钢板的表面温度达到 500 ，钢板中心温度为：Cuout,duout
9、dx=pdeval(0,x,sol(37,:,:),0*1e-2)371.2850 C有关传热量问题：（钢板的表面温度达到 500 时，所需的总传热量）C，0001tcVQBFoA21exp对平板：（参见传热学）， ，681.17329.1cosin21A 9273.0sin15800.exp532.210 eBF4176.8.1000 tcVQJ08.02.5.17450 tcQ J1564.186.20MatLab 解法：x=0:1:10*1e-2;t=0:60:2160; （表面温度达到 500 时所需时间为 2160s）Csol=pdepe(0,pdefun1,icfun1,bcfun
10、1,x,t);size(sol)length(t)q=0;for i=1:36q=q+174*(sol(i,11,1)+sol(i+1,11,1)/2-1000)*60;（所需的热量均是从表面传递进入钢板的）endq结果： J10479.28Q有一直径为 40cm 钢锭温度为 20 ，将其置于 900 的炉中加热，平壁导热系数CC为 ，对流传热系数 ，导温系数为 ；CW/m8.34 W/m1742h /sm10695.2a试分析温度分布随时间的变化及钢锭表面温度加热到 750 时所需的时间；钢锭可近似为无限长的圆柱。 xtat初始条件为： 20,t边界条件为：（平壁中心坐标为 0，绝热） ，0
11、,xt xttth,2.0,.函数：pdefun2.micbun2.mbcfun2.m命令：x=0:1:20*1e-2;t=0:120:5520; (结合图 5520s 时，钢锭表面温度达到 750 左右)Csol=pdepe(1,pdefun2,icfun2,bcfun2,x,t); m=1(圆柱)mesh(x,t,sol)PDETOOL 工具求解二维稳态与动态 PDE：如图偏心环形空间内表面温度为 100 ，外表面温度为 20 ；试给出其温度分布。CC如图导热物体，下表面温度为 20 ，上部三角形截面处温度为 100 ，其余各面C C绝热，试给出其温度分布。有一砖砌的烟气通道，其截面形状如
12、附图所示，边长为 1m，内管道直径为 0.5m。已知内外壁温分别为 80 、25 ，砖的导热系数为 1.5 ，试确定该通道的温度C CW/m分布、距离任意相邻两直角边各 0.1m 处的温度及每米长烟道上的散热量。解：距离相邻两直角边 0.1m 处的温度，即以上图中坐标（ 0.4,0.4） 、 （-0.4,0.4） 、 （0.4,-0.4）或（-0.4,-0.4）处；命令：x=-0.5:0.05:0.5;y=-0.5:0.05:0.5;uxy=tri2grid(p,t,u,x,y);interp2(x,y,uxy,-0.4,-0.4) 结果：29.9010interp2(x,y,uxy,-0.4
13、, 0.4) 结果：29.9004interp2(x,y,uxy, 0.4,-0.4) 结果：29.9015interp2(x,y,uxy, 0.4, 0.4) 结果：29.9063故其温度为 29.9 C ntAQ根据本题其四面对称，所以计算其中一面；传递的热量在内表面处难以计算（圆形的表面） ，但传递的热量必然通过外表面，其传热量只需将表面处各节点的一阶导数与节点间的面积、导热系数相乘即可得到传热量。命令：x=-0.5:0.05:0.5;y=-0.5:0.05:0.5;uxy=tri2grid(p,t,u,x,y);dudx=( uxy(:,2)-uxy(:,1)/(x(2)-x(1) x
14、 方向一阶导数q=1.5*0.05*sum(dudx)*4 结果：675.0801dudx=( uxy(2,:)-uxy(1,:)/(y(2)-y(1) y 方向一阶导数q=1.5*0.05*sum(dudx)*4 结果：675.0683由以上命令得x 方向一阶导数：0 26.1727 47.8501 74.5754 99.8289 120.8720 141.5859 160.9858 175.3053 184.6430 188.6017 184.9701 174.9577 160.6397 141.7009 120.5195 97.6419 72.9719 50.2694 26.1751 0
15、y 方向一阶导数0 26.1727 50.2649 72.9663 97.6365 120.5146 141.6968 160.6367 174.9557 184.9691 188.6015 184.6433 175.3060 160.9868 141.5868 120.8722 99.8279 74.5720 47.8485 26.1696 0由以上边界处的一阶导数结合图可知在边界中心处的热通量最大；每米烟道的传热量为675W/m（根据传热学形状因子的计算方法得 672.8W/m）同上题，烟气通道原温度为 25 ，零时刻开始内外壁温分别维持在为 80 、25C C，砖的导温系数为 0.000
16、01 ，试确定该通道温度分布随时间的变化及温度分CW/m布稳定所需的时间。 2ytxat解：pdetool 工具各菜单中的参数设置略，其中 solve 菜单 parameters 选项的时间设置为 0:100:8000。首先在稳态情况下（椭圆形）求解（三角形网格经过两次细化） ，导出解 u；而后解以上非稳态（抛物形）并导出解 u1（三角形网格细化与稳态相同）和三角形网格参数p，e，t；命令：pdemesh(p,e,t,u1(:,1) 0spdemesh(p,e,t,u1(:,2) 100spdemesh(p,e,t,u1(:,4) 300spdemesh(p,e,t,u1(:,7) 600sp
17、demesh(p,e,t,u1(:,11) 1000spdemesh(p,e,t,u1(:,21) 2000spdemesh(p,e,t,u1(:,31) 3000spdemesh(p,e,t,u1(:,51) 5000spdemesh(p,e,t,u1(:,81) 8000s稳定时间（以某时刻的温度分布与稳态温度分布的相对误差来判断）norm(u-u1(:,21)./u) 结果：2.4511norm(u-u1(:,31)./u) 0.8079norm(u-u1(:,51)./u) 0.0770norm(u-u1(:,81)./u) 0.0047可见在 8000s 以后的温度分布已基本没有变化
18、。1u0h=1 r=0 c=1 a=0 f=1误差 plot Height3D 中选择 user entry : u-(1-x.2-y.2)/4-1 -0.50 0.51-1-0.500.5100.050.10.150.20.25-1 -0.50 0.51-1-0.500.51-1.5-1-0.500.511.5x 10-4-x-y1u0h=1 r=0 c=1 a=0 f=1-x-y需采用 solve parameters 中 use nonlinear solver 方可收敛-1 -0.50 0.51 -1 -0.50 0.51-0.08-0.06-0.04-0.0200.020.040.060.08
