颗粒材料微观统计力学研究综述

姚绪坤，秦建敏

（大连大学 土木工程技术研究与开发中心，辽宁 大连 116622）

颗粒材料微观统计力学研究综述

姚绪坤，秦建敏*

（大连大学 土木工程技术研究与开发中心，辽宁 大连 116622）

基于颗粒材料的多尺度结构，探讨和搭建微观结构与宏观力学行为之间的联系，研究新的关系式来解决跨尺度耦合的过程是颗粒材料研究的挑战之一。由于微观结构层次的非平衡性和强相互作用，对微观物理量和几何参数的统计描述可以提供连接微细观运动与宏观力学行为之间的联系。本文介绍了颗粒材料不同结构层次的统计参数，总结了研究者关于统计参数与力学概念之间关系的成果，得出对于颗粒材料，正确选取微观参数并予以采用科学的统计方法，是准确预测宏观力学行为的关键。

颗粒材料；多尺度；统计参数

0 引言

自然界中岩土、积雪、浮冰，工业生产活动中的药品、煤炭、矿石、化工产品，日常生活中的糖盐、谷物等都是颗粒材料，颗粒材料在人们生活中普遍存在并与社会生产活动密不可分[1]1。颗粒材料是以接触力为主要作用的多体系统。从材料结构角度看，颗粒材料比常见的水、空气等普通液体，或晶体、金属等普通固体多出一个以邻近颗粒大小为特征长度的细观层次，是典型的多尺度结构并且不同尺度上有不同的结构层次及不同的演化物理和速率，因此颗粒材料表现出比普通固体和普通液体等一般材料更为丰富的力学行为[2]5。

颗粒材料与连续介质材料相比，具有许多特殊的力学特性：（1）颗粒材料的变形是由颗粒移动所产生的结构变形和颗粒的自身变形所共同引起的，结构变形是不可恢复的，且在加载和卸载时都会产生；（2）颗粒材料在常规静水压力作用下会产生不可恢复的体积变形；（3）颗粒材料具有压硬性和剪胀性，静水压力与剪应变以及剪应力与体应变之间具有耦合关系；（4）颗粒材料在不同初始密度情况下可表现出不同的硬化或软化特性，并且其软化特性是一种结构软化；（5）颗粒材料没有明显的弹塑性阶段[3]。

对于颗粒材料的本构关系的研究，在我国有近20年，国外已经有30多年历史。传统方法中，利用弹塑性基本理论，通过修正土力学中的一些本构模型，使之符合颗粒材料的应力应变特征。但颗粒材料的变形，除包括瞬时变形外，还存在着明显的体积流变和剪切流变。颗粒材料本构关系的数学模型很多，但无法真正为工程界所接受，解决工程实际问题。

1 颗粒材料的多尺度结构

1.1 多尺度框架

基于颗粒材料的多尺度结构特征，研究者提出了不同的尺度分割方法。要根据问题的重点将颗粒材料分为不同的尺度。孙其诚和金峰[4]等将颗粒材料分为颗粒（微观）－力链（细观）－颗粒体系（宏观）这三个尺度，对以力链为核心的多尺度结构特征和其演变规律进行研究，在颗粒固体应力传播模式、颗粒流体流动本构关系等方面做出很大贡献。秦建敏和张洪武[5-7]等对颗粒体系的微、细和宏观的多尺度概念进行改进，提出以孔隙胞元作为中间尺度的“颗粒接触对－孔隙胞元－颗粒体系的微－细－宏观多尺度体系，用基于孔隙胞元的离散元方法对颗粒体进行双轴加载数值试验，模拟了以滑动变形表征的变形局部化现象，研究了颗粒材料剪胀性的微观机制。另外一些学者[8-12]对颗粒材料的接触角、接触力、接触时间数、组构张量、配位数等物理和几何参数进行统计，并采用粗粒化方法和平均场理论得到颗粒材料的弹性模量、应力张量、应变张量、能量密度、热流等，从而构建颗粒材料的本构关系，描述其宏观力学性质。

1.2 各尺度变量之间的联系

