|
|
发表于 2008-10-26 19:07:26
|
显示全部楼层
来自: 中国四川攀枝花
在网上搜得一篇论文,相信对你有所帮助
8 F j# i2 @" F) Y* `橡胶复合材料结构非线性有限元分析
I8 \6 m: l/ z" J# [: P
. _$ |1 \) y2 ~- }2 e危银涛 杜星文 王友善* 吴宝国5 d( N; s0 f% x
(哈尔滨工业大学复合材料研究室, 哈尔滨 150001)8 w% h. X0 q/ d; t. k; t
(*桦林集团总公司, 牡丹江 157032)
4 k$ }% p4 _: X4 z6 o: @
7 Z8 A, r6 S3 \: H摘 要 用非线性有限元法对一种典型的橡胶复合材料结构──轮胎的力学性能进行了分析。给出了轮胎滚动边界值问题的数学描述和有限元列式及相应的求解策略。取旋转刚体构形为参考构形,得到消去时间变量的惯性场和控制方程。用增量约束方法处理轮胎滚动接触问题可以达到收敛快、精度高。研究了滚动参数对轮胎总体和局部变形和受力的影响。计算结果与已有数据和试验相吻合,对轮胎设计和汽车动力学有指导意义。
: _* r/ e) M* _6 M; ^% T+ h关键词 橡胶复合材料, 轮胎, 滚动接触, 约束增量, 非线性有限元7 N; V0 ~5 J& E3 I5 @
中图分类号 TB332, O242.21, V255.3' h/ {1 D8 S" R
; Y4 K' M8 Q) \ 轮胎是一种典型的橡胶复合材料结构,它是汽车最重要的部件之一,支撑汽车的重量,传递汽车的牵引制动力矩,并起缓冲隔振等作用。轮胎力学性能的好坏直接关系到汽车各种性能的好坏。由于轮胎结构和材料的非均质各向异性,载荷的复杂性,以及大变形,轮胎力学性能的分析一直是十分困难的工作,轮胎力学经历了从简单模型分析、网络分析到薄壳、厚壳分析,三维分析的过程。尽管壳单元可以大大减少分析的自由度,而且各种高阶层合理论及其相应的有限元法也可以提供较精确的横向变形及应力分析,但在边界条件的处理及提高精度方面壳单元有难以克服的困难﹝2﹞。目前对于轮胎静态分析已有了较成熟的方法和结果﹝1,6,7﹞,然而对于滚动接触问题分析尚缺乏有效分析方法﹝2,3﹞,很少有文献发表。
+ y: t, }' _$ O) O, S 本文作者根据几种新采用的技朮,考虑轮胎的稳态惯性场,接地面内的粘着-滑动效应,用三维非线性有限元方法分析了轮胎滚动接触问题,主要内容包括:
: ?* V/ i8 m4 q& b! | (1) 基于中介构形的滚动轮胎惯性场和控制方程的表述﹔5 K1 }- ^. `+ I/ f2 ^* A
(2) 轮胎滚动接触问题的增量约束处理方法﹔: h1 l; R) B9 u1 n$ }
(3) 发展了三维等参元分析轮胎滚动力学性能﹔
/ w: [, S$ N7 X6 J% ?0 | (4) 非对称有限元方程组的波前解法﹔
* m% z9 R9 B. h X0 ? (5) 轮胎有限元网格的自动划分及优化。! f, d8 s' T. o; c, W w! G
基于本文的方法,对某轮胎厂生产的轿车用子午胎给出了数值算例,与已有数据和试验相比表明本文的方法是可靠有效的。
/ O4 g/ c5 {' n6 N; [, M, l2 t" d8 k% C0 v; O
1 轮胎结构与载荷
& d" ?) K; a& e' y+ u4 }. R9 J, a$ p( M4 L* w1 s! U/ ]; M
图1 半钢丝子午胎断面结构及材料& ~( U S" z" k( H0 z2 a/ h
Fig.1 Section structure and materials' |) N4 u; p9 W3 w& p+ n V* |) p6 e
of semi-radial tire
1 w! \2 M* z$ S- w7 s7 X+ _/ g7 j 轮胎是由橡胶及帘线增强橡胶复合材料构成的复合结构体(见图1),即使在一般情况下,轮胎材料也要发生较大的变形及转动。轮胎载荷包括内部空气压力,与地面的接触载荷。接触问题在本质上是非线性的,接触区域大小和形状以及接触载荷大小事先是不知道的,轮胎结构分析是一个各向异性,非均质,几何非线性,材料非线性及边界非线性问题,分析起来困难较大。随着有限元单元法,复合材料力学,计算方法以及计算器科学等相关科学发展,在轮胎结构分析已经可以比较精确地进行,带动了轮胎设计理论的发展。本文中将帘线-橡胶复合材料看成是横观各向同性弹性材料,材料特性均由试验确定,橡胶则被认为是各向同性弹性材料。
- _0 J1 ~: h0 C H$ V/ ]& i4 Y
5 ^1 s" H5 ?( F/ [% t, m7 C2 滚动轮胎结构分析边界值问题 F: n: `6 q1 r- g4 m/ l
2.1 大变形情况下的平衡方程和虚功方程
' a5 @" u6 D% ~5 \2 q; e Green应变和Second-Kirchhoff应力定义如下:
4 v4 c k; ^3 }- N0 D; `1 W; U' F: I) H/ ~8 @1 N- O
大变形情况下的平衡方程和虚功方程为
! q$ f/ M: p3 i, x5 n! D
& f+ |$ Q8 s; f1 w# u) @2.2 牵引滚动轮胎结构分析的边界值问题描述
2 P* ^( h% N. A$ z8 U( n3 u! Z2 V
0 |7 }( s# I" y4 Y7 c' V, N2 ]2.3 中介参考构形──移动的Lagrangian构形
9 }" b7 T& }/ Y0 c9 { 轮胎是在地面上不断滚动的,以初始构形为参考构形分析要求必须追踪轮胎的所有的滚动过程,在计算量上是不可接受的。为此我们引入轮胎刚体滚动构形作为参考构形(如图2所示),从而将轮胎的运动分解为刚体滚动和纯变形﹝2,6﹞。基于此,滚动轮胎的速度场和加速度场可以表述如下:
% m b* t V s# Q- \, V; I$ H
4 Q( U! P6 L! Z% g0 v) w4 X" z 参考构形y随时间的变化关系可以表示为
: H- v) D( ?" o6 h p
/ t* G' X+ g5 ~ 将式(4)、(5)代入式(3)即得到消除时间变量的以参考构形为基础的稳态速度场和加速度场,从而我们可以对轮胎稳态问题进行拟静态分析。) y! O2 k2 X$ E* g' z
. \0 a: m! J) J% Y
3 有限元列式
5 p7 g0 r/ c9 S1 B% Z 将单元应变、应力、节点位移、形状函数表示如下:9 Y6 M5 x# {1 u# j9 B. Y
8 m5 g% X6 w4 A E% K$ r; B 其中形状函数上标1、2、3表示其行数,下标表示对应的节点,则有
7 z! R( m# G. }1 k* e# I
1 y4 x- P/ @% ?% y+ w 稳态加速度的分量形式为
; e5 F0 {( b) [ v0 h2 u0 h: S+ z. _3 s7 i( v% H
上式中第一项为相对加速度,第二项为向心加速度,第三项为科氏加速度。离散形式的虚功方程为* o. n% Y3 a, {& X$ _$ k8 S1 i
& x+ y. I/ N3 \6 v+ E% d; B将式(8)代入(9)并整理得到
& f6 ?0 Q, ]7 _3 [式(10)是非线性的,用Newton-Raphson方法对其线性化8 g$ S0 ~/ ]+ ^# \% c! M8 Z
: A% z5 |0 r! _9 P; {, A; R其中各矩阵定义如下:+ C9 o% W+ E2 v/ `, ?& s- [
& B' @' A5 ^& ? 矩阵KL、KN、KS的定义见文献﹝6﹞,由KI的定义,显然它是非对称的,从而导致总体切线刚度矩阵的非对称,需要特别的求解器。; ]4 B7 F& e8 e# d
* t3 S, q2 z9 _4 _* B9 Z- \, V4 接触问题求解方法
; A# l2 U3 G& o. a, O: }/ v' O+ @1 ]9 j9 r! J8 c7 T
图2 初始约束赋值2 V# u8 }0 v4 G
Fig.2 Schematic of contact algorithm by' Y0 T6 \" E2 S
variable constraint method0 a% X/ z8 Y, I' j7 C) G+ U
采用由本文作者之一所提出的可变约束法处理接触问题,其详细内容参见文献﹝1,6,7﹞,在此不再赘述。图2是赋初约束示意图,其中Uc2=Yc2-Y2是纵向(y2方向)的约束,先预设一接触区,并赋予落入接触区内的所有点上述约束,在每步迭代后对各约束点进行约束力和位移校核,如果约束是合理的,则其约束增量为零,对自由点进行位移校核,如果其落入接触区,则赋予其增量约束﹝6,7﹞。将接触点分成两类:驻定的和滑动的,如果约束切反力超过最大摩擦力,则将驻定接触点释放成滑动点,并将摩擦力化成主动力进行计算,否则保持其粘连。对整个过程进行迭代直到满足收敛条件,本文中所用收敛判据是最大位移判据。& O; h" O$ r, K; G
& g x9 N' [% H3 J% ~6 p; z
5 滚动轮胎三维有限元模拟2 A8 W. S0 q% \% ]
将轮胎结构用六面体八节点等参元(五面体等参元作为过渡单元)仿真。注意到在微分算子Φ()中包含形状函数对参考坐标的一次及二次导数,在等参元的列式中,必须将其转化成对自然坐标的导数,一次导数的转化借助雅可比矩阵很容易实行,下面导出二次导数的转化公式1 x# J: j) C* ^$ w
设
/ p7 V" k8 G. |+ y$ U9 g" C, m$ g
}5 J, N! B1 D3 |3 \' g其中:﹝J1﹞、﹝J2﹞可由形状函数对自然坐标的导数及节点坐标求得,由式(16)有2 m# p2 x# L0 C3 e* M# |
2 @* z1 f* e$ Q% x0 k& x9 |3 ?+ z- S
式(17)右端的一阶导数和二阶导数均是对单元自然坐标,从而得到形状函数对参考坐标的二阶导数,该变换式对位移也适用# q4 w( }; ]9 I& }& F9 K5 f8 I
- ?2 F% v+ A* g2 u9 {六面体等参元的形状函数及其一次和二次导数的值列于表1。经过以上的变换,就可以建立起式(10)所表示的有限元总体列式。5 v6 B7 r, {( W. r
表1 形状函数的一次和二次导数
, ?& e) M- h0 e/ |0 f A- R4 L. f8 ?. m/ a) \. d% A4 D
6 非对称有限元方程组的波前解法及网格优化
3 I/ n# t( d4 @- K/ i' M) N! u6.1 非对称有限元方程组的波前解法# }; |5 w7 J" R: e, ?
0 c( `2 p! |- K2 x5 B+ L
(a) 对称 (b) 非对称
1 [8 B( ]4 _9 j# R图3 对称与非对称矩阵波前区, p" k3 i8 P; v2 T9 F% M
Fig.3 Wave-front area for symmetric
; t/ v# C# L2 ?7 j(a) and nonsymmetric (b) matrix! f) k- R! A0 H/ i
如前所述,本文中所建立的有限元总体矩阵是非对称非线性的,采用Newton-Raphson方法可以将非线性矩阵线线化得到(11)式从而迭代求解,在每一迭代步内采用波前法﹝10﹞1。波前法是适用于中小型微机求解大型有限元方程组的常用解法,是对高斯消去法的改进,主要特点是集成和消元交替进行,从而有效地节省了内存。非对称刚阵波前解法与对称刚阵波前解法并无原则的不同,仅在波前区的存储,方程的消元和主元的选择上稍有差异。如图3所示,对于对称矩阵仅需存储波前三角形中的元素,并可用一维存储的方法,而对非对称矩阵消元所及的行和列均需存储。2 _/ n" t$ q) ^/ n
6.2 网格自动生成及优化9 A' V* m- E5 F+ I1 @( D- k1 h3 J' z; ~
波前法不受节点编号的影响,但受单元编号影响很大。本文中在自动网格生成之后,根据波宽极小的原则将单元编号优化,消元次序按元素新编号控制。轮胎断面网格自动生成是先沿经向,再沿径向,然后再周向(见图4、5)展开。优化过程则是先沿周向,再沿径向,然后是经向方向形成。
( \9 Y1 `4 T( u" v; h1 P
" K+ ^' c7 p2 o) g2 k2 V. z图4 轮胎断面有限元网格 图5 轮胎总体网格变形图
: C) _8 X0 d: A- u: p/ o9 H; HFig.4 Tire section FEM mesh Fig.5 Tire total deformed mesh ; m% V' e6 _- y' l
' O7 B0 v) l- f4 j
! J/ ?. f& {) t! i+ l7 数值结果与结论0 \/ a2 A' ^# N! }& ~
对145R12LT子午胎进行了静载和自由滚动分析,其断面和总体见图4、5。每个断面划分140个单元,167个节点,总共20个断面,2800个单元,3340个节点。计算结果分为三类:(a)变形,(b)应力应变,(c)接地面信息,包括接地面积和接地压力分布。
! B4 z. M/ @7 n9 E: Q! @; e 图5是下沉量为30mm时的轮胎总体变形图,为节省篇幅没有给出各断面的变形情况。该图直观地显示了轮胎在受到接地载荷及稳态滚动情况下的变形。每一点的变形,最大变形等信息均可从计算结果中获得。这些信息可为轮胎轮廓设计、结构设计与改进提供依据。
# P+ {1 q* i: }% U 图6是第一层带束单元8(胎肩)周向正应力随滚动速度和周向位置变化曲线,可以看出随着轮胎滚动一周带束胎肩经历了受压至接触面内的受拉再受压的过程,这种循环变化的应力是导致疲劳破坏的主要原因,而且该应力幅随着速度的增加而增加,轮胎的设计应尽量减小这种应力幅。本文中首次用有限元方法显示了带束胎肩应力的这种变化,为轮胎设计提供了科学依据。图7是胎体单元15(胎侧)经向正应力随滚动速度和周向位置变化曲线,可以看出,胎侧胎体在接触面内受力减小,所受应力被松弛,即胎体由于充气所收的张力在接地面内被松弛,这与一般的结论相吻合﹝4﹞,这种张紧-松弛的应力变化是导致胎侧疲劳的主要原因,本文中用有限元方法首次定量显示了这种应力变化。4 d R6 N3 T. Y( M9 ~5 _
图8显示接触摩擦力变化,该图表明自由滚动摩擦力呈反对称分布﹝5,9﹞,胎肩处的摩擦力要大于胎冠中心。这是导致磨胎肩的原因。图9显示接地压力变化。从图9可以看出,在下沉量较大时胎肩处的接地压力大于接地中心﹝5,9﹞,这是导致胎肩摩擦力大的原因,所以轮胎设计应使接地压力分布尽量平坦。
) V$ U, r; K$ T, _" i: m+ ?( Z' r$ w 图10是第一层带束周向正应力在横截面内随经向位置变化曲线,图10表明带束在胎冠受拉力最大,胎肩最小。图11是胎体经向正应力在横截面内随经向位置变化曲线,图11表明胎体受力在胎圈处有显著变化﹝1,6,7﹞,在接触面内胎体受力变化复杂。
; K+ K7 R, a5 v# g) H: o8 j' d 值得指出的是本文作者所编制的程序可以提供每一点的变形、应力应变和接地面信息及这些量随操作条件而变化的特性。由于考虑了稳态惯性场和可变约束条件,加上采用三维等参元、Newton-Raphson方法和波前法,使得本文的分析结果更加全面、精确。
; i& V k: E& i$ b6 a8 i 轮胎的试验特别是动态试验是件特别困难的工作。因此,精确的有限元分析技朮才显得更加重要。本文的分析结果大部分是和其它分析结果对比,对轮胎的变形轮廓进行了试验校核。由于动态试验的困难,对转速为零即静态时变形轮廓进行了校核,具体的试验过程和方法详见文献﹝6,7﹞。计算与试验相比,在轮胎的最大位移处胎侧的最大位移误差为3/100,表明本文的方法是有效的合理的。图6、7中横坐标为周向角度,纵坐标为应力(Mpa),图10、11中横坐标为经向位置,纵坐标为应力(MPa),NS为断面号,图8、9横坐标为周向角度位置,纵坐标为接触力(N),下沉量为30mm。
9 M, ]8 P* `( A; G0 _/ z. g( g! _
图6 带束单元8周向正应力随速度、位置变化 图7 胎体单元15经向正应力随速度、位置变化
7 U: Y8 Z9 t0 m) ^Fig.6 Variation of circumferential normal stress of belt element 8 with velocity and positions
% ]3 x2 F" r7 U, G7 t, P( ]# iFig.7 Variation of meridian normal stress of body element 15 with velocity and positions
$ E' k7 g3 J) l% w
/ F- a# u# R% @; Q2 K图8 接触摩擦力变化图(自由滚动,下沉量30mm) 图9 接地压力变化(自由滚动,下沉量30mm)& o4 `+ a& t. I# p& P: O
Fig.8 Variation of contact friction forces (freely rolling, Δu=30mm)
+ `8 n7 w3 N7 }- O* P; P" E. G# v* DFig.9 Variation of contact pressure (freely rolling, Δu=30mm)
* k3 a- X" }' X1 f" C2 E3 l3 @! T! V% |0 T
图10 第一层带束周向正应力随经向变化 图11 胎体经向正应力随经向变化- Z$ V1 n0 y4 m E8 d# h
Fig.10 Variation of circumferential normal stress of first belt. f/ Q2 I. J' i5 G' E
Fig.11 Variation of body meridian normal stress' C @/ f9 y& i. b& V# ?
" a( u' l/ x) C' c0 @1 j
8 总结与展望% Y1 w0 k7 O3 \# B6 q
本文中简洁地描述了橡胶复合材料轮胎结构滚动接触问题,使得复杂问题可以分析求解。轮胎的非线性有限元分析是现代轮胎力学的主要内容,它对轮胎的优化设计,及汽车动力学来说都非常重要。现代汽车工业对轮胎的要求越来越高从而对轮胎力学分析提出了更高的要求,尽管国内外学者在这方面作出了很多努力,仍然有许多工作要做,包括: 更精确地仿真轮胎的构形和结构, 更真实的轮胎材料模型, 更强的模拟轮胎所受各种载荷的能力, 统计模拟、损伤传播与失效预测以及车-胎-地一体化分析。本文的工作是这种努力的一部分,利用本文提出的方法和编制的程序可以较为精确地分析滚动轮胎内的受力、变形、接地面内的压力分布和摩擦力分布。部分结果已应用于工程实际,更精确的分析需要更强大的计算器支持和对轮胎材料(橡胶及其复合材料)力学性能的更深入的了解,以及更有效的计算方法。 |
评分
-
查看全部评分
|