薄板非线性传热研究

打开实时脚本

本例展示了如何对薄板进行传热分析。

板是方形的，温度沿底部边缘固定。没有热量从其他三个边传递(即它们是绝缘的)。热量通过对流和辐射从板的顶部和底部传递。因为辐射是包含在内的，所以这个问题是非线性的。本例的目的之一是展示如何处理PDE问题中的非线性问题。

同时进行了稳态和瞬态分析。在稳态分析中，我们感兴趣的是在板达到平衡状态后不同位置的最终温度。在瞬态分析中，我们感兴趣的是板中的温度作为时间的函数。通过这种瞬态分析可以回答的一个问题是，板块达到平衡温度需要多长时间。

平板的传热方程

平板的平面尺寸为1米乘1米，厚度为1厘米。由于平板相对于平面尺寸较薄，因此在厚度方向上可以假定温度恒定;得到的问题是二维的。

对流和辐射换热被假定发生在板的两个面和指定的环境温度之间。

每单位面积上每个板面由于对流而传递的热量定义为

$问_{c} ＝ h_{c} （ T - T_{一个} ）$

在哪里 $T_{一个}$ 是环境温度， $T$ 在平板表面某一特定x和y位置的温度，和 $h_{c}$ 是指定的对流系数。

每单位面积上由于辐射而从每个板面传递的热量定义为

$问_{r} ＝ ϵ σ （ T^{4} - T_{一个}^{4} ）$

在哪里 $ϵ$ 面和的发射率 $σ$ 为斯特凡-玻尔兹曼常数。由于辐射传递的热量与表面温度的四次幂成正比，所以这个问题是非线性的。

描述薄板内温度的偏微分方程为

$ρ C_{p} t_{z} \frac{\partial T}{\partial t} - k t_{z} \nabla^{2} T + 2 问_{c} + 2 问_{r} ＝ 0$

在哪里 $ρ$ 是物质密度， $C_{p}$ 是比热， $t_{z}$ 一是板厚，二是板面传热的因素。

用PDE工具箱所期望的形式重写这个方程是很方便的

$ρ C_{p} t_{z} \frac{\partial T}{\partial t} - k t_{z} \nabla^{2} T + 2 h_{c} T + 2 ϵ σ T^{4} ＝ 2 h_{c} T_{一个} + 2 ϵ σ T_{一个}^{4}$

问题的设置

该板由铜组成，具有以下特性:

K = 400;铜的热导率%，W/(m-K)Rho = 8960;铜的密度%，kg/m^3specificHeat = 386;铜比热% J/(kg-K)厚= .01;%板厚，单位为米stefanBoltz = 5.670373e-8;% Stefan-Boltzmann常数，W/(m^2-K^4)hCoeff = 1;%对流系数W/(m^2-K)假定环境温度为300开尔文。Ta = 300;发射热量= .5;板表面的%发射率

创建带有单个因变量的PDE模型。

numberOfPDE = 1;model = createpde(numberOfPDE);

对于一个正方形，几何图形和网格很容易定义如下所示。

宽度= 1;高度= 1;

通过给出4个角的x-位置和4个角的y-位置来定义正方形。

GDM =[3 4 0宽度宽度0 0 0高度高度]';G = decsg(gdm，“S1 ', (“S1 ') ');

将DECSG几何图形转换为几何对象，这样它就被附加到PDEModel中

geometryFromEdges(模型中,g);

绘制几何图形并显示用于边界条件定义的边缘标签。

