力学学报  2018 , 50 (2): 315-328 https://doi.org/10.6052/0459-1879-17-420

Orginal Article

基于微形态模型的颗粒材料中波的频散现象研究

修晨曦, 楚锡华

武汉大学土木建筑工程学院工程力学系,武汉 430072

STUDY ON DISPERSION BEHAVIOR AND BAND GAP IN GRANULAR MATERIALS BASED ON A MICROMORPHIC MODEL

Xiu Chenxi, Chu Xihua

Department of Engineering Mechanics,School of Civil Engineering, Wuhan University,Wuhan 430072,China

中图分类号:  O34

文献标识码:  A

收稿日期: 2017-12-20

接受日期:  2017-12-20

网络出版日期:  2018-04-17

版权声明:  2018 《力学学报》编辑部 《力学学报》编辑部 所有

基金资助:  国家自然科学基金资助项目(11772237,11472196).

作者简介:

作者简介:楚锡华,教授,主要研究方向:计算固体力学,颗粒材料力学. E-mail: Chuxh@whu.edu.cn

展开

摘要

基于颗粒材料冲击与波动响应特性的调控波传播行为的超材料设计受到广泛关注,设计这类材料需要对颗粒材料的波传播机制及调控机理有深入认识. 波在颗粒材料中传播的频散现象及频率带隙等行为与材料的非均匀性密切相关,通常讨论频散现象是基于弹性理论框架建立微结构连续体或高阶梯度连续体等广义连续体模型来进行. 本研究基于细观力学给出了一个颗粒材料的微形态连续体模型. 在该模型中,考虑了颗粒的平动和转动,且颗粒间的相对运动分解为两部分:即宏观平均运动和细观真实运动. 基于此分解,提出了一个完备的变形模式,得到了对应于不同应变及颗粒间运动的宏细观本构关系. 结合宏观变形能的细观变形能求和表达式,获得了基于细观量表示的宏观本构模量. 应用所建议模型考察了波在弹性颗粒介质的传播行为,给出了不同形式的波的频散曲线,结果显示此模型具有预测频率带隙的能力.

关键词: 颗粒材料 ; 微形态模型 ; 波传播 ; 频散 ; 频率带隙

Abstract

The design of metamaterials is paid more attention to control the behaviors of the wave propagation based on response characteristics of shock and wave in granular materials, and it requires in-depth understanding of the propagation mechanism and control mechanism of waves for granular materials. The dispersion behavior and frequency band gap of granular materials are closely related to the heterogeneity. Generally, the dispersion behavior and frequency band gap are based on the elastic theory framework to establish a generalized continuum model including the microstructural continuum or the high order gradient continuum. This study proposes a micromorphic continuum model based on micromechanics for granular materials. In this model, the translation and the rotation of particles are taken into consideration, and the relative motion between particles is decomposed into two parts: the macroscopic mean motion and the microscopic actual motion. Based on this decomposition, a complete pattern of deformation is obtained. The macroscopic deformation energy is defined by a summation of the microscopic deformation energy at each contact. As a result, the micromorphic constitutive relation can be derived, and the corresponding constitutive modulus can be derived by microscopic parameters in terms of contact stiffness parameters and microscopic geometric parameters. The proposed model investigates the propagation of waves in an elastic granular medium, give dispersion curves for different waves such as longitudinal, transverse and rotational waves and predict the frequency band gap. It proves that the proposed model has the ability to describe dispersion behaviors and predict the frequency band gap in granular materials.

Keywords: granular materials ; micromorphic model ; wave propagation ; dispersion ; frequency band gap

0

PDF (2987KB) 元数据 多维度评价 相关文章 收藏文章

本文引用格式 导出 EndNote Ris Bibtex

修晨曦, 楚锡华. 基于微形态模型的颗粒材料中波的频散现象研究[J]. 力学学报, 2018, 50(2): 315-328 https://doi.org/10.6052/0459-1879-17-420

Xiu Chenxi, Chu Xihua. STUDY ON DISPERSION BEHAVIOR AND BAND GAP IN GRANULAR MATERIALS BASED ON A MICROMORPHIC MODEL[J]. Chinese Journal of Theoretical and Applied Mechanics, 2018, 50(2): 315-328 https://doi.org/10.6052/0459-1879-17-420

引 言

颗粒材料是指由大量离散固体颗粒及颗粒间孔隙构成的集合体,其具有强烈的非均匀性与复杂的力学性质,诸如多尺度行为、可 破碎性、各向异性 等 [1,2,3]. 深入研究颗粒材料的基本力学行为不仅有助于减少其生产运输的过程中不必要的浪费,也有助于提高预测和控制诸如泥石流、山 体滑坡等由于颗粒材料变形、破坏而导致的自然灾害[4]. 考虑到颗粒材料的离散特性,离散单元法是模拟其力学行为的更加自然与合理的一个方法. 但是在实际应用中,颗粒数量往往十分庞大,其运动规律也很复杂,基于离散单元法模拟和预测工程问题仍然是一个巨大的挑战. 因而,连续体方法或者多尺度方法模拟颗粒材料仍是比较合适与便捷的途径,建立连续体框架下的模型研究其力学行为具有十分重 要的意义[5,6].

颗粒材料的宏观性质与细观尺度上的颗粒排布、接触刚度等信息密切相关 [3]. Chang等 [7]认为通过微结构连续体途径(microstructural continuum approach)可以考虑细观结构的影响,基于能够反映细观离散特性的宏观变形量,得到宏观尺度上的本构关系. 这种方法已被 大量用于研究颗粒材料的力学行为中[7,8,9,10,11,12,13,14,15]. 其中,Cosserat理论[8-13,16-18]十分具有代表性且得到了广泛应用. 该理论通过引入旋转自由度,推导出与其对应的偶应力及表征细观结构信息的特征长度,从而适合模拟颗粒材料的微结构, 且更加符合颗粒材料的物理实际. 但是,经典Cosserat理论只能部分表征微结构信息,不能与颗粒材料的离散结构特性建立显式联系[14]. Nemat-Nasser等 [14,15]为表征微结构特性,引入了结构张量,并建立起其与宏观应力的联系,使本构关系能 够反映微结构的演化,从而发展出基于细观力学的连续介质描述. 该理论通过应力、应变概念描述颗粒材料的宏观变形行为,基于相互接触的颗粒间的局部运动及力学行为建立细观本构关系, 其关键是如何建立颗粒材料离散的细观力学响应与宏观上连续行为之间的联 系[19]. 考虑颗粒材料的细观离散特性,基于经典连续体的应力应变定义不再适用,大量研究[19,20,21,22,23,24,25]致力于发展颗粒材料细观结构的描述、细观结构应力应变的定义以及采用均一化方法建立宏细观的联系. 如Chang等[20]将颗粒材料的本构关系划分为3个层次:表征元、微单元及颗粒接触,以表征元体现宏观应力应变关系,以微单元体现局部应力应变关系,以颗粒接触体现力与位移的关系. 并基于Vogit均一化理论应用最优拟合假设、运动假设、静态假设等建立颗粒材料宏细观联系及本构关系[26,27].

