近年来,随着纳米尺度力学研究的兴起,传统连续介质力学遇到了新挑战. 比如,根据传统连续力学计算结果,微裂纹尖端的应力 场是奇异的,但这个结论在物理上难以解释,且已有实验证明,微裂纹尖端的应力分布是非奇异的. 再如,声子散射实验证明了高频弹性波是弥散的,其传播速度与频率有关;然而在传统连续力学中波速是一个不依赖于频率的常数. 所以,寻找和发展适应于纳米力学研究的新途径和新方法,是当前连续力学的研究热点之一. 目前,针对纳米力学开展的实验研究、原子模拟及新连续模型已日趋成熟. 其中,在不同于传统连续介质力学的新连续理论中,由艾林根等[1, 2]提出的非局部理论在纳米力学研究中获得广泛关注.
非局部连续力学将分子间作用力是长程力的思想,引入到了传统连续介质力学中. 该理论认为,连续体内某一点的应力不仅取决于该点的应变,还是连续体内所有点的应变及应变历史的函数[1]. 已有大量文献报道非局部理论在位错力学、断裂力学、波传播理论、流体力学等领域的应用[3, 4, 5]. 特别地,非局部理论能够很好地解决微尺度下传统连续力学难以处理的问题,比如根据非局部理论,微裂纹尖端的应力是非奇异的, 这是一个合理的结果. 因此,非局部理论已成为分析纳米力学的主要手段之一[6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21],在研究碳纳米管、石墨烯等热门纳米材料的力学性能上 扮演着重要角色.
令人奇怪的是,利用非局部理论研究纳米结构力学行为时,一些学者的工作表明[6, 7, 8, 9, 10, 11, 12, 13],考虑非局部效应的纳米结构 的弯曲刚度降低;另一些学者[14, 15, 16, 17, 18, 19, 20, 21]的工作却预示了与之相反的结果,即纳米结构的非局部效应显著增强其弯曲刚度. 这些彼此矛盾的研究结论引起了工程技术人员的困惑,也在一定程度上不利于纳米技术的快速发展. 这两种观点均有分子动力学模拟的支持[9, 15],特别地,通过对比[15],发现刚度降低与刚度提高情形分别位于分子动 力学模拟的上下方或者左右两侧,且接近程度相当,因此仅以分子动力学模拟无法解决这个争议. 本文基于梯度型非局部微分本构,结合迭代求解法与泰勒级数展开法,以相对简洁的途径求得了非局部高阶应力表达,讨论了非局 部弹性应力场的物理意义,并进一步推导并发展了非局部高阶梁理论的挠曲轴微分方程.} 然后,以非局部高阶梁模型考察了纳米梁的非局部弯曲问题,利用级数法求得了弯曲挠度的精确非局部表达式,并在数值算例中分析 了挠度的非局部效应,并由此判断纳米结构在考虑非局部效应情况下弯曲刚度的变化. 结果证明:在梯度型非局部微分本构框架下,纳米结构的非局部刚度可能降低也可能提高也可能不变,刚度的变化取决于纳米梁的边 界支撑条件以及外载荷的具体形式.
1 梯度型非局部梁的高阶弯曲模型针对各向同性的均匀弹性连续体,艾林根提出的非局部弹性理论的核心思想可分别归纳描述 为[1, 2]
${\sigma _{ij}}\left( r \right) = \int {_V} \alpha \left( {\left| {r' - r} \right|,\tau } \right){{\sigma '}_{ij}}\left( {r'} \right)\;{\rm{d}}V\left( {r'} \right)$ | (1) |
$\sigma _{ij}-\left( {e_0 a} \right)^2\nabla ^2\sigma _{ij} = {\sigma }'_{ij} $ | (2) |
其中,$\sigma_{ij}$为非局部应力张量,${ r}$为所考察连续体内某参考点,$\alpha $为非局部权函数,它与距离$\left| {{r^\prime } - r} \right|$以及非局部尺度因子$\tau $有关,$\sigma'_{ij}$为经典弹性应力张量,$V$为弹性体体积,$e_{0}$为材料常数,$a$为内特征尺度,$\nabla ^2 = \dfrac{\partial ^2}{\partial r_k \partial r_k }$是拉普拉斯算子. 式(1)通常称为非局部积分本构,其在数学求解上不方便. 相比之下,转换之后的非局部微分本构式(2)则容易求解,它将非局部应力与经典应力以梯度关系相关联,因此也称为梯度型非局部本构.
基于式(2)进一步化简可得非局部梁理论下的应力本构关系为
$\sigma - {\mkern 1mu} {\left( {{e_0}a} \right)^2}\frac{{{{\rm{d}}^2}\sigma }}{{{\rm{ d}}{x^2}}} = \sigma '$ | (3) |
其中,$x$为轴坐标. 在建立梁的弯曲变形方程时,横向坐标的方向需要事先规定,否则容易造成符号混乱甚至错误,有些文 献(比如[7, 14])就存在这个问题. 本文在梁截面的纵向对称轴上建立横向坐标$y$,并取其向上为正,则梁的经典弯曲正应力为
$\sigma ' = - Ey\frac{{{{\rm{d}}^2}w}}{{{\rm{d}}{x^2}}}$ | (4) |
其中,$E$为弹性模量,$w$为弯曲挠度,其方向同$y$向上为正. 采用迭代法求解方程(3),其中经典应力 $\sigma'$作为非局部应力 $\sigma $ 的零阶近似 $\sigma_{0}$,得
${\sigma _n} - {\left( {{e_0}a} \right)^2}\frac{{{d^2}{\sigma _{n - 1}}}}{{{\rm{d}}{x^2}}} = \sigma '$ | (5) |
计算求得第$n$阶非局部应力 $\sigma_{n}$为
${\sigma _n} = - Ey\sum\limits_{m = 0}^n {{{\left( {{e_0}a} \right)}^{2m}}} \frac{{{d^{2m + 2}}w}}{{d{x^{2m + 2}}}}$ | (6) |
对式(6)取极限,即可得精确非局部应力的表达为
$\begin{array}{l} \sigma = - Ey\sum\limits_{m = 0}^\infty {{{\left( {{e_0}a} \right)}^{2m}}} \frac{{{{\rm{d}}^{2m + 2}}w}}{{{\rm{d}}{x^{2m + 2}}}} = \\ \qquad \sigma ' - Ey\sum\limits_{m = 1}^\infty {{{\left( {{e_0}a} \right)}^{2m}}} \frac{{{{\rm{d}}^{2m + 2}}w}}{{{\rm{d}}{x^{2m + 2}}}} \end{array}$ | (7) |
梯度型非局部高阶应力可以视为由经典弯曲应力与非局部无穷级数项组成,其中无穷级数通项为非局部挠度对轴坐标的导数,并且与 非局部参数$e_{0}a$有关. 因此,梯度型非局部高阶应力恰为非局部挠度的逐阶梯度之和,其第一项即为经典应力值. 事实上,非局部高阶应力式(7)还可以相对简洁地通过直接泰勒展开求出,即
$\begin{array}{l} \sigma = \frac{1}{{\left[{1 - {{\left( {{e_0}a} \right)}^2}\frac{{{{\rm{d}}^2}}}{{{\rm{d}}{x^2}}}} \right]\frac{1}{{\sigma '}}}} = \\ \qquad \sigma ' + {\left( {{e_0}a} \right)^2}\frac{{{{\rm{d}}^2}\sigma '}}{{{\rm{d}}{x^2}}} + {\left( {{e_0}a} \right)^4}\frac{{{{\rm{d}}^4}\sigma '}}{{{\rm{d}}{x^4}}} + \cdots = \\ \qquad \sum\limits_{m = 0}^\infty {{{\left( {{e_0}a} \right)}^{2m}}} \frac{{{{\rm{d}}^{2m}}\sigma '}}{{{\rm{d}}{x^{2m}}}} = - Ey\sum\limits_{m = 0}^\infty {{{\left( {{e_0}a} \right)}^{2m}}} \frac{{{{\rm{d}}^{2m + 2}}w}}{{{\rm{d}}{x^{2m + 2}}}} \end{array}$ | (8) |
以往文献对非局部应力的处理方式有两类,一类是回避非局部应力,直接根据非局部应力本构式(3)推导非局部弯矩关系,再由非局部 弯矩关系结合纳米梁受力平衡方程推得问题的非局部控制方 程[7, 8];另一类是通过求解 式(3)的微分方程得到非局部应力 $\sigma $ 的通解,但是在确定积分常数时往往缺乏直接说服力,且程序较为繁杂[14, 15]. 相比之下,本文提出的求解思路和方法显得更加简洁. 进一步地,非局部梁模型的高阶挠曲轴微分方程为
$\frac{M}{{E{I_z}}} = \frac{{{{\rm{d}}^2}w}}{{{\rm{d}}{x^2}}} + {\rm{ }}\sum\limits_{m = 1}^\infty {{{\left( {{e_0}a} \right)}^{2m}}} \frac{{{{\rm{d}}^{2m + 2}}w}}{{{\rm{d}}{x^{2m + 2}}}}$ | (9) |
其中${I_z} = {\rm{ }}\int {_A} {y^2}{\rm{ d}}A$为惯性矩,$A$为梁的横截面积. 弯矩方向以使得梁发生上凹变形为正,当$y$坐标向上为正时,$\sigma $ 与$y$总是异号的,因此非局部弯矩表达为$M = - {\rm{ }}\int {_A} \sigma {\rm{ d}}A$. 显然,当$e_{0}a=0$即忽略结构的非局部效应,高阶挠曲轴微分方程退化为经典挠曲轴近似微分方程,其中弯矩与弯曲变形曲率的正 负号总是一致的.
2 高阶弯曲方程的正则摄动解法考虑欧拉-伯努利梁模型的非局部弯曲问题,假设梁受到集度为$q_{0}$的横向均布载荷以及右端面力 偶矩$M_{\rm e}$的作用,如图 1所示,此时高阶挠曲轴微分方程为
$\begin{array}{l} \frac{{{q_0}L}}{{2E{I_z}}}x - \frac{{{q_0}}}{{2E{I_z}}}{x^2} + \frac{{{M_{\rm{e}}}x}}{{E{I_z}L}} = \\ \qquad \frac{{{{\rm{d}}^2}w}}{{{\rm{d}}{x^2}}} + \sum\limits_{m = 1}^\infty {{{\left( {{e_0}a} \right)}^{2m}}} \frac{{{{\rm{d}}^{2m + 2}}w}}{{{\rm{d}}{x^{2m + 2}}}} \end{array}$ | (10) |
其中,$L$为梁的长度. 文献中要么没有直接推导出非局部高阶梁的弯曲方程(其控制方程阶数与相应的经典理论一致)[7, 8],要 么采用较复杂的过程推得类似的但形式不同的高阶方程[14, 15, 16],而我们这里通过相对简单的途径,得出了 非局部高阶微分方程.
方程(10)的无量纲形式为
$\begin{array}{l} \frac{{{{\bar q}_0}}}{2}\bar x - \frac{{{{\bar q}_0}}}{2}{{\bar x}^2} + {{\bar M}_{\rm{e}}}\bar x = \\ \qquad \frac{{{{\rm{d}}^2}\bar w}}{{{\rm{d}}{{\bar x}^2}}} + \sum\limits_{m = 1}^\infty {{\tau ^{2m}}} \frac{{{{\rm{d}}^{2m + 2}}\bar w}}{{{\rm{d}}{{\bar x}^{2m + 2}}}} \end{array}$ | (11) |
其中,无量纲参数 $\bar {q}_0 = \dfrac{q_0 L^3}{EI_z }$,$\bar {M}_{\rm e} = \dfrac{M_{\rm e} L}{EI_z }$,$\bar {x} = \dfrac{x}{L}$,$\bar {w} = \dfrac{w}{L}$,$\tau = \dfrac{e_0 a}{L}$.
考虑简支梁两端铰支条件,有
$\bar {w}\left( 0 \right) = \bar {w}\left( 1 \right) = 0$ | (12) |
以往文献针对非局部高阶梁模型,往往采用截断法,仅取无穷级数的前几项即所谓非局部效应最强的几项做近似求 解[14, 15, 16, 17]. 理由是 $\tau $ 是一个小量,随着级数项$m$的增大,$\tau^{2m}$是一个可以忽略的高阶小量. 然而,尽管$\tau^{2m}$所带来的非局部尺度因子越来越小,但是其乘数${{\rm{d}}^{2m + 2}}\bar w/{\rm{ d}}{{\bar x}^{2m + 2}}$却不一定很小,并且无穷多项叠加之后,可能会产生预料不到的结果,其非局部效应未必可以忽略. 因此,本文将不会对非局部高阶梁控制方程做任何截断近似,而是直接求解无穷级数高阶方程. 根据扰动理论,假设方程(11)的解为级数形式$\bar w = {\rm{ }}\sum\limits_{n = 0}^\infty {{{\bar w}_n}{\tau ^n}} $, 采用逐阶求解法,将级数解代入式(11)可得
$\begin{array}{l} \frac{{{{\bar q}_0}}}{2}\bar x - \frac{{{{\bar q}_0}}}{2}{{\bar x}^{{\kern 1pt} 2}} + {{\bar M}_{\rm{e}}}\bar x = \\ \qquad {\tau ^0}\left( {\frac{{{{\rm{d}}^2}{{\bar w}_0}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 2}}}} + \tau \frac{{{{\rm{d}}^2}{{\bar w}_1}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 2}}}} + {\tau ^2}\frac{{{{\rm{d}}^2}{{\bar w}_2}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 2}}}} + {\tau ^3}\frac{{{{\rm{d}}^2}{{\bar w}_3}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 2}}}} + \cdots } \right) + \\ \qquad {\tau ^2}\left( {\frac{{{{\rm{d}}^4}{{\bar w}_0}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 4}}}} + \tau \frac{{{{\rm{d}}^4}{{\bar w}_1}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 4}}}} + {\tau ^2}\frac{{{{\rm{d}}^4}{{\bar w}_2}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 4}}}} + {\tau ^3}\frac{{{{\rm{d}}^4}{{\bar w}_3}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 4}}}} + \cdots } \right) + \\ \qquad {\tau ^4}\left( {\frac{{{{\rm{d}}^6}{{\bar w}_0}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 6}}}} + \tau \frac{{{{\rm{d}}^6}{{\bar w}_1}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 6}}}} + {\tau ^2}\frac{{{{\rm{d}}^6}{{\bar w}_2}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 6}}}} + {\tau ^3}\frac{{{{\rm{d}}^6}{{\bar w}_3}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 6}}}} + \cdots } \right) + \\ \qquad {\tau ^6}\left( {\frac{{{{\rm{d}}^8}{{\bar w}_0}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 8}}}} + \tau \frac{{{{\rm{d}}^8}{{\bar w}_1}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 8}}}} + {\tau ^2}\frac{{{{\rm{d}}^8}{{\bar w}_2}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 8}}}} + {\tau ^3}\frac{{{{\rm{d}}^8}{{\bar w}_3}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 8}}}} + \cdots } \right) + \cdots = \\ \qquad \sum\limits_{m = 0}^\infty {{\tau ^{2m}}} \left( {\sum\limits_{n = 0}^\infty {{\tau ^n}} \frac{{{{\rm{d}}^{2m + 2}}{{\bar w}_n}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 2m + 2}}}}} \right) \end{array}$ | (13) |
令 $\tau $ 的各次幂系数均为零,可得
$\left. \begin{array}{l} \frac{{{{\bar q}_0}}}{2}\bar x - \frac{{{{\bar q}_0}}}{2}{{\bar x}^{{\kern 1pt} 2}} + {{\bar M}_{\rm{e}}}\bar x = \frac{{{{\rm{d}}^2}{{\bar w}_0}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 2}}}}\\ \frac{{{{\rm{d}}^2}{{\bar w}_1}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 2}}}} = 0\\ \frac{{{{\rm{d}}^2}{{\bar w}_2}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 2}}}} + \frac{{{{\rm{d}}^4}{{\bar w}_0}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 4}}}} = 0\\ \frac{{{{\rm{d}}^2}{{\bar w}_3}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 2}}}} + \frac{{{{\rm{d}}^4}{{\bar w}_1}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 4}}}} = 0\\ \frac{{{{\rm{d}}^2}{{\bar w}_4}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 2}}}} + \frac{{{{\rm{d}}^4}{{\bar w}_2}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 4}}}} + \frac{{{{\rm{d}}^6}{{\bar w}_0}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 6}}}} = 0\\ \frac{{{{\rm{d}}^2}{{\bar w}_5}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 2}}}} + \frac{{{{\rm{d}}^4}{{\bar w}_3}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 4}}}} + \frac{{{{\rm{d}}^6}{{\bar w}_1}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 6}}}} = 0\\ \quad \quad \quad \quad \quad \vdots \\ \frac{{{{\rm{d}}^2}{{\bar w}_{2j - 1}}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 2}}}} + \frac{{{{\rm{d}}^4}{{\bar w}_{2j - 3}}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 4}}}} + \frac{{{{\rm{d}}^6}{{\bar w}_{2j - 5}}}}{{{\rm{d}}{{\bar x}^6}}} + \cdots + \\ \qquad \frac{{{{\rm{d}}^{2j}}{{\bar w}_1}}}{{{\rm{d}}{{\bar x}^{2j}}}} = 0\frac{{{{\rm{d}}^2}{{\bar w}_{2j}}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 2}}}} + \frac{{{{\rm{d}}^4}{{\bar w}_{2j - 2}}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 4}}}} + \frac{{{{\rm{d}}^6}{{\bar w}_{2j - 4}}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 6}}}} + \cdots + \\ \qquad \frac{{{{\rm{d}}^{2j + 2}}{{\bar w}_0}}}{{{\rm{d}}{{\bar x}^{{\kern 1pt} 2j + 2}}}} = 0 \end{array} \right\}$ | (14) |
其中$j =1$,2,3,$\cdots$.
从而解得
$\bar {w}_0 = \dfrac{\bar {q}_0 \bar {x}^{3}}{12} - \dfrac{\bar {q}_0 \bar {x}^{4}}{24} + \dfrac{\bar {M}_{\rm e} \bar {x}^{3}}{6} + C_0 \bar {x} + D_0$ | (15a) |
$\bar {w}_2 = \dfrac{\bar {q}_0 \bar {x}^{2}}{2} + C_2 \bar {x} + D_2 $ | (15b) |
$\bar {w}_j = C_j \bar {x} + D_j ,\ \ \left( {j = 1,3,4,5,\cdots } \right)$ | (15c) |
考虑边界条件(12),对级数试解中的各项应该同时满足$\bar {w}_n ( 0 ) = \bar {w}_n ( 1 ) = 0$,代入式(15)可得
$\bar {w}_0 = \dfrac{\bar {q}_0 \bar {x}^{3}}{12} - \dfrac{\bar {q}_0 \bar {x}^{4}}{24} + \dfrac{\bar {M}_{\rm e} \bar {x}^{3}}{6} - \dfrac{\bar {q}_0 + 4\bar {M}_{\rm e} }{24}\bar {x} $ | (16a) |
$\bar {w}_2 = \dfrac{\bar {q}_0 \bar {x}^{2}}{2} - \dfrac{\bar {q}_0 \bar {x}}{2} $ | (16b) |
$\bar {w}_j = 0 \ \ \ \left( {j = 1,3,4,5,\cdots } \right) $ | (16c) |
因此,非局部挠度表达式可推得为
$\begin{array}{l} \bar w = \frac{{{{\bar q}_0}{{\bar x}^{{\kern 1pt} 3}}}}{{12}} - \frac{{{{\bar q}_0}{{\bar x}^{{\kern 1pt} 4}}}}{{24}} + \frac{{{{\bar M}_{\rm{e}}}{{\bar x}^{{\kern 1pt} 3}}}}{6} - \\ \qquad \frac{{{{\bar q}_0} + 4{{\bar M}_{\rm{e}}}}}{{24}}\bar x + {\tau ^2}\left( {\frac{{{{\bar q}_0}{{\bar x}^{{\kern 1pt} 2}}}}{2} - \frac{{{{\bar q}_0}\bar x}}{2}} \right) \end{array}$ | (17) |
事实上,非局部挠度可仅由边界条件(12)确定,而不使用其演化的边界条件$\bar {w}_n ( 0 ) = \bar {w}_n ( 1 ) = 0$. 首先,直接写出带有未定系数的挠度表达式为
$\begin{array}{l} \bar w = \frac{{{{\bar q}_0}{{\bar x}^3}}}{{12}} - \frac{{{{\bar q}_0}{{\bar x}^4}}}{{24}} + \frac{{{{\bar M}_{\rm{e}}}{{\bar x}^3}}}{6} + {C_0}\bar x + \\ \qquad {D_0} + {\tau ^2}\left( {\frac{{{{\bar q}_0}{{\bar x}^2}}}{2} + {C_2}\bar x + {D_2}} \right) + \\ \qquad \tau \left( {{C_1}\bar x + {D_1}} \right) + \sum\limits_{n = 3}^\infty {{\tau ^n}} \left( {{C_n}\bar x + {D_n}} \right) \end{array}$ | (18) |
将边界条件(12)代入式(18),有
${D_0} + \sum\limits_{n = 1}^\infty {{\tau ^n}} \;{D_n} = 0$ | (19a) |
$\begin{array}{l} \frac{{{{\bar q}_0} + 4{{\bar M}_{\rm{e}}}}}{{24}} + {C_0} + {D_0} + \\ \qquad \frac{{{\tau ^2}{{\bar q}_0}}}{2} + {\rm{ }}\sum\limits_{n = 1}^\infty {{\tau ^n}} \;\left( {{C_n} + {D_n}} \right) = 0 \end{array}$ | (19b) |
以上两式相减可得
$\frac{{{{\bar q}_0} + 4{{\bar M}_{\rm{e}}}}}{{24}} + {C_0} + \frac{{{\tau ^2}{{\bar q}_0}}}{2} + {\rm{ }}\sum\limits_{n = 1}^\infty {{\tau ^n}} \;{C_n} = 0$ | (20) |
综合考虑式(19a)及式(20)可得
$\left. \begin{array}{l} {C_0} = - \frac{{{{\bar q}_0} + 4{{\bar M}_{\rm{e}}}}}{{24}}\\ {C_2} = - \frac{{{{\bar q}_0}}}{2}\\ {C_j} = 0\;\;\;\left( {j = 1,3,4,5,\cdots } \right)\\ {D_i} = 0\;\;\;\left( {i = 0,1,2,3,\cdots } \right) \end{array} \right\}$ | (21) |
因此,由式(18)及式(21)两式,同样可得挠度的非局部表达式,与式(17)完全一致.
以上方法适应于非局部高阶梁弯曲的任意形式,比如纳米悬臂梁在线性分布载荷且最大集度为$p_{0}$的作用下,如图 2所示, 此时高阶挠曲轴微分方程为
$\begin{array}{l} \frac{{{p_0}L}}{{2E{I_z}}}x - \frac{{{p_0}}}{{6E{I_z}L}}{x^3} - \frac{{{p_0}{L^2}}}{{3E{I_z}}} = \\ \qquad \frac{{{{\rm{d}}^2}w}}{{{\rm{d}}{x^2}}} + \sum\limits_{m = 1}^\infty {{{\left( {{e_0}a} \right)}^{2m}}} \;\frac{{{{\rm{d}}^{2m + 2}}w}}{{{\rm{d}}{x^{2m + 2}}}} \end{array}$ | (22) |
根据同样的求解方法,结合悬臂梁的边界条件,得到此时非局部挠度的无量纲表达式为
$\bar w = \frac{{10{{\bar p}_0}{{\bar x}^{{\kern 1pt} 3}} - {{\bar p}_0}{{\bar x}^{{\kern 1pt} 5}} - 20{{\bar p}_0}{{\bar x}^{{\kern 1pt} 2}} + 20{\tau ^2}{{\bar p}_0}{{\bar x}^{{\kern 1pt} 3}}}}{{120}}{\rm{ }}$ | (23) |
再如,图 3所示的固支梁承受均布载荷,其高阶挠曲轴微分方程为
$\begin{array}{l} - \frac{{{q_0}}}{{2E{I_z}}}{x^2} + \frac{{{q_0}L}}{{2E{I_z}}}x - \frac{{{q_0}{L^2}}}{{12E{I_z}}} = \\ \qquad \frac{{{{\rm{d}}^2}w}}{{{\rm{d}}{x^2}}} + \sum\limits_{m = 1}^\infty {\;{{\left( {{e_0}a} \right)}^{2m}}} \frac{{{{\rm{d}}^{2m + 2}}w}}{{{\rm{d}}{x^{2m + 2}}}} \end{array}$ | (24) |
经过正则摄动法计算发现,若要满足全部的4个边界条件(左右固定端挠度和转角均为0),必须使得$\bar {w}_j$ ($j =0$,1,2,3,$\cdots$)表达式中的全部待定系数$C_{j}=D_{j} =0$且 $\tau =0$,即
$\bar w = - \frac{{{{\bar q}_0}{{\bar x}^{{\kern 1pt} 4}}}}{{24}} + \frac{{{{\bar q}_0}{{\bar x}^{{\kern 1pt} 3}}}}{{12}} - \frac{{{{\bar q}_0}{{\bar x}^{{\kern 1pt} 2}}}}{{24}}{\rm{ }}$ | (25) |
此时无量纲的非局部挠度与经典挠度表达式一致,非局部效应消失.
3 数值算例与分析 3.1 算例非局部效应仅在纳米尺寸材料和结构中才可能显著体现,对于宏观尺度非局部效应可以忽略不计. 因此在数值算例中, 考虑矩形截面纳米梁,长为$L = 1\mu $m,宽和高分别为50nm与10nm,弹性模量$E=1$TPa. 以简支梁为例,设外载$q_{0}=0.01$N/m,$M_{\rm e}=2\times 10^{ - 14}$Nm,据此可求得无量纲参数$\bar {q}_0 = 2.4$,$\bar {M}_{\rm e} = 4.8$. 将上述无量纲参数代入式(17),可得非局部挠度的变化关系如图 4所示.
从图 4可见,非局部挠度随着无量纲化的轴坐标及非局部尺度因子的变化而变化,其中随着轴坐标的变化趋势与经典连续力学结果完全一致. 为了更加清晰地展示非局部尺度因子的作用效果,考察中点挠度,即在式(17)中令$\bar {x} = 1 /2$得
${\bar w_{{\rm{mid}}}} = - \frac{{5{{\bar q}_0}}}{{384}} - \frac{{{{\bar M}_{\rm{e}}}}}{{16}} - \frac{{{{\bar q}_0}}}{8}{\tau ^2}$ | (26) |
当$\tau =0$即忽略非局部效应时,中点挠度回归到经典解. 考虑非局部解与经典解之比,即
$\frac{{{{\bar w}_{{\rm{mid}}}}}}{{{{\bar w}_{{\rm{c - mid}}}}}} = 1 + \frac{{48{{\bar q}_0}{\tau ^2}}}{{5{{\bar q}_0} + 24{{\mathop {\bar M}\limits^\; }_{\rm{e}}}}}$ | (27) |
图 5中三角形数据点曲线给出了梁中点的非局部与经典挠度之比与尺度因子的关系,其中$\bar {q}_0 = 2.4$,$\bar {M}_{\rm e} = 4.8$. 可见随着非局部尺度因子的增大,即非局部效应的增强,挠度之比也得以提高. 因此,考虑非局部效应之后,此时非局部弯曲挠度高于相应的经典挠度,非局部挠度的提高意味着非局部效应降低了结构的弯曲刚度. 另一方面,从图 5可见,当非局部尺度因子从0增加到0.2时,非局部中点挠度比经典挠度约提高了3.62%;当尺度因子从0增加 到0.05时,非局部中点挠度仅提高了0.23%. 考虑到对具体材料而言,非局部参数即内特征尺度因子$e_{0}a$一般不变,若材料的外特征尺度$L$远大于内特征尺度 $e_{0}a$,即当$\tau =e_{0}a/L$较小时,非局部效应可以忽略不计. 当$\tau =0$即内特征尺度相比外特征尺度完全可以忽略时,非局部结果回归到宏观经典理论结果.
类似地,图 2所示悬臂纳米梁的自由端非局部挠度与其经典值之比,可由式(23)计算为
$\frac{{{{\bar w}_{{\rm{free}}}}}}{{{{\bar w}_{{\rm{c - free}}}}}} = 1 - \frac{{20{\tau ^2}}}{{11}}$ | (28) |
可见自由端挠度之比与所受的单一载荷$p_{0}$无关. 图 5中矩形数据点曲线给出了挠度之比与尺度因子的变化关系,这里显示了与上述三角形数据点曲线相反的结果. 考虑非局部效应的悬臂梁的自由端挠度减小,因此纳米结构的弯曲刚度提高.
以上数值算例表明,考虑非局部效应之后,纳米结构弯曲刚度或降低或提高或不变都有可能,这主要依赖于外载荷的具体形式 以及结构的边界约束条件,不能一概而论. 因此,两类刚度预测模型[6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21]不存在孰是孰非,关于该问题的学术争论也无需再进一步展开,形成两类研究 结论彼此矛盾的非局部梁模型的根本原因在于不同边界和载荷的自然属性. 物理上看,这与非局部理论建立在原子间长程力作用的基础上有关,即非局部理论认为整体不等于各部分之和,若将整体分割 成若干部分,各部分不能再保持整体所具备的原始性能,即各部分是``非局部的''. 所以,在此框架下边界约束与外载类型将会产生显著影响.
3.2 分析与讨论首先,考察梁的最大挠度所在的位置. 以图 2所示的纳米悬臂梁受线性分布载荷为例,其最大挠度可由下式决定
$6{\bar x_m} - \bar x_m^3 + 12{\tau ^2}{\bar x_m} - 8 = 0{\rm{ }}$ | (29) |
由式(29)可见,与经典连续理论不同,最大挠度所处的位置$\bar {x}_m $与非局部尺度因子$\tau $有关. 图 6给出了最大挠度所处位置随着非局部尺度因子的变化关系,可见当非局部作用十分强烈时(即外特征尺度极其微小或非局部尺度因子$\tau $很大时),最大挠度发生位置随着非局部效应的增大而逐步偏离自由端. 然而,对于一般纳米结构的非局部问题,$\tau $的取值不会很大,因此可以忽略非局部尺度因子对最大挠度发生位置的影响. 在本例中,当$0 < \tau \le 0.5$时,最大挠度的发生位置仍然在悬臂梁自由端,与经典连续理论即$\tau =0$时预测结果一致.
其次,分析边界条件相同、载荷不同时的结果. 相比于图 2所示的悬臂梁受线性分布载荷,这里再分析悬臂梁自由端受集中载荷的情形,也是文献中颇具争议的问题[6, 14]. 此时根据高阶梁的挠曲轴微分方程可得
$\frac{{{F_0}}}{{E{I_z}}}x + \frac{{{F_0}L}}{{E{I_z}}} = \frac{{{{\rm{d}}^2}w}}{{{\rm{d}}{x^2}}} + \sum\limits_{m = 1}^\infty {{{\left( {{e_0}a} \right)}^{2m}}} \;\frac{{{{\rm{d}}^{2m + 2}}w}}{{{\rm{d}}{x^{2m + 2}}}}$ | (30) |
其中$F_{0}$表示悬臂梁自由端所受的集中载荷. 根据同样的求解过程,可得非局部挠度中仅摄动解的第一项$\bar {w}_0 $为非零解,即无量纲的非局部挠度表达式为
$\bar w = \frac{{{{\bar F}_0}{{\bar x}^3}}}{6} - \frac{{{{\bar F}_0}{{\bar x}^2}}}{2}$ | (31) |
其中$\bar {F}_0 = {F_0 L^2}/ {EI_z }$. 可见悬臂梁的自由端承受集中载荷时,其非局部挠度解与经典解一致. 因此,在边界条件相同的情况下,不同的载荷形式对非局部效应有着不同的影响.
再次,分析载荷相同、边界条件不同时的结果. 同样相比于图 2所示的悬臂梁受线性分布载荷,这里考虑图 7所示的固支梁承受线性分布载荷,其高阶挠曲轴微分方程为
$\begin{array}{l} \frac{{3{p_0}L}}{{20E{I_z}}}x - \frac{{{p_0}}}{{6E{I_z}L}}{x^3} - \frac{{{p_0}{L^2}}}{{30E{I_z}}} = \\ \qquad \frac{{{{\rm{d}}^2}w}}{{{\rm{d}}{x^2}}} + \sum\limits_{m = 1}^\infty {{{\left( {{e_0}a} \right)}^{2m}}} \;\frac{{{{\rm{d}}^{2m + 2}}w}}{{{\rm{d}}{x^{2m + 2}}}} \end{array}$ | (32) |
类似于图 3的分析过程,此时欲满足全部的边界条件,非局部尺度因子必须取零,那么非局部挠度也退化至其经典解,即
$ \bar w = - \frac{{{{\bar p}_0}{{\bar x}^5}}}{{120}} + \frac{{{{\bar p}_0}{{\bar x}^3}}}{{40}} - \frac{{{{\bar p}_0}{{\bar x}^2}}}{{60}}$ | (33) |
此时挠度之比随非局部尺度因子的变化如图 5中圆形数据点曲线所示. 因此,在载荷形式相同的情况下,不同的边界条件对非 局部效应也有着不同的影响.
以上分析可见,纳米梁结构的非局部效应除了增强或减弱两种情况外,还可能不存在,该结论也解决了一类争议课题[6, 14]. 非局部效应不存在的情形包括悬臂梁自由端承受集中载荷,以及一些超静定结构,比如图 3及图 7所示的固支梁. 此外,其他超静定结构比如一端固定一端简支梁承受分布载荷等情形,经过计算发现其挠度的非局部效应亦消失.
值得一提的是,本文的所有讨论都是从非局部微分本构关系式(2)出发的. 针对原始的积分本构式(1),通过研究非局部权 函数$\alpha $的性质,艾林根提出了不同维度下的权函数表达式,基于这些权函数表达,非局部理论下的平面波的散射曲线预测结果与原子晶 格动力学的预测结果非常吻合,最大误差不超过1.2%. 特别地,当权函数取为线性微分算子的格林函数时,便可推得非局部微分本构表达[2]. 因为模型求解上的方便,目前在纳米力学的研究中基本都采用非局部微分本构[3, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21];而在带有两平行裂纹的功能 梯度材料的非局部理论解[4],以及能量的非局部新模型及基于此的新的应变梯度理论[5]等研究中,采用了积分本构 表达式. 因此,微分本构式(2)的建立是基于积分本构 式(1),尽管它们形式上不同并且存在一定程度上的简化和近似,但是可以预测,分别从积分本构与微分本构出发,所得出的定性研 究结论是一致的.
4 结论本文以非局部弹性理论应用于纳米结构性能预测上的不一致性为出发点,采用相对简洁的思路,利用迭代法求解了非局部梁的 高阶应力,并进一步推导了梯度型非局部高阶梁模型的控制方程. 梁的挠曲轴微分方程表达为无穷级数形式,其中包含了非局部挠度对轴坐标的无穷高阶导数. 在梁的高阶弯曲方程中,采用正则摄动法求解了非局部弯曲挠度的级数解. 结果表明:随着非局部尺度因子的增大亦即非局部效应的增强,纳米梁的弯曲挠度相比于其经典解表现出增大、减小、不变3种 趋势,具体结果依赖于纳米梁的边界条件以及外载荷的形式. 因此,基于梯度型非局部微分本构模型,计及非局部效应对纳米结构弯曲刚度的影响程度为:可能降低也可能提高还可能维持不变. 在非局部效应影响弯曲挠度的算例中,非局部挠度最大值发生位置在一定程度上也受到非局部效应的影响. 本文研究结果解决了非局部弹性理论应用于纳米力学的现存争议,以相对简单的思路和途径,从理论上解释了现有不同研究结论的成因,为纳米技术的快速进步奠定理论基础.
1 | Eringen AC, Edelen DGB. On nonlocal elasticity. International Journal of Engineering Science, 1972, 10: 233-248 |
2 | Eringen AC. On di erential equations of nonlocal elasticity and solutions of screw dislocation and surface waves. Journal of Applied Physics, 1983, 54: 4703-4710 |
3 | 郑长良. 非局部弹性直杆振动特性及Eringen 常数的一个上限. 力学学报, 2005, 37(6): 796-798 (Zheng Changliang. The free vibration characteristics of nonlocal continuum bar and an upper bound of material constant in Eringen's nonlocal model. Chinese Journal of Theoretical and Applied Mechanics, 2005, 37(6): 796-798 (in Chinese)) |
4 | Liang J,Wu SP, Du SY. The nonlocal solution of two parallel cracks in functionally graded materials subjected to harmonic anti-plane shear waves. Acta Mechanica Sinica, 2007, 23(4): 427-435 |
5 | 易大可, 王自强. 能量非局部模型和新的应变梯度理论. 力学学 报, 2009, 41(1): 60-66 (Yi Dake, Wang Tzuchiang. Energy nonlocal model and new strain gradient theory. Chinese Journal of Theoretical and Applied Mechanics, 2009, 41(1): 60-66 (in Chinese)) |
6 | Peddieson J, Buchanan GG, McNitt RP. Application of nonlocal continuum models to nanotechnology. International Journal of Engineering Science, 2003, 41: 305-312 |
7 | Wang Q, Varadan VK. Vibration of carbon nanotubes studied using nonlocal continuum mechanics. Smart Materials & Structures,2006, 15: 659-666 |
8 | Lu P, Lee HP, Lu C, et al. Application of nonlocal beam models for carbon nanotubes. International Journal of Solids and Structures,2007, 44: 5289-5300 |
9 | Hu YG, Liew KM, Wang Q, et al. Nonlocal shell model for elastic wave propagation in single- and double-walled carbon nanotubes. Journal of the Mechanics and Physics of Solids, 2008, 56: 3475-3485 |
10 | Kiani K. Vibration analysis of elastically restrained double-walled carbon nanotubes on elastic foundation subjected to axial load using nonlocal shear deformable beam theories. International Journal of Mechanical Sciences, 2013, 68: 16-34 |
11 | Li YS, Cai ZY, Shi SY. Buckling and free vibration of magnetoelectroelastic nanoplate based on nonlocal theory. Composite Structures,2014, 111: 522-529 |
12 | Ke LL, Wang YS, Yang J, et al. Free vibration of size-dependent magneto-electro-elastic nanoplates based on the nonlocal theory. Acta Mechanica Sinica, 2014, 30(4): 516-525 |
13 | Hosseini-Hashemi S, Kermajani M, Nazemnezhad R. An analytical study on the buckling and free vibration of rectangular nanoplates using nonlocal third-order shear deformation plate theory. European Journal of Mechanics—A/Solids, 2015, 51: 29-34 |
14 | Lim CW. On the truth of nanoscale for nanobeams based on nonlocal elastic stress field theory: equilibrium, governing equation and static deflection. Applied Mathematics and Mechanics, 2010, 31:37-54 |
15 | Lim CW, Yang Y. Wave propagation in carbon nanotubes: Nonlocal elasticity-induced sti ness and velocity enhancement e ects. Journal of Mechanics of Materials and Structures, 2010, 5: 459-476 |
16 | Lim CW, Yang Q. Nonlocal thermal-elasticity for nanobeam deformation: exact solutions with sti ness enhancement e ects. Journal of Applied Physics, 2011, 110: 013514 |
17 | Wang L. A modified nonlocal beam model for vibration and stability of nanotubes conveying fluid. Physica E, 2011, 44: 25-28 |
18 | Yang Q, Lim CW. Thermal e ects on buckling of shear deformable nanocolumns with von Kármán nonlinearity based on nonlocal stress theory. Nonlinear Analysis: Real World Applications, 2012,13: 905-922 |
19 | Lim CW, Xu R. Analytical solutions for coupled tension-bending of nanobeam-columns considering nonlocal size e ects. Acta Mechanica,2012, 223: 789-809 |
20 | Yu YM, Lim CW. Nonlinear constitutive model for axisymmetric bending of annular graphene-like nanoplate with gradient elasticity enhancement e ects. Journal of Engineering Mechanics-ASCE,2013, 139: 1025-1035 |
21 | Islam ZM, Jia P, Lim CW. Torsional wave propagation and vibration of circular nanostructures based on nonlocal elasticity theory. International Journal of Applied Mechanics, 2014, 6: 1450011 |