非稳态导热微分方程的分析解

掌握从非稳态过渡到稳态的感性认识

传热学-第三章-一维平壁非稳态导热的分析解

温度分布和热流量分布随时间和空间的变化规律的函数原型：

$t=f(x,y,z,\tau); \Phi=f(\tau)$

导热微分方程的一般形式(第二章推导过)

$\rho c_p \frac{\partial t}{\partial \tau}=\frac{\partial}{\partial x}(\lambda \frac{\partial t}{\partial x})+\frac{\partial}{\partial y}(\lambda \frac{\partial t}{\partial y})+\frac{\partial}{\partial z}(\lambda \frac{\partial t}{\partial z})+q_v\tag{1}$

求解方法有三种，分析解、数值解、实验解

一.下面为非稳态导热，有限厚度无限大平板非稳态导热的解析解（分离变量法）

1.物理模型为 $2\delta$ 厚度的平板， $\tau=0$ 初始温度分布为 $t_0(x)$ ，放入 $t_\infty$ 的流体中，求温度分布 $t(x,\tau)$

该物理模型简化为一维问题，非稳态，常物性（ $\lambda , a$ 为常数），两侧温度分布对称，中心为原点，无內热源，边界条件为在传热系数为 $h$ 的流体中的第三类边界条件,微分方程 $(1)$ 简化为：

$\frac{\partial t}{\partial \tau}=a\frac{\partial^2t}{\partial x^2}\tag{2}$

①时间条件为初始条件,该条件只有瞬态条件才有，稳态热传导没有 $\tau=0$ 的初始条件

②温度分布对称,中心温度最大或者最低，因此中心为绝热边界条件，外边界条件以第三类边界条件（举例）:

$\begin{cases} \frac{\partial t}{\partial \tau}=a\frac{\partial^2t}{\partial x^2}\\ -\lambda\frac{\partial t}{\partial x}|_{x=\delta}=h(t|_\delta-t_\infty)\\ \frac{\partial t}{\partial x}|_{x=0}=0\\ \tau=0, t(x,0)=t_0\\ \tag{3} \end{cases}$

③分析条件的齐次与非齐次问题，方程（2）为齐次， $x=\delta$ 处的边界条件为非齐次，x=0的边界条件为齐次，全部齐次化： $\theta(x,\tau)=t（x,\tau）-t_\infty$ ,因此 $t=\theta+t_\infty$ ,并且注意到由于 $t_\infty$ 为常数，在求一阶偏微分，二阶偏微分的时候不影响，其他方程仍为齐次( $\frac{\partial t}{\partial \tau}=\frac{\partial \theta}{\partial \tau}$ , $a\frac{\partial^2t}{\partial x^2}=a\frac{\partial^2\theta}{\partial x^2}$ )

方程组引入过余温度后， $(2)(3)$ 转为为如下方程组:

$\begin{cases} \frac{\partial \theta}{\partial \tau}=a\frac{\partial^2\theta}{\partial x^2}\\ -\lambda\frac{\partial \theta}{\partial x}|_{x=\delta}=h\theta\\ \frac{\partial \theta}{\partial x}|_{x=0}=0\\ \tau=0,\theta=\theta_0=t_0-t_\infty \tag{4} \end{cases}$

2.分离变量法,假设 $\theta=f(x，\tau)=X(x)\cdot\phi(\tau)$

并带入到 $(4)$ 中第一个式子
$X\frac{\partial \phi}{\partial \tau}=a\phi \cdot \frac{\partial^2X}{\partial x^2}\tag{5}$

$X(x)与\tau无关，\phi(\tau)与x无关$

偏微分退化为常微分，微分量可以自由移动，交换变量变化如下
$\frac{1}{a\phi}\frac{d \phi}{d \tau}=\frac{1}{X} \frac{d^2X}{d x^2}\tag{6}$

分析上述公式，左边是关于时间的函数，右边是关于空间的函数, 秒=米？

要让秒=米是不可能的，等式两边必须同等于一常数

$\frac{1}{a\phi}\frac{d \phi}{d \tau}=\frac{1}{X} \frac{d^2X}{d x^2}=D\tag{6}$

转换为两个独立方程，如下

$\begin{cases} \frac{1}{a\phi}\frac{d \phi}{d \tau}=D\\ \frac{1}{X} \frac{d^2X}{d x^2}=D\\ \end{cases} \tag{7}$

解常微分方程, $(7)的$ 第一个方程很好解, $\frac{1}{\phi}\frac{d \phi}{d \tau}=aD$ , $ln\phi=D+c_1$

$\phi=c_1 e^{aD\tau}\tag{8}$

分析以下上面的数据, $a$ 为热扩散率，大于0，时间恒为正，如果D为正数，那么 $\phi随时间增大会不会趋近0$ ，而实际情况是时间越长温度越接近流体温度，过余温度接近0，因此D一定为负数

设 $D=-\beta^2$ ，那么方程组(7)变化为如下：

$\begin{cases} \frac{1}{a\phi}\frac{d \phi}{d \tau}=-\beta^2\\ \frac{1}{X} \frac{d^2X}{d x^2}=-\beta^2\\ \end{cases} \tag{9}$

对于 $关于时间\tau的函数\phi$ 来说，解微分方程如下为

$\phi=c_1 e^{-a\beta^2\tau}\tag{10}$

对于关于x位置的函数，实际上为 $X^{''}+\beta^2X=0$ ,是二阶导数与函数值的线性组合，这个解在高数(下)中查到结果

