刚塑性有限元解题法的改进及在轧制中的应用

D0I:10.13374/j.issn1001-053x.1987.04.009 北京钢铁学院学报 第9卷第4期 Journal of Beijing University Vol,9 No.4 1987年10月 of Iron and Steel Technology 0ct.1987 刚塑性有限元解题法的改进 及在轧制中的应用 硕卓 贺毓辛 (压力加工系) 摘 要 本文对刚塑性有限元的初速度场及收敛性进行了专门的研究和改进。用于解决 轧制工程问题，计算精度较高、CPU时间铰少，它是一种可靠的理论分析方法. 关键词：刚塑性有限元法，轧制，初速度场，收敛性 Improvement on the Rigid-Plastic Finite Element Method and Its Application in Rolling Engineering Gu Zhuo He Yuxin Abstract Rigid-plastic Finite Element Method RFEM)is non-lincar,so it is necessary to give a good initial velocity field and to ensure for evaluation. Ori K.I.proposed "function G"to determine the initial velocity field, but he didn't make a through analysis,Newton-Raphson method is gene- rally used in RFEM,but it has some drawbacks: 1 it is possible to deverge during itterating;and 2 )the step is difficult to determine,so that the calculating time is increased. In order to overcome the drawbacks above mentioned,we advanced a new 1986-06一02收 57

冲 第 卷第 期 年 月 北 京 钢 铁 学 院 学 报 。 。 卜 刚塑性有限元解题法的改进 及在轧制中的应用 硕卓 贺毓辛 压力加工系 户 摘 要 本文对刚塑性有限元的初速度场及收敛性进行了专门的研究和改进 用于解决 轧制工程问题 ， 计算精度较高 ， 时间较少 ， 它是一种可靠的理论分析方法 关镇词 刚塑性有限元法 ， 轧制 ， 初速度场 ， 收敛性 一 目 “ “ 万 ‘劣 ” 一 一 ， 。 。 “ ， ， 一 ， 多 ， ， ， 一 一 一玫稿 DOI ：10．13374／j ．issn1001－053x．1987．04．009

method-penality method of function G and carried on a detail analysis, As iteration method the Gill-Murrary mothod (improved Newton method) is used instead of Newton-Raphson method and the linear search is made by 0.618 method, The RFEM program is worked out for calculating parameters of pla- stic working processes.The calculating results obtained in studying roll- ing process coincide with the:expefimentar data.The'research shows cle- arly that the improved RFEM method is satisfactory for analysing plastic working problems. Key words:rigid-plastic finite element method,rolling,initial velo- city field,convergence. 前 言 近年来，随着电子计算机的普及，许多计算方法已引入到塑性加工中来，1973年小 林史郎(1)推出了刚塑性有限元的矩阵表达式，这个方法已在塑性加工中有了许多运用， 如计算轧制压力、力矩〔2,3)，以简化的三维单元模拟钢锭的动态轧制4)等。 刚塑性有限元法以虚功原理为基础，能量做泛函，速度“为未知量，即： 中()=Svoedy-∫srTi4:ds (1) 式中σ、e分别为等效应力和等效应变速率，V为体积，T为已知边界的应力，5r与“分 别为已知应力边界的表面与速度。 第(1)式中各项的物理意义为：第一项表示变形能，第二项是外力所做的功。以 其为基本式，引入体积不变约束条件求最小值而得到问题的真实解。因引入体积不变条 件的方法不同有三种不同的刚塑性有限元法：Lagrange乘子法、罚函数法、可压缩法。 这几种方法在使用中都存在着一些问题和困难：泛函④与引入的体积不变约束都是非线 性的，进行优化求最小值时必须迭代求解。减少迭代次数和时间需要给出一个较好的初 速度场，为了保证计算的可靠性，必须采用有严格数学证明的方法。虽然上述两个问题 已有一些解决方法，但在使用中，由于塑性加工问题大多是3维问题，用刚塑性有限元 法求解单元多、节点多。这些方法求解难以收敛，机时长且不可靠，特别是轧制问题， 边界较复杂，在实际中很少用于？维问题。为了将有限元法广泛应用于塑性加工，必须 解决这些问题。 本文专门对刚塑性有限元的初速度场及收敛性进行了研究和改进，用改进后的方法 编制了程序对高件轧制的3维变形进行了模拟，对平辊薄件轧制（平面问题）进行了计 算。 1初速度场的设定 森谦一郎(5)等对可压缩材料的刚塑性有限元法提出了以G函数求初速度场的方法， 58

