两相流：颗粒群受到的曳力

颗粒和流体之间的作用力是两相流领域最重要的本构关系。任何科学问题的研究都遵循先易后难的原则，对颗粒流体间作用力，同样是先从单颗粒、低雷诺数的情况研究，这篇文章介绍了处在爬流状态下单颗粒的受力。本文介绍颗粒群与流体间的作用力，尤其是颗粒群所受曳力。

一、颗粒群曳力关联式研究思路

符号说明：\(u\)无穷远处来流速度；\(u_0\)床层内流体真实速度；\(U\)表观速度，即床层进口流体速度，也是化工原理中的“空塔速度”，\(U=(1-\phi)u_0\)；颗粒群床层固含率\(\phi\)；粒径\(d\)

单颗粒曳力

1. \(\phi \rightarrow 0\)，\(Re\rightarrow 0\)，Stokes定律：\(C_D=24/Re\)

\[F=3\pi \mu d u=\frac{1}{2}\rho_g u^2 \frac{\pi}{4}d^2 C_D \]

2. \(\phi \rightarrow 0\)，\(Re<1\)，Oseen公式，\(C_D=\frac{24}{Re}(1+\frac{3}{16}Re)\)

3. \(\phi \rightarrow 0\)，\(500<Re<2\times 10^5\)，Newton公式，\(C_D=0.44\)

4. \(\phi \rightarrow 0\)，标准阻力公式，\(C_D=\frac{24}{Re}(1+\frac{1}{6}Re^{2/3})\)

颗粒群曳力

流体流过静止的颗粒床层，床层固含率为\(\phi=NV_p/V\)，底部入口流速为\(U\)，\(U=(1-\phi)u_0\)

颗粒受到两种力：由流体运动产生的曳力\(F_d\)，压力梯度力\(F_p=-V_p\nabla P\)，因此颗粒受到的合力为\(F_{gs}=F_d+F_p\)

流体流过床层的压降损失大小应等于颗粒对其的作用力：

\[(P_{bottom}-P_{top})S_{bed}=V\nabla P=NF_{gs} \]

因此，\(F_d=(1-\phi)F_{gs}\)，定义颗粒群中单颗粒曳力与Stokes曳力之比为无量纲曳力，

\[F(\phi,Re)=F_d/3\pi \mu d U , \quad F(0,0)=1, \quad Re=\frac{dU\rho_g}{\mu} \]

工程中常用床层压降来关联曳力，转化成床层压降为

\[F(\phi,Re)=-\frac{1-\phi}{\phi} \frac{\nabla P d^2}{18\mu U} \]

目前有了单颗粒无量纲曳力的表达式，研究者进一步从以下两种途径来从简单情况推向复杂情况

方式一：从低雷诺数出发

\[F(\phi, Re)=F(\phi,0)+\alpha (\phi)Re \tag{1} \]

通过第二项来添加\(Re\)的影响。
研究者提出了以上两个方程的形式：

\[F(\phi,0)=\frac{a\phi}{18(1-\phi)^2},\quad \alpha(\phi)=\frac{b}{18(1-\phi)^2} \tag{2} \]

Ergun总结了大量实验，得出上式中\(a=150,b=1.75\)，其他研究者也提出了不同的\(ab\)值，例如著名的Koch-Hill-Ladd曳力。Van der Hoef给出了目前最精确的\(F(\phi,0)\)表达式：

\[F(\phi,0)=\frac{180\phi}{18(1-\phi)^2}+(1-\phi)^2 (1+1.5\sqrt{\phi}) \tag{3} \]

方式二：从低固含率出发

\[F(\phi, Re)=F(0,Re)(1-\phi)^{-\beta} \tag{4} \]

第一项\(F(0,Re)\)在前文《单颗粒曳力》已介绍过。通过第二项来添加\(\phi\)的影响，其确定方式常通过不同\(\phi\)的颗粒群沉降实验测量，测量沉降速度\(v_t(0)/v_t(\phi)\)的比值，然后将其关联为\(\phi\)的函数：

\[v_t(0)/v_t(\phi)=(1-\phi)^{-n} \]

注意，上式中的\(n\)与第二项中的\(\beta\)含义不同。对于低\(Re\)，\(\beta=n-1\)；对于高\(Re\)，\(\beta=2n-1\)。Richardson提出，

