前文提到，经过求解 pEqn，并修正速度Ua和Ub以后，总体的一连性便得到了保证。为了得到分散相的体积分率alpha，还需要使用得到的速度场来求解分散相的一连性方程，即alphaEqn。分散相一连性方程可以表达如下：
∂ α a ∂ t + ∇ ⋅ ( α a U a ) = 0 \frac{\partial \alpha_a}{\partial t}+\nabla \cdot (\alpha_a U_a)=0 ∂t∂αa+∇⋅(αaUa)=0
为了让每一项都写成守恒情势，而且保证 α a \alpha_a αa的有界性，Weller 将分散相一连性方程写成如下情势
∂ α a ∂ t + ∇ ⋅ ( α a U ) + ∇ ⋅ ( U r α a ( 1 − α a ) ) = 0 \frac{\partial \alpha_a}{\partial t}+\nabla \cdot (\alpha_a U)+\nabla \cdot(U_r\alpha_a(1-\alpha_a))=0 ∂t∂αa+∇⋅(αaU)+∇⋅(Urαa(1−αa))=0
其中
U = α a U a + α b U b U r = U a − U b U=\alpha_a U_a + \alpha_b U_b \\ U_r=U_a-U_b U=αaUa+αbUbUr=Ua−Ub
OpenFOAM-2.1.1 的alphaEqn求解的正是 Weller 提出的情势。主要的代码如下：
有关 OpenFOAM 中使用的 KTGF 模型的理论可以参考Derivation, Implementation, and Validation of Computer Simulation Models for Gas-Solid Fluidized Beds，以及B.G.M. van Wachem 的其他相干论文。
g0
这一段的结果，相当于在本系列第一篇最后的分散相动量方程的等式右边额外增加一项：
− g 0 ∗ m i n ( e p r e A l p h a E x p ∗ ( α a − α M a x ) , e x p M a x ) ∇ α a / ( α a ρ a ) -g0*min(e^{preAlphaExp*(\alpha_a - \alpha_{Max})},expMax)\nabla \alpha_a / (\alpha_a \rho_a) −g0∗min(epreAlphaExp∗(αa−αMax),expMax)∇αa/(αaρa)
可见，这一项的作用是给固相额外施加了一个力，可以认为是固相压力，这个力与 α a \alpha_a αa的梯度有关，且与梯度的方向相反。
留意这一项的特点：g0, preAlphaExp, alphaMax, expMax 都是需要用户指定的参数，一般将g0, preAlphaExp和expMax设置为一个比力大的正整数（ 1 0 3 10^3 103的量级），当alpha - alphaMax < 0时， e p r e A l p h a E x p ∗ ( α a − α M a x ) e^{\,preAlphaExp*(\alpha_a - \alpha_{Max})} epreAlphaExp∗(αa−αMax) 将是一个很小的数，此时整个增加的一项也将是一个较小的数，以是，alpha - alphaMax < 0时额外增加的那一项对动量方程的贡献很小。但是，当alpha - alphaMax > 0时，随着 alpha 偏离 alphaMax 越来越远，exp(preAlphaExp*(alpha - alphaMax)将敏捷增大，min(exp(preAlphaExp*(alpha - alphaMax)), expMax)的值也将敏捷增大，直到设定值expMax。可见，g0项的作用可以明白为为了防止固相过度堆积。气固系统中，固相的体积分率是有上限的，可以将alphaMax设置为这个上限值，当某个网格里的alpha超过设定的alphaMax时，就让固相敏捷弹开，以防止固相体积分率过大。
而 fvm::div(-fvc::flux(-phir, beta, schemer), alpha, schemer) 则比 $ \nabla \cdot (\alpha_a U)+\nabla \cdot(U_r\alpha_a(1-\alpha_a) $ 多出了 fvm::div(-fvc::flux(-phipp, beta, schemer), alpha, schemer) 一项，写成公式就是 $ \nabla \cdot [\alpha (\beta *\mathrm{ppMagf}) \nabla \alpha] $ 。
多出的这两项加起来，为
∇ ⋅ [ α ∗ p p M a g f ∇ α ] \nabla \cdot [\alpha *\mathrm{ppMagf} \ \nabla \alpha] ∇⋅[α∗ppMagf ∇α]
对比上面的 alphaEqn 减去的 laplacian 项，会发现这一项跟谁人 laplacian 是一样的。
以是，g0 > 0 时，先将固相压力带来的通量代入到 alphaEqn 的构建中，然后再减去对应的 laplacian 项，这么做的目的，应该是为了计算稳固性以及保证 alpha 的有界性（保证 0 < a l p h a < a l p h a M a x 0<alpha<alphaMax 0<alpha<alphaMax）。但是这背后的数值原理，我也没有完全明白。
最后，有必要说明一下，g0相干的项着实就是在动量方程中增加了一项颗粒压力，想要达到的结果是防止固相过度堆积。kineticTheory开启后将计算颗粒相的应力，其中包括了颗粒相压力。以是 g0可以当作 KTGF 的某种简朴的替换。从原理上讲，kineticTheory和g0不应该同时开启。
packingLimiter
Henrik Rusche， PHD Thesis， Computational Fluid Dynamics of Dispersed Two-Phase Flows at High Phase Fractions， Imperial College of Science, Technology & Medicine, Department of Mechanical Engineering, 2002