非稳态热传导

发表于 更新于 分类于

本文参考杨世铭、陶文铨所著《传热学》第四版第三章。物体的温度随时间而变化的导热过程称为非稳态导热(unsteady heat conduction)。根据物体温度随时间的推移而变化的特性，非稳态导热分为两类：(1) 物体的温度随时间的推移逐渐趋近于恒定值和 (2) 物体的温度随时间而做周期性的变化。在周期性的非稳态导热过程中，物体中各点的温度及热流密度都随时间做周期性的变化。例如，由于太阳辅射的周期性变化而引起的房屋的墙壁、屋顶等的温度场随时间的变化（以 24 h 为周期），地球表面层的温度由于季节更替而引起的周期性变化（以 1 年为周期）等等。

非稳态导热基本概念

非稳态导热特点

工程上几种典型的非稳态导热过程温度变化率的数量级如下图所示。在该图极高速非稳态导热区域（例如短脉冲、高强度激光处理）应当考虑非傅里叶导热的影响。

几种典型的非稳态导热过程温度变化率

非稳态导热过程中在热量传递方向上不同位置处的导热量是处处不同的，不同位置间导热量的差别用于（或来自）该两个位置间的物体内能随时间的变化，这是区别于稳态导热的一个特点。因此，对非稳态导热一般不能用热阻的方法来作问题的定量分析。

导热微分方程解的唯一性定律

假定物体热物性参数均为常数，则导热微分方程可写作 \[ \rho C_p \frac{\partial T}{\partial t} = \lambda \nabla^2 T + \dot{\Phi}, \] 其中\(~\nabla^2~\)是拉普拉斯（Laplace）算子。如果用热扩散率\(~a=\lambda/(\rho C_p)~\)进行变量替换，于是有 \[ \begin{align}\label{ht1} \frac{\partial T}{\partial t} = a \nabla^2 T + \frac{\dot{\Phi}}{\rho C_p}. \end{align} \] 一个实用上经常遇到的简单特例是初始温度均匀，即\(~T(x,y,z;0) = T_0\)。鉴于第三类边界条件比较常见，本文将着重讨论物体处于恒温介质中的第三类边界条件的非稳态导热，即： \[ -\lambda \left( \frac{\partial T}{\partial n} \right)_\text{w} = h(T_\text{w}-T_\text{f}). \] 如果某一函数\(~T(x,y,z;t)~\)同时满足以上微分方程、边界条件以及一定初始值，则该函数就是这一特定导热问题的唯一解。该结论称为解的唯一性定律。

第三类边界条件下\(~Bi~\)数对平板中温度分布的影响

为了说明第三类边界条件下非稳态导热时物体中的温度变化特性与边界条件参数的关系，下面来分析以下简单情形。设有一块厚为\(~2\delta~\)的金属平板，初始温度为\(~T_0\)，突然将它置于温度为\(~T_\infty~\)的流体中进行冷却，表面传热系数为\(~h\)，平板的导热系数为\(~\lambda\)。根据平板导热热阻\(~\delta/\lambda~\)与表面对流热阻\(~1/h~\)的相对大小，平板中温度场变换有以下三种情况

Bi~数对平板中温度分布的影响
  1. \(1/h \ll \delta/\lambda\)

这时，由于表面对流换热热阻\(~1/h~\)几乎可以忽略，因而过程一开始平板的表面温度就被冷却到\(~T_\infty\)。随着时间的推移，平板内部各点的温度逐渐下降而趋近于\(~T_\infty\)，如上图（a）。

  1. \(1/h \gg \delta/\lambda\)

这时，平板内部导热热阻\(~\delta/\lambda~\)几乎可以忽略，因而任一时刻平板中各点的温度接近均匀，并随着时间的推移整体地下降，逐渐趋近于\(~T_\infty\)，如上图（b）。

  1. \(1/h \sim \delta/\lambda\)

这时，平板中不同时刻的温度分布介于上述两种极端情况之间，如上图（c）。

由此可见，上述两个热阻的相对大小对于物体中非稳态导热的温度场的变化具有重要影响。我们知道，表征这两个热阻比值的量纲为一的量（习惯上称为无量纲量或无量纲数）就是毕渥（Biot）数 \[ Bi = \frac{\delta/\lambda}{1/h} = \frac{\delta h}{\lambda}. \] 像毕渥数、雷诺数这一类表征某一类物理现象或物理过程特征的无量纲数称为特征数（characteristic number），习惯上又称准则数。出现在特征数定义式中的几何尺度称为特征长度（characteristic length）,一般用符号\(~l~\)表示。在这里以平板的半厚作为特征长度，即取\(~l=δ\)。

