CN113297781A - 一种基于内作用力的颗粒簇建模方法 - Google Patents

Abstract

本发明提供一种基于内作用力的颗粒簇建模方法，包括以下步骤：步骤10)在边界内初始化生成n个位置随机的元素；n表示大于等于2的整数；步骤20)依次计算每个元素所受的其它元素的作用力以及所受的边界的作用力，调整元素的位置；步骤30)判断是否满足预设条件，若不满足，则执行步骤20)；步骤40)得到颗粒簇模型。本发明基于内作用力的颗粒建模方法，计算元素间以及元素与边界间的作用力，通过迭代的方式不断调整元素的位置，以尽可能减小簇内的总作用力，从而使簇内的元素拟合任意边界形状时尽可能消除内部的重叠量。

Description

一种基于内作用力的颗粒簇建模方法

技术领域

本发明涉及元素颗粒簇仿真建模技术领域，具体涉及一种基于内作用力的颗粒簇建模方法。

背景技术

颗粒系统是一定数量颗粒的集合，在自然界、工业生产和日常生活中都十分常见。这些颗粒系统的种类十分多样，从宏观尺度到微观尺度处处都有着颗粒的存在。颗粒的物理性质非常复杂，不仅具有固体和液体的双重特性，而且有许多其他特殊性质。这些特殊性决定了颗粒系统在宏观上表现的复杂又特别，如筒仓效应、振动分离特性、颗粒的自分离现象等。

随着计算机技术的高速发展，越来越多的模型已通过仿真模拟实现。离散元素法是一种能够描述颗粒系统力学行为的数值计算方法。近几年离散元素法被广泛地应用于颗粒系统的仿真分析与研究，展现出了其较好的模拟精确度与计算性能。在使用离散元进行仿真的过程中，对物料进行建模的步骤是不可或缺的。其中，使用粘结颗粒簇建模方法生成的颗粒簇可以合适地仿真可破碎或弹塑性材料。因此，如何建模颗粒簇就成了一个至关重要的问题。颗粒簇建模的目的是确定颗粒簇内部中元素颗粒的位置分布，从而用大量元素颗粒拟合任意形状的颗粒簇。

目前，主要采用圆盘(二维中)或圆球(三维中)的元素颗粒来建立任意形状的颗粒簇模型。对于任意形状的颗粒建模，通常有两类方法：颗粒填充和颗粒挤压。而在这些方法中，基本元素之间是可重叠的。当仿真过程中发生破碎情况时，该方法无法避免元素间由于初始存在重叠而产生的内应力，该内应力的存在会夸大破碎时应有的使碎片分离的力，而对仿真结果造成的显著误差。因此，需要提出一种颗粒簇建模方法，使得簇内的元素颗粒在拟合任意边界形状时尽可能消除内部的重叠量。

发明内容

本发明要解决的技术问题是：提出一种基于内作用力的颗粒簇建模方法，可以使得簇内的元素颗粒在拟合任意边界形状时减小颗粒间的重叠量。

为解决上述技术问题，本发明实施例提供一种基于内作用力的颗粒簇建模方法，包括以下步骤：

步骤10)在边界内初始化生成n个位置随机的元素；n表示大于等于2的整数；

步骤20)依次计算每个元素所受的其它元素的作用力以及所受的边界的作用力，调整元素的位置；

步骤30)判断是否满足预设条件，若不满足，则执行步骤20)；

步骤40)得到颗粒簇模型。

作为本发明实施例的进一步改进，所述元素为圆球或圆盘。

作为本发明实施例的进一步改进，所述步骤20)具体包括：

步骤201)根据元素与其它元素之间的位置，得到元素所受的其它元素的作用力；

步骤202)根据元素与边界之间的位置，得到元素所受的边界的作用力；

步骤203)结合元素所受的其它元素的作用力和元素所受的边界的作用力，得到元素的调整偏移量；

步骤204)根据所述元素的调整偏移量，调整元素的位置。

作为本发明实施例的进一步改进，所述步骤201)中，利用式(1)计算元素所受的其它元素的作用力：

式中，

表示元素i所受的元素j的作用力，r_i表示元素i的半径，r_j表示元素j的半径，O_i表示元素i的圆心或球心，O_j表示元素j的圆心或球心，d表示元素i和元素j之间的圆心距离或球心距离。

