编辑: hyszqmzc | 2019-09-21 |
6 No.
1 M ar.
2006 收稿日期 : 2005-
08 -
10 基金项目 : 国家自然科学基金项目(50175110);
高等学校博士学科点专项科研基金项目(2002533007) 作者简介 : 姚松(1975 - ), 男 ,湖北公安人 , 中南大学工学博士研究生 , 从事车辆结构分析和撞击动力学研究 . 文章编号 : 1671- 1637(2006)01 -
0016 -
04 DKQ 弯曲单元的构造及应用 姚松,田红旗 (中南大学 轨道交通安全教育部重点实验室, 湖南 长沙 410075) 摘要: 为了研究复杂载荷作用下薄板结构的受力 , 基于离散 Kirchhoff 原理, 推导了三角形弯曲 薄板单元的应变矩阵 、 刚度矩阵中的显式面积坐标积分方法 ,根据变分原理构造了任意四边形弯曲 单元 DKQ , 引入对称及反对称畸变模式进行网格畸变敏感性分析 ,将DKQ 弯曲单元与平面应力 单元组合, 得到了用于薄板分析的四边形平板单元, 并编制了有限元分析程序 .计算结果表明 : 相 对偏移量从 - 0.
18 变化至
0 .
18 ,反对称畸变模式下挠度最大误差为
1 .
83 %, 而对称畸变模式下挠 度最大误差为
0 .
99 %;
对某地铁车体结构 ,计算结果与 ANS YS 结果误差在
3 .
5 % 之内 ,这说明构 造的 DKQ 弯曲单元对网格畸变不敏感 ,具备良好的位移解收敛性和计算精度 . 关键词: 车辆工程;
结构分析;
有限元;
四边形弯曲单元 中图分类号 : U270.
12 文献标识码 : A Construction and Application of DKQ Bending Element YAO Song , TIAN Hong- qi (Key Laboratory for Traffic Safety on T rack of Ministry of Education , Central South University , Changsha 410075, Hunan, China) Abstract: In order to analyze the states of thin- walled structures under complicated loads by discrete Kirchhoff principle, the strain matrix of triangle bending element and explicit area coordinate integral method were discussed in detail. According to variant principle , DKQ(Discrete Kirchhoff Quadrilateral) bending element w as constructed. By introducing symmetric and antisymmetric mesh distortion modes , mesh distortion sensitivity analysis w as carried out. Through combining DKQ element w ith planar stress element , quadrilateral plate element used for general thin - w alled structures analysis w as implemented by finite element program. Numerical test result indicates that w hen relative offset value varies from -
0 .
18 to 0.18 , the deflection maximum errors of antisymmetric and symmetric modes are 1.
83 %and 0.
99 % respectively. For a car- body of some subway vehicle, the differences between the analysis values and those obtained by ANSYS are less than
3 .
5 % . The result shows that DKQ bending element is insensitive to mesh distortion and has good convergence capability and accuracy for displacement field.
2 tabs ,
6 figs ,
10 refs. Key words: vehicle engineering ;
structure analysis ;
finite element ;
qudrilateral bending element Author resume: YAO Song (1975- ), male , doctoral student , 86- 731-
2655294 , dynacn @ 126. com.
0 引言从力学角度分析 ,轨道车辆 [ 1- 2] 薄板结构所受到 的载荷可以分解为平行于薄板中面的载荷和垂直于 薄板中面的载荷, 应力分析由平面应力分析和薄板 弯曲应力分析叠加而成 .早期的弯曲板单元大多基 于经典的薄板理论 ,在以该理论为基础的板单元的 能量泛函中,包含了挠度的二阶偏导数 ,根据收敛准 则,要求位移在单元交界面上必须保持挠度为 C
1 连续,这给构造板单元带来了极大的困难.由此研究 人员将注意力转向了基于 M indlin- Reissner 理论的 Mindlin 板 ,对挠度 w 、法线转角 θ x 、θ y 采用各自独 立的场函数 ,利用罚函数的方法将 w 、θ x 和θy之间 的约束条件引入能量泛函 ,从而将 C
1 类连续转化为 C0 连续性问题 .该板单元对中厚板很有效 ,当板逐 渐变薄时 ,单元刚度矩阵中的剪切项占主导地位, 计 算出的弯曲变形远小于实际变形 ;
当板非常薄时求 得的位移趋向于零 , 从而产生了 剪切闭锁 现象 .
1980 年Batoz[ 3] 通过挠度和转角分别独立插值 , 然 后在单元的若干个离散点上强迫挠度 w 与法线转 角θ x 、θ y 满足经典薄板理论中的 Kirchhoff 约束, 提 出了基于离散 Kirchhoff 假设的三角形 DKT [ 4- 7] (Discrete Kirchhoff Triangle)弯曲板单元 .该单元 为C0 型单元 ,克服了 剪切闭锁 现象 .多个文献 表明 DK T 弯曲单元在求解薄板弯曲问题时都显示 出良好的性能, 具有较高的精度 .但在求解平面问 题时 ,由于三角形单元为线性位移模式 ,为常应力单 元,而四边形单元为双线性位移模式,在同样网格密 度下 ,四边形单元能够获得精度更高的数值解 ,因此 有必要构造任意四边形弯曲单元 ,便于与四边形平 面应力单元组合成四边形平板单元进行结构分析 .
1 DKT 单元 图1为DKT 三角形弯曲板单元, 其单元节点 位移向量 ae 由
1、
2、3 三个角点的挠度 wi 、 转角 Χ xi 和Χyi 组成 ae =(a T
1 , a T
2 ,a T
3 ) T ai =(wi , Χ xi , Χ yi )T Χ xi = w / x | i Χ yi = w / y | i (i =1 ,
2 ,3) 弯曲单元挠度场的构造直接通过三角点的挠度构造 w = ∑
3 i =1 Niw i
4 、
5 、6 三个中间节点仅有转角自由度, 并引入如下 假设 : 单元每边的法向转角 Χ n 为线性分布 ,例如横 边上 Χ n4 =(Χ n1 +Χ n2 )/2 Χ n = w / n 单元边界上切向转角为 Χ s = w / s 由角点挠度和角点切向转角所确定的三次式决定, 例如横边 Χ s4 = - 3(w1 - w2 )/(2l12 )- (Χ s1 +Χ s2)/4 图1DKT 弯曲单元 Fig.
1 DK T Bending Element 基于以上假设可 以得到 : 由于单元各边 界上切向转角完全由 所在边角节点的切向 转角和挠度的三次式 插值得到 ,法向转角由 所在边角节点的法向 转角线性插值得到 , 所 以相邻单元公共边界 上的切向转角和法向转角是完全协调的.单元转角 场构造时利用了单元
3 个中间节点
4 、
5 、
6 的转角 Χ (x , y) = ∑
3 i=1 Ni Χ (x, y)i + ∑
6 m=4 Nm Χ (x , y)m 根据函数梯度的定义可以得到 w / n w / s = cos(n, x) cos(n, y) cos(s , x) cos(s , y) w / x w / y 通过复杂的推导可以将单元
3 个中节点的转角 Χ (x , y)m 表达成
3 个角节点挠度 、转角 wi 、Χ xi 、Χ yi 的 函数.这样中节点的转角自由度被 凝聚 掉 ,由此 转角自由度的引入不仅进一步提高了单元位移场的 插值精度,同时还保证了单元间转角的协调性. 在得到 Χ x 、Χ y 的插值表达式以后 ,可以根据下 式得到单元广义应变 κ= - Χ x / x Χ y / y Χ x / y + Χ y / x =B ae 式中: B 为关于面积坐标 L
1 、L2 、L3 的矩阵 ,单元弯 曲刚度矩阵为 K =(t3 /12) BT DBdA 根据面积坐标幂函数在三角形全面积上的精确 积分公式 L k
1 1 L k
2 2 L k
3 3 d A = 2A(k1 ) (k2 )(k3 ) (k1 +k2 +k3 +2) 由于B 中关于面积坐标Li 的最高次数为
1 次, 可以 将B矩阵进行如下分解 B =C +B1 L1 +B2 L2 +B3 L3 经过这样的分解以后得到 B1 、B2 、 B3 均为常量 系数矩阵,将其代入单元弯曲刚度矩阵中 ,然后利用 幂函数精确积分公式即可完成单元弯曲刚度矩阵的 显式数值积分,该数值积分是完全精确的 .
17 第1期姚松,等: DKQ 弯曲单元的构造及应用
2 任意四边形弯曲单元 DKQ 本文基于三角形弯曲板单元 DKT 构造了任意 四边形 DKQ(Discrete Kirchhoff Quadrilateral)弯 曲板单元 ,并对这种弯曲板单元的收敛性和对网格 畸变的敏感性进行了详细探讨 . 图2说明了任意四边形弯曲单元的构造过程 , 图2DKQ 的构造 Fig.
2 Construction of DKQ 沿四边形的两条对角 线将其划分为 ①、 ②、 ③、④四个 DK T 三角 形弯曲板单元 , 因此四 边形单元的内力虚功等 于4 个三角形单元的内 力虚功之和的一半.根 据虚功原理 ,DKQ 的单元刚度矩阵可以通过如下方 式构造 K - a - e =
1 2 (K1 ae1 +K2 ae2 + K3 ae3 +K4 ae4 ) 式中 : a - e 为DKQ 弯曲单元的节点位移向量 ;
aei 为4个DK T 弯曲单元的局部节点位移向量, Ki 为4个DKT 弯曲单元的局部弯曲刚度矩阵 ;
K - 为DKQ 弯 曲单元的刚度矩阵.从上式可以看出, 在编制有限 元程序时 ,由三角形单刚组集成四边形单刚时根据 其所对应的节点位移编号来进行, 容易实现模块化 .
2 .
1 DKQ单元的收敛性分析 通过上述方式构造的四边形单元能否具备好的 收敛性是决定其能否应用于车辆结构分析的关键 , 采用 Fortran[ 8] 语言编制了 DKQ 弯曲单元的有限 元程序 ,分别计算了边长为 L 的简支和固支方板在 中心集中载荷 P 作用下的最大挠度;
采用四分之一 模型与规则网格计算, N 表示二分之一的边长上的 单元数;
为了便于与参考文献[ 9] 中的理论解进行比 较,文中引入挠度系数 C ,定义如下 C = D0 wmax PL
2 D0 = Et
3 12(1 - u
2 ) D0 为板的弯曲刚度 ,表1列出了模型厚跨比 t /L 为0.05 时的挠度系数的数值解 . 当采用规则网格, 沿四分之一结构边长划分为
5 个单元时 , 简支条件下最大挠度与解析解仅相差
0 .
84 %,固支条件下最大挠度与解析解也仅相差
1 .
67 %.从图
3 的曲线可以看出, 随着网格的加密 , 采用 DKQ 单元的有限元数值解与理论解之间的误 图3C误差与 N 的关系曲线 Fig.
3 Curves of C Difference and N 差呈负幂指数下降 , 说 明该单元在规则网格 条件下具有较好的收 敛性, 在相同网格条件 下 ,简支边界条件的误 差要小于固支边界条 件的误差. 表1C的数值解 Tab.
1 Numeric Solution of C 理论解[ 9] 简支 :
0 .
011 60 固支 :
0 .
005 61 N 数值解 数值解
2 0 .
012 063
0 .
005 979
3 0 .
011 833
0 .
005 822
4 0 .
011 744
0 .
005 745
5 0 .
011 698
0 .
005 704
2 .
2 DKQ 单元网格畸变敏感性分析 在实际进行车辆结构有限元分析时, 由于结构 几何复杂性 ,划分的网格往往是不太规则的, DKQ 弯曲单元能否在网格畸变程度不太高的情况下保证 计算的精度是其中的关键 ,因此在本文中通过引入 两种简单的畸变模式(图4),探讨了 DKQ 弯曲单元 对网格畸变的敏感性.计算时所取的网格参数 N 为2,Δ表示内部节点与中心偏离的距离, 图4中所 示的偏离方向定义 Δ为正值. 图4反对称与对称畸变模式 Fig.
4 Antisymmetric and Symm etric Distortion Modes 根据分析,从图
5 看到 ,当内部节点的相对偏移 量Δ/(0. 5L)从-0.18 变化至 0.
18 时 ,与规则网格 相比,反对称畸变模式下挠度误差没有超过
2 %,而 对称畸变模式下挠度误差在
1 %以内, 说明构造的 DKQ 弯曲单元对网格畸变不敏感 .
3 平板单元 在进行车辆结构分析时, 结构的受力状态不再 是纯弯曲,同时可能还有面内载荷,因此要将四边形 弯曲板单元和四边形平面应力单元组合起来形成平 板单元 [ 10] .在局部坐标系中 ,节点位移参数本不包 含θzi , 但是为了进行总刚集成时避免共面单元组集
18 交通运输工程学报2006 年 后出现 θ zi 方向刚度为零的情况,在组合过程中给 θ zi 方向赋予单位刚度, K(m) 为单元平面应力刚度矩阵 , K - 为单元弯曲刚度矩阵 ,组集过程为 Kp = K (m)
0 0
0 0
0 0........