$X(x)=c_2cos(\beta x)+c_3sin(\beta x)\tag{11}$

由于之前分离变量的假设为 $\theta(x,\tau)=X(x)\cdot\phi(\tau)$ ,那么我们把 $(10)\cdot (11)$ 就得到了过余温度关于 $x,\tau$ 函数的一般解形式:

$\theta(x,\tau)=[Acos(\beta x)+Bsin(\beta x)]\cdot e^{-a\beta^2\tau}\tag{12}$

其中未知系数c1,c2,c3的换为新的未知数 $A=c_1\cdot c_2,B=c_1\cdot c_3$

3.边界条件及初始条件确定未知系数得到定解

$A，B，\beta$ （β是我们引入的）三个未知数，三个定解条件。

方程组 $(4)$ 中有一个初始条件，二个定解条件。三个方程决定三个未知数，此方程有定解

先看x=0处的绝热边界条件 $x=0,\frac{\partial \theta}{\partial x}=0$ ,将通解形式(12)对x求偏导带入，如下：

$\frac{\partial \theta}{\partial x}=[-A\beta sin(\beta x)+B\beta cos(\beta x)]\cdot e^{-a\beta^2\tau}\tag{13}$

x=0带入，求得 $B\beta\cdot e^{-a\beta^2\tau}=0$ （指数函数永远不为0）,而 $\beta \neq 0$ ,那么只有B=0。

一般解简化为如下，1个定解条件解得一个B=0

$\theta(x,\tau)=Acos(\beta x)\cdot e^{-a\beta^2\tau}\tag{14}$
第二个定解条件，为第三类边界条件 $x=\delta，-\lambda\frac{\partial \theta}{\partial x}=h\theta|_{x=\delta}$ ,将 $x=\delta以及\theta的表达式（14）带入(13)$
$-\lambda\cdot [-A\beta sin(\beta \delta)]\cdot e^{-a\beta^2\tau}=h\cdot Acos(\beta \delta)\cdot e^{-a\beta^2\tau}\tag{15}$
上式中
$\require{cancel} -\lambda\cdot [-\cancel{A}\beta sin(\beta \delta)]\cdot \cancel {e^{-a\beta^2\tau}}=h\cdot \cancel{A}cos(\beta \delta)\cdot \cancel{e^{-a\beta^2\tau}}$
简化为 $\lambda\cdot\beta sin(\beta \delta)=h\cdot cos(\beta \delta)$
进一步
$tan(\beta \delta)=\frac{h}{\lambda\cdot \beta}\tag{16}$

我们学过了Bi数的定义，为内阻：外阻, $Bi=\frac{\delta}{\lambda}/\frac{1}{h}=\frac{h\delta}{\lambda}$ ,将等式（16）右边分子分母同时乘以 $\delta$ 凑出一个Bi数

（16）转换为 $tan(\beta \delta)=\frac{Bi}{\beta\delta}$ ,换元 $\mu=\beta\delta$

$tan\mu=\frac{Bi}{\mu}或者 cot\mu=\frac{\mu}{Bi}\tag{17}$

这个方程为超越方程，左边是tan x 三角函数，右边是1/x 函数，并且tan x为周期函数，两个曲线有无穷个交点，就有无穷个解.

000.png

那么，这个超越方程的解集合称之为特征方程，无穷多个特征解为 $\mu_1,\mu_2……\mu_n$

这些特征解的无穷级数求和，组成新的解形式,一般情况下用6-8项就够了，下面我们简单分析一下。

1.当 $Bi\to\infty,也就是函数tan\mu\to\infty，特征值\mu_n为\frac{1}{2}\pi,\frac{3}{2}\pi,\frac{5}{2}\pi……$

2.当 $Bi\to 0,tan\mu\to 0,，特征值\mu_n为0,\pi,2\pi……$

3.在给定Bi数的情况下，对应每一个特征值 $\mu_n=\beta_n\delta\to\beta_n=\frac{\mu_n}{\delta}$ ，温度分布的特解分别为如下：

$\begin{cases} \theta_1(x,\tau)=A_1cos(\beta_1x)\cdot e^{-a\beta_1^2\tau}\\ \theta_2(x,\tau)=A_2cos(\beta_2x)\cdot e^{-a\beta_2^2\tau}\\ ………\\ \theta_2(x,\tau)=A_ncos(\beta_nx)\cdot e^{-a\beta_n^2\tau}\\ \end{cases} \tag{18}$

反直觉之——有多个温度分布的解析式:

公式集合(18)为 $x=\delta，-\lambda\frac{\partial \theta}{\partial x}=h\theta|_{x=\delta}$ 及 $x=0,\frac{\partial \theta}{\partial x}=0$ 得到的特解，常数 $A_1,A_2……A_n$ 为任何值，都满足导热微分方程式和这关于x位置的两个定解条件，但是一个都不满足关于时间 $\tau=0,\theta=\theta_0$ 时刻的初始条件的温度分布——因为 $A_i$ 只要有一个不相等，就会有多个初始温度分度，产生悖论。

知识点补充：导热微分方程为线性方程 $\frac{\partial \theta}{\partial \tau}=a\frac{\partial^2\theta}{\partial x^2}$

温度对时间一阶导数的系数式1，温度对位置x的二阶导数为 $a$ ，都与温度无关

该导热问题的通解为各个特解的线性求和叠加：导热微分方程式和边界条件都是线性的——温度和温度的各阶导数项的系数都与温度无关。