作为本发明实施例的进一步改进，所述步骤202)中，利用式(2)计算元素所受的边界的作用力：

式中，

表示元素i所受的来自边界点P_k的作用力，k_s为弹性系数，k_c为调节系数，O_i表示元素i的圆心或球心，r_i表示元素i的半径，O₀表示边界的形心，Q表示元素i的球心与边界点P_k的连线与元素i的远离O₀的交点。

作为本发明实施例的进一步改进，所述步骤203)中，利用式(3)计算元素的调整偏移量：

式中，

表示元素i的调整偏移量，

表示元素i所受的其它所有元素的作用合力，

表示元素i所受的边界点P_k的作用力。

作为本发明实施例的进一步改进，所述步骤203)中，利用式(4)计算元素的调整偏移量：

式中，

表示元素i的调整偏移量，

表示元素i所受的其它所有元素的作用合力，

表示元素i所受的边界点P_k的作用力,

表示元素i上一次的调整偏移量。

作为本发明实施例的进一步改进，所述步骤30)中，所述预设条件为：

执行步骤20)的次数等于预设次数或重叠度因数小于预设阈值；

所述重叠度因数利用式(5)计算：

式中，n表示元素总数，R_eq表示边界的等效半径，

表示所有元素的总合力。

作为本发明实施例的进一步改进，所有元素的半径均相等。

作为本发明实施例的进一步改进，所述元素的半径

φ表示填充率，V_b表示边界内空间的体积。

与现有技术相比，本发明的技术方案具有以下有益效果：本发明实施例提供的基于内作用力的颗粒建模方法，计算元素间以及元素与边界间的作用力，通过迭代的方式不断调整元素的位置，以尽可能减小簇内的总作用力，从而使簇内的元素拟合任意边界形状时尽可能消除内部的重叠量。

附图说明

图1是元素之间的位置关系图；

图2是元素与边界之间的位置关系图；

图3是本发明具体实例的方法流程图；

图4是二维中元素建模的迭代过程示意图；

图5是三维中元素建模的迭代过程示意图；

图6是二维中使用不同数量元素拟合边界图；

图7是二维中采用不同数量元素拟合的重叠度因数下降曲线图；

图8是三维中使用不同数量元素拟合边界图；

图9是三维中采用不同数量元素拟合的重叠度因数下降曲线图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行更清楚、完整地描述。

本发明实施例提供一种基于内作用力的颗粒簇建模方法，包括以下步骤：

步骤10)在边界内初始化生成n个位置随机的元素；n表示大于等于2的整数；

步骤20)计算元素之间的作用力以及元素和边界之间的作用力，调整元素的位置；

步骤30)判断是否满足预设条件，若不满足，则执行步骤20)；

步骤40)得到颗粒簇模型。

上述实施例的基于内作用力的颗粒建模方法，计算元素间以及元素与边界间的作用力，以尽可能减小簇内的总作用力为目标，通过迭代的方式不断调整元素的位置，从而使簇内的元素拟合任意边界形状时尽可能消除内部的重叠量。

步骤10)中，在边界内初始化生成n个元素，它们具有随机的位置。为了简化建模方法，本实施例中，n个元素的半径均相等。它们的半径为r₀，

φ表示填充率，V_b表示边界内空间的体积。

优选的，步骤20)中，计算一个元素所述的其它元素的作用力以及所受的边界的作用力，调整该元素的位置，具体包括：

步骤201)根据元素与其它元素之间的位置，得到元素所受的其它元素的作用力。

具体的，利用式(1)计算元素所受的其它元素的作用力：

式中，

表示元素i所受的元素j的作用力，r_i表示元素i的半径，r_j表示元素j的半径，O_i表示元素i的圆心或球心，O_j表示元素j的圆心或球心，d表示元素i和元素j之间的圆心距离或球心距离。

本实施例中，元素之间的位置关系分为三种情况，如图2所示，包含、相交和相离。本实施例考虑到元素位置关系的三种情况，利用式(1)计算元素所受的其它元素的作用力，通过它们之间的距离来确定元素之间作用力。当两个元素发生重叠(即包含或相交)时将会分别受到相互之间的排斥力，而当它们相离时该作用力为零。因此它们相互之间需要一个合适的作用力，元素之间应当尽可能靠近又避免存在重叠量。通过迭代优化元素之间的作用力，从而来减少元素之间的重叠。

