如果用DPMFoam求解稀相流会怎么样？误差大么

[Alvin]

@dzw05 参照 多相流数学模型 http://dyfluid.com/docs/multiphase.html 3.3 连续相模型之欧拉模型（宏观模型）：

“从可压缩动量方程出发，假设密度是常数，那可以将密度提出来，就是现在采用的不可压缩动量方程形式啊，为什么会默认密度为1呢？”
假设密度是常数，那么即使将密度提出来，DPMFoam求解器中的动量方程中压力项还少除以了一个密度（除非密度是1）

“具体来说在Openfoam中，particleForces中的pressureGradient这项才是浮力项。”
颗粒所受的浮力在颗粒所受重力项中已经体现（）： pressureGradient就是单纯的压力梯度力，但是大小和浮力相反（可参考颗粒在流体中的受力分析 https://wenku.baidu.com/view/171125dcd15abe23482f4d93.html、董长银-牛顿流体中的固体颗粒运动模型分析及应用https://wenku.baidu.com/view/8fafedfef705cc175527097a.html）

“DPMFoam动量方程中的重力项有一部分是流体的重力，另一部分和pressureGradient这一项组合起来作为浮力。” 您是说连续相的动量方程中的重力项要被拆分吗？据我调查，颗粒所受重力中已经考虑了流体的浮力影响，那么，颗粒对流体应该有相反的作用力（体积力的源项），那么连续相动量方程里的重力项似乎可以解释了，是作为颗粒的耦合项加入的。

@东岳 “pressureGradient是粒子上下表面引致的压力差项，粒子若非常小，则可忽略。但若重力gravity若被忽略，则完全没有浮力影响。” 颗粒浮力的影响通过连续相的密度以及重力加速度体现出来，关于https://wenku.baidu.com/view/8fafedfef705cc175527097a.html文献中提到的考虑了压力梯度力后，浮力项不起作用了，因为两者大小相等，方向相反，在openfoam中还没有验证，应该根据具体案例和压力分布相关

[Alvin]

