恒容容器放气的瞬时流量的计算

有时候，你会遇到一个问题，该问题的描述如下：

你有一个已知体积的容器，设容器体积为V，里面装有一定压力(初始压力)的气体，如空气或氢气等，设初始压力为1MPa，容器出口连接着一个阀门开关，开关后面接1/4in.的钢管，钢管出口即为气体出口。当阀门瞬间全开时，气体出口的瞬时流量值随时间变化到底是怎么样的呢？

该问题相当于在不考虑管壁与管长对气体产生粘滞阻力的影响下，已知气体管道直径d，即管道横截面积 $A=\frac{\pi d^{2}}{4}$ ，已知管子进口静压为 $P_0=1MPa$ ，已知管子出口静压为 $P_{b}=P_{atm}$ ，即一个大气压，同时知道进口气体总温 $T_{0}$ 为323K，求出口瞬时流量随时间 $t$ 的变化关系和曲线。

1. 第一种方法：根据哈根泊谡叶方程

利用理想气体方程： $PV=\frac{m}{M}RT$ 和哈根泊谡叶关系式： $q=\frac{\pi r^{4}}{8\mu L}\Delta P,\Delta P=P-P_{atm}$ ，两个方程联立，可得 $ln(\frac{P-P_{atm}}{P_0-P_{atm}})=-ABt,A=\frac{\pi r^{4}}{8\mu L},B=\frac{\rho RT}{MV}$ ,利用该关系式，得到 $P$ 随时间 $t$ 的关系，然后再根据质量流量的关系，或者直接套回哈根泊谡叶关系式得到瞬时流量： $\frac{dm}{dt}=-\frac{MV}{R}\cdot \frac{1}{T_0}\frac{dP}{dt}, T_0=323K, V=0.001m^{3}, R=8.314J/(mol\cdot K),M=0.002 kg/mol$

可以知道，瞬时流量是一指数函数曲线形式，见图2。通过积分，可得知积分总流量为 $7.45L$ ，这个跟理想气体方程 $\Delta PV=\Delta mRT,\Delta m=\frac{\Delta PV}{RT}\Leftrightarrow \Delta q=\frac{\Delta m}{\rho }=\frac{\Delta PV}{\rho RT}=7.46L$ 差分得出的总流量非常接近。不过根据推导，一开始的瞬时流量值非常离谱，可以去到 $409094SLM(L/min)$ ，根据 $\frac{dm}{dt}=\rho \cdot u\cdot A\Leftrightarrow \frac{dq}{dt}=\frac{d(m/\rho )}{dt}=u\cdot A$ ，可以知道出口流体平均速度 $u=409094\cdot 10^{-3}/60/(\pi \cdot 0.00635^{2}/4)=2320388.88m/s$ ，光速是 $u_c=299792458m/s$ ，出口速度已经达到 $0.007$ 倍的光速，也超过空气声速 $\approx 314m/s$ ，岂不是一个超音速流体？

关于哈根泊谡叶关系式的推导，见下图1。

图1

图2

2. 第二种方法：根据气体动力学推算

假设排气过程与气体管道壁面的换热忽略不计，即壁面是绝热的，气体流体是一个准稳态问题，排气口相当于是收缩，没有扩张，根据气体动力学可知，出口气体流速只能加速到马赫数1，即 $Ma=1$ 。根据总静温关系式 $\frac{T_{0}}{T_{b}}=1+\frac{k-1}{2}Ma^{2},T_0=323K$ ，得知 $T_b\approx 269.17K$ 。再根据马赫数定义式 $Ma^{2}=\frac{v^{2}}{c^{2}}=\frac{v^{2}}{kR_gT}$ ，这里 $k$ 是气体比热比，定义为定压比热 $C_p$ 与定容比热 $C_v$ 之比，即 $v_b=Ma\cdot \sqrt{kR_gT_b},k=1.4,R_g=\frac{R=8.314J/(mol\cdot K)}{0.002kg/mol}=4157J/(kg\cdot K),T_b=269.17K$

得到氢气气体流速 $v_b\approx 1295.22m/s$ 。根据 $\frac{dm}{dt}=\rho \cdot u\cdot A$ ， $P=\rho R_gT$ , $\frac{T_{0}}{T_{b}}=1+\frac{k-1}{2}Ma^{2},T_0=323K$ ， $\frac{P}{P_b}=(1+\frac{k-1}{2}Ma^{2})^{\frac{k}{k-1}}$ , $u=Ma\cdot \sqrt{kR_gT_b}$ ，在 $Ma=1$ 的壅塞流阶段，可解得 $ln\frac{P}{P_0}=\frac{\sqrt{1.2kR_gT_0}}{V}MaAt\Leftrightarrow P=f(t)=P_0e^{-47.734t}$ ，这阶段，理解为 $u$ 不变， $P$ 变化导致的 $\rho$ 变化，瞬时流量也会随之变化，得到下图图3中从左到右的第二张图的蓝色部分的曲线，蓝色部分为积分流量，通过积分得到壅塞流下的总流量为 $6.69L$ 。后面非壅塞流状态下的亚声速流，原则上也是利用 $\frac{dm}{dt}=\rho \cdot u\cdot A$ ， $P=\rho R_gT$ ， $u=Ma\cdot \sqrt{kR_gT_b}$ ， $\frac{T_{0}}{T_{b}}=1+\frac{k-1}{2}Ma^{2},T_0=323K$ ， $\frac{P}{P_b}=(1+\frac{k-1}{2}Ma^{2})^{\frac{k}{k-1}}$ ，这5个式子得到 $\frac{dm}{dt}=f(P)$ 的关系，我用欧拉法获得解析解的近似值，得到后续的流量曲线，如下图图3中从左到右的第三张图的后续蓝色部分的曲线，通过积分面积算得亚声速流下总流量为 $0.747L$ ，因此总流量为 $6.69L+0.747L\approx 7.44L$ 。