有限元4薄板弯曲问题.docx

有限元4薄板弯曲问题.docx

《有限元4薄板弯曲问题.docx》由会员分享，可在线阅读，更多相关《有限元4薄板弯曲问题.docx（22页珍藏版）》请在冰豆网上搜索。

有限元4薄板弯曲问题

第4章弹性薄板弯曲问题的有限元法

薄板弯曲问题在理论上和应用上都具有重要意义,并有专门著作加以论述（如耀乾《平板理论》）。

象其它弹性力学问题一样，用微分方程、差分法等经典方法所能求解的薄板问题很有限，一般只能解决等厚、小孔口、支承情况较简单的单跨板。

故工程设计中以往多采用简化、近似、图表等方法来解决板的设计问题。

在板的分析中，常取板的中面为xoy平面（如图）。

平板结构按其厚度t与短边a的比值大小而分为:

厚板（Thickplate）和

薄板（Thinplate）两种。

当

时称为薄板

平板上所承受的荷载通常有两种:

1.面拉压荷载。

由面拉压刚度承担,属平面应力问题。

2.垂直于板的法向荷载,弯扭变形为主,具有梁的受力特征,即常说的弯曲问题。

平板在垂直于板面的荷载作用下产生挠度W。

当最大挠度w远小于t时,称为小挠度问题（or刚性板）（stiffnessplate）

当最大挠度w与t相差不大时,称为大挠度问题（or柔性板）（flexureplate）

（工程定义:

为刚性板；

为柔性板；

为绝对柔性板。

）

4.1基本理论

一、基本假定

１、略去垂直于中面的法向应力。

（

），即以中面上沿Z方向的挠度W代表板的挠度）

２、变形前垂直中面的任意直线，变形后仍保持为垂直中面的直线。

（─法向假定

）

３、板弯曲时，中面不产生应力。

（─中面中性层假定）

上述假定常称为薄板小挠度问题假定（or柯克霍夫假定）。

符合上述假定的平板即为刚性板。

二、基本方法

以上述假定为基础，板分析中常用挠度

作为基本未知量,下面介绍以

为基本未知量所导出的有关方程。

１、几何方程（应变─挠度关系）

①弹性曲面沿x,y方向的倾角

从中面取出一微小矩形ABCD，如图所示，设其边长为dx,dy，变形后弯曲成曲面A'B'C'D'

设A点挠度

则沿x方向倾角（绕y轴）

（B’点绕度

）

沿y方向倾角（绕x轴）

（D’点绕度

）

②沿x,y方向位移

作平行于

平面,设中面上点A到A1的距离为Z，变形后,A点有挠度W,同时发生弯曲，

曲面沿x方向的倾角为

根据法线假定,则A1点沿x方向的位移:

（负号为方向与x相反）

同理取

平面得：

（4-1-1）

③Z平面的应变分量和曲、扭率

基本假定，由于

故板任意点的应变与平面问题相同:

（4-1-2）

此为Z平面的应变─挠度度几何方程。

上式中的

为曲面在X,Y方向的曲、扭率，记为：

（4-1-3）

所以,

２、物理方程（应力─挠度关系）

由于忽略σz对变形的影响,因此z平面的应力─应变关系具有与平面问题相同的形式:

将（4-1-2）代入得:

或简写为：

（4-1-4）

式中弹性矩阵:

３、力方程（力─挠度关系）

从板取微元体

由其上正应力

和剪应力

可在截面上合成合力矩:

（

面上由

产生的绕Y轴弯矩）

（

面上由

产生的绕X轴弯矩）

扭矩:

（由剪应力产生，如图）

假定

分别表示单位宽度上的力矩。

如是，力矩阵：

简写成

（4-1-5）

比较（4-1-4）和（4-1-5）可得用力矩表示的平板应力:

由此可见，平板上、下表面处的应力最大：

