杆系结构有限元分析1

1、第一章杆系结构有限元分析1-1 概述 杆系是工程中常见的结构体系，比较简单，其中每一个杆件都可以看作是一个单元，单元受力与位移的关系很容易求得而且物理概念清晰、直观。结构力学中介绍的矩阵位移法是采用经典的方法讲述的，它是利用转角位移方程来建立单元特性公式，所以只适用于杆系。有限元方法是在结构矩阵分析的矩阵位移法基础上发展起来的，在建立位移场的过程中采用的是最有普遍意义的方法，即建立单元位移场函数，通过最小势能原理进行单元和整体分析。杆系结构分类按变形分类轴向变形杆件桁架结构扭转变形杆件传动轴系弯曲变形杆件刚架结构按轴线分布分类平面结构体系空间结构体系杆系结构有限元分析方法局部坐标系下单元刚度矩

2、阵整体坐标系下单元刚度矩阵集成结构整体刚度矩阵单元分析整体分析平衡求解引入边界条件、求解线性方程组，得未知量解答一般规定 杆单元ij，单元局部坐标系为oxyz，i点为原点，x轴沿着杆轴线，其正方向为由i指向j，其余各轴按右手螺旋规则确定。设ui，vi，wi，uj，vj，wj为杆元结点位移分量；Ui，Vi，Wi，Uj，Vj，Wj为杆单元结点力分量，一律规定和坐标轴正向一致时为正。设杆的长度为l，弹性模量为E，横截面积为A。ijxyzui( Ui )vi( Vi )wi( Wi )uj( Uj)vj( Vj )wj( Wj )1-1 拉(压)杆单元位移函数 对于铰接杆单元，在小变形假设的前提下，与

3、杆垂直方向的位移并不使杆产生应变和应力。对每一个结点只需考虑一个结点位移和结点力，因而只需研究一维杆单元。ijxui( Ui )uj( Uj)lx单元在结点力作用下各点的位移叫内位移描绘内位移的函数叫位移函数设位移函数：u(x)=a1+a2xa1,a2是两个待定常数，可由i，j两结点的位移唯一确定。 1-1 拉(压)杆单元x=0，u(0)=uix=l，u(l)=uj代入位移函数：u(x)=a1+a2xa1=ui , a2=(uj-ui)/l则:u(x)=(1-x/l)ui +(x/l) uj ui或写成: u(x)=1-x/l x/l uj ui=Ni Nj uj=Nue在有限元法中，Ni、N

4、j分别称为i、j点的形状函数N称为形状函数矩阵形状函数矩阵把单元的结点位移和单元的内位移连接起来，其每一个元素都是坐标的函数。1-1 拉(压)杆单元当ui=1，uj=0时，杆单元的位移u(x)就是Ni当ui=0，uj=1时，杆单元的位移分布就是Nj位移分布规律：N=1-x/l x/l则:u(x)=(1-x/l)ui +(x/l) ujNi1Nj11-1 拉(压)杆单元uiuj1-1 拉(压)杆单元形状函数的力学含义：当单元的一个结点位移为单位值，其他结点的位移为零时，单元内位移的分布规律。形状函数的两个重要性质为： 1、一结点为1，其他结点为0； 2、任意一点总和为1。自然结点离散化为有限元的

5、集合，实现了结构模型离散化，那么，形状函数完成了数学模型离散化，这两个离散化的步骤构成了有限元法的理论基础。形状函数把两孤立的常值位移，化为连续函数几何关系和物理关系或写成单元的应变和应力，根据应变定义 ui带入位移函数： u(x)=1-x/l x/l uj应变矩阵 1-1 拉(压)杆单元对于拉(压)杆，应力与应变之间的关系有用矩阵表示为D为弹性矩阵 其中 S=DB 称为应力矩阵对于拉（压）杆单元有1-1 拉(压)杆单元平衡关系杆单元结点力向量单元在外力和内力作用下处于平衡状态，反映单元平衡状态的关系式就是刚度方程。下面利用最小势能原理推导单元的刚度方程。最小势能原理：在满足连续条件和边界条件

