基于简化的应变梯度理论下Kirchhoff板模型边值问题的提法及其应用*

2022-04-27徐晓建邓子辰

应用数学和力学 2022年4期

徐晓建，邓子辰

（1.长安大学公路学院特殊地区公路工程教育部重点实验室，西安 710064；2.西北工业大学工程力学系，西安 710072）

引言

板结构广泛应用于航空航天、土木、机械等领域，在几乎所有工程科学中广泛存在.因此，对板结构的经典理论力学的研究一直是前沿课题.然而，随着科技的快速发展，尤其是纳米技术的革新，人们发现，越来越多工程结构的力学行为用经典结构力学理论预测会得到很大的误差[1-2].例如，Lam 等[3]的微梁弯曲实验结果表明，梁的归一化刚度与其厚度有密切联系；而传统梁模型的归一化刚度却与梁厚度无关.Chen 等[4]研究表明，ZnO 结构的弹性模量具有尺寸效应.随着研究的进一步深入，人们发现，材料所表现出的尺寸效应随着材料/结构尺寸的减小逐渐明显[5-6].

因此，如何建立能够准确表征结构/材料尺寸效应的力学模型，近年来一直是力学工作者研究的关键科学问题和热点问题[7].基于此，Wang 等[8-9]采用分子动力学和力学模型研究发现，碳纳米管具有尺寸效应的波动行为，可由应变梯度模型较好地进行预测.Ansari 等[10]研究了薄板振动问题，并给出了非经典板模型与分子动力学结果吻合得很好的结论.Khakalo 等[11]研究了具有三维三角形结构的蜂窝板结构的力学模型，并利用均匀化方法给出了非经典材料参数.王碧蓉等[12]通过修正剪力非局部参数，研究了双臂碳纳米管的波传播特性.

近年来，研究者们广泛关注板挠度的解析解及数值解的研究[13-14].然而，基于简化的应变梯度理论的薄板模型的研究还较少.为弥补这一不足，本文提出一种新型薄板理论，并构建其边值问题.以弹性地基板为例，研究了尺寸效应参数对薄板挠度解及自由振动频率的影响.

1 应变梯度理论及Kirchhoff 板模型

1.1 应变梯度理论

根据简化的应变梯度理论[15-16]，体积为V、平面域为 Ω的线弹性体的应变能为

式中

分别为应变张量和应变梯度张量（不失一般性，本文中所用符号均符合张量的符号约定）；高阶应力张量 τi jk为

其中，l2为表征由于应变梯度引起的材料尺度效应的参数.

根据文献[17]，外力功为

系统的总动能包括经典动能和速度梯度引起的动能，即

式中，ρ，l1和分别为体密度、表征速度梯度的材料参数和速度梯度.

系统的总势能由动能、应变能和外力功三个部分组成，其表达式为

1.2 薄板方程及边界条件

考虑一各向同性的线弹性均质材料，厚度为h，在xOy平面内所围成的平面域为Ω，该域边界为分段光滑曲线Γ，切线方向为s，外法线方向为n，如图1所示.规定n到s的转向与x轴到y轴的转向相同时为正，并且n和坐标x的夹角为θ.假定薄板承受一横向载荷p3、边界弯矩、高阶边界弯矩和边界力的作用，横向位移为w(x,y).

图1 板边界及荷载Fig.1 Boundary conditions and loadings

由薄板的假定，得以下几何关系：

本文中，下标中的希腊字母仅取x,y.

在平面应力状态下的各向同性薄板的本构方程为

式中，E，µ和δ 分别为板的弹性模量、Poisson 比和δ 函数.

令Rαβ和Mαβ分别为厚度引起的单位宽度广义力和经典弯矩，定义如下物理量：

则以位移w表示的弯矩为

由式(9)和(10)可知，由厚度引起的广义力与弯矩存在如下关系：

在平面应力状态下，可以忽略所有带有z指标的应力，如 σαz=0.但由式(3)和(8)可知，ηαβz=εαβ,z=−w,αβ≠0.因此，文献[18]简单地忽略ηαβz对应变能的影响，会产生较大的误差.

根据式(1)、(3)和(4)，应变能的一阶变分为

结合式(8)和式(13)，并考虑到式(10)中弯矩的定义，可以改写式(13).随后对所得方程的指标γ 进行一次分部积分后，式(13)最终可改写为如下三个积分之和：

将指标依次代入，式(14)中I1可展开为如下分量形式：

如无特殊说明，本文公式中的偏微分算子 L2仅作用于弯矩，且L2(·)=1−(·),γγ.

式(15)可通过使用两次Gauss 散度定理，最后结果为

令任意函数f1=f1(x,y)和f2=f2(x,y)，将Stokes 公式

代入式(16)，可将面积分转化为线积分.随后，利用附录A 并经过繁琐的推导，得以下线积分表达式：

其中

式中，Qn1，Mnn1和Mns分别为仅考虑x,y影响的非经典剪力、法向弯矩和切向弯矩.

利用式(11)，式(18)～(20)可以改写为直角坐标系下xOy的结果.随后，利用附录A，式(18)～(20)可进一步改写成局部坐标系(n, s)的形式，结果为

式中，[Mns]k为第k个角点的弯矩.

相似地，通过繁琐的推导，式(14)中I3为

式中

式（23）的第一个等号利用了Stokes 公式，第二个等号利用了附录A 中的转换关系和两次对s的分部积分.

结合式(17)、(23)和I2，应变能式(14)的一阶变分为

式中

分别为非经典总等效剪力和总法向弯矩.

利用式(21)和(24)，式(26)和(27)可以改写为如下位移形式：

由式(25)，可知外力虚功的形式为

