中厚板理论的适用范围和精确程度的研究＊

龙述尧，姜 琛

（湖南大学 汽车车身先进设计制造国家重点实验室，湖南 长沙 410082）

板和壳在生产建设中被广泛使用，随着科学技术的发展，在某些领域中，除了传统的薄板（壳）外，还使用了各种各样的特殊板壳结构.例如，气冷核反应堆的预应力混凝土压力容器，航空航天中的各种夹层板（壳）和复合材料板（壳）等.在板壳问题中，工程上广泛使用的是Kirchhoff－Love薄板理论，Re－issner中厚板理论[1]和 Mindlin中厚板理论[2].在工程实际问题中，由于结构的几何形状、边界条件和荷载的复杂性，使用解析方法求解几乎不可能，大多使用各种数值方法进行分析，其中广泛使用的是有限元法.但是有限元法需域内划分单元、应力解精度差、计算成本高，所以在20世纪70年代，计算机发展未能满足有限元法的需求时，许多研究者开始研究边界积分方程法 （Boundary Integration Equation Method，BIEM），提出了一种只需在求解域边界上划分单元的新方法，由 C.A.Brebbia[3]首先命名为边界元法 （Boundary Element Method，BEM）.之后，边界元法被应用到各个领域，其中F.Vander Weeёn[4]和 S.Y.Long，C.A.Brebbia，J.C.F.Telles[5]研究了边界元法在Reissner理论中厚板中的应用.

本文采用Reissner中厚板理论边界元程序、Mindlin中厚板理论有限元程序和三维弹性力学有限元程序计算各种厚跨比下四边固支和四边简支矩形板，并进行对比，来说明三维理论、中厚板理论和薄板理论分析板时厚跨比的适用范围以及精确程度，为采用何种理论来分析板提供理论依据和参考；同时也验证了采用缩减积分的Mindlin中厚板有限元法消除了剪切锁死现象.

1 Reissner中厚板理论

Reissner中厚板理论是考虑了横向剪切变形和挤压变形的理论.考虑笛卡尔坐标系（x1，x2，x3），其中x1，x2轴在板中面内，x3则垂直于板中面，如图1所示.在板的上下表面切应力为零而正应力为分布载荷集度，即x3＝±h／2，σ33＝±q／2.

图1 中厚板几何形状和坐标系Fig.1 Geometry shape and coordinate of moderately thick plate

各应力在厚度方向上的变化关系[1]如下：

式中：i，j＝1，2；Mij为单位长度弯矩；Qi为剪力.为了求解Reissner中厚板控制方程的基本解，在此假定体力在厚度方向上也按应力一样变化，即

由以上假定，Reissner中厚板的控制方程变为：

其中的弯矩、扭矩和剪力用转角和位移表示为：

式中：D＝Eh3／[12（1－ν2）]为板的弯曲刚度；ν为泊松比；是厚度参数.式（4）中转角和位移为沿板厚的平均值，即：

式中：vi和v3为坐标轴方向的位移分量.

为了简便，将βi（i＝1，2）和w用ui（i＝1，2，3）表示，则对于区域Ω的边界Γ，位移和广义面力的边界条件可以表示为：

其中S＝Su∪St，ti为广义面力.假定边界上某点外法线方向余弦为ni＝cos（n，xi），则广义面力为

2 Reissner中厚板边界元基本解和边界积分方程

基本解可由式（3）中体力项定义为一块无限大中厚板作用单位横向集中力或单位弯矩求得，具体推导见文献[3].仅列出关于位移和面力的基本解：

上两式中：

其中K0（z）和K1（z）为修正贝塞尔函数.

使用Betti功互等定理、格林函数或加权残值法，可以得到Reissner中厚板的边界积分方程：

式中：cij是一个与边界形状有关的系数，对于光滑边界cij＝δij／2，对于有角点和尖点的边界法线不连续的边界，可以解析计算cij，但更方便的是间接计算cij的值.当整个边界无外力和载荷作用时，式（9）变为：

上式的非零解为板的3种独立刚体位移的任意组合，把这些位移代入式（10）有

当均布荷载情况下，可以应用散度定理把式（9）中的域积分简化为边界积分.

数值方法求解式（9）后，未知边界位移和边界力可以全部求出.令cij＝δij可以求出内部点的位移，而内部点的内力解，需要对式（9）关于源点ξ求微分，有如下形式：

式中：i，j≠3.P的分量为Pij＝Mij，P3j＝Qj.其中wij，Dijk和Sijk见文献[6]中附录Ⅰ.

