对称结构有限元分析----3节点三角形单元的分析一问题分析(对称框架线弹性实体的静力平衡问题)图是一个方形弹性实体,单位边长、单位厚度、承受等效竖向压力21m,其中边界条KN件暗示着存在两组相对称的平面,因此现考虑的仅是问题的。
每个节点上的自由度号码代表了各自在x和y方向上可能的位移。
结构和单元信息NELS NCE NN NIP8 2 9 1AA BB E V.5 .55 1.E6 .3约束节点自由度信息NR5K , NF(:,K), I=1,NR10 1 4 0 1 7 0 0 8 1 9 1 0 载荷信息LOADED_NODES3(K, LOADS(NF(:,K)), I=1 , LOADED_NODES)1 .0 -.252 .0 -.53 .0 -.253333节点三角形单元网络的总体节点和单元编号3节三角形单元局部坐标系中节点和自由度编号二理论基础(有限元方法原理)通过弹性力学变分原理建立弹性力学问题有限元方法表达格式的基本步骤。
最小位能原理的未知场变量是位移,以结点位移为基本未知量,并以最小位能原理为基础建立的有限元为位移元。
它是有限元方法中应用最为普遍的单元,也是本书主要讨论的单元。
对于一个力学或无力问题,在建立其数学模型以后,用有限元方法对它进行分析的首要步骤是选择单元形式。
平面问题3结点三角形单元是有限元方法最早采用,而且至今仍经常采用的单元形式。
我们将以它作为典型,讨论如何应用广义坐标建立单元位移模式与位移插值函数,以及如何根据最小位能原理建立有限元求解方程的原理、方法与步骤,并进而引出弹性力学问题有限元方法的一般表达格式。
对于前一问题,着重讨论选择广义坐标和有限元位移模式的一般原则和建立其位移插值函数的一般步骤。
对于后一问题,着重讨论单元刚度矩阵和单元载荷向量的形式,总体刚度矩阵和总体载荷向量集成的原理和方法,以及它们各自的特性。
作为一种数值方法,有限元解的收敛性无疑是十分重要的问题,以后将讨论解的收敛准则及其物理意义,所阐明的原则在以后还将得到进一步的应用和具体化。
在建立了有限元的一般表达格式以后,原则上可以将它推广到平面问题以外的其他弹性力学问题和采用任何形式的单元。
轴对称问题具有很广泛的应用领域,轴对称问题3结点三角形 单元的表达格式可以看作是平面问题此种单元表达格式的直接推广。
一)弹性力学平面问题的有限元格式结点三角形单元是有限元方法中最早提出,并且至今仍广泛应用的单元,由于三角形单元对复杂边界有较强的适应能力,因此很容易将一个二维离散成有限个三角形单元,如图1所示。
在边界上以若干段直线近似原来的曲线边界,随着单元增多,这种拟合将趋于精确。
我们在讨论如何应用有限元方法分析各类具体问题的开始,将以平面问题3结点三角形单元为例来阐明弹性力学问题有限元分析的表达格式和一般步 1.1)单元位移模式及插值函数的构造典型的3节点三角形单元节点编码i,j,m ,以逆时针方向编码为正向。
每个节点有位移分量如图所示。
⎥⎦⎤⎢⎣⎡=i i v u i a (i,j,m)每个单元有6个节点位移即6个节点自由度,亦即[]Tmm j j i im j i ev u v u v u a a a =⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=a1.2) 单元的位移模式和广义坐标在有限元方法中单元的位移模式或称位移函数一般采用多项式作为近似函数,因为多项式运算简便,并且随着项数的增多,可以逼近任何一段光滑的函数曲线。
多项式的选取由低次到高次。
3结点三角形单元位移模式选取一次多项式 y x u 321βββ++=(1)y x v 654βββ++=它的矩阵表达式是φβ=u (2) 其中⎥⎦⎤⎢⎣⎡=v u u ⎥⎦⎤⎢⎣⎡=ϕϕφ00 []y x1=ϕ []T621ββββ=φ称为位移模式,它表示位移作为坐标x ,y 的函数中所包含的项次。
对于现在的情况,单元内的位移是坐标x ,y 的线性函数;61~ββ 是待定系数,称之为广义坐标。
6个广义坐标可由单元的6个结点位移来表示。
在(1)的1式中代入结点i 的坐标(x ,y )可得到结点i 在x 方向的位移i μ,同理可得j μ和 m μ。
它们表示为i i i y x u 321βββ++=j j j y x u 321βββ++= (3) m m m y x u 321βββ++=解(2)式可以得到广义坐标有结点位移表示的表达式。
上式的系数行列式是A x x y x y x D mmj j ii 2111== (4) 其实A 是三角形单元的面积。
广义坐标61~ββ为()m m j j iimmj j j i i i u a u a ua Ay x u y x u y x u Dm++==2111β)(2111112m m j j i i m m j j i i u b u b u b Ay u y u y u D ++==β (5))(2111113m m j j i i mmj j i i u c u c u c Au x u x u x D++==β同理,利用3个结点y 方向的位移,即(a )式的第2式可求得 )(214m m j j i i v a v a v a A ++=β)(215m m j j i i v b v b v b A ++=β (6) )(216m m j j i i v c v c v c A++=β在(c )式和(d )式中 j m m j mmj j i y x y x y x y x a -==m j m j i y y y y b -=-=11 ),,m j i ( (7)m j mj i x x x x c +-==11上式(i ,j ,m )表示下标轮换,如i →j ,j →m ,m →i 。
以下同此。
1.3) 位移插值函数将求得的广义坐标 代入(1)式,可将位移函数表示成结点位移的函数,即 m m j j i i u N u N u N u ++=m m j j i i v N v N v N v ++= (8) 其中)(21y c x b a AN i i i i ++=(i,j,m) (9)m j i N N N ,, 称为单元的插值函数或形函数,对于当前情况,他的坐标x 、y 的一次函数。
其中的 m i i i c c b a ,,,.....,是常熟,取决于单元的3个结点坐标。
g )式中的单元面积A 可通过(7)的系数表示为 )(21)(2121i j j i m j i c b c b a a a D A -=++== (10)(f )式的矩阵形式是[][]eem jim j i mjim m jj i i m jim jiNaa N NN a a a IN ININ v u v u v u N NN N N N v u u ==⎪⎪⎪⎭⎫ ⎝⎛=⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎥⎦⎤⎢⎣⎡=⎥⎦⎤⎢⎣⎡=0000 (11)N 称为插值函数矩阵或形函数矩阵,e a 称为单元结点位移列阵。
(1) 插值函数具有如下性质 在结点上插值函数的值有⎩⎨⎧≠===ij i j y x N ij j j i 当当01),(δ (i,j,m)(12)即有0),(),(,1),(===m m j j j i j j i y x N y x N y x N 。
也就是说在i 结点上1=i N ,在j ,m 结点上0=i N 。
由(8)式可见,当i i y y x x ==,即在结点i ,应有i u u =,因此也必然要求0,1===m ji N N N 。
其他两个形函数也具有同样的性质。
此性质称为kronecker delta性质。
(2)在单元中任一点各插值函数之和应等于1,即1=++m j i N N N (13) 因为若单元发生刚体位移,如x 方向有刚体位移0μ,则单元内(包括结点上)到处应有位移0μ,即0u uu u mj i ===,又由(g )式得到00)(u u N N N u N u N u N u m j i m m j j i i =++=++=因此必然要求1=++m j i N N N 。
若插值函数不满足此要求,则不能反映单元的刚体位移,用以求解必然得不到的正确的结果。
单元的各个结点位移插值函数之和等于1的性质称为规一性。
(3)对于现在的单元,插值函数是线性的,在单元内部及单元的边界上位移也是线性的,可由结点上的位移值唯一地确定。
由于相邻单元公共结点的结点位移是相等的,因此保证了相邻单元在公共边界上位移的连续性。
1.4)应变矩阵和应力矩阵确定了单元位移后,可以很方便地利用几何方程和物理方程求得单元的应变和应力。
在几何方程中,位移用(11)式代入,得到单元应变为 [][]ee m j iem jiexyyx Baa B B B aN NN L LNa LU =====⎪⎪⎪⎭⎫⎝⎛=γεεε (14)B 称为应变矩阵,L 是平面问题的微分算子。
应变矩阵B 的分块子矩阵是⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎣⎡∂∂∂∂∂∂∂∂=⎥⎦⎤⎢⎣⎡⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎣⎡∂∂∂∂∂∂∂∂==xN yNyN x NN N x yy x LNB i ii ii iii 000000 (i,j,m) (15) 对(9)式求导可得 i i b AxN 21=∂∂i i c AyN 21=∂∂ (16)代入(15)式得到⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=i ii ii b c c b 00A 21B (i,j,m) (17) 3结点单元的应变矩阵是[]⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡==m mjjii m j i m i mj i b c b c b c c c c b b b AB B B B j 00000021 (18)式中m j i m j i c c c b b b ,,,,, 由(7)式确定,它们是单元形状的参数。
当单元的结点坐标确定后,这些参数都是常量(与坐标变量x ,y 无关),因此B 是常量阵。
当单元的结点位移 e a 确定后,由B 转换求得的单元应变都是常量,也就是说在载荷作用下单元中各点具有同样的 x ε 值、y ε值及xyγ 值。
因此3结点三角形单元称为常应变单元。
在应变梯度较大(也即应力梯度较大)的部位,单元的划分应适当密集,否则将不能反映应变的真实变化而导致较大的误差。
单元应力可以根据物理方程求得,即在物理方程中代入(14)式可以得到eexyyx Sa DBa D ===⎪⎪⎪⎭⎫⎝⎛=ετσσσ (19)其中[][]mjimj i S SS B B B D DB S === (20)S 称为应力矩阵。
将平面应力或平面应变的弹性矩阵及(18)式代入(20)式,可以得到计算平面应力或平面应变问题的单元应力矩阵。