6、的位移中，满足平衡条件的位移其总势能最小，反之亦然。单元总势能其中U e为单元的应变能，V e为单元的外力势能。 1-1 拉(压)杆单元令则外力势能总势能1-1 拉(压)杆单元根据最小势能原理，势能泛函取驻值的必要条件即：杆单元的平衡方程 1-1 拉(压)杆单元整体坐标系的刚度矩阵设oxyz为总体坐标系，oxyz为局部坐标系。同理对于j 结点用矩阵表示为1-1 拉(压)杆单元对两结点杆单元，当用总体坐标系位移ue表示局部坐标系中位移(u)e 时有转换关系T 矩阵称为坐标变换矩阵，是以为子矩阵的对角方阵，为正交矩阵。当用局部坐标系位移表示总体坐标系中的位移时有1-1 拉(压)杆单元令 类似方法，

7、总体坐标系与局部坐标系间结点力的关系式将式 代入上式将式 代入上式单元总体坐标系的平衡方程1-1 拉(压)杆单元1-1 拉(压)杆单元整体坐标系下的单元刚度矩阵1-1 拉(压)杆单元 设局部坐标的x，y，z 轴与整体坐标系x，y，z轴之间夹角方向余弦为：cos(x,x) cos(x,y) cos(x,z) cos(y,x) cos(y,y) cos(y,z) cos(z,x) cos(z,y) cos(z,z)则相应的位移关系为空间杆单元1-1 拉(压)杆单元方向余弦矩阵空间杆单元坐标变换矩阵单元在两个坐标系中刚度矩阵转换关系同样有矩阵中仅仅包含有坐标的倾角，仅平行移动坐标轴，刚度矩阵中元素值

8、不变，矩阵的阶数也不改变。1-2 扭转杆单元在局部坐标系中，每一个点将具有一个基本未知位移，最简单的单元位移函数可以设为杆件发生自由扭转时，待求位移是截面的扭转角结点位移向量结点力向量其中的待定常数可以用两端节点的扭角 i 、j表示，从而任意截面的扭转角可由结点位移和形函数表示为其中由材料力学扭转可知其中 为形函数。1-2 扭转杆单元单元刚度方程扭转杆的势能为单元刚度方程为局部坐标扭转杆单元刚度矩阵由泛函取驻值的必要条件有其显式为1-2 扭转杆单元结点位移矢量：ue=ui vi i uj vj j T结点力矢量：f e=Ni Qi Mi Nj Qj Mj T正负规定：轴力、剪力与局部坐标系坐标

9、轴正向一致为正； 弯矩以顺时针转动为正。梁单元结点位移梁单元结点力1-3 平面直梁单元ijuiviiujvjjlijNiQiMiNjQjMjl 先将轴力与剪力、弯矩分开考虑。直梁弯曲时满足平截面假设，原来垂直轴线的平面，变形后仍垂直于轴线。若梁中面挠度为v，则因弯曲而引起的轴向位移为其中，y 是所讨论点距中性轴的距离，位移函数应变根据合成关系有（其中，I 是截面惯性矩 ）根据胡克定律1-3 平面直梁单元又根据荷载集度与弯矩的微分关系：现讨论的梁EI为常数若梁上无分布荷载，即q=0由此可以判断v是x的三次函数则梁的转角1-3 平面直梁单元将i、j结点位移代入上两式用矩阵表示从而得1-3 平面直梁

10、单元单元的形状函数矩阵Nv为其中1-3 平面直梁单元设则 以上四个函数即是二结点梁单元的形状函数。由于每结点有两个位移参数，则每结点有两个形状函数。 H0i是指i结点零阶导数（即i点位移）对应的形函数，H1i指i结点一阶导数（即i结点转角）对应的形函数，位移函数用内插多项式表示为1-3 平面直梁单元当x=0，即=0， ；当x=l，即=1, H0i(1)=0；当x=0，x=l时H0i曲线在i、j点切线平行于x轴，即H0i不引起结点转角改变。H0i表示当i点垂直位移vi=1，其他位移等于零时梁元的挠曲形状。形状函数的性质:xyij1H0i(x)xyij1H0j(x)xyiji=1H1i(x)xyi