例如，F(x)=a,F(y)=b, 叠加起来F(x+y)=a+b, 给一个信号x和y的叠加信号，等于单独给信号以后的结果之和

线性叠加的组合得到一个与初始温度分布相等的表达式：

$\theta(x,\tau)=\sum_{n=1}^{\infty}A_ncos(\beta_n x)e^{-a\beta^2_n \tau}\tag{19}$

公式（19）中 $\beta_n$ 为已知数， $A_n$ 为未知数，带入初始条件 $\tau=0,\theta=\theta_0=t_0-t_\infty$ 求解得到 $A_n$ 即可得到唯一特解

$\beta_n=\frac{\mu_n}{\delta}$

$\theta_0=\sum_{n=1}^{\infty}A_ncos(\beta_n x)=\sum_{n=1}^{\infty}A_ncos(\mu_n\frac{x}{\delta}) \tag{20}$

A_n还是不知道，级数方程没有复变的基础下不好解，下面还是利用特征函数的特性-正交性来解决,当下标不等（非对角元乘积的积分为0）

正交性问题来源于常微分方程的施图姆-刘维尔理论 $

这里不展开，二阶常微分方程为亥姆霍兹方程，利用该结论即可。施图姆-刘维尔理论 (Sturm–Liouville theory ) - 知乎 (zhihu.com)

002.png

需要利用的结论是

$当m\neq n$ 时， $\int_0^\delta A_ncos(\mu_n\frac{x}{\delta})\cdot cos(\mu_m \frac{x}{\delta})dx=0$

方程(20)两边同时乘以 $cos(\mu_m \frac{x}{\delta})$ ,并在 $0\leq x\leq \delta$ 范围内积分,这个m可以等于也可以不等于n

$\theta_0\cdot \int_0^\delta cos(\mu_m \frac{x}{\delta})d_x=\int_0^\delta\sum_{n=1}^{\infty}A_ncos(\mu_n\frac{x}{\delta})\cdot cos(\mu_m \frac{x}{\delta})dx \tag{21}$

那么公式中无穷求和只剩下 m等于n这一项，去掉级数符号

$\theta_0\cdot \int_0^\delta cos(\mu_n \frac{x}{\delta})d_x=\int_0^\delta A_ncos^2(\mu_n\frac{x}{\delta})dx\tag{22}$

这个积分就好解了，可以求得 $A_n$ 等于两个积分相除

$\require{cancel} A_n=\theta_0\frac{\int_0^\delta cos(\mu_n\frac{x}{\delta})dx}{\int_0^\delta cos^2(\mu_n\frac{x}{\delta})dx}=\theta_0\frac{\int_0^\delta cos(\mu_n\frac{x}{\delta})d(\mu_n\frac{x}{\delta})\cdot\frac{\delta}{\mu_n}}{\int_0^\delta\frac{1}{2}(1+cos(2\mu_n\frac{x}{\delta}))d(2\mu_n\frac{x}{\delta})\cdot \frac{\delta}{2\mu_n}}\\ A_n=\theta_0\frac{sin(\mu_n\frac{x}{\delta})|_0^\delta\cdot \frac{\delta}{\mu_n}}{\frac{1}{2}(2\mu_n\frac{x}{\delta}+sin(2\mu_n \frac{x}{\delta}))|_0^\delta \cdot \frac{\delta}{2\mu_n}}=\theta_0\frac{sin(\mu_n\frac{x}{\delta})|_0^\delta\cdot \cancel{\frac{\delta}{\mu_n}}}{\frac{1}{2}(\cancel{2}\mu_n\frac{x}{\delta}+\cancel{2}sin(\mu_n \frac{x}{\delta})cos(\mu_n \frac{x}{\delta}))|_0^\delta \cdot \cancel{\frac{\delta}{2\mu_n}}}=\theta_0\cdot\frac{2sin\mu_n}{\mu_n+sin\mu_n cos\mu_n}\tag{23}$

最后将 $A_n$ 的表达式带入公式（19）中

$\theta(x,\tau)=\sum_{n=1}^{\infty} A_n cos(\beta_nx)\cdot e^{-a\beta^2_n\tau}\tag{19}$

得到

$\theta(x,\tau)=\theta_0 \sum_{n=1}^{\infty}\frac{2sin\mu_n}{\mu_n+sin\mu_n cos\mu_n}cos(\beta_nx)\cdot e^{-a\beta^2_n\tau}\tag{24}$
其中 $\mu_n$ 为超越方程 $tan\mu=\frac{Bi}{\mu}$ 的系列特征根并且 $\beta_n=\frac{\mu_n}{\delta}$ ,上式（24）为无未知数的无穷级数组合

唯一决定解析解的公式已经推导完毕，我们进行一下变形然后进行对时间和位置变化的趋势行为的分析

已知毕沃准则数 $Bi=\frac{h\delta}{\lambda}$ ,为内阻与外阻之比；傅里叶准则数 $F_o=\frac{a\tau}{L^2}$ ，L为特征长度，在 $2\delta$ 的平板中, $L=\delta，F_o=\frac{a\tau}{\delta^2}$ ,为无量纲时间

将特解(24)转换为这两个准则数的自变量的形式,利用 $\mu=\beta\delta;\beta_n=\frac{\mu_n}{\delta}$

