《环件轧制三维有限元模拟中质量缩放方法的运用.docx》由会员分享,可在线阅读,更多相关《环件轧制三维有限元模拟中质量缩放方法的运用.docx(10页珍藏版)》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、环件轧制三维有限元模拟中质量缩放方法的运用 钱东升华林左治江袁银良摘要:阐述了环件轧制三维有限元模拟中运用质量缩放方法的意义及其理论依据、并对显式 动力学有限元模拟中各种质量缩放方法进行了分类,给出了环件轧制模拟中质量缩放方法的 选取原那么。最后基于ABAQUS/EXPLICIT操作平台,建立了铅环件的环件轧制三维有限元 模型,通过具体模拟及分析比拟计算结果对选取原那么进行了验证,从而总结出了适合于环件 轧制模拟的有效的质量缩放方法。关键词:环件轧制有限元模拟质量缩放中图分类号:TG1刖百大变形动态非线性问题(nonlinear problem of dynamic large deforma
2、tion) 1的计算和分 析可以说代表目前有限元分析的最高水平,也是目前还在处于迅速开展的领域,其包括的方 面主要有:板料成型工艺的数值模拟(汽车覆盖件),汽车平安性碰撞的全过程计算,材料 的体积成形工艺(挤压、锻造、轧制)的全过程数值模拟等,环件轧制的三维有限元模拟就 属于其中的一种材料体积成形工艺的全过程数值模拟。对于此类问题中存在的几何非线性和 材料非线性问题,用于计算求解的方法主要有两种:静力隐式算法和动力显式算法。静力隐 式算法的特征是迭代计算,在每一时间步内需反复迭代,迭代收敛性会受许多因素的影响, 在计算中须调整迭代以满足收敛,而且其计算时间随模型单元数量呈指数增长,所以其计算
3、时间长。而动力显式算法是递推计算,无需迭代,就时间步长直接进行递推计算,计算时间 随单元数量呈线性增长,计算时间短,而且可以通过质量缩放来缩短计算时间,所以对于大 规模计算和高度非线性问题2,动力显式算法总体效率要比静力隐式算法高,优势明显。近年来,随着在板料成型及碰撞冲击领域得到口益广泛的应用,对于材料的体枳成形工艺领 域,人们也都开始采用动力显式有限元法对其进行过程模拟研究3,而且取得了较好的成 效。动力显式算法在有限元分析中通常被用来解决两类典型的问题:结构动力学响应分析和准静 态分析(其中包括复杂的非线性效应和接触条件),而环件轧制的模拟那么是其中的一种准静 态分析。由于动力显式算法采
4、用的中心差分法是通过对时间求导进行计算,引入了质量矩阵 4,所以在准静态分析中用于平衡方程式中的离散质量矩阵对计算效率和精度都起着至关 2 Belyschko T et al. Nonlinear Finite Elements for Continua and Stuctures.3刘立忠.隐式静力和显式动力有限元在轧制过程模拟中的应用.塑性工程学报,2001, 8 (4) : 81834 Guyan R J. Reduction in stiffness and mass matrices. AIAA journal, 1965, 3(2):3805刘立忠,刘相华.利用显式动力学有限元法模
5、拟平板轧制过程.塑性工程学报,2001, 8(1): 51 546张国名,肖宏,用显式动力学有限元法分析板带轧制压力分布.燕山大学学报,2004, 2 (24) : 1461497庄茁等译.ABAQUS/Explicit有限元软件入门指南.北京:清华大学出版社,1998网罗洲.环件轧制过程的显式有限元模拟分析.塑性工程学报,2004, 1 (11) : 68709华林,黄兴高,朱春东,环件轧制理论和技术.北京:机械工业出版社,200110江雄心,万平荣,游步东.三维有限元模拟中的网格重划.金属成形工艺,2002, 2 (20) : 21 2411 Lee N K, Yoon J H, Yang
6、 D Y . Finite element analysis of large deformation by automatic renoding as a weak remeshing technique . Int J Mesh Sci, 1992, 34(4):255 273Abaqus Analysis Users Manual. ABAQUS Corp.,200313华林.环件轧制成形原理和技术设计方法:博士学位论文.西安:西安交通大学,2000(end)重要的作用。从而在这类问题的模拟中,如何正确适当的运用质量缩放方法,使其能够起到 在保持计算精度和稳定结果的同时提高计算效率的作用
7、,也成为了一个至关重要的问题。目 前,针对这个问题,人们也开始对其进行了相应的分析研究,东北大学在利用显式动力学有 限元法模拟平板轧制过程中,通过计算实例研究了采用提高模拟轧制速度和质量缩放技术来 缩短计算时间等方法的有效性,在忽略惯性惯性力影响的前提下得出了与实验结果相符合的 计算结果5。燕山大学在用显式动力学有限元法分析板带轧制压力分布时提出了为了加快 计算速度而采用的虚拟速度和虚拟质量带来的虚拟惯性效应会影响计算的结果6。而在环 件轧制成形模拟方面,并没有人通过具体的实例分析来验证质量缩放技术在环轧模拟中运用 的有效性和合理性。本文运用Abaqus软件建立环件轧制的三维有限元模型,并采用
8、动力 显式算法对铅环件的轧制过程模拟进行了实例计算,验证了质量缩放技术在环轧模拟中运用 的有效性,并通过对计算结果的分析比拟,给出了适合环件轧制模拟的有效的质量缩放方法。2质量缩放在环轧模拟中的意义及理论依据在Abaqus中利用Explict7来进行准静态模拟分析的求解,为了在较快的时间内获得可以 接受的静态解,以提高求解的效率,通常通过两种方法来到达目的:一种是人为的提高成型 模拟过程的速度;另一种那么是通过质量缩放来提高模拟速度。对于环件轧制,模型的材料性 能随应变率的改变而变化8,人为的提高进给速度就会人为的改变整个变形的过程,使得 与实际不相符,而且在实际的环轧生产工艺中,进给速度最大
9、值受环件咬入条件的限制9, 人为的提高进给速度,如果超过了最大值,即使模拟能够得到稳定解,在实际轧制中由于环 件不能满足咬入条件,轧制过程也是不可能正常进行的。所以人为的提高进给速度会使得模 拟与实际脱离,失去了模拟用于验证和指导实践的作用,是不可取的。所以对于环轧的模拟, 只能通过质量缩放来提高模拟速度。在Abaqus中,质量缩放是通过质量缩放因子来完成的, 质量缩放因子取的太小达不到提高计算速度的效果,取的太大那么会使得随着质量的增大而随 之产生的虚拟惯性力的增大影响到计算的结果,甚至使得结果不能收敛。所以质量缩放因子 的正确选取也就成为了环轧模拟中,质量缩放方法运用的关键。动力显式算法的
10、中心差分法是条件收敛的,其稳定的收敛条件是:其中at为时间步长,小加为临界时间步长,Tn是有限元系统的最小固有振动周期。在利用 中心差分法求解具体问题时,时间步长必须小于由该问题求解方程性质所决定的临界时间步长ter,否那么算法将是不收敛的。临界时间步长通常可以近似认为是应力波通过任一单元网格所消耗时间的最小值,即:式中Lmin为模型中单元的最小长度,与网格划分的疏密相关,cd为应力波在单元中传播的 速度,对于金属材料可以表示为:其中P为材料密度,入和口为拉梅常数,可以通过杨氏模量E和泊松比v来表示:(1 + 欧1-2廿)E 12(1 +7)代入上式那么可以表示为成min(i)+ v)(l -
11、 2廿)在Abaqus中用Explict进行求解计算时,所消耗的时间长短与增量数n有直接联系,n可 以表示为n=T/atcr ,其中T为整个模拟过程完成所需要的时间周期(对于环件轧制模拟, T即为轧制时间)。,因为入和U与材料的物理性能相关,是常数,由此可以看出n值大小与时间周期T成正 比,与单元长度L及p值成反比。n值越小,Explict在计算时处理这些增量所消耗的时间 也就越短,计算效率也就越高。所以为了减小n值,可以通过减小时间周期和增大单元长 度和密度来完成。在环轧模拟中,减小时间周期即为减小轧制时间,而轧制时间与进给速度 成反比,减小轧制时间即意味着增大进给速度,这在前面已经说过是不
12、可取的,而单元长度 的大小代表着网格的疏密,增大单元长度即为疏划网格,网格划的太疏那么会降低计算精度, 达不到进行精确分析的效果,一般来说也是不可取的。而对于增大密度,根据塑性变形的体 积不变定理,增大密度也就是增大质量,质量增大。倍,根据上面公式临界时间步长值Mcr 就会增大而倍,相应的n值就会减小而倍。所以在Abaqus中,通过质量放大来到达等效 于增大密度从而减小n值的作用,即能保证模拟的稳定性和精确度,又能有效的提高计算 效率,是最可取的一种方法。3质量缩放方法的分类及选取原那么前面已经提过,动力显式算法主要用于解决结构动力学响应分析和准静态分析两种类型的问 题,虽然都是用同一种算法,
13、但这两类问题仍然有很大差异,Abaqus/Explict中用于动力显 式算法的质量缩放方法有很多种,但这些方法对于两种类型问题的分析并不是通用的,大多 数适合于准静态分析的质量缩放方法并不适合于结构动力学响应分析。由于环件轧制模拟属 于一种准静态分析。所以这里只介绍适合于准静态分析的质量缩放方法。Abaqus/Explict所提供的用于准静态分析的质量缩放方法主要有两种类型:1 .通过用户给定的质量缩放因子进行缩放。这是最简单直接的一种质量缩放方法,用户输入 质量想要放大的倍数即缩放因子值,那么Abaqus/Explict会将临界时间步长值放大倍,而且 在计算过程的每一时间步中,都会将此时间步
14、的DT值(Abaqus/Explict中以DT表示时间 步长,以下均用DT表示)放大、回倍。2 .通过用户给定的DT值进行缩放。这种方法就是用户直接输入所需求的DT值,而不需要 通过给出质量缩放因子进行缩放,但前提条件是所给出的DT值不能小于Abaqus/Explict 根据用户所给出的材料性能参数计算出的原始DT值,否那么就起不到质量放大的作用。这里 需要说明的是Abaqus/Explict在进行求解计算时会对用户建模所生成的求解文件进行检查, 并根据文件中给出的材料性能参数p,A和p从而估算给出一个原始DT值,也就是未进行 质量缩放前的DT值,这个值其实也就是通常所说的临界步长值,但是Ab
15、aqus/Explict为了 保证计算的稳定和收敛,所给出的原始DT值一般都要比实际临界步长值要小,对于三维模 型,一般为实际值的布/31倍。然后Abaqus/Explict再根据用户给出的质量缩放方法对 DT值进行放大,如果给出的是缩放因子a,那么Abaqus/Explict将DT值增大以倍,如果用 户直接给出了 DT值,Abaqus/Explict那么直接将其替换成用户给出的DT值,两个DT值的 比值也就是用户运用这种缩放方法进行放大后,计算时间减少的倍数。所以用户在运用这种 方法时必须保证给定的DT值不能小于原始DT值。对于这种方法,根据Abaqus/Explict 对DT值进行不同控制
16、的操作,还可以将其分为三种方法。a直接根据给出的DT值进行缩放,对DT值不进行控制。这种方法其实和通过给出缩放因 子进行缩放是同一原理。因为给出缩放因子,在表现方法上也是通过给出的因子值,将DT 值进行相应倍数的放大。如果用户给出的DT值是原始DT值的4倍,效果就相当于给出的 缩放因子值为16,那么Abaqus/Explic会将每一时间步的DT值放大4倍。b用户给出DT值,在计算的每一时间步中,Abaqus/Explict自动对每一步的DT值进行控 制,使得每一步的DT值在进行缩放后不超过用户最初给定的DT值。因为在用Abaqus建 立模型时,网格划分采用的是自适应网格划分,那么在模拟过程中,
17、随着模型形状变形, Abaqus会自动对模型进行自适应网格重划由于模型各个区域的变形程度不同,所 以重划后的网格与最初网格相比,其各个单元的单元长度也发生了变化,有的单元的长度会 减小,有的那么会增大,随着模型变形的加剧,进行屡次网格重划后的模型,其各个单元的单 元长度变化也会越来越大,有的单元的长度会逐渐减小,有的那么会逐渐增大,Abaqus/Explict 在确定每一时间步的DT值时,为了保证所有单元的计算收敛,是以这一时刻的应力波穿过 长度最大的一个单元所消耗的时间为这一时间步的DT值进行计算的。这样,如果这个DT 值不超过临界时间步长那么就能保证这个时间步计算的收敛性。所以在具体确定某
18、一时间步的 DT值时,随着有些单元的长度的逐渐增大,应力波穿过最大一个单元所消耗的时间也就越 大,那么这一时间步的DT值就会增大,可能会超过原始DT值,甚至有可能超过了临界时 间步长值,那么这一时间步的计算就不收敛了。为了保证计算的精确度,Abaqus/Explict 中的自适应网格划分在每一时间步都会对模型网格进行重划,也就是说整个计算过程的每一 时间步的DT值都是在变化的。所以这种缩放方法通过对每一时间步的DT值进行控制,使 其不超过原始DT值,以防止出现DT值超出临界时间步长值,使得计算结果不收敛的现象, 对于大变形的准静态分析是非常合适的。c用户给出DT值,在计算的每一时间步中,Aba
19、qus/Explict自动对每一步的DT值进行捽 制,使得每一步的DT值在进行缩放后都等于用户最初给定的DT值。这种方法通过缩放 DT值使每一时间步的DT值都恒等于最初的DT值,这样也就防止出现上面所说的随着网 格重划而出现的逐渐增大的DT值超过临界时间步长值使得计算不收敛的现象。所以这种方 法对于保证准静态分析过程的收敛性也是很有效的。以上几种缩放方法都适合于准静态分析。但是,无论对于通过质量缩放因子缩放,还是给出 DT值缩放,如果缩放系数太大,那么随之产生的虚拟惯性力会影响到计算的精度和收敛性, 缩放系数太小那么达不到提高计算效率的目的。在准静态分析中12,为了保证模拟过程的稳 定性,一般
20、要求变形的动能值较小,且变化平稳,势能值那么较大,动能与势能的比值较小, 一般不能超过0.1,这样使得能量大局部用于模型的变形。所以一般在准静态分析中,质量 缩放方法的基本选取原那么就是:无论采用何种方法进行缩放,必须保证在整个过程中的动能 值要较小,动能与势能的比值不能超过0.1。在满足这个前提条件的情况下进行缩放,一般 都是可行的。4质量缩放方法的实例比拟本文以Abaqus软件为操作平台,以工业纯铅为材料,建立铅环件轧制的三维有限元模型。 并用以动力显式算法为基础的求解器Abaqus/Explict进行求解计算。有限元模型如下:驱动辐驱动辐图1环件轧制三维有限元模型图导向银芯辐环件工业纯铅
21、密度p=11340kg/m3,弹性摸量E=17GPa,泊松比v=0.3,最初网格划分的单元长 度为0.01m,根据上文所给出Zktcr公式,将上述值代入求得临界时间步长理论值为tcr=4.955x10%,而Abaqus/Explict在检查求解文件时所给出的原始DT值为2.511 x10 6s由此可以看出为了保证计算过程的稳定性,Abaqus/Explict所给出原始DT值要比理论值小。首先采用给出缩放因子进行缩放的方法,来比拟不同数量级的缩放因子,对模拟结果的影响, 本文分别取缩放因子为4、100、1000对环件轧制过程进行了模拟,环件最终的成形状况如 图2所示:图2不同数量级的缩放因子时的
22、环件最终的成形状况图2中从左至右依次分别为质量缩放因子为4、100、1000时环件最终的成形网格图。通过 比照可以看出,最左边的图效果最好,为圆环形;中间的图基本为圆环状,但左中局部呈明 显的椭圆形,效果比要差;最右边的图效果最差,环件基本没有变形左中局部就被压扁了。 由此,单从形状上看,采用1000的缩放因子是不合理的。缩放因子分别为4和100时,环 件端面宽展变形,环件变形的动能历史,及动、势能比值比拟如下:6555453325150(山陵花也纹0 5 10 15 20 25 30环件在向IB J)1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.11.41.31.2
23、01234时间/S图3环件端面的宽展变形图4环件变形的动能历史图5动、势能比值图3横轴是环件的径向尺寸,从环件的内环面开始为0,纵轴是环件端面的宽展距离。从图 中可以看出,缩放因子为4时,环件沿径向的宽展变形较平缓,基本符合实际轧制中两端 宽展大,中间宽展小的情况,成形后的端面凹陷不会很明显。而缩放因子为100时,环件 沿径向由内向外呈明显的上升趋势,外端面的宽展要比内端面大很多,那么成形后就会形成 很明显的端面凹陷13。再比拟环件变形的动能历史,从图4可以看出,缩放因子为4时, 环件在整个变形过程的动能值较小,且变化平缓,符合上文所述的准静态分析中,要求有较 小且变化平稳的动能值的原那么,这
24、样环件变形过程就会较稳定。缩放因了为100时,在同 一时刻环件变形的动能值儿乎都要比前者大很多,而且整个变形过程中动能变化很不平稳, 升降幅度很大,不符和上述原那么,这样环件变形过程就会不稳定,产生振动或摇摆的现象。 最后比拟环件的动、势能比值,从图四中可以看出,缩放因子为4时,在变形的初始阶段, 环件的动、势能比值超过了 0.1,最大值为0.3,且不稳定,有波动现象,到大约0.4秒时比值下降为小于0.1,随后趋于平稳。而缩放因子为100时,动、势能比值不仅明显超过了 0.1,而且比值较大,最大值到达了 1.3,比同一时刻的缩放因子为4时的比值耍大很多。 在接近1秒时比值才小于0.1趋于稳定,
25、不稳定时间持续较长,几乎为整个过程的1/3,而 在稳定阶段,比值也比缩放因子为4时高。由此,从以上3幅图的比拟可以看出:无论从 工艺需求的宽展变形和端面凹陷方面,还是从Abaqus/Explict本身进行准静态分析时对动 能输出及动、势能比值的要求方面来看,缩放因子取4比取100要更合理。以上通过比照分析确定了缩放因子取4时符合本文所建铅环件模型的模拟需求。下面比拟 直接给出缩放因子进行缩放即上文所述方法1和通过给出DT值进行缩放,并对DT值进行 不同控制的缩放方法即上文所述的方法2中的方法b和方法c(前面已经讲过方法2中的方 法a和方法一其实是一样的,所以在下面的比拟中以方法a代替方法1)这
26、三种方法对模拟 成形的影响。模拟过程中仍以缩放因子4为前提,转化为DT值即为放大2倍。采用3种不 同方法时,环件端面宽展变形,环件变形的动能历史,及动、势能比值比拟如下:5 3525 15 O 一 3 2 1030252015105 03052205 O 5 O2S5 MLWI533图6环件端面宽展变形图7环件变形的动能历史从图6可以看出,三种方法中环件端面的宽展变形并无明显差异,由此可以看出这三中缩 放方法对宽展变形的影响不大。从图7可以看出,虽然三种方法的动能变化都较平缓,但 同一时刻方法b的动能值比方法a和c明显要小,方法a和c的动能值变化那么基本一致。 图8为三种方法的环件变形动、势能
27、比值比拟,从图中可以看出,方法a和c刚开始阶段 比值都超过了 0.1,但c相对a要好,最大比值为0.17,而“最大值为0.31, c的比值波动 的幅度要比a小,比值降为小于0.1时所经历的时间也要比a短。方法b的效果那么最好, 环件变形的整个过程中动、势能比值都要小于0.1,最大值为0.095,而在随后的稳定阶段, 同一时刻方法b的动、势能比值也要比方法a和c的小。图8为计算过程中每个时间步的 DT值输出。由于Abaqus/Explict给出原始DT值为2.511x10%,而最初都是将DT值放大 2倍即为5.022x10-6s所以从图中可以看出三种方法刚开始的DT值都为5.022x10%。方
28、法a不对DT值进行任何控制,从图中可以看出有些时间步的DT值明显大于最初值,虽然 在本模型的分析中,这些较大的DT值并没有超过实际的临界时间步长而导致计算结果不收 敛,但在其它模型的分析中,如果Abaqus/Explict所给出的原始DT值与实际临界时间步长 值很接近时,那么运用方法a进行缩放时,较大的DT值就可能会超过实际临界时间步长值而使得计算结果不收敛。方法b是控制DT值,使每步的DT值都不超过最初值,从图中也 可以看出,整个过程的DT值处在一个下降的趋势中,旦都不超过5.022x10上。方法c是 控制DT值,使每步的DT值都等于最初值,从图中明显可以看出,整个DT值变化为一条 直线都等
29、于5.022x10%。另外,对于使用三种缩放方法进行模拟时,记录CPU计算消耗 口寸间分别为1小时39分,1小时40分和1小时42分。可见三种缩放方法计算效率儿乎一 样。 be 法法法方方方654321001234时向/S图8环件变形的动、势能比值图9环件变形的DT值比拟从以上四幅图的分析比拟中可以看出方法b既符合准静态分析对动能输出及动、势能比值 的要求,又能保证有效的计算效率,所以是三种方法中最适合本模型分析的质量缩放方法。总结以上的分析比拟可以得出:准静态分析中质量缩放方法的选取原那么对环件轧制的模拟是 适用的。对于环轧模拟中缩放方法的运用,首先应该从环件成形的形状、宽展变形、动能变 化
30、及动、势能输出等方面进行比拟,确定一个较合适的缩放因子,然后再运用方法b,通过 对每一时间步的时间步长值进行控制,保证其不超过临界时间步长值。这样,即提高了计算 效率又保证了计算过程的稳定性和收敛性。5结论本文采用不同的质量缩放方法及放大系数对铅环件三维有限元模型进行了实例分析,通过分 析比拟不同缩放系数时的模拟结果,验证了在准静态分析中,为了保证模拟过程的稳定性, 质量缩放方法的选取必须遵循:通过缩放后的整个模拟过程其变形的动能值应较小,且变化 平稳,势能值那么较大,动能与势能的比值应较小,一般不能超过0.1,这一选取原那么对环件 轧制模拟是适用的;通过分析比拟采用不同缩放方法时的模拟结果,总结了一套适用于环件 轧制模拟的最有效的质量缩放方法。参考文献1曾攀.有限元分析及应用.北京:清华大学出版社,2004