\[n=4.65,\quad Re<0.2 \\ n=2.4, \quad Re>500 \]

Wen-Yu曳力认为\(\beta=3.7\)，著名的Di Felice曳力认为\(\beta\)与\(Re\)之间存在弱相关。

两种方式比较

第一种方式从低雷诺数出发，其适用范围是稠密体系即\(\phi\)较大的体系。第二种方式则适用于低固含率即稀疏体系。因此Gidaspow综合以上两种方式的优缺点，提出在\(\phi<0.2\)应用第二种，\(\phi>0.2\)应用第一种。Gibilaro同样综合以上两种方式，提出了新的关联式

\[F(\phi,Re)=f(Re)(1-\phi)^{-3.8} \]

使用模拟数据构建曳力关联式的研究流程

前文介绍了通过理论分析和实验测量获取颗粒群中单颗粒曳力的流程。本节介绍一种从模拟数据出发构建曳力关联式的流程，该曳力关联式为著名的Beetstra曳力（参考文献1）。

模拟方法为直接数值模拟，流体采用LBM求解，LBM网格为粒径的\(1/20\)左右，颗粒和流体之间无滑移。模拟体系为长方体，体系中包含54个颗粒，通过改变长方体大小来获得不同体系固含率。颗粒群按固定速度\(v\)运动，相对来看，流体速度为\(U=-v\)，\(Re=dU\rho_g/\mu\)

直接数值模拟可以得到每一个颗粒所受的流体作用力\(F_{gs,i}\)，则平均每个颗粒所受的流体作用力为

\[\overline{F}_{gs}=\frac{\sum F_{gs,i}}{N} \]

该力与曳力之间的关系在上文中已经提到，\(F_d=(1-\phi)F_{gs}\)，因此，无量纲曳力关联式的形式为：

\[F(\phi,Re)=-\frac{(1-\phi)\overline{F}_{gs}}{3\pi \mu d U} \]

前文已经提到构建曳力关联式的两种方式，作者在这里同样从两种方式出发来构建。

依据方式1

从低雷诺数出发即公式1。

由公式1，作者计算了\(\alpha\)，其中\(F(\phi,0)\)由公式3所得。

\[\alpha(\phi,Re)=\frac{F(\phi,Re)-F(\phi,0)}{Re} \]

若\(\alpha\)符合公式2，那么\(\alpha\cdot (1-\phi)^2\)应为常数，而作者文中图3的结果显然不是这样，说明\(\alpha\)应是\((\phi,Re)\)两者的函数。因此，由模拟数据拟合\(\alpha(\phi,Re)\)即可得到曳力关联式。

依据方式2

从低固含率出发即公式4。其中待定系数为\(\beta\)，待定项为\(F(0,Re)\)，因为目前还没有十分统一的\(F(0,Re)\)的表达形式。作者提出White总结的表达式能较好的拟合\(F(0,Re)\)，因此使用White的表达式。

公式4表达的含义是曳力可以表达为两个函数乘积的形式且两个变量各自独立：\(F(\phi,Re)=f_1(\phi) f_2(Re)\)。不过模拟结果表明，该假设并不成立。如果一定要使用这种形式，则需要给出\(\beta\)，即

\[\beta(\phi,Re)=-log\left( \frac{F(\phi,Re)}{F(0,Re)} \right)/log(1-\phi) \]

使用模拟数据拟合\(\beta\)，最后可以得到依据方式2建立的曳力关联式

两种方式比较

先说结论：作者推荐用方式1构建曳力关联式，并且作者给出了\(\alpha\)表达式的拟合结果。关联式预测结果与模拟结果间的误差在3%以内。不用方式2构建主要是因为表达式中含有\(F(0,Re)\)，用不同研究者提出的\(F(0,Re)\)会得到差异很大的\(\beta\)

参考文献

1. Beetstra R., et al. AIChE J. 2007, 53(2).
2. Di Felice R. Chem Eng Sci. 1995, 50.
3. Van der Hoef MA, et al. J Fluid Mech. 2005, 528.
4. White FM. Viscous fluid flow. pp 206–210. New York: McGraw-Hill; 1974.