CN104915471A - 基于离散元法的结构表面磨损仿真方法 - Google Patents

Abstract

本发明提供了一种基于离散元法的结构表面磨损仿真方法，包括如下步骤：S1：获取结构表面形状，并将目标表面划分为由若干三角单元组成的三角网格；S2：对三角单元进行离散元仿真，仿真其被散体颗粒碰撞的过程，从而计算并记录各三角单元的碰撞能量数据；S3：依据步骤S2得到的所述碰撞能量数据计算各三角单元的表面上剥落的材料体积，从而依据该材料体积得到各三角单元质心的高度变化量；S4：依据步骤S3得到的各三角单元质心的高度变化量重新拟合所述目标表面。本发明基于离散元法的结构表面磨损仿真方法为散体物料处理设备的结构设计、功能设计与寿命设计等提供了参考。

Description

基于离散元法的结构表面磨损仿真方法

技术领域

本发明涉及散体无聊处理机械结构仿真设计领域，尤其涉及一种基于离散元法的结构表面磨损仿真方法。

背景技术

散体物料处理机械的主要失效形式之一就是结构表面磨损。表面磨损专指由大量砂砾、岩土、尘埃或类似散体颗粒频繁碰撞摩擦机械结构，并造成其表面材料体积剥落的现象。

结构部件的表面磨损机理复杂，不仅与物料自身的材料和工况等因素有关，而且机械设备内散体颗粒与结构的碰撞过程与碰撞能量有紧密的联系，所以目前对结构表面磨损的理论研究和实验研究都存在较高的技术难度。

发明内容

本发明要解决的技术问题是如何实现结构表面磨损的仿真。

为了解决这一技术问题，本发明提供了一种基于离散元法的结构表面磨损仿真方法，包括如下步骤：

S1：获取结构表面形状，并将目标表面划分为由若干三角单元组成的三角网格；

S2：对三角单元进行离散元仿真，仿真其被散体颗粒碰撞的过程，从而计算并记录各三角单元的碰撞能量数据；

S3：依据步骤S2得到的所述碰撞能量数据计算各三角单元的表面上剥落的材料体积，从而依据该材料体积得到各三角单元质心的高度变化量；

S4：依据步骤S3得到的各三角单元质心的高度变化量重新拟合所述目标表面。

可选的，在所述步骤S1中，所述三角单元为平直不可弯曲的区域面。

可选的，在所述步骤S2中，在离散元仿真时，将被散体颗粒碰撞的过程区分为若干时间步长进行仿真，记录下各个仿真时步下散体颗粒在三角单元上产生的接触信息，并根据该接触信息计算并记录各三角单元的碰撞能量数据。

可选的，在所述步骤S2中，一个时间步长内一个三角单元的碰撞能量数据包括接触阻尼能量消耗E_damp和摩擦能量消耗E_slip，分别通过以下公式计算：

