炉子的热计算

**图14:** 炉子描述
$\begin{figure}\epsfig{file=furngeo.eps,width=11cm}\end{figure}$

本问题涉及图14所示炉子的热计算。炉子由一块底板组成，底板温度是预设的。它在极短的时间内从300 K线性变化到1000 K，之后保持恒定。炉子的侧壁与外界隔离，但通过辐射与其他炉壁进行热交换。侧壁和底板的发射率为 $\epsilon=1$ 。炉顶通过辐射与其他炉壁以及固定温度为300 K的环境温度进行热交换。炉顶的发射率为 $\epsilon=0.8$ 。此外，炉顶还通过对流与以恒定速率0.001 kg/s流动的流体（空气）进行热交换。流体在右上角的温度为300 K。炉壁由10厘米厚的钢制成。钢的材料常数为：热导率 $\kappa=50 \mathrm{W/m K}$ ，比热 $c=446 \mathrm{W/kg K}$ 和密度 $\rho=7800 \mathrm{kg/m^3}$ 。空气的材料常数为：比热 $c_p=1000 \mathrm{W/kg K}$ 和密度 $\rho=1 \mathrm{kg/m^3}$ 。对流系数为 $h=25 \mathrm{W/m^2 K}$ 。炉子的尺寸为 $0.3 \times 0.3 \times 0.3 \mathrm{m^3}$ （立方体）。在时，所有部分温度均为 $T=300 \mathrm{K}$ 。我们想知道位置A、B、C、D和E的温度随时间的变化。

**
**   结构：炉子。
**   测试目标：带对流和辐射的壳单元。
**
*NODE, NSET=Nall
       1,  3.00000e-01,  3.72529e-09,  3.72529e-09 
       ...
*ELEMENT, TYPE=S6, ELSET=furnace
     1,      1,      2,      3,      4,      5,      6
     ...
*ELEMENT,TYPE=D,ELSET=EGAS
301,603,609,604
...
*NSET,NSET=NGAS,GENERATE
603,608
*NSET,NSET=Ndown 
1, 
...
*PHYSICAL CONSTANTS,ABSOLUTE ZERO=0.,STEFAN BOLTZMANN=5.669E-8
*MATERIAL,NAME=STEEL
*DENSITY
7800.
*CONDUCTIVITY
50.
*SPECIFIC HEAT
446.
*SHELL SECTION,ELSET=furnace,MATERIAL=STEEL
0.01
*MATERIAL,NAME=GAS
*DENSITY
1.
*SPECIFIC HEAT
1000.
*FLUID SECTION,ELSET=EGAS,MATERIAL=GAS
*INITIAL CONDITIONS,TYPE=TEMPERATURE
Nall,300.
*AMPLITUDE,NAME=A1
0.,.3,1.,1.
*STEP,INC=100
*HEAT TRANSFER
0.1,1.
*VIEWFACTOR,WRITE
*BOUNDARY,AMPLITUDE=A1
Ndown,11,11,1000.
*BOUNDARY
603,11,11,300.
*BOUNDARY,MASS FLOW
609,1,1,0.001
...
*RADIATE
** Radiate based on down
1, R1CR,1000., 1.000000e+00
...
** Radiate based on top
51, R1CR,1000., 8.000000e-01
...
** Radiate based on side
101, R1CR,1000., 1.000000e+00
...
** Radiate based on top
51, R2,300., 8.000000e-01
...
*FILM
51, F2FC, 604, 2.500000e+01
...
*NODE FILE
NT
*NODE PRINT,NSET=NGAS
NT
*END STEP

上面列出了输入文件。它以节点定义开始。结构中最高节点号为602。节点603到608是流体节点，即在流体中定义了额外节点（z=0.3对应炉顶，z=0对应底部）。流体节点603对应流体温度为300 K（"入口"）的位置，节点608对应"出口"，其他节点位于两者之间。实际上，流体节点的坐标不进入计算。只有带有*FILM关键字对流定义支配炉子与流体之间的交换。使用*ELEMENT卡定义了构成炉壁的6节点壳单元。此外，流体节点也被分配到单元（单元类型D），即所谓的网络单元。这些单元需要将材料属性分配给流体。实际上，传统上材料属性分配给单元而不是节点。每个网络单元由两个端节点（温度未知）和一个中节点组成，中节点用于定义通过单元的质量流量。流体节点603到613分配到网络单元301到305。

