压电粘弹性纳米梁的振动

The Vibration of Piezoelectric Viscoelastic Nanobeam

摘要: 压电材料因其具有独特的力电耦合特性而受到研究者的关注。当材料发生变形时可能表现出既有粘性又有弹性的特性，分数阶微积分理论可以很好地刻画材料所具备的这种特性。本文考虑了非局部弹性理论，利用Hamilton原理推导了粘弹性压电纳米材料的分数阶控制方程，采用Galerkin方法和预估校正法求解控制方程，讨论了外加电压、非局部参数、分数阶等参数对粘弹性压电纳米梁的影响。结果表明：分数阶对粘弹性材料产生影响；系统阻尼随着分数阶的增加而增加。此外，非局部因子和压电场对其非线性振动行为有显著影响。

Abstract: Piezoelectric materials have attracted the attention of researchers because of their unique electromechanical coupling characteristics. When the material is deformed, it may show both viscous and elastic characteristics, and the fractional calculus theory can well describe the characteristics of the material. In this paper, considering the theory of non-local elasticity, the fractional governing equations of viscoelastic piezoelectric nanomaterials are derived from Hamilton’s principle. The Galerkin method and the predictive correction method are used to solve the governing equations. The effects of applied voltage, non-local parameters, fractional order and other parameters on the viscoelastic piezoelectric nanobeams are discussed. The results show that the fractional order affects viscoelastic materials. The damping of the system increases as the fractional order increases. In addition, non-local factors and piezoelectric fields have a significant effect on the nonlinear vibration behavior.

1. 引言

随着科技的发展，对于材料的要求也越来越高，压电材料作为一种新型的智能材料，广泛用于医学、声学、航空航天等诸多领域。分数阶微积分理论已经广泛被用于处理粘弹性问题，由于分数阶具有更好的时间记忆性，且可利用较少的参数更精确地模拟粘弹性材料的本构模型。因此，分数阶微积分理论与粘弹性理论结合起来分析。近年来，研究者将粘弹性材料与尺寸相关的非局部理论和表面效应理论等结合起来考虑粘弹性结构的力学[1]-[4]。Lei等[5]推导了考虑非局部弹性理论的分数阶粘弹性欧拉–伯努利纳米梁的自由振动方程，研究了非局部参数、粘弹性常数和外部阻尼参数对其固有频率的影响。他们还研究了非局部粘弹性Timoshenko梁的相同问题[6]。Qiu等[7]分析了初始位移、分数阶、粘弹性系数、长厚比、阻尼系数和力幅值对欧拉–伯努利纳米梁线性和非线性振动时间响应的影响。

基于已有的分数阶微积分理论，研究者们将非局部弹性理论、表面弹性理论、非局部应变梯度理论等扩展到压电材料中，分析在压电智能材料中所表现出的力学现象。Ansari等人[8]利用非局部弹性理论，通过数值方法研究了非局部参数、温度变化和外加电压对不同端部和不同边界条件下的压电纳米梁尺寸相关振动特性的影响。Gheshlaghi等[9]推导了简支条件下压电纳米线的固有频率和屈曲电压表达式，分析了压电纳米线的振动特性，强调了表面效应和小尺度效应的影响。Wang [10]研究了在表面效应下，表面能密度和表面松弛参数对压电纳米线自由横向振动和屈曲的影响。Ghorbanpour-Arani等[11]将非局部粘弹性理论与Euler-Bernoulli梁模型相结合，分析了双粘弹性压电纳米梁系统的自由振动和强迫振动，并考虑了非局部参数、粘–巴斯特尔纳克常数和电压Maxwell系数对双粘弹性压电纳米梁振动特性的影响。Ke和Wang [12]讨论了非局部理论对压电纳米梁非线性振动特性的影响，分析了非局部参数、温度变化和外加电压对压电纳米梁量纲相关非线性振动特性的影响。Zenkour等人[13]基于非局部板理论讨论了压电型Kelvin-Voigt粘弹性纳米片在粘性帕斯捷尔纳克介质中的振动特性。分析了不同边界条件下粘弹性压电纳米板的长径比、边厚比、温度变化和湿度浓度对粘弹性压电纳米板特征频率的影响。Alireza等人[14]建立了基于Timoshenko非局部应变梯度理论(NGST)的纳米梁的耦合连续粘弹性模型，使用有限分析方法进行研究，讨论了Timoshenko纳米梁中粘度是负功模拟能量耗散的主要原因。Mohamed等人[15]利用Bernstein多项式分析，研究了具有粘弹性边界条件的分数阶粘弹性梁的非线性受迫振动。探讨了Caputo分数阶导数阶数、地基参数、初曲率幅值、粘弹性参数和激振幅值对粘弹性梁的非线性受迫振动的影响。Shao等[16]研究了热载荷作用下张力时变的粘弹性运动薄膜的参数振动，利用Kelvin粘弹性本构关系和Hamilton原理推导出了非线性振动方程，采用多尺度方法和Routh-Hurwitz准则求解了薄膜系统的失稳响应。Makkad等人[17]研究了非局部分数阶粘弹性纳米梁在斜坡型加热下的弯曲分析，采用了积分变换方法导出了温度、弯矩、挠度和热应力的封闭解。考察了松弛时间、斜坡时间参数和分数阶阶数在各个方面的影响。