一 。 一 一 以 飞旋 一 圣 ， ‘ 誉 ‘ 奋 几 补 ’ ‘ 电 ， ‘ 一 ， ， ， 月 舌 近年来 ， 随着 电子 计算机的 普及 ， 许多计算方法 己 引入 到塑 性加工 中来 ， 年小 林史郎巾 推出 了刚 塑 性有限元 的矩 阵表达 式 ， 这 个方法 已在 塑 性加工 中有 了许多运用 ， 如计算轧制压 力 、 力矩〔 “ ， 〕 ， 以 简化的三维单元模拟钢锭 的动态轧制 印等 。 刚塑性有限元法以 虚功 原理 为基础 ， 能 量做泛 函 ， 速 度 为未知 量 ， 即 巾 厂 一 ， 式 中石 、 秘别为等效应 力 和等效应 变速率 ， 为体积 ， 云为 已知 边 界的应力 ， 与。 分 别为 已知应力边界的表 面与速 度 。 第 式中各项的物理意义 为 第一项表示 变形 能 ， 第二项是 外力所 做的功 。 以 其为基本式 ， 引人体积 不变约束条件求 最小值而得到 问题 的真实解 。 因引入体积 不 变条 件的方法不 同有三种 不 同的刚塑 性 有限元法 乘子法 、 罚 函数法 、 可压缩法 。 这几种方法在使用 中都存在着 一些 问题 和 困难 泛 函中与引人的体积 不 变约 束都是非线 性的 ， 进 行优 化求 最小值时必须迭 代求 解 。 减 少迭代次 数和时 间需 要给 出一 个较好 的初 速度场， 为 了保证 计算的可靠性 ， 必须采用 有严格数学证 明的方法 。 虽然上述 两个问题 已 有一些解决方法 ， 但在使用 中 ， 由于塑性加工 问题大 多是 维问题 ， 用刚塑性有限元 法 求解单元 多 、 节 点多 。 这些方法 求解难以收敛 ， 机时长且不可靠 ， 特别是轧制问题 ， 边界较复杂 ， 在实际 中很少用于 维问题 。 为 了将有限元法 广泛应用 于塑 性加工 ， 必须 解决这 些问题 。 本文专门对刚 塑性 有限元的 初速 度场及 收 敛性进 行 了研究 和改进 ， 用改进后的方法 编制 了程序对高件轧制 的 维变形进 行 了模拟 ， 对平辊 薄件轧制 平面 问题 进 行 了计 算 。 初速度场的设定 森谦一 郎〔 〕 等对可压 缩材料 的刚塑 性有限 元 法提 出 了以 函数求初速 度场 的方法