国外研究者 Love[13]，Goddard[14]，Rothenburg 和Selvadurai[15]，Rothenburg 和 Kruyt[16]，Bagi[17]和 Li等[18]用材料应力均匀的假设，基于颗粒的受力方程或者虚功原理对无体积力和有体积力的的静力平衡状态进行了推导。Li等[19]考虑了动力问题中颗粒旋转问题，De Saxce[20]对动力问题中的时间效应进行了考虑。秦建敏等[2]8-15引用Li等中的相关内容，基于牛顿第二运动定律和均匀应力场的假设，考虑了颗粒旋转的影响，推导了在静力和动力情况下都适合的柯西应力细观表达式：

图1 颗粒-颗粒、颗粒-边界的相互作用示意图

对应变张量的细观定义多种多样。Bagi[21]认为应变张量的细观定义可以分为两类，第一类是将颗粒材料看作为等效连续体，假设颗粒中心的位移集中发生在接触点处，基于均匀场的假设，建立相关几何量和颗粒位移与位移梯度张量的关系，这是基于等效几何图形的应变定义。第二类是基于平均位移梯度(应变)拟合出颗粒材料体系的位移场，找出位移场与颗粒材料体系的特征位移偏差最小的平均位移梯度，这是基于最优拟合的应变定义，特征位移可以选取接触点处的相对平动位移或者颗粒中心的平动位移[22]。

基于等效几何图形和颗粒平动位移的应变张量的细观定义，Rothenburg和Kruyt提出适用于二维颗粒材料基于 Satake的图形系统的定义，Bagi[23]给出适用于由单纯形为单元的几何系统的定义，Li等提出适合于多面体/多变形单元系统的定义。Li等二维情况下的表达式与Kruyt和Rothenburg的表达式都不要求对应的几何系统由单纯形组成，这里给出Li等的应变张量细观定义：

图2 面积矢量χ的定义

如果相对位移有清晰定义，那么上述应变张量细观表达式的推导适用于对任何由多面体/多边形构成的几何系统，包括孔隙胞元系统，但并不只限于孔隙胞元系统。

应力张量的细观表达式基于牛顿第二定律，首先，对颗粒材料进行体积平均得到平均应力，用高斯公式转换成边界上的力，接着表示成离散元的边界力，在每个接触点上，作用力满足平衡条件，大小相等方向相反，最后用牛顿第二定律将局部几何变量（接触矢量）耦合到表达式。同样，应变张量也是对颗粒材料的体积平均，接着基于高斯公式转换成边界上的量，表示成离散元的相对位移，在每个接触点上，相对位移满足相容性条件，然后用位移相容性条件将局部几何变量（互补面矢量）耦合到表达式。比较应力、应变张量的细观结构表达式的推导过程，两者最后都是用与微元联系的细观量来表示以边界量表示的表达式来得到细观定义，推导过程很相似。在这里推导宏细观联系的方法并不局限于相关量的物理背景。

应力和应变张量的表达式表明应力与单位体积的接触力和接触矢量有关，应变与单位体积接触点的相对位移和局部几何矢量有关，然而研究表明单位体积的接触以力及相对位移不同，颗粒材料还是可以表现出完全相同的宏观力学性质，这说明颗粒材料宏观表现出来的应力与应变只是与微观量（接触力，相对位移等）的统计量有关。何旭珍[24]在接触方向、接触力、接触向量、孔隙向量满足统计规律的条件下，得到了积分形式的应力关系：

2 统计力学

由于颗粒材料的宏观行为不是由某个特定接触点的性质决定，而是由大量的微观量的统计规律决定，统计方法可以摒弃大量无用信息，常常被用来连接不同的尺度，因此统计方法可以是一种研究大量颗粒组成的宏观物体行为行之有效的方法。白以龙和夏蒙棼[25-27]详细地介绍了固体统计细观损伤力学和损伤演化诱致突变理论，尝试将力学及统计物理学和非线性科学结合起来研究固体材料损伤和破坏问题。

2.1 统计参数

在研究颗粒材料动力特性时，需要研究基本物理量和几何参数如局部体积分数、配位数、接触力、接触角、力链在空间的真实分布与变化，但它们都会在较小的范围内发生很大的变化。为了较好的描述颗粒材料系统的宏观状态变化特征，必须对这些基本物理量和几何参数进行有效而准确的统计与平均。