以上是薄板弯曲问题中的基本公式,从中可见其挠度W是弯曲问题中的基本未知函数。

且由于忽略了z方向的变化,因此它只是x，y的函数:

w=w（x,y）。

若w已知，则位移，力、应力均可按上述相应公式求出。

在经典解析法中，W（x,y）常设为三角级数形式。

例如，四边简支矩形板的W（x,y）设为:

（纳维尔解）

式中

为待定系数。

假定荷载

则可得位移函数:

4.2有限元分析方法

一、矩形单元的典型形式

将图示矩形薄板沿x,y方向划分成若干小矩形（常取等分）

从中取出一小矩形（单元），共有四个结点，此时不能象在平面问题中一样，将结点视为“铰”，而是“刚性的”，即每个结点有三个位移分量：

挠度

绕x、y轴转角

即结点i的位移

同理,相应的结点力

符号重新定义是为了有限元表示的方便，由此得单元结位移向量

节点力

二、位移模式（函数）

１、位移模式的选取

插值多项式取为：

（4-2-1）

在上式中,前10项取到了三次项的全部,最后两项则是从五个四次项

中选用了两个。

没选

是因为它没有多一项与其配对，没选

它们在边界上结出的挠度函数是四次的，比

和

要高一次,较之更难满足边界的协调和条件。

２、位移模式的检验

（三个基本要求：

刚体位移，常应变，尽可能的边界协调）

①前三项含单元的刚体位移状态：

第一项

与坐标x,y无关,表示z方向的挠度是─常量,刚体移动

②二次项代表均匀变形状态:

曲率

，

，

③能保证相邻单元在公共边界上挠度的连续性。

④不能保证相邻单元在公共边界上法线转角的连续性。

以单元１～２边界为例，在此边界上

=常量，代入位移模式4-2-1，可知边界上的挠度W是x的三次函数，合并整理后可得：

两个端点共有4个边界条件,（结点1,2的挠度W1,W2,和转角

。

利用他们可唯一确定四个常数C1～C4。

因为相邻单元在结点1,2的W,θy对应相同，则两个单元依据四个条件得到的C1～C4亦相同，即两单元在边界具有同一挠度函数W。

⑤法线转角

仍以1-2边界为例，将y=-b代入后，此时

但对θx来讲，1,2结点只能提供2个已知条件，不能完全确定上式，故边界的法线转角不能保证连续性。

因此，这种单元是非协调元，但可以验证这种非协调远是能通过分片试验的。

（即当单元划分不断缩小时，计算结果仍能收敛于精确解。

）

三、形函数和形函数矩阵。

分别将单元结点1,2,3,4的坐标值代入（4-2-1），并事先求出

，

，便可得到各结点的位移值。

一共可得12个关于

的方程组，联立求解可得：

（4-2-2）

形函数矩阵:

式中形函数:

（4-2-3）（i=1234）

在上面的推导中,我们仍然选用了局部坐标（无因次坐标）。

局部坐标与整体坐标的关系为:

四、单元的几何矩阵[B]和力矩阵[S]

1．几何矩阵[B]

由前可知

将（4－2－3）代入（4-2-4）得到几何矩阵：

（4-2-5）

或以子块形式表示：

[B]=[B1B2B3B4]。

式中：

2.力矩阵[S]

由基本方程（4-2-5）可得到：

（4-2-6）

称为力矩阵，把单元的四个结点坐标分别代入4-2-4，求得

后，即可获得

，各节点力矩阵

的显式：

五、单元刚度矩阵

由一般公式得:

。

将几何矩阵[B]和弹性矩阵[D]的表达式代入，积分可得薄板弯曲问题矩形单元的单元刚度矩阵的显示：

六、荷载等效变换

由荷载等效变换的一般公式可得

1．法向均布荷载q

代入上述公式得：

2．单元中心点受法向集中力P

代入上述公式可得：

七、位移边界条件

对称、固定边和简支边上支点的已知位移条件如下：

对称轴:

法线转角=0

固定边:

挠度=0（或已知值）

边线转角=0（或已知值）

法线转角=0（或已知值）

简支边:

挠度=0（或已知值）

边线转角=0（或已知值）

自由边上节点的挠度、边线和法线转角均为特定参数，同部节点一样。

与板铰接的固定立柱，其节点挠度Ｗ=0，也可以是已知值。

八、计算例题

例题1:

计算图示四边固定方板

方板的边长为l，厚度为t，弹性模型量为E，波松比μ=0.3，全板承受均布法向荷载ｑ，求薄板中的挠度和力。

单元划分：

为了说明解题方法，采用最简单的网络2×2，

即把方板分成四个矩形单元。

由于对称性，只需计

算一个单元，例如，计算图中有阴影的单元，单元

的节点编号为１，２，３，４。

此时，单元的a,b是

计算节点荷载:

由前面的均布荷载计算公式得：

边界条件：

边界23和34为固定边，因此节点2,3,4的挠度、边线和法线转角均为零。

边界12和14为对称轴，因此θx1=0、θy1=0。

于是，在4个节点和12个位移分量中,只有一个待求的未知量

。

结构的代数方程组：

这是一个单元的计算题目，单元刚度矩阵在此处即为总刚度矩阵。

引入支承条件后，在总刚度矩阵中只取第一行、列元素，在方程组右端项中只保留第一个元素。

于是结构的代数方程为：

同此解出

。

其中

力:

利用式（4-2-6）可求得方板中点力矩为：

由表看出，网格越密，计算结果越接近于精确答案。

还可看出，位移的精度一般比力的精度高，这是因为在位移法中，位移是由基本方程直接求出的，而力则是根据位移间接求出的。

4.3薄板有限元程序设计

一、总框图

根据弯曲板有限元分析方法的解题过程，可写出其总框图如下：

┌───────┐

│输入原始数据│

│orCAI│

└───┳───┘┌──────┐

↓┌──┤算等效结点力│

┌───┻───┐│└──────┘

│形成荷载列阵├←┘┌─────┐

│├←───┤形成单元│

└───┳───┘┌─┤定位向量├─┐

↓│└─────┘│

┌───┻───┐││

│形成总刚├←─┘┌─────┐│

│├←───┤单刚││

└───┳───┘└─────┘│

↓│

┌───┻────┐│

│解方程输出位移││

└───┳────┘│

↓┌──────┐│

┌───┻───┐│几何矩阵[B]││

│├←───┤弹性矩阵[D]││

│计算单元力等│└──────┘│

│├←───────────┘

└───┳───┘

↓

┌──┻──┐

│结束│

└─────┘

下面结合程序对框图中的容加以说明。

二、子框图

１、单元坐标结点编号及单刚形式。

为了取挠度向下为正，又能与前述坐标系统统一，特将前述坐标前翻180°（如图）

为了能适用板的弹型性分析,程序采用了应力元和弯曲元的组合形式，即每个结点考虑5个位移分量:

U,V,W,θx,θy,前2个为平面应力问题的未知量,后3个为弯曲板的结点未知量。

当只作弹性分析时，平面应力元和弯曲元是非藕连的，即单刚的两个副块垣为0,单刚的形式为：

u1v1…u4v4w1θx1θy4…w4θx4θy4

┌┐

│平面应力元│0│

[K]e=│（8×8）││

├──────┼──────────────┤

││弯曲元│

│0│（12×12）│

└┘

程序中单刚数组为DK（20,20）,子程序:

SubroutineDG（A,B,E,T,U）为其形成单刚的子程序。

２、自动形成单元编号信息（单元信息数组：

[IB]）。

３、结点定位向量。

４、形成荷载列向量。

（a.结点力；b.非结点力（只考虑均布力））

５、总刚，SubroutineZG（M,N,LD,A,B