一种离散单元法的弹性可变形颗粒模型

第32卷第7期重庆大学学报

Vol.32No.72009年7月

Journal of Chongqing University

J ul.2009

文章编号:10002582X (2009)0720743204

一种离散单元法的弹性可变形颗粒模型

温　彤,雷　杰,

裴春雷

(重庆大学材料科学与工程学院,重庆400030)

摘　要:基于弹性变形的Hoo ke 定律,提出了一种考虑颗粒变形以及不同材料特性的离散单

元法(discrete element met hod ,D EM )的多边形颗粒模型,根据该模型开发了相应的DEM 程序。应用有限元方法和D EM 模拟了弹性颗粒的碰撞过程。通过与有限元计算结果比较,证明在处理颗粒的接触问题时,该弹性可变形颗粒模型比传统刚性模型能够准确地反映颗粒介质的实际变形和接触力的变化,从而能够提高DEM 分析的精度。

关键词:离散单元法;颗粒;变形;碰撞 中图分类号:TF124文献标志码:A

E lastic deform able particle model in discrete element method

WE N Tong ,LEI Jie ,PEI Chun 2lei

(College of Material Science and Engineering ,Chongqing U niversity ,Chongqing 400030,P.R.China )Abstract :A polygon particle model wit h discrete element met hod (DEM )is developed based on t he Hooke ’s law ,in which t he geomet ry change caused by t he elastic deformation and feat ures of material can be taken into account.A DEM p rogramme is developed based on t he p roposed model and collision processes of elastic particles is st udied wit h finite element met hod (FEM )and paring t he result of D EM wit h t hat of FEM.When dealing wit h t he problem of particles contact ,t he real deformation and t he contact force variation of t he particles can be presented more accurately wit h elastic deformable particle model ,compared wit h t hat f rom t raditional rigid particle model.

K ey w ords :discrete element met hod ;particle ;deformation ;collision

离散单元法(discrete element met hod ,DEM )是由Cundall 等人在20世纪70年代提出的一种分析离散体力学问题的数值方法[1]。该方法通过跟踪每一个颗粒的运动以及颗粒与周围环境的相互作用来认识整个颗粒系统,可以提供每个时间步中颗粒的位置、位移增量、速度以及角速度等重要信息。该方法有效弥补了连续介质力学在处理离散颗粒系统方面的局限,经过30多年的发展,成为了模拟非连续体的代表性方法,近年来在岩土工程、粉末冶金以及粉体工程等领域的研究中越来越得到重视[227]。

但现有的DEM 分析中,大多把颗粒假设为刚

性体,不能直接考虑实际颗粒受到外力作用时产生的弹性甚至塑性变形,同时通过颗粒间的几何叠加来处理和近似计算颗粒的接触、体积变化等,与实际情况有较大出入。笔者对传统的刚性模型进行了改进,提出了一种考虑颗粒弹性变形引起几何形状改变的颗粒模型,并开发了相应的DEM 程序。

1　离散单元法简介

常用的DEM 颗粒模型有圆形颗粒、椭圆形颗

粒以及多边形颗粒等类型。圆形颗粒建模简单、计算容易,使用较普遍,但与颗粒实际形状常常有较大差别。椭圆形颗粒和多边形颗粒能考虑实际形状的不规则性,从统计意义上讲,其计算结果的误差更小一些,但计算量更大。

笔者采用图1所示的多边形颗粒定义方式[224,8210],每个颗粒的几何形状由参数n 、r 和θ控制。其中:n 为多边形颗粒的边数;r 为各顶点的半径;θ为各边与颗粒中心所夹的角度。

控制颗粒运动的方程为牛顿第二定律

m i d 2x i d t 2=m i g +f c ,I i d 2θi

d t

2

=M i ,(1)其中:m i 、I i 分别为颗粒i 的质量和转动惯量;x i 、θi

分别为颗粒的质心位置和角度方向;g 为重力加速度;t 为时间;f c 、M i 分别为接触力和作用在颗粒i 上的力矩

。

图1　DEM 的颗粒构造模型

颗粒间的接触处理是DEM 模拟的重要一环。

首先需要从几何上判断哪些颗粒之间有接触,然后再计算颗粒间的接触力,具体的判断和计算过程与颗粒的定义方式有直接的联系。笔者采用“两步探

测法”[8210]

,即:首先分别以2个颗粒的最大半径画圆(见图2),如果2个颗粒中心之间的连线小于2个最大半径之和,则认为颗粒“相邻”;记录所有相邻点的信息后,在“相邻”颗粒之间判断颗粒顶点间的“刺入”情况,如果有,就认为存在接触,然后就可以计算它们之间的接触力

。

图2　颗粒的接触判断

目前有2种DEM 的接触处理模型[11213]:一种

称为硬球模型,另一种称为软球模型。前者采用的基本方程是积分形式的牛顿方程,它假设2个颗粒的碰撞不需要时间,一接触马上弹开,因此不能处理几个颗粒同时碰撞的情况。后者采用的是微分形式的牛顿方程,可以处理多颗粒的同时碰撞,求出碰撞过程中的力速度位移的演化过程,但计算时间更长。文中采用的是软球模型。

2　可变形的DEM 颗粒模型

颗粒受到外力作用后,会产生弹性及塑性变形。当应力不大时,可以仅考虑弹性变形。通常D EM 模型中假定颗粒形状不变,利用颗粒间的相互叠加来计算接触力,叠加量越大,则相互作用力越大。但实际上2个颗粒接触并相互挤压时,都会有一定的变形,因此只有考虑颗粒的变形才能得到合理的结果和更高的分析精度。2.1　颗粒间的接触力

将颗粒间相互作用的接触力分为法向力f n 和切向力f t 两部分[4],如图3所示。在假定颗粒材料为单一的条件下,法向力用式(2)计算。

f n =K n α-χv n ,

(2)其中:K n 为法向接触刚度,由颗粒材料特性决定;v n 为法向相对速度;χ为一个粘性参量,用来计算粘性力(不考虑粘性力时可以不计);α为2个颗粒之间的叠加量,它决定了接触力的大小

。

图3　颗粒的接触力模型

实际中,相互作用的颗粒可能由2种以上的不同材料组成,力学特性的差别主要体现在密度、杨氏模量等参数的不同。当不同材料的2个颗粒接触时,笔者根据分别赋予的不同密度和接触刚度系数将法向力f n 定义为

f n =(K 1α1+K 2α2)-χ

v n ,(3)

其中:K 1和K 2分别为2种材料的法向接触刚度;α1

和α2分别为2个颗粒相互刺入的深度。

切向力有2种机制:接触处有相对运动的“滑

4

47重庆大学学报 第32卷