所设G函数的形式为 M M2 M3 G= ∑(oeV)2+∑（x·Vu·s)2±∑(T::)2 (2) n11 n=1 m=1 式中：x为摩擦边界的摩擦力，△v为工件与工具的相对速度。v:为表面节点速度。M1 为总单元数，M2为有摩擦力的单元数，M3为有外力的表面单元数。正负号视T::值而 定，若T::>0为负，反之为正。 认为G函数中各项是泛函： Φ=vvedV+ .r…adg-Jw. (8) 中各项的平方和，使G为最小的速度场接近于使泛函中最小的速度场。设立的G函数是 一个二次函数，用Newton法求解，只须一次计算即可求得G的最小值。白光润r6)等用 有限元计算三维问题时也应用了G函数。但是，从森等所设的G函数可看出，在每个单 元中，将o、ε视为常量时才能运用。这与大多数实际问题不相符。另外，对G函数在实 际中运用时，什么条件下能得到接近真实解的初速度场，设有进行深入的讨论。 本文采用的是罚函数法，设泛函为 o-∫oay+∫s.r…ads-∫s,Tids+%(∫vdr）° 式中M是一个大的正数，称为罚函数因子，εv为体积应变速率。 在森谦一郎G函数的基础上，我们提出了罚函数的G函数。设： o=g小a)'+三(…a)广±三不ao） (6) 式中的σ、ε等均为速度的函数，这样G函数就带有普遍性。为了使(5)式对任意节点 的等参单元均能使用，本文以矩阵推导使用公式。 设}=〔B)《a,〔K)=〔B)〔B,{R}=∫T(N)d, {Q}=∫va〔B)rds(C}d,则e=(号{u}〔K){u})古 式中：〔B〕为应变矩阵；{“}为节点速度列阵，{C}为常数矩阵，在三维问题时为 (1,1,1,0,0,0)T,〔N)为形函数。 为求解轧制问题，设座标如图1所示。山图可知： {△v}=(wx-vcos0,vy,vz+vin0)T。 式中：vx、,、vz为轧件与轧辊接触面上某点的三个速度分量，v为轧辊线速度，为该 点的夹角。 59

所设 函数的形式为 斌 一 － 一 － 一 艺 厂 “ 艺 · · 士 艺 一 一 卜 式 中 为摩擦边 界的摩擦力 ， △ 为工 件与工具 的相 对速度 。 。 ，为表面 节 点速 度 。 为总单元 数 为有摩擦力的单元数， 为有外力的表面单元 数 。 正负 号视 ，值而 定 ， 若 为负 ， 反 之为正 。 认 为 函数中各项是 泛 函 ， 声 犷 丘二 △·‘一 。 了一 ‘ · 中各项的平方和 ， 使 为最小 的速度场接近 于使泛 函中 最小的速度场 。 设立 的 函数 是 一个二次 函数 ， 用 法求 解 ， 只须一次计算即可求得 的最小值 。 白光润卿 等 用 有限元计算三维问题 时也应用 了 函数 。 但是 ， 从森等所设的 函数可看 出 ， 在每 个 单 元 中 ， 将 、 。 视为常量 时才能运 用 。 这与大 多数实际问题 不相符 。 另外 ， 对 函 数在实 际 中运用 时 ， 什 么 条件下能得到接近 真实解的初速度场 ， 没 有进行深 入的讨论 。 本文采 用 的是 罚函数法 ， 设泛 函 为 ， 丁 ‘ 昆 丁二 △· ‘一 。 厂一 ‘ · 一养 一 二 · ‘ 厂 式 中 是 一 个大的正 数 ， 称 为 罚函数因子 。 为体积应变速率 。 在森谦 一 郎 函数的基础上 ， 我们提 出了罚 函数的 函数 。 、， 、 ‘ 、少 崖丁 几 ￡ 厂 、 “ 瞥 ， 少 十 ， 雇 ‘ 二 · 。 艺 ” 留二二 “ … ‘别卜七、 一 、 艺 诊留二二 汽 、 二 · 厂 式 中的 、 石等均 为速度的 函数 ， 这 样 函数就带有普遍性 。 为 了使 式对任意节点 的 等参单元 均能 使用 ， 本文以 矩阵推导 使用公式 。 设 八卜 〔 。 〕 、 · ， ， 〔 、 卜 〔 〕 〔 。 〕 ， 、 卜 丁 。 式 〔 、 〕 ‘ 一 ， 丁 〔 〕 · · 一 贝。奋已 ， 、 ， 、 气 百 戈“ ， ‘ 八 ’‘ 了 ， 式 中 〔 〕 为应 变矩 阵， 。 为 节点速度 列阵 为常 数矩阵 ， 在三 维问题时为 ， ， ， ， ， 〔 〕 为形 函数 。 为求 解轧制 问题 ， 设 座标如 图 所示 。 山图可知 一 ， 。 ， ， 。 式 中 、 。 ， 、 为轧件与轧辊接触 面 上某 点的 三 个速 度分 量， 。 为轧辊线速度， 为该 点的夹 角