步骤202)根据元素与边界之间的位置，得到元素所受的边界的作用力。

具体的，利用式(2)计算元素所受的边界的作用力：

式中，

表示元素i所受的来自边界点P_k的作用力，k_s为弹性系数，k_c为调节系数，O_i表示元素i的圆心或球心，r_i表示元素i的半径，O₀表示边界的形心，Q表示元素i的球心与边界点P_k的连线与元素i的远离O₀的交点。k_c一般取值为2，当元素数量较大时，增大k_c的值。

本实施例中，将具有不规则形状的边界表示为一个离散点集合，称这些离散的点为边界点。对于任意的不规则边界，将其形心记为O₀。元素与边界之间的位置分为四种情况，如图3所示，包含、相交且元素圆心或球心在边界内、相交且元素圆心或球心在边界外和相离。本实施例考虑到元素与边界位置关系的四种情况，利用式(2)计算元素所受的边界的作用力，通过它们之间的距离来确定元素之间作用力。当元素位于边界内时收到恒定的力将其向边界吸引，这个力不应当过大以避免所有的元素都集中在了边界上。当元素与边界相交或相离时，其受到的作用力应与P_k点和Q点之间的距离成正比。通过迭代优化元素与边界之间的作用力，从而增加元素与边界的拟合程度。

步骤203)结合元素所受的其它元素的作用力和元素所受的边界的作用力，得到元素的调整偏移量。

具体的，利用式(3)计算元素的调整偏移量：

式中，

表示元素i的调整偏移量，

表示元素i所受的其它所有元素的作用合力，

表示元素i所受的边界点P_k的作用力。

为了减少动态优化过程中元素球坐标的震荡情况并加快优化速度，本发明优选实施例利用式(4)计算元素的调整偏移量：

式中，

表示元素i的调整偏移量，

表示元素i所受的其它所有元素的作用合力，

表示元素i所受的边界点P_k的作用力,

表示元素i上一次的调整偏移量。

步骤204)根据所述元素的调整偏移量，调整元素的位置。

具体的，调整元素i的位置：

其中，

表示第t步调整前的坐标，

表示第t步调整后的位置坐标，l_r表示学习率。

本发明实施例，综合元素之间的作用力以及元素与边界间的作用力作为内作用力，根据内作用力来调整元素的位置，从而通过迭代优化元素的位置来减小内作用力，进而减少元素间的重叠量。

本发明实施例中，根据不同元素间的位置关系，提出了元素间的作用力以及式(1)以衡量元素间的重叠程度。通过进行多次迭代过程，使元素间的重叠量较快减少，同时内作用力相应的减小。根据元素与边界之间的位置关系，提出元素与边界间的作用力以及式(2)以衡量元素与边界间的重叠程度，通过进行多次迭代过程，可减少元素与边界的重叠量。能够比较快速且稳定地减少元素间的重叠量，同时使各元素能够尽可能的靠近边界且在边界内。

优选的，步骤30)中的预设条件为：

执行步骤20)的次数等于预设次数或重叠度因数小于预设阈值。

其中，利用式(5)计算重叠度因数：

式中，n表示元素总数，R_eq表示边界的等效半径，

表示所有元素的总合力。

本发明实施例中，利用式(5)计算的重叠量因数，反映当前情况下元素与元素间的重叠程度及与元素与边界轮廓的拟合程度。重叠量因数的数值由所有元素间的作用力总和以及元素与边界间的作用力总和决定的，因此当重叠量因数的数值越小，反应元素之间的重叠程度和元素与边界的拟合程度越好。将迭代次数和重叠量因数作为迭代终止条件，可以使得当迭代终止时，元素之间的重叠程度和元素与边界的拟合程度达到较好的稳定状态，即继续迭代也很难使得拟合程度大幅度优化。

优选的，元素为圆球或圆盘。本发明实施例方法既可用于二维建模，也可用于三维建模。

下面提供14个实例来验证本发明方法具有良好性能。实例均按照本发明方法建立颗粒簇模型，如图3所示，包括：