2. 基础知识

2.1. 非局部压电弹性理论

压电材料是机械能和电能之间的转换。压电材料的本构方程为[18]

$\begin{array}{l}{\sigma }_{ij}-{\left(e_{0}a\right)}^2{\nabla }^2{\sigma }_{ij}=C_{ijkl}{\epsilon }_{kl}-e_{kij}{E}_k,\\ D_i-{\left(e_{0}a\right)}^2{\nabla }^2{D}_i=e_{ikl}{\epsilon }_{kl}+{\zeta }_{ik}{E}_k,\end{array}$ (1)

其中，${\sigma }_{ij}$ ，${\epsilon }_{kl}$ ，$E_i$ ，$D_i$ 分别是弹性应力、弹性应变、电场强度、电位移，$c_{ij}$ ，$e_{kij}$ ，${\zeta }_{ik}$ 为弹性常数、压电常数和介电常数。

2.2. 分数阶微积分

当$0<\alpha <1$ 时，Riemann-Liouville型分数阶微分可定义为[19]

$$ $\begin{array}{c}{}_a{}^{RL}D{}_t^{\alpha }f\left(t\right)=\frac{1}{\Gamma \left(1-\alpha \right)}\frac{\text{d}}{\text{d}t}{\displaystyle {\int }_a^t\frac{f\left(\tau \right)}{{\left(t-\tau \right)}^{\alpha }}\text{d}\tau }\\ =\frac{f\left(0\right)}{\Gamma \left(1-\alpha \right)t^{\alpha }}+\frac{1}{\Gamma \left(1-\alpha \right)}{\displaystyle {\int }_0^t\frac{f^{\prime }\left(\tau \right)}{{\left(t-\tau \right)}^{\alpha }}\text{d}\tau }\end{array}$ (2)

其中，$\Gamma \left(z\right)={\displaystyle {\int }_0^{\infty }{\text{e}}^{-t}{t}^{z-1}\text{d}t}$ 为Gamma函数。当$0<\alpha <1$ 时，Caputo型分数阶导数为

${}_a{}^{C}D{}_t^{\alpha }f\left(t\right)=\frac{1}{\Gamma \left(1-\alpha \right)}{\displaystyle {\int }_a^t\frac{f^{\prime }\left(\tau \right)}{{\left(t-\tau \right)}^{\alpha }}\text{d}\tau }$ (3)

因此，对于函数$f\left(t\right)$ ，Riemann-Liouville型和Caputo型之间的关系为

${}_a{}^{RL}D{}_t^{\alpha }f\left(t\right)={}_a{}^{C}D{}_t^{\alpha }f\left(t\right)+\frac{f\left(0\right)}{\Gamma \left(1-\alpha \right)t^{\alpha }}$ (4)

2.3. 分数阶Kelvin-Voigt模型

如图1，分数阶Kelvin-Voigt模型[19]用于描述粘弹性固体，其模型由一个弹簧元件和Abel粘壶并联组成。

由于弹簧元件与Abel粘壶是并联的，因此有