微形态理论[28,29,30,31,32,33]同样具有描述复杂微结构变形的能力,可以被用于模拟颗粒材料. Misra [28]根据Mindlin[29]和Eringen[30]的微结构线弹性理论建议了基于细观力学方法的颗粒材料微形态模型. 微形态理论假设物体是由宏观物质与细观物质组成的变形协调体,其每个物质点都是一个微结构胞元,并且引入9个变形自由度用以描述微结构的变形与旋转,适用于复杂微结构的变形模拟. 针对颗粒材料Misra等 [28]认为一个物质点即为一个颗粒集合体体积单元,将颗粒运动分解为宏观运动与细观运动两个部分,并基于细观变形能求和得到宏观变形能函数,再通过对宏观变形能求导构建基于细观结构信息表示的宏观本构关系.

近年来,基于颗粒材料冲击与波动响应特性的调控波传播行为的超材料设计引起了研究者的兴趣,为此研究者已开展了大量基于广义连续 体模型(高阶梯度模型或Cosserat连续体模型)的波动分析 [34-36]. 正如上文所述,Cosserat连续体模型因包含了可与细观上颗粒转动相关联的转动自由度以及能够反映一定细观结构的内禀特征长度,且相对 简单而受到颗粒材料研究领域的广泛关注[17,30,34-35,37]. 但是Cosserat连续体仅能反映细观结构的刚体转 动 [35,38],考虑到颗粒材料细观结构或者中尺度结构上诸如拉压与剪切等其他变形 模式的重要性[39,40], 因而,能够反映细观结构完整运动模式的微形态连续体模型 [29,38]引起了研究者的重视[28,41-43]. 一阶微形态连续体模型除了描述宏观运动的自由度外,还包含描述微/细观结构变形的自由度,因此在其框架下基于细观途径发展颗粒 材料连续体模型能在一阶梯度情况下为更多细观信息的均匀化提供对应的描述变量.

本文基于Mindlin-Eringen的微结构理论[29,30]和Misra的微形态理论 [28],给出了一个颗粒材料的微形态连续体模型. 不仅对颗粒的平动进行分解,同时对颗粒转动也进行了分解:颗粒运动分解为宏观平均运动(平动与转动)和细观真实运动(平动与转动)两 部分,并认为颗粒的细观真实运动(平动与转动)是由宏观平均运动(平动与转动)与相对运动(平动与转动)组成. 基于此分解,求得相应的细观本构关系和细观变形能函数. 结合宏观变形能的细观变形能求和表达式与细观本构关系的表达式,可获得基于细观量表示的宏观本构关系与相应的本构模量,其中包括了偶应力(相对偶应力)与微曲率(相对微曲率)的关系及其本构模量. 通过推导的微形态模型探讨波在弹性颗粒介质中传播的现象,考察纵波、横波、横向剪切波、横向旋转波的频散关系 [44],并对其频率带隙[44]进行预测.

1 颗粒的运动描述

根据Mindlin的微结构理论 [29],建立一个全局坐标系 x,并考虑一个体积单元为 V,此时, V代表了一个物质点 P,如图1所示. 然后在 P的重心建立一个局部坐标系 x',保持 x'的坐标轴平行于 x的坐标轴,见图1. 考虑 V中的两个颗粒分别为 p和它的临近颗粒 q,那么颗粒 p的位移可由颗粒 q的位移用泰勒级数展开[28]

ϕip=ϕiq+ϕi,jqlj+(1)

式中, ϕi为颗粒中心的位移, lj为连接颗粒 p中心与颗粒 q中心的分支向量. 值得注意的是, ϕi应假 定为 x'的一个特定函数与 xt的任意函数之和29]. 由此,可以求得细观的位移梯度

ψijψij=ϕi,j=ϕix'j(2)

此时, ψij只是 xt的函数. 那么,通过宏细观的位移梯度,可以定义一个相对变形量

γijγij=ψij-ϕ̅i,j(3)

式中, ϕ̅i,j为宏观位移 ϕ̅i的梯度,并且,假定 ϕ̅ixt的函数,即 ϕ̅i=ϕ̅ix,t. 那么,宏观应变可以求得为

εij=ϕ̅(i,j)12ϕ̅i,j+ϕ̅j,i(4)

式中, ϕ̅(i,j)表示 ϕ̅i,j的对称部分.

将式(2)代入式(1)中,可以得到颗粒 pq的相对位移为

δipq=ϕip-ϕiq+eijkχjprkp-χjqrkq=ϕ̅i,jlj+γijlj+eijkχj,lllrk(5)

式中, eijk为置换张量, χi表示颗粒转动, rk为颗粒的半径. 此时,将相对位移分解,分别定义为

δiM=ϕ̅i,jlj,δim=γijlj,δiR=eijkχj,lllrk(6)

其中, δiM表示宏观平均应变引起的位移, δim表示颗粒位移波动的梯度引起的位移, δiR表示可以转动引起的位移,所以,颗粒的相对位移可以看作由宏观上的平均位移,相对波动的位移,以及颗粒相对转动的位移3个部分组成. 同理,颗粒 p的转动也可以用颗粒 q的转动泰勒级数展开

χip=χiq+χi,jqlj+(7)

那么,这两个颗粒间的相对转动为

θipq=χi,jlj(8)

与位移的分解相似,颗粒的转动可以分解为整个体积元的转动与相对于体积元的转动两个部分,并定义这个相对的转动的梯度

αijαij=(χi-Γi),j(9)

其中, Γi是整个体积元的转动,并且有

Γi=eijkϕk,j(10)

那么,式(8)中的相对转动可以表示为