步骤1)根据预设定的轮廓形状创建出所对应目标边界的坐标位置。

步骤2)根据给定的填充率φ和边界内空间的体积V_b初始化生成n个具有半径r₀和坐标随机的元素。

步骤3)计算每个元素所受的来自相接触的所有元素和边界的作用力，计算当前迭代步的调整偏移量，根据计算得到的调整偏移量调整每个元素坐标，计算当前重叠度因数。

步骤4)如果重叠度因数OF不小于预设阈值且迭代次数没有达到最大迭代次数，则执行步骤3)，否则执行步骤5)。

步骤5)导出当前元素分布坐标。

实例1，n为10。

建模过程如图4所示。实例1是在二维中进行建模，在初始化后生成了10个大小一致、位置随机的圆盘。这些圆盘之间存在大量的重叠，因而为动态优化提供了很大的调整偏移量。在Step2中，原本大量重叠在一起的圆盘散开并超出了边界，在后面的迭代步中，每个圆盘都不断调整自己的位置，直到最后重叠度因数变化率低于设定阈值停止迭代。

从图4中可以发现，在元素的迭代过程中不可避免出现元素间存在重叠，甚至存在大部分都在边界外的个别元素，如Step2所示。其原因是当小区域内存在一定数量的高重叠量元素群体时，如Step1所示，该元素群体存在的高重叠量会在所在小范围产生相当大的作用力。从而导致有个别元素在某个时刻暂时被挤出边界，形成大部分面积都在边界外的元素。在这之后，由于这些大部分面积都在边界外的元素仍然受到来自边界的作用力，且此时该作用力被设定为远大于那些大部分面积或者全部面积都在边界内的元素，因而这些短暂出现在边界外的元素会迅速回到边界内。结果如Step12所示，元素之间重叠少且元素预边界的拟合程度高。

实例2，元素为圆球，n为10。

建模过程如图5所示。实例2是在三维中进行建模，在初始化后生成了10个大小一致、位置随机的圆球。结果如step12所示，元素之间重叠少且元素预边界的拟合程度高。结合图4和图5所示，元素在不同维度下表现出的迭代的规律相似。

实例3-8，元素为圆盘，n分别为6、10、25、50、75和100。得到的建模结果如图6所示，重叠度因数下降曲线如图7所示。其中，为了更加准确地获得实施例的重叠度因数，对每一种元素数目在固定迭代步数相同的情况下计算100次，并取它们在不同迭代步数时的重叠度因数平均值作为该时刻下的重叠度因数值。

从图6和图7中不难发现，当减小元素的粒径并且增加元素球数量对边界尽可能地进行拟合时，应采用更高的填充率，从而使边界内的空白区域更少。即在同样边界的情况下，采用越高的填充率(粒径越小，数量越大)能够使得元素可以在不增加重叠量的前提下更加充分的拟合边界。这表现在重叠度因数OF上就是，随着元素数量的增加，用元素拟合边界的最终重叠度因数也越来越小。这表明了，使用越多数量的元素便可以获得更好的建模效果，更好的缩小颗粒簇与实际被拟合颗粒间的误差。

观察图7中的重叠度因数下降曲线可以发现，重叠度因数在迭代最开始的阶段一般呈现出较大的波动，在某些阶段会出现“超调”的现象，即重叠度因数不减反增。然而，随着迭代的进行，重叠度因数曲线会逐渐平滑，并在下降到一定的值以后趋于稳定。这表明了本实施例方法的结果是收敛的。

实例9-14，元素为圆球，n分别为6、10、25、50、75和100。得到的建模结果如图8所示，重叠度因数下降曲线如图9所示。同样的，为了更加准确地获得实施例的重叠度因数值，我们对每一种元素数目在固定迭代步数相同的情况下计算100次，并取它们在不同迭代步数时的重叠度因数平均值作为该时刻下的重叠度因数值。

观察图8和图9中可以发现，与二维中类似，随着元素数量的增加，可以采用更高的填充率，所组成的颗粒簇的形状与目标边界越来越接近，同时迭代的最终重叠度因数也越来越小。这表明了，使用越多数量的元素便可以获得更好的建模效果，更好的缩小重叠度因数。同时，随着迭代的进行，重叠度因数曲线会逐渐平滑，并在下降到一定的值以后趋于稳定。这表明了本实施例方法的结果是收敛的。