$\{\begin{array}{c}{\sigma }_e=E{\epsilon }_e,{\sigma }_v=\eta \frac{{\text{d}}^{\alpha }{\epsilon }_v}{\text{d}{t}^{\alpha }}\\ \epsilon ={\epsilon }_e={\epsilon }_v,\sigma ={\sigma }_e+{\sigma }_v\end{array}$ (5)

式中，下标$e$ 表示弹性元件，下标$v$ 表示粘性元件。分数阶Kelvin-Voigt模型的本构方程为

$\begin{array}{c}\sigma =E\epsilon \left(x,t\right)+E_{t_{\alpha }}{D}_t^{\alpha }\epsilon \left(x,t\right)\\ =E\left(1+\overline{g}\frac{{\partial }^{\alpha }}{\partial t^{\alpha }}\right)\epsilon \left(x,t\right)\end{array}$ (6)

Figure 1. Fractional Kelvin-Voigt model

图1. 分数阶Kelvin-Voigt模型

其中，$E$ 为材料的杨氏模量，$E_{t_{\alpha }}$ 为阻尼系数，$\overline{g}$ 为粘弹性系数且$\overline{g}=E_{t_{\alpha }}/E$ ，$\epsilon \left(x,t\right)$ 为轴向应变，$D_t^{\alpha }$ 为分数阶导数。

3. 模型的建立和控制方程

为了研究粘弹性压电纳米梁的力学行为，如图2展示了电压$\varphi \left(x,z,t\right)$ 下的长度为$L$ 、宽度$b$ 和厚度$h$ 的压电纳米梁模型。

Figure 2. Schematic diagram of viscoelastic piezoelectric nanobeam

图2. 粘弹性压电纳米梁示意图

3.1. 控制方程得推导

基于Euler-Bernoulli梁理论，由粘弹性纳米梁的位移场及根据von-Kármán的几何非线性假设，可以得到应变与位移之间的关系[5]

${\epsilon }_{xx}=\frac{\partial u_x}{\partial x}+\frac{1}{2}{\left(\frac{\partial u_z}{\partial x}\right)}^2=\frac{\partial U}{\partial x}-z\frac{{\partial }^{2}W}{\partial x^2}+\frac{1}{2}{\left(\frac{\partial W}{\partial x}\right)}^2$ (7)

根据分数阶Kelvin-Voigt模型，及压电材料的非局部弹性理论就可以得到

$\begin{array}{l}{\sigma }_{xx}-{\left(e_{0}a\right)}^2{\nabla }^2{\sigma }_{xx}=\left(1+\overline{g}\frac{{\partial }^{\alpha }}{\partial t^{\alpha }}\right)c_{11}{\epsilon }_{xx}-e_{31}{E}_z\hfill \\ D_x-{\left(e_{0}a\right)}^2{\nabla }^2{D}_x={\zeta }_{11}{E}_x\hfill \\ D_z-{\left(e_{0}a\right)}^2{\nabla }^2{D}_z=e_{31}{\epsilon }_{xx}+{\zeta }_{33}{E}_z\hfill \end{array}$ (8)

3.2. Hamilton原理

关于Hamilton原理，它是以变分为基础的建模方法，通过计算出势能和动能，就可以推导出分数阶粘弹性压电纳米梁的控制方程

${\displaystyle {\int }_0^t\left(\delta {\Pi }_k-\delta {\Pi }_s\right)\text{d}t}=0$ (9)

应变能的变分为

$\begin{array}{c}\delta {\Pi }_s={\displaystyle {\int }_0^L\left[N_x^{\ast }\left(\frac{\partial }{\partial x}\delta U+\frac{\partial W}{\partial x}\frac{\partial }{\partial x}\delta W\right)-M_x^{\ast }\frac{{\partial }^2}{\partial x^2}\delta W\right]\text{d}x}\\ -{\displaystyle {\int }_0^L{\displaystyle {\int }_{-\frac{h}{2}}^{\frac{h}{2}}\left[D_z\beta \sin \left(\beta z\right)\delta \phi -D_x\cos \left(\beta z\right)\frac{\partial }{\partial x}\delta \phi \right]\text{d}z\text{d}x}}\end{array}$ (10)

其中，并且法向合力和弯矩为

$N_x^{\ast }={\displaystyle {\int }_{-\frac{h}{2}}^{\frac{h}{2}}{\sigma }_{xx}\text{d}z},M_x^{\ast }={\displaystyle {\int }_{-\frac{h}{2}}^{\frac{h}{2}}z\sigma {}_{xx}\text{d}z}$ (11)

动能的变分为

$\delta {\Pi }_k={\displaystyle {\int }_0^{L}I\left(\frac{\partial U}{\partial t}\frac{\partial }{\partial t}\delta U+\frac{\partial W}{\partial t}\frac{\partial }{\partial t}\delta W\right)\text{d}x}$ (12)

其中，$I=\rho h$ 。

因此，可以根据变分原理得到压电纳米梁的控制方程

$\begin{array}{l}\frac{\partial N_x^{*}}{\partial x}=I\frac{{\partial }^{2}U}{\partial t^2}\\ \frac{{\partial }^2{M}_x^{*}}{\partial x^2}+\frac{\partial }{\partial x}\left(N_x^{*}\frac{\partial W}{\partial x}\right)=I\frac{{\partial }^{2}W}{\partial x^2}\\ {\displaystyle {\int }_{-\frac{h}{2}}^{\frac{h}{2}}\left[D_z\beta \sin \left(\beta z\right)+\frac{\partial D_x}{\partial x}\cos \left(\beta z\right)\right]\text{d}z}=0\end{array}$ (13)

3.3. 控制方程

通过把(10)和(12)式代入(9)式，应用Hamilton原理，得到控制方程，对该方程进行无量纲化处理

$\begin{array}{l}w=\frac{W}{L},\xi =\frac{x}{L},\mu =\frac{e_{0}a}{L},O=L^2\sqrt{\frac{I_1}{D_{11}}},\left\{\tau ,g\right\}=\left\{\frac{t}{O},\frac{\overline{g}}{O^{\alpha }}\right\},\\ \beta =\frac{L^2{A}_{11}}{D_{11}},m=\frac{X_{11}}{X_{33}{L}^2},\gamma =\frac{F_{31}^2}{X_{33}{D}_{11}},p=\frac{2e_{31}{V}_0{L}^2}{D_{11}}\end{array}$(14)

从而得到无量纲下的分数阶粘弹性压电纳米梁的控制方程

$\begin{array}{l}\left[m+\frac{{\mu }^{2}m\beta }{2}{{\displaystyle {\int }_0^1\left(\frac{\partial w}{\partial \xi }\right)}}^2\text{d}\xi \right]\left(1+g\frac{{\partial }^{\alpha }}{\partial {\tau }^{\alpha }}\right)\frac{{\partial }^{6}w}{\partial {\xi }^6}+\frac{\beta }{2}{{\displaystyle {\int }_0^1\left(\frac{\partial w}{\partial \xi }\right)}}^2\text{d}\xi \left(1+g\frac{{\partial }^{\alpha }}{\partial {\tau }^{\alpha }}\right)\frac{{\partial }^{2}w}{\partial {\xi }^2}\\ -\left\{\left[1+\left({\mu }^2+m\right)\frac{\beta }{2}{{\displaystyle {\int }_0^1\left(\frac{\partial w}{\partial \xi }\right)}}^2\text{d}\xi \right]\left(1+g\frac{{\partial }^{\alpha }}{\partial {\tau }^{\alpha }}\right)+\gamma +pm\left(1-{\mu }^2\frac{{\partial }^2}{\partial {\xi }^2}\right)\right\}\frac{{\partial }^{4}w}{\partial {\xi }^4}+p\left(1-{\mu }^2\frac{{\partial }^2}{\partial {\xi }^2}\right)\frac{{\partial }^{2}w}{\partial {\xi }^2}\\ =\frac{{\partial }^2}{\partial {\tau }^2}\left[\left(w-\left({\mu }^2+m\right)\frac{{\partial }^{2}w}{\partial {\xi }^2}\right)+{\mu }^2\frac{{\partial }^{4}w}{\partial {\xi }^4}\right]\end{array}$ (15)

3.4. 问题求解

本篇文章采取Galerkin方法进行求解问题，它是指采用微分方程对应的弱形式，其原理为通过选取有限多项式函数(又称基函数)，将他们叠加，再要求结果在求解域内及边界上的加权积分(权函数为试函数本身)满足原方程，便可以得到一组易于求解的线性代数方程，且自然边界条件能够自动满足。取梁的瞬时模态函数作为试函数，对振动的偏微分方程进行离散求解，简支梁(SS)的阵型函数为

$p_n\left(x\right)=A_n\sin \frac{n\pi x}{L}\left(n=1,2,\cdots ,\infty \right)$ (16)

应用Galerkin方法将积分–分数阶微分方程转化为常微分方程，以状态空间的形式表示，再通过预估校正法进行数值求解。

对于简支压电纳米梁，式(15)的解可以表示为

$w\left(\xi ,t\right)={\displaystyle \sum _{n=1}^N{\phi }_n\left(\tau \right)p_n\left(\xi \right)}$ (17)

其中，${\phi }_n\left(\tau \right)$ 表示随时间变化的函数，$P_n\left(\xi \right)=\sin \left(n\pi \xi \right)$ 为纳米梁所对应的模态函数。

把式(17)代入(15)式中，应用Galerkin过程，即有

$\begin{array}{l}m{\Gamma }^{\left(6\right)}\left(\frac{{\mu }^2\beta }{2}{\rm X}\Lambda +\Phi \right)-\frac{\beta }{2}\left[\left({\mu }^2+m\right){\Gamma }^{\left(4\right)}-{\Gamma }^{\left(2\right)}\right]{\rm X}\Lambda -{\Gamma }^{\left(4\right)}\Phi -\gamma {\Gamma }^{\left(4\right)}\left({\phi }_j\right)\\ +p\left[m{\mu }^2{\Gamma }^{\left(6\right)}+{\Gamma }^{\left(2\right)}-\left({\mu }^2+m\right){\Gamma }^{\left(4\right)}\right]\left({\phi }_j\right)\\ =\left\{{\Gamma }^{\left(0\right)}-{\mu }^2{\Gamma }^{\left(2\right)}-m\left[{\Gamma }^{\left(2\right)}+{\mu }^2{\Gamma }^{\left(4\right)}\right]\right\}\left({\ddot{\phi }}_j\right),\left(i=1,\cdots ,n\right)\end{array}$(18)

其中

$\begin{array}{l}{\Gamma }^{\left(\delta \right)}={\displaystyle \sum _{j=1}^N{\displaystyle {\int }_0^1{p}_i{p}_j^{\left(\delta \right)}\text{d}\xi }},\left(\delta =0,2,4,6\right)\\ \Lambda =\left[{\phi }_j{\phi }_k{\phi }_l+g{\phi }_j{\left({\phi }_k{\phi }_l\right)}^{\left(\alpha \right)}\right]\\ {\rm X}={\displaystyle \sum _{k=1}^N{\displaystyle \sum _{l=1}^N{\displaystyle {\int }_0^1{p^{\prime }}_k{p^{\prime }}_l\text{d}\xi },}}\Phi =\left({\phi }_j+g{\phi }_j^{\left(\alpha \right)}\right)\end{array}$ (19)

从而利用预估矫正法求解非线性方程。

4. 数值模拟

压电纳米梁是由PZT-4 [12]制成，PZT-4是一种经过改性的压电陶瓷材料，它能够快速在机械能与电能之间转换，它具有低机械损耗和介电损耗、较大的交流退极化场，以及较高的介电常数、机电耦合系数和压电常数，其材料性能如下所示：$c_{11}=132\text{Gpa}$ ，$e_{31}=-4.1C/{\text{m}}^2$ ，${\zeta }_{11}=5.841\times {10}^{-9}C/\text{V}\cdot \text{m}$ ，${\zeta }_{33}=7.124\times {10}^{-9}N/{\text{m}}^2\cdot K$ ，${\lambda }_1=4.738\times {10}^{5}N/{\text{m}}^2\cdot K$ ，$\rho =7500\text{kg}/{\text{m}}^3$ ，$L=12\text{nm}$ ，$h=2\text{nm}$ 。下面是基于此材料特性，对该模型的讨论。

将本文的结果退化为参考文献[12]中进行比较，如图3所示，展示了在考虑非局部理论时，本文研究的分数阶粘弹性压电纳米梁的非线性振动频率高于Ke等人，振幅比较低，这主要是由于分数阶粘弹性对系统的抑制作用。Ke等人仅考虑了非局部理论，在弹性情况下的纳米梁的非线性振动，本文不考虑温度对粘弹性纳米梁的影响，引入了分数阶粘弹性，从而使梁既具备粘性特性，也具备弹性特性。观察图中的两条曲线，两者之间存在差距，因此，在分析考虑问题的相关影响时，要综合进行考量。下面展示了几类参数对分数阶粘弹性压电纳米梁非线性振动时间相应的影响。

如图4(a)所示，从图中可以看出，分数阶的阶数影响着粘弹性压电纳米梁的非线性振动，当分数阶的阶数增加时，系统的非线性振动频率和振幅在减小。分数阶的阶数从0.2增加到0.9，系统引入分数阶

Figure 3. Comparison of results considering nonlocal theory with those of Ke et al.

图3. 考虑非局部理论的结果与Ke等人结果的比较

粘弹性，阻尼增强。结果表明，随着分数阶的阶数的增加，系统的阻尼增大，刚度减小，导致梁模型振动减弱。考虑非局部理论时，研究整个系统的非线性时间响应，如图4(b)所示，从图中可以明显看出，随着非局部参数的不断增加，纳米梁的非线性频率和振幅呈现出逐渐减小的趋势。这一现象的出现，可能与材料本身的特性密切相关。具体而言，尺寸依赖性在一定程度上降低了分数阶粘弹性压电梁的刚度。而刚度的降低使得纳米梁在振动过程中需要更长的时间来完成一个循环，从而导致其非线性频率显著降低。

(a) (b)

Figure 4. Influence of the fractional order and nonlocal parameter on the time-response of fractional viscoelastic piezoelectric nanobeam

图4. 分数阶阶数和非局部参数对分数阶粘弹性压电纳米梁时间响应的影响

图5(a)显示了粘弹性系数对分数阶粘弹性压电纳米梁的非线性振动时间的响应。从图中可以观察出，当粘弹性系数的增大时，分数阶粘弹性压电纳米梁系统的振动幅值减小。当我们取定粘弹性系数g为0时，此时，不考虑粘性情况，该模型退化为弹性模型。因此，增加粘弹性系数可以使系统的阻尼增加，

(a) (b)

Figure 5. Influence of viscoelastic coefficient and external voltage on the time-response of fractional viscoelastic piezoelectric nanobeam

图5. 粘弹性系数和外加电压对分数阶粘弹性压电纳米梁时间响应的影响

系统的振动减缓。外加电压改变如图5(b)所示，显示了施加正电压降低了纳米梁的非线性频率，这是因为外加的正向电压在纳米梁中产生了轴向压缩力，而负向电压增加了纳米梁的非线性频率，因为外加负电压在纳米梁中产生了张力。给定初始振幅，外加电压从负变正时，非线性频率增大。在正电压条件下，随着电压的减小，系统的非线性振动频率也随之降低；而在负电压条件下，情况则完全相反。从图中可以清晰地推断出，负电压下的非线性振动频率显著高于正电压下的频率。因此，在进行建模和求解问题时，必须综合考虑各种因素对材料的影响，以确保模型的准确性和可靠性。

5. 结论

随着分数阶阶数增大，系统的振幅减小，阻尼增大。

当非局部参数$\mu$ 的增加时，非线性频率减小。

粘弹性系数的增加，降低了分数阶粘弹性压电纳米梁的非线性振动幅值，增加了系统的阻尼。

施加电压为正时，增加电压时系统的非线性频率和振幅都降低；当外加电压为负时，则相反。

综上所述，非局部参数的增加会导致非线性频率的降低；相反，从正电压到负电压的外加电压会导致非线性频率的增加。值得注意的是，由于粘性效应，振动振幅随着时间的推移而减小；然而，随着分数阶与非局部参数和粘弹性系数的增加，系统的阻尼增强，更快地抑制系统的振动。

基金项目

甘肃省自然基金重点项目(25JRRA802)。

NOTES

^*通讯作者。