$\require{cancel} \theta(x,\tau)=\theta_0 \sum_{n=1}^{\infty}\frac{2sin\mu_n}{\mu_n+sin\mu_n cos\mu_n}cos(\cancel{\beta_n}\mu_n \frac{x}{\delta})\cdot e^{-a\cancel{\beta^2_n}\frac{\mu_n^2}{\delta^2}\tau}=\theta_0 \sum_{n=1}^{\infty}\frac{2sin\mu_n}{\mu_n+sin\mu_n cos\mu_n}cos(\mu_n\frac{x}{\delta})e^{-\mu_n^2 F_o}\tag{25}$

在（25）中， $\mu$ 为Bi数相关， $\frac{x}{\delta}$ 为位置相关， $F_0$ 为时间相关，因此该解析式包含了①内部导热与外部条件的因素②位置自变量③时间自变量这三个因变量

$\theta(x,\tau)=f(Bi,F_o,\frac{x}{\delta})$

先对傅里叶数进行分析，当 $F_0\geq 0.2$ 时（根据热扩散率和厚度估算一般为20秒），该指数函数衰减的速度超过了前面函数(三角函数)增加的速度，无穷级数在第1项就几乎收敛！那么我们用第一项 $\mu_1$ 就可以得到精度很高的近似解。

在上述特殊情况，可以用上式取级数①第一项进行计算，②可以用诺莫图计算，③也可以用campo近似拟合公式计算,将过余温度表达式与初始时刻的过于温度相除，并且分子分母同时乘以中心处任意时间的过余温度

$\frac{\theta(x,\tau)}{\theta_0}=\frac{\theta(x,\tau)}{\theta _m(\tau)}\cdot\frac{\theta _m(\tau)}{\theta_0}\\ f(Bi,\frac{x}{\delta})\cdot f(Bi,F_O)\tag{26}$

有了温度分布，根据 $Q=cm\Delta t$ 计算换热量

$Q\tau=\rho c_p\int_{-\delta}^{+\delta}(t_0-t)dx=2\rho c_p\delta\theta_0[1-\sum_{n=1}^{\infty}\frac{2sin^2\beta_n}{\beta^2_n+\beta_nsin\beta_ncos\beta_n}e^{-\beta_n^2F_o}] [J/m^2]\tag{27}$
其中持续到稳态放完的总热量为 $Q_0=2\delta\rho c_p(t_0-t_\infty)=2\rho c_p\delta\theta_0$ ,因此 $\tau$ 时刻下释放的热量与 $Q_0$ 比值: $\frac{Q_\tau}{Q_0}=f(F_o,Bi)$ ,为傅里叶数和毕渥数的函数，与位置无关

二、准则数对温度分布的影响

1. $F_o$ 数对温度分布的影响

$F_o=a\tau/\delta^2>0.2$ 时，无量纲温度表达式取无穷级数第一项,有一个根 $\mu_1$ 即可

$\theta(x,\tau)=\theta_0 \frac{2sin\mu_1}{\mu_1+sin\mu_1 cos\mu_1}cos(\mu_1\frac{x}{\delta})e^{-\mu_1^2 F_o}\tag{28}$
两边取对数,并注意到超越方程的特征根 $\mu_1$ 为Bi数和位置 $\frac{x}{\delta}$ 的函数：
$\require{cancel}\\ ln\theta=-(\beta_1^2\frac{a}{\delta^2})\tau+ln[\theta_0\frac{2sin\mu_1}{\mu_1+sin\mu_1 cos\mu_1}cos(\mu_1\frac{x}{\delta})]=-\cancel{(\beta_1^2\frac{a}{\delta^2})}m\tau+\cancel{ln[\theta_0\frac{2sin\mu_1}{\mu_1+sin\mu_1 cos\mu_1}cos(\mu_1\frac{x}{\delta})]}K(Bi,\frac{x}{\delta})$
简化为 $ln\theta=-m\tau+K(Bi,\frac{x}{\delta})\tag{29}$
给定边界条件确定了对流换热系数h，就给定了Bi数，以及位置给定，m为常数，与时间和位置都无关，称之为冷却率或加热率。所以采用第一项级数 $\mu_1$ 获得的温度与时间的变化曲线，不同位置下的斜率都相等。
假设我们给定中心位置，那么过余温度为 $\theta_m, 那么我们在一个确定的Bi数下，\frac{\theta_m}{\theta_0}对时间的斜率就是常数，参考课本P112页的诺莫图$

① $F_o>0.2$ ，过余温度的对数随时间变化的斜率为常数，与时间和位置无关，后称之为正规状况阶段：初始温度分布的影响已经消失

② $F_o<0.2$ 对应的非正规状况阶段:中心处收到初始温度的影响变化速度比表面上更慢。

除开大平壁以外其他模型也有类似的结论

非稳态导热分为

三个阶段：①非正规状况阶段；②正规状况阶段；③新的稳态

000.png

2.Bi准则数对温度分布的影响:

$Bi=\frac{h\delta}{\lambda}=\frac{\delta/\lambda}{1/h}=\frac{物体内部导热热阻}{物体表面对流换热热阻}$

无限大平板在冷却时，其第三类边界条件: $x=\delta,-\lambda \frac{\partial t}{\partial x}|_{x=\delta}=h(t|_{x=\delta}-t_\infty)$ ,将边界条件进行处理

回顾一下温度分布的曲线,中心处冷却的慢，边界上冷却的快

001.png

$-\frac{\partial t}{\partial x}|_{x=\delta}=\frac{t|_{x=\delta}-t_\infty}{\lambda/h}=\frac{t|_{x=\delta}-t_\infty}{\delta/Bi}\tag{30}$
将边界条件处理后，其物理意义为边界处的温度分布随x变化的斜率等于温差除以 $\delta/Bi$