零维问题的分析法——集总参数法

当固体内部的导热热阻远小于其表面的换热热阻时，任何时刻固体内部的温度都趋于一致，以致可以认为整个固体在同一瞬间均处于同一温度下。这时所要求解的温度仅是时间\(~t~\)的一元函数而与空间坐标无关，好像该固体原来连续分布的质量与热容量汇总到一点上，而只有一个温度值那样。这种忽略物体内部导热热阻的简化分析方法称为集中参数法（lumped parameter method）。显然，如果物体的导热系数相当大，或者几何尺寸很小，或表面传热系数极低，则其非稳态导热都可能属于这一类型的问题。例如，测量变化着的温度的热电偶就是个典型的实例。

集总参数法温度场的分析解

设有一任意形状的固体，其体积为\(~V\)，表面积为\(~A\)，并具有均匀的初始温度\(~t_0\)。在初始时刻，突然将它置于温度恒为\(~T_\infty~\)的流体中，设\(~T_0>T_\infty\)，固体与流体间的表面传热系数\(~h~\)及固体的物性参数均保持常数，试求物体温度随时间的依变关系。此问题可应用集总参数法分析。

由于物体的内部热阻可以忽略，温度与空间坐标无关，所以式\(~\eqref{ht1}~\)可简化为 \[ \frac{\mathrm{d}T}{\mathrm{d}t} = \frac{\dot{\Phi}}{\rho C_p}. \] 其中\(~\dot{\Phi}~\)是广义热源，界面上的热量边界条件为 \[ -\dot{\Phi}V = hA(T-T_\infty). \] 因为物体被冷却，\(T>T_\infty\)，故\(~\dot{\Phi}~\)为负值。联立以上两式，则有 \[ \rho C_p V\frac{\mathrm{d}T}{\mathrm{d}t} = -hA(T-T_\infty). \] 引入过余温度\(~\theta = T - T_\infty\)，则上式变成 \[ \begin{align}\label{eq1} \rho C_p V\frac{\mathrm{d}\theta}{\mathrm{d}t} = -hA\theta,\ \end{align} \] 初始条件为\(~\theta(0) = T_0 - T_\infty = \theta_0\)。将式\(~\eqref{eq1}~\)分离变量，得 \[ \frac{\mathrm{d}\theta}{\theta} = - \frac{hA}{\rho C_p V}\mathrm{d}t. \] 对上式进行积分 \[ \int_{\theta_0}^\theta \frac{\mathrm{d}\theta}{\theta} = - \int_0^t \frac{hA}{\rho C_p V}\mathrm{d}t, \] 求解得 \[ \begin{align}\label{eq2} \frac{\theta}{\theta_0} = \frac{T - T_\infty}{T_0 - T_\infty} = \exp\left(- \frac{hA}{\rho C_p V}t\right). \end{align} \] 令\(~l_\text{c} = V/A\)，则 \[ \frac{hA}{\rho C_p V}t = \frac{hl_\text{c}}{\lambda} \frac{\lambda}{\rho C_p}\frac{t}{l_\text{c}^2} = \frac{hl_\text{c}}{\lambda} \frac{at}{l_\text{c}^2} = Bi\cdot Fo, \] 其中\(~Fo~\)是傅里叶数。式\(~\eqref{eq2}~\)可写作
\[ \frac{\theta}{\theta_0} = \exp(-Bi\cdot Fo). \]

导热量计算式、时间常数与傅里叶数

  1. 导热量计算式

采用集总参数法分析时，从初始时刻到某一瞬间为止的时间间隔内，导热物体与流体间所交换的热量可由瞬时热流量对时间做积分而得。导热物体的瞬时热流量为 \[ \Phi = -\rho C_p V \frac{\mathrm{d}T}{\mathrm{d}t} = (T_0-T_\infty)hA\exp\left(- \frac{hA}{\rho C_p V}t\right) \] 从\(~t=0~\)时刻到\(~t~\)时刻之间所交换的总热量为 \[ Q_t = \int_0^t \Phi \mathrm{d}t = (T_0-T_\infty)\rho C_p V\left[1-\exp\left(- \frac{hA}{\rho C_p V}t\right)\right]. \]

  1. 时间常数