3 数值计算

把所考虑的板的边界划分为M个单元，以每个节点为源点并在M个单元上积分，离散式（9），并将其中域积分转换为边界积分，则可得线性方程组

式中：H，G为影响矩阵；D为对应于载荷q的边界积分列向量.代入已知的边界条件后，将其中的已知项移到方程等号右边与D向量合并，未知项移到等号左边，得到标准代数方程组

采用标准高斯消去法求解上述方程组，即可求得未知量.之后可以利用式（9）计算内部任意一点的位移，利用式（12）求内部点的内力.

由于式中的奇异性，在计算奇异积分时，采用了三次多项式非线性变换[7].此种变换可以在奇异点附近集中较多的高斯点，因此对奇异积分采用标准高斯积分能达到相当高的精度.

对于存在角点和尖点等边界法线不连续的情况，采用双节点或部分不连续单元.双节点单元中2个节点有不能在同一方向上给定位移的限制，而不连续单元没有这种限制.

根据前述离散化的公式，编制了Reissner中厚板边界元的FORTRAN程序[8].

4 算 例

采用了2个具有代表性的算例来研究薄板和中厚板理论适用范围和精确程度.算例为：1）受均布荷载作用的四边固支矩形板，如图2所示；2）受均布荷载作用的四边简支矩形板.两算例的长、宽取a＝b＝1，横向荷载集度q＝1，板的弯曲刚度D＝1，泊松比ν＝0.3.

4.1 受均布荷载作用的四边固支矩形板

利用本文编制的Reissner中厚板边界元程序计算了在不同厚跨比（h／a＝0.001，0.005，0.01，0.05，0.1，0.2，0.5）下的板中心点的挠度w，板中心点下表面应力σx和板固支边中点的弯矩Mx.并与利用有限元软件ABAQUS中对应的三维有限元解和板壳有限元解，还有部分薄板解析解进行对比.

利用ABAQUS的三维8节点单元、二维4节点板壳单元和2节点边界单元进行分析计算，其中以三维8节点有限元分析的结果作为基准，以说明二维4节点板壳单元和2节点边界单元计算结果的准确性.

图2 四边固支矩形板Fig.2 Four edges clamped rectangular plate

在三维有限元、板壳有限元和中厚板边界元计算中，取中面中心点的挠度和中心点板的下表面点的应力σx作为比较对象.因为三维单元只能得到节点的位移和应力，要想得到截面上的内力时需要进行合成计算，为此取中心点板的下表面点的应力σx作为比较对象.板壳有限元分析结果既可输出内力，也可输出应力，但由于边界元只能输出内力，所以需要利用弯矩Mx与应力σx的关系求得：

式中：z是应力计算点到中面的距离.

为了更好、更精确地对比各方法的计算结果，对Reissner中厚板边界元模型每边划分16个线性单元，整个板共64个单元；对板壳有限元模型同样每边划分16个线性单元，共16×16个单元；三维有限元模型则在板面内每边划分16个线性单元，厚度方向上网格的划分是依照三维单元尽量要成正六面体的准则来进行的，这样能保证作为基准的三维实体有限元结果的精度.

图3画出了各种厚跨比下板中点的挠度.从图3中可以看出在所有厚跨比情况下，3种方法的计算结果十分吻合，说明了考虑剪切变形的Reissner中厚板理论是对三维弹性力学的较精确的简化.从图中还可以看出，当厚跨比大于0.2时，与三维单元的结果误差明显增大.

图3 各厚跨比情况下四边固支板中点挠度Fig.3 Central deflections of clamped plate for several thickness－to－span ratios

从图4可看到，在厚跨比较小（h／a＜0.2）时，边界元、板壳有限元和三维有限元计算出的σx十分接近，但是在h／a＞0.2后，误差越来越大.但是边界元解精度比板壳有限元的解精度高，证明了边界元法在理论上比有限元法应力（内力）的精度高.综合各种厚跨比挠度和应力的误差，中厚板与厚板的厚跨比分界线一般取0.2为宜，厚跨比超过了0.2，无论是挠度还是应力的误差都会急剧增加.同时，从图3和图4还可以看出：利用缩减积分技术的Mindlin板壳有限元，在计算非常薄的板时大大缓解了剪切锁死现象.

图4 各厚跨比情况下四边固支板中点下表面的σx误差Fig.4 Central errors ofσxof clamped plate for several thickness－to－span ratios

