[发明专利]一种基于CFD-DEM耦合模型的磨粒流场分析方法

摘要

本发明公开了一种基于耦合模型的磨粒流场分析方法，包括以下步骤：(1)建立并求解面向磨粒流加工的流体相控制方程组；(2)根据步骤(1)的求解结果，通过DEM理论求解磨粒运动方程；(3)根据磨粒运动方程，更新磨粒在流道中的位置，并求解新的体积分数以及磨粒受到的作用力；(4)结合新的流体体积分数以及磨粒受到的作用力，将其代入到步骤(1)中开始新一轮的计算；依照上述计算过程，不断循环更新每个磨粒的速度和位置，进而演化整个磨粒流加工系统，得到磨粒流运动规律。本发明将计算流体力学模型与离散元理论模型相结合的方法引入到磨粒流加工领域，能够实现对磨粒运动的精确分析，进而实现对磨粒流加工工艺参数的准确调控。

主权项

一种基于CFD‑DEM耦合模型的磨粒流场分析方法，其特征在于：主要包括如下步骤：(1)建立并求解面向磨粒流加工的流体相控制方程组：在CFD‑DEM耦合策略分析中，磨粒采用DEM理论建模并视为真正的离散磨粒存在，流体相和磨粒之间不存在质量传递；假设由磨粒占据后，流体网格单元中流体相剩余的体积分数为αl，磨粒受到的作用力为磨粒受到的平均作用力为S，由此得到CFD‑DEM耦合框架下的流体相控制方程组为：其中，ρl与分别为流体相的密度及速度；为液相应力‑应变张量，表达式如下：其中，μl以及λl为流体相剪切粘度和体积粘度；为流体相湍动能强度；磨粒受到的平均作用力S表达式如下：其中，Vl为流体单元网格体积；(2)根据步骤(1)的求解结果，通过DEM理论求解磨粒运动方程：磨粒运动求解采用DEM理论求解，在DEM理论求解框架下，每一个磨粒运动在Lagrange坐标系下通过牛顿第二运动定律求解，磨粒运动求解方程为：其中，与分别代表碰撞接触力以及重力，为流体对磨粒的作用力，为磨粒的重力，表达式为：其中mp为磨粒质量；磨粒运动方程中磨粒各项受力求解如下：2.1)流体对磨粒的作用力求解：当磨粒处于高速流体当中，较高的速度梯度会造成磨粒表面的压力差，因此，必须考虑Saffman升力对磨粒的作用；流体对磨粒的驱动作用可用下式表示：式中，为流体对磨粒的曳力，为流体对磨粒的压力梯度力，为流体对磨粒的浮力，为流体对磨粒的Saffman升力；在磨粒流加工中，存在高浓度磨粒区域，因此必须考虑空隙率对曳力的影响，曳力由下式计算：其中，χ＝3.7‑0.65exp[‑(1.5‑logRep)2/2]，式中，ρp为磨粒密度，为磨粒速度以及dp为磨粒直径，Cd为曳力系数，Rep为磨粒雷诺数；曳力系数Cd通过下式计算：Cd=(0.63+4.8Rep0.5)2,]]>磨粒雷诺数Rep通过下式计算：其中μl为流体相动力粘性系数；流体对磨粒的压力梯度力由下式计算：流体对磨粒的浮力由下式计算：流体对磨粒的Saffman升力由下式计算：Vp为磨粒的体积；2.2)磨粒的碰撞接触力求解：由于在流体驱动下，磨粒不断冲击工件，依靠磨粒与工件的碰撞接触力实现工件的光整加工，因此采用软球模型精确计算碰撞接触力，并采用弹性系数及阻尼系数来量化弹簧、阻尼器、滑动器的作用；接触力可以分解为法向接触力和切向接触力，假设弹性系数及阻尼系数分别为k及η，则法向接触力可由下式计算：式中，α为磨粒接触时的法向重叠量；为磨粒i相对于磨粒j的速度，为磨粒i球心位置到磨粒j球心位置的单位矢量；和分别为磨粒i和磨粒j的速度；kn是磨粒的法向弹性系数，其计算公式为：kn=43(1-υi2Ei+1-υj2Ej)-1(ai+ajaiaj)-1/2,]]>ηn是磨粒的法向阻尼系数，其计算公式为：ηn=2mkn,]]>式中，E和υ是磨粒材料的弹性模量和泊松比，a是磨粒半径；切向接触力可由下式计算：式中，kt和ηt是切向弹性系数和切向阻尼系数；δ是接触点的切向位移；是接触点的滑移速度，可由下式计算：式中，ai和aj分别是磨粒i和磨粒j的半径；和分别是磨粒i和磨粒j的角速度；kt是磨粒的切向弹性系数，由下式计算：kt=8α1/2(1-υi2Gi+1-υj2Gj)-1(ai+ajaiaj)-1/2]]>ηt是磨粒的切向阻尼系数，由下式计算：ηt=2mkt,]]>综合考虑法向力和切向力的作用，磨粒i所受到的合力为：磨粒i所受到的力矩为：磨粒流加工过程中，多个磨粒可能同时发生接触，此时，磨粒i所受到的合力为：磨粒i所受到的合力矩为：(3)根据步骤(2)中求得的磨粒运动方程，更新磨粒在流道中的位置，并求解新的体积分数αp以及磨粒受到的作用力αp计算过程如下：首先在每一个磨粒周围采用边界盒将其包围起来；之后在每一个边界盒内规则的抽取样本点，如果样本点位于磨粒表面内，则将其保存；最后，对每一个样本点进行检验，确认样本点是否位于网格单元内；因此，某个网格单元内的固体体积分数就是网格单元内的样本点数量在所有样本点数量中所占的比例，由下式计算：αp=1-αl=ΣparticlesncNVparticle,]]>其中，Vparticle为磨粒体积，nc为网格单元中样本点数量，N为边界盒内样本点总数量；(4)结合步骤(3)中新的流体体积分数αl以及磨粒受到的作用力将流体体积分数αl以及磨粒受到的作用力代入到步骤(1)中开始新一轮的计算；依照步骤(1)、步骤(2)和步骤(3)的计算过程，不断的循环更新每个磨粒的速度和位置，进而确定性的演化整个磨粒流加工系统，得到磨粒流运动规律。