基于二维最简颗粒排列的线性与赫兹模型参数研究

2019-07-09张科强陈能远马帅帅

陕西理工大学学报(自然科学版) 2019年3期

郭鸿，曹龙，张科强，陈能远，马帅帅，王普

(1.陕西理工大学土木工程与建筑学院，陕西汉中 723000；2.信息产业部电子综合勘察研究院，陕西西安 710054)

离散元方法(Discrete Element Method，DEM)在当今的科研领域应用愈来愈广泛，如岩土工程[1-2]、地质构造、机械工程等领域。离散元方法着眼于对固体和松散介质力学特性的深入研究，适合描述固体材料中细观、宏观裂纹扩展、破坏累积并断裂、破坏冲击和微震响应等问题。PFC2D/3D(二维/三维颗粒流程序)是离散元方法的主流软件之一，其中线性接触模型和赫兹接触模型是最常用的两个模型。相对于线性接触模型，赫兹接触模型能更真实地模拟球体与墙体、球体与球体的接触效果。刘畅等人[3]通过PFC2D建立了岩石材料的单轴和双轴离散元[4]数值模拟试验，进行了平行黏结有关参数[5-7]标定过程的研究，总结了细观参数与宏观参数之间的关系，方便平行黏结有关参数的确定。徐小敏等人[8]通过室内三轴试验的PFC3D模拟和结果的回归分析，基于线性接触模型建立了颗粒材料初始杨氏模量、初始泊松比等宏观弹性常数与颗粒法向刚度、颗粒刚度比等细观弹性常数间的经验公式。然而，就常用的两种接触模型(赫兹模型和线性模型)而言，其参数之间的相互联系和标定目前尚未见文献报道。

另外，在PFC2D运用过程中[9]，首先需要确定颗粒的相关参数，包括颗粒的尺寸，颗粒的刚度、强度，以及剪切模量、泊松比[10-11]等，计算模型参数的正确取值，有利于减少程序调试的次数。

鉴于此，本文以二维情况下规则排列的颗粒为研究对象，以赫兹模量和线性模型为理论基础，在目前国内外学者普遍运用的“试错法”[12]的基础上，深入研究线性接触模型和赫兹接触模型相关参数的关系，推导两种模型关键参数之间的关系。

1 数学模型推导

PFC2D是通过离散单元方法来模拟圆形颗粒介质的运动及其相互作用[13]。赫兹接触和线性接触可承受有限压力，两种模型用不同的计算方法模拟颗粒接触，可模拟材料的本构行为。

线性模型：线性接触模型由常用的剪切刚度Kn和Ks定义相关模型参数(墙体和球体、球体和球体)。

赫兹(Hertz)模型：Hertz-Mindlin 接触是基于Mindlin和德雷谢维奇(1953)理论的近似计算模型，在Cundall(1988)中描述了非线性接触公式。它仅适用于接触球面的情况。赫兹模型的基本参数为剪切模量G和泊松比v。

在PFC2D手册的基础上，对公式进行处理，得出线性模型与赫兹模型在形变量与受力相同的情况下，剪切模量G的计算方法。

线性模型中关于球体的法向刚度：

(1)

切向刚度：

(2)

赫兹模型中球体的法向刚度：

(3)

切向刚度：

(4)

在球体与球体、球体与墙体接触的部分参数计算如下：

球体与球体接触：

球体与墙体接触：

(5)

在赫兹模型中力与形变量的关系为曲线，对形变量积分得出形变量与力的关系式：

(6)

(7)

由式(7)可知：

(8)

由式(5)可知：

(9)

(10)

由式(10)得：

(11)

进一步，可得出：

(12)

由于赫兹模型的参数标定程序相对复杂，而线性模型参数标定容易，因此，在已知最终应力和变形等效的情况下，用线性模型代替赫兹模型进行数值模拟，可以提高计算效率。

2 PFC2D模型建立及数据分析

2.1 规则模型建立

图1 颗粒模型

2.2 线性模型和赫兹模型的数值压缩试验对比

4种规则排列颗粒体系的线性接触模型参数和赫兹模型参数如表1所示，相应的压缩试验结果如图2所示。不失一般性，对于2×1排列颗粒，泊松比取值为0.4，其余排列均为0.1，法向刚度均取为同一数值1.0×107N/m。

表1 颗粒模型参数表

(a) 2×1排列颗粒 (b) 10×1、10×3、10×5排列颗粒图2 线性接触模型和赫兹接触模型接触力与总变形量关系曲线

2.3 相同力与变形量条件下剪切模量与泊松比关系

以10×1颗粒体系为例，分别建立线性接触模型和赫兹接触模型与上节例子相同，线性模型的参数法向刚度取值为1.0×107N/m。以10×1颗粒体系为数值试样进行压缩试验，设定一定压缩量3.79 cm，则线性接触模型的最大接触力为21.1 kN。赫兹模型的参数泊松比取5个值，分别是0.1、0.2、0.3、0.4、0.5，则可以根据两个模型的最终等效确定剪切模量G。赫兹接触模型的数值压缩试验和用式(11)计算所得的剪切模量G和泊松比ν的关系如图3所示。可以看出，模拟值和理论值非常接近，从而验证了式(11)的正确性。