中国有色金属学报 2003,(01),65-70 DOI:10.19476/j.ysxb.1004.0609.2003.01.011
有限体积法模拟铝型材挤压成形过程
周飞 苏丹 彭颖红 阮雪榆
上海交通大学塑性成形工程系,上海交通大学塑性成形工程系,上海交通大学塑性成形工程系,上海交通大学塑性成形工程系 上海200030 ,上海200030 ,上海200030 ,上海200030
摘 要:
研究了有限体积法模拟金属塑性成形的基本理论和关键技术 ,得到了适用于塑性成形的有限体积控制方程 ,给出了模拟过程中有限体积网格体系的建立、成形过程中金属流动的跟踪描述和时间增量步长的确定等技术处理方法。然后采用该数值模拟方法 ,对铝合金门窗型材的挤压成形过程进行了仿真。详细地分析了该零件在挤压成形过程中金属的流动情况 ,给出了成形各阶段等效应变、温度和速度等物理场量的分布情况以及整个成形过程中模具载荷 /行程曲线的变化情况。研究结果证明有限体积法是一种行之有效的铝型材挤压数值模拟方法 ,它可以为铝型材新产品的开发设计提供理论指导
关键词:
铝型材挤压 ;有限体积法 ;数值模拟 ;变形跟踪 ;
中图分类号: TG376.2
作者简介: 周 飞(1971),男,博士;
收稿日期: 2002-04-04
Simulation of aluminum material extrusion process with finite volume method
Abstract:
The basic theory and key technologies of simulation of me ta l forming processes by finite volume method (FVM) were studied firstly. The gove rning equations of FVM for metal forming were given. The simulation technologies such as mesh building, material deformation tracking and time increment computa tion were also provided. Then the extrusion process of an aluminum decoration pa rt in windows was simulated by FVM. The metal flow process was analyzed in detai ls. The distributions of effective strain, temperature and velocity during the a luminum extrusion process were discussed. The curve of load versus stroke was al so given. So FVM is proved to be an especially useful numerical simulation metho d for aluminum extrusion processes. The FVM simulation technology can help the r esearch and development of new aluminum extrusion products.
Keyword:
aluminum extrusion; finite volume method; numerical simula tion; deformation tracking;
Received: 2002-04-04
铝型材挤压有限元数值模拟方法在国内外已有不少研究
[1 ,2 ,3 ]
。 但是, 由于有限元采用Lagrange方法, 在模拟塑性成形过程时常常会碰到网格再划分问题。 多数民用建筑型材的挤压比为40~90, 铝锭坯料直径为8~300 mm, 型材厚度为0.6~2 mm, 因此, 铝型材挤压属于极度大变形塑性成形问题。 在有限元模拟过程中, 当材料流经工作带时, 不可避免地进行着一次又一次的网格再划分。 这种网格再划分不但耗时, 损失计算精度, 而且过度畸变的网格在再划分时往往出现失败, 使得模拟无法完成, 得不到流出工作带后的制品情况, 从而达不到铝型材挤压数值模拟的真正目的。
有限体积法以前多用于模拟流体的流动过程
[4 ,5 ]
。 近年来, 部分学者也逐渐将有限体积法用于模拟金属的塑性成形问题
[6 ]
。 基于Euler的有限体积法是将网格固定在空间, 材料在流动过程中Euler网格不发生变化。 因此, 用有限体积法模拟大变形塑性成形问题可以很好地避免网格再划分问题。 本文作者将对有限体积法用于铝型材挤压的数值模拟技术进行详细的研究。
1 塑性材料有限体积法基本理论
有限体积法将材料流动所要经历的空间采用Euler网格进行离散。 Euler网格是一个固定在空间的参考框架, 通常采用六面体单元。 Euler单元的节点在空间固定不动, 因此, Euler单元的体积是不变的。 材料只是从一个单元流到另一个单元, 并且材料的质量、 动量和能量也随之从一个单元流到另一个单元。 在材料流动的过程中, 必须满足质量守恒、 动量守恒、 能量守恒、 本构方程、 状态方程和热平衡方程等控制方程
[7 ]
。
1)质量守恒方程
? ρ ? t + ? ( ρ v i ) ? x i = 0 ? ? ? ( 1 )
?
ρ
?
t
+
?
(
ρ
v
i
)
?
x
i
=
0
?
?
?
(
1
)
2)动量守恒方程
? ( ρ v i ) ? t + ? ( ρ v i v j + p δ i j ? S i j ) ? x j = 0 ? ? ? ( 2 )
?
(
ρ
v
i
)
?
t
+
?
(
ρ
v
i
v
j
+
p
δ
i
j
-
S
i
j
)
?
x
j
=
0
?
?
?
(
2
)
3)能量守恒方程
? ( ρ E ) ? t + ? ( ρ v i E ) ? x i = ? [ v i ( S i j ? p δ i j ) ] ? x j ? ? ? ( 3 )
?
(
ρ
E
)
?
t
+
?
(
ρ
v
i
E
)
?
x
i
=
?
[
v
i
(
S
i
j
-
p
δ
i
j
)
]
?
x
j
?
?
?
(
3
)
上述3式中, ρ , v i , S ij , p , E 分别为材料的密度、 速度矢量、 应力张量、 静水压力和内能。
4)本构方程
? σ i j ? t = G [ ? v i ? x j + ? v j ? x i ] + λ ? v k ? x k δ i j ? ? ? ( 4 ) ε ˙ p i j = Λ ? f ? σ i j ? ? ? ( 5 )
?
σ
i
j
?
t
=
G
[
?
v
i
?
x
j
+
?
v
j
?
x
i
]
+
λ
?
v
k
?
x
k
δ
i
j
?
?
?
(
4
)
ε
˙
i
j
p
=
Λ
?
f
?
σ
i
j
?
?
?
(
5
)
式中
ε ˙ p i j
ε
˙
i
j
p
为塑性应变速率张量, Λ 为与应力、 应变以及应变速率有关的系数, f 为各向同性材料的Von Mises屈服函数。 因此, 应力张量可由下式得到:
σ ij =-pδ ij +S ij (6)
假设材料为致密体材料, 材料的密度保持不变, 则从式(1)和式(2)可得到速度变化率在有限体积单元上的积分式:
∫ V ? v i ? t d V = 1 ρ ∫ V ? σ i j ? x j d V ? ? ? ( 7 )
∫
V
?
v
i
?
t
d
V
=
1
ρ
∫
V
?
σ
i
j
?
x
j
d
V
?
?
?
(
7
)
利用高斯定理, 可以将体积积分转化为在表面积上的积分, 于是式(7)变为
∫ V ? σ i j ? x j d V = ∫ S σ i j n j d S ? ? ? ( 8 )
∫
V
?
σ
i
j
?
x
j
d
V
=
∫
S
σ
i
j
n
j
d
S
?
?
?
(
8
)
式中 V 为单元体积, S 为单元的表面积, n i 为表面法向的单位矢量。
5)状态方程
? ( ρ Φ ) ? t + ? ( ρ Φ v i ) ? x i = 0 ? ? ? ( 9 )
?
(
ρ
Φ
)
?
t
+
?
(
ρ
Φ
v
i
)
?
x
i
=
0
?
?
?
(
9
)
式中 Φ 为任意状态变量, 如体积、 应力、 压力等。 式(9)在有限体积上相应的数值积分方程为
d d t ∫ V ρ Φ d V = ? ∫ ? V ρ Φ v i n i d S ? ? ? ( 1 0 )
d
d
t
∫
V
ρ
Φ
d
V
=
-
∫
?
V
ρ
Φ
v
i
n
i
d
S
?
?
?
(
1
0
)
6)热平衡方程
根据热力学第一定律, 单元体内能量的变化率等于外力对单元体所做的功率和单元与外界进行热交换率的总和。 因此, 材料单元的能量平衡方程为
ρ c p T ˙ ? ? ? x i ( k ? T ? x i ) ? σ i j ε ˙ p i j = 0 ? ? ? ( 1 1 )
ρ
c
p
Τ
˙
-
?
?
x
i
(
k
?
Τ
?
x
i
)
-
σ
i
j
ε
˙
i
j
p
=
0
?
?
?
(
1
1
)
式中 k 为材料的导热系数, T 为工件温度, c p 为比热容。
若考虑工件与模具及环境之间的热交换, 运用Newton冷却法则, 有
k ? T ? x i n i = ? h ( T ? T w ) ? ? ? ( 1 2 )
k
?
Τ
?
x
i
n
i
=
-
h
(
Τ
-
Τ
w
)
?
?
?
(
1
2
)
式中 h 为工件与模具之间的热交换系数, T w 为模具或环境的温度。
上述微分方程可以用对时间的龙格-库塔法(Runge-Kutta time-integration scheme)来进行求解
[8 ]
, 从而得到有限体积内变形体的各种物理场量。
2 有限体积数值模拟的关键技术
2.1 有限体积网格的建立
对于三维塑性成形问题, 采用六面体单元建立有限体积的Euler网格体系。 由于Euler网格在成形过程中不随金属一起流动, 也不发生变形, 因此, Euler网格可以采用长方体或正方体单元。
在有限体积法模拟金属塑性成形中, Euler网格的建立应遵循以下两条原则:
1)Euler网格需覆盖金属流动所经历的整个区域;
2)为了保证模拟精度, 在成形工件的最小壁厚处至少应包含两个Euler单元。
上述两条原则决定了Euler网格在空间的分布及网格单元的大小。
2.2 变形体边界的跟踪描述
有限体积法中, 可以用MAC(marker and cell)法
[9 ]
、 VOF(volume of fluid)法
[10 ]
和TS(triangle segments)法
[11 ]
来跟踪材料流动的自由表面。 MAC法和VOF法都是非常粗略和近似地表达流体的自由表面, 无法得到准确的描述。 因此, 这两种方法常用来描述液体和气体的流动。 在描述金属的塑性流动过程中, 多采用TS法。
TS法是在材料的外围包裹一层类似于有限元网格的三角形几何网格——表面网格, 这个表面网格随着材料一起流动, 从而跟踪材料流动的自由表面。 为了精确起见, 每个表面网格结点之间的距离应小于有限体积网格的最小长度。 表面网格的运动位置可用类似于有限元法中表面节点的运动来确定
[12 ]
。 求得表面各网格结点的空间位置后, 即得到变形体在每一成形加载步上的边界形状。 因此, 表面网格可以很好地对成形过程中变形体的几何边界进行跟踪描述。
与有限单元法相比, TS法更加灵活, 特别是在为了更好地拟合模具表面而优化变形体表面网格质量时, 表面网格的任意一个结点可以随意地增加或删减。
2.3 时间增量步长的确定
有限体积法在控制体积上采用显式时间积分, 根据Courant准则
[7 ]
, 时间步长由应力波通过单元的最短时间来决定, 并且与单元的最小网格尺寸有关。 有限体积法模拟金属塑性成形中, 时间增量步长Δt 由下式决定:
Δt =L /(C +U ) (13)
式中 L 为有限体积网格单元的最小长度, C 为声速, U 为物质流动速度。 为了保证模拟精度, 在每一个加载步上, U 取有限体积网格中心格点中的最大速度。
3 铝型材挤压模拟实例
图1所示为一门窗锁扣铝型材产品, 其断面图如图2所示。 该型材采用平模挤压成形, 模具的工作带长度分布如图3所示。
3.1 模拟参数
材料:6061铝合金;
坯料成形温度:480 ℃;
模具预热温度:460 ℃;
坯料尺寸:d 78 mm×13 mm;
摩擦因数:μ =0.45;
挤压速度:10 mm/s;
挤压比:λ =38.19;
图1 门窗锁扣零件
Fig.1 A window lock part
图2 门窗锁扣型材断面图 (mm)
Fig.2 Profile of window lock part (mm)
图3 工作带长度分布 (mm)
Fig.3 Different length of die land (mm)
压下量:7.5 mm;
Euler网格单元大小:0.6 mm;
Euler网格数量:392万。
3.2 模拟结果分析
由于铝型材挤压成形的挤压比非常大, 难以一次成形, 所以常常采用导流模或分流模进行预成形, 以分解挤压比。 铝型材的挤压成形过程一般分为3个阶段:镦粗阶段、 预成形阶段和终成形阶段。
图4所示为铝合金门窗锁扣型材挤压过程中金属的流动情况。 成形之初, 金属在挤压头的作用下首先充填挤压筒, 该阶段的挤压成形相当于镦粗。 镦粗结束后, 金属开始流入预成形模具, 将坯料金属根据制品的形状进行预先合理的分配, 如图4(b)所示。 当金属充满预成形模具型腔后, 材料开始流向凹模工作带, 金属发生剧烈的塑性变形, 如图4(c)所示。 此时, 挤压力急剧上升。 当制品挤出一定的长度后, 成形进入稳定挤压阶段, 如图4(d)所示, 图中挤压头的压下量为7.5 mm。
图4 挤压过程中金属的流动情况
Fig.4 Metal flow during extrusion process
(a)—Stroke 0 mm;(b)—Stroke 2.67 mm;(c)—Stroke 5.00 mm;(d)—Stroke 7.50 mm;(e)—Top view;(f)—Front view
图4(e)所示为7.5 mm压下量时型材制品的俯视图, 由图中可以看出, 有限体积数值模拟所得到的制品形状与图1中产品的实际形状一致。 由此可以说明, 数值模拟可以真实地反映金属在成形过程中的流动情况。 图4(f)所示为7.5 mm压下量时型材制品的正视图。 图中显示, 制品前沿金属的流速不均匀, 这说明模具工作带长度的设计不够合理, 有待改善。
图5所示为模具载荷随挤压行程的变化曲线, 它清楚地表示出了铝型材成形挤压的3个不同阶段。 图中, 最大的挤压力可以为挤压机的吨位选择提供依据。图6、 图7所示分别为变形体在挤压成形过程中等效应变和温度场的分布及变化情况。 从中可以看出, 温度场和等效应变场的分布情况非常类似。 它表明, 变形体的温度升高主要来自于成形过程中所产生的塑性变形能。
图8所示为沿挤压方向金属流速的分布和变化情况。 由图8(c)可知, 进入稳定挤压阶段后, 制品前沿的流速为371.8 mm/s, 而根据挤压速度和挤压比计算得到的理论速度为381.9 mm/s, 二者相差2.6%, 说明有限体积数值模拟的结果与理论值吻合较好。
图5 载荷—行程曲线
Fig.5 Curve of load vs stroke
图6 等效应变分布
Fig.6 Distribution of effective strain
(a)—Stroke 2.67 mm;(b)—Stroke 5.00 mm;(c)—Stroke 7.50 mm
图7 温度场分布 (单位:℃)
Fig.7 Distribution of temperature(Unit: ℃)
(a)—Stroke 2.67 mm;(b)—Stroke 5.00 mm;(c)—Stroke 7.50 mm
图8 金属流速沿挤压方向分布 (单位:mm/s)
Fig.8 Distribution of velocity along extrusion direction (Unit: mm/s)
(a)—Stroke 2.67 mm;(b)—Stroke 5.00 mm;(c)—Stroke 7.50 mm
4 结论
研究表明, 有限体积法不但可以实时地跟踪描述塑性成形过程中金属的流动情况, 而且可以给出变形体在成形过程中物理场量的分布和模具载荷的变化情况。 有限体积法模拟计算的结果与理论计算以及实际生产的结果吻合较好。 因此, 有限体积法是继有限元法之后的又一种有效的金属塑性成形数值模拟方法。
对于铝型材挤压这类极度大变形塑性成形问题, 采用有限体积法进行数值模拟可以解决有限元法无法避免的网格再划分问题。 它一方面减少了由于网格再划分带来的精度损失; 另一方面, 由于有限体积法无需网格再划分, 它还可以保障模拟过程能够顺利地完成。 因此, 有限体积法特别适用于模拟复杂薄壁铝型材的挤压成形问题。
参考文献
[1] ChandaT ,ZhouJ ,KowalskiL ,etal.3DFEMsimula tionofaluminumextrusionintheinitialnon steadystate[A].AdvancedTechnologyofPlasticity,Proceedingsof6thICTP[C].1999,3:18131818.
[2] 周 飞,彭颖红,阮雪榆.铝型材挤压过程有限元数值模拟[J].中国有色金属学报,1998,8(4):637643. ZHOUFei,PENGYing hong,RUANXue yu.Numer icalsimulationofaluminummaterialextrusionprocessbyfiniteelementmethod[J].TheChineseJournalofNon ferrousMetals,1998,8(4):637643.
[3] HAONan hai,LIKe zhi.Numericaldesignofdielandforshapeextrusion[J].ChineseJournalofMech anicalEngineering,2001,14(1):9193.
[4] AndersonDA ,Tannehill.ComputationalFluidMechan icsandHeatTransfer[M ].NewYork:HemispherePub lishingCorporation,Taylor&FrancisGroup,1984.321508.
[5] NaotD ,Radi.Numericalsimulationofsecondarycur rentsinchannelflow[J].JournalofHydraulicDivision,ASCE ,1982,108:948968.
[6] deVriesE ,DingP .Simulationof3Dforgingandextru sionproblemsusingafinitevolumemethod[A].Proceed ingof17thMSCJAPANUsersConference[C].Japan,1999.155161.
[7] VersteegHK ,MalalasekeraW .AnIntroductiontoComputationalFluidDynamics-TheFiniteVolumeMethod[M ].LongmanGroupLtd,1995.5137.
[8] 沈永欢,梁在中,许履瑚,等.实用数学手册[M ].北京:科学出版社,2001.726732. SHENYong huan,LIANGZai zhong,XULu hu,etal.AppliedMathematicManual[M ].Beijing:SciencePress,2001.726732.
[9] SwaminathanCR ,VollerVR .Atime implicitfillingal gorithm[J].ApplMathModelling,1994,18:101109.
[10] Hirt,NicholsCW .Volumeoffluid(VOF)methodforthedynamicsoffreeboundaries[J].JournalofCompu tationalPhysics,1981,39:201225.
[11] SlagterWJ,FlorieCJL ,VenisACJ .Advancesin3Dforgingprocessmodeling[J].JournalofComputation alPhysics,1991,92:82105.
[12] 彭颖红.金属塑性成形仿真技术[M].上海:上海交通大学出版社,1999.5497. PENGYing hong.SimulationTechnologiesofMetalForming[M ].Shanghai:ShanghaiJiaotongUniversityPress,1999.5497.