冲击载荷作用下柱壳链中的弹性波传播简化模型及其解析解

彭克锋; 崔世堂; 潘昊; 郑志军; 虞吉林

doi:10.11883/bzycj-2020-0246

冲击载荷作用下柱壳链中的弹性波传播简化模型及其解析解

doi: 10.11883/bzycj-2020-0246

1.
中国科学技术大学近代力学系中国科学院材料力学行为和设计重点实验室，安徽合肥 230027
2.
北京应用物理与计算数学研究所，北京 100088

基金项目: 国家自然科学基金(11772330, 11872360)；中央高校基本科研业务费专项(WK2480000003, WK2090050043)

详细信息

作者简介:
彭克锋（1994－　），男，博士研究生，pkf@mail.ustc.edu.cn

通讯作者:
郑志军（1979－　），男，副教授，zjzheng@ustc.edu.cn

中图分类号: O347.4
计量
- 文章访问数: 1636
- HTML全文浏览量: 835
- PDF下载量: 166
- 被引次数: 8
出版历程
- 收稿日期: 2020-07-17
- 修回日期: 2020-08-31
- 刊出日期: 2021-01-05

Simplified model of elastic wave propagation in cylindrical shell chain under impact load and its analytical solution

1.
CAS Key Laboratory of Mechanical Behavior and Design of Materials, Department of Modern Mechanics, University of Science and Technology of China, Hefei 230027, Anhui, China
2.
Institute of Applied Physics and Computational Mathematics, Beijing 100088, China

摘要

摘要: 柱壳链能引起波形的弥散，具备操控波形的潜力。建立了柱壳链结构的等效连续介质模型和细观有限元模型，研究了质量块冲击作用下柱壳链中的弹性应力波传播过程及其几何弥散特性。基于考虑横向惯性修正的Rayleigh-Love波动方程，建立了柱壳链在质量块冲击下的控制方程，采用Laplace变换及其逆变换获得了位移场、速度场和应变场的解析解，所得结果与细观有限元模拟结果较好吻合。结果表明，在冲击过程中应变和速度峰值均逐渐减小，应变峰值、振荡幅度和波形前沿宽度与泊松比和惯性半径相关，泊松比和惯性半径越大，应变峰值越小，应变分布振荡越剧烈，波形前沿宽度越宽。
- 柱壳链 /
- 应力波 /
- 横向惯性修正 /
- 解析解 /
- 几何弥散
Abstract: Cylindrical shell chain can cause dispersion of waveform and has the potential to manipulate waveform. An equivalent continuum model and a mesoscopic finite element model of cylindrical shell chain structure were established, and the stress wave propagation process and its geometric dispersion characteristics in cylindrical shell chains under mass impact were studied. For the in-plane compression of a single cylindrical shell, the deformation in the out-of-plane direction is very small and the in-plane deformation perpendicular to the loading direction is relatively large. Thus, the out-of-plane Poisson’s ratio can be taken as 0, but the in-plane one cannot be ignored. For the simplification of analysis, the cylindrical shell chain is considered as a rod composed of anisotropic homogeneous continuum. Based on the Rayleigh-Love wave equation with the transverse inertia correction, the governing equation of elastic wave propagation in a cylindrical shell chain under mass impact was obtained and rewritten in a dimensionless form. The analytical solutions of displacement, velocity and strain fields were obtained by using Laplace transform and its inverse transform and expressed in the form of infinite series. A chain with 30 cylindrical shells was constructed numerically and its dynamic impact behavior was simulated with finite element code ABAQUS/Explicit. The theoretical predictions of mechanical responses are in good agreement with the results of mesoscopic finite element simulation. During the impact process, the peak values of strain and velocity decrease gradually. The peak strain, the oscillation amplitude of waveform and the width of waveform front are related to the in-plane Poisson's ratio and the radius of gyration of the cylindrical shell chain. The larger the in-plane Poisson's ratio and the radius of gyration of the cylindrical shell chain, the smaller the peak strain, the stronger the oscillation of the strain waveform and the wider the width of the waveform front.
- cylindrical shell chain /
- stress wave /
- transverse inertia correction /
- analytical solution /
- geometric dispersion

HTML全文

由颗粒、柱壳或球壳组成的链结构能够实现应力波的传输，并具备较好的波形操控能力，如波形的放大和衰减^[1-3]，因此可作为声学透镜^[4]、减震材料^[5-6]、声学开关^[7]等，在工程领域有广泛的应用前景。链结构中的波传播过程因离散特性而显得十分复杂，一维链系统的波传播行为可能呈现线性、弱/强非线性的特征^[8]，不同的波传播特征主要受组成链的单个胞元基体材料性质和胞元结构（如胞元形状和胞元尺寸）的影响^[9-12]。其中，由柱壳组成的链结构，因胞元中空且壁厚可控，在操控波形方面表现出优异的性能^{[9-10, 13-14]}。因此，了解柱壳链中应力波的传播规律及弥散机理，对设计新型高效的波形控制器具有重要意义。