考虑到式(6)，可得动能在时间间隔[t0,t1]的一阶变分为

式中，L1(·)=1(·),γγ.值得一提的是，当薄板的位形在初始时刻t=t0和最终时刻t=t1给定时，式(31)最后结果中等号右端第二项为零.

利用Hamilton 原理，有

将式(25)、(30)和(31)代入到式(32)并考虑到式(11)，得承受横向分布载荷下的薄板控制方程为

边界条件为

显然，当应变梯度参数和速度梯度参数都取零时，式(33)退化为经典薄板的控制方程，式(34)和(35)退化为经典薄板的边界条件，式(36)将自动满足.对于一般梯度薄板来讲，式(34)～(36)给出了光滑平面域内的非经典边界条件.当不考虑z方向的微分（即令I2=0）和速度梯度影响时，式(33)退化成文献[18]的结果.

本小节将讨论工程中常见非经典薄板的边界条件：

(ⅰ)对于周边固支圆板，其周向挠度和法向转角为零.因此，式(34)给出的经典边界条件为w=0，式(35)给出的经典边界条件为w,n=0.由于边界高阶弯矩不为零，式(36)给出的非经典高阶边界条件为w,nn=0.

(ⅱ)对于周边简支圆板，其周向挠度和弯矩为零.因此，式(34)给出的经典边界条件为w=0，式(35)给出的经典边界条件为==0.由于边界存在高阶弯矩或法向曲率不为零的情况，式(36)给出的非经典高阶边界条件为Mnnn=0或w,nn=0.

对于矩形薄板，以边界线法向为x轴正向为例，此时有 θ =0，且(n,s)与(x,y)坐标系重合.这表明，可以简单地把边界条件(34)～(36)中的n,s分别替换为x,y.因此，3 种常见薄板的边界条件列于表1中.

表1 矩形薄板3 种常见的边界条件Table 1 Three common boundary conditions (BCs)for a rectangular plate

1.3 角点条件

本小节将考虑板边界为分段光滑边界的角点问题.在推导边界条件时，角点条件由以下五部分组成.

第一部分由式(17)等号右边最后一项产生，为

第二部分由I2产生，为

第三部分由式(23)第一个等号右边第一项产生，为

第四部分由式(23)第一个等号右边第二项产生，为

第五部分由式(23)第一个等号右边第三项产生，为

叠加以上五式，得以位移表示的角点条件为

值得一提的是，转化为本文符号后，文献[18]通过变分原理导出的角点条件为

显然，文献[18]中把 δw,s项遗漏了.与文献[18]不同，本文的创新性体现在以下几个方面：1)考虑了速度梯度影响，即在总动能式(6)中考虑了含有l1的微分项；2)考虑了沿厚度方向的微分对板有效抗弯刚度的贡献，即考虑了由厚度引起的广义力方程(12)对板刚度的影响.

2 弹性地基上的周边简支矩形板

考虑一放置于弹性地基上的均质等厚度各向同性薄板，其长宽高分别为a,b和h，如图2所示.弹性地基符合Winkler 地基模型，则地基反力为-Kw(x,y)，K为地基刚度系数.

图2 弹性地基上的周边简支矩形板Fig.2 A fully simply supported rectangular plate resting on an elastic foundation

2.1 静位移分析

首先考虑一承受横向均布荷载为p的四边简支薄板的静位移问题.此时板的控制微分方程(33)改写为

挠度可用双三角级数表示为

式中，wmn为待求参数，λm=mπ/a，λn=nπ/b.

均布荷载p可用双三角级数表示为

将式(45)和(46)同时代入式(44)，可得wmn的表达式.再将该表达式回代入式(45)，得挠度解为

当l2不存在时，上式可退化成经典文献[19]的结果.

2.2 自由振动分析

其次考虑一周边简支薄板的自由振动问题.此时式(33)可改写为以下频域方程：

将式(45)代入式(48)，可得板固有频率解为

如令l1=l2=0，即不考虑应变梯度和速度梯度参数的影响时，式(49)可退化为经典结果[19].

3 结果和讨论

为研究梯度参数和地基参数对板结构力学行为的影响，本节分别选取板的静位移和固有频率进行研究.根据经典板理论，引入以下无量纲化位移及无量纲化固有频率参数：

如无特殊说明，板为正方形板，边长为a，且计算参数为

3.1 有效性验证

当不考虑地基影响时，Ansari 等[10]利用分子动力学（MD）的方法研究了zigzag 型方形石墨烯的自由振动频率问题.为验证本文模型的有效性，本小节将对比他们的结果.文献[10]中板的计算参数为

式中，D，ρ，h分别为板的抗弯刚度、体密度和厚度.

将以上材料参数代入到式(49)，并令m=n=1,K=0，得到拟合函数.利用拟合函数拟合文献[10]中的数据，得尺寸效应参数为

需要说明的是，式(52)中的材料参数仅适用于文献[10]中拟合石墨烯MD 的结果.对于碳纳米管等结构来说，基于应变梯度和速度梯度理论的材料参数l1和l2也可通过类似的方法通过拟合得到，如文献[20]给出的拟合结果为l2=0.035 5 nm，l1/l2=2.5.

不失一般性，后文采用式(52)的材料参数计算.

图3给出了基频随方形板边长的变化关系.可以看到，本文结果与MD 方法得到的结果吻合得很好，而经典薄板解却不能很好地预测板振动频率随其边长的变化，尤其当边长a<15 nm 时.同时，当不考虑速度梯度影响时（即令l1=0），得到的基频频率与文献[10]中的结果差异较大.这些间接地证明了本文模型的精确性及在工程实践中预测微纳米板结构力学行为的必要性.