@dzw05
“没错，在离散相的gravity中，表达式是massg_(1.0 - p.rhoc()/p.rho())，即考虑了连续相和离散相的密度比，这个是浮力的%(#c43d3d)[主要来源]。” 浮力还有哪些其它来源呢？
“如果单纯的将连续相动量方程中的重力g删除，而不是设为0，那么其实并没有影响粒子浮力的计算，因为粒子运动方程中根本没有出现连续相的压力。” 粒子运动方程中出现的连续相的压力，不是可以通过压力梯度力来体现出来了吗？这样一来，连续相动量方程中的耦合项或者压力项需做相应调整了。

@东岳 您提出存在下面的两套方程，是等效的吗？

[李东岳]

这俩套方程都是文献中看到的，可以参考本帖24楼。我自己的求解器中采用这种形式

$$\begin{array}{}\text{(1)}& m_{\mathrm{d}\mathrm{p}\mathrm{m}}\frac{\mathrm{\partial }{\mathbf{U}}_{\mathrm{d}\mathrm{p}\mathrm{m}}}{\mathrm{\partial }t}={\mathbf{F}}_{\mathrm{d}\mathrm{r}\mathrm{a}\mathrm{g}}+m_{\mathrm{d}\mathrm{p}\mathrm{m}}\mathbf{g}+\frac{m_{\mathrm{d}\mathrm{p}\mathrm{m}}}{{\rho }_{\mathrm{d}\mathrm{p}\mathrm{m}}}\mathrm{\nabla }{p}_{\mathrm{c}}+\mathrm{O}\mathrm{t}\mathrm{h}\mathrm{e}\mathrm{r}\mathrm{s}\end{array}$$

俩种形式可以互相转换。依据$p=p_0+{\rho }_{\mathrm{c}}\mathbf{g}\mathbf{h}$即可相互推导

[dzw05]

@alvin @东岳 贴出来的两套方程才是问题的根源所在。目前OpenFOAM中的压力梯度项并没有出现连续相的压力，而是所谓的rhoc*Du/Dt。

@alvin 不可压缩动量方程中的p，其实是p/rho，这个和可压缩方程中是不同的，simpleFoam中也是这样处理的啊。难道咱们讨论的方程不是同一个吗？至少采用DPMFoam计算液-固流动时，连续相的压力是没有问题的。

[李东岳]

$$\begin{array}{}\text{(2)}& m_{\mathrm{d}\mathrm{p}\mathrm{m}}\frac{\mathrm{d}{\mathbf{U}}_{\mathrm{d}\mathrm{p}\mathrm{m}}}{\mathrm{d}t}=-\frac{m_{\mathrm{d}\mathrm{p}\mathrm{m}}}{{\rho }_{\mathrm{d}\mathrm{p}\mathrm{m}}}\mathrm{\nabla }p+m_{\mathrm{d}\mathrm{p}\mathrm{m}}\mathbf{g}+\frac{m_{\mathrm{d}\mathrm{p}\mathrm{m}}}{{\rho }_{\mathrm{d}\mathrm{p}\mathrm{m}}}\frac{3}{4}\frac{C_{\mathrm{D}}{\rho }_{\mathrm{c}}}{d_{\mathrm{d}}}|{\mathbf{U}}_{\mathrm{c}}-{\mathbf{U}}_{\mathrm{d}\mathrm{p}\mathrm{m}}|({\mathbf{U}}_{\mathrm{c}}-{\mathbf{U}}_{\mathrm{d}\mathrm{p}\mathrm{m}})+\mathrm{O}\mathrm{t}\mathrm{h}\mathrm{e}\mathrm{r}\mathrm{F}\mathrm{o}\mathrm{r}\mathrm{c}\mathrm{e}\mathrm{s}.\end{array}$$
伯努利：
$$\begin{array}{}\text{(3)}& \mathrm{\nabla }p={\rho }_{\mathrm{c}}\mathbf{g}\end{array}$$
代入：
$$\begin{array}{}\text{(4)}& m_{\mathrm{d}\mathrm{p}\mathrm{m}}\frac{\mathrm{d}{\mathbf{U}}_{\mathrm{d}\mathrm{p}\mathrm{m}}}{\mathrm{d}t}=m_{\mathrm{d}\mathrm{p}\mathrm{m}}\mathbf{g}(1-\frac{{\rho }_{\mathrm{c}}}{{\rho }_{\mathrm{d}\mathrm{p}\mathrm{m}}})+\frac{m_{\mathrm{d}\mathrm{p}\mathrm{m}}}{{\rho }_{\mathrm{d}\mathrm{p}\mathrm{m}}}\frac{3}{4}\frac{C_{\mathrm{D}}{\rho }_{\mathrm{c}}}{d_{\mathrm{d}}}|{\mathbf{U}}_{\mathrm{c}}-{\mathbf{U}}_{\mathrm{d}\mathrm{p}\mathrm{m}}|({\mathbf{U}}_{\mathrm{c}}-{\mathbf{U}}_{\mathrm{d}\mathrm{p}\mathrm{m}})+\mathrm{O}\mathrm{t}\mathrm{h}\mathrm{e}\mathrm{r}\mathrm{F}\mathrm{o}\mathrm{r}\mathrm{c}\mathrm{e}\mathrm{s}.\end{array}$$

[Alvin]

@dzw05 抱歉，理解有限，“ 目前OpenFOAM中的压力梯度项并没有出现连续相的压力，而是所谓的rhoc*Du/Dt。” 这句话不是很理解.
“不可压缩动量方程中的p，其实是p/rho，这个和可压缩方程中是不同的，simpleFoam中也是这样处理的啊。”从量纲上确实是这样，查看过simpleFoam的方程，代码中没发现有p/rho这样的处理，您能从源码角度讲一下它的实现过程吗（源码位置？），希望不吝赐教？

@东岳 再次谢谢老师的指导

[dzw05]

@alvin （1）rhoc*Du/Dt，这个是在src\lagrangian\intermediate\submodels\Kinematic \ParticleForces \PressureGradient\PressureGradientForce.T.C代码中出现的，应该对应于东岳写的第一套方程。
（2）至于压力这个，源码确实没有出现，但推导的时候就是这样出来的。simpleFoam计算输出的p也不是真正的压力，而是p/rho，你要看真实的压力的话，应该乘以连续相的密度。

[Alvin]

@东岳 这样说来，利用现有的DPMFoam（3.0.1版本）计算气液两相流，若气泡同时考虑压力梯度力和重力（massg_(1.0 - p.rhoc()/p.rho())），可能会导致气泡浮力失效？

@dzw05 谢谢您，我再看一下。“至于压力这个，源码确实没有出现，但推导的时候就是这样出来的。simpleFoam计算输出的p也不是真正的压力，而是p/rho，你要看真实的压力的话，应该乘以连续相的密度。”
谢谢您的解答。

@东岳 @dzw05 离散相方程中的压力梯度力确实是rhocDu / Dt这种形式，和浮力（rhocg）是两个概念

[李东岳]

@Alvin 不知道你的问题解决没有。我在我自己写的求解器中，通过

1. 将p方程换为p_rgh方程
2. 将g去掉

解决了。因为轻颗粒的浮力在方程$\text{(???)}$中包含了，重颗粒也如此。就像 @dzw05 说的，连续相可以认为是单向不可压缩流，类似simpleFoam。 @Alvin 你试试将DPMFoam中的phig去掉，运算看看？

[Alvin]

@东岳 问题解决了，如您所说，去掉动量方程中的重力加速度项，案例结果就合理了。我在考虑，动力方程中的重力项的物理含义是什么，颗粒受到周围流体对它施加的浮力，反之，根据牛顿第三定律，流体也应该受到对应颗粒浮力的反作用力，方向相反，这个重力项是否就是这个反作用力，尚未查到相关文献，只是一种猜测。

[李东岳]

@alvin DPMFoam中还有个小bug，不知道会不会影响你的结果

http://www.cfd-china.com/topic/969/mppicinterfoam中的bug/18

[Alvin]

@ 东岳 您方便说明一下“......我验证了双向耦合中界面交换相、相分数的耦合作用大小。在我测试的算例中，相分数的耦合作用较小......”测试过程，修改了哪些参数？ 相分数的耦合作用较小是通过和哪个对比的？

[李东岳]

多相的耦合主要在于相分数耦合${\alpha }_{\mathrm{d}}$和动量界面交互力$\mathbf{M}$的耦合。要验证相互作用的大小，要么把$\mathbf{M}$设为0计算，要么把${\alpha }_{\mathrm{d}}$设为1计算。看二者的相对重要程度。

[Alvin]

@东岳 再次感谢您的耐心讲解，受教

[yfclark]

@dzw05 @东岳 从大家的讨论看是不是把重力设为0之后，颗粒相控制方程中的压力梯度和重力项都变成了0，还有那个ducdt我也看了一下好像和公式（6）也对不上啊

[李东岳]

我写的那个公式6？openfoam里面压力梯度项和我那个不一样，用的你说的DUc/Dt

[myjimp]

@李东岳 李老师,这样修改DPMFoam可以计算稀相的吗？

[李东岳]

@myjimp 可以，没问题

[dxl]

@alvin 请问这个重力加速度项是 phicForces里面的 rAUcf*(g & mesh.Sf()) 吗

[dxl]

@李东岳 在 如果用DPMFoam求解稀相流会怎么样？误差大么 中说：

@Alvin 不知道你的问题解决没有。我在我自己写的求解器中，通过

1. 将p方程换为p_rgh方程
2. 将g去掉

解决了。因为轻颗粒的浮力在方程$\text{(???)}$中包含了，重颗粒也如此。就像 @dzw05 说的，连续相可以认为是单向不可压缩流，类似simpleFoam。 @Alvin 你试试将DPMFoam中的phig去掉，运算看看？

这个1,2两项是二选一即可吗，phig指的是不是phicForces里面的rAUcf*（g&mesh.Sf()）