$\left\{\begin{array}{c}{E}_{\mathrm{damp}}=\underset{c}{\mathrm{\Σ}}[\left({\mathrm{\η}}_n\right|{\stackrel{\·}{u}}_n\left|{\mathrm{\Δx}}_n\right)+\left({\mathrm{\η}}_t\right|{\stackrel{\·}{u}}_t\left|{\mathrm{\Δx}}_t\right)]\\ E_{\mathrm{slip}}=\underset{c}{\mathrm{\Σ}}\left[\right|F_t|\×{\mathrm{\Δx}}_t]\end{array}\right.$

其中，c是某个接触信息的索引号，η_n和η_t分别是法向和切向阻尼系数，和分别是颗粒在接触过程中的法向速度和切向速度，Δx_n和Δx_t分别是颗粒在接触过程中的法向和切向位移量，F_t是接触面上的切向力。

可选的，在所述步骤S3中，通过以下公式计算得到第i个三角单元的表面上剥落的材料体积：

V_i＝δ(α×E_slip,i+β×E_damp,i)

其中，δ是破坏率，单位为m³/J，α是摩擦、胶合作用的影响系数，β是颗粒碰撞过程中阻尼作用的影响系数；E_damp,i是该三角单元的接触阻尼能量消耗和E_{slip， i}是该三角单元的摩擦能量消耗。

可选的，在所述步骤S3中，三角单元质心的坐标P_i在磨损方向n’上新的高度值h_i由以下公式获得：

h_i＝h_i0-ΔH_i

其中，h_i0为初始高度；ΔH_i由以下公式获得：

${\mathrm{\&Delta;H}}_i=\frac{V_i}{S_i\&CenterDot;\cos (<P_i,n^{\&prime;}>)}$

其中，S_i为三角单元的面积，<P_i,n’>表示坐标P_i的向量和向量n’间的夹角。

可选的，在所述步骤S4中，重新拟合所述目标表面后，还包括对重新拟合后的目标表面进行平滑处理的过程。

可选的，在进行平滑处理时，使得：

(1)目标表面上两个连续坡度的斜率之差最小；

(2)两个相邻三角单元的质心坐标在磨损方向上的高度变化量之差最小；

(3)各三角单元的质心坐标在磨损方向上的高度随时间递减。

本发明将离散元法应用于散体物料处理设备的表面磨损仿真问题。具体来说，依据利用离散元法将散体颗粒与结构间的碰撞能量和磨损量建立定量关系，并用可视化的流程展现出结构CAD模型的表面磨损特征。在离散元法中，结构的CAD模型和散体颗粒均假设是刚体，所以结构的弹塑性形变不作考虑，结构的表面磨损问题在于如何定量地计算刚性表面的磨损量，并以可视化的方式展现出受磨损结构表面的几何特征。基于离散元法的结构表面磨损仿真方法为散体物料处理设备的结构设计、功能设计与寿命设计等提供了参考。

附图说明

图1是本发明一实施例中工业水泥磨机的仿真模型示意图；

图2和图3是本发明一实施例中衬板的仿真模型示意图；

图4是本发明一实施例中三角网格示意图；

图5是未磨损表面效果示意图；

图6是磨损后表面效果示意图。

具体实施方式

以下将结合图1至图6，通过一个实施例对本发明提供的基于离散元法的结构表面磨损仿真方法进行详细的描述，其为本发明可选的实施例，可以认为，本领域的技术人员在不改变本发明精神和内容的范围内，能够对其进行修改和润色。

本发明提供了一种基于离散元法的结构表面磨损仿真方法，包括如下步骤：S1：获取结构表面形状，并将目标表面划分为由若干三角单元组成的三角网格；

在所述步骤1中，在结构CAD模型上选择目标表面，作为磨损仿真的主要研究区域，目标表面将被划分成由三角单元组成的网格表面，其中三角单元基于以下定义：

所述三角单元是组成网格表面的最小单元，其来源包括标准的基于三角网格的三维模型文件；三角单元是平直不可弯曲的区域面，其数量和尺寸决定了网格表面的几何精度和磨损表面的平滑程度。

S2：对三角单元进行离散元仿真，仿真其被散体颗粒碰撞的过程，从而计算并记录各三角单元的碰撞能量数据；

在所述步骤S2中，在离散元仿真时，将被散体颗粒碰撞的过程区分为若干时间步长进行仿真，记录下各个仿真时步下散体颗粒在三角单元上产生的接触信息，并根据该接触信息计算并记录各三角单元的碰撞能量数据。具体来说，在所述步骤2中，离散元仿真程序计算并记录各个仿真时步下，散体颗粒在三角单元上产生的所有接触信息，包括颗粒-三角单元的接触点位置、接触力、接触过程中颗粒的运动速度和位移等。

在步骤S2中，一个时间步长内结构表面的总体碰撞能量E_impact由下式计算：

E_impact＝E_damp+E_slip

在所述步骤S2中，一个时间步长内一个三角单元的碰撞能量数据包括接触阻尼能量消耗E_damp和摩擦能量消耗E_slip，分别通过以下公式计算：

$\left\{\begin{array}{c}{E}_{\mathrm{damp}}=\underset{c}{\mathrm{\&Sigma;}}[\left({\mathrm{\&eta;}}_n\right|{\stackrel{\&CenterDot;}{u}}_n\left|{\mathrm{\&Delta;x}}_n\right)+\left({\mathrm{\&eta;}}_t\right|{\stackrel{\&CenterDot;}{u}}_t\left|{\mathrm{\&Delta;x}}_t\right)]\\ E_{\mathrm{slip}}=\underset{c}{\mathrm{\&Sigma;}}\left[\right|F_t|\&times;{\mathrm{\&Delta;x}}_t]\end{array}\right.$

其中，c是某个接触信息的索引号，η_n和η_t分别是法向和切向阻尼系数，它们是仿真的输入参数，和分别是颗粒在接触过程中的法向速度和切向速度，Δx_n和Δx_t分别是颗粒在接触过程中的法向和切向位移量，F_t是接触面上的切向力，这些数据均可由离散元仿真程序直接得到。

S3：依据步骤S2得到的所述碰撞能量数据计算各三角单元的表面上剥落的材料体积，从而依据该材料体积得到各三角单元质心的高度变化量；

在所述步骤S3中，通过以下公式计算得到第i个三角单元的表面上剥落的材料体积：

V_i＝δ(α×E_slip,i+β×E_damp,i)

其中，δ是破坏率，单位为m³/J，α是摩擦、胶合作用的影响系数，β是颗粒碰撞过程中阻尼作用的影响系数，这两个参数假设与材料的特性相关，通常取β＝10α；E_damp,i是该三角单元的接触阻尼能量消耗和E_slip,i是该三角单元的摩擦能量消耗。

在所述步骤S3中，三角单元质心的坐标P_i在磨损方向n’上新的高度值h_i由以下公式获得：

h_i＝h_i0-ΔH_i

其中，h_i0为初始高度；ΔH_i由以下公式获得：

${\mathrm{\&Delta;H}}_i=\frac{V_i}{S_i\&CenterDot;\cos (<P_i,n^{\&prime;}>)}$

其中，S_i为三角单元的面积，<P_i,n’>表示坐标P_i的向量和向量n’间的夹角。

换言之，在所述步骤3中，目标表面上各三角单元的质心坐标P_i在磨损方向n’上的初始高度h_i0＝P_i·n’，受磨损后的结构表面，P_i在磨损方向上新的高度值h_i＝h_i0-△H_i，△H_i由下式计算：

${\mathrm{\&Delta;H}}_i=\frac{V_i}{S_i\&CenterDot;\cos (<P_i,n^{\&prime;}>)}$

其中S_i为三角单元的面积，<P_i,n’>表示向量P_i和向量n’间的夹角。

S4：依据步骤S3得到的各三角单元质心的高度变化量重新拟合所述目标表面。

在所述步骤S4中，重新拟合所述目标表面后，还包括对重新拟合后的目标表面进行平滑处理的过程。

在进行平滑处理时，使得：

(1)目标表面上两个连续坡度的斜率之差最小；

(2)两个相邻三角单元的质心坐标在磨损方向上的高度变化量之差最小；

(3)各三角单元的质心坐标在磨损方向上的高度随时间递减。

以下以工业水泥磨机进行进一步地举例，应明确，下文所阐述的内容包含上文所阐述方案所有的特征。只是以更具体的方式进行阐述。

本实施例中，图1为的工业水泥磨机离散元仿真模型，该模型将磨机的筒体1抽象为具有材料特性的圆柱体，其内壁上固定有一定数量的平面阶梯衬板2，平面阶梯衬板2的尺寸如图3所示。为了节省仿真计算时间，截取其中长度为0.5m的一段圆柱体作为求解边界，转速2.59rad/s。筒体1内的研磨介质3抽象为具有材料特性与尺寸参数的球形颗粒单元，球径随机分布于50～90mm之间，介质填充率为18％。

考虑衬板布置的对称性，取其中一块衬板作为目标衬板进行表面磨损仿真分析，而其它衬板将具有与之相似的表面磨损特征，如图4所示。

首先，将目标衬板的两个工作平面作为目标表面，并将其划分成286片三角单元，工作平面磨损方向n’＝-n，n是目标平面的法向量，划分的各三角单元的质心坐标P_i在磨损方向上的初始高度h_i0＝P_i·n’。

然后，开始执行离散元仿真计算。当磨机筒体旋转一周后，筒内的研磨介质达到稳定运动状态，在此过程中，离散元仿真程序计算工作平面的三角单元上发生的所有碰撞能量，并作记录。一个时间步长，工作平面的总体碰撞能量E_impact采用如下公式计算：

E_impact＝E_damp+E_slip

其中E_damp是颗粒接触结构表面过程中的阻尼能量消耗，E_slip是颗粒相对结构表面滑移的摩擦能量消耗，它们分别可以由下式计算得到：

$\left\{\begin{array}{c}{E}_{\mathrm{damp}}=\underset{c}{\mathrm{\&Sigma;}}[\left({\mathrm{\&eta;}}_n\right|{\stackrel{\&CenterDot;}{u}}_n\left|{\mathrm{\&Delta;x}}_n\right)+\left({\mathrm{\&eta;}}_t\right|{\stackrel{\&CenterDot;}{u}}_t\left|{\mathrm{\&Delta;x}}_t\right)]\\ E_{\mathrm{slip}}=\underset{c}{\mathrm{\&Sigma;}}\left[\right|F_t|\&times;{\mathrm{\&Delta;x}}_t]\end{array}\right.$

其中c是某个接触信息的索引号，η_n和η_t分别是法向和切向阻尼系数，它们是仿真的输入参数，和分别是颗粒在接触过程中的法向速度和切向速度，△x_n和△x_t分别是颗粒在接触过程中的法向和切向位移量，F_t是接触面上的切向力，这些数据均可由离散元仿真程序直接得到。

将E_damp和E_slip的数值代入下式：

V_i＝δ(α×E_slip,i+β×E_damp,i)

可以求得从三角单元i表面上剥落材料的体积，此处磨损系数δ＝1.5×10^-7m³/J，阻尼能量权重系数α＝1，摩擦能量权重系数β＝0.1。

三角单元的质心坐标在磨损方向上新的高度值h_i＝h_i0-△H_i，△H_i可以由下式计算：

${\mathrm{\&Delta;H}}_i=\frac{V_i}{S_i\&CenterDot;\cos (<P_i,n^{\&prime;}>)}$

其中S_i为三角单元的面积，<P_i,n’>表示向量P_i和向量n’间的夹角。

最后，根据各三角单元的h_i值，在n’方向上拟合磨损表面，设定以下平滑目标：(1)目标平面上两个连续坡度的斜率之差最小；(2)两个相邻三角单元的质心坐标在磨损方向上的高度变化量之差最小；(3)各三角单元的质心坐标在磨损方向上的高度随时间递减。最终，经平滑处理后的磨损表面轮廓，如图5和图6所示。

衬板在研磨介质的冲击和磨损的作用下，在工作表面上形成了沟槽，这些沟槽是球磨机衬板的典型失效形式。参考图5和图6显示衬板的工作表面形成了沟槽，若将其与实际的残体照片对比，可以说明本发明提出的一种基于离散元法的结构表面磨损仿真方法能够准确有效地展现由大量散体颗粒在结构表面上形成的主要磨损特征。

综上所述，本发明将离散元法应用于散体物料处理设备的表面磨损仿真问题。具体来说，依据利用离散元法将散体颗粒与结构间的碰撞能量和磨损量建立定量关系，并用可视化的流程展现出结构CAD模型的表面磨损特征。在离散元法中，结构的CAD模型和散体颗粒均假设是刚体，所以结构的弹塑性形变不作考虑，结构的表面磨损问题在于如何定量地计算刚性表面的磨损量，并以可视化的方式展现出受磨损结构表面的几何特征。基于离散元法的结构表面磨损仿真方法为散体物料处理设备的结构设计、功能设计与寿命设计等提供了参考。

Claims (8)

1.一种基于离散元法的结构表面磨损仿真方法，其特征在于：包括如下步骤：

S1：获取结构表面形状，并将目标表面划分为由若干三角单元组成的三角网格；

S2：对三角单元进行离散元仿真，仿真其被散体颗粒碰撞的过程，从而计算并记录各三角单元的碰撞能量数据；

S3：依据步骤S2得到的所述碰撞能量数据计算各三角单元的表面上剥落的材料体积，从而依据该材料体积得到各三角单元质心的高度变化量；

S4：依据步骤S3得到的各三角单元质心的高度变化量重新拟合所述目标表面。

2.如权利要求1所述的基于离散元法的结构表面磨损仿真方法，其特征在于：在所述步骤S1中，所述三角单元为平直不可弯曲的区域面。

3.如权利要求1所述的基于离散元法的结构表面磨损仿真方法，其特征在于：在所述步骤S2中，在离散元仿真时，将被散体颗粒碰撞的过程区分为若干时间步长进行仿真，记录下各个仿真时步下散体颗粒在三角单元上产生的接触信息，并根据该接触信息计算并记录各三角单元的碰撞能量数据。

4.如权利要求3所述的基于离散元法的结构表面磨损仿真方法，其特征在于：在所述步骤S2中，一个时间步长内一个三角单元的碰撞能量数据包括接触阻尼能量消耗E_damp和摩擦能量消耗E_slip，分别通过以下公式计算：

$\left\{\begin{array}{c}{E}_{\mathrm{damp}}=\underset{c}{\mathrm{\&Sigma;}}[\left({\mathrm{\&eta;}}_n\right|{\stackrel{\&CenterDot;}{u}}_n\left|\mathrm{\&Delta;}{x}_n\right)+\left({\mathrm{\&eta;}}_t\right|{\stackrel{\&CenterDot;}{u}}_t\left|\mathrm{\&Delta;}{x}_t\right)]\\ E_{\mathrm{slip}}=\underset{c}{\mathrm{\&Sigma;}}\left[\right|F_t|\&times;\mathrm{\&Delta;}{x}_t]\end{array}\right.$

其中，c是某个接触信息的索引号，η_n和η_t分别是法向和切向阻尼系数，和分别是颗粒在接触过程中的法向速度和切向速度，Δx_n和Δx_t分别是颗粒在接触过程中的法向和切向位移量，F_t是接触面上的切向力。

5.如权利要求1所述的基于离散元法的结构表面磨损仿真方法，其特征在于：在所述步骤S3中，通过以下公式计算得到第i个三角单元的表面上剥落的材料体积：

V_i＝δ(α×E_slip,i+β×E_damp,i)

其中，δ是破坏率，单位为m³/J，α是摩擦、胶合作用的影响系数，β是颗粒碰撞过程中阻尼作用的影响系数；E_damp,i是该三角单元的接触阻尼能量消耗和E_slip,i是该三角单元的摩擦能量消耗。

6.如权利要求1所述的基于离散元法的结构表面磨损仿真方法，其特征在于：在所述步骤S3中，三角单元质心的坐标P_i在磨损方向n′上新的高度值h_i由以下公式获得：

h_i＝h_i0-ΔH_i

其中，h_i0为初始高度；ΔH_i由以下公式获得：

${\mathrm{\&Delta;H}}_i=\frac{V_i}{S_i\&CenterDot;\cos (<P_i,n^{\&prime;}>)}$

其中，S_i为三角单元的面积，<P_i,n′>表示坐标P_i的向量和向量n′间的夹角。

7.如权利要求1所述的基于离散元法的结构表面磨损仿真方法，其特征在于：在所述步骤S4中，重新拟合所述目标表面后，还包括对重新拟合后的目标表面进行平滑处理的过程。

8.如权利要求7所述的基于离散元法的结构表面磨损仿真方法，其特征在于：在进行平滑处理时，使得：

(1)目标表面上两个连续坡度的斜率之差最小；

(2)两个相邻三角单元的质心坐标在磨损方向上的高度变化量之差最小；

(3)各三角单元的质心坐标在磨损方向上的高度随时间递减。