11、jH1j(x)j=1(a) (b) (c) (d)1-3 平面直梁单元由此得其中轴力N引起的位移u(x)仍设为线性，u(x)=a0+a1x，将结点位移引入可求得1-3 平面直梁单元位移函数改写为则形状函数N为现将结点位移列阵合并为1-3 平面直梁单元位移函数求得后，可得到应变和应力的表达式：忽略剪切影响，设N是轴力引起的应变，b是弯曲引起的应变应力为梁的应力应是N 和b 代数和1-3 平面直梁单元梁上的结点力f e=Ni Qi Mi Nj Qj Mj T，并有分布力q(x)作用。在选位移函数时虽然假设了q=0，若作用有分布载荷q(x)，位移函数仍可用三次幂函数近似，分析过程完全同前。只是外力势

12、中增加梁元的刚度矩阵现根据最小势能原理求梁单元刚度矩阵，梁的应变能1-3 平面直梁单元总势能式中， 取驻值时有：是分布载荷的等效结点力。梁单元刚度方程刚度矩阵1-3 平面直梁单元将B带入刚度矩阵并积分，得： 1-3 平面直梁单元对短梁（hl/5），应计及剪切影响，对刚度矩阵作如下修正式中， 是剪切影响系数，AS是有效抗剪面积。1-3 平面直梁单元整体坐标系的刚度矩阵设oxyz为总体坐标系，oxyz为局部坐标系。规定由总体坐标系x轴到局部坐标系x轴的夹角逆时针为正。杆单元总体坐标系下的结点位移分量用u,v,表示，局部坐标下的位移分量用u,v,表示。则平面杆单元结点i在总体坐标系和局部坐标系下的位

13、移分量关系有1-3 平面直梁单元同理对于j结点有用矩阵表示为1-3 平面直梁单元总体坐标系位移ue和(u)e的转换关系为 T 矩阵称为坐标变换矩阵，是以为子矩阵的对角方阵，为正交矩阵。当用局部坐标系位移表示总体坐标系中的位移时有类似于拉(压)杆单元，可得直梁元刚度矩阵矩阵变换关系仍为1-3 平面直梁单元1-3 平面直梁单元等效结点载荷所做的功等效结点力 等效结点力是根据功互等原理，将分布载荷转移到结点上所得到的载荷。这样的变更，对全结构的计算不会带来明显的误差。但对载荷区域单元的应力分布，将有较大的影响。设单元发生虚位移r*时，分布力q(x)所做的虚功为根据虚功原理，单元分布力所作虚功与等效荷

14、载虚功相等r*=N(u*)e1-3 平面直梁单元从而得到分布载荷的等效结点力计算公式Nu是轴向位移形函数，线位移时当p(x)为均布载荷时,p(x)=p，则，即两结点各半。（1）分布轴力p(x)的等效结点力p(x)ij1-3 平面直梁单元(2) 分布横向荷载q(x)的等效结点力q(x)ij1-3 平面直梁单元1-3 平面直梁单元(3) 分布力矩mz(x)的等效结点力mz(x)ij ，因此相应的形状函数矩阵为对应mz的转角当挠度为x的三次式时1-3 平面直梁单元式中，当 为均匀分布情况时1-4 总体刚度矩阵直接刚度法（定位向量） 结构总体刚度矩阵可以通过直接处理整个结构的方法（虚位移原理，势能原理

15、）而得到。通过这种途径得到的结构总体刚度矩阵物理意义清晰，但它不是一种容易编制程序的系统化方法。1-4 总体刚度矩阵 通过组装单元刚度矩阵以形成总体刚度矩阵。 结构总体刚度矩阵的元素是由单元刚度矩阵的元素组成的，只要确定了单元刚度矩阵各元素在结构总体刚度矩阵中的位置，就可以由单元刚度矩阵直接集成结构总体刚度矩阵。 把单元杆端位移分量（局部码）所对应的结构结点位移向量的序号（整体位移码）组成一向量，它成为单元的定位向量。利用单元定位向量可以完全确定单元刚度矩阵的每个元素在结构（原始）刚度矩阵中的行码和列码。直接刚度法（定位向量）1-4 总体刚度矩阵由定位向量对号入座集装结构总体刚度矩阵的具体做法