当采用集总参数法分析时，物体中的过余温度随时间成指数曲线关系变化。在过程的开始阶段温度变化很快，随后逐渐减慢，如下图所示

过余温度变化

如果\(~t=\rho C_p V/(hA)\)，则有 \[ \frac{\theta}{\theta_0} = \frac{T - T_\infty}{T_0 - T_\infty} = \exp(-1) = 36.8\%. \] \(\rho C_p V/(hA)~\)称为时间常数（time constant），记为\(~t_\text{c}\)。当时间\(~t=t_\text{c}~\)时，物体的过余温度已经降低到了初始过余温度值的 36.8%。在用热电偶测定流体温度的场合，热电偶的时间常数是说明热电偶对流体温度变动响应快慢的指标。显然，时间常数越小，热电偶越能迅速反映出流体温度的变动。时间常数不仅取决于热电偶的几何参数\(~V/A\)、物理性质\(~\rho\)、\(C_p\)，还同换热条件（\(h\)）有关。从物理意义上来说，热电偶对流体温度变化反应的快慢取决于其自身的热容量\(~\rho C_p V~\)及表面换热条件（\(hA\)）。热容量越大，温度变化得越慢；表面换热条件越好（\(hA~\)越大），单位时间内传递的热量越多，则越能使热电偶的温度迅速接近被测流体的温度。\(\rho C_p V~\)与\(~hA~\)的比值反映了这两种影响的综合结果。

  1. 傅里叶数的物理意义

\(Bi~\)数具有固体内部单位导热面积上的导热热阻与单位表面积上的换热热阻（即外部热阻）之比的意义\(~Bi=(l/\lambda)/(1/h)\)。\(Bi~\)数越小，意味着内热阻越小或外热阻越大，这时采用集中参数法分析的结果就越接近实际情况。例如，对于用热电偶测定流体温度的场合，\(Bi~\)数的值大概只有 0.001（或更小）这样一个数量级。试验证实，这时式\(\eqref{eq2}~\)同实测结果符合得很好。

傅里叶数的物理意义可以理解为两个时间间隔相除所得的无量纲时间，即\(~Fo=t/(l_\text{c}^2/a)\)，分子\(~t~\)是从边界上开始发生热扰动的时刻起到所计算时刻为止的时间间隔，分母\(~l_\text{c}^2/a~\)可以视为使边界上发生的有限大小的热扰动穿过一定厚度的固体层扩散到的面积上所需的时间。因此，\(Fo~\)数可以看成是表征非稳态过程进行深度的无量纲时间。在非稳态导热过程中，这一无量纲时间越大，热扰动就越深人地传播到物体内部，因而物体内各点的温度越接近周围介质的温度。

集总参数法的适用范围及应用举例

对于平板、圆柱与球中的一维非稳态第三类边界条件下的导热问题，当按特征长度 \[ \left . \begin{aligned} l &= \delta,\text{厚度为}~2\delta~的平板 \\ l &= R,\text{圆柱} \\ l &= R, 球 \end{aligned} \right \}. \] \(Bi~\)满足 \[ Bi = \frac{hl}{\lambda} \leq 0.1 \] 时，物体中最大与最小过余温度之差小于 5%，对于工程计算，此时已经足够精确地可以认为整个物体温度均匀。对于平板、圆柱与球，\(Bi~\)应该分别小于 0.1、0.05 和 0.033。

典型一维物体非稳态导热的分析解

三种几何形状物体的温度场分析解

  1. 平板
平板第三类边界条件

设有一块厚为\(~2\delta~\)的无限大平板，初始温度为\(~T_0\)。在初始瞬间将它放置于温度为\(~T_\infty~\)的流体中，设平板两边对称受热，板内温度分布必以其中心截面为对称面。因此，只要研究厚为\(~δ~\)的半块平板的情况即可。将\(~x~\)轴的原点置于板的中心截面上，对于\(~x\geq 0~\)地半平板，可以导出如下微分方程和边界条件 \[ \begin{aligned} \frac{\partial T}{\partial t} = a \frac{\partial^2 T}{\partial x^2},\ (0<x<\delta,t>0) \\ h[T(\delta,t)-T_\infty] = \left .\lambda\frac{\partial T(x,t)}{\partial x}\right|_{x=\delta}. \end{aligned} \] 引入过余温度\(~\theta = T - T_\infty\)，采用分离变量法求解得 \[ \frac{\theta(\eta,t)}{\theta_0} = \sum_{n=1}^\infty C_n\exp(-\mu_n^2 Fo)\cos(\mu_n \eta), \] 其中\(~Fo = at/\delta^2,~\eta = x/\delta\)。系数\(~C_n~\)满足 \[ C_n = \frac{2\sin \mu_n}{\mu_n + \cos\mu_n\sin\mu_n}. \] \(\mu_n~\)是下列超越方程得根，称为特征值（eigenvalue）： \[ \tan \mu_n = \frac{Bi}{\mu_n},\ n = 1,2,\cdots \] 其中\(~Bi=h\delta/\lambda\)。

  1. 圆柱