建立链结构中相邻胞元之间的接触关系，是研究柱壳链中波传播过程的基础。在一维均匀颗粒链中，Hertz定律常用来表征相邻颗粒之间接触力F和相对位移δ的关系，如F∝δ^1.5^[15]。对于空心壳，接触力和相对位移的关系为F∝δⁿ，其中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以初速度v₀沿X方向撞击一个静止的柱壳链的一端，柱壳链的另一端固定，如图1(a)所示。所有柱壳排布成一列，相邻柱壳相互绑定，链总长度为L，含有N个柱壳。柱壳的壁厚为h，外径为D=L/N，垂直方向厚度为a。柱壳的基体材料为弹性材料，密度为ρ_s，弹性模量为E_s，泊松比为ν_s。在X方向压缩时，柱壳在Y方向将产生不可忽略的变形，即Y方向的泊松比较大，而在面外Z方向的变形量很小，可以忽略，即Z方向的泊松比为0。

图 1 柱壳链简化模型

Figure 1. Simplified model of a cylindrical shell chain

下载: 全尺寸图片幻灯片

在质量块的冲击下，柱壳链中的应力波从近端传向远端。柱壳链的几何非线性将导致复杂的非线性应力波传播问题。为简化分析，假设柱壳链为由均匀的连续介质组成的杆，杆在X方向的长度为L，在Y方向的高度为D，在Z方向的宽度为a，如图1(b)所示。杆的横截面面积为A=aD，等效密度为ρ₀=πρ_s(D−h)h/D²，等效弹性模量为E，在Y方向的等效泊松比为ν，而在Z方向的等效泊松比为0。

1.2 理论分析

1.2.1 计及横向惯性的波传播模型

Rayleigh-Love波动方程可用来描述弹性杆中的应力波传播过程和由横向惯性引起的弥散效应^[22]。Rayleigh-Love波动方程描述了弹性杆中质点沿X方向的位移u(X, t)的演化过程：

$\frac{{{{\text {∂}} ^2}u}}{{{\text {∂}} {t^2}}}=c_0^2\frac{{{{\text {∂}} ^2}u}}{{{\text {∂}} {X^2}}} + {\nu ^2}{\kappa ^2}\frac{{{{\text {∂}} ^4}u}}{{{\text {∂}} {X^2}{\text {∂}} {t^2}}}$

(1)

式中：t为时间， ${c_0}=\sqrt {E/{\rho _0}}$ 为初等理论（不考虑横向惯性效应）下杆中的纵波波速，к为截面对X轴的回转半径（惯性半径）。

在本文中，因为Z方向等效泊松比为0，所以惯性半径к为：

$\kappa =\sqrt {\frac{1}{A}\int_A {{Y^2}{\rm{d}} Y{\rm{d}} Z} } =\sqrt {\frac{1}{D}\int_{ - D/2}^{D/2} {{Y^2}{\rm{d}} Y} } =\frac{{\sqrt 3 }}{6}D$

(2)

初始条件为：

$u(X,0)=0,\quad\quad\quad \frac{{{\text {∂}} u}}{{{\text {∂}} t}}(X \ne 0,0)=0,\quad\quad\quad \frac{{{\text {∂}} u}}{{{\text {∂}} t}}(0,0)={v_0}$

(3)

边界条件为：

$\frac{M}{A}\frac{{{{\text {∂}} ^2}u}}{{{\text {∂}} {t^2}}}(0,t)=E\frac{{{\text {∂}} u}}{{{\text {∂}} X}}(0,t),\quad u(L,t)=0$

(4)

以长度L、波速c₀和密度ρ₀为基本量构造特征参数，定义无量纲量：

$\bar u=\frac{u}{L},\quad {\kern 1pt} \bar X=\frac{X}{L},\quad {\kern 1pt} {\kern 1pt} \bar t=\frac{{{c_0}t}}{L},\quad {\kern 1pt} {\kern 1pt} \bar m=\frac{M}{{A{\rho _0}L}},\quad {\bar v_0}=\frac{{{v_0}}}{{{c_0}}},\quad {\kern 1pt} {\kern 1pt} \bar \kappa =\frac{\kappa }{L}$

(5)

于是，可以得到质量块冲击下弹性杆中应力波传播的无量纲化控制方程：