接下来定义两个节点集合：GAS包含所有流体节点，Ndown包含炉底的所有节点。

*PHYSICAL CONSTANTS卡在辐射起作用的分析中是必需的。它定义绝对零度（此处为0，因为我们使用开尔文）和斯特藩-玻尔兹曼常数。在本输入文件中，全程使用SI单位。

接着定义了钢（STEEL）的材料常数。热分析必须定义导热率、比热和密度。*SHELL SECTION卡将STEEL材料分配给之前*ELEMENT语句定义的单元集合FURNACE。它包含属于炉子的所有单元。此外，还指定了0.01 m的厚度。

材料GAS的材料常数由密度和比热组成。这些是流体的常数。不考虑流体中的传导。材料GAS分配给包含所有网络单元的单元集EGAS。

*INITIAL CONDITIONS卡为所有节点（炉子节点和流体节点）定义了300 K的初始温度。*AMPLITUDE卡定义了从0.0的0.3线性增加到1.0的1.0的斜坡函数。它将用于定义炉底的温度边界条件。这结束了模型定义。

第一步描述了温度边界条件从到的线性增加。*STEP卡上的INC=100参数允许此步有100个增量。过程是*HEAT TRANSFER，即我们希望执行纯热分析：唯一未知量是温度，没有机械未知量（如位移）。步时间为1，初始增量大小为0.1。两者都出现在*HEAT TRANSFER卡下的行上。*HEAT TRANSFER卡上缺少STEADY STATE参数表明这是瞬态分析。

接下来是温度边界条件：炉底板保持在1000 K，但由幅值A1调制。结果是温度边界条件从0.3 x 1000 = 300K开始，线性增加到t=1 s时的1000 K。第二个边界条件指定（流体）节点603的温度保持在300 K。这是入口温度。注意"11"是温度自由度。

流体中的质量流量通过应用于网络单元中节点中节点第一自由度的*BOUNDARY卡定义。第一行告诉我们（流体）节点609的质量流量为0.001。节点609是网络单元301的中节点。由于该速率为正，流体从节点603流向节点604，即从网络单元301的第一节点到第三节点。用户必须确保质量守恒（程序也会对此进行检查）。

第一组辐射边界条件指定炉底顶部表面以腔辐射形式辐射，发射率为1，环境温度为1000 K。对于腔辐射，当某位置处的视角因子不等于1时，使用环境温度。不足1的部分向环境辐射。每行中的第一个数字是单元，标签中的数字（每行第二个条目）是暴露于辐射的单元面。通常，这些行在cgx（CalculiX GraphiX）中自动生成。

第二和第三块定义了炉内顶部和侧面的内部腔辐射。第四块定义了炉顶板顶部向环境（保持300 K）的辐射。顶板的发射率为0.8。

接下来是膜条件。对炉顶板顶部表面定义了受迫对流，对流系数为 $h=25 \mathrm{W/mK}$ 。*FILM关键字下的第一行表明单元51的第二面通过受迫对流与（流体）节点604相互作用。该行中的最后一个条目是对流系数。因此，对于每个与流体相互作用的面，必须指定适当的流体节点进行相互作用。

最后，*NODE FILE卡确保温度存储在.frd文件中，*NODE PRINT卡确保流体温度存储在.dat文件中。

完整输入文件是CalculiX测试示例的一部分（furnace.inp）。对于当前分析，添加了第二步，将底部温度额外保持3000秒。

计算过程中会发生什么？炉壁和炉顶由于壁中的传导和来自底部的辐射而加热。然而，炉顶也通过与环境辐射和对流与流体失去热量。由于与流体的相互作用，温度不对称：在入口处流体较冷，炉子将比出口处失去更多热量，而出口处流体温度较高，与炉子的温差较小。因此，由于对流，我们期望温度从入口到出口增加。由于传导，我们期望在顶部中间出现温度最小值。这两种效应叠加。图15显示了t=3001 s时的温度分布。从炉底到炉顶存在温度梯度。在顶部温度确实不对称。图16也显示了这一点，其中绘制了位置A、B、C、D和E的温度随时间的变化。

注意，尚未达到稳态。还请注意，2D单元（如壳单元）自动扩展为具有正确厚度的3D单元。因此，从CalculiX GraphiX绘制的图片显示的是3D单元。