设半径为\(~R~\)的一实心圆柱，初始温度为\(~T_0\)，在初始瞬间将它置于温度为\(~T_\infty~\)的流体中，流体与圆柱表面间的表面传热系数\(~h~\)为常数。圆柱中无量纲温度的分析解如下： \[ \frac{\theta(\eta,t)}{\theta_0} = \sum_{n=1}^\infty C_n\exp(-\mu_n^2 Fo) J_0(\mu_n \eta), \] 其中\(~Fo = at/\delta^2,~\eta = r/R\)。系数\(~C_n~\)满足 \[ C_n = \frac{2}{\mu_n}\frac{J_1 (\mu_n)}{J_0^2(\mu_n)+J_1^2(\mu_n)}. \] \(\mu_n~\)是下列超越方程得根： \[ \mu_n \frac{J_1(\mu_n)}{J_0(\mu_n)} = Bi,\ n=1,2,\cdots \] 其中\(~Bi=hR/\lambda\)。

  1. 球

设半径为\(~R~\)的一实心球，初始温度为\(~T_0\) ，在初始瞬间将它置于温度为\(~T_\infty~\)的流体中，流体与球表面间的表面传热系数\(~h~\)为常数。球中无量纲温度的分析解如下： \[ \frac{\theta(\eta,t)}{\theta_0} = \sum_{n=1}^\infty C_n\exp(-\mu_n^2 Fo) \frac{\sin(\mu_n \eta)}{\mu_n\eta}, \] 其中\(~Fo = at/\delta^2,~\eta = r/R\)。系数\(~C_n~\)满足 \[ C_n = 2\frac{\sin(\mu_n)-\mu_n\cos(\mu_n)}{\mu_n-\sin(\mu_n)\cos(\mu_n)}. \] \(\mu_n~\)是下列超越方程得根： \[ 1-\mu_n\cos(\mu_n) = Bi,\ n=1,2,\cdots \] 可见，平板、圆柱和球中无量纲过余温度\(~\theta/\theta_0~\)与\(~Fo~\)数、\(Bi~\)数及无量纲距离\(~\eta~\)有关，即： \[ \frac{\theta(\eta,t)}{\theta_0} = \frac{T(\eta,t)-T_\infty}{T_0-T_\infty} = f(Fo,Bi,\eta). \]

非稳态导热正规状况阶段分析解的简化

  1. 非稳态导热正规状况阶段的物理概念与数学含义

非周期性的非稳态导热过程在进行到一定深度后，初始条件对物体中无量纲温度分布的影响基本消失,温度分布主要取决于边界条件的影响。非稳态导热的这一阶段称为正规状况阶段。

现在从分析解的数学表达式来揭示正规状况阶段的数学含义。三个解的特征值都是\(~Bi~\)数的函数，在一定的\(~Bi~\)数下\(~\mu_n~\)之值随\(~n~\)的增加而迅速增加。例如，对于平板，在\(~Bi=1.0~\)时前 4 个根分别为 0.8603、3.4256、6. 4373 和 9.5293。由三个分析解中反映时间影响的部分\(~\exp(-\mu_n^2 Fo)~\)可见，无穷级数第一项以后的各项，会随着\(~Fo~\)数的增加而迅速衰减。数值计算表明，当\(~Fo~\)数大于 0.2 以后，略去无穷级数中第二项及以后各项所得的计算结果与按完整级数计算结果的偏差小于 1%。这相当于将无穷级数解中的系数\(~C_n\ (n\geq 2)~\)取为零。因为\(~C_n~\)的无穷系列值是为了使分析解满足初始条件而引人的，这样的处理就意味着初始条件的影响已经消失，所以三个分析解无穷级数的第一项就是正规状况阶段温度场的解。对于非周期性的非稳态导热过程，从过程的开始到温度分布趋近于稳态分布的时间间隔中，初始条件影响基本消失的阶段占了极大部分的比例，故称这一阶段为“正规状况”。

  1. 正规状况阶段三个分析解的简化表达式

