Dispersion and Attenuation Characteristics of Fractional-Order Thermoelastic Guided Waves in Functionally Graded Piezoelectric Hollow Cylinders
-
摘要: 基于分数阶热电弹性理论和Legendre多项式方法,构建了功能梯度空心圆柱中导波传播的数学模型. 讨论了分数阶次、压电效应、径厚比等对导波传播,特别是对其衰减的影响规律. 数值结果表明,压电效应对衰减的影响主要集中在截止频率和突变频率附近,并使得突变频率发生前移;分数阶对热波模态相速度和衰减的影响较大,且热波相速度存在模态交叉,在交叉频率点附近分数阶对相速度的影响相反;热波衰减随着分数阶增大而逐渐减小;第一阶纵向模态衰减受到了压电效应的抑制,其余模态衰减都显著增大,并且电开路受到的影响要比电短路状态大.Abstract: Based on the fractional-order thermo-electric-elastic theory and the Legendre polynomial series method, a mathematical model for guided wave propagation in functionally graded hollow cylinders was established. The effects of the fractional order, the piezoelectric effect, and the radius-thickness ratio on the wave propagation, especially on its attenuation, were discussed. The numerical analysis results indicate that, the piezoelectric effect on attenuation mainly concentrates near the cutoff frequency and the mutation frequency, and causes the mutation frequency to shift forward. The fractional order has a great impact on the phase velocity and attenuation of the thermal wave mode, and has an opposite impact on the phase velocity around the crossover frequency point where the crossover mode occurs with the thermal wave velocity. But the thermal wave attenuation gradually decreases with the fractional order. Meanwhile, the 1st longitudinal mode attenuation is suppressed by the piezoelectric effect. However, the attenuation of other modes significantly increases, and the impact of the electrical open circuit is greater than that of the electrical short circuit.
-
0. 引言
功能梯度压电材料(FGPM)通过对不同材料的适当组合,可以制造出满足各种需求和改变电弹性性能的期望方向的独特压电材料,同时材料的组合亦由最初的金属/陶瓷扩展到金属/合金,非金属/非金属,非金属/陶瓷,陶瓷/陶瓷等多种组合形式[1]. 同时,工作环境的温度对FGPM力学性能和寿命的影响非常大,研究功能梯度结构的热力学性能是具有重要现实意义的活跃课题.
经典热弹性理论的提出为热力学研究提供了基础,然而其暗示热波是以无限速度传播的,这与实际情况不符. 为弥补经典理论的不足,学者们建立了各种广义热弹理论,如带有单时间因子的L-S(Lord-Shulman)广义热弹动力学理论[2]、将两个不同的热松弛时间因子分别引入了本构方程和能量方程的G-L(Green-Lindsay)理论[3]、无耗散的G-N(Green-Naghdi)热弹理论[4]以及三相滞广义热弹理论. 基于这些广义热弹理论,学者们对热弹固体中的动力学特性进行了深入分析和研究[5-10]. 然而,热电弹结构波动特性的研究,特别是FGPM结构中波动特性的研究不足. 需要指出的是,FGPM结构的波特性与它们在传感器、谐振器和振动滤波器中的性能直接相关,因而引起了广泛的关注. Heydarpour等[11]考虑了带有压电层的功能梯度圆柱形面板中的热弹性响应. Khoshgoftar等[12]研究了功能梯度材料厚壁圆筒的热压电行为,结果发现,通过在FGPM固体结构中应用适当的力和热分布,可以控制FGPM中的应力、电位和电场分布. Dai和Jiang[13]对由FGPM组成的圆柱体在均匀磁场中,在电、热和机械载荷作用下的电磁-热弹性行为进行了分析研究. Ootao等[14]研究了由均匀供热引起的功能梯度热电空心圆柱体瞬态压热弹性问题,材料的非均匀性对温度变化、位移、应力和电势的影响. Arefi和Rahimi[15]推导了基于能量法的FGPM任意结构热弹性分析的一般公式,使用能量法分析了功能梯度压电旋转压力容器的热弹性问题. El-Naggar等[16]基于L-S热弹性理论研究了理想导电平面中的电磁-热弹性体中磁场、旋转、热场、初始应力以及空隙对P波反射的影响. Zhu等[17]基于射线矩阵公式分析了流体填充的功能梯度压电空心圆柱体中的波传播,研究表明只有在小波数和高模的情况下,材料的不均匀性才对波的传播起着重要作用.
由上述文献可知,尽管热弹性波是波动力学领域的一个重要研究方向,但热电弹耦合FGPM管中导波的研究还鲜见报道. 压电空心圆柱作为力电设备中最常见的结构形式之一,被广泛应用于声波器件和机电系统. 然而,在压电结构工作时,其热力电三种物理场互相耦合,结构温度升高,三种能量形式也相互影响. 因此,研究管中的热电弹导波传播对于压电管的结构设计和无损检测具有重要意义. 另外,分数阶微积分在热弹性力学中应用广泛[5],取得了一系列成果,但对热电弹管中的导波传播还未有研究. 本文基于分数阶L-S型广义热电弹耦合理论,构建了FGPM空心圆柱中导波传播的数学模型. 然后采用Legendre多项式级数方法(Legendre polynomial series approach,LPSA)[18-20]求解耦合波动方程,将微分方程转化为线性特征值问题,可直接求得代表导波传播和衰减的复波数解. 讨论了分数阶次、压电效应、径厚比等对导波传播,特别是对其衰减的影响规律.
1. FGPM空心圆柱中热弹波传播特性
1.1 FGPM空心圆柱热弹波传播模型
假设空心圆柱足够长,且材料只在径向表现出非均匀性;a和b分别为内半径和外半径,厚度h=b-a. 根据L-S广义热电弹性耦合理论[5, 8],FGPM空心圆柱中热弹性导波传播问题的控制方程可以表达如下:
$$ \left\{\begin{array}{l} \frac{\partial T_{r r}}{\partial r}+\frac{1}{r} \frac{\partial T_{r \theta}}{\partial \theta}+\frac{\partial T_{r z}}{\partial z}+\frac{T_{r r}-T_{\theta \theta}}{r}=\rho \frac{\partial^{2} u_{r}}{\partial t^{2}}, \\ \frac{\partial T_{r \theta}}{\partial r}+\frac{1}{r} \frac{\partial T_{\theta \theta}}{\partial \theta}+\frac{\partial T_{\theta z}}{\partial z}+\frac{2 T_{r \theta}}{r}=\rho \frac{\partial^{2} u_{\theta}}{\partial t^{2}}, \\ \frac{\partial T_{r z}}{\partial r}+\frac{1}{r} \frac{\partial T_{\theta z}}{\partial \theta}+\frac{\partial T_{z z}}{\partial z}+\frac{T_{r z}}{r}=\rho \frac{\partial^{2} u_{z}}{\partial t^{2}}, \\ \frac{\partial D_{r}}{\partial r}+\frac{1}{r} \frac{\partial D_{\theta}}{\partial \theta}+\frac{\partial D_{z}}{\partial z}+\frac{D_{r}}{r}=0, \\ K_{1} \frac{\partial^{2} T}{\partial r^{2}}+K_{2}\left(\frac{1}{r^{2}} \frac{\partial^{2} T}{\partial \theta^{2}}+\frac{1}{r} \frac{\partial T}{\partial r}\right)+K_{3} \frac{\partial^{2} T}{\partial z^{2}}= \\ \quad\left(1+\frac{\tau_{0}^{\alpha}}{\Gamma(\alpha+1)} \frac{\partial^{\alpha}}{\partial t^{\alpha}}\right)\left(\rho C_{\mathrm{e}} \frac{\partial T}{\partial t}+T_{0} \beta_{1} \frac{\partial \varepsilon_{r r}}{\partial t}+T_{0} \beta_{2} \frac{\partial \varepsilon_{\theta \theta}}{\partial t}+T_{0} \beta_{3} \frac{\partial \varepsilon_{z z}}{\partial t}-T_{0} P_{3} \frac{\partial \varPhi}{\partial z \partial t}\right). \end{array}\right. $$ (1) 几何方程为
$$ \left\{\begin{array}{l} \varepsilon_{\theta \theta}=\frac{1}{r} \frac{\partial u_{\theta}}{\partial \theta}+\frac{u_{r}}{r}, \varepsilon_{r r}=\frac{\partial u_{r}}{\partial r}, \varepsilon_{z z}=\frac{\partial u_{z}}{\partial z}, \\ \varepsilon_{r \theta}=\frac{1}{2}\left(\frac{1}{r} \frac{\partial u_{r}}{\partial \theta}+\frac{\partial u_{\theta}}{\partial r}-\frac{u_{\theta}}{r}\right), \\ \varepsilon_{\theta z}=\frac{1}{2}\left(\frac{\partial u_{\theta}}{\partial z}+\frac{1}{r} \frac{\partial u_{z}}{\partial \theta}\right), \varepsilon_{r z}=\frac{1}{2}\left(\frac{\partial u_{r}}{\partial z}+\frac{\partial u_{z}}{\partial r}\right) ,\\ E_{\theta}=-\frac{1}{r} \frac{\partial \varPhi}{\partial \theta}, E_{r}=-\frac{\partial \varPhi}{\partial r}, E_{z}=-\frac{\partial \varPhi}{\partial z}. \end{array}\right. $$ (2) 考虑无应力边界条件和电学开路条件,正交各向异性空心圆柱的本构关系为
$$ \left\{\begin{array}{l} T_{r r}=\left(C_{11} \varepsilon_{r r}+C_{12} \varepsilon_{\theta \theta}+C_{13} \varepsilon_{z z}-\beta_{1} T-e_{31} E_{z}\right) g(r), \\ T_{\theta \theta}=C_{12} \varepsilon_{r r}+C_{22} \varepsilon_{\theta \theta}+C_{23} \varepsilon_{z z}-\beta_{2} T-e_{32} E_{z}, \\ T_{z z}=C_{13} \varepsilon_{r r}+C_{23} \varepsilon_{\theta \theta}+C_{33} \varepsilon_{z z}-\beta_{3} T-e_{33} E_{z}, \\ T_{\theta z}=2 C_{44} \varepsilon_{\theta z}-e_{24} E_{\theta}, \\ T_{r z}=\left(2 C_{55} \varepsilon_{r z}-e_{15} E_{r}\right) g(r), \\ T_{r \theta}=2 C_{66} \varepsilon_{r \theta} g(r), \\ D_{r}=\left(2 e_{15} \varepsilon_{r z}+\epsilon_{11} E_{r}+P_{1} T\right) g(r), \\ D_{\theta}=2 e_{24} \varepsilon_{\theta z}+\epsilon_{22} E_{\theta}+P_{2} T, \\ D_{z}=e_{31} \varepsilon_{r r}+e_{32} \varepsilon_{\theta \theta}+e_{33} \varepsilon_{z z}+\epsilon_{33} E_{z}+P_{3} T . \end{array}\right. $$ (3) 式(1)—(3)中, Tij,εij,ui分别为应力、应变和位移分量,$C_{i j}, e_{i j}, \epsilon_{i j}$分别为材料的弹性系数、压电常数与介电常数,Di表示电位移分量,Ej为电场强度,Φ为电势,Pj为热释电常数,t和τ0表示时间和松弛时间,βi和Ce分别为恒定应变下的体积膨胀系数和比热,Kj为材料常数,T和T0分别为温度变化和参考温度,其中环境温度T0 =296 K,ρ为密度,α为分数阶次,Γ是Gamma函数,
$$g(r)= \begin{cases}1, & a \leqslant r \leqslant b ,\\ 0, & \text { others }.\end{cases} $$ (4) 在式(3)中,边界条件通过矩形窗函数引入至本构关系,进而在式(8)中引入至控制方程,免去了单独加载边界条件的过程,并且可使最终的特征值问题为线性特征值问题.
绝热边界条件和电学闭路条件后文再叙. 为便于求解,进行无量纲化:
$$ \left\{\begin{array}{l} \hat{C}_{i j}=\frac{C_{i j}}{\bar{C}_{11}}, \hat{\rho}=\frac{\rho}{\bar{\rho}}, \hat{x}_{i}=\frac{v_{x}}{k_{x}} x_{i}, \hat{u}_{i}=\frac{v_{x}^{3} \bar{\rho}}{k_{x} \bar{\beta}_{1} T_{0}} u_{i}, \tau_{0}=\frac{v_{x}^{2}}{k_{x}} t_{0}, \\ \hat{T}=\frac{T}{T_{0}}, \hat{\beta}_{i}=\frac{\beta_{i}}{\bar{\beta}_{1}}, \eta=\frac{\left(\bar{\beta}_{1}\right)^{2} T_{0}}{\bar{\rho}^{2} \bar{C}_{\mathrm{e}} v_{x}^{2}}, \hat{\varPhi}=\frac{v_{x} \bar{e}_{33}}{k_{x} \bar{\beta}_{1} T_{0}} \varPhi, \hat{P}_{i}=\frac{P_{i} \bar{C}_{11}}{\bar{\beta}_{1} \bar{e}_{33}}, \\ \hat{\epsilon}_{i j}=\frac{\epsilon_{i j} \bar{C}_{11}}{\left(\bar{e}_{33}\right)^{2}}, \hat{e}_{i j}=\frac{e_{i j}}{\bar{e}_{33}}, \hat{K}_{i}=\frac{K_{i}}{\bar{K}_{1}}, \hat{C}_{\mathrm{e}}=\frac{C_{\mathrm{e}}}{\bar{C}_{\mathrm{e}}}, \end{array}\right. $$ (5) 其中,$v_{x}=\sqrt{C_{11} / \rho}$为纵波波速,$k_{x}=K_{1} /\left(\rho C_{\mathrm{e}}\right)$为$x$方向上的热扩散系数,$\bar{C}_{11}, \bar{\rho}, \bar{\beta}_{1}, \bar{C}_{\mathrm{e}}, \bar{e}_{33}, \bar{K}_{1}$分别是指定某一材料参数,本文为PZT-5A的相关材料参数.
空心圆柱中轴向导波为z方向传播的自由谐波,则可以假设位移、温度和电势表达式设为
$$\left\{u_{r}, u_{\theta}, u_{z}, T, \varPhi\right\}=\{U(r), V(r), W(r), X(r), Y(r)\} \exp (\mathrm{i} k z+\mathrm{i} N \theta-\mathrm{i} \omega t), $$ (6) 其中,U,V,W分别为r,θ,z方向的位移幅值,X为温度幅值, Y为电势幅值;k是波数,k=Re(k)+i · Im(k),k的实部表示波的传播,虚部表示衰减;ω为无量纲角频率.
类似于文献[21],功能梯度材料的属性可做如下多项式拟合:
$$ \left\{\begin{array}{l} \hat{C}_{i j}(r)=C_{i j}^{(m)}\left(\frac{2 r-(b+a)}{b-a}\right)^{m}, \hat{\rho}_{i}(r)=\rho_{i}^{(m)}\left(\frac{2 r-(b+a)}{b-a}\right)^{m}, \\ \hat{e}_{i j}(r)=e_{i j}^{(m)}\left(\frac{2 r-(b+a)}{b-a}\right)^{m}, \hat{C}_{\mathrm{e}}(r)=C_{\mathrm{e}}^{(m)}\left(\frac{2 r-(b+a)}{b-a}\right)^{m}, \\ \hat{\epsilon}_{i j}(r)=\epsilon_{i j}^{(m)}\left(\frac{2 r-(b+a)}{b-a}\right)^{m}, \hat{K}_{i}(r)=K_{i}^{(m)}\left(\frac{2 r-(b+a)}{b-a}\right)^{m}, \\ \hat{\beta}_{i}(r)=\beta_{i}^{(m)}\left(\frac{2 r-(b+a)}{b-a}\right)^{m}, \hat{P}_{i}(r)=P_{i}^{(m)}\left(\frac{2 r-(b+a)}{b-a}\right)^{m}, \end{array}\right. ~~~ m=0,1,2, \cdots, L \text {. } $$ (7) 式(7)中,L表明材料参数的梯度表达式进行了L阶次的多项式近似. 将式(2)—(7)代入式(1), 简化推导过程,并略去上标“^”,可得
$$ \begin{aligned} & \frac{2 m}{b-a}\left[\frac{2 r-(b+a)}{b-a}\right]^{m-1}\left[C_{11}^{(m)} r^{2} U^{\prime}+C_{12}^{(m)}(\mathrm{i} N r V+r U)+C_{13}^{(m)} \mathrm{i} k r^{2} W-\beta_{1}^{(m)} r^{2} X+e_{31}^{(m)} \mathrm{i} k r^{2} Y\right] g(r)+ \\ & \quad\left[\frac{2 r-(b+a)}{b-a}\right]^{m}\left\{C_{11}^{(m)} r^{2} U^{\prime \prime}+C_{12}^{(m)}\left(\mathrm{i} N r V^{\prime}+r U^{\prime}-\mathrm{i} N V-U\right)+\right. \\ & C_{13}^{(m)} \mathrm{i} k r^{2} W^{\prime}-\beta_{1}^{(m)} r^{2} X^{\prime}+e_{31}^{(m)} \mathrm{i} k r^{2} Y^{\prime}+e_{15}^{(m)} \mathrm{i} k r^{2} Y^{\prime}+ \\ & C_{60}^{(m)}\left(-N^{2} U+\mathrm{i} N r V^{\prime}-\mathrm{i} N V\right)+C_{55}^{(m)}\left(-k^{2} r^{2} U+\mathrm{i} k r^{2} W^{\prime}\right)+C_{11}^{(m)} r U^{\prime}+C_{12}^{(m)}(\mathrm{i} N V+U)+ \\ & C_{13}^{(m)} \mathrm{i} k r W-\beta_{1}^{(m)} r X+e_{31}^{(m)} \mathrm{i} k r Y-\left(C_{12}^{(m)} r U^{\prime}+C_{22}^{(m)}(\mathrm{i} N V+U)+\right. \\ & \left.C_{23}^{(m)} \mathrm{i} k r W-\beta_{2}^{(m)} r X+e_{32}^{(m)} \mathrm{i} k r Y\right)+\left(C_{11}^{(m)} r^{2} U^{\prime}+C_{12}^{(m)}(\mathrm{i} N r V+r U)+\right. \\ & \left.\left.\quad C_{13}^{(m)} \mathrm{i} k r^{2} W-\beta_{1}^{(m)} r^{2} X+e_{31}^{(m)} \mathrm{i} k r^{2} Y\right) \frac{\partial g(r)}{\partial r}+\rho^{(m)} \omega^{2} r^{2} U\right\}=0, \end{aligned} $$ (8a) $$ \begin{aligned} & \frac{2 m}{b-a}\left[\frac{2 r-(b+a)}{b-a}\right]^{m-1} C_{66}^{(m)}\left(\mathrm{i} N r U+r^{2} V^{\prime}-r V\right) g(r)+ \\ & \quad\left[\frac{2 r-(b+a)}{b-a}\right]^{m}\left[C_{66}^{(m)}\left(\mathrm{i} N r U+r^{2} V^{\prime}-r V\right) \frac{\partial \pi(r)}{\partial r}+\rho^{(m)} \omega^{2} r^{2} V\right]+ \\ & \quad\left[\frac{2 r-(b+a)}{b-a}\right]^{m}\left\{C_{66}^{(m)}\left(\mathrm{i} N r U^{\prime}+r^{2} V^{\prime \prime}-r V^{\prime}-\mathrm{i} N U+V\right)+C_{12}^{(m)} \mathrm{i} N r U^{\prime}+\right. \\ & C_{22}^{(m)}\left(-N^{2} V+\mathrm{i} N U\right)-C_{23}^{(m)} k N r W-\beta_{2}^{(m)} \mathrm{i} N r X-e_{32}^{(m)} k N r Y+ \\ & \left.C_{44}^{(m)}\left(-k^{2} r^{2} V-k N r W\right)-e_{24}^{(m)} k N r Y+2 C_{66}^{(m)}\left(\mathrm{i} N U+r V^{\prime}-V\right)\right\} g(r)=0, \end{aligned} $$ (8b) $$ \begin{array}{l} \frac{2 m}{b-a}\left[\frac{2 r-(b+a)}{b-a}\right]^{m-1}\left[C_{55}^{(m)}\left(\mathrm{i} k r^{2} U+r^{2} W^{\prime}\right)+e_{15}^{(m)} r^{2} Y^{\prime}\right] \pi(r)+\\ \quad\left[\frac{2 r-(b+a)}{b-a}\right]^{m}\left[C_{55}^{(m)}\left(\mathrm{i} k r^{2} U+r^{2} W^{\prime}\right)+e_{15}^{(m)} r^{2} Y^{\prime}\right] \frac{\partial g(r)}{\partial r}+ \\ \quad\left[\frac{2 r-(b+a)}{b-a}\right]^{m}\left\{C_{55}^{(m)}\left(\mathrm{i} k r^{2} U^{\prime}+r^{2} W^{\prime \prime}\right)+e_{15}^{(m)} r^{2} Y^{\prime \prime}+C_{44}^{(m)}\left(-k N r V-N^{2} W\right)-\right. \\ \quad e_{24}^{(m)} N^{2} Y+C_{13}^{(m)} \mathrm{i} k r^{2} U^{\prime}+\rho^{(m)} \omega^{2} r^{2} W+C_{23}^{(m)}(-k N r V+\mathrm{i} k r U)-C_{33}^{(m)} k^{2} r^{2} W- \\ \left.\quad \beta_{3}^{(m)} \mathrm{i} k r^{2} X-e_{33}^{(m)} k^{2} r^{2} Y+C_{55}^{(m)}\left(\mathrm{i} k r U+r W^{\prime}\right)+e_{15}^{(m)} r Y^{\prime}\right\}=0, \end{array} $$ (8c) $$ \begin{aligned} & {\left[\frac{2 r-(b+a)}{b-a}\right]^{m}\left\{K_{1}^{(m)} \hat{r}^{2} X^{\prime \prime}+K_{2}^{(m)}\left(-N^{2} X+\hat{r} X^{\prime}\right)-K_{3}^{(m)} k^{2} \hat{r}^{2} X+c \eta \beta_{1}^{(m)} \mathrm{i} \omega \hat{r}^{2} U^{\prime}+\right.} \\ & \left.\quad c \eta\left[\left(-\beta_{2}^{(m)} \omega N \hat{r} V+\beta_{2}^{(m)} \mathrm{i} \omega \hat{r} U\right)-\beta_{3}^{(m)} k \omega \hat{r}^{2} W+P_{3}^{(m)} k \omega \hat{r}^{2} Y\right]\right\}= \\ & \quad-\left[\frac{2 r-(b+a)}{b-a}\right]^{m+n} c \rho^{(m)} C_{\mathrm{e}}^{(n)} \mathrm{i} \omega r^{2} X, \\ \end{aligned} $$ (8d) $$ \begin{array}{l} \frac{2 m}{b-a}\left[\frac{2 r-(b+a)}{b-a}\right]^{m-1}\left[e_{15}^{(m)}\left(\mathrm{i} k r^{2} U+r^{2} W^{\prime}\right)-\epsilon_{11}^{(m)} r^{2} Y^{\prime}+P_{1}^{(m)} r^{2} X\right] g(r)= \\ \quad-\left[\frac{2 r-(b+a)}{b-a}\right]^{m}\left\{\left[e_{15}^{(m)}\left(\mathrm{i} k r^{2} U^{\prime}+r^{2} W^{\prime \prime}\right)-\epsilon_{11}^{(m)} r^{2} Y^{\prime \prime}+\right.\right. \\ P_{1}^{(m)} r^{2} X^{\prime}+e_{24}^{(m)}\left(-k N r V-N^{2} W\right)+\epsilon_{22}^{(m)} N^{2} Y+ \\ P_{2}^{(m)} \mathrm{i} N r X+e_{31}^{(m)} \mathrm{i} k r^{2} U^{\prime}+e_{32}^{(m)}(-k N r V+\mathrm{i} k r U)-e_{33}^{(m)} k^{2} r^{2} W+\epsilon_{33}^{(m)} k^{2} r^{2} Y+ \\ \left.P_{3}^{(m)} \mathrm{i} k r^{2} X+e_{15}^{(m)}\left(\mathrm{i} k r U+r W^{\prime}\right)-\epsilon_{11}^{(m)} r Y^{\prime}+P_{1}^{(m)} r^{2} X\right]+ \\ \left.\left[e_{15}^{(m)}\left(\mathrm{i} k r^{2} U+r^{2} W^{\prime}\right)-\epsilon_{11}^{(m)} r^{2} Y^{\prime}+P_{1}^{(m)} r^{2} X\right] \frac{\partial g(r)}{\partial r}\right\} . \end{array} $$ (8e) 上述等式中$c=1+\left(-\mathrm{i} \omega \tau_{0}\right)^{\alpha} / \Gamma(\alpha+1), ()^{\prime}$和()″分别表示对r的一阶和二阶导数,为求解方程组(8),获得热电弹性类Lamb波的特性,可进一步将位移、温度、电势等物理场振幅变量进行Legendre多项式展开
$$ \{U, V, W, X, Y\}=\sum\limits_{m=0}^{\infty}\left\{p_{m}^{(1)}, p_{m}^{(2)}, p_{m}^{(3)}, p_{m}^{(4)}, p_{m}^{(5)}\right\} Q_{m}(r), $$ (9) 其中
$$Q_{m}(r)=\sqrt{\frac{2 m+1}{b-a}} P_{m}\left(\frac{2 r-(b+a)}{b-a}\right), $$ (10) pmj是Legendre多项式及其展开系数. 随着展开阶m的增加,多项式的高阶多项式项对总体结果的影响降低,因此,对式(9)的求和可以停止在某个值M.
绝热边界需用引入一个新的Legendre多项式级数来表示:
$$ Z(r)=\pi(r) \sum\limits_{m=0}^{M} p_{m}^{(6)} Q_{m}(r). $$ (11) 式(9)中温度表达式X(r)对r求导后结果等于式(11)中的温度梯度表达式,即
$$ \sum\limits_{m=0}^{M} p_{m}^{(4)} \frac{\mathrm{d}}{\mathrm{d} r} Q_{m}(r)=\pi(r) \sum\limits_{m=0}^{M} p_{m}^{(6)} Q_{m}(r). $$ (12) 然后该式两边同乘以Qn(r),n=0, 1, …, M,并进行积分,可得p(4)和p(6)的关系为
$$ \boldsymbol{H} \boldsymbol{p}^{(4)}=\boldsymbol{G} \boldsymbol{p}^{(6)} \Rightarrow \boldsymbol{p}^{(6)}=\boldsymbol{H} \boldsymbol{p}^{(4)}, $$ (13) 式中,G为单位矩阵. 通过温度及其梯度的表达式及其关系,即可利用扩展Legendre多项式法求解绝热边界条件下的热弹性导波传播问题.
类似地,式(8)两边同乘以Qn(r),n=0, 1, …, M,并进行积分,可将耦合方程转为以下的矩阵形式:
$$ k^{2} \boldsymbol{A} \boldsymbol{p}+k \boldsymbol{B} \boldsymbol{p}+\boldsymbol{C} \boldsymbol{p}=\mathbf{0}, $$ (14) 其中
$$ \begin{aligned} \boldsymbol{A} & =\left[\begin{array}{ccccc} \boldsymbol{A}_{11}^{n, m} & \mathbf{0} & \mathbf{0} & \mathbf{0} & \mathbf{0} \\ \mathbf{0} & \boldsymbol{A}_{22}^{n, m} & \mathbf{0} & \mathbf{0} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \boldsymbol{A}_{33}^{n, m} & \mathbf{0} & \boldsymbol{A}_{35}^{n, m} \\ \mathbf{0} & \mathbf{0} & \mathbf{0} & \boldsymbol{A}_{44}^{n, m} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \boldsymbol{A}_{53}^{n, m} & \mathbf{0} & \boldsymbol{A}_{55}^{n, m} \end{array}\right], \boldsymbol{B}=\left[\begin{array}{ccccc} \mathbf{0} & \mathbf{0} & B_{13}^{n, m} & \mathbf{0} & \boldsymbol{B}_{15}^{n, m} \\ \mathbf{0} & \mathbf{0} & \boldsymbol{B}_{23}^{n, m} & \mathbf{0} & \boldsymbol{B}_{25}^{n, m} \\ \boldsymbol{B}_{31}^{n, m} & \boldsymbol{B}_{32}^{n, m} & \mathbf{0} & \boldsymbol{B}_{34}^{n, m} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \boldsymbol{B}_{43}^{n, m} & \mathbf{0} & \boldsymbol{B}_{45}^{n, m} \\ \boldsymbol{B}_{51}^{n, m} & \boldsymbol{B}_{52}^{n, m} & \mathbf{0} & \boldsymbol{B}_{54}^{n, m} & \mathbf{0} \end{array}\right], \\ \boldsymbol{C} & =\left[\begin{array}{ccccc} \boldsymbol{C}_{11}^{n, m} & \boldsymbol{C}_{12}^{n, m} & \mathbf{0} & \boldsymbol{C}_{14}^{n, m} & \mathbf{0} \\ \boldsymbol{C}_{21}^{n, m} & \boldsymbol{C}_{22}^{n, m} & \mathbf{0} & \boldsymbol{C}_{24}^{n, m} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \boldsymbol{C}_{33}^{n, m} & \mathbf{0} & \boldsymbol{C}_{35}^{n, m} \\ \boldsymbol{C}_{41}^{n, m} & \boldsymbol{C}_{42}^{n, m} & \mathbf{0} & \boldsymbol{C}_{44}^{n, m} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \boldsymbol{C}_{53}^{n, m} & \boldsymbol{C}_{54}^{n, m} & \boldsymbol{C}_{55}^{n, m} \end{array}\right], \\ \boldsymbol{p} & =[\begin{array}{rl} \boldsymbol{p}^{(1)}, \boldsymbol{p}^{(2)}, \boldsymbol{p}^{(3)}, \boldsymbol{p}^{(4)}, \boldsymbol{p}^{(5)} ]& { }^{\mathrm{T}} \end{array} . \end{aligned} $$ 矩阵Aij,Bij,Cij的维数分别为(M+1)×(M+1). 为了降低求解难度,引入一个新的波数相关向量:
$$ \boldsymbol{q}=k \boldsymbol{p}, $$ (15) 则可转化为如下特征方程:
$$ \left[\begin{array}{cc} \mathbf{0} & \boldsymbol{I} \\ -\boldsymbol{A}^{-1} \boldsymbol{C} & -\boldsymbol{A}^{-1} \boldsymbol{B} \end{array}\right]\left[\begin{array}{l} \boldsymbol{p} \\ \boldsymbol{q} \end{array}\right]=k\left[\begin{array}{l} \boldsymbol{p} \\ \boldsymbol{q} \end{array}\right], $$ (16) I是单位矩阵,通过式(15),可以快速求解热弹波传播的频散和衰减关系以及物理场分布.
1.2 电短路基本方程与求解
电短路边界条件限制电势在边界处为零,式(8e)可以进行以下更改:
$$ \varPhi(r, \theta, z, t)=(r-a)(r-b) \exp (\mathrm{i} k z+\mathrm{i} N \theta-\mathrm{i} \omega t) Y(r) . $$ (17) 另外,式(8e)中Dr需要删除矩形窗函数. 用式(17)替换式(8e)中相应的表达式,其余等式关系保持不变,可得
$$ \begin{gathered} \frac{2 m}{b-a}\left[\frac{2 r-(b+a)}{b-a}\right]^{m-1}\left[e_{15}^{(m)}\left(\mathrm{i} k r^{2} U+r^{2} W^{\prime}\right)-\epsilon_{11}^{(m)} r^{2}[(r-a)(r-b) Y]^{\prime}+P_{1}^{(m)} r^{2} X\right]+ \\ {\left[\frac{2 r-(b+a)}{b-a}\right]^{m}\left\{e_{15}^{(m)}\left(\mathrm{i} k r^{2} U^{\prime}+r^{2} W^{\prime \prime}\right)-\epsilon_{11}^{(m)} r^{2}[(r-a)(r-b) Y]^{\prime \prime}+\right.} \\ P_{1}^{(m)} r^{2} X^{\prime}+e_{24}^{(m)}\left(-k N r V-N^{2} W\right)+\epsilon_{22}^{(m)} N^{2}[(r-a)(r-b) Y]+P_{2}^{(m)} \mathrm{i} N r X+ \\ e_{31}^{(m)} \mathrm{i} k r^{2} U^{\prime}+e_{32}^{(m)}(-k N r V+\mathrm{i} k r U)-e_{33}^{(m)} k^{2} r^{2} W+\epsilon_{33}^{(m)} k^{2} r^{2}[(r-a)(r-b) Y]+ \\ \left.P_{3}^{(m)} \mathrm{i} k r^{2} X+e_{15}^{(m)}\left(\mathrm{i} k r U+r W^{\prime}\right)-\epsilon_{11}^{(m)} r[(r-a)(r-b) Y]^{\prime}+P_{1}^{(m)} r X\right\}=0 . \end{gathered} $$ (18) 上述公式的主要变化在于电势Y(r),后续过程与电学开路相同,可得到电短路边界条件的波数复数解.
2. 数值验证和收敛性判断
2.1 数值验证
为了验证LPSA的有效性,图 1对比了本文结果和已有结果. 图 1(a)为相速度, 图 1(b)为虚部. 另外,如前文所述,由于缺乏可用的比较对象,本文进行如下2个退化比较研究:图 1(a)计算了PZT-4圆柱体中的纵波传播,并与文献[22]中的结果进行对比,周向阶数为1,材料参数如表 1所示;图 1(b)采用本文方法和GMM[23]计算了cobalt/steel功能梯度空心圆柱中热弹波的衰减曲线,分数阶次α=0.5,材料参数如表 2所示. 两项对比表明:本文计算结果与已知文献结果均一致,表明了本文所提出的LPSA的有效性.
property C11/(N·m-2) C12/(N·m-2) C13/(N·m-2) C22/(N·m-2) C23/(N·m-2) C33/(N·m-2) PZT-4 1.39×1011 7.78×1010 7.43×1010 1.39×1011 7.43×1010 1.15×1011 property C44/(N·m-2) C55/(N·m-2) C66/(N·m-2) e15/(C·m-2) e24/(C·m-2) e31/(C·m-2) PZT-4 2.56×1010 2.56×1010 3.06×1010 12.7 12.7 -5.2 property e32/(C·m-2) e33/(C·m-2) ε11/ (C2·N-1·m-2) ε22/(C2·N-1·m-2) ε33/(C2·N-1·m-2) ρ/(kg·m-3) PZT-4 -5.2 15.1 6.46×10-9 6.46×10-9 5.62×10-9 7.5×103 material C11/(N·m-2) C12/(N·m-2) C13/(N·m-2) C33/(N·m-2) C44/(N·m-2) C66/(N·m-2) steel 2.692 3×1011 1.153 8×1011 1.153 8×1011 2.692 3×1011 7.692×1010 7.692×1010 cobalt 3.071×1011 1.65×1011 1.027×1011 3.581×1011 7.55×1010 7.105×1010 material ρ/(kg·m-3) Ce/(J·kg-1·K-1) β1/(N·K-1·m-2) β3/(N·K-1·m-2) K1/(W·m-1·K-1) K3/(W·m-1·K-1) steel 7.85×103 477 6.0×106 6.0×106 52 52 cobalt 8.836×103 427 7.04×106 6.9×106 69 69 2.2 收敛性判断
本小节及后文继续分析讨论了方法的收敛性,以及压电效应对准弹性波和热波的影响. FGPM材料由PZT-5A和CdSe组成,材料参数[24]见表 3. 功能梯度材料参数可以表示为P(r)=V1P1(r)+V2P2(r),Vi, Pi(r)分别表示上表面和下表面的体积分数和相应的性质,∑Vi=1.
property CdSe PZT-5A property CdSe PZT-5A C11/(N·m-2) 7.41×1010 1.39×1011 ε11/(C2·N-1·m-2) 8.26×10-11 6.00×10-9 C12/(N·m-2) 4.52×1010 7.78×1010 ε22/(C2·N-1·m-2) 8.26×10-11 6.00×10-9 C22/(N·m-2) 7.41×1010 1.39×1011 ε33/(C2·N-1·m-2) 9.03×10-11 5.47×10-9 C13/(N·m-2) 3.93×1010 7.54×1010 K1/(W·m-1·K-1) 9 1.5 C23/(N·m-2) 3.93×1010 7.54×1010 K2/(W·m-1·K-1) 9 1.5 C33/(N·m-2) 8.36×1010 1.13×1011 K3/(W·m-1·K-1) 9 1.5 C44/(N·m-2) 1.32×1010 2.56×1010 β1/(N·K-1·m-2) 6.21×105 1.52×106 C55/(N·m-2) 1.32×1010 2.56×1010 β2/(N·K-1·m-2) 6.21×105 1.52×106 C66/(N·m-2) 1.445×1010 3.06×1010 β3/(N·K-1·m-2) 5.51×105 1.53×106 e31/(C·m-2) -0.16 -6.98 P1/(C·K-1·m-2) 0 0 e32/(C·m-2) -0.16 -6.98 P2/(C·K-1·m-2) 0 0 e33/(C·m-2) 0.347 13.8 P3/(C·K-1·m-2) -2.94×10-6 -4.52×10-4 e15/(C·m-2) -0.138 13.4 Ce/(J·kg-1·K-1) 260 420 e24/(C·m-2) -0.138 13.4 ρ/(kg·m-3) 5.504×103 7.75×103 图 2和3为周向阶数N=3时的非轴对称弯曲模态实部与虚部收敛图及其局部放大图. 图 2(a)为准弹性波频散曲线,图 3(a)为准弹性波衰减曲线, 右图均为局部放大图. 可以明显看出,当M=6和7时,图 2(a)中的前四阶色散曲线已经完全一致,随着M的增加,后续更多的模态曲线也越来越一致,图 2(b)局部放大图对比显示得更加明显,在M=8时,实部已经完全收敛. 在图 3所示的衰减曲线上也可以发现相同的情况,这表明前四阶模态在M=6时就已经是基本收敛的,由此可以判断,图中所示模态无论是频散曲线还是衰减曲线在M=8时便已经完全收敛. 该算例表明LSPA具有良好的收敛性,且是逐阶收敛的. 在后文分析过程中采用M=10.
3. 结果分析
3.1 压电效应对波的传播和衰减的影响
研究应力波在空心圆柱壳中的传播时,应分别研究3种不同的传播模态,即纵向模态L(0, m)、扭转模态T(0, m)和弯曲模态F(N, m),考虑到所有沿z轴方向传播的波模态,其中纵向模态和扭转模态是轴对称模态,而弯曲模态F(N, m)是非轴对称模态. 周向阶次N=0, 1, 2, 3, …,模数m=1, 2, 3, …,其中,当N=0时,表示解耦的轴对称模态;当N≥1时,表示非轴对称弯曲模态. 为充分研究压电效应对空心圆柱中纵向导波的影响,在这里分别取两个不同的N值进行研究,N=0时的轴对称模态和N=3时的非轴对称弯曲模态.
图 4—6所示为周向阶数N=0时,轴对称模态的热波与准弹性波模态的相速度与衰减曲线. 图 4(a)为相速度曲线, 图 4(b)为衰减曲线; 图 5(a)为纵向模态, 图 5(b)为扭转模态; 图 6(a)为衰减图, 图 6(b)为局部放大图. 由图 4可以看出,在轴对称模态下,热波的相速度与衰减受到压电效应的影响均比较微小. 但对于准弹性波来说,则截然不同,图 5(a)中纵向模态受到压电效应影响较大,压电效应增强了各纵向模态的相速度,且电开路的相速度要大于电短路的情况,而且模态越高,则影响越明显,而图 5(b)中扭转模态则表现完全不同,其相速度几乎不受压电效应影响. 对于衰减来说,从图 6中的衰减曲线中可以看到,仅第一阶的纵向模态衰减受到了压电效应的抑制,其余模态衰减都显著增大,并且电开路受到的影响要比电短路状态大.
周向阶数N=3时,非对称模态的准弹性波模态相速度与衰减曲线如图 7、8所示. 图 8(a)为衰减图, 图 8(b)为局部放大图. 由于压电效应对热波影响微弱,在此不再考虑. 从图 8可以看出,纵向弯曲模态在Ω=1.5时发生模态突变,使得相速度与衰减都发生突变,这与文献[23]情况一致,压电效应使得突变频率发生略微前移. 其原因在于非轴对称模态中(N≠0), 弯曲扭转模态和弯曲纵向模态的接近,导致能量耗散的突然变化. 模态衰减突变时,衰减的极速下降表明在该频率波的能量耗散较小,有利于进行无损检测. 对比几幅图可以发现,与轴对称模态不同的是,在非轴对称模态中,热波的相速度与衰减在电开路情况下与无压电情况相差不大,但在电短路情况下,热波的相速度与衰减都有明显减小. 准弹性波模态的相速度和衰减总体上也受压电影响而增大,此现象随模态升高,愈发明显,另外压电效应对衰减的影响主要集中在截止频率和突变频率附近. 压电效应的分析结果表明,电弹耦合对不同模态传播和衰减的影响强弱不同,这些分析对于声波器件的设计具有指导意义.
3.2 分数阶的影响
不同分数阶α对线性梯度空心圆柱中纵向导波频散和衰减曲线的影响如图 9—11所示,此处选取周向阶数N=3. 图 9(a)为相速度曲线, 图 9(b)为衰减曲线;图 11(a)为衰减图, 图 11(b)为局部放大图. 从图 9中可以看出,热波模态的相速度和衰减受到分数阶影响较大,热波相速度存在模态交叉的频率点,在频率点前后α对相速度的影响相反,但是热波衰减则是随着α的增大而逐渐减小. 图 10、11显示,准弹性波模态衰减随着α增大而减小,但其相速度不会受到分数阶变化的影响. 时间分数阶表征扩散行为的历史记忆性,因此调整分数阶次可以更加匹配实际的传播行为.
3.3 径厚比的影响
图 12—14给出了不同半径厚度比μ下的纵向弯曲模态内第一热波与准弹性波的频散和衰减情况,在这里,周向阶数N=3. 图 12(a)为相速度曲线, 图 12(b)为衰减曲线;图 14(a)为衰减图, 图 14(b)为局部放大图. 随着μ的增大,热波相速度与衰减均降低,同时准弹性波的相速度减小,衰减增大. 在截止频率和突变频率处表现尤为明显,同样有在高频下逐渐收敛的趋势. 这表明径厚比对热波和弹波模态的相速度与衰减具有显著影响,可为圆柱壳的超声无损检测提供理论依据.
3.4 梯度的影响
梯度阶次L代表两种材料组合而成的FGPM参数因配比修改后发生改变,引起的空心圆柱内轴向热弹性导波频散衰减变化情况,如图 15—17所示. 图 15(a)为相速度曲线, 图 15(b)为衰减曲线;图 17(a)为衰减图, 图 17(b)为局部放大图. 由图 15可以明显看出,梯度阶次改变对第一热波的影响相当之大,在梯度阶次呈线性变化由L=1变为L=2时,第一热波模态的相速度明显变小,而衰减变化则完全相反,并且随频率增大,L对衰减的影响越来越大. 而从图 16、17可观察得出,准弹性波相速度随L变大而增大,而衰减曲线则发生较大改变,但是总体也呈现出增大趋势,并且梯度阶次使波的模态突变频率与截止频率延后,在高阶模态上表现尤为突出. 梯度的研究结果表明,通过梯度函数的设计,可以满足不同波速、衰减和截止频率的声波器件设计要求.
4. 结论
1) 压电效应抑制第一阶纵向模态的衰减,但增强了其余模态衰减;并且电开路边界受到的影响要比电短路状态大.
2) 压电效应对衰减的影响主要集中在截止频率和突变频率附近;同时,压电效应使得突变频率发生前移.
3) 分数阶对热波模态相速度和衰减的影响较大;同时,热波相速度存在模态交叉,在交叉频率附近α对相速度的影响相反,但是热波衰减则是随着α增大而逐渐减小.
4) 热波相速度与衰减随着径厚比增大而减小;而随着径厚比减小,准弹性波的相速度增大,衰减减小.
5) 随着梯度阶次增大,第一热波模态的相速度减少,而其衰减变化则完全相反;同时,准弹性波相速度和衰减增大. 并且梯度阶次使波的模态突变频率与截止频率延后,在高阶模态上表现尤为突出.
-
property C11/(N·m-2) C12/(N·m-2) C13/(N·m-2) C22/(N·m-2) C23/(N·m-2) C33/(N·m-2) PZT-4 1.39×1011 7.78×1010 7.43×1010 1.39×1011 7.43×1010 1.15×1011 property C44/(N·m-2) C55/(N·m-2) C66/(N·m-2) e15/(C·m-2) e24/(C·m-2) e31/(C·m-2) PZT-4 2.56×1010 2.56×1010 3.06×1010 12.7 12.7 -5.2 property e32/(C·m-2) e33/(C·m-2) ε11/ (C2·N-1·m-2) ε22/(C2·N-1·m-2) ε33/(C2·N-1·m-2) ρ/(kg·m-3) PZT-4 -5.2 15.1 6.46×10-9 6.46×10-9 5.62×10-9 7.5×103 material C11/(N·m-2) C12/(N·m-2) C13/(N·m-2) C33/(N·m-2) C44/(N·m-2) C66/(N·m-2) steel 2.692 3×1011 1.153 8×1011 1.153 8×1011 2.692 3×1011 7.692×1010 7.692×1010 cobalt 3.071×1011 1.65×1011 1.027×1011 3.581×1011 7.55×1010 7.105×1010 material ρ/(kg·m-3) Ce/(J·kg-1·K-1) β1/(N·K-1·m-2) β3/(N·K-1·m-2) K1/(W·m-1·K-1) K3/(W·m-1·K-1) steel 7.85×103 477 6.0×106 6.0×106 52 52 cobalt 8.836×103 427 7.04×106 6.9×106 69 69 property CdSe PZT-5A property CdSe PZT-5A C11/(N·m-2) 7.41×1010 1.39×1011 ε11/(C2·N-1·m-2) 8.26×10-11 6.00×10-9 C12/(N·m-2) 4.52×1010 7.78×1010 ε22/(C2·N-1·m-2) 8.26×10-11 6.00×10-9 C22/(N·m-2) 7.41×1010 1.39×1011 ε33/(C2·N-1·m-2) 9.03×10-11 5.47×10-9 C13/(N·m-2) 3.93×1010 7.54×1010 K1/(W·m-1·K-1) 9 1.5 C23/(N·m-2) 3.93×1010 7.54×1010 K2/(W·m-1·K-1) 9 1.5 C33/(N·m-2) 8.36×1010 1.13×1011 K3/(W·m-1·K-1) 9 1.5 C44/(N·m-2) 1.32×1010 2.56×1010 β1/(N·K-1·m-2) 6.21×105 1.52×106 C55/(N·m-2) 1.32×1010 2.56×1010 β2/(N·K-1·m-2) 6.21×105 1.52×106 C66/(N·m-2) 1.445×1010 3.06×1010 β3/(N·K-1·m-2) 5.51×105 1.53×106 e31/(C·m-2) -0.16 -6.98 P1/(C·K-1·m-2) 0 0 e32/(C·m-2) -0.16 -6.98 P2/(C·K-1·m-2) 0 0 e33/(C·m-2) 0.347 13.8 P3/(C·K-1·m-2) -2.94×10-6 -4.52×10-4 e15/(C·m-2) -0.138 13.4 Ce/(J·kg-1·K-1) 260 420 e24/(C·m-2) -0.138 13.4 ρ/(kg·m-3) 5.504×103 7.75×103 -
[1] 沈璐璐, 蔡方圆, 杨博. 功能梯度压电板柱面弯曲的弹性力学解[J]. 应用数学和力学, 2023, 44(3): 272-281. doi: 10.21656/1000-0887.430224SHEN Lulu, CAI Fangyuan, YANG Bo. Elasticity solutions for cylindrical bending of functionally graded piezoelectric material plates[J]. Applied Mathematics and Mechanics, 2023, 44(3): 272-281. (in Chinese) doi: 10.21656/1000-0887.430224 [2] LORD H W, SHULMAN Y. A generalized dynamical theory of thermoelasticity[J]. Journal of the Mechanics and Physics of Solids, 1967, 15(5): 299-309. doi: 10.1016/0022-5096(67)90024-5 [3] GREEN A E, LINDSAY K. Thermoelasticity[J]. Journal of elasticity, 1972, 2(1): 1-7. doi: 10.1007/BF00045689 [4] GREEN A, NAGHDI P. On undamped heat waves in an elastic solid[J]. Journal of Thermal Stresses, 1992, 15(2): 253-264. doi: 10.1080/01495739208946136 [5] 段晓宇, 马永斌. 分数阶热弹理论下重力场对二维纤维增强介质的影响[J]. 应用数学和力学, 2021, 42(5): 452-459. doi: 10.21656/1000-0887.410125DUAN Xiaoyu, MA Yongbin. Effects of the gravity field on 2D fiber-reinforced media under the fractional order theory of thermoelasticity[J]. Applied Mathematics and Mechanics, 2021, 42(5): 452-459. (in Chinese) doi: 10.21656/1000-0887.410125 [6] 刘旭, 姚林泉. 热环境中旋转功能梯度纳米环板的振动分析[J]. 应用数学和力学, 2020, 41(11): 1224-1236. doi: 10.21656/1000-0887.410090LIU Xu, YAO Linquan. Vibration analysis of rotating functionally gradient nano annular plates in thermal environment[J]. Applied Mathematics and Mechanics, 2020, 41(11): 1224-1236. (in Chinese) doi: 10.21656/1000-0887.410090 [7] ABD-ALLA A, ABO-DAHAB S, AHMED S, et al. Rayleigh surface wave propagation in an orthotropic rotating magneto-thermoelastic medium subjected to gravity and initial stress[J]. Mechanics of Advanced Materials and Structures, 2020, 27(16): 1400-1411. doi: 10.1080/15376494.2018.1512019 [8] SHARMA J N, PAL M. Propagation of Lamb waves in a transversely isotropic piezothermoelastic plate[J]. Journal of Sound & Vibration, 2004, 270(4/5): 587-610. [9] 王现辉, 李方琳, 刘宇建, 等. 板中热弹波传播: 一种改进的勒让德多项式方法[J]. 力学学报, 2020, 52(5): 1277-1285. https://www.cnki.com.cn/Article/CJFDTOTAL-LXXB202005007.htmWANG Xianhui, LI Fanglin, LIU Yujian, et al. Thermoelastic wave propagation in plates: an improved Legendre polynomial approach[J]. Chinese Journal of Theoretical and Applied Mechanics, 2020, 52(5): 1277-1285. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-LXXB202005007.htm [10] PENG W, CHEN L, HE T. Nonlocal thermoelastic analysis of a functionally graded material microbeam[J]. Applied Mathematics and Mechanics, 2021, 42(6): 855-870. doi: 10.1007/s10483-021-2742-9 [11] HEYDARPOUR Y, MALEKZADEH P, DIMITRI R, et al. Thermoelastic analysis of functionally graded cylindrical panels with piezoelectric layers[J]. Applied Sciences, 2020, 10(4): 1397. doi: 10.3390/app10041397 [12] KHOSHGOFTAR M, ARANI A G, AREFI M. Thermoelastic analysis of a thick walled cylinder made of functionally graded piezoelectric material[J]. Smart Materials and Structures, 2009, 18(11): 115007. doi: 10.1088/0964-1726/18/11/115007 [13] DAI H L, JIANG H J. Analytical study for electromagnetothermoelastic behavior of a functionally graded piezoelectric solid cylinder[J]. Mechanics of Advanced Materials and Structures, 2013, 20(10): 811-818. doi: 10.1080/15376494.2012.676715 [14] OOTAO Y, AKAI T, TANIGAWA Y. Transient piezothermoelastic analysis for a functionally graded thermopiezoelectric hollow cylinder[J]. Journal of Thermal Stresses, 2008, 31(10): 935-955. doi: 10.1080/01495730802250508 [15] AREFI M, RAHIMI G. General formulation for the thermoelastic analysis of an arbitrary structure made of functionally graded piezoelectric materials, based on the energy method[J]. Mechanical Engineering, 2011, 62: 221-235. [16] EL-NAGGAR A, KISHKA Z, ABD-ALLA A, et al. On the initial stress, magnetic field, voids and rotation effects on plane waves in generalized thermoelasticity[J]. Journal of Computational and Theoretical Nanoscience, 2013, 10(6): 1408-1417. doi: 10.1166/jctn.2013.2862 [17] ZHU J, CHEN W, YE G, et al. Waves in fluid-filled functionally graded piezoelectric hollow cylinders: a restudy based on the reverberation-ray matrix formulation[J]. Wave Motion, 2013, 50(3): 415-427. doi: 10.1016/j.wavemoti.2012.10.006 [18] LIU C, YU J, XU W, et al. Theoretical study of elastic wave propagation through a functionally graded micro-structured plate base on the modified couple-stress theory[J]. Meccanica, 2020, 55: 1153-1167. doi: 10.1007/s11012-020-01156-8 [19] OTHMANI C, ZHANG H, LV C, et al. Orthogonal polynomial methods for modeling elastodynamic wave propagation in elastic, piezoelectric and magneto-electro-elastic composites: a review[J]. Composite Structures, 2022, 286: 115245. [20] ZHENG M, MA H, LYU Y, et al. Derivation of circumferential guided waves equations for a multilayered laminate composite hollow cylinder by state-vector and Legendre polynomial hybrid formalism[J]. Composite Structures, 2021, 255: 112950. [21] CAO X, JIN F, JEON I. Calculation of propagation properties of Lamb waves in a functionally graded material (FGM) plate by power series technique[J]. NDT & E International, 2011, 44(1): 84-92. [22] SHATALOV M Y, EVERY A G, YENWONG-FAI A S. Analysis of non-axisymmetric wave propagation in a homogeneous piezoelectric solid circular cylinder of transversely isotropic material[J]. International Journal of Solids & Structures, 2010, 46(3/4): 837-850. [23] WANG X H, LI F L, ZHANG B, et al. Wave propagation in thermoelastic inhomogeneous hollow cylinders by analytical integration orthogonal polynomial approach[J]. Applied Mathematical Modelling, 2021, 99: 57-80. [24] GUHA S, SINGH A K. Plane wave reflection/transmission in imperfectly bonded initially stressed rotating piezothermoelastic fiber-reinforced composite half-spaces[J]. European Journal of Mechanics A: Solids, 2021, 88: 104242. -