2.1.1 局部体积分数

体积分数指的是颗粒实体面积（或体积）占整个统计域的面积（或体积）的百分比，体积分数的大小表示的是颗粒材料内部颗粒分布的疏密程度。在二维情况下基本的计算公式如下：

2.1.2 配位数

配位数指的是颗粒集合内单个颗粒体与周围颗粒的平均接触数，其计算公式如下：

公式中Nb是颗粒体中心处于统计域内的颗粒个数；是第b 个颗粒周围接触的颗粒数在对颗粒材料加载的过程中，局部配位数会发生明显的变化，因此需要计算颗粒材料的不同区域内相应颗粒的配位数，才能描述其在空间的变化特征。局部配位数的统计方法与局部体积分数类似。

2.1.3 接触力

当外部荷载或约束边界作用在颗粒材料时，颗粒材料内部颗粒之间会形成分布不均匀的接触力，一小部分颗粒参与了接触力的传递并承担了大部分的外部荷载，而其它颗粒只会承受微弱的接触力或直接不受力。假设在计算区域内满足条件的颗粒接触点个数为由于每个接触点上两个颗粒的接触力的方向刚好是相反的，故统计每个颗粒的接触力每个接触点上的接触力实际计算了两次。如果是颗粒和边界之间发生接触，接触力只计算一次。在统计满足条件的接触力之前用下面的公式对均值进行归一化：

Fij表示第个j颗粒对i个颗粒产生的接触力大小。归一化以后，就可以得到接触力F的分布区间为如果将颗粒间的每次相互接触都考虑成等概率事件，那么每一个接触力的分布概率都为1 2Nc。将接触力的取值区间等分成n个小区间，在任意一个小的接触力区间内，接触力F 分布的概率密度函数为：

2.1.4 接触方位角

乡村产业兴旺，意味着乡村经济发展充满活力。尽管反映某个乡村的产业兴旺可能是一个或两三个产业，但总体上衡量乡村产业兴旺的角度和标准应是多维的，因而乡村产业兴旺的内涵特征展现得比较丰富，至少有以下方面。

2.1.5 力链

在重力或外荷载作用下，密集排布的颗粒材料中，颗粒自由活动空间小，颗粒间相互挤压变形。变形较大且连接成准直线形的颗粒，承担较大份额的重力或外荷载，形成强力链；其它颗粒间接触变形较小，传递的外力较小，形成弱力链[1]73。孙其诚[28-29]等通过离散动力学模拟，对静态堆积颗粒不同接触力的力链长度和颗粒间接触力与重力夹角进行统计分析得到计算不同接触力下的力链长度和夹角的公式：

并对二维颗粒体系进行单轴压缩，计算发现了强力链长度的幂率分布规律，设力链长度为其中是体系的平均粒径，表示该力链上颗粒数目和颗粒粒径。对于发现强力链长度的分布按照幂率分布

2.2 统计参数与力学概念之间的关系

王等明[11]57-60通过对剪切颗粒流中的局部体积分数和配位数的变化进行统计，研究分析了剪切颗粒流的剪胀性，得到剪切速率会影响系统的剪胀特征。随着剪切速率的增大，局部体积分数变小，平均配位数降低，剪胀效果更加明显。

朱纪跃[30]统计分析了微观尺度下颗粒间的力学行为，得到不同应变阶段接触力大小的分布特征，发现在摩擦系数和围压一定的前提下，开始平均接触力随轴应变近似于线性增大，当轴向应变大于某个值，平均接触力不再随着轴向应变增大，而是随着轴向应变先减小，接着在一定范围内波动。并通过不同加载阶段接触力方向的分布特征，发现在颗粒样本的加载过程中，偏应力峰值与接触力方向的分布特征有着某种内在的关联。