\[ \begin{align} &\text{平板}\quad \frac{\theta(\eta,t)}{\theta_0} = \frac{2\sin \mu_1}{\mu_1 + \cos\mu_1\sin\mu_1}\exp(-\mu_1^2 Fo)\cos(\mu_1 \eta), \label{eq3} \\ &\text{圆柱}\quad \frac{\theta(\eta,t)}{\theta_0} = \frac{2}{\mu_1}\frac{J_1 (\mu_1)}{J_0^2(\mu_1)+J_1^2(\mu_1)}\exp(-\mu_1^2 Fo) J_0(\mu_1 \eta),\\ &\text{球}\quad \frac{\theta(\eta,t)}{\theta_0} = 2\frac{\sin(\mu_1)-\mu_1\cos(\mu_1)}{\mu_1-\sin(\mu_1)\cos(\mu_1)}\exp(-\mu_1^2 Fo) \frac{\sin(\mu_1 \eta)}{\mu_1\eta}. \end{align} \]

以平板的解为例，正规状况阶段的任何时刻，平板中任意处（\(\eta\)）与平板中心（\(\eta=0\)）处的过余温度之比为 \[ \begin{align} \frac{\theta(\eta,t)}{\theta(0,t)} = \cos(\mu_1\eta),\label{eq4} \end{align} \] 可见这一比值与时间无关，只取决于特征值\(~\mu_1\)，即取决于边界条件，这是与“正规状况”4 个字的含义相一致的。

  1. 一段时间间隔内所传导的热量计算式

从初始时刻到平板与周围介质处于热平衡这一过程中所传递的热量为 \[ Q_0 = \rho C_p V(T_0-T_\infty). \] 这是非稳态导热过程中所能传递的最大热量，从初始时刻到某一时刻\(~t~\)这一阶段中所传递的热量\(~Q~\)与\(~Q_0~\)之比为 \[ \begin{aligned} \frac{Q}{Q_0} &= \frac{\rho C_p\int_V [T_0-T(x,t)]\mathrm{d} V}{\rho C_p V(T_0 - T_\infty)} \\ &= \frac{1}{V}\int_V\frac{(T_0-T_\infty)-(T-T_\infty)}{T_0-T_\infty}\mathrm{d}V \\ &= 1-\frac{1}{V}\int_V\frac{T-T_\infty}{T_0-T_\infty}\mathrm{d}V = 1-\frac{\bar{\theta}}{\theta_0}. \end{aligned} \] 平板、圆柱和球的正规状况阶段温度场与导热量的计算式可以统一表示为： \[ \begin{aligned} \frac{\theta(\eta,t)}{\theta_0} = A\exp(-\mu_1^2 Fo)f(\mu_1 \eta), \\ \frac{Q}{Q_0} = 1 - A\exp(-\mu_1^2 Fo)B. \end{aligned} \] 三种形状物体的\(~A\)、\(B\) 、\(f(\mu_1 η)~\)的表达式如下表

~A、B 、f(\mu_1 η)~的表达式

非稳态导热正规状况阶段的工程计算方法

利用上述公式计算时，需要根据不同的\(~Bi~\)数查出相应的特征值，并且要涉及贝塞尔（Bessel）函数等的计算，不甚方便。在传热学的发展史上先后提出了两种简化计算手续的方法，即由海斯勒（Heisler）等提出的诺谟图（nomogram） 方法；及由 Campo 提出的近似拟合公式方法，现分别介绍如下。

  1. 图线法

历史上曾广泛采用按分析解的级数第一项而绘制的用以确定温度分布的曲线,称为海斯勒图。以无限大平板为例，它首先根据式\(~\eqref{eq3}~\)给出\(~\theta_\text{m}/\theta_0~\)随\(~Fo~\)及\(~Bi~\)变化的曲线（此时\(~x/\delta=0\)），随后再根据式\(~\eqref{eq4}~\)确定\(~\theta/\theta_\text{m}~\)的值。于是平板中任意一点值为 \[ \frac{\theta}{\theta_0} = \frac{\theta_\text{m}}{\theta_0}\frac{\theta}{\theta_\text{m}}. \]

  1. 近似拟合公式法

