深部岩体结构面动力特性与效应专刊简介
doi: 10.11883/bzycj-2025-0158
-
-
由颗粒、柱壳或球壳组成的链结构能够实现应力波的传输,并具备较好的波形操控能力,如波形的放大和衰减[1-3],因此可作为声学透镜[4]、减震材料[5-6]、声学开关[7]等,在工程领域有广泛的应用前景。链结构中的波传播过程因离散特性而显得十分复杂,一维链系统的波传播行为可能呈现线性、弱/强非线性的特征[8],不同的波传播特征主要受组成链的单个胞元基体材料性质和胞元结构(如胞元形状和胞元尺寸)的影响[9-12]。其中,由柱壳组成的链结构,因胞元中空且壁厚可控,在操控波形方面表现出优异的性能[9-10, 13-14]。因此,了解柱壳链中应力波的传播规律及弥散机理,对设计新型高效的波形控制器具有重要意义。
建立链结构中相邻胞元之间的接触关系,是研究柱壳链中波传播过程的基础。在一维均匀颗粒链中,Hertz定律常用来表征相邻颗粒之间接触力F和相对位移δ的关系,如F∝δ1.5 [15]。对于空心壳,接触力和相对位移的关系为F∝δn,其中n为表征非线性程度的接触指数,取值与空心壳壁厚及形状相关。对于球壳,n通常为1.2~1.5[10, 16]。Kim等[9]指出,对于椭圆柱壳,可能n<1,从而表现出应变软化的特征。
基于链结构中相邻胞元的接触关系,可以用离散元方法求解链结构中的波传播过程,即对每个胞元的力学守恒方程进行单独求解[9, 17]。为了方便分析以获得解析解,也可将离散方程连续化后采用应力波理论来研究[18-22]。Pochhammer-Chree精确理论可给出谐波激励下弹性圆杆中应力波传播的精确解,但其控制方程形式复杂,尤其对于瞬态冲击载荷,方程解析求解往往十分困难[23]。Rayleigh-Love近似理论中,假设杆在受到加载时,其纵向运动会带来杆的横向涨缩,并且假设杆中质点的横向位移沿径向线性分布,从而推导出计及横向惯性效应的修正的波动方程[22]。Rayleigh-Love理论具有较高精度,并成功地用于研究瞬态加载边界条件的波传播问题,如阶跃应力脉冲[24]、梯形应力脉冲[25]和指数衰减的边界速度载荷[26]。质量块初速冲击是一种生活中常见的冲击加载情形,但因其边界载荷随时间变化且需要根据牛顿运动定律进行反解,导致方程求解困难,尚未见相应的解析解。
本文中,在理论上分析冲击载荷下柱壳链中的弹性波传播过程及其几何弥散特性,将柱壳链等效为弹性均质方杆,并采用细观有限元方法构建柱壳链的细观结构对理论模型进行验证。首先,基于Rayleigh-Love横向惯性修正理论和量纲分析方法,获得柱壳链在质量块初速冲击下的无量纲控制方程。再进一步通过Laplace积分变换方法,对该控制方程进行求解。然后,将理论预测结果与细观有限元结果进行比较,并讨论该理论结果对于不同壁厚柱壳链的适用性。最后,基于理论模型的解析解,讨论柱壳链在初速冲击载荷下的波传播过程,和泊松比、惯性半径、冲击质量和速度对波形弥散的影响。
1. 理论分析和数值模拟
1.1 问题描述和模型简化
考虑一个刚性质量块M以初速度v0沿X方向撞击一个静止的柱壳链的一端,柱壳链的另一端固定,如图1(a)所示。所有柱壳排布成一列,相邻柱壳相互绑定,链总长度为L,含有N个柱壳。柱壳的壁厚为h,外径为D=L/N,垂直方向厚度为a。柱壳的基体材料为弹性材料,密度为ρs,弹性模量为Es,泊松比为νs。在X方向压缩时,柱壳在Y方向将产生不可忽略的变形,即Y方向的泊松比较大,而在面外Z方向的变形量很小,可以忽略,即Z方向的泊松比为0。
在质量块的冲击下,柱壳链中的应力波从近端传向远端。柱壳链的几何非线性将导致复杂的非线性应力波传播问题。为简化分析,假设柱壳链为由均匀的连续介质组成的杆,杆在X方向的长度为L,在Y方向的高度为D,在Z方向的宽度为a,如图1(b)所示。杆的横截面面积为A=aD,等效密度为ρ0=πρs(D−h)h/D2,等效弹性模量为E,在Y方向的等效泊松比为ν,而在Z方向的等效泊松比为0。
1.2 理论分析
1.2.1 计及横向惯性的波传播模型
Rayleigh-Love波动方程可用来描述弹性杆中的应力波传播过程和由横向惯性引起的弥散效应[22]。Rayleigh-Love波动方程描述了弹性杆中质点沿X方向的位移u(X, t)的演化过程:
∂2u∂t2=c20∂2u∂X2+ν2κ2∂4u∂X2∂t2 (1) 式中:t为时间,
c0=√E/ρ0 为初等理论(不考虑横向惯性效应)下杆中的纵波波速,к为截面对X轴的回转半径(惯性半径)。在本文中,因为Z方向等效泊松比为0,所以惯性半径к为:
κ=√1A∫AY2dYdZ=√1D∫D/2−D/2Y2dY=√36D (2) 初始条件为:
u(X,0)=0,∂u∂t(X≠0,0)=0,∂u∂t(0,0)=v0 (3) 边界条件为:
MA∂2u∂t2(0,t)=E∂u∂X(0,t),u(L,t)=0 (4) 以长度L、波速c0和密度ρ0为基本量构造特征参数,定义无量纲量:
ˉu=uL,ˉX=XL,ˉt=c0tL,ˉm=MAρ0L,ˉv0=v0c0,ˉκ=κL (5) 于是,可以得到质量块冲击下弹性杆中应力波传播的无量纲化控制方程:
{∂2ˉu∂ˉt2=∂2ˉu∂ˉX2+α∂4ˉu∂ˉX2∂ˉt20<ˉX<1,ˉt>0ˉu(ˉX,0)=0,∂ˉu∂ˉt(ˉX≠0,0)=0,∂ˉu∂ˉt(0,0)=ˉv0,ˉm∂2ˉu∂ˉt2(0,ˉt)=∂ˉu∂ˉX(0,ˉt),ˉu(1,ˉt)=0 (6) 式中:
α=ν2ˉκ2 。该控制方程由一个四阶偏微分方程和一些定解条件组成,其中含有2个无量纲自变量ˉX 和ˉt 以及3个无量纲参数α 、ˉv0 和ˉm 。尽管计及横向惯性修正的波动方程的形式一致[20, 25],但对于不同的定解条件,其解大为不同[24-26]。此外,文献中大都将Rayleigh-Love波动方程应用于均质圆杆中应力波传播行为的研究,而应用到Hopkinson杆实验以对波形进行修正[25]时,则未见用于链结构中的波传播过程及规律的研究。1.2.2 波传播模型的解析解
先对
ˉu(ˉX,ˉt) 进行Laplace变换,时域ˉt 变为像域s,即:U(ˉX,s)=∫∞0ˉu(ˉX,ˉt)e−sˉtdˉtRes>0 (7) 此时,方程(6)可以变作:
{s2U=d2UdˉX2+αs2d2UdˉX20<ˉX<1ˉm(−ˉv0+s2U(0,s))=dUdˉX(0,s)U(1,s)=0 (8) 式(8)的微分方程是关于
ˉX 的二阶常系数齐次方程,结合定解条件,得到:U(ˉX,s)=ˉmˉv0[exp(−sˉX√1+αs2)−exp(sˉX−2s√1+αs2)]√1+αs2s[1+ˉms√1+αs2+(1−ˉms√1+αs2)exp(−2s√1+αs2)] (9) 作Laplace逆变换:
ˉu(ˉX,ˉt)=12πi∫β+i∞β−i∞U(ˉX,s)esˉtds (10) 式中:β±i∞为复变量,β为实变数。将式(9)代入式(10),可得:
ˉu(ˉX,ˉt)=ˉmˉv02πi∫β+i∞β−i∞[exp(sˉt−sˉX√1+αs2)−exp(sˉt+sˉX−2s√1+αs2)]√1+αs2s[1+ˉms√1+αs2+(1−ˉms√1+αs2)exp(−2s√1+αs2)]ds (11) 式(11)的积分可以采用留数定理来求解。令被积函数的分母为零,可得式(11)有奇点0和简单奇点
sk=iηk(k=±1,±2,±3,⋯) ,其中实数ηk由此计算:ˉmηk√1−αη2k=cot(ηk/√1−αη2k) (12) 式中:
η2k≤1/α 。式(12)存在无穷多个根,且其正根和负根一一对应。值得注意的是,ηk只依赖于α和ˉm ,而与ˉv0 无关。由留数定理,有:ˉu(ˉX,ˉt)=ˉmˉv0Res[U(ˉX,ˉs)esˉt,0]+ˉmˉv0∑kRes[U(ˉX,ˉs)esˉt,sk]=ˉmˉv0lims→0U(ˉX,ˉs)esˉts+ˉmˉv0∑klims→ηkiU(ˉX,ˉs)esˉt(s−iηk)=ˉmˉv0∑k[1−αη2k)2(exp(iηkˉt−iηkˉX√1−αη2k)−exp(iηkˉt+iηk(ˉX−2)√1−αη2k)](ˉmηk√1−αη2k+i)−2ηk{ˉm(1−αη2k)[η2k(ˉm−2α)+1]+1} (13) 可见奇点0对上述结果没有贡献,而所有贡献来自于简单奇点
sk=iηk(k=±1,±2,±3,⋯) 。考虑式(12)中正根与负根的对应关系,式(13)的虚部可以相互抵消。因此,式(13)可以改写为:ˉu(ˉX,ˉt)=2ˉmˉv0∞∑k=1ηk>0(1−αη2k)2sin(ηkˉt)sin[ηk(1−ˉX)/√1−αη2k]ηk{1+ˉm(1−αη2k)[η2k(ˉm−2α)+1]}sin(ηk/√1−αη2k) (14) 进一步得出应变和无量纲速度:
ε=−∂ˉu(ˉX,ˉt)∂ˉX=2ˉmˉv0∞∑k=1ηk>0(1−αη2k)3/2sin(ηkˉt)cos[ηk(1−ˉX)/√1−αη2k]{1+ˉm(1−αη2k)[η2k(ˉm−2α)+1]}sin(ηk/√1−αη2k) (15) ˉv=∂ˉu(ˉX,ˉt)∂ˉt=2ˉmˉv0∞∑k=1ηk>0(1−αη2k)2cos(ηkˉt)sin[ηk(1−ˉX)/√1−αη2k]{1+ˉm(1−αη2k)[η2k(ˉm−2α)+1]}sin(ηk/√1−αη2k) (16) 不考虑横向惯性影响(α=0)时,式(14)~(16)可分别简化为:
ˉu(ˉX,ˉt)=2ˉmˉv0∞∑k=1ηk>0sin(ηkˉt)sin[ηk(1−ˉX)]ηk(1+ˉm+ˉm2η2k)sinηk (17) ε=−∂ˉu(ˉX,ˉt)∂ˉX=2ˉmˉv0∞∑k=1ηk>0sin(ηkˉt)cos[ηk(1−ˉX)](1+ˉm+ˉm2η2k)sinηk (18) ˉv=∂ˉu(ˉX,ˉt)∂ˉt=2ˉmˉv0∞∑k=1ηk>0cos(ηkˉt)sin[ηk(1−ˉX)](1+ˉm+ˉm2η2k)sinηk (19) 其中,ηk满足的方程(12)可以改写为:
ˉmηk=cotηk (20) 值得注意,式(14)~(19)均为无穷级数,其数值计算的结果受到求和项数目k的影响,k越大,结果越精确。因此,后续的分析为保证结果的准确性,k需足够大,本文中均取k>50。
1.3 数值模拟
1.3.1 细观有限元模型
采用细观有限元方法,建立质量块初速冲击柱壳链的有限元模型。对每个柱壳单独划分网格,再导入有限元软件ABAQUS/Explicit进行求解。本文中,柱壳链几何参数取值分别为:单个柱壳的壁厚h=0.56 mm,外径D=10 mm,厚度a=4 mm,柱壳链中共含30个柱壳,链总长L=300 mm。通过网格收敛性分析,网格尺寸取为0.3 mm,单个柱壳含1 920个六面体单元(C3D8R)。模拟中,所有的面接触方式均采用硬接触,设置摩擦因数为0.2。柱壳基体材料密度ρs=966 kg/m3,弹性模量Es=1 600 MPa,泊松比νs=0.3。柱壳链近端受到刚性质量块(M=1 g)初始速度v0=1 m/s的冲击,远端固定。不同时刻柱壳链中的应力云图如图2所示。结果表明,应力波从冲击端往支撑端沿着柱壳逐个传播,且应力主要集中在柱壳之间的接触区域。
1.3.2 柱壳链材料参数的确定
等效连续介质模型中的材料参数,可通过对柱壳链进行单轴准静态压缩测试获得。壁厚0.56 mm的柱壳链在单轴准静态恒速压缩(v=0.05 m/s)下的名义应力应变曲线如图3(a)。观察发现,其应力应变曲线呈现轻微软化。采用名义应力应变关系式
σN=BεnN 进行拟合,可得B=2.6 MPa,n=0.93。指数n接近于1,为了分析方便,本文中采用线性关系σN=EεN 进行拟合,可得等效弹性模量E=3.2 MPa,其中考虑名义应变的拟合范围为0<εN<0.1,当名义应变超过0.1时,应变软化较明显,线性拟合偏差较大。泊松比定义为弹性段横向名义应变εY的负值与轴向名义应变εX的比值,即ν=−εY/εX。观察发现,随着压缩的进行,泊松比并不能保持为一个恒定的值,如图3(b),为分析方便,取平均测得等效泊松比ν约为0.85。采用类似的方法,对壁厚为0.33、1.07 mm的柱壳链进行测试,测得对应的等效杨氏模量分别为0.3、10.8 MPa,而等效泊松比均约为0.85。2. 结 果
2.1 位移和载荷曲线
质量块(M=1 g)以初速v0=1 m/s冲击壁厚0.56 mm的柱壳链,其冲击端位移随着加载时间的变化,如图4(a)所示。结果表明,冲击端位移随着时间增大,然后逐渐趋于平缓,理论分析与数值模拟结果基本一致。理论分析与细观有限元模拟所得的支撑端和冲击端的载荷历史曲线如图4(b)所示。观察发现,理论分析能够较好地捕捉载荷曲线的整体特征,但曲线振荡存在一些差异,理论计算的载荷曲线呈现低频振荡,但在数值模拟中,载荷曲线呈现高频振荡,这是由于罚函数接触算法所致。由应力波理论,对于不考虑几何弥散的均质弹性杆,支撑端与冲击端峰值载荷的比为2。而该算例中,由理论分析得到的支撑端和冲击端的峰值载荷分别为2.1、1.3 N,支撑端与冲击端峰值载荷的比小于2。因此,柱壳链可作为一种有效的冲击缓冲器。
2.2 应变和速度分布
对于第2.1节的算例,理论分析和数值模拟给出的应变场和速度场如图5所示。需要说明的是,在数值模拟中,难以计算沿冲击方向连续变化的宏观应变场,因此定义一些局部点上的应变以检验理论结果。在过柱壳中心的YZ横截面上,取位于其表面的所有点作为该柱壳的参考点,提取这些点在冲击过程中的平均位移和速度,进一步通过相邻两个柱壳参考点的相对位移ΔX计算应变,即
ε=ΔX/D 。为方便比较,理论上也可以计算与数值模拟中对应参考点处的速度和应变。图5(a)为不同时刻的应变分布,结果表明,随着波往远端传播,应变峰值略有减小,且随着时间的推移,波后振荡的区域越来越大。理论计算与细观有限元模拟所得结果整体趋势一致,波后的波形振荡和波形前沿宽度也大致相同。不同时刻的速度分布如图5(b)所示。从速度分布可以发现,峰值速度随着波的传播逐渐降低,理论计算的粒子速度与细观有限元结果较好吻合。当质量块(M=1 g)冲击速度为1 m/s时,柱壳链中的应变较小,增大冲击速度可以检验应变较大时理论模型的适用性。图6为冲击速度为5、10 m/s时柱壳链中的应变分布。结果表明,理论分析与数值模拟的结果较好吻合,冲击速度越大,峰值应变越大,表明理论能够很好地预测峰值较大应变(约0.1)时柱壳链中不同时刻的应变分布。
进一步比较不同壁厚柱壳链的理论分析与数值模拟的结果。通过理论分析与细观有限元模拟计算,壁厚0.33、1.07 mm的柱壳链在受到质量块冲击(M=1 g, v0=1 m/s)时,链中的应变分布如图7所示。结果表明,理论与数值模拟很好吻合,理论能够很好地预测不同壁厚柱壳链中不同时刻的应变分布。
3. 讨 论
3.1 泊松比和无量纲惯性半径对波传播过程的影响
在理论分析中,可以发现泊松比ν和无量纲惯性半径
ˉκ 总是以组合的形式出现,因此可以直接讨论无量纲参数α=ν2ˉκ2 的影响。但是,为了更直观地考察泊松比ν和无量纲惯性半径ˉκ 的影响,这里分别改变这两个参数。图8(a)为不同泊松比某个时刻杆中的应变分布,其中无量纲惯性半径为0.04,泊松比为0、0.25、0.50、0.75和1.00。结果表明,泊松比为0时,应变分布曲线光滑无振荡,应变波形前沿宽度很窄。当泊松比不为0时,应变分布存在明显的振荡现象,且随着泊松比的增大,振荡幅值增大,应变波形前沿宽度增宽。此外,应变峰值也随着泊松比的增大而逐渐降低。以应变峰值处作为波阵面位置,泊松比越大,波速越小。类似地,固定泊松比为某个值,如ν=0.5,讨论不同惯性半径ˉκ 对变波形的影响。ˉκ 受柱壳半径的影响,图8(b)中ˉκ 由小到大分别为0.02、0.04和0.06。与泊松比的影响类似,惯性半径越大,应变峰值越小,波速越小,振荡幅值越大。综上可见,一维杆中横向惯性效应的影响越大,应变分布振荡越剧烈,并且波形前沿宽度增宽、应变峰值减小。3.2 无量纲速度和质量对波传播过程的影响
无量纲速度
ˉv0 和质量ˉm 对应变波形的影响如图9所示。可以看出,随着无量纲速度的增大,应变波形的振动波长不变,应变最大位置处对应的波速不变,波形前沿宽度不变,变化的仅是应变幅值。其实,由式(14)可以发现,无量纲速度仅出现在影响应变幅值的项。不同于无量纲速度,无量纲质量不仅影响应变幅值,也影响波形的波长。无量纲质量越大,应变幅值越大,但由图9(b)可见,质量在一定范围内变化时,对波形前沿宽度、应变最大位置处对应波速的影响均较小。3.3 弥散效应的分析
进一步分析Rayleigh-Love波动方程的弥散关系。在计及横向惯性的弹性杆中,纵波波速不再恒定[20, 25]。将行波解[20]:
u(X,t)=u0exp[iξ(X−ct)] (21) cc0=1√1+ν2κ2ξ2≈1−2π2ν2κ2λ2 (22) 式中:ξ=2π/λ,λ为波长。在圆杆的弥散分析[20, 25]中,其惯性半径
κ=√2D/4 。在本文中,长方形截面杆在垂直于冲击方向的面内各向异性,即面内不同方向上具有不同的泊松比,其惯性半径由式(2)为κ=√3D/6 。可见,波速与泊松比ν和D/λ相关,随着ν或D/λ的增大,波速越来越小,这与上述参数分析的结论一致。4. 结 论
建立柱壳链结构的等效连续介质模型,对冲击载荷作用下柱壳链中的应力波传播过程及其几何弥散特性进行了分析,并用细观有限元模型对理论模型进行验证。
(1)基于考虑横向惯性修正的Rayleigh-Love波动方程,分析了几何弥散对应力波传播的影响,给出了柱壳链在质量块初速冲击下的无量纲控制方程。基于Laplace变换及逆变换,获得了该控制方程的解析解。
(2)解析解的载荷和位移曲线、应变和速度分布均与细观有限元结果较好吻合,且适用于不同壁厚的柱壳链。
(3)讨论了泊松比、惯性半径、冲击质量和速度对冲击载荷作用下柱壳链中波形弥散的影响。结果表明,应变峰值、波形振荡幅度、波形前沿宽度和应变最大位置处对应的波速主要与泊松比和惯性半径相关,泊松比和惯性半径越大,应变峰值越小,波形振荡幅度越大,波形前沿宽度越宽,应变最大位置处对应的波速越小。
本文中建立的冲击载荷作用下柱壳链中的弹性波传播简化模型及解析解,可为揭示柱壳链的弥散机理和波形控制器的优化设计提供参考。
期刊类型引用(6)
1. 彭克锋,常白雪,虞吉林,郑志军. 椭圆柱壳链的波形弥散特性研究. 固体力学学报. 2023(06): 805-814 . 百度学术
2. 杨洪升,郑宇轩,周风华,李玉龙. 弹性SHPB试件内部纵向应力波传播及所产生的横向约束应力. 固体力学学报. 2023(06): 815-830 . 百度学术
3. 刘颢,白震,李志强,李世强. 冲击载荷下轻质夹芯拱最大刚度拓扑优化及动力响应. 爆炸与冲击. 2022(06): 100-114 . 本站查看
4. 杨智程,刘龙飞,刘炼煌,殷鹏志,吴志强. 外部爆炸载荷下表面粗糙度对45钢柱壳剪切带行为的影响. 高压物理学报. 2022(04): 114-126 . 百度学术
5. 彭克锋,郑志军,周风华,虞吉林. 密度梯度柱壳链的弹性波传播特性研究. 力学学报. 2022(08): 2131-2139 . 百度学术
6. 余同希,朱凌,许骏. 结构冲击动力学进展(2010-2020). 爆炸与冲击. 2021(12): 4-64 . 本站查看
其他类型引用(2)
-
计量
- 文章访问数: 302
- HTML全文浏览量: 31
- PDF下载量: 241
- 被引次数: 8