θiR=θipq=χi,jlj=αijlj+eijkφk,jplp(11)

对其进行分解,令

θis=αijlj,θiu=eijkϕk,jplp(12)

以求得微曲率

κij=χi,j(13)

图1   物质点P与其微体积元V

Fig. 1   Material point P and its microstructural volume element V

2 颗粒的宏细观本构关系

颗粒的宏观变形能密度 W可以表示为连续变形量的函数: W=Wεij,γij,κij,αij. 从而可以求得与这些变形量共轭的宏观应力量

τij=Wεij,σij=Wγij,μij=Wκij,νij=Wαij(14)

其中, τij为柯西应力, σij为相对应力, μij为偶应力, νij为相对偶应力. Tordesillas和Walsh等[45]认为颗粒材料的宏观变形能密度为颗粒接触间的细观变形能密度 w之和,而细观变形能密度 w是颗粒间相对运动的函数,即 w=wδiM,δim,δiR,θiR,θis,θiu. 那么,整个体积元的 V的变形能密度为

W=1VwδiM,δim,δiR,θiR,θis,θiu(15)

对细观变形能密度 w求偏导,可以得到对应细观运动量的颗粒接触力 fi与接触力矩 mi

wδiξ=fiξ,wθiη=miη(16)

其中, ξ=M,m,R, η=R,s,u. 式(15)和式(16)代入式(14)中,宏观应力可以用细观尺度上的颗粒接触力来表示

τij=Wεij=1VwδkMδiklj=1VfiMlj(17)

σij=Wγij=1Vwδkmδiklj=1Vfimlj(18)

μij=Wκij=1VwδlRδlRκij+wθlRθlRκij=1VflRelikljrk+miRlj(19)

νij=Wαij=1Vwθlsθlsαij=1Vmislj(20)

为建立颗粒的细观本构关系,首先对一个颗粒接触对建立局部坐标系,如图2所示[46]. 向量 n是接触面的法向向量,另外两个正交向量 st在接触面切向平面上

n=cosαe1+sinαcosβe2+sinαsinβe3(21)

s=dndα=-sinαe1+cosαcosβe2+cosαsinβe3(22)

t=n×s=-sinβe2+cosβe3(23)

其中, e1, e2, e3分别为平行于坐标轴的3个单位向量.

忽略法向与切向的交叉项,细观变形能密度 w可以简单写为

w=ξfnξδnξ+fsξδsξ+ftξδtξ+mnηθnη+msηθsη+mtηθtη(24)

图2   颗粒接触的局部坐标[46]

Fig. 2   Local coordinate system at a contact of particles[46]

其中, ξ=M,m,R; η=R,s,u. 下标 n, s, t分别为上述局部坐标系的分量. 则,进一步对于线弹性各向同性材料,细观变形能密度为

w=12ξKnξδnξ2+Ksξδsξ2+Ktξδtξ2+Gnηθnη2+Gsηθsη2+Gtηθtη2(25)

其中, Knξ, Ksξ, KtξGnη, Gsη, Gtη分别为对应于上文分解后的位移与转动的接触刚度. 对于一个接触对的两个圆形颗粒来说,接触面是一个圆形区域,那么 st方向的接触刚度应该一致,则令

Ksξ=Ktξ=Kwξ,Gsη=Gtη=Gwη(26)

通过旋转张量将局部坐标转换为全局坐标,颗粒的细观本构关系可以写为

fiξ=Kijξδjξ,Kijξ=Knξninj+Kwξsisj+titj(27)

miη=Gijηθjη,Gijη=Gnηninj+Gwηsisj+titj(28)

式中, i,j=1,2,3. 将式(27)与式(28)代入式(17) ~式(20)中,可以求出基于细观信息表示的宏观本构关系为

τij=1VfiMlj=1VKikMδkMlj=1VKikMεkllllj=1VKikMllljεkl=CijklMεkl(29)

σij=1Vfimlj=1VKikmδkmlj=1VKikmγkllllj=1VKikmllljγkl=Cijklmγkl(30)

μij=1VflRelikljrk+miRlj=1VKlmRδmRelikljrk+GikRθjRlj=1VKnmRemkqχk,lllrqenipljrp+GikRχk,llllj=

1VKnmRemkqllrqenipljrp+GikRllljχk,l=BijklR+CijklRχk,l(31)

νij=1Vmislj=1VGiksθkslj=1VGiksαkllllj=1VGiksllljαkl=Cijklsαkl(32)

3 颗粒的宏观本构参数识别

从式(29) ~式(32)可以看到,宏观本构参数是接触刚度和分支向量的函数. 对于一个体积元 V来说,其中的每个接触对的接触刚度和分支向量通常是不同的,为了简化计算,假设材料为各向同性,Chang等 [47]提出了一个方向分布密度函数 ξ,并令 ξα,β=1/4π,这样可以保证下式成立

0π02πξsinαdβdα=1(33)

由此能够将式(29) ~式(32)的离散求和转化为积分形式

CijklM=1VKikMlllj=l2NV0π02πKikMnlnjξsinαdβdα(34)

Cijklm=1VKikmlllj=l2NV0π02πKikmnlnjξsinαdβdα(35)

BijklR=1VKnmRemkqllrqenipljrp=l2r2NV0π02πKnmRemkqnlnqenipnjnpξsinαdβdα(36)

CijklR=1VGikRlllj=l2NV0π02πGikRnlnjξsinαdβdα(37)

Cijkls=1VGikslllj=l2NV0π02πGiksnlnjξsinαdβdα(38)

式中, NV表示体积元的体积密度. 那么,求解上述积分,就可以得到宏观本构模量

Ciiiiξ=l2NV153Knξ+2Kwξ,Ciijjξ=l2NV15Knξ-KwξCijijξ=l2NV15Knξ+4Kwξ,Cijjiξ=l2NV15Knξ-KwξCijklξ=0,otherwise(39)

Ciiiiη=l2NV153Gnη+2Gwη,Ciijjη=l2NV15Gnη-GwηCijijη=l2NV15Gnη+4Gwη,Cijjiη=l2NV15Gnη-GwηCijklη=0,otherwise(40)

BiiiiR=l2r2NV152KwR,BiijjR=l2r2NV15-KwRBijijR=l2r2NV154KwR,BijjiR=l2r2NV15-KwRBijklR=0,otherwise(41)

式中, ξ=M,m,R; η=R,s; i,j=1,2,3ij.