刘斯宏[12]通过双轴压缩试验得到颗粒接点数按接触角的统计分布数等向压缩时，锯齿状分布可用虚线圆来拟合，颗粒结构为各向同性；剪切时锯齿状分布变成花生状或椭圆形，颗粒结构的各向异性增大。王子健[31]为更好的考虑颗粒骨架的局部力学性能在 N（α）的基础上引入 F（α），表示沿着颗粒接触角 所统计出来的颗粒接触挤压力平均值，并分析了分布随着加载时间变化的规律，从细观机理上解释了宏观试验中破坏滑动面一般呈现左右的事实。

陈立平[9]等在细观颗粒层面研究砂土破坏和摩擦特性的研究中引入了数理统计方法，将组成砂土的单个颗粒简化为具有一定长细比的空间椭球体，假设方向角满足正态分布规律，应用细观统计的方法，从颗粒层面对砂土的摩擦特性和破坏机制进行研究，建立了内摩擦角与外荷载之间的对应关系。

张惠新[33]通过对颗粒间强力链长度分布规律的统计分析，得到了颗粒间强力链长度分布的概率密度函数；从Matlab绘制颗粒间接触力箭头分布图中得到了刚性压条压入颗粒体系时强弱力链的出现、发展、演变以及力链网络的形成过程，对力链的形成和发展进行了定量分析，对获得颗粒材料宏观参数的微观表达，搭接从微观到宏观的多尺度研究起到推动作用。

3 结语

颗粒的空间排列及相互作用力决定了颗粒材料的宏观变形和强度特性，建立微细观结构与宏观特性的联系，是颗粒材料力学的研究思路之一，也是不同于传统连续介质力学的本质所在。统计方法常常用来连接两个不同尺度，基于微观结构、微观运动量、微观接触力的统计来描述颗粒材料的宏观行为，但是低尺度的物理本质如何传递到更高尺度，与更高尺度强耦合，以及如何用公式表示对更高尺度的影响，它应该建立在实验和对力学规律的结果的深入的分析的基础上。到目前为止，颗粒材料的微观统计力学的研究还比较欠缺。

[1]孙其诚, 王光谦. 颗粒物质物质力学导论[M]. 北京∶ 科学出版社, 2009.

[2]孙其诚, 厚美瑛. 颗粒物质物理与力学[M]. 北京∶ 科学出版社, 2010.

[3]秦建敏. 基于离散元模拟的岩土力学性能研究及应变局部化理论分析[D]. 大连∶ 大连理工大学, 2007∶ 3-5.

[4]孙其诚, 金峰. 颗粒物质的多尺度结构及其研究框架[J].物理, 2009(4)∶ 225-232.

[5]张洪武, 秦建敏. 颗粒材料内部滑动变形局部化的离散元模拟[J]. 计算力学学报, 2011, 28∶ 6-9.

[6]秦建敏, 迟璐璐. 颗粒材料剪胀性的微观力学分析[J]. 岩土力学, 2013, 34(5)∶ 1508-1514.

[7]秦建敏, 迟璐璐. 颗粒材料多尺度力学行为的离散元模拟[C]//颗粒材料计算力学研究进展∶ 颗粒材料计算力学会议论文, 2012.

[8]郑晓静. 散体物质多尺度力学行为的若干研究[C]//颗粒材料计算力学研究进展∶ 颗粒材料计算力学会议论文,2012.

[9]陈立平. 基于细观统计的砂土摩擦特性与破坏机制研究[J]. 岩土力学, 2013, 34(12)∶ 3471-3478.

[10]叶晓燕. 冲击作用下颗粒介质的动力学特性及尺度关系的理论研究[D]. 兰州∶ 兰州大学, 2013.

[11]王等明. 密集颗粒系统的离散单元模型及其宏观力学行为特征的理论研究[D]. 兰州∶ 兰州大学, 2009.

[12]刘斯宏. 基于细观结构的颗粒介质应力应应变关系研究[J].科学通报, 2009, 54(11)∶ 1496-1503.

[13]Love A E H. A Treatise of Mathematical Theory of Elasticity[M]. Cambridge University Press, 1927.

[14]Goddard J. An elastohydrodynamics theory for the rheology of concentrated suspensions of deformable particles [J].Journal of Non-Newton i an Fluid Mechanics, 1977(2)∶169-189.

[15]Rothenburg L, Selvadurai A P S. A micromechanical definition of the Cauchy stress tensorfor particular media [J].Mechanics of Structured Media. Elsevier, Amsterdam, 1981,469-486.

[16]Kruyt N P, Rothenburg L. Micromechanical definition of the strain tensor for granular materials [J]. Journal of Applied Mechanics, 1996, 118∶ 706-711.

[17]Bag K. Microstructural stress tensor of granular assemblies with volume forces [J]. Journal of Appl Mech. 1999, 66(4)∶934-936.

[18]Li X H, Yu H S, Li X S. Macromicro relations in granular mechanics [J]. International Journal of Solids and Structures,2009, 46∶4331-4341.

[19]Li X, Li X S. Micromacro quantification of the internal structure of granular materials [J]. Journal of EngMech, 2009,135(7)∶641-656.

[20]De Saxce G, Fortin J, Millet O. About the numerical simulation of the dynamics of granular media and the definition of the mean stress tensor [J]. Mech Mater, 2004, 36∶1175-1184.

[21]Bagi K. Analysis of microstructural strain tensor for granular assemblies [J]. InternationalJournal of Solids and Structures,2006, 43(10)∶ 3166-3184.

[22]孙其诚, 程晓辉, 季顺迎, 等. 岩土类颗粒物质宏—细观力学研究进展[J]. 力学进展， 2011,41(3)∶ 351-371.

[23]Bagi K. Stress and strain in granular assemblies [J].Mechanics of Materials, 1996, 22∶ 165-177.

[24]何旭珍. 基于微观力学的颗粒材料宏观微观力学性质研究[D]. 北京∶ 北京交通大学, 2012.

[25]白以龙, 汪海英, 夏蒙棼, 等. 固体的统计细观力学-连接多个耦合的时空尺度[J]. 力学进展, 2006, 36(2)∶286-305.

[26]白以龙, 柯孚久, 夏蒙棼. 固体中微裂纹系统统计演化的基本描述[J]. 力学学报, 1991, 23(3)∶ 291-298.

[27]夏蒙棼, 韩闻生, 柯孚久, 等. 统计细观损伤力学和损伤演化诱致突变[J]. 力学进展, 1995, 25(2)∶ 145-159.

[28]孙其诚, 王光谦. 静态堆积颗粒中的力链分布[J]. 物理学报, 2008, 57(8)∶ 4668-4674.

[29]孙其诚, 金峰, 王光谦, 等. 二维颗粒体系单轴压缩形成的力链结构[J]. 物理学报, 2010, 59(1)∶ 31-37.

[30]朱纪跃. 基于离散单元法的颗粒物质静动力学行为研究[D]. 兰州∶ 兰州大学, 2013.

[31]王子健, 刘斯宏, 王怡舒, 等. 不同应力路径下颗粒材料细观结构分析[J]. 水力学报, 2013,44(7)∶ 772-778.

[32]张慧新. 钢珠颗粒体系受集中力时的力链研究[D]. 青岛∶青岛理工大学, 2014.

The Research Review of Granular Material Microscopic Statistical Mechanics

YAO Xu-kun，QIN Jian-min*
(Research and Development Center of Civil Engineering Technology, Dalian University, Dalian 116622, China)

∶Based on Multi-scale structure of granular materials, it is one of the challenges of granular materials research to set up the connection between the microscopic structure and macroscopic mechanical behavior and explore new formula to solve the coupled across scales process. Due to nonequilibrium and strong interaction of the microstructure, a statistical description of microcosmic physical quantities and geometrical parameters can provide the connection between the micro motion and macroscopic mechanical behavior. This paper introduces the statistical parameters of different structure levels of the granular material and summarizes achievements of the researchers about the relationship between the statistical parameters and the mechanical concept. It is concluded that for granular materials, correct selection of microscopic parameters and using scientific statistical methods are key to predict the macroscopic behavior.

∶granular materials; multi-scale; statistical parameters

O347.7

A

1008-2395(2015)06-0011-06

2015-05-04

国家自然科学基金资助项目（11272069）。

姚绪坤（1989-），男，硕士研究生，研究方向：防震减灾与土木环境工程。

秦建敏（1979-），女，博士，副教授，研究方向：岩土介质力学等。