$\left\{ \begin{aligned} &\frac{{{{\text {∂}} ^2}\bar u}}{{{\text {∂}} {{\bar t}^2}}}=\frac{{{{\text {∂}} ^2}\bar u}}{{{\text {∂}} {{\bar X}^2}}} + \alpha \frac{{{{\text {∂}} ^4}\bar u}}{{{\text {∂}} {{\bar X}^2}{\text {∂}} {{\bar t}^2}}}\quad \quad\quad 0 {\text{＜}} \bar X {\text{＜}} 1,\quad \bar t {\text{＞}} 0 \\ &\bar u(\bar X,0)=0,\quad \frac{{{\text {∂}} \bar u}}{{{\text {∂}} \bar t}}(\bar X \ne 0,0)=0,\quad \frac{{{\text {∂}} \bar u}}{{{\text {∂}} \bar t}}(0,0)={{\bar v}_0},\quad \bar m\frac{{{{\text {∂}} ^2}\bar u}}{{{\text {∂}} {{\bar t}^2}}}(0,\bar t)=\frac{{{\text {∂}} \bar u}}{{{\text {∂}} \bar X}}(0,\bar t),\quad \bar u(1,\bar t)=0 \end{aligned} \right.$

(6)

式中： $\alpha ={\nu ^2}{\bar \kappa ^2}$ 。该控制方程由一个四阶偏微分方程和一些定解条件组成，其中含有2个无量纲自变量 $\bar X$ 和 $\bar t$ 以及3个无量纲参数 $\alpha$ 、 ${\bar v_0}$ 和 $\bar m$ 。尽管计及横向惯性修正的波动方程的形式一致^{[20, 25]}，但对于不同的定解条件，其解大为不同^[24-26]。此外，文献中大都将Rayleigh-Love波动方程应用于均质圆杆中应力波传播行为的研究，而应用到Hopkinson杆实验以对波形进行修正^[25]时，则未见用于链结构中的波传播过程及规律的研究。

1.2.2 波传播模型的解析解

先对 $\bar u(\bar X,\bar t)$ 进行Laplace变换，时域 $\bar t$ 变为像域s，即：

$U(\bar X,s)=\int_0^\infty {\bar u(\bar X,\bar t){{\rm{e}} ^{ - s\bar t}}{\rm{d}} \bar t} \quad\quad\quad {\rm{Re}} \;s {\text{＞}} 0$

(7)

此时，方程(6)可以变作：

$\left\{ \begin{aligned} & {s^2}U=\frac{{{{\rm{d}} ^2}U}}{{{\rm{d}} {{\bar X}^2}}} + \alpha {s^2}\frac{{{{\rm{d}} ^2}U}}{{{\rm{d}} {{\bar X}^2}}}\quad \quad \quad 0 {\text{＜}} \bar X {\text{＜}} 1 \\ & \bar m\left( { - {{\bar v}_0} + {s^2}U(0,s)} \right)=\frac{{{\rm{d}} U}}{{{\rm{d}} \bar X}}(0,s)\quad \quad \quad U(1,s)=0 \end{aligned} \right.$

(8)

式(8)的微分方程是关于 $\bar X$ 的二阶常系数齐次方程，结合定解条件，得到：

$U(\bar X,s)=\frac{{\bar m{{\bar v}_0}\left[ {{{\rm{exp}}\left({ - \dfrac{{s\bar X}}{{\sqrt {1 + \alpha {s^2}} }}}\right)} - {{\rm{exp}}\left({\dfrac{{s\bar X - 2s}}{{\sqrt {1 + \alpha {s^2}} }}}\right)}} \right]\sqrt {1 + \alpha {s^2}} }}{{s\left[ {1 + \bar ms\sqrt {1 + \alpha {s^2}} + \left( {1 - \bar ms\sqrt {1 + \alpha {s^2}} } \right){{\rm{exp}}\left({\dfrac{{ - 2s}}{{\sqrt {1 + \alpha {s^2}} }}}\right)}} \right]}}$

(9)

作Laplace逆变换：

$\bar u(\bar X,\bar t)=\frac{1}{{2{\rm{{\text{π}} }}{\rm{i}} }}\int_{\beta - {\rm{i}} \infty }^{\beta + {\rm{i}} \infty } {U(\bar X,s){{\rm{e}} ^{s\bar t}}{\rm{d}} s}$

(10)

式中：β±i∞为复变量，β为实变数。将式(9)代入式(10)，可得：

$\bar u(\bar X,\bar t)=\frac{{\bar m{{\bar v}_0}}}{{2{\rm{{\text{π}} i}}}}\displaystyle\int_{\beta - {\rm{i}} \infty }^{\beta + {\rm{i}} \infty } {\frac{{\left[ {{{\rm{exp}}\left({s\bar t - \dfrac{{s\bar X}}{{\sqrt {1 + \alpha {s^2}} }}}\right)} - {{\rm{exp}}\left({s\bar t + \dfrac{{s\bar X - 2s}}{{\sqrt {1 + \alpha {s^2}} }}}\right)}} \right]\sqrt {1 + \alpha {s^2}} }}{{s\left[ {1 + \bar ms\sqrt {1 + \alpha {s^2}} + \left( {1 - \bar ms\sqrt {1 + \alpha {s^2}} } \right){{\rm{exp}}\left({\dfrac{{ - 2s}}{{\sqrt {1 + \alpha {s^2}} }}}\right)}} \right]}}{\rm{d}} s}$

(11)

式(11)的积分可以采用留数定理来求解。令被积函数的分母为零，可得式(11)有奇点0和简单奇点 ${s_k} = i{\eta _k}(k = \pm 1,{\rm{ }} \pm 2,{\rm{ }} \pm 3,{\rm{ }} \cdots )$ ，其中实数η_k由此计算：

$\bar m{\eta _k}\sqrt {1 - \alpha \eta _k^2} =\cot \left( {{\eta _k}/\sqrt {1 - \alpha \eta _k^2} } \right)$

(12)

式中： $\eta _k^2 {\text{≤}} 1/\alpha$ 。式(12)存在无穷多个根，且其正根和负根一一对应。值得注意的是，η_k只依赖于α和 $\bar m$ ，而与 ${\bar v_0}$ 无关。由留数定理，有：

$\begin{split} \bar u(\bar X,\bar t)=&\bar m{{\bar v}_0}{\rm{Res}} \left[ {U(\bar X,\bar s){{\rm{e}} ^{s\bar t}},0} \right] + \bar m{{\bar v}_0}\sum\limits_k^{} {{\rm{Res}} \left[ {U(\bar X,\bar s){{\rm{e}} ^{s\bar t}},{s_k}} \right]} \\ =&\bar m{{\bar v}_0}\mathop {\lim}\limits_{s \to 0} U(\bar X,\bar s){{\rm{e}} ^{s\bar t}}s + \bar m{{\bar v}_0}\sum\limits_k^{} {\mathop {\lim}\limits_{s \to {\eta _k}{\rm{i}}} U(\bar X,\bar s){{\rm{e}} ^{s\bar t}}\left( {s - {\rm{i}}{\eta _k}} \right)} \\ =&\bar m{{\bar v}_0}\sum\limits_k^{} {\frac{{{{\left[ {1 - \alpha \eta _k^2} \right)}^2}\left( {{{\rm{exp}}\left({{\rm{i}}{\eta _k}\bar t - \dfrac{{{\rm{i}}{\eta _k}\bar X}}{{\sqrt {1 - \alpha \eta _k^2} }}}\right)} - {{\rm{exp}}\left({{\rm{i}}{\eta _k}\bar t + \dfrac{{{\rm{i}}{\eta _k}\left( {\bar X - 2} \right)}}{{\sqrt {1 - \alpha \eta _k^2} }}}\right)}} \right]\left( {\bar m{\eta _k}\sqrt {1 - \alpha \eta _k^2} + {\rm{i}}} \right)}}{{ - 2{\eta _k}\left\{ {\bar m\left( {1 - \alpha \eta _k^2} \right)\left[ {\eta _k^2\left( {\bar m - 2\alpha } \right) + 1} \right] + 1} \right\}}}} \end{split}$

(13)

可见奇点0对上述结果没有贡献，而所有贡献来自于简单奇点 ${s_k} = i{\eta _k}(k = \pm 1,{\rm{ }} \pm 2,{\rm{ }} \pm 3,{\rm{ }} \cdots )$ 。考虑式(12)中正根与负根的对应关系，式(13)的虚部可以相互抵消。因此，式(13)可以改写为：

$\bar u(\bar X,\bar t)=2\bar m{\bar v_0}\sum\limits_{\begin{subarray}{l} k=1 \\ {\eta _k} {\text{＞}} 0 \end{subarray}} ^\infty {\frac{{{{\left( {1 - \alpha \eta _k^2} \right)}^2}\sin \left( {{\eta _k}\bar t} \right)\sin \left[ {{\eta _k}\left( {1 - \bar X} \right)/\sqrt {1 - \alpha \eta _k^2} } \right]}}{{{\eta _k}\left\{ {1 + \bar m\left( {1 - \alpha \eta _k^2} \right)\left[ {\eta _k^2\left( {\bar m - 2\alpha } \right) + 1} \right]} \right\}\sin \left( {{\eta _k}/\sqrt {1 - \alpha \eta _k^2} } \right)}}}$

(14)

进一步得出应变和无量纲速度：

$\varepsilon = - \frac{{\partial \bar u(\bar X,\bar t)}}{{\partial \bar X}}=2\bar m{\bar v_0}\sum\limits_{\begin{subarray}{l} k=1 \\ {\eta _k} {\text{＞}} 0 \end{subarray}} ^\infty {\frac{{{{\left( {1 - \alpha \eta _k^2} \right)}^{3/2}}\sin \left( {{\eta _k}\bar t} \right)\cos \left[ {{\eta _k}\left( {1 - \bar X} \right)/\sqrt {1 - \alpha \eta _k^2} } \right]}}{{\left\{ {1 + \bar m\left( {1 - \alpha \eta _k^2} \right)\left[ {\eta _k^2\left( {\bar m - 2\alpha } \right) + 1} \right]} \right\}\sin \left( {{\eta _k}/\sqrt {1 - \alpha \eta _k^2} } \right)}}}$

(15)

$\bar v=\frac{{\partial \bar u(\bar X,\bar t)}}{{\partial \bar t}}=2\bar m{\bar v_0}\sum\limits_{\begin{subarray}{l} k=1 \\ {\eta _k} {\text{＞}} 0 \end{subarray}} ^\infty {\frac{{{{\left( {1 - \alpha \eta _k^2} \right)}^2}\cos \left( {{\eta _k}\bar t} \right)\sin \left[ {{\eta _k}\left( {1 - \bar X} \right)/\sqrt {1 - \alpha \eta _k^2} } \right]}}{{\left\{ {1 + \bar m\left( {1 - \alpha \eta _k^2} \right)\left[ {\eta _k^2\left( {\bar m - 2\alpha } \right) + 1} \right]} \right\}\sin \left( {{\eta _k}/\sqrt {1 - \alpha \eta _k^2} } \right)}}}$

(16)

不考虑横向惯性影响（α=0）时，式(14)～(16)可分别简化为：

$\bar u(\bar X,\bar t)=2\bar m{\bar v_0}\sum\limits_{\begin{subarray}{l} k=1 \\ {\eta _k} {\text{＞}} 0 \end{subarray}} ^\infty {\frac{{\sin \left( {{\eta _k}\bar t} \right)\sin \left[ {{\eta _k}\left( {1 - \bar X} \right)} \right]}}{{{\eta _k}\left( {1 + \bar m + {{\bar m}^2}\eta _k^2} \right)\sin {\eta _k}}}}$

(17)

$\varepsilon = - \frac{{\partial \bar u(\bar X,\bar t)}}{{\partial \bar X}}=2\bar m{\bar v_0}\sum\limits_{\begin{subarray}{l} k=1 \\ {\eta _k} {\text{＞}} 0 \end{subarray}} ^\infty {\frac{{\sin \left( {{\eta _k}\bar t} \right)\cos \left[ {{\eta _k}\left( {1 - \bar X} \right)} \right]}}{{\left( {1 + \bar m + {{\bar m}^2}\eta _k^2} \right)\sin {\eta _k}}}}$

(18)

$\bar v=\frac{{\partial \bar u(\bar X,\bar t)}}{{\partial \bar t}}=2\bar m{\bar v_0}\sum\limits_{\begin{subarray}{l} k=1 \\ {\eta _k} {\text{＞}} 0 \end{subarray}} ^\infty {\frac{{\cos \left( {{\eta _k}\bar t} \right)\sin \left[ {{\eta _k}\left( {1 - \bar X} \right)} \right]}}{{\left( {1 + \bar m + {{\bar m}^2}\eta _k^2} \right)\sin {\eta _k}}}}$

(19)

其中，η_k满足的方程(12)可以改写为：

$\bar m{\eta _k}=\cot {\eta _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/m³，弹性模量E_s=1 600 MPa，泊松比ν_s=0.3。柱壳链近端受到刚性质量块（M=1 g）初始速度v₀=1 m/s的冲击，远端固定。不同时刻柱壳链中的应力云图如图2所示。结果表明，应力波从冲击端往支撑端沿着柱壳逐个传播，且应力主要集中在柱壳之间的接触区域。

图 2 不同时刻柱壳链中的应力云图

Figure 2. Von Mises stress distributions in the cylindrical shell chain at different times

下载: 全尺寸图片幻灯片

1.3.2 柱壳链材料参数的确定

等效连续介质模型中的材料参数，可通过对柱壳链进行单轴准静态压缩测试获得。壁厚0.56 mm的柱壳链在单轴准静态恒速压缩(v=0.05 m/s)下的名义应力应变曲线如图3(a)。观察发现，其应力应变曲线呈现轻微软化。采用名义应力应变关系式 ${\sigma _{\rm{N}}}=B\varepsilon _{\rm{N}}^n$ 进行拟合，可得B=2.6 MPa，n=0.93。指数n接近于1，为了分析方便，本文中采用线性关系 ${\sigma _{\rm{N}}}=E{\varepsilon _{\rm{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。

图 3 名义应力应变曲线和等效泊松比

Figure 3. Nominal stress-strain curves and equivalent Poissons ratios

下载: 全尺寸图片幻灯片

2. 结　果

2.1 位移和载荷曲线

质量块（M=1 g）以初速v₀=1 m/s冲击壁厚0.56 mm的柱壳链，其冲击端位移随着加载时间的变化，如图4(a)所示。结果表明，冲击端位移随着时间增大，然后逐渐趋于平缓，理论分析与数值模拟结果基本一致。理论分析与细观有限元模拟所得的支撑端和冲击端的载荷历史曲线如图4(b)所示。观察发现，理论分析能够较好地捕捉载荷曲线的整体特征，但曲线振荡存在一些差异，理论计算的载荷曲线呈现低频振荡，但在数值模拟中，载荷曲线呈现高频振荡，这是由于罚函数接触算法所致。由应力波理论，对于不考虑几何弥散的均质弹性杆，支撑端与冲击端峰值载荷的比为2。而该算例中，由理论分析得到的支撑端和冲击端的峰值载荷分别为2.1、1.3 N，支撑端与冲击端峰值载荷的比小于2。因此，柱壳链可作为一种有效的冲击缓冲器。

图 4 冲击端的位移、冲击端和支撑端的载荷

Figure 4. Displacements at the support end, forces at the impact and support ends

下载: 全尺寸图片幻灯片

2.2 应变和速度分布

对于第2.1节的算例，理论分析和数值模拟给出的应变场和速度场如图5所示。需要说明的是，在数值模拟中，难以计算沿冲击方向连续变化的宏观应变场，因此定义一些局部点上的应变以检验理论结果。在过柱壳中心的YZ横截面上，取位于其表面的所有点作为该柱壳的参考点，提取这些点在冲击过程中的平均位移和速度，进一步通过相邻两个柱壳参考点的相对位移ΔX计算应变，即 $\varepsilon =\Delta X/D$ 。为方便比较，理论上也可以计算与数值模拟中对应参考点处的速度和应变。图5(a)为不同时刻的应变分布，结果表明，随着波往远端传播，应变峰值略有减小，且随着时间的推移，波后振荡的区域越来越大。理论计算与细观有限元模拟所得结果整体趋势一致，波后的波形振荡和波形前沿宽度也大致相同。不同时刻的速度分布如图5(b)所示。从速度分布可以发现，峰值速度随着波的传播逐渐降低，理论计算的粒子速度与细观有限元结果较好吻合。

图 5 不同时刻杆中的应变和速度分布

Figure 5. Strain and velocity distributions in rod at different times

下载: 全尺寸图片幻灯片

当质量块（M=1 g）冲击速度为1 m/s时，柱壳链中的应变较小，增大冲击速度可以检验应变较大时理论模型的适用性。图6为冲击速度为5、10 m/s时柱壳链中的应变分布。结果表明，理论分析与数值模拟的结果较好吻合，冲击速度越大，峰值应变越大，表明理论能够很好地预测峰值较大应变(约0.1)时柱壳链中不同时刻的应变分布。

图 6 不同冲击速度下杆中的应变分布

Figure 6. Strain distributions in rod under different impact velocities

下载: 全尺寸图片幻灯片

进一步比较不同壁厚柱壳链的理论分析与数值模拟的结果。通过理论分析与细观有限元模拟计算，壁厚0.33、1.07 mm的柱壳链在受到质量块冲击(M=1 g, v₀=1 m/s)时，链中的应变分布如图7所示。结果表明，理论与数值模拟很好吻合，理论能够很好地预测不同壁厚柱壳链中不同时刻的应变分布。

图 7 不同壁厚时杆中的应变分布

Figure 7. Strain distributions in rod with different wall thicknesses

下载: 全尺寸图片幻灯片

3. 讨　论

3.1 泊松比和无量纲惯性半径对波传播过程的影响

在理论分析中，可以发现泊松比ν和无量纲惯性半径 $\bar \kappa$ 总是以组合的形式出现，因此可以直接讨论无量纲参数 $\alpha ={\nu ^2}{\bar \kappa ^2}$ 的影响。但是，为了更直观地考察泊松比ν和无量纲惯性半径 $\bar \kappa$ 的影响，这里分别改变这两个参数。图8(a)为不同泊松比某个时刻杆中的应变分布，其中无量纲惯性半径为0.04，泊松比为0、0.25、0.50、0.75和1.00。结果表明，泊松比为0时，应变分布曲线光滑无振荡，应变波形前沿宽度很窄。当泊松比不为0时，应变分布存在明显的振荡现象，且随着泊松比的增大，振荡幅值增大，应变波形前沿宽度增宽。此外，应变峰值也随着泊松比的增大而逐渐降低。以应变峰值处作为波阵面位置，泊松比越大，波速越小。类似地，固定泊松比为某个值，如ν=0.5，讨论不同惯性半径 $\bar \kappa$ 对变波形的影响。 $\bar \kappa$ 受柱壳半径的影响，图8(b)中 $\bar \kappa$ 由小到大分别为0.02、0.04和0.06。与泊松比的影响类似，惯性半径越大，应变峰值越小，波速越小，振荡幅值越大。综上可见，一维杆中横向惯性效应的影响越大，应变分布振荡越剧烈，并且波形前沿宽度增宽、应变峰值减小。

图 8 泊松比和惯性半径对应变分布的影响

Figure 8. Influences of Poisson’s ratio and inertia radius on strain distributions

下载: 全尺寸图片幻灯片

3.2 无量纲速度和质量对波传播过程的影响

无量纲速度 ${\bar v_0}$ 和质量 $\bar m$ 对应变波形的影响如图9所示。可以看出，随着无量纲速度的增大，应变波形的振动波长不变，应变最大位置处对应的波速不变，波形前沿宽度不变，变化的仅是应变幅值。其实，由式(14)可以发现，无量纲速度仅出现在影响应变幅值的项。不同于无量纲速度，无量纲质量不仅影响应变幅值，也影响波形的波长。无量纲质量越大，应变幅值越大，但由图9(b)可见，质量在一定范围内变化时，对波形前沿宽度、应变最大位置处对应波速的影响均较小。

图 9 无量纲速度和无量纲质量对应变分布的影响

Figure 9. Influences of dimensionless velocity and mass on strain distributions

下载: 全尺寸图片幻灯片

3.3 弥散效应的分析

进一步分析Rayleigh-Love波动方程的弥散关系。在计及横向惯性的弹性杆中，纵波波速不再恒定^{[20, 25]}。将行波解^[20]：

$u(X,t)={u_0}\exp [{\rm{i}}\xi \left( {X - ct} \right)]$

(21)

代入式(1)，可以得到^{[20, 25]}：

$\frac{c}{{{c_0}}}=\frac{1}{{\sqrt {1 + {\nu ^2}{\kappa ^2}{\xi ^2}} }} \approx 1 - \frac{{2{{\rm{{\text{π}} }}^2}{\nu ^2}{\kappa ^2}}}{{{\lambda ^2}}}$

(22)

式中：ξ=2π/λ，λ为波长。在圆杆的弥散分析^{[20, 25]}中，其惯性半径 $\kappa ={{\sqrt 2 D} / 4}$ 。在本文中，长方形截面杆在垂直于冲击方向的面内各向异性，即面内不同方向上具有不同的泊松比，其惯性半径由式(2)为 $\kappa ={{\sqrt 3 D} / 6}$ 。可见，波速与泊松比ν和D/λ相关，随着ν或D/λ的增大，波速越来越小，这与上述参数分析的结论一致。

4. 结　论

建立柱壳链结构的等效连续介质模型，对冲击载荷作用下柱壳链中的应力波传播过程及其几何弥散特性进行了分析，并用细观有限元模型对理论模型进行验证。

（1）基于考虑横向惯性修正的Rayleigh-Love波动方程，分析了几何弥散对应力波传播的影响，给出了柱壳链在质量块初速冲击下的无量纲控制方程。基于Laplace变换及逆变换，获得了该控制方程的解析解。

（2）解析解的载荷和位移曲线、应变和速度分布均与细观有限元结果较好吻合，且适用于不同壁厚的柱壳链。

（3）讨论了泊松比、惯性半径、冲击质量和速度对冲击载荷作用下柱壳链中波形弥散的影响。结果表明，应变峰值、波形振荡幅度、波形前沿宽度和应变最大位置处对应的波速主要与泊松比和惯性半径相关，泊松比和惯性半径越大，应变峰值越小，波形振荡幅度越大，波形前沿宽度越宽，应变最大位置处对应的波速越小。

本文中建立的冲击载荷作用下柱壳链中的弹性波传播简化模型及解析解，可为揭示柱壳链的弥散机理和波形控制器的优化设计提供参考。

图 1 柱壳链简化模型

Figure 1. Simplified model of a cylindrical shell chain

下载: 全尺寸图片幻灯片

图 2 不同时刻柱壳链中的应力云图

Figure 2. Von Mises stress distributions in the cylindrical shell chain at different times

下载: 全尺寸图片幻灯片

图 3 名义应力应变曲线和等效泊松比

Figure 3. Nominal stress-strain curves and equivalent Poissons ratios

下载: 全尺寸图片幻灯片

图 4 冲击端的位移、冲击端和支撑端的载荷

Figure 4. Displacements at the support end, forces at the impact and support ends

下载: 全尺寸图片幻灯片

图 5 不同时刻杆中的应变和速度分布

Figure 5. Strain and velocity distributions in rod at different times

下载: 全尺寸图片幻灯片

图 6 不同冲击速度下杆中的应变分布

Figure 6. Strain distributions in rod under different impact velocities

下载: 全尺寸图片幻灯片

图 7 不同壁厚时杆中的应变分布

Figure 7. Strain distributions in rod with different wall thicknesses

下载: 全尺寸图片幻灯片

图 8 泊松比和惯性半径对应变分布的影响

Figure 8. Influences of Poisson’s ratio and inertia radius on strain distributions

下载: 全尺寸图片幻灯片

图 9 无量纲速度和无量纲质量对应变分布的影响

Figure 9. Influences of dimensionless velocity and mass on strain distributions

下载: 全尺寸图片幻灯片

参考文献(26)

[1]	JIANG S, SHEN L M, GUILLARD F, et al. Energy dissipation from two-glass-bead chains under impact [J]. International Journal of Impact Engineering, 2018, 114: 160–168. DOI: 10.1016/j.ijimpeng.2018.01.002.
[2]	KIM E, YANG J, HWANG H, et al. Impact and blast mitigation using locally resonant woodpile metamaterials [J]. International Journal of Impact Engineering, 2017, 101: 24–31. DOI: 10.1016/j.ijimpeng.2016.09.006.
[3]	PARK C M, PARK J J, LEE S H, et al. Amplification of acoustic evanescent waves using metamaterial slabs [J]. Physical Review Letters, 2011, 107(19): 194301. DOI: 10.1103/physrevlett.107.194301.
[4]	DONAHUE C M, ANZEL P W J, BONANOMI L, et al. Experimental realization of a nonlinear acoustic lens with a tunable focus [J]. Applied Physics Letters, 2014, 104(1): 014103. DOI: 10.1063/1.4857635.
[5]	TAN K T, HUANG H H, SUN C T. Blast-wave impact mitigation using negative effective mass density concept of elastic metamaterials [J]. International Journal of Impact Engineering, 2014, 64: 20–29. DOI: 10.1016/j.ijimpeng.2013.09.003.
[6]	DARAIO C, NESTERENKO V F, HERBOLD E B, et al. Energy trapping and shock disintegration in a composite granular medium [J]. Physical Review Letters, 2006, 96(5): 058002. DOI: 10.1103/physrevlett.96.058002.
[7]	LI F, ANZEL P, YANG J, et al. Granular acoustic switches and logic elements [J]. Nature Communications, 2014, 5(1): 5311. DOI: 10.1038/ncomms6311.
[8]	NESTERENKO V F. Propagation of nonlinear compression pulses in granular media [J]. Journal of Applied Mechanics and Technical Physics, 1983, 24(5): 733–743. DOI: 10.1007/bf00905892.
[9]	KIM H, KIM E, CHONG C, et al. Demonstration of dispersive rarefaction shocks in hollow elliptical cylinder chains [J]. Physical Review Letters, 2018, 120(19): 194101. DOI: 10.1103/physrevlett.120.194101.
[10]	NGO D, GRIFFITHS S, KHATRI D, et al. Highly nonlinear solitary waves in chains of hollow spherical particles [J]. Granular Matter, 2013, 15(2): 149–155. DOI: 10.1007/s10035-012-0377-5.
[11]	ON T, WANG E H, LAMBROS J. Plastic waves in one-dimensional heterogeneous granular chains under impact loading: Single intruders and dimer chains [J]. International Journal of Solids and Structures, 2015, 62: 81–90. DOI: 10.1016/j.ijsolstr.2015.02.006.
[12]	GANESH R, GONELLA S. Nonlinear waves in lattice materials: adaptively augmented directivity and functionality enhancement by modal mixing [J]. Journal of the Mechanics and Physics of Solids, 2017, 99: 272–288. DOI: 10.1016/j.jmps.2016.11.001.
[13]	YIN S, CHEN D H, XU J. Novel propagation behavior of impact stress wave in one-dimensional hollow spherical structures [J]. International Journal of Impact Engineering, 2019, 134: 103368. DOI: 10.1016/j.ijimpeng.2019.103368.
[14]	KIM H, KIM E, YANG J. Nonlinear wave propagation in 3D-printed graded lattices of hollow elliptical cylinders [J]. Journal of the Mechanics and Physics of Solids, 2019, 125: 774–784. DOI: 10.1016/j.jmps.2019.02.001.
[15]	JOHNSON K L. Contact mechanics [M]. Cambridge: Cambridge University Press, 1985.
[16]	PAUCHARD L, RICA S. Contact and compression of elastic spherical shells: the physics of a“ping-pong”ball [J]. Philosophical Magazine B, 1998, 78(2): 225–233. DOI: 10.1080/13642819808202945.
[17]	KIM E, YANG J. Review: wave propagation in granular metamaterials [J]. Functional Composites and Structures, 2019, 1(1): 012002. DOI: 10.1088/2631-6331/ab0c7e.
[18]	DENG B L, MO C Y, TOURNAT V, et al. Focusing and mode separation of elastic vector solitons in a 2D soft mechanical metamaterial [J]. Physical Review Letters, 2019, 123(2): 024101. DOI: 10.1103/physrevlett.123.024101.
[19]	HERBOLD E B, NESTERENKO V F. Propagation of rarefaction pulses in discrete materials with strain-softening behavior [J]. Physical Review Letters, 2013, 110(14): 144101. DOI: 10.1103/physrevlett.110.144101.
[20]	WANG L L. Foundations of stress waves [M]. New York: Elsevier, 2007.
[21]	CHREE C. The equations of an isotropic elastic solid in polar and cylindrical co-ordinates their solution and application [J]. Transactions of the Cambridge Philosophical Society, 1889, 14: 250.
[22]	LOVE A E H. A treatise on the mathematical theory of elasticity [M]. Cambridge: Cambridge University Press, 2013.
[23]	BRIZARD D, JACQUELIN E, RONEL S. Polynomial mode approximation for longitudinal wave dispersion in circular rods [J]. Journal of Sound and Vibration, 2019, 439: 388–397. DOI: 10.1016/j.jsv.2018.09.062.
[24]	DAVIES R M. A critical study of the Hopkinson pressure bar [J]. Philosophical Transactions of the Royal Society A: Mathematical Physical Sciences, 1948, 240(821): 375–457. DOI: 10.1098/rsta.1948.0001.
[25]	杨洪升, 李玉龙, 周风华. 梯形应力脉冲在弹性杆中的传播过程和几何弥散 [J]. 力学学报, 2019, 51(6): 1820–1829. DOI: 10.6052/0459-1879-19-183. YANG H S, LI Y L, ZHOU F H. The propagation process and the geometric dispersion of a trapezoidal stress pulse in an elastic rod [J]. Chinese Journal of Theoretical and Applied Mechanics, 2019, 51(6): 1820–1829. DOI: 10.6052/0459-1879-19-183.
[26]	NAITOH M, DAIMARUYA M. The influence of a rise time of longitudinal impact on the propagation of elastic waves in a bar [J]. Bulletin of JSME, 1985, 28(235): 20–25. DOI: 10.1299/jsme1958.28.20.

施引文献

期刊类型引用(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 . 本站查看