000.png

上图为定性描述分析，将边界处的温度分布的曲线做切线，与 $t_\infty$ 流体温度相交于 $0'$ ，与横坐标夹角为 $\alpha$ ,因此这条切线的斜率为 $tan\alpha,\alpha+\phi=180°,tan\alpha=-tan\phi$ ，即：

$\frac{\partial t}{\partial x}|_{x=\delta}=tan\alpha=-tan\phi\\ -\frac{\partial t}{\partial x}|_{x=\delta}=tan\phi$

而根据平面几何: $tan\phi=\frac{t|_{x=\delta}-t_\infty}{x'}$

$-\frac{\partial t}{\partial x}|_{x=\delta}=\frac{t|_{x=\delta}-t_\infty}{\delta/Bi}\tag{30}$
对比公式(30) $tan\phi$ 等于其右边，那么我们可以得到
$x'=\frac{\lambda}{h}=\frac{\delta}{Bi}\tag{31}$

思考，为什么所有不同时间的温度分布曲线的切线与 $t=t_\infty$ 这条温度轴交于同一点？

$\because x'$ 与时间无关，只与板的半个厚度 $\delta$ 和毕渥数Bi有关

$\therefore$ 一旦板的物性确定（ $\lambda$ ），边界流体的条件确定( $h$ )，那么交点都汇聚一点，也可以认为板厚确定，Bi数确定之后交点就确定.

点O'距离壁面的距离为 $x'=\frac{\lambda}{h}=\frac{\delta}{Bi}$

任何时刻，壁面温度分布的切线都通过坐标为 $(\delta+\delta/Bi,t_\infty)$ 的O'点——第三类边界条件的定向点

定性分析

①. $Bi\to\infty$ ,相当于h无穷大，从纯物理意义上看，对流换热热阻趋于0，平壁的表面温度从冷却过程一开始就立即降温（垂直下降,不经过空间距离）到流体温度 $t_\infty$ ，从几何上看，x'=0,定向点O'就在平壁表面上

003.png

②. $Bi\to 0$ ,物理上分析意味着导热热阻非常小,内部导热性能非常好，看成是等温分布，任意时刻温度分布曲线为一水平线；几何上分析，交点x'在无限远处，所以温度分布必须为水平线，才会有定向点为无穷远的情况。

004.png

③.Bi为中间状况是，交点在 $(-\delta-\delta/Bi,t_\infty，\delta+\delta/Bi,t_\infty)$ 两点，温度分布为如下:

005.png

④.当Bi<0.1的时候，工程上就可以把他作为接近极限的判据，中心温度与表面温度之差<5%,接近均匀一致，采用集中参数法（集总参数法），当做0维物体处理.

最后编辑于：2024.11.14 09:44:13

人面猴
序言：七十年代末，一起剥皮案震惊了整个滨河市，随后出现的几起案子，更是在滨河造成了极大的恐慌，老刑警刘岩，带你破解...
沈念sama阅读 204,793评论 6赞 478
死咒
序言：滨河连续发生了三起死亡事件，死亡现场离奇诡异，居然都是意外死亡，警方通过查阅死者的电脑和手机，发现死者居然都...
沈念sama阅读 87,567评论 2赞 381
救了他两次的神仙让他今天三更去死
文/潘晓璐我一进店门，熙熙楼的掌柜王于贵愁眉苦脸地迎上来，“玉大人，你说我怎么就摊上这事。” “怎么了？”我有些...
开封第一讲书人阅读 151,342评论 0赞 338
道士缉凶录：失踪的卖姜人
文/不坏的土叔我叫张陵，是天一观的道长。经常有香客问我，道长，这世上最难降的妖魔是什么？我笑而不...
开封第一讲书人阅读 54,825评论 1赞 277
港岛之恋（遗憾婚礼）
正文为了忘掉前任，我火速办了婚礼，结果婚礼上，老公的妹妹穿的比我还像新娘。我一直安慰自己，他们只是感情好，可当我...
茶点故事阅读 63,814评论 5赞 368
恶毒庶女顶嫁案：这布局不是一般人想出来的
文/花漫我一把揭开白布。她就那样静静地躺着，像睡着了一般。火红的嫁衣衬着肌肤如雪。梳的纹丝不乱的头发上，一...
开封第一讲书人阅读 48,680评论 1赞 281
城市分裂传说
那天，我揣着相机与录音，去河边找鬼。笑死，一个胖子当着我的面吹牛，可吹牛的内容都是我干的。我是一名探鬼主播，决...
沈念sama阅读 38,033评论 3赞 399
双鸳鸯连环套：你想象不到人心有多黑
文/苍兰香墨我猛地睁开眼，长吁一口气：“原来是场噩梦啊……” “哼！你这毒妇竟也来了？” 一声冷哼从身侧响起，我...
开封第一讲书人阅读 36,687评论 0赞 258
万荣杀人案实录
序言：老挝万荣一对情侣失踪，失踪者是张志新（化名）和其女友刘颖，没想到半个月后，有当地人在树林里发现了一具尸体，经...
沈念sama阅读 42,175评论 1赞 300
护林员之死
正文独居荒郊野岭守林人离奇死亡，尸身上长有42处带血的脓包…… 初始之章·张勋以下内容为张勋视角年9月15日...
茶点故事阅读 35,668评论 2赞 321
白月光启示录
正文我和宋清朗相恋三年，在试婚纱的时候发现自己被绿了。大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
茶点故事阅读 37,775评论 1赞 332
活死人
序言：一个原本活蹦乱跳的男人离奇死亡，死状恐怖，灵堂内的尸体忽然破棺而出，到底是诈尸还是另有隐情，我是刑警宁泽，带...
沈念sama阅读 33,419评论 4赞 321
日本核电站爆炸内幕
正文年R本政府宣布，位于F岛的核电站，受9级特大地震影响，放射性物质发生泄漏。R本人自食恶果不足惜，却给世界环境...
茶点故事阅读 39,020评论 3赞 307
男人毒药：我在死后第九天来索命
文/蒙蒙一、第九天我趴在偏房一处隐蔽的房顶上张望。院中可真热闹，春花似锦、人声如沸。这庄子的主人今日做“春日...
开封第一讲书人阅读 29,978评论 0赞 19
一桩弑父案，背后竟有这般阴谋
文/苍兰香墨我抬头看了看天上的太阳。三九已至，却和暖如春，着一层夹袄步出监牢的瞬间，已是汗流浃背。一阵脚步声响...
开封第一讲书人阅读 31,206评论 1赞 260
情欲美人皮
我被黑心中介骗来泰国打工，没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留，地道东北人。一个月前我还...
沈念sama阅读 45,092评论 2赞 351
代替公主和亲
正文我出身青楼，却偏偏与公主长得像，于是被迫代替她去往敌国和亲。传闻我的和亲对象是个残疾皇子，可洞房花烛夜当晚...
茶点故事阅读 42,510评论 2赞 343

非稳态导热微分方程的分析解

掌握从非稳态过渡到稳态的感性认识

温度分布和热流量分布随时间和空间的变化规律的函数原型：

导热微分方程的一般形式(第二章推导过)

求解方法有三种，分析解、数值解、实验解

一.下面为非稳态导热，有限厚度无限大平板非稳态导热的解析解（分离变量法）

1.物理模型为厚度的平板，初始温度分布为，放入的流体中，求温度分布

该物理模型简化为一维问题，非稳态，常物性（为常数），两侧温度分布对称，中心为原点，无內热源，边界条件为在传热系数为的流体中的第三类边界条件,微分方程简化为：

①时间条件为初始条件,该条件只有瞬态条件才有，稳态热传导没有的初始条件

②温度分布对称,中心温度最大或者最低，因此中心为绝热边界条件，外边界条件以第三类边界条件（举例）:

③分析条件的齐次与非齐次问题，方程（2）为齐次，处的边界条件为非齐次，x=0的边界条件为齐次，全部齐次化：,因此,并且注意到由于为常数，在求一阶偏微分，二阶偏微分的时候不影响，其他方程仍为齐次(,)

方程组引入过余温度后，转为为如下方程组:

2.分离变量法,假设

分析上述公式，左边是关于时间的函数，右边是关于空间的函数, 秒=米？

要让秒=米是不可能的，等式两边必须同等于一常数

转换为两个独立方程，如下

解常微分方程,第一个方程很好解,,

分析以下上面的数据,为热扩散率，大于0，时间恒为正，如果D为正数，那么，而实际情况是时间越长温度越接近流体温度，过余温度接近0，因此D一定为负数

设，那么方程组(7)变化为如下：

对于来说，解微分方程如下为

对于关于x位置的函数，实际上为,是二阶导数与函数值的线性组合，这个解在高数(下)中查到结果

由于之前分离变量的假设为,那么我们把就得到了过余温度关于函数的一般解形式:

其中未知系数c1,c2,c3的换为新的未知数

3.边界条件及初始条件确定未知系数得到定解

（β是我们引入的） 三个未知数，三个定解条件。

方程组中有一个初始条件，二个定解条件。三个方程决定三个未知数，此方程有定解

先看x=0处的绝热边界条件 ,将通解形式(12)对x求偏导带入，如下：

x=0带入，求得（指数函数永远不为0）,而,那么只有B=0。

一般解简化为如下，1个定解条件解得一个B=0

我们学过了Bi数的定义，为内阻：外阻, ,将等式（16）右边分子分母同时乘以凑出一个Bi数

（16）转换为,换元

这个方程为超越方程，左边是tan x 三角函数，右边是1/x 函数，并且tan x为周期函数，两个曲线有无穷个交点，就有无穷个解.

那么，这个超越方程的解集合称之为特征方程，无穷多个特征解为

这些特征解的无穷级数求和，组成新的解形式,一般情况下用6-8项就够了，下面我们简单分析一下。

1.当

2.当

3.在给定Bi数的情况下，对应每一个特征值，温度分布的特解分别为如下：

反直觉之——有多个温度分布的解析式:

公式集合(18)为及得到的特解，常数为任何值，都满足导热微分方程式和这关于x位置的两个定解条件，但是一个都不满足关于时间时刻的初始条件的温度分布——因为只要有一个不相等，就会有多个初始温度分度，产生悖论。

知识点补充：导热微分方程为线性方程

温度对时间一阶导数的系数式1，温度对位置x的二阶导数为，都与温度无关

该导热问题的通解为各个特解的线性求和叠加：导热微分方程式和边界条件都是线性的——温度和温度的各阶导数项的系数都与温度无关。

线性叠加的组合得到一个与初始温度分布相等的表达式：

公式（19）中为已知数，为未知数，带入初始条件求解得到即可得到唯一特解

A_n还是不知道，级数方程没有复变的基础下不好解，下面还是利用特征函数的特性-正交性来解决,当下标不等（非对角元乘积的积分为0）

正交性问题来源于常微分方程的施图姆-刘维尔理论 $

这里不展开，二阶常微分方程为亥姆霍兹方程，利用该结论即可。施图姆-刘维尔理论 (Sturm–Liouville theory ) - 知乎 (zhihu.com)

需要利用的结论是

方程(20)两边同时乘以,并在范围内积分,这个m可以等于也可以不等于n

那么公式中无穷求和只剩下 m等于n这一项，去掉级数符号

这个积分就好解了，可以求得等于两个积分相除

最后将的表达式带入公式（19）中

得到

唯一决定解析解的公式已经推导完毕，我们进行一下变形然后进行对时间和位置变化的趋势行为的分析

已知毕沃准则数,为内阻与外阻之比；傅里叶准则数，L为特征长度，在的平板中,,为无量纲时间

将特解(24)转换为这两个准则数的自变量的形式,利用

在（25）中，为Bi数相关，为位置相关，为时间相关，因此该解析式包含了①内部导热与外部条件的因素②位置自变量③时间自变量这三个因变量

先对傅里叶数进行分析，当时（根据热扩散率和厚度估算一般为20秒），该指数函数衰减的速度超过了前面函数(三角函数)增加的速度，无穷级数在第1项就几乎收敛！那么我们用第一项就可以得到精度很高的近似解。

在上述特殊情况，可以用上式取级数①第一项进行计算，②可以用诺莫图计算，③也可以用campo近似拟合公式计算,将过余温度表达式与初始时刻的过于温度相除，并且分子分母同时乘以中心处任意时间的过余温度

有了温度分布，根据计算换热量

二、准则数对温度分布的影响

1.数对温度分布的影响

时，无量纲温度表达式取无穷级数第一项,有一个根即可

①，过余温度的对数随时间变化的斜率为常数，与时间和位置无关，后称之为正规状况阶段：初始温度分布的影响已经消失

②对应的非正规状况阶段:中心处收到初始温度的影响变化速度比表面上更慢。

除开大平壁以外其他模型也有类似的结论

非稳态导热分为

2.Bi准则数对温度分布的影响:

回顾一下温度分布的曲线,中心处冷却的慢，边界上冷却的快

上图为定性描述分析，将边界处的温度分布的曲线做切线，与流体温度相交于，与横坐标夹角为,因此这条切线的斜率为，即：

而根据平面几何:

思考，为什么所有不同时间的温度分布曲线的切线与这条温度轴交于同一点？

与时间无关，只与板的半个厚度和毕渥数Bi有关

一旦板的物性确定（），边界流体的条件确定()，那么交点都汇聚一点，也可以认为板厚确定，Bi数确定之后交点就确定.

点O'距离壁面的距离为

任何时刻，壁面温度分布的切线都通过坐标为的O'点——第三类边界条件的定向点

定性分析

①.,相当于h无穷大，从纯物理意义上看，对流换热热阻趋于0，平壁的表面温度从冷却过程一开始就立即降温（垂直下降,不经过空间距离）到流体温度，从几何上看，x'=0,定向点O'就在平壁表面上

②.,物理上分析意味着导热热阻非常小,内部导热性能非常好，看成是等温分布，任意时刻温度分布曲线为一水平线；几何上分析，交点x'在无限远处，所以温度分布必须为水平线，才会有定向点为无穷远的情况。

③.Bi为中间状况是，交点在两点，温度分布为如下:

1.物理模型为 $2\delta$ 厚度的平板， $\tau=0$ 初始温度分布为 $t_0(x)$ ，放入 $t_\infty$ 的流体中，求温度分布 $t(x,\tau)$

该物理模型简化为一维问题，非稳态，常物性（ $\lambda , a$ 为常数），两侧温度分布对称，中心为原点，无內热源，边界条件为在传热系数为 $h$ 的流体中的第三类边界条件,微分方程 $(1)$ 简化为：

①时间条件为初始条件,该条件只有瞬态条件才有，稳态热传导没有 $\tau=0$ 的初始条件

方程组引入过余温度后， $(2)(3)$ 转为为如下方程组:

2.分离变量法,假设 $\theta=f(x，\tau)=X(x)\cdot\phi(\tau)$

解常微分方程, $(7)的$ 第一个方程很好解, $\frac{1}{\phi}\frac{d \phi}{d \tau}=aD$ , $ln\phi=D+c_1$

分析以下上面的数据, $a$ 为热扩散率，大于0，时间恒为正，如果D为正数，那么 $\phi随时间增大会不会趋近0$ ，而实际情况是时间越长温度越接近流体温度，过余温度接近0，因此D一定为负数

设 $D=-\beta^2$ ，那么方程组(7)变化为如下：

对于 $关于时间\tau的函数\phi$ 来说，解微分方程如下为

对于关于x位置的函数，实际上为 $X^{''}+\beta^2X=0$ ,是二阶导数与函数值的线性组合，这个解在高数(下)中查到结果

由于之前分离变量的假设为 $\theta(x,\tau)=X(x)\cdot\phi(\tau)$ ,那么我们把 $(10)\cdot (11)$ 就得到了过余温度关于 $x,\tau$ 函数的一般解形式:

其中未知系数c1,c2,c3的换为新的未知数 $A=c_1\cdot c_2,B=c_1\cdot c_3$

$A，B，\beta$ （β是我们引入的）三个未知数，三个定解条件。

方程组 $(4)$ 中有一个初始条件，二个定解条件。三个方程决定三个未知数，此方程有定解

先看x=0处的绝热边界条件 $x=0,\frac{\partial \theta}{\partial x}=0$ ,将通解形式(12)对x求偏导带入，如下：

x=0带入，求得 $B\beta\cdot e^{-a\beta^2\tau}=0$ （指数函数永远不为0）,而 $\beta \neq 0$ ,那么只有B=0。

我们学过了Bi数的定义，为内阻：外阻, $Bi=\frac{\delta}{\lambda}/\frac{1}{h}=\frac{h\delta}{\lambda}$ ,将等式（16）右边分子分母同时乘以 $\delta$ 凑出一个Bi数

（16）转换为 $tan(\beta \delta)=\frac{Bi}{\beta\delta}$ ,换元 $\mu=\beta\delta$

那么，这个超越方程的解集合称之为特征方程，无穷多个特征解为 $\mu_1,\mu_2……\mu_n$

1.当 $Bi\to\infty,也就是函数tan\mu\to\infty，特征值\mu_n为\frac{1}{2}\pi,\frac{3}{2}\pi,\frac{5}{2}\pi……$

2.当 $Bi\to 0,tan\mu\to 0,，特征值\mu_n为0,\pi,2\pi……$

3.在给定Bi数的情况下，对应每一个特征值 $\mu_n=\beta_n\delta\to\beta_n=\frac{\mu_n}{\delta}$ ，温度分布的特解分别为如下：

知识点补充：导热微分方程为线性方程 $\frac{\partial \theta}{\partial \tau}=a\frac{\partial^2\theta}{\partial x^2}$

温度对时间一阶导数的系数式1，温度对位置x的二阶导数为 $a$ ，都与温度无关

公式（19）中 $\beta_n$ 为已知数， $A_n$ 为未知数，带入初始条件 $\tau=0,\theta=\theta_0=t_0-t_\infty$ 求解得到 $A_n$ 即可得到唯一特解

方程(20)两边同时乘以 $cos(\mu_m \frac{x}{\delta})$ ,并在 $0\leq x\leq \delta$ 范围内积分,这个m可以等于也可以不等于n

这个积分就好解了，可以求得 $A_n$ 等于两个积分相除

最后将 $A_n$ 的表达式带入公式（19）中

已知毕沃准则数 $Bi=\frac{h\delta}{\lambda}$ ,为内阻与外阻之比；傅里叶准则数 $F_o=\frac{a\tau}{L^2}$ ，L为特征长度，在 $2\delta$ 的平板中, $L=\delta，F_o=\frac{a\tau}{\delta^2}$ ,为无量纲时间

将特解(24)转换为这两个准则数的自变量的形式,利用 $\mu=\beta\delta;\beta_n=\frac{\mu_n}{\delta}$

在（25）中， $\mu$ 为Bi数相关， $\frac{x}{\delta}$ 为位置相关， $F_0$ 为时间相关，因此该解析式包含了①内部导热与外部条件的因素②位置自变量③时间自变量这三个因变量

有了温度分布，根据 $Q=cm\Delta t$ 计算换热量

1. $F_o$ 数对温度分布的影响

$F_o=a\tau/\delta^2>0.2$ 时，无量纲温度表达式取无穷级数第一项,有一个根 $\mu_1$ 即可

① $F_o>0.2$ ，过余温度的对数随时间变化的斜率为常数，与时间和位置无关，后称之为正规状况阶段：初始温度分布的影响已经消失

② $F_o<0.2$ 对应的非正规状况阶段:中心处收到初始温度的影响变化速度比表面上更慢。

上图为定性描述分析，将边界处的温度分布的曲线做切线，与 $t_\infty$ 流体温度相交于 $0'$ ，与横坐标夹角为 $\alpha$ ,因此这条切线的斜率为 $tan\alpha,\alpha+\phi=180°,tan\alpha=-tan\phi$ ，即：

而根据平面几何: $tan\phi=\frac{t|_{x=\delta}-t_\infty}{x'}$

思考，为什么所有不同时间的温度分布曲线的切线与 $t=t_\infty$ 这条温度轴交于同一点？

$\because x'$ 与时间无关，只与板的半个厚度 $\delta$ 和毕渥数Bi有关

$\therefore$ 一旦板的物性确定（ $\lambda$ ），边界流体的条件确定( $h$ )，那么交点都汇聚一点，也可以认为板厚确定，Bi数确定之后交点就确定.

点O'距离壁面的距离为 $x'=\frac{\lambda}{h}=\frac{\delta}{Bi}$

任何时刻，壁面温度分布的切线都通过坐标为 $(\delta+\delta/Bi,t_\infty)$ 的O'点——第三类边界条件的定向点

①. $Bi\to\infty$ ,相当于h无穷大，从纯物理意义上看，对流换热热阻趋于0，平壁的表面温度从冷却过程一开始就立即降温（垂直下降,不经过空间距离）到流体温度 $t_\infty$ ，从几何上看，x'=0,定向点O'就在平壁表面上

②. $Bi\to 0$ ,物理上分析意味着导热热阻非常小,内部导热性能非常好，看成是等温分布，任意时刻温度分布曲线为一水平线；几何上分析，交点x'在无限远处，所以温度分布必须为水平线，才会有定向点为无穷远的情况。

③.Bi为中间状况是，交点在 $(-\delta-\delta/Bi,t_\infty，\delta+\delta/Bi,t_\infty)$ 两点，温度分布为如下: