（基础数学专业论文）改进循环中点求积newton法在板材轧制中的应用.pdf

- 7 一 at h e s i si nf u n d a m e n t a lm a t h e m a t i c s i i l lll iri l l 1l llj illl y 18 4 12 9 4 l 一- l-1-。 a p o l i c a t i o no li m p r o v e dr o t a t l v em i d - p o i n t m e n s u r a t i o na n dn e w t o nm e t h o d i nf l a tr o l l i n g b yw a n gh a i p i n g s u p e r v i s o r ：p r o f e s s o rs o n gs h u n i n or t h e a s t e r nu n i v e r s i t y o c t o b e r2 0 0 7 ， 一| - 卜l i， 、e 0+iiii 卜，。1 ， 。，_一 。 ， 独创性声明 本人声明，所呈交的学位论文是在导师的指导下完成的。论文中 取得的研究成果除加以标注和致谢的地方外，不包含其他人已经发表 或撰写过的研究成果，也不包括本人为获得其他学位而使用过的材料。 与我一同工作的同志对本研究所做的任何贡献均已在论文中作了明确 的说明并表示谢意。 学位论文作者签名：兰沤孕 1 日期：聊年，i ，月2 自 学位论文版权使用授权书 本学位论文作者和指导教师完全了解东北大学有关保留、使用学 位论文的规定：即学校有权保留并向国家有关部门或机构送交论文的 复印件和磁盘，允许论文被查阅和借阅。本人同意东北大学可以将学 位论文的全部或部分内容编入有关数据库进行检索、交流。 ( 如作者和导师不同意网上交流，请在下方签名；否则视为同意。) 学位论文作者签名： 签字日期： 导师签名： 签字日期： 东北大学硕士学位论文摘要 改进循环中点求积n e w t o n 法在板材轧制中的应用 摘要 有限元法( f e m ) 是一种广泛应用于板带轧制过程的有效数值计算方法目前在 保证精度的前提下，提高计算效率是有限元法在线应用的核心问题大量的数值 方法( 如n e w t o n 法) 被应用于求解有限元中形成的非线性方程组，形成f e m 高速在 线计算的理论本文针对求解有限元法中形成的非线性方程组，提高f e m 的计算 效率展开研究，主要内容如下： ( 1 ) 提出了改进循环中点求积法，它可为n e w t o n r a p h s o n ( n r ) 法寻找一个较 好的初始值改进循环中点求积法和n r 法相结合，称作改进循环中点求积n e w t o n 法( 简称为m p 法) 可用来求解非线性方程组 ( 2 ) 给出了m p 法求解非线性方程组的具体计算步骤和程序框图，并对改进 循环中点求积法的稳定性、收敛性和局部截断误差进行了分析，证明了它是稳定 的、收敛的二阶方法 ( 3 ) 为了检验改进循环中点求积n e w t o n 法的可行性，将其与实际应用相结合 求解刚塑性有限元中形成的速度增量方程，并选取某钢厂的现场数据进行模 拟试验 ( 4 ) 给出步长参数对算法的影响，并在计算时间、迭代次数、轧制力等方面将 m p 法与n r 法进行了比较n r 法的计算时间比m p 法的计算时间多且m p 法 计算的轧制力与现场实测值吻合较好数值实验验证了m p 法在板材轧制中的有 效性和稳定性 关键词：中点求积法；n e w t o n 。r a p h s o n 法：刚塑性有限元法；轧制力 6 l i i 东北大学硕士学位论文 a b s t r a c t a p p l i c a t i o no f 1 r o t a t i v em i dp o i n t i m o r o v e or o t a t i v em i d - o o i n t m e n s u r a t i o na n dn e w t o nm e t h o d i nf l a tr o l l i n g a b s t r a c t f i n i t ee l e m e n tm e t h o d ( f e m ) h a sb e e nw i d e l ya p p l i e dt os t r i pr o l l i n gp r o c e s s a t p r e s e n t ，h e i g h t e n i n g t h e c a l c u l a t i n g e f f i c i e n co ff e mw i t he n s u r i n gc a l c u l a t i n g p r e c i s i o nh a sb e e nt h ek e yi s s u ei ns o l v i n gp r o b l e m so fo n l i n ea p p l i c a t i o n m a n y n u m e r i c a lm e t h o d ss u c ha sn e w t o nm e t h o dh a v eb e e ne m p l o y e dt os o l v et h es y s t e mo f s i m u l t a n e o u se q u a t i o n so b t a i n e di nf ea n df o r mt h e o r yr e s e a r c ho ff e mi ns t r i po n l i n e r o l l i n g i nt h i sp a p e r ，t h ei m p r o v e dr o t a t i v em i d p o i n tm e n s u r a t i o ni sp r e s e n t e d t os o l v e t h es y s t e m so fs i m u l t a n e o u se q u a t i o n so b t a i n e di nf ea n dt h es t u d yi sp r o g r e s s e dt o h e i g h t e nt h ec a l c u l a t i n ge f f i c i e n co ff e m t h ef o l l o w i n g a r eo b t a i n e d ： ( 1 ) t h ei m p r o v e dr o t a t i v em i d - p o i n tm e n s u r a t i o ni sp r e s e n t e da n di tc a np r o v i d ea g o o dv a l u ef o rn e w t o n r a p h s o n ( n r ) m e t h o d t h ei m p r o v e dr o t a t i v em i d - p o i n t m e n s u r a t i o nc o m b i n i n gw i t ht h en e w t o n - r a p h s o nm e t h o dd e s i g n a t e dt h em - pm e t h o d ， w h i c hi se m p l o y e dt os o l v et h en o n l i n e a re q u a t i o n s ( 2 ) t h ei d i o g r a p h i cc o m p u t i n gp r o c e s sa n dp r o c e d u r et h a tt h em - p m e t h o ds o l v e s n o n l i n e a re q u a t i o n sa r eg i v e d t h ea s t r i n g e n c y , s t a b i l i t ya n de r r o ro ft h ei m p r o v e d r o t a t i v e m i d p o i n t m e n s u r a t i o na r ea n a l y z e d t h e i m p r o v e d r o t a t i v e m i d 。p o i n t m e n s u r a t i o ni ss t e a d ya n dc o n v e r g e n t ( 3 ) t h em pm e t h o di s u s e dt os o l v et h ee q u a t i o n so fv e l o c i t yi n c r e m e n ti n r i g i d p l a s t i cf e m f o rp r o v i n gi t s f e a s i b i l i t y w eu s et h em - pm e t h o dt oc a l c u l a t es p o t d a t ao fs o m es t e e l w o r k s ( 4 ) i nt h i sp a p e r , i th a sb e e na n a l y z e dt h a tt h ep e d o m e t e rp a r a m e t e ri n f l u e n c e st h e m pm e t h o d t h em pm e t h o da n dt h en - r m e t h o da r ec o m p a r e di nt e r m so fc p ut i m e ， t h en u m b e ro fi t e r a t i o n sa n dr o l l i n gf o r c ea n ds oo n t h ec p ut i m er e q u i r e df o r c a l c u l a t i o no ft h en rm e t h o di sm o r et h a nt h em pm e t h o d t h ec o m p u t e dr o l l i n g i i i o jl ir、 东北大学硕士学位论文 f o r c e sb yt h em pm e t h o dc o r r e l a t e sb e t t e rw i t ht h em e a s u r e dv a l u e n u m e r i c a lt e s t s s h o wt h a tt h em pm e t h o di sf e a s i b l ea n ds t e a d y k e y w o r d s ：m i d - - p o i n tm e n s u r a t i o n ；n e w t o n - r a p h s o nm e t h o d ；r i g i d - p l a s t i cf e m ； r o l l i n gf o r c e 参 l 东北大学硕士学位论文 目录 目录 独创性声明i 摘要i i a b s t r a c t i i i 目录v 第1 章绪论1 1 1 轧制理论中刚塑性有限元法的概述1 1 1 1l a g r a n g e 乘数法2 1 1 2 罚函数法2 1 1 3 可压缩法2 1 2 有限元软件在轧制中的应用3 1 3 有限元法求解存在的问题3 1 4 循环中点求积法的发展4 1 5 本文的主要研究内容4 第2 章刚塑性有限元理论5 2 1 刚塑性有限元的理论问题5 2 2 刚塑性有限元的求解途径5 2 3 有限元公式6 2 4 总能耗率泛函的最小化9 2 5 收敛准则1 0 第3 章改进循环中点求积n e w t o n 法1 1 3 1 准备知识1 1 3 2 改进循环中点求积法1 6 3 - 3 性态分析1 9 3 3 1 局部截断误差1 9 3 3 2 收敛性2 0 3 3 3 稳定性2 3 第4 章数值分析2 7 v 目录东北大学硕士学位论文 4 1 网格划分。2 7 4 2 速度边界条件2 7 4 3 变形抗力模型2 8 4 4 步长参数对m p 法的影响2 9 4 4 1 步长参数对m p 法的迭代次数和计算时间的影响3 0 4 4 2 步长参数对m p 法的轧制力的影响3l 4 5 划分不同的网格数对m p 法的影响3 2 4 6 改进循环中点求积n e w t o n 法的有效性3 5 第5 章结论3 7 参考文献3 9 致谢4 3 v i 东北大学硕士学位论文第1 章绪论 第1 章绪论 随着现代科学技术和经济的发展，产品精度要求越来越高，控制参数越来越 多，各种高新轧制技术不断出现在计算机技术和数值分析方法支持下发展起来 的有限元法( f e m ，f i n i t ee l e m e n tm e t h o d ) 为解决这些复杂的工程分析计算问题提 供了有效的途径随着大型有限元模拟软件的日益成熟和完善，有限元法已成为 一种被广泛应用于板带轧制过程的有效的数值计算方法对于板材轧制问题，提 高求解速度是实现f e m 在线应用的核心问题，为此应加快有限元法的求解速度， 可从轧制过程参数积分量快速计算方法、有限元模拟简化算法、加快h e s s i a n 矩阵 的形成、大型方程组或非线性方程组的快速解法、计算机硬件条件等方面入手， 对求解板材轧制过程的f e m 计算速度进行系统研究，形成一套高效、准确的新型 有限元算法模式，解决相关基础理论和应用的关键技术问题，为f e m 在$ l n 过程 控制中的在线应用奠定基础 1 1 , l n 理论中刚塑性有限元法的概述 研究金属塑性加工问题时，弹性变形与塑性变形相比在总变形量中所占的比 例很小，例如当压下率大于1 0 时，冷轧钢的弹性变形一般不大于总变形量的5 ， 热轧钢的弹性变形一般不大于总变形量的1 刚塑性有限元法忽略了金属变形中 的弹性效应，以速度场为基本量，形成有限元列式这种类型主要有刚塑性有限 元法和刚粘性有限元法刚塑性有限元法虽无法考虑弹性变形问题和残余应力问 题，但计算程序大大简化在弹性变形较小甚至可以忽略时，采用刚塑性有限元 法可以达到较高的计算效率i l 圳 刚塑性有限元法一般是从刚塑性材料的变分原理或上界定理出发，按有限元 模式把能耗率泛函表示为节点速度的非线性函数，利用最优化理论得出满足极值 条件的最优解，即总能耗率取最小值的运动许可速度场，从中进一步利用塑性力 学的基本关系得出变形速度场、应力场以及各种变形参数和力能参数 刚塑性有限元法主要有l a g r a n g e 乘数法、罚函数法、可压缩法等计算方法这 些方法的出现使刚塑性有限元法在理论体系上更为完善，大大地促进了刚塑性有 限元法求解各类金属成形过程中参数的应用 第1 章绪论东北大学硕士学位论文 1 1 1l a g r a n g e 乘数法 刚塑性有限元法中的l a g r a n g e 数乘法是求总能耗率泛函： 在体积不可压缩条件： 下得极小值 1 1 2 罚函数法 驴= 肛q 一，p 。v f l s ， q 屯= 宅x + 宅y + ez = 0 ， 为了处理体积不可压缩条件的约束，z i e n k i e w i c z 提出了有限元分析中的罚函 数法1 9 8 2 年，g u o j il i ( 李国基) 等用罚函数法对带宽展的轧制过程进行了有限 元分析，其后，顾卓、贺毓辛用罚函数法进行了轧制方坯时端部的变形和平辊轧 制等问题的刚塑性有限元分析，卢于逑、汪时融用这种方法求解了 p p m ( p r e s s p i e r c i n gm i l l ) 穿孔过程，n a r i m 等求解了板带变形问题【5 】，成田健次郎 等用罚函数法计算了轧制过程中薄板的厚度分布 1 1 3 可压缩法 塑性力学中常采用体积不可压缩这个假设条件，由此给刚塑性有限元法求解 带来了两方面困难： ( 1 ) 体积不可压缩导致屈服于静水压力无关的结论，因而同一种变形状态可由 多种迭加上不同静水压力的应力状态所对应，反映在刚塑性有限元法的求解中， 不能由变形速度场直接求出应力场 ( 2 ) 所设的运动许可速度场需要满足体积不变约束条件，从而增加了设定初始 速度场的难度 1 9 7 3 年，大矢根在研究粉末冶金烧结材料的塑性理论时，提出了屈服条件不 仅与偏差应力的二次不变量有关，也与应力的一次不变量( 即静水压力) 有关由 此发展起来的刚塑性有限元可压缩法，使得体积不可压缩不再成为运动许可速度 场的约束条件，同时可由变形速度场直接求出应力场，从而巧妙地克服了前面提 到的两个困难 _ 东北大学硕士学位论文第1 章绪论 刘相华、吴迪、白光润等利用这一方法在国内首先求解了三维平板轧制问题 和三维高件的轧制问题，其后又进一步求解了万能孔型中轧制h 型钢和h 型钢轧 边端等问题刘全用类似的方法研究了轧件在椭圆一圆孔型系统中的变形过程之 后，可压缩法在轧制问题的求解过程中得到了广泛的应用【6 。0 1 刚塑性有限元中的这三种方法从不同角度来处理体积不变条件，都能够成功 地用于求解轧制问题 1 2 有限元软件在轧制中的应用 目前国际上在板带轧制方面有一定影响的有限元软件约有6 0 套，如a b q u e s 、 n a s t r a n 、m a r c 、a n s y s 、f o r g e 、d e f o r m 和q f o r m 等 孙铁铠等【l l 】应用大型有限元软件a n s y s l s d y n a 对双金属复合板的轧制过 程进行了数值模拟，获得了单位轧制力的分布情况，并由此计算出双金属复合板 的轧制力，计算所得轧制力与试验结果吻合较好李世芸等【1 2 】也进行了类似工作 康永林等【13 j 采用多孔材料的几何模型，利用m a r c 软件对弹簧钢6 0 s i 2 m n 半 固态轧制进行了有限元模拟 有限元软件的专业化和智能化必将成为另一个研究热点目前标准的有限元 分析软件对使用者的有限元建模和分析能力的要求相当高，并且当前有限元软件 的计算效率仍达不到在线应用所需要的速度因此，提高有限元法的求解速度是 实现f e m 在线应用的核心问题 1 3 有限元法求解存在的问题 有限元法虽然计算精度高，但是还存在计算效率低、计算时间长等缺点主 要有以下几点： ( 1 ) 追求精度，过分细化单元，加大求解规模、牺牲计算时间求解场变量思想 占据主导地位； ( 2 ) 迭代求解占据大量计算时间，导致有限元法求解时间长； 刚塑性有限元法中，采用n e w t o n 法求解非线性方程组时，需要进行反复迭代 来逼近真实解，迭代次数的多少直接影响计算时间的长短，迭代求解过程占据大 量计算时间，导致有限元法求解时间长 第1 章绪论 东北大学硕士学位论文 1 4 循环中点求积法的发展 延拓法作为数值工具用于方程求根最早是l a h a y e 于1 9 3 4 年提出的，1 9 4 8 年 他又将延拓法用于求方程组的解，a v i l a 给出了有关数值延拓法的收敛定理把延 拓法与c a u c h y 问题联系起来研究参数微分法最早是d a v i d e n k o 于1 9 5 3 年给出的， 以后又有很多作者在这方面做了大量研究m e y e r 给出了有关d a v i d e n k o 微分方程 解的存在性定理，文献 1 4 d p 给出了怎样选择求解d a v i d e n k o 微分方程的最优公式 及由此给出的数值方法 循环中点求积法就是文献 1 4 中提到的最优方法之一，循环中点求积法是解同 伦方程的参数微分法，是一种大范围收敛算法文献 1 5 】将循环中点求积法与 n e w t o n 法相结合形成了循环中点求积n e w t o n 法，循环中点求积法为n e w t o n 法寻 找一个比较好的初始值，循环中点求积n e w t o n 法可用于求解非线性方程组2 0 0 6 年黄辉【l6 j 将该算法运用于求解机械多体系统中 1 5 本文的主要研究内容 对于板材轧制问题，目前研究的重点之一是在保证计算精度的前提下，提高 计算效率使有限元法能够在线应用用有限元法求解板材轧制问题时，大型非线 性方程组的迭代求解占了一定的计算时间为加快收敛速度，本文将从以下几个 方面展开研究： ( 1 ) 提出了求解大型非线性方程组的新算法一一m p 法 ( 2 ) 计算出新算法的收敛阶，并证明新算法的收敛性及其稳定性 ( 3 ) 用f o r t r a n 语言把新算法编写成程序 ( 4 ) 用新算法求解板材轧制中形成的速度增量方程，并与n r 法进行比较 东北大学硕士学位论文第2 章刚塑性有限元理论 第2 章刚塑性有限元理论 2 1 刚塑性有限元的理论问题 应用有限元法时，一般要考虑如下几个问趔1 7 】： ( 1 ) 将微分方程或一些物理问题等转化成变分问题 ( 2 ) 变分问题解的存在及唯一性问题 ( 3 ) 对区域作剖分，构造基函数或单元形函数 ( 4 ) 形成有限元函数( 若是线性问题，该方程为线性代数方程组；若是非线性 问题，该方程一般为非线性代数方程组) ( 5 ) 提供有限元方程的有效解法 ( 6 ) 理论上论证上述每一步骤的可行性以及互相关联的问题 对于问题( 1 ) ，文献 1 8 】建立了许多力学问题中的变分原理文献 1 9 1 对问题( 2 ) ， ( 6 ) 给出了证明自提出刚塑性有限元的可压缩法后，在应用于求解各种具体金属 变形问题时，对处理好上述问题中( 3 ) ，( 4 ) ，( 5 ) 的作了大量的工作但对于问题( 5 ) 仍有许多工作需要做，对解法的改善有利于提高解的精确性、收敛速度等 2 2 刚塑性有限元的求解途径 变分原理指出了一条塑性力h - r _ 问题的求解途径，即在运动许可速度场中设法 找出满足泛函驻值或最小值条件的速度场得到了真解具体求解过程可以分为以 下几个步骤： ( 1 ) 设定初始速度场，当使用马可夫变分原理时，运动许可速度场须满足速度 边界条件和体积不变条件，使用可压塑材料变分原理或不完全广义变分原理时， 只须满足速度边界条件 ( 2 ) 找出泛函与速度场之间的关系，建立总能耗率泛函，把泛函表述为速度的 函数 ( 3 ) 利用极值理论求泛函的驻值或最小值，通常是令泛函的一阶变分为零，从 中得出速度场的真解或最佳近似解 ( 4 ) 利用塑性力学关系由位移速度场进步求出变形速度场、应力场、应变场， 第2 章刚塑性有限元理论东北大学硕士学位论文 进而对接触表面应力积分得出总轧制力和平均单位压力，由位移速度场计算出宽 展、侧边形状、前后滑、中性面等，由总能耗率换算出轧制力矩和功率 设定一个支配整个变形物体质点运动规律的速度场是非常困难的，利用有限 元中化整为零的基本思想，着眼于处理局部区域的一个小单元的时候，问题往往 会变得简单一些如果单元的尺寸足够小，数目足够多，即使在一个单元内选择 简单的速度插值函数，当把这些单元组合起来时，就能够描述非常复杂的速度分 布刚塑性有限元法就是按照这种思路将轧制变形区和总能耗率泛函进行离散化 2 3 有限元公式 可压缩刚塑性有限元法已经成功用于解决许多轧制问题真实的速度场要使 如下泛函耿极小值： 伊= 熹i l k o - , 方d a + l 等( 孵一岛) 舭l t v d l ( 2 。) = 9 p + 9 j + 9 j ， 其中：矽p 是塑性变形功率； 厅是等效应力； 享是等效变形速度； k 是一个依赖材料状态的参数，在塑性区取l ，在刚性区取0 5 ； 够是摩擦变形功率： 舰是摩擦因子； 是轧件表面和轧辊之间相对滑动速度； 岛是一个小的正常数； 7 是张力变形功率； 丁是张力，后张力取负值，前张力取正值； v 是轧件在张力方向的线速度 对可压缩材料，屈服准则表示如下： 孑= g 吒+ 吾吒嘭= 吒， ( 2 2 ) 东北大学硕士学位论文 第2 章刚塑性有限元理论 其中：吒是偏差应力张量； 是平均应力； g 是一个小的正常数( 0 0 1 o 0 0 0 1 ) 等效应变速度享表示如下： 方：， 2 + 1 6 , 2 ， ( 2 3 ) 占2 、| ，6 ( 2 3 ) 其中：乇是应变速率张量； 童。是体积应变速率 假设材料是刚塑性材料且遵守如下的本构关系： 吒= d g ”，( 2 4 ) 其中：d 是一个应变和速度函数； 聊是应变速度敏感指数 轧辊与轧件相对滑动速度： = ( u v 尺c 。s 口) 2 + 嘭 ， ( 2 5 ) 其中：匕，0 分别表示接触面上任一点x ，y 方向的速度分量； 是接触面上轧辊任一点的线速度； 口接触面上任一点对应的圆心角；，是变形区长度； 尺是轧辊半径： 日和h 分别是轧件的轧前厚度和轧后厚度； 矽咬入角 应用有限元离散口0 1 ，变形区通过连接, 个节点分离成e 个单元，每个单元的速 度场“可近似成 “= n 7 1 ，( 2 6 ) 其中：1 ，是节点速度矢量； 是形函数矩阵 应变速率张量可写成如下的矩阵形式： 第2 章刚塑性有限元理论 东北大学硕士学位论文 占= b v ， 其中：b 矩阵中的元素是对形函数求偏导 用节点速率表示言，矩阵形式为 享= 厢= f 一，其中：【c 】_ 【1 1 o 】。； z 是常矩阵： z = 42 99 24 99 0o 一个单元的平均矩阵表示如下： 否= 瓦1 善u + e ， 其中：m 是每个单元中高斯点的个数； 巨是b 矩阵在高斯点i 的值 式( 2 8 ) 刺o n 下形式： 享= f 孓 能耗率泛函( 2 1 ) 有如下离散形式： 噶归喜鼎k 叹乒霹面赢 + 乩等 k l ( 2 7 ) ( 2 8 ) ( 2 9 ) ( 2 1 0 ) ( 2 1 1 ) 窆1 = 1 瓢比 用刚塑性有限元法求解节点速度矢量的解时，根据变分原理应求总能耗率泛 函( 2 1 2 ) 的最小值 8 i、0l， 东北大学硕士学位论文 第2 章刚塑性有限元理论 2 4 总能耗率泛函的最小化 根据多元函数求极值的方法，求总能耗率泛函q , x 寸n 个未知速度分量v j 的偏导 数并令其偏导数为零则有 即 _ o ) ，( ，刊， ( 2 1 3 ) 啡昏抖+ e 南叫 v 斗尚圈叫砒卜 ( 2 1 4 ) 可得到1 个与v 有关的非线性方程，求解此非线性方程组，可得速度场( v l ，v 2 ，) 采用n r 法对非线性方程组( 2 1 4 ) 进行求解把未知数序列( h ，屹，屹) 记为矢 量矿，将泛函伊= 厂( 哥) 在；= 巩的邻域内按泰勒级数展开并取前三项得： 在此 妒= 厂( ；) q ( ；) = ( ；t ) + 夥( ；) ( ；一；t ) + 三( 函t ) rv 2 f ( v t ) ( 函t ) ， 1 ，- - 1 , k2a 1 ，t ( 2 1 5 ) ( 2 1 6 ) 诜为速度修正量，w 是( ；) 的梯度，v 2 厂是厂( ；) 的二阶导数，即h e s s i a n 矩阵： q ( ；) 是啪二次函数，令其一阶偏导数等于零，则 v q ( 哥) = v 2 ( 琉) v 破+ 夥( 碗) = o ， ( 2 1 7 ) ( 2 1 8 ) 、卜，铲一 鼢播 ，”u = 耵 俨 第2 章刚塑性有限元理论 东北大学硕士学位论文 或 v 2 厂( 诜) v 诜= 一耵( 讣 ( 2 1 9 ) 由上式解出的v ；t 可使q ( ；) 达到极值且使泛函妒接近其极值由式( 2 1 6 ) 取尼+ 1 迭 代步中： v t + 1 = v t + 肚v ，( 2 2 0 ) 其中：a 是阻尼因子 求得的新速度场玩+ 。更接近真实解，直到经过m 个迭代步之后满足收敛条件， ；一0 ，此时；。即为最终解 2 5 收敛准则 刚塑性有限元法得到的只能是满足某种精度要求的近似解在刚塑性有限元 计算中通常采取以下两种收敛准则： ( i ) 相对速度增量准则 钳嗨 其中：6 1 是预先给定的一个小的正常数( 1 0 - 4 ) ，玩是速度增量， o 碗i i = ( 喜砰) v 2 ，o 瓦l i = ( 喜砰) 啦 ( 2 ) 能耗率泛函增量准则 鲤 岛，( 2 2 2 ) 纯 其中：g ：是预先给定的一个小的正常数( 1 0 。4 ) 仇= i 仇一纯一l i ， 这里吼是计算能耗率泛函第七次的值，仇一。是计算能耗率泛函第七一1 次的值 当迭代过程同时满足式( 2 2 1 ) 、( 2 2 2 ) 式，可以保证求得的速度场和总能耗率泛 函都能达到所要求的精度 1ljjlliiilj 东北大学硕士学位论文笫3 章改进循环中点求积n e w t o n 法 第3 章改进循环中点求积n e w t o n 法 现代计算数学和优化理论提供了很多求解大型非线性方程组 f ( x ) = 0 ， ( 3 1 ) ( x r ”，f ：dcr ”寸尺”为一阶可微函数) 的方法【2 1 1 目前已经证明当目标函数的 二阶偏导数矩阵( h e s s i a n 矩阵) 可求，且能给出较好初值时，用n e w t o n 法求解非线 性问题的收敛速度很快【2 2 1 n e w t o n 法虽有很快的收敛速度，但它对初始点有一定的要求，即当初始点 接近真实解_ 时才能快速收敛，由于一是待求的，很难检验维否靠近t 初值 选择不好容易导致n e w t o n 法不能快速收敛，为了解决这一问题，延拓法提供了 一种大范围收敛的迭代算法目前延拓法已被广泛应用于求解非线性问题【2 4 之6 1 对 任意给定的初值，用延拓法常可求得问题的一个解 3 1 准备知识 延拓法的思想是对所考虑的i 司题引入参数f ，构造一簇映象日，使当，取某一 特定值时，日就是映象f ，而当，= 0 时得出方程( x ) = o 的解 确切的说，就 是构造一簇映象：d 【o ，1 】( 2 2r ”1 一r ”，代替单个映象f ，使日满足条件 h ( x ，1 ) = f ( x ) ，h ( x ，o ) = f o ( x ) ，v x ed ， ( 3 2 ) 其中：f o ( x ) = 0 的解讷已知，而方程h ( x ，1 ) = o 就是原来的非线性方程组( 3 1 ) 求 问题( 3 1 ) 的解就转化成求同伦方程 日( x ，f ) = o ，f 【o ，1 】，x d ， ( 3 3 ) 的解x = x ( ，) ，这里x ：【o ，1 】寸尺”连续依赖于f ，当扣。时，方程( 3 3 ) 的解纳已知， 当t = 1 时，方程( 3 3 ) 的解x = x ( 1 ) 就是方程组( 3 1 ) 的解也就是说，方程 日( x ( ，) ，) = 0 的解x = x ( f ) 为r ”中的一条连续曲线，或称之为r ”中的一条道路， ；e - n n x o = x ( o ) ，另一端为五= x ( 1 ) 构造满足条件( 3 2 ) 的同伦映象日，有很 1 第3 章改进循环中点求积n e w t o n 法东北大学硕士学位论文 多种方法，例如司令 h ( x ，f ) = ( x ) + ( 1 一f ) 磊( x ) ，x d ，【o ，1 1 ， ( 3 4 ) 这种同伦称为线性同伦，它可认为是圪与f 的线性组合f o 为已知映象，f o ( x ) = 0 的解堀知 在f 的定义域d 中任取知令磊( x ) = f ( x ) 一f ( ) 则 h ( x ，f ) = f ( x ) - ( 1 一t ) f ( x o ) ，x d ， o ，1 】 ( 3 5 ) 这种同伦称为牛顿同伦 令f o ( x ) = x - - x 0 ，则 日( x ，f ) = ( 1 一t ) ( x - x o ) + w ( x ) ， ( 3 6 ) 这种同伦称为不动点同伦，e ( x ) 有解确即x ( o ) = x o ，坍任意选取 也可构造同伦 m ( x ，f ) = ( 1 一f ) 彳( x 一) + 旷( x ) ， ( 3 7 ) 这里4 为，z 阶非奇异矩阵，它相当于在( 3 4 ) 中取f o ( x ) = a ( x - x o ) ，这时矗( x ) 有唯 一解x 0坍任意选取同伦( 3 6 ) 是( 3 7 ) 中取a = i ( i 为n 阶单位矩阵) 的特例 对同伦参数，作变换构造新的同伦以同伦( 3 5 ) 为例，若令f = 1 一e ，则有 h ( x ，s ) = f ( x ) 一g f ( x o ) ，s 【o ，佃) ， ( 3 8 ) 当s = o 时，h ( x ，o ) 有解粕当j 一+ o o 时，h ( x ，s ) 专f ( x ) ，即同伦方程h ( x ，s ) = o 的解x = x ( s ) _ x ( s 哼4 - 0 0 ) ，这里x 为v ( x ) = 0 的解 用延拓法求非线性方程组( 3 1 ) 的解，可构造不同的同伦h ，并通过解同伦方 程( 3 3 ) 得到但需要解决如下几个问题： ( 1 ) 映象日在什么条件下，同伦方程( 3 3 ) 有光滑的解曲线c ( 或表示为x = x ( f ) ) 存在? ( 2 ) 当，从0 到1 变化时，如何保证曲线c 不停止在任何t d ，满足p ( o ) = ，且f ( p ( ，) ) = g ( ，) ，对所有的f o ，1 】成立 引理3 2 设f ：dcr ”专r ”在开集d 的每一点是局部同胚映象，如果f 对 任意线性函数g ( ，) - - ( 1 一t ) y o + 砒具有延拓性，其中f 【o ，1 】，y o ，y l 尺”是任意的， 则f ( d ) = r ” 引理3 3 设映象f ：dcr ”一r ”在d 上连续可微，并假定存在一个开球 s = s ( ，) cd ，使得j | f ( x ) 。10 ，蜥s 成立，其中，- f ( ) lj ，则方程 f ( x ) = ( 1 - t ) f ( x o ) ，【o ，1 】，x s ， ( 3 9 ) 存在唯一解x ：【o ，l 】专scr ”，且x ( ，) 连续可微并满足微分方程c a u c h y 问题 x ( f ) = 一 乞擘( ，) ) f ( ) ， v ， o ，1 】 ( 3 1 0 ) 【x ( o ) = ， 1 引理3 3 说明同伦方程h ( x ，) = f ( x ) 一( 1 - t ) f ( x o ) = o 存在唯一解x = x ( f ) ，且 x ( t ) 是微分方程( 3 1 0 ) 的解，因此，通过求解微分方程初值问题( 3 1 0 ) 的数值解，可 第3 章改进循环中点求积n e w t o n 法 东北大学硕士学位论文 得剑非线性万程组( 3 1 ) 的解，这种方法称为参数微分法 引理3 4 1 5 1 假定： f ：dc r ”一r n 在d 上g 可导，若初值d 满足下列条 件： ( f ) i | f ( ) ，( ) 0 玎； ( i i ) 对任何x ，y s ( x o ，) 满足 i l f 。( 而) 。1e f ( x ) 一f ( y ) i l l i i x y o ； ( i i i ) h o = 切，且s ( ，r o ) cd ，其中，。1 - x 压了- 一2 h o ；则c a u c h y 问题 ( 3 1 0 ) 在【o ，1 】上有解x = x ( ，) s ( j c 0 ，r o ) ，且f ( x ( 1 ) ) = o 可用数值方法【2 8 1 求c a u c h y 问题( 3 10 ) 的数值解若令 气2 专，拈0 ，1 ， ( 3 1 1 ) ( 3 1 0 ) 的解x = x ( f ) 在点t k 的数值解讫，若用p 阶的数值求积公式计算，可表示为 讫+ l = ( 坼，一，；h ，七= o ，1 ，n - 1 ， o 为多步法计算由开始，用( 3 1 2 ) 由_ 逐 步计算到h 由于每步计算均有误差，故h 只是x ( 1 ) 的近似，它们误差可表示为 忙一x ( 1 ) 忙o h p 其中c p 是不依赖于乃地常数，p 1 只要足够大，由( 3 1 2 ) 数值方法求得的h 与x ( 1 ) 就足够近似，故可用h 作为求解非线性方程组( 3 1 ) 的迭 代初始近似，若迭代法用n e w t o n 法，则得参数微分法与n e w t o n 法组合程序： f+ 。= 矽( 致，吒一，；厅) ， l 晚“= 坼一 f ( ) - 1f ( 坼) ， ( 3 1 3 ) k = 0 ，1 ，n 一1 ；， k ，k = n ，n + 1 ， 文献 1 5 已证明了这个组合是大范围收敛的 在用参数微分法求解非线性方程q t ( 3 1 ) 时，如何选择数值求积公式( 3 1 2 ) 是实 际计算中的一个重要问题，也是本文研究的一个重点若从计算简单考虑通常可 用e u l e r 法，这时p = l ， 容易证明l i x ( 1 ) - x 1 q h ， 其中 东北大学硕士学位论文笫3 章改进循环中点求积n e w t o n 法 c l = 丢护( 扎) o p 矿巾一1 于是只要选择使忙( 1 ) 一h l l 6 且 办：专掣，这时由计剿h 共需计算个右端函黼由于捌、碱 n 很大，计算量也就大若选另外求积公式可能由于步长办增大，使较小，从而 到计算右端函数次数也可能减少当i x ( 1 ) - x v 忙艿时用何种求积公式能使参数 微分法( 3 1 3 ) 计算的工作量最省，这是方法的效率问题在文献【1 4 】中给出了计算 工作量较省的1 次2 阶求积公式： 中点求积公式 _ = 一如，( ) f ( x o ) ， t + ；2 + l ( x k - x t _ 1 ) ， 七= l ，一1 ( 3 1 4 ) 屯“= 一专 f 。( t + ； f c ， 公式( 3 1 4 ) 每步只计算一次右端函数值，其总体误差为= 1 1 h x ( 1 ) 忙c ：局2 ，1 龇一黼俐删鼽专 ( 象卜( 3 1 4 ) 叭z o j - 算次右端函数，并得到h ，把它做初始近似，再用n e w t o n 法 + 。= - i f ( 稚) - 1f ( 吒) ，k = ，n + i ， ( 3 1 5 ) 进行精确化，就可求得非线性方程组( 3 1 ) 的解五= x ( 1 ) 组合程序( 3 1 4 ) ，( 3 1 5 ) 称为中点求积n e w t o n 法在用这个程序求解非线性方程组( 3 1 ) 时要确定 厅