三种几何形状的第一特征值\(~μ_1\)，以及\(~A\)、\(B~\)和零阶贝塞尔函数\(~J_0(x)~\)通过以下拟合公式：

\[ \begin{aligned} \mu_1^2 &= \left ( a + \frac{b}{Bi} \right)^{-1}, \\ A &= a + b(1-d^{-eBi}), \\ B &= \frac{a + cBi}{1 + bBi}, \\ J_0 &= a + bx + cx^2 + dx^3. \end{aligned} \]

上式中常数
计算~J_0(x)~的常数

上述分析解无论对物体被加热或冷却都是适用的。对于一维平板，在热边界条件方面还可以应用于以下两种情形：(1) 平板一侧绝热，另一侧为第三类边界条件；(2) 平板两侧面均为第一类边界条件且维持在相同的温度。

半无限大物体的非稳态导热

半无限大物体（semi-infinite body）可以看成是一维平板的一种特殊情况。所谓半无限大物体，是指几何上如下图所示 的那样的物体，其特点是从\(~x = 0~\)的界面开始可以向正向以及，上、下方向上无限延伸，而在每一个与\(~x~\)坐标垂直的截面，上物体的温度都相等。现实世界中不存在这样的半无限大物体，但是在研究物体中非稳态导热的初始阶段时，则有可能把实际物体当作半无限大的物体来处理。例如，假设有一块几何上为有限厚度的平板，起初具有均匀的温度，然后其一侧表面突然受到热扰动，如壁温突然升高到一定值并保持不变，或者突然受到恒定的热流密度加热，或者受到温度恒定的流体的加热或冷却。当扰动的影响还局限在表面附近而尚未深入到平板内部中去时，就可有条件地把该平板视为一“半无限大物体”。工程导热问题中有不少情形可按半无限大物体处理。

半无限大物体

三种边界条件下半无限大物体温度场的分析解

有一半无限大物体，初始温度均匀为\(~T_0\)。在\(~t = 0~\)时刻，\(x=0~\)的侧面突然受到热扰动，这种情况可以归纳为以下三种边界条件：(1) 表面温度突然变化到\(~T_\text{w}\)，并保持恒定（第一类）；(2) 受到恒定的热流密度加热（第二类）；(3) 与温度为\(~T_\infty~\)的流体进行热交换（第三类）。这三类边界条件如下图所示。假定物体的热物性为常数，没有内热源。

三种边界条件

上述条件下物体中温度的控制方程和定解条件为 \[ \begin{aligned} \frac{\partial T}{\partial t} &= a \frac{\partial^2T}{\partial x^2}, \ 0<x<\infty; \\ t &= 0,\ T(x,t) = T_0; \\ x &= 0. \end{aligned} \] 温度场的分析解为

第一类边界条件 \[ \frac{T(x,t)-T_\text{w}}{T_0-T_\text{w}}=\mathrm{erf}\left(\frac{x}{2\sqrt{at}}\right), \] 第二类边界条件 \[ T(x,t)-T_0 = \frac{2q_0\sqrt{at/\pi}}{\lambda}\exp\left(-\frac{x^2}{4at}\right)-\frac{q_0 x}{\lambda}\mathrm{erfc}\left(\frac{x}{2\sqrt{at}}\right), \] 第三类边界条件 \[ \begin{align}\label{eq5} \frac{T(x,t)-T_0}{T_\infty-T_0}=\mathrm{erf}\left(\frac{x}{2\sqrt{at}}\right)-\exp\left(\frac{hx}{\lambda}+\frac{h^2at}{\lambda^2}\right)\mathrm{erfc}\left(\frac{x}{2\sqrt{at}}+\frac{h\sqrt{at}}{\lambda}\right) \end{align} \] 其中\(~\mathrm{erf}(x)~\)称为误差函数，\(\mathrm{erfc}(x)=1-\mathrm{erf}(x)~\)称为余误差函数。

导热量计算式

下面以上述第一种边界条件为例，来导出从初始时刻到某一指定时刻\(~t~\)之间，即在时间间隔\(~[0,t]~\)内半无限大物体表面与外界的换热量（亦即半无限大物体内的导热量）。

通过任意截面\(~x~\)处的热流密度 \[ \begin{aligned} q_x &= -\lambda\frac{\partial T}{\partial x} \\ &= -\lambda(T_0-T_\text{w})\frac{\partial\mathrm{erf}(\eta)}{\partial x} \\ &= \lambda\frac{T_\text{w}-T_0}{\sqrt{\pi at}}\exp\left(\frac{x^2}{4at}\right), \end{aligned} \] 在表面上的导热量为 \[ \begin{aligned} Q &= A\int_0^t q_\text{w}\mathrm{d}t \\ &= A\int_0^t\frac{\lambda(T_\text{w}-T_0)}{\sqrt{\pi at}}\mathrm{d}t \\ &= 2A\sqrt{t/\pi}\sqrt{\rho C_p \lambda}(T_\text{w}-T_0). \end{aligned} \] 由以上两式可见，表面上的瞬时热流密度与时间的平方根成反比，而总的导热量则与时间的平方根成正比。此外，\(Q~\)还与物体的\(~\sqrt{\rho C_p \lambda}~\)成正比（注意：在稳态导热中，导热量与\(~\rho C_p~\)无关，而只与\(~\lambda~\)成正比）。在材料成形工业中称\(~\sqrt{\rho C_p \lambda}~\)为吸热系数，它的大小代表了物体向与其接触的高温物体吸热的能力。在选择造型材料与冷铁时，吸热系数是一个重要指标，它关系到物体（如铸件）的冷却速度。

简单几何形状物体多维非稳态导热的热分解

在多维导热问题中，几种简单儿何形状物体的非稳态导热问题的分析解，可以用几个相应的一维非稳态导热分析解相乘得出，称为乘积解法（production solution method）。方形截面的二维柱体、短圆柱体以及立方体是这类几何形状物体的典型例子。本节中先介绍获得无量纲温度场的乘积解法的基本思想，并对二维问题给出数学证明，以加深对乘积解法应用条件的理解。然后介绍确定导热量的方法。通过计算例题还可以进一步体会到采用拟合公式方法计算的优点。

获得无量纲温度场的乘积解法

如下图所示的三种固体。其中：图 (a) 所示为方形截面的二维导热物体，沿着\(~z~\)方向物体温度没有变化；图 (b) 所示为短圆柱体，物体温度在半径及轴向发生变化；图 (e) 所示是立方体，温度在三个坐标方向均改变。这三种固体在几何上可以分别看成是由两块平板、一块平板与一个圆柱以及三块平板相贯而成,如图 (d)、(e) 和 (f) 所示。假设三个物体的初始温度都是均匀的，记为\(~T_0\)，然后与周围介质之间发生对流传热，流体温度为\(~T_\infty\)，表面传热系数为\(~h\)。在这样的条件下，三个固体中的无量纲温度场可以由其几何上的相贯体的一维分析解相乘而得：

乘积解法 \[ \begin{aligned} &\text{二维柱体}\quad \Theta = \frac{\theta(x,y,t)}{\theta_0} = \Theta_{\rho 1}(x,t)\cdot \Theta_{\rho 2}(y,t), \\ &\text{短圆柱}\quad \Theta = \frac{\theta(x,r,t)}{\theta_0} = \Theta_{\rho 1}(x,t)\cdot \Theta_{o}(r,t), \\ &\text{立方体}\quad \Theta = \frac{\theta(x,y,z,t)}{\theta_0} = \Theta_{\rho 1}(x,t)\cdot \Theta_{\rho 2}(y,t)\cdot \Theta_{\rho 3}(z,t). \end{aligned} \] 其中\(~\Theta_\rho~\)、\(\Theta_o~\)分别表示一维平板及圆柱在第三类边界体条件下无量纲温度的分析解。

二维问题的证明

下面以无限长方柱体（即其截面为长方形的柱体）的非稳态导热问题为例来作分析。这是一个二维问题。截面尺寸为\(~2\delta_1\times 2\delta_2~\)的方柱体可以看成是两块厚度分别为\(~2\delta_1~\)及\(~2\delta_2~\)的无限大平板垂直相交所截出的物体。

方柱体

如上图所示，设方柱体初始温度为\(~T_0\)，开始时被置于温度为\(~T_\infty~\)的流体中，表面与流体间的表面传热系数为\(~h\)。上图阴影截面上的无量纲温度分布\(~\Theta(x,y,t)~\)由下列微分方程和定解条件决定： \[ \begin{aligned} &\frac{\partial \Theta}{\partial t} = a\left( \frac{\partial^2\Theta}{\partial x^2} + \frac{\partial^2\Theta}{\partial y^2} \right), \\ &\Theta(x,y,0) = 1, \\ &\Theta(\delta_1,y,t) + \frac{\lambda}{h}\left.\frac{\partial \Theta(x,y,t)}{\partial x}\right|_{x=\delta_1} = 0, \\ &\Theta(x,\delta_2,t) + \frac{\lambda}{h}\left.\frac{\partial \Theta(x,y,t)}{\partial y}\right|_{y=\delta_2} = 0, \\ &\left.\frac{\partial \Theta(x,y,t)}{\partial x}\right|_{x=0} = 0, \\ &\left.\frac{\partial \Theta(x,y,t)}{\partial y}\right|_{y=0} = 0. \end{aligned} \] 其中 \[ \Theta = \frac{T(x,y,t)-T_\infty}{T_0-T_\infty} = \frac{\theta}{\theta_0} \] 为无量纲过余温度。

这两块无限大平板分析解的乘积就是上述无限长方柱体的解，即 \[ \Theta(x,y,t) = \Theta_x(x,t)\Theta_y(y,t). \] 乘积解法只适用于：(1) 物体初始温度均匀；(2) 周围介质温度均匀；(3) 表面传热系数均匀；(4) 常物性、没有内热源。对于表面传热系数为无限大的情形，即第一类边界条件，当然也是适用的。

导热量的计算

为了获得非稳态导热过程从初始时刻到时刻\(~t~\)的导热量，可以首先计算该热量占非稳态导热过程总导热量的百分数，然后与从热平衡角度得出的总导热量相乘即可。所以,计算导热量的关键是要得出这个导热量的百分数。朗司顿（Langston）已证明，对于多维非稳态导热，导热量百分数也可以通过类似乘积解法的模式得出。假设\(~(Q/Q_0)_1\)、\((Q/Q_0)_2\)、\((Q/Q_0)_3~\)分别是构成一个二维与三维非稳态导热物体的一维几何体的导热量百分数，则可有：

二维问题 \[ \frac{Q}{Q_0}=\left(\frac{Q}{Q_0}\right)_1+\left(\frac{Q}{Q_0}\right)_2\left[1-\left(\frac{Q}{Q_0}\right)_1\right], \] 三维问题 \[ \begin{aligned} \frac{Q}{Q_0} &= \left(\frac{Q}{Q_0}\right)_1+\left(\frac{Q}{Q_0}\right)_2\left[1-\left(\frac{Q}{Q_0}\right)_1\right] \\ &+ \left(\frac{Q}{Q_0}\right)_3\left[1-\left(\frac{Q}{Q_0}\right)_1-\left(\frac{Q}{Q_0}\right)_2\right]. \end{aligned} \]

小结

当遇到一个工.程非稳态导热问题时，建议按照以下步骤进行求解：

  1. 如果物体的形状复杂而且需要获得物体中的温度随时间变化的详细的信息，例如要确定发动机转子和气缸壁中的温度在启动过程中的分布及变化，以计算金属中由于温度不均匀而引起的热应力，则应求助于数值计算的方法。
  2. 如果物体的形状相对比较简单，或者允许作一定的近似处理，则可考虑采用分析解的方法。对此类问题，首先要尽可能确切地决定非稳态导热物体的热边界条件。在许多情形下，导热物体的边界热作用可以由第三类边界条件近似描述，这时获取比较准确的表面传热系数的数值对于解的准确性具有重要意义。
  3. 对于在整个非稳态导热过程中表面传热系数有相当大变化的情形，作为计算的第一步可以取其平均值来计算。进一步的方法是将非稳态过程分为几个时段，每一个时段中表面传热系数取为常数。这样的处理往往带有迭代的性质。
  4. 对第三类边界条件的问题,可以采用以下解题策略：
    1. 首先尝试集中参数法，计算\(~Bi~\)数（其中可以用\(~V/F~\)作为特征尺度）。如果\(~Bi<0.1\)，则就采用集中参数法。
    2. 如果\(~Bi>0.1\)，则计算\(~Fo~\)数。如果\(~Fo<0.05\sim0.06\)，则可将导热物体看成是半无限大的物体，采用式\(~\eqref{eq5}~\)计算物体中的温度。
    3. 如果\(~Bi>0.1,0.06<Fo<0.2\)，则对可以作为一维问题处理的导热物体，需采用完全的级数解。注意，求解多维问题的乘积解法对于非稳态导热的初始阶段也是适用的。
    4. 如果\(~Bi>0.1,Fo>0.2\)，则可采用正规状况阶段的简化解法。