4 平衡方程

根据能量守恒,外力做功转化为动能与势能,则,由哈密顿原理可以得到

δ0tT-Wdt+0tδWextdt=0(42)

其中, W为总势能, T为总动能, Wext为外力功.

上文中假设总势能密度 W是宏观应变量的函数,根据式(14),总势能密度 W的变分应为

δW=τijδεij+σijδγij+μijδκij+νijδαij=

τijδϕ̅i,j+σijδϕi,j-δϕ̅i,j+μijδχi,j+νijδχi,j-eipqδϕq,pj=

τij-σijδϕ̅i,j+σijδψij+(μij+νij)δχi,j-νijeipqδψpq,j=

τij-σijδϕ̅i,j-τij-σij,jδϕ̅i+σijδψij+

μij+νijδχi,j-(μij+νij),jδχiνijeipqδψpq,j+νkl,lekijδψij(43)

那么,可以求得总势能 W的变分

δW=VδWdV=-Vτij-σij,jδϕ̅idV-

V(μij+νij),jδχidV+Vσij+νkl,lekijδψijdV+

Sτij+σijnjδϕ̅idS+Sμij+νijnjδχidS-SνijeipqnjδψpqdS(44)

同时,根据Mindlin的微结构理论[29],求得动能的变分为

δ0tTdt=-0tdtV(ρϕ̅̈iδϕ̅i+13ρ'd2ψ̈ijδψij+23ρ'd2χ̈iδχi)dV (45)

其中, d为体积元的边长, ρ为宏观质量密度, ρ'为细观质量密度. 式(42)中的外力功的变分可以定义为[29]

δWext=Vfiδϕ̅idV+VmiδχidV+VΦijδψijdV+Stiδϕ̅idS+SMiδχidS+STijδψijdS(46)

其中, fi, mi, Φij分别为单位体积上的体力、力矩、二阶体力, ti, Mi, Tij分别为单位面积上的面力、力矩、二阶面力.

从而,将式(43) ~式(46)代入式(42)中,可以得到

Vτij-σij,j+fi-ρϕ̅̈iδϕ̅idV+V(μij+νij),j+mi-23ρ'd2χ̈iδχidV+V-σij-νkl,lekij+Φij-13ρ'd2ψ̈ijδψijdV+Sti-τij-σijnjδϕ̅idS+SMi-μij+νijnjδχidS+STi+νklekijnlδψijdS(47)

为了使式(47)成立,要求式中每一项均为0,那么,就可以得到3组平衡方程的变分

τij-σij,j+fi=ρϕ̅i(μij+νij),j+mi=23ρ'd2χ̈iΦij-σij-νkl,lekij=13ρ'd2ψ̈ij(48)

以及3组边界条件

τij-σijnj=tiμij+νijnj=MiTij+νklekijnl=0(49)

5 运动方程与波的传播

将式(29) ~式(32)代入上文的平衡方程(48)中,得到3组15个运动方程 CijklM+Cijklmϕ̅k,lj-Cijklmψkl,j=ρϕ̅iBijklR+CijklR+Cijklsχk,lj-Cijklsekmnψmn,lj=23ρ'd2χ̈iCijklmϕ̅k,l-Cijklmψkl-Cklmnsekijχm,nl+Cklmnsekijempqψpq,nl=13ρ'd2ψ̈ij(50)

其中,忽略 fi, mi, Φij项. 并且,假定波沿 x1方向传播,那么为了求解平面波,要求运动量只是 x1的函数

ϕ̅i=ϕ̅ix1,tχi=χix1,tψij=ψijx1,t(51)

将式(51)代入式(50)中,得到这15个运动方程具体的分量形式

ρϕ̅1=C1111M+C1111mϕ̅1,11-C1111mψ11,1+C1122mψ22,1+C1133mψ33,1ρϕ̅2=C2121M+C2121mϕ̅2,11-C2112mψ12,1+C2121mψ21,1ρϕ̅3=C3131M+C3131mϕ̅3,11-C3113mψ13,1+C3131mψ31,1(52)

23ρ'd2χ̈1=B1111R+C1111R+C1111sχ1,11-C1111sψ23,11-ψ32,1123ρ'd2χ̈2=B2121R+C2121R+C2121sχ2,11-C2121sψ31,11-ψ13,1123ρ'd2χ̈3=B3131R+C3131R+C3131sχ3,11-C3131sψ12,11-ψ21,11(53)

13ρ'd2ψ̈11=C1111mϕ̅1,1-C1111mψ11+C1122mψ22+C1133mψ33(54a)

13ρ'd2ψ̈22=C2211mϕ̅1,1-C2211mψ11+C2222mψ22+C2233mψ33(54b)

13ρ'd2ψ̈33=C3311mϕ̅1,1-C3311mψ11+C3322mψ22+C3333mψ33(54c)

13ρ'd2ψ̈12=C1221mϕ̅2,1-C1212mψ12+C1221mψ21-C3131sχ3,11+C3131sψ12,11-ψ21,11(54d)

13ρ'd2ψ̈13=C1331mϕ̅3,1-C1313mψ13+C1331mψ31+C2121sχ2,11-C2121sψ31,11-ψ13,11(54e)

13ρ'd2ψ̈21=C2121mϕ̅2,1-C2112mψ12+C2121mψ21+C3131sχ3,11-C3131sψ12,11-ψ21,11(54f)

13ρ'd2ψ̈23=-C2323mψ23+C2332mψ32-C1111sχ1,11+C1111sψ23,11-ψ32,11(54g)

13ρ'd2ψ̈31=C3113mϕ̅3,1-C3113mψ13+C3131mψ31-C2121sχ2,11+C2121sψ31,11-ψ13,11(54h)

13ρ'd2ψ̈32=-C3223mψ23+C3232mψ32+C1111sχ1,11-C1111sψ23,11-ψ32,11(54i)

,式(52) ~式(54)的解是下列简谐波函数形式

ϕ̅i=Aiiexpiξx1-ωtχi=Biexpiξx1-ωtψij=Cijexpiξx1-ωt(55)

其中, ξ为单位长度内的波数,即波长的倒数, ω为角频率, Ai, Bi, Cij为振幅. 然后,将式(52) ~式(54)的15个方程进行分组,具体分为5个方程组

(A) 式(54b)与式(54c)

13ρ'd2ψ̈22=C2211mϕ̅1,1-C2211mψ11+C2222mψ22+C2233mψ3313ρ'd2ψ̈33=C3311mϕ̅1,1-C3311mψ11+C3322mψ22+C3333mψ33

(B)式(54g)与式(54i)

13ρ'd2ψ̈23=-C2323mψ23+C2332mψ32-C1111sχ1,11+C1111sψ23,11-ψ32,1113ρ'd2ψ̈32=-C3223mψ23+C3232mψ32+C1111sχ1,11-C1111sψ23,11-ψ32,11

(C) 式(53a), 式(54g)与式(54i)

23ρ'd2χ̈1=B1111R+C1111R+C1111sχ1,11-C1111sψ23,11-ψ32,11

13ρ'd2ψ̈23=-C2323mψ23+C2332mψ32-C1111sχ1,11+C1111sψ23,11-ψ32,11

13ρ'd2ψ̈32=-C3223mψ23+C3232mψ32+C1111sχ1,11-C1111sψ23,11-ψ32,11(D) 式(52b), 式(53c), (54d)与(54f)

ρϕ̅̈2=C2121M+C2121mϕ̅2,11-C2112mψ12,1+C2121mψ21,1

23ρ'd2χ̈3=B3131R+C3131R+C3131sχ3,11-C3131sψ12,11-ψ21,11

13ρ'd2ψ̈12=C1221mϕ̅2,1-C1212mψ12+C1221mψ21-C3131sχ3,11+C3131sψ12,11-ψ21,11

13ρ'd2ψ̈21=C2121mϕ̅2,1-C2112mψ12+C2121mψ21+C3131sχ3,11-C3131sψ12,11-ψ21,11

或者式(52c), (53b), (54e)与式(54h)

ρϕ̅̈3=C3131M+C3131mϕ̅3,11-C3113mψ13,1+C3131mψ31,1

23ρ'd2χ̈2=B2121R+C2121R+C2121sχ2,11-C2121sψ31,11-ψ13,11

13ρ'd2ψ̈13=C1331mϕ̅3,1-C1313mψ13+C1331mψ31+C2121sχ2,11-C2121sψ31,11-ψ13,11

13ρ'd2ψ̈31=C3113mϕ̅3,1-C3113mψ13+C3131mψ31-C2121sχ2,11+C2121sψ31,11-ψ13,11

(E) 式(52a), 式(54a), 式(54b)与式(54c)

ρϕ̅̈1=C1111M+C1111mϕ̅1,11-C1111mψ11,1+C1122mψ22,1+C1133mψ33,1

13ρ'd2ψ̈11=C1111mϕ̅1,1-C1111mψ11+C1122mψ22+C1133mψ33

13ρ'd2ψ̈22=C2211mϕ̅1,1-C2211mψ11+C2222mψ22+C2233mψ33

13ρ'd2ψ̈33=C3311mϕ̅1,1-C3311mψ11+C3322mψ22+C3333mψ33

再将式(55)代入这5个方程组,就可以求得关于波数 ξ与频率 ω的频散方程. (1)对于方程组(A)与(B),波沿 x1方向传播,形变方向为 x2x3方向,且其形变由剪切应变组成,即可称为横向剪切波;(2)对于方程组(C),形变方向同样 为 x2x3方向( χ1表示以 x1为轴进行旋转的转动自由度),且形变由剪切应变与旋转组成,可称为横向旋转波;(3)对于方程组(D),波沿 x1方向传播,形变方向沿 x2或者 x3方向传播,可视为经典波动理论的横波概念;(4)对于方程组E,波沿 x1方向传播,形变方向沿 x1方向传播,即为经典波动理论的纵波:

(A)横向剪切波

13ρd2ω2=C2222m-C2233m(56)

(B)横向剪切波

13ρd2ω2=C2323m-C2332m(57)

(C)横向旋转波

B1111r+C1111r+C1111s-23ρ'd2ω2C2323m-C2332m-13ρ'd2ω2+2C1111sξ2=2C1111sξ22(58)

(D)横波

C3131M+C3131mξ2-ρω2C3113m-C3131m2ξC3113m+C3131m2ξ0C1331m-C3131mξ2C2121sξ2+C1313m-C1331m-13ρ'd2ω202C2121sξ2C1331m+C1313mξ0C1331m+C1313m-13ρ'd2ω200C2121sξ20B2121R+C2121R+C2121sξ2-23ρ'd2ω2=0(59)

(E)纵波

C1111M+C1111mξ2-ρω2C1111m-C1122mξC1122mξC1111mξC1111m-C1122m-13ρ'd2ω2C1122mC1111m+2C1122mξ0C1111m+2C1122m-13ρ'd2ω2=0(60)

其中,式(58) ~式(60)中波数 ξ与频率 ω的关系比较复杂,具体的解并未给出. 由于假设材料各向同性,式(56)与式(57)求得的频率是相等的. 令 ξ=0,上述各种波的截止频率为:

(A)/(B) 横向剪切波 ωξ=0=ω2(61)

(C) 横向旋转波 ωξ=0=ω3ω4(62)

(D) 横波 ωξ=0=ω2ω30(63)

(E) 纵波 ωξ=0=ω1ω20(64)

其中,将截止频率从大到小排列有 ω1>ω2>ω3>ω4>0,且

ω1=l2NVρ'd2Knm(65a)

ω2=l2NVρ'd22Knm+3Kwm5(65b)

ω3=l2NVρ'd2Kwm(65c)

ω4=l2NVρ'd22r2KwR+3GnR+2GwR+3Gns+2Gws10(65d)

对于行波来说,可根据截止频率是否为零分为两类:截止频率为零的声波(acoustic wave)与截止频率非零的光波(optical wave). 从式(61) ~式(64)可以看到,横波与纵波均有一个声波分量和两个光波分量,横向旋转波具有两个光波分量,横向剪切波的形式为一条光波. 此外,截止频率与颗粒接触尺度上的细观参数相关,当接触刚度取为零时,所有波型均变为声波.

6 结果与讨论

上一部分中推导了横向剪切波、横向旋转波、横波、纵波的频散方程,根据此频散方程,可以求得波数与频率的频散关系. 从频散方 程可以看到,频散关系与细观参数相关. 通常来讲,基于细观尺度的颗粒材料研究一般取法向接触刚度为10 ~1 000 MN/m量级,切向接触刚度取为法向接触刚度的0.5倍,滚动刚度取为100 Nm量级. 理论研究中,颗粒的尺度常取为1 mm量级. 本文选取的参数来自Misra等的研究[28],如表1所示,那么,上文中非零的截止频率分别为: ω1=8.19×106rad/s, ω2=6.86×106rad/s, ω3=5.74×106rad/s, ω4=77.5rad/s.

表1   细观参数取值

Table 1   Values of micro parameters

ParametersValuesParametersValuesParametersValues
l,d10-3 mKnM200 MN/mKwM100 MN/m
r0.5×10-3 mKnm200 MN/mKwm100 MN/m
NV10-3 mKnR0.02 MN/mKwR0.01 MN/m
ρ1570 kg/m3GnR2×10-8 MNmGwR10-8 MNm
ρ'3000 kg/m3Gns2×10-8 MNmGns10-8 MNm

新窗口打开

图3分别给出了表示横向剪切波、横向旋转波、横波与纵波的频散曲线. 从图3(a)可以看到,对于横向剪切波来讲,只有一条呈平直的光波型的频散曲线,在本文参数取值下, 无论波数为多大,频率始终为 6.86×106rad/s,保持截止频率不变化,即此介质无论只能传递频率为 6.86×106rad/s 的横向剪切波. 同样地,在图3(b)中,横向旋转波也具有一条呈平直的光波型的频散曲线,保持频率为 5.74×106rad/s,此外,当波数较小时,有一条频率较小(小于80 rad/s)且迅速衰减的双曲线型光波的频散曲线,表明频率为 5.74×106rad/s 与0 ~77.5 rad/s的横向旋转波能够传递. 从图3(c)可以看到,横波具有四条频散曲线:第1条为低频线性增长的声波,最后达到 7.26×105rad/s;第2条频率从0开始随波数增大线性增长,后增速减小,频率逐渐趋向平稳达到 4.46×106rad/s;第3条起始频率(即截止频率)为 5.74×106rad/s,随后略微增大达到稳定的频率 6.34×106rad/s;第4条起始频率为 6.86×106rad/s,随后快速地呈双曲线式增长. 可以看到,第2条与第3条曲线之间、第3条与第4条曲线之间分别存在一个狭长带,期间没有任何频率的波通过,即为横波的频率带隙: 4.46 ×10 6~5.74×106rad/s 与6.34 ×10 6~6.86×106rad/s. 对于纵波,如图3(d),具有3条频散曲线,除不存在低频的曲线外,其变化规律与横波基本相同:第1条频散曲线达到稳定频率 5.13×106rad/s,第2条频散曲线的频率从 6.86×106rad/s增长为 7.75×106rad/s, 第3条从截止频率为 8.19×106rad/s呈双曲线式快速增大,且有存在两条频率带隙:5.13 ×10 6~6.86×106rad/s与7.75 ×10 6~8.19×106rad/s. 此频散曲线的产生与式(56) ~式(60)的频散方程直接相关. 其中,横向剪切波只受与细观的平动项相关的参数的影响,纵波的频散由与宏观和细观平动项相关的参数共同影响,而横向 旋转波由细观的平动与转动项控制频散现象,横波的频散则受宏--细观的平动与转动项的共同影响. 所以,各种形式的波的频散现象会存在差异. 值得注意的是,已有物理实验与离散元模拟[34,48-49]证实六方密排结构中横向旋转波与频率带隙的 存在.

图3   横向剪切波、横向旋转波、横波与纵波的频散曲线

Fig. 3   Dispersion curves for transverse shear wave, transverse rotational wave, transverse wave and longitudinal wave

图4   所建议模型与Misra模型[28]的横向剪切波、横向旋转波、横波与纵波的频散曲线对比图

Fig. 4   Dispersion curves of the proposed model vs. Misra model [28] for transverse shear wave, transverse rotational wave, transverse wave and longitudinal wave

图4所示为所建议模型与Misra模型 [28]的横向剪切波、横向旋转波、横波与纵波的频散曲线对比图. 对于横向剪切波,频 散曲线基本重合. 对于横向旋转波,本文所建议模型多具有1条低频的频散曲线,而另1条基本重合. 对于横波,本文所建议模型多具有1条较低频的直线型频散曲线,而另外3条与Misra模型 [28]的频散曲线的变化规律相近,其中1条频散曲线基本重合,所建议模型低频的频率带隙比Misra模型[28]更窄,高频的频率带隙基本一致. 对于纵波,两个模型均具有3条频散曲线,且相互重合,频率带隙基本一致. 横向旋转波与横波同Misra模型所预测存在差异,而其他形式的波基本一致,主要因为本文所建议模型多考虑了一个独立的转动自由度,因而具备了一个考虑转动的变形模式,通过式(52)~式(54)可以反映出来,转动自由度主要影响横向剪切波与横波的频散方程.

图5所示为所建议模型预测频率带隙的宽度. 从图中可以看到,综合考虑横向剪切波、横向旋转波、横波与纵波的频散曲线,存在两个狭长带(图中红色与蓝色阴影区域),无论波数如何变化,期间均无任何频率的波能够在此介质传播,即此介质的频率带隙为5.13 ×106~5.74×106rad/s 与6.34 ×10 6~6.86×106rad/s. 对比于Misra模型[28]预测的频率带隙,总体上,频率带隙的宽度相一致,但横波的频率带隙的宽度,本文建议模型所预测要窄于Misra模型,这还是由于所建议模型具有多考虑转动的变形模式的存在所引起的.

图5   所建议模型预测频率带隙

Fig. 5   The frequency band gaps predicted by the proposed model

7 结论

颗粒材料的细观力学方法已经被应用于发展颗粒材料的微形态模型. 基于此方法发展的颗粒材料宏观连续体模型可以反映细观结构信息. 本文针对颗粒材料提出了一个基于细观力学的微形态模型. 首先采用Mindlin-Eringen的微结构理论 [29-30]进行运动分析,对颗粒的平动与转动均进行分解,分解为宏观平均运动(平动与转动)和细观真实运动(平动与转动)两部分,并认为颗粒的细观真实运动(平动与转动)是由宏观平均运动(平动与转动)与相对运动(平动与转动)组成. 基于此分解,求得相应的细观本构关系和细观变形能函数,并且宏观变形能由细观变形能求和得到. 由此可获得基于细观量表示的宏观本构关系与相应的本构模量.

基于所建议的微形态模型预测和探讨波在弹性颗粒介质中传播的现象,考察了纵波、横波、横向剪切波、横向旋转波的频散关系,并对其频率带隙进行预测,主要结论如下:

(1)所建议模型能够模拟频散现象. 由于纵波、横波、横向剪切波、横向旋转波的频散关系分别受宏-细观尺度上的与平动项与转动项相关的参数的影响,其频散曲线各不相同. 横向剪切波具有1条平直的频散曲线;横向旋转波具有1条平直的频散曲线与1条呈双曲线衰减型频散曲线;纵波具有3条频散曲线,且与Misra [28]所预测的频散曲线基本重合;由于所建议模型具有考虑了转动自由度的变形模式,横波具有四条频散曲线,比Misra [28]所预测的多出一条低频的呈线性增长的频散曲线,另外3条变化规律相似且其中一条频散曲线基本重合.

(2)所建议模型预测了弹性颗粒介质的频率带隙. 无论波数如何变化,频率范围为5.13 ×10 6~5.74×106rad/s 与6.34 ×10 6~6.86×106rad/s 的波均不能在此介质中传播,即为所建议模型预测的频率带隙.

The authors have declared that no competing interests exist.


参考文献

[1] Li XK, Chu XH, Feng YT.

A discrete particle model and numerical modeling of the failure modes of granular materials

.Engineering Computations, 2005, 22(8): 894-920

[本文引用: 1]     

[2] Chu XH, Yu C, Xiu CX, et al.

Two scale modeling of behaviors of granular structure: Size effects and displacement fluctuations of discrete particle assembly

.Structural Engineering & Mechanics, 2015, 55(2): 315-334

[本文引用: 1]     

[3] 冯春, 李世海, 刘晓宇.

基于颗粒离散元法的连接键应变软化模型及其应用

. 力学学报, 2016, 48(1): 76-85

[本文引用: 2]     

(Feng Chun, Li Shihai, Liu Xiaoyu.

Particle-dem based linked bar strain softening model and its application

.Chinese Journal of Theoretical and Applied Mechanics, 2016, 48(1): 76-85(in Chinese))

[本文引用: 2]     

[4] 季顺迎, 孙珊珊, 陈晓东.

颗粒材料剪切流动状态转变的环剪试验研究

. 力学学报, 2016, 48(5): 1061-1072

[本文引用: 1]     

(Ji Shunying, Sun Shanshan, Chen Xiaodong.

Shear cell test on transition of shear flow states of granular materials

.Chinese Journal of Theoretical and Applied Mechanics, 2016, 48(5): 1061-1072(in Chinese))

[本文引用: 1]     

[5] 白以龙, 汪海英, 夏蒙棼.

固体的统计细观力学—连接多个耦合的时空尺度

. 力学进展, 2006, 36(2): 286-305

[本文引用: 1]     

(Bai Yilong, Wang Haiying, Xia Mengfen, et al.

Statistical mesomechanics of solid, linking coupled multiple space and time scales

.Advances in Mechanics, 2006, 58(6): 286-305 (in Chinese))

[本文引用: 1]     

[6] Torquato S, Haslach H.

Random heterogeneous materials: Microstructure and macroscopic properties

.Applied Mechanics Reviews, 2002, 55(4): B62

[本文引用: 1]     

[7] Chang CS, Ma L.

Modeling of discrete granulates as micropolar continua

.Journal of Engineering Mechanics, 1990, 116(12): 2703-2721

[本文引用: 2]     

[8] Eringen AC, Suhubi ES.

Nonlinear theory of simple micro-elastic solids—I

.International Journal of Engineering Science, 1964, 2(2): 189-203

[本文引用: 2]     

[9] Eringen AC, Suhubi ES.

Nonlinear theory of simple micro-elastic solids—II

.International Journal of Engineering Science, 1964, 2(4): 389-404

[本文引用: 1]     

[10] Stojanovi R.

On the mechanics of materials with microstructure

.Acta Mechanica, 1972, 15(3-4): 261-273

[本文引用: 1]     

[11] Mühlhaus HB.

Application of Cosserat theory in numerical solutions of limit load problems

.Archive of Applied Mechanics, 1989, 59(2): 124-137

[本文引用: 1]     

[12] Li XK, Tang HX.

A consistent return mapping algorithm for pressure-dependent elastoplastic Cosserat continua and modelling of strain localisation

.Computers & Structures, 2005, 83(1): 1-10

[本文引用: 1]     

[13] Granik VT, Ferrari M.

Microstructural mechanics of granular media

.Mechanics of Materials, 1993, 15(4): 301-322

[本文引用: 2]     

[14] Christoffersen J, Mehrabadi MM, Nematnasser S.

A micromechanical description of granular material behavior

.Journal of Applied Mechanics, 1981, 48(2): 339

[本文引用: 3]     

[15] Mehrabadi MM, Nemat-Nasser S.

Stress, dilatancy and fabric in granular materials

.Mechanics of Materials, 1983, 2(2): 155-161

[本文引用: 2]     

[16] 唐洪祥, 李锡夔.

Cosserat连续体模型中本构参数对应变局部化模拟结果影响的数值分析

. 计算力学学报, 2008, 25(5): 676-681

[本文引用: 1]     

(Tang Hongxiang, Li Xikui.

Numerical analysis for the effects of constitutive parameters in Cosserat continuum model on the simulation results of the strain localization

.Chinese Journal of Computational Mechanics, 2008, 25(5): 676-681 (in Chinese))

[本文引用: 1]     

[17] Li XK, Liu QP, Zhang JB.

A micro-macro homogenization approach for discrete particle assembly-Cosserat continuum modeling of granular materials

.International Journal of Solids and Structures, 2010, 47(2): 291-303

[本文引用: 1]     

[18] Chang JF, Chu XH, Xu YJ.

Finite-element analysis of failure in transversely isotropic geomaterials

. International Journal of Geomechanics, 2014, 15(6): 04014096

[本文引用: 1]     

[19] Gennes PG.

Granular matter: A tentative view

.Reviews of Modern Physics, 1999, 71(2): S374

[本文引用: 2]     

[20] Chang CS, Chang Y, Kabir MG.

Micromechanics modeling for stress-strain behavior of granular soils I: Theory

.Journal of Geotechnical Engineering, 1992, 118(12): 1959-1974

[本文引用: 2]     

[21] Chang CS, Kuhn MR.

On virtual work and stress in granular media

.International Journal of Solids & Structures, 2005, 42(13): 3773-3793

[本文引用: 1]     

[22] Saxcé GD, Fortin J, Millet O.

About the numerical simulation of the dynamics of granular media and the definition of the mean stress tensor

.Mechanics of Materials, 2004, 36(12): 1175-1184

[本文引用: 1]     

[23] Kruyt NP.

Statics and kinematics of discrete Cosserat-type granular materials

.International Journal of Solids & Structures, 2003, 40(3): 511-534

[本文引用: 1]     

[24] Bonelli S, Millet O, Nicot F, et al.

On the definition of an average strain tensor for two-dimensional granular material assemblies

.International Journal of Solids & Structures, 2012, 49(7-8): 947-958

[本文引用: 1]     

[25] Kruyt NP, Millet O, Nicot F.

Macroscopic strains in granular materials accounting for grain rotations

.Granular Matter, 2014, 16(6): 933-944

[本文引用: 1]     

[26] Chang CS, Liao CL.

Constitutive relation for a particulate medium with the effect of particle rotation

.International Journal of Solids & Structures, 1990, 26(4): 437-453

[本文引用: 1]     

[27] Chang CS, Gao J.

Second-gradient constitutive theory for granular material with random packing structure

.International Journal of Solids & Structures, 1995, 32(16): 2279-2293

[本文引用: 1]     

[28] Misra A, Poorsolhjouy P.

Granular micromechanics based micromorphic model predicts frequency band gaps

.Continuum Mechanics & Thermodynamics, 2016, 28(1-2): 215-234

[本文引用: 15]     

[29] Mindlin RD.

Micro-structure in linear elasticity

.Archive for Rational Mechanics and Analysis, 1964, 16(1): 51-78

[本文引用: 8]     

[30] Eringen AC.

Microcontinuum Field Theories: Foundations and Solids

. New York: Springer, 1999

[本文引用: 4]     

[31] Neff P, Ghiba ID, Madeo A, et al.

A unifying perspective: the relaxed linear micromorphic continuum

.Continuum Mechanics & Thermodynamics, 2014, 26(5): 639-681

[本文引用: 1]     

[32] Ghiba ID, Neff P, Madeo A, et al.

The relaxed linear micromorphic continuum: Existence, uniqueness and continuous dependence in dynamics

.Mathematics & Mechanics of Solids, 2015, 20(10): 1171-1197

[本文引用: 1]     

[33] Neff P, Ghiba ID, Lazar M, et al.

The relaxed linear micromorphic continuum: Well-posedness of the static problem and relations to the gauge theory of dislocations

.The Quarterly Journal of Mechanics and Applied Mathematics, 2015, 68(1): 53-84

[本文引用: 1]     

[34] Merkel A, Luding S.

Enhanced micropolar model for wave propagation in ordered granular materials

.International Journal of Solids & Structures, 2017, 106: 91-105

[本文引用: 2]     

[35] Maugin GA, Metrikine AV.

Mechanics of Generalized Continua: One Hundred Years After the Cosserats. Berlin, Heidelerg:

Springer-Verlag, 2010

[本文引用: 2]     

[36] Askes H, Aifantis EC.

Gradient elasticity in statics and dynamics: An overview of formulations, length scale identification procedures, finite element implementations and new results

.International Journal of Solids & Structures, 2011, 48: 1962-1990

[37] Muhlhaus HB, Oka F.

Dispersion and wave propagation in discrete and continuous models for granular materials

.International Journal of Solids & Structures, 1996, 33: 2841-2858

[本文引用: 1]     

[38] Huang WX, Sloan SW, Sheng DC.

Analysis of plane Couette shear test of granular media in a Cosserat continuum approach

.Mechanics of Materials, 2014, 69: 106-115

[本文引用: 2]     

[39] Tordesillas A, Muthuswamy M, Walsh SDC.

Mesoscale measures of nonaffine deformation in dense granular assemblies

.Journal of Engineering Mechanics, 2008, 134: 1095-1113

[本文引用: 1]     

[40] Jiang MJ, Yu HS, Harris D.

Kinematic variables bridging discrete and continuum granular mechanics

.Mechanics Research Communications, 2006, 33: 651-666

[本文引用: 1]     

[41] Lee JD, Wang XQ.

Generalized micromorphic solids and fluids

.International Journal of Solids & Structures, 2011, 49: 1378-1387

[本文引用: 1]     

[42] Hutter G.

Homogenization of a Cauchy continuum towards a micromorphic continuum

.Journal of the Mechanics and Physics, 2017, 99: 394-408

[43] Forest S.

Nonlinear regularization operators as a derived from the micromorphic approach to gradient elasticity, viscoplasticity and damage

.Proceedings of the Royal Society A, 2016, 472: 20150755

[本文引用: 1]     

[44] 马天雪, 苏晓星, 董浩文.

声光子晶体带隙特性与声光耦合作用研究综述

. 力学学报, 2017, 49(4): 743-757

[本文引用: 2]     

(Ma Tianxue, Su Xiaoxing, Dong Haowen, et al.

Review of bandgap characteristics and acousto-optical coupling in phoxonic crystals

.Chinese Journal of Theoretical and Applied Mechanics, 2017, 49(4): 743-757(in Chinese))

[本文引用: 2]     

[45] Tordesillas A, Walsh DCS.

Incorporating rolling resistance and contact anisotropy in micromechanical models of granular media

.Powder Technology, 2002, 124(1-2): 106-111

[本文引用: 1]     

[46] Hicher PY, Chang CS.

Elastic model for partially saturated granular materials

.Journal of Engineering Mechanics, 2008, 134(6): 505-513

[本文引用: 3]     

[47] Chang CS, Lun M.

Elastic material constants for isotropic granular solids with particle rotation

.International Journal of Solids & Structures, 1992, 29(8): 1001-10

[本文引用: 1]     

[48] Merkel A, Tournat V, Gusev V.

Dispersion of elastic waves in three-dimensional noncohesive granular phononic crystals: properties of rotational modes

.Physical Review E, 2010, 82(3 Pt 1): 031305

[本文引用: 1]     

[49] Merkel A, Tournat V, Gusev V.

Experimental evidence of rotational elastic waves in granular phononic crystals

.Physical Review Letters, 2011, 107(22): 225502

[本文引用: 1]     

/