16、：（1）求单元(e)在整体坐标系中的刚度矩阵k(e);（2）将单元(e)的定位向量分别写在单元刚度矩阵k(e)的上方和右侧（或左侧）,形成k(e)的每一行或每一列与单元定位向量的一个分量相对应关系，这个分量即为k(e)中相应的行和列在结构整体刚度矩阵k中的行码或列码;（3）按照由单元定位向量中分量的行码和列码，将单元刚度矩阵k(e)的元素正确地叠加到结构刚度矩阵k中。1-4 总体刚度矩阵1-4 总体刚度矩阵下图所示的平面桁架结构1(1,2)2(3,4)3(5,6)x,uy,v 1(1,2)、2(3,4)、3(5,6)中括号前面的数字代表整体结点编号，括号中的数字即为位移向量的序号。 拉(压)杆

17、单元刚度矩阵为一个四阶矩阵，设每个杆单元的单刚矩阵分别如下：1-4 总体刚度矩阵利用单元定位向量叠加结构刚度矩阵为刚度矩阵的物理意义1-4 总体刚度矩阵对应的节点力记为其刚度方程记为设某单元或结构有n个节点位移或者1-4 总体刚度矩阵在上式中，令也就是说，如果节点位移ui为单位值而其他为零(固定不动) 则 刚度矩阵的诸元素(刚度系数)在数值上等于造成某种特定变形状态时所需要施加的各节点上的力。刚度矩阵的第i列，就是为使第i个节点位移为“1”，而其余节点位移均为零时，所施加于节点上的全部力系。1-4 总体刚度矩阵1K1iK2iKiiUi=1KjiKmi 令节点k的垂直方向位移Ui=1，而该点的水

18、平位移以及其他节点的全部位移均固定不动，则k点施加的垂直力的大小就等于Kii，而其他支反力的大小分别等于K1i，K2i，KjiKni，这些系数的值愈大，说明这个结构(或单元)愈难以发生变形，反映了单元或结构抵抗变形的能力，因而得名为“刚度系数”。另外，这些系数所体现的刚度特性又反映了各节点位移和节点力相互之间的影响，所以又称之为“刚度影响系数”。1-4 总体刚度矩阵 其中Kij(i=j)反映节点位移与本身对应的节点力之间的影响，称为“自身影响系数”。在刚度矩阵中总处在主对角线上，故又称为“主系数”；Kij(ij)反映不对应的节点力同节点位移之间的影响，称为“交叉影响系数”。它们在刚度矩阵中总是

19、占据主对角线两旁的位置，故又称为“副系数”。 总之，任一刚度系数在数值上等于使结构(或单元)产生单位位移而其余节点位移为零时第 i 个节点力的值。尽管它们在数值上可能相等，但因次却不一定相同。刚度系数的因次是力的因次/位移的因次（力和位移都是指的广义力和广义位移）刚度矩阵的性质1-4 总体刚度矩阵将刚度方程代入上式对于线性弹性体体系，应变能表达式刚度矩阵有以下性质:1. 主系数必为正值1-4 总体刚度矩阵2. 刚度矩阵是正定的将应变能表达写成展开的形式 W是关于变量u的二次齐次多项式。在线性代数学里，这种多项式称为“二次型”。而由多项式系数所组成方阵K称为这二次型的矩阵。不论位移列阵u取何种数

20、值，除非u=0，应变能W总是正值。二次型的矩阵也就称作是“正定矩阵”。1-4 总体刚度矩阵3. 刚度矩阵是对称矩阵 这一性质是由功的互等定理决定的。即：对于线性弹性体，第一种加载状态下的诸力在第二种加载状态下移动相应位移时所作的功等于第二种加载状态下的诸力在第一种加载状态下移动相应位移时所做的功。根据前面所述刚度系数的物理意义，系数Kij和Kji分别包含在这样两种加载状态的诸力中，其对应的外力功是1Kij=Kji1，故 Kij=Kji。1-4 总体刚度矩阵4. 刚度矩阵的任一行(或列)代表一个平衡力系。 当节点位移列阵的单元全部为线位移时，任一行(或列)的代数和应为零。 由前所述，刚度矩阵的任一列在数值上等于某种特定位移状态下的全部外力和支反力，它们当然构成一个平衡力系。而由对称性可知，任一行也就具有同样