高阶剪切和法向变形板理论在功能梯度板中的应用

高阶剪切和法向变形板理论在功能梯度板中的应用
摘要：将高阶剪切变形和法向变形板理论应用于有限元理论中，位移沿板厚度方向并展开为勒让德多项式级数，将其系数作为广义位移自由度，并带入控制方程和能量方程，导出有限元计算格式。

采用Fortran语言将其编程实现，通过不同材料梯度参数功能梯度板数值计算算例与精确解的对比，两者都可以较好的吻合。

沿界面的挠度和应力计算结果的分布，说明了一般的板理论的假设并不适合功能梯度板的分析，这种高阶板单元相对更适用于功能梯度板的分析与计算。

关键词：功能梯度材料板，高阶剪切变形和法向变形板理论，有限单元法，Mori-Tanaka模型
引言：
功能梯度材料是一种新型复合材料，是由两相或两相以上颗粒复合材料合成且沿着一个或多个方向每种组分不断变化的预制构件，具有变化的微观结构和连续变化的力学或热学性能。

传统的复合材料是分层离散的材料组合，并伴随着界面的产生和材料热力、机械和物理特性的突变。

由于非完全粘结、残余应力等因素，通常这些界面都包含一定的缺陷。

即使没有，这种材料性质的不匹配会导致应力集中，使得界面变成一个在正常工作状态下容易开裂、脱粘和层裂的温床。

而功能梯度材料则根据具体要求，选择使用两种或两种以上具有不同性能的材料，通过连续的改变材料的组成和结构，使其内部界面消失，从而得到功能相应于组成和结构的变化而渐变的非均质材料，以减小和克服结合部位的性能不匹配因素[1]。

目前，国内外已有大量对于功能梯度材料结构的研究工作。

由于功能梯度材料的材料参数变化的，所以问题的控制方程是变系数偏微分方程，即使对于简单问题也难以求得解析解。

现有的研究主要的分析方法有直接将考虑剪切变形的的板壳理论套用到功能梯度材料结构当中[2-4]；层合模型法[5-7]，将功能梯度材料结构沿材料变化方向分为若干层，每一层近似看做均质材料，然后加上界面处的连续条件；渐进解法[8,9]等。

而对于数值计算方法而言，如有限元法，由于材料不均匀，利用常规实体元对功能梯度材料进行三维分析，需要划分大量的单元。

所以针对这种材料特点，需要新的有限单元与计算格式，有人提出了梯度有限元方法[10,11]。

本文作者将采用同时考虑了剪切变形和层间变形的高阶板理论(a higher order shear and normal deformable plate theory)[12]构造了一种新的有限板单元，并该单元对功能梯度板进行分析，并通过算例说明了这种单元的可行性。

HOSNDPT板理论与控制方程的建立
如图1所示，考虑一材料参数沿板厚Z方向变化的功能梯度板。

其中未变形的板的在直角坐标系中的域为。

板的中面为坐标平面OXY记为S，沿x, y, z 轴的位移分别记为u, v, w。

图1
将位移场函数展开为的级数：
其中为定义在上单位正交化的次勒让德多项式，即满足：
为狄拉克函数。

因为勒让德多项式为正交完备化的函数系，所以当时，该级数收敛于位移场函数。

在小变形的情况下，根据应变与位移的关系可得：
是关于的次多项式。

可以写成：
有限元计算格式的推导
功能梯度板所占的域记为，四边边界面记为，上下表面分别记为与。

板所受的体力为，在上的位移边界条件为：，在上的应力边界条件为，在上下表面受面力。

问题的静力平衡方程可以写成：
由板的受力条件与自然边界条件，可求的板的总位能：
设为形函数，其中为单元节点数，令
其中为单元广义节点位移，则单元内任意一点位移可表示为：
单元节点基本未知量：
单元节点位移向量：
为单元的阶形函数矩阵。

将式代入几何方程中，可得单元应变与节点基本未知量的关系式：
为单元的阶应变矩阵。

计算算例与结果
4.1两相复合材料的有效模量
功能梯度材料是由两种不同材料混合而成，通常其材料颗粒的精确分布则无法获得，所以其有效材料模量只能基于材料的体积分量和分散相的大致形状来估算。

目前，已有几种以从微观力学观点得出计算复合材料宏观模量的模型，本文选择比较流行的Mori-Tanaka模型[13]。

Mori-Tanaka模型适用于计算微观结构为一种材料构成连续的基质和另一种材料以离散的颗粒包裹于其中的复合材料有效模量。

基质相材料1的体积模
量、剪切模量和体积分量记为、和，颗粒相材料2的体积模量、剪切模量和体积分量记为、和，显然有。

4.2模型参数
假定是由陶瓷(SiC)和铝(Al)组成的功能梯度板，其板厚。

对于Al：弹性模量，泊松比，对于SiC：弹性模量，泊松比。

体积模量和剪切模量与弹性模量和泊松比的关系式为和，其材料组分沿厚度方向的变化函数为
和分别是陶瓷相在板的上下表面的体积分量。

并令金属相(Al)为基质相，陶瓷相(SiC)为颗粒相，则有和。

4.3计算结果对比与分析
假定板四边简支，即在上；，在上。

则在有限元模型中，对于广义自由度有，在上；，在上。

在板的上表面施加大小为的分布正应力，对于本例，取。

将该单元划分为个四边形版单元，则可令表示单元划分密度。

将计算结果无量纲化，令
图3、图4、图5分别为单元密度为，材料参数，，时用本文方法求得的对于不同的长厚比的功能梯度板板中心挠度、板中心上表面正应力和板边缘中点处的剪应力的数值解与Vel和Batra[13]求得的相应的精确解的对比。

从2到40，可以看出本文方法不论对于厚板还是薄板，都取得令人满意的结果。

图3 图4
图5图六
图6为，，，长厚比时，本文方法求得的功能梯度板板中心沿截面的挠度、正应力分布和板边缘中点沿截面的剪应力分布的数值解与Vel和Batra[13]求得的精确解的误差随网格密度增加，即计算自由度的增加的误差的变化趋势。

可以看出随着自由度增加数值解趋向于精确解，当误差不大于，当时，数值计算结果已趋于稳定，且误差不大于。

结论
本文在分析功能梯度板时，摒弃了一般板理论的简化性假设，而直接从三维弹性力学方程出发，将高阶剪切变形和法向变形板理论的位移模式代入，并推导出在新的位移模式下的有限元计算格式。

然后采用一种平面4节点单元，对四边简支的功能梯度材料板进行计算分析，其中假设梯度混合材料遵循Mori-Tanaka 模型，材料的体积组分函数沿板厚方向以多项式函数形式分布。

用本文方法得到的不同材料梯度参数的功能梯度板数值计算结果与精确解能够较好的吻合。

并且通过增加模型的计算自由度数，验证了其数值解结果收敛于精确解的趋势。

说明本文提出的有限板单元适用于不同梯度分布的薄板与厚板的静力分析计算，因为在计算刚度矩阵的过程与具体材料分布函数和采用的复合材料本构无关，克服了一般解析方法在计算功能梯度板时只能针对某种特定材料分布或本构模型的局限性。

而且这种有限板单元易于扩展，可进一步将其应用于功能梯度材料板的热动力分析中。

注：文章内所有公式及图表请用PDF形式查看。