图;pdegplot(模型,“EdgeLabels”，“上”）;轴([-。1 1.1 -。1 1.1]);标题(显示边缘标签的几何图形）

图中包含一个轴对象。标题为Geometry with Edge Labels的axis对象包含5个类型为line, text的对象。

指定系数。通过将上述方程与PDE工具箱文档中的标量抛物方程进行比较，PDE工具箱所需的系数表达式可以很容易地识别出来。

C =厚*k;

由于辐射边界条件的存在，“a”系数是温度u的函数。它被定义为MATLAB表达式，因此在分析过程中可以对不同的u值进行计算。

a = @(~，状态)2*hCoeff + 2*emiss*stefanBoltz*state.u.^3;f = 2*hCoeff*ta + 2*emiss*stefanBoltz*ta^4;d = thick*rho*specificHeat;specifyCoefficients(模型,“m”0,“d”0,“c”c“一个”一个,“f”f);

板的底边设置为1000开尔文。

应用边界条件。三个板边是绝缘的。由于诺伊曼边界条件等于零是有限元公式的默认值，这些边的边界条件不需要明确设置。在底边，边1上的所有节点上设置狄利克雷条件，

applyBoundaryCondition(模型,“边界条件”，“边缘”,1，“u”, 1000);

指定初始猜测。

setInitialConditions(模型中,0);

在正方形上创建三角形网格，每个方向上大约有十个元素。

Hmax = .1;元素尺寸百分比msh = generateMesh(模型，“Hmax”, hmax);图;pdeplot(模型);轴平等的标题(“三角网格板”)包含(“x坐标,米”) ylabel (“坐标,米”）

图中包含一个轴对象。标题为Plate with triangle Element Mesh的axis对象包含2个类型为line的对象。

稳态解

由于a和f系数是温度的函数(由于辐射边界条件)，solvepde自动选取非线性解算器得到解。

R = solvepde(模型);u = r . nodesolution;图;pdeplot(模型,“XYData”u“轮廓”，“上”，“ColorMap”，“喷气机”）;标题(板内温度，稳态解)包含(“x坐标,米”) ylabel (“坐标,米”)轴平等的

图中包含一个轴对象。标题为Temperature In The Plate, Steady State Solution的坐标轴对象包含12个类型为patch, line的对象。

p = msh.Nodes;plotAlongY (p, u, 0);标题(温度作为y坐标的函数)包含(“坐标,米”) ylabel (“温度,开尔文”）

图中包含一个轴对象。标题为Temperature As a Function of Y-Coordinate的axis对象包含一个类型为line的对象。

流(['板上边缘的温度='.．.' %5.1f degrees-K\n'), u (4));

板顶部边缘的温度= 449.8度- k

临时的解决方案

包括d系数。

specifyCoefficients(模型,“m”0,“d”d“c”c“一个”一个,“f”f);endTime = 5000;tlist = 0:50:endTime;numNodes = size(p,2);

设置所有节点的初始温度为ambient, 300k。

u0(1:numNodes) = 300;

设置底边E1的初始温度为常数BC, 1000 K。

1000年setInitialConditions(模型,“边缘”1);

设置以下解算器选项。

model.SolverOptions.RelativeTolerance = 1.0e-3;model.SolverOptions.AbsoluteTolerance = 1.0e-4;

通过使用来解决问题solvepde．求解器自动选取抛物线求解器来获得解。

R = solvepde(model,tlist);u = r . nodesolution;图;:情节(tlist u (3));网格在标题([沿顶部边缘的温度.．.“板块的时间函数”])包含(“时间,秒”) ylabel (“温度,开尔文”）

图中包含一个轴对象。标题为Temperature Along The Top Edge of The Plate as a Function of Time的坐标轴对象包含一个line类型的对象。

图;pdeplot(模型,“XYData”u(:,结束),“轮廓”，“上”，“ColorMap”，“喷气机”）;标题(sprintf ([《盘子里的温度》.．.'瞬态解(%d秒)\n'), tlist (1)));包含(“x坐标,米”) ylabel (“坐标,米”)轴平等的；

图中包含一个轴对象。标题为Temperature In The Plate,Transient Solution(5000秒)的坐标轴对象包含12个类型为patch, line的对象。

流(['\ ntop edge(t = %5.1f secs) = '.．.' % 5.1 f degrees-K \ n”), tlist(结束),u (4));

上边缘的温度(t = 5000.0秒)= 441.8度- k

总结

从稳态解和结束时刻的瞬态解得到的板内温度图非常接近。也就是说，经过5000秒左右，暂态解已经达到稳态值。两种溶液在平板顶部边缘的温度一致在1%以内。金宝搏官方网站