图1计算的座标示意图 Fig.1 Schematie drawing of coordinate used in calculating 又设：{vo}=(-ucos9,0,si0)T,从有限元的基本公式可知：{△v}= 〔N){“}-{vo},这样，（5)式可表示为： G-((号(K)业) +三(s《a如ao)) M3 ±2{u}T{R}{R}T{w} 1T世压 +之，g{u)r(Q}{Q)r{u} + (6) 而(6)式的平方项在积分号内，不是二次函数，必须进一步简化。先以高斯积分法公 式引人(6)式再简化得： G=空(o空空支{u}tK){u}1J1HHiH) m=1 i=1j=1R=1 +竖2空名(N){u}-{})r n=1i=1j-1 M3 (〔N){u}-{oo})lJ2I2HH;±Σ{u}r{R}{R}T{u} 有12g】 M{ur1QiQ'w} +Σ m=1 2Vm (7) 式中：|J21,IJaI分别表示2阶与3阶雅可比行列式的值，H,H,Hx为高斯 积分的加权系数，1为高斯点数。 这样，G为{}的一个2次函数，在给定边界条件下一次即可求得最小值，得到初 速度场。 60

川月 图 计算的座标示意图 。 匕 一 又 设 。 。 一 。 ， ， ， 从 有限元 的基本公 式 可知 〔 〕 一 二 ， 。 。 ， 这 样 ， 式可 表示 为 丁沪 令 ‘ · ’ · 崖 丁 · 〔 、 〕 厦 一艺 、 一牙 犷 、 △。 △。 毛 艺 一 粤二 二一犷芬一 艺 厂 夏 夏 而 式的平方项在积分号 内 ， 不是 二次 函数 ， 必须进一步简化 。 先以 高斯积分 法公 式引人 式再简化得 艺 艺 艺 艺 蚤 〔 〕 。 食 贾 畏 万 艺 〔 〕 一 。 · 〔 〕 “ 一 。 。 贾 爹士 艺 “ 卜 “ 里、 艺 几 脚 丁 艺 式 中 ， 分别表示 阶与 阶雅可 比 行列式的值， 。 ， ， 为高 斯 积分 的加权系数 ， 为 高斯点数 。 这样 ， 为 。 的 一 不 次 函数 ， 在给定边界条件下 一次 即 可求得最小 值 ， 得 到 初 速度场

如果就这样给出初速度场，并不一定能够得到理想的接近真实解的结果。本文对G 函数进行了分析：G函数在公式中引人了体积不变条件作为约束，计算中令速度场满足 应变速率和速度的关系及边界条件，基本上得到的是一个运动许可的速度场：采用了有 限元的基本方法分割单元，变形体内各点连续。用G函数的目的是求初速度场，得到的 结果与给出的边界条件有关。如果给出速度边界条件，得到的速度场一般较好，如果给 出的应力边界或混合边界条件，虽然初速度场满足边界应力，但不满足应力与应变速率 的关系，这时的初速度场不一定满足边界的真实速度。又由于连续性，内部的速度场也 不一定接近真实速度场。实际的计算中证实了这个问题。 因此，在使用G函数时，尽可能给出速度边界条件。如果遇到应力边界条件和混合 边界条件，可根据实际问题，假设一个近似真实速度场的速度边界条件，代人求初速度 场，在随后的迭代中再精确满足真实应力边界或混合边界。 2收敛性的讨论 刚塑性有限元法求最小值的过程属于最优化技术中的非线性规划问题一等式约束 问题，即求解： minΦ(a) WEE 满足约束条件：ev1(a)=0 j=1,2,,M1。将约束条件引入泛函后， 变为无约束优化问题。在求最小值的过程中，常用Newton法进行迭代。这种方法的思 路是将中(4)用一个二次函数甲（山）来逼近。如果4：为中(4)极小点的一个近似点， 将中(4)在4K处做Taylor展开，且略去高于二次的项，得： 西()=p()=中(4k)+(#-)△(.) +女（4-)△2中(）(-) (8) 式中：△中(ux)为中(u)在4x点的梯度函数，V20(4:)为中(4)在u点的二阶导 数矩阵一Hessian矩阵，可以得知V2中(x)是对称矩阵。 取p(4)的极小点作为中(4)极小点的一个近似点4K+1,对(8)式必有： 7pk（4k)=0 (9) 从(8)式得： Vpk(.)=V中(k)+V2中(.）(a-) 代入（日）式： -=-〔v2中(4）〕-1v西() (10) 设： P.=-〔V2Φ(k)〕17φ(4.) 1) 61

卜 如果就这样给出初速度场 ， 并不一 定能够得到理想的接近真实解的结果 。 本文对 函数进行 了分析 函 数在公式 中引人 了体积不 变条件作为约束 ， 计算中令速度场满足 应 变速率 和速度的关 系及边 界条件 ， 基本上得到的是 一 个运动许可的速度场， 采用 了有 限元的基本方法分割单元 ， 变形 体内各点连续 。 用 函 数的 目的是 求初速度场 ， 得到 的 结果与给 出的边界条件有关 。 如果给 出速度边 界条件 ， 得到的速度场一般较好 如果给 出的应力边界或混合边界条件 ， 虽然初速度场满 足 边 界应力 ， 但不满 足 应力与应变速率 的关 系 ， 这时的初速 度场不 一 定满 足边 界的真实速 度 。 又 由于连续性 ， 内部的速度场也 不一 定接近 真实速度场 。 实际 的计算 中证实 了这个问题 。 因此 ， 在使用 函数时 ， 尽 可能给 出速度边 界条件 。 如果遇 到应 力边界条件和混 合 边界条件 ， 可根 据实际问题 ， 假设 一 个近 似真实速度场 的速度边界条件 ， 代 人求初速度 场 ， 在 随后 的迭 代 中再精确满 足 真实应力边界或混合边 界 。 收敛性的讨论 刚 塑性有限元法求 最小 值的过 程属于 最优 化技术 中的非线性规划 问题－ 等式约 束 问题 ， 即求 解 中 “ 〔 满 足 约 束条件 “ ， ， ” 一 ， 。 将约 束条件引 人泛 函 后 ， 变为无约 束优 化问题 。 在求最小 值的过 程 中 ， 常用 法进 行迭代 。 这种 方法 的 思 路是 将中 “ 用一 个二次 函数甲 来逼近 。 如果“ 为巾 极小 点的一个近 似 点 ， 将中 在“ 处 做 展开 ， 且略去 高于二次 的项 ， 得 中 自 甲 巾 “ 、 一 ， △少 、 于 一 。 △“ 中 、 一 式中 △中 “ 为巾 “ 在“ 点的梯度 函数， “ 少 七 为必 。 在 “ ， 点的 二 阶 导 数矩阵－ 矩阵 ， 可以得知 “ 巾 是 对称矩阵 。 取甲、 “ 的极小 点作为巾 “ 极小 点的一 个近似 点。 ‘ ， ， 对 式必 有 、 、声矛、 从 式得 甲 、 代 人 式 甲 。 “ ， 中 ， “ 中 一 ， “ 一 “ 七 一 〔 “ 巾 “ 〕 一 中 “ 。 设 尸 〔 少 ‘ 〕 一 巾 “ ‘