对比图7和图9可以发现，与二维平面的重叠度因数曲线相比，在三维空间中的重叠度因数在迭代开始时呈现出的波动明显较小，出现“超调”的现象数量较少以及最终的稳定重叠度因数较大，这是由于在三维空间中单个元素直接接触的元素数量较多，所受来自不同方向的作用力的数量增加，从而使元素受到的限制作用较大，每次迭代发生的位移较小，重叠度因数曲线比较平缓。

以上显示和描述了本发明的基本原理、主要特征和优点。本领域的技术人员应该了解，本发明不受上述具体实施例的限制，上述具体实施例和说明书中的描述只是为了进一步说明本发明的原理，在不脱离本发明精神和范围的前提下，本发明还会有各种变化和改进，这些变化和改进都落入要求保护的本发明范围内。本发明要求保护的范围由权利要求书及其等效物界定。

Claims (10)

1.一种基于内作用力的颗粒簇建模方法，其特征在于，包括以下步骤：

步骤10)在边界内初始化生成n个位置随机的元素；n表示大于等于2的整数；

步骤20)依次计算每个元素所受的其它元素的作用力以及所受的边界的作用力，调整元素的位置；

步骤30)判断是否满足预设条件，若不满足，则执行步骤20)；

步骤40)得到颗粒簇模型。

2.根据权利要求1所述的基于内作用力的颗粒簇建模方法，其特征在于，所述元素为圆球或圆盘。

3.根据权利要求1所述的基于内作用力的颗粒簇建模方法，其特征在于，所述步骤20)具体包括：

步骤201)根据元素与其它元素之间的位置，得到元素所受的其它元素的作用力；

步骤202)根据元素与边界之间的位置，得到元素所受的边界的作用力；

步骤203)结合元素所受的其它元素的作用力和元素所受的边界的作用力，得到元素的调整偏移量；

步骤204)根据所述元素的调整偏移量，调整元素的位置。

4.根据权利要求3所述的基于内作用力的颗粒簇建模方法，其特征在于，所述步骤201)中，利用式(1)计算元素所受的其它元素的作用力：

式中，

表示元素i所受的元素j的作用力，r_i表示元素i的半径，r_j表示元素j的半径，O_i表示元素i的圆心或球心，O_j表示元素j的圆心或球心，d表示元素i和元素j之间的圆心距离或球心距离。

5.根据权利要求3所述的基于内作用力的颗粒簇建模方法，其特征在于，所述步骤202)中，利用式(2)计算元素所受的边界的作用力：

式中，

表示元素i所受的来自边界点P_k的作用力，k_s为弹性系数，k_c为调节系数，O_i表示元素i的圆心或球心，r_i表示元素i的半径，O₀表示边界的形心，Q表示元素i的球心与边界点P_k的连线与元素i的远离O₀的交点。

6.根据权利要求3所述的基于内作用力的颗粒簇建模方法，其特征在于，所述步骤203)中，利用式(3)计算元素的调整偏移量：

式中，

表示元素i的调整偏移量，

表示元素i所受的其它所有元素的作用合力，

表示元素i所受的边界点P_k的作用力。

7.根据权利要求3所述的基于内作用力的颗粒簇建模方法，其特征在于，所述步骤203)中，利用式(4)计算元素的调整偏移量：

式中，

表示元素i的调整偏移量，

表示元素i所受的其它所有元素的作用合力，

表示元素i所受的边界点P_k的作用力,

表示元素i上一次的调整偏移量。

8.根据权利要求1所述的基于内作用力的颗粒簇建模方法，其特征在于，所述步骤30)中，所述预设条件为：

执行步骤20)的次数等于预设次数或重叠度因数小于预设阈值；

所述重叠度因数利用式(5)计算：

式中，n表示元素总数，R_eq表示边界的等效半径，

表示所有元素的总合力。

9.根据权利要求1所述的基于内作用力的颗粒簇建模方法，其特征在于，所有元素的半径均相等。

10.根据权利要求9所述的基于内作用力的颗粒簇建模方法，其特征在于，所述元素的半径

φ表示填充率，V_b表示边界内空间的体积。

Images