图5还给出了固支边中点处的弯矩Mx.在计算弯矩时，三维实体单元选取20节点二次单元，因为二次单元的位移是二次分布的，但是应力是线性分布的，所以二次单元积分出的弯矩应该更精确.从图5可以看出边界元解比板壳有限元的解精确得多，而且在薄板阶段和解析解十分靠近，充分说明了中厚板边界元解对薄板有很好的精度.边界元解和有限元三维实体单元解有差异，应当是来源于弯矩计算时积分的误差，因为三维单元只能输出应力，想获得弯矩还需对应力进行积分，而节点并不是在高斯积分点上，只能使用梯形法近似积分，致使精度较差.

图5 各厚跨比下四边固支板边中点处的MxFig.5 CentralMxof clamped plate for several thikness－to－span ratios

4.2 受均布载荷作用的四边简支矩形板

如图6所示四边简支矩形板，所有数据同例4.1.

图6 四边简支矩形板Fig.6 Simple－supported rectangular plate

板壳有限元仍然采用双线性板壳单元，而中厚板边界元采用2节点线性单元.另外单元划分准则也跟算例4.1相同，边界元在每边上划分16个单元，整板共64个单元；板壳有限元也在每边划分16个单元.

图7 各厚跨比下四边简支板中点挠度Fig.7 Central deflection of simple－supported plate for several thickness－to－span ratios

从图7中可以看出，3种单元的四边简支板中心的挠度的解十分接近.在厚板阶段（h／a＞0.1），图中清楚地表明了本文所用中厚板边界元法比板壳有限元法更接近三维有限元的解，说明了边界元法精度比有限元高.

图8 各厚跨比下四边简支板中心处弯矩的误差Fig.8 Central errors ofMxof simple－supported plate for several thickness－to－span ratios

从图8中可以看出，无论是边界元还是板壳有限元的解都十分接近三维有限元的解.从图中可以看到，在薄板阶段（h／a＜0.1），边界元解和薄板理论十分接近，和三维有限元解、板壳有限元解有些许误差，因为板壳有限元解的内力精度要低于边界元解.

5 结 论

本文将Reissner中厚板理论运用到边界元法中，并通过数值方法将其编制成FORTRAN程序，通过几个典型算例的计算，证明了边界元法分析中厚板的程序的正确性.另外，通过与三维有限元、板壳有限元的对比，证明了Reissner理论是对三维弹性力学的正确简化，同时，中厚板边界元法比板壳有限元法在挠度、内力、应力上具有更高的精度，更高的计算效率.因此，使用边界元方法对中厚板进行数值计算，是十分合适和有效的.

另外，从两个算例可以得到结论，薄板和中厚板的分界线一般取厚跨比等于0.1，此后薄板理论解和中厚板理论解的误差越来越大；而中厚板和厚板的分界线应当在厚跨比等于0.2处，因为此后中厚板理论解均和三维弹性力学解有相当大的误差.同时说明了中厚板理论适用于薄板情况，并验证了利用缩减积分技术的Mindlin板壳有限元，在计算非常薄的板时大大缓解了剪切锁死现象.

[1] REISSNER E.The effect of transverse shear deformation on the bending of elastic plate[J].J Appl Mech，1945，12：69－77.

[2] MINDLIN R D.Influence of rotatory inertia and shear on flexural motions of isotropic elastic plates[J].Journal of Applied Mechanics，1951，18：31－38.

[3] BREBBIA C A.The boundary element method for engineers[M].London：Pentech Press，1978：159－211.

[4] VANDER WEEËN F.Application of the boundary integral equation method to Reissner's plate model[J].International Journal for Numerical Methods in Engineering，1982，18（1）：1－10.

[5] LONG S Y，BREBBIA C A，TELLS J C F.Boundary element bending analysis of moderately thick plates[J].Engineering Analysis，1988，5（2）：64－74.

[6] KARAM V J，TELLES J C F.On boundary elements for Reissner’s plate theory[J].Engineering Analysis，1988，5（1）：21－28.

[7] TELLES J C F.A self－adaptive coordinate transformation for efficient numerical evaluation of general boundary element integral[J].Int J Numer Engineering，1987，24：959－973.

[8] 龙述尧.边界单元法概论[M].北京：中国科学文化出版社，2002：214－253.LONG Shu－yao.Introduction of boundary elements method[M].Beijing：China Scientific ＆ Cultural Publishing House，2002：214－253.（In Chinese）