一步成形模擬中一種新的松弛因子搜索算法

張亞奇，王 舒，劉偉杰，張向奎，胡 平，鄭國君

（1. 大連理工大學 運載工程與力學學部 汽車工程學院，大連 116024；2. 大連理工大學 工業裝備結構分析國家重點實驗室，大連 116024；3. 大連迪克斯科技有限公司，大連 116024）

一步成形模擬中一種新的松弛因子搜索算法

張亞奇1,2，王 舒2，劉偉杰1,2，張向奎1,2，胡 平1,3，鄭國君3

（1. 大連理工大學 運載工程與力學學部 汽車工程學院，大連 116024；2. 大連理工大學 工業裝備結構分析國家重點實驗室，大連 116024；3. 大連迪克斯科技有限公司，大連 116024）

0 引言

車身覆蓋件沖壓成形是汽車生產中一個很重要的環節，制約著汽車質量的提高和改型換代周期的縮短，傳統的依靠經驗進行反復調試的模具設計方法已經不能滿足高質量、高效率和低成本的要求。隨著CAD和CAE技術的發展，板料沖壓成形模擬技術越來越多的被用于預測變形路徑和成形缺陷如裂紋和起皺。

板料成形模擬技術是基于大變形彈塑性有限元理論的一種數值模擬算法，主要可分為增量理論和全量理論[1]。增量理論的思想是將整個變形分為若干個小增量步，第i步的求解是建立在第i-1步所得的結果之上的，迭代到最后可以得到最終的模擬結果，增量法最主要的優點是模擬精度高，但其計算復雜，數據準備繁瑣，計算時間長，不能滿足在模具最初設計階段快速估計產品可制造性的要求。

一步成形模擬方法是一種基于全量理論的有限元逆算法[1]，該方法最早由Kobayashi等人[2]提出，J. L. Batoz，Y. Q. Guo 和H. Hub都進行過研究[3～5]。這種方法假設變形過程是比例加載的，只考慮初始毛胚板和變形終了構型兩個狀態，忽略中間狀態，建立有限元方程。其突出優點是計算效率高并且擁有較好的精確性。既可以在產品的初始設計階段快速預測沖壓件的成形性，還可以在制造階段快速估計毛胚板料的尺寸和形狀，應用前景廣泛。

一步成形模擬最終需要求解一個大規模非線性方程組，一般采用Newton-Raphson（N-R）算法進行迭代求解[5]，但是經典的Newton-Raphson迭代對初始解的精度要求比較苛刻，初始解必須在真實解附近，迭代才能收斂。對于一些復雜問題尤其是高度非線性的問題，很難找到讓經典的N-R迭代收斂的初始解。因此為了放松對初始解的要求，提高迭代穩定性和收斂速度，通常需要引入松弛因子。松弛因子的取值對求解的穩定性和計算效率有很大的影響。如果松弛因子取值太大，則容易發生“指數溢出[6]（Exponential Overflow）”，即計算將以無結果中止；如果松弛因子太小，會使計算效率變差，增加程序的運行時間。

在一步成形模擬求解過程中，選取合適的松弛因子比較困難，而且到目前為止國內外對這方面的研究還比較少，國外還沒有文章對其進行專門的研究，即便是涉及到松弛因子的文章，大部分也沒有給出算法的具體步驟。Sven K. Esche等[7]比較了多種松弛因子線性搜索方法；L. Armijo[8]首先采用“折半搜索”的方式確定松弛因子，該方法也被稱為Armijo Rule，在一步成形模擬中有廣泛應用。國內吉林大學的那景新和陸善斌等人[9]提出了松弛因子的黃金分割算法，用一維搜索0.618法確定松弛因子的取值；那景新，高華等人[1]采用插值拋物線近似替代目標函數，優化松弛因子的取值。上述兩種算法的迭代收斂性和計算效率都比固定松弛因子的算法要好，但都比較復雜，需要的附加計算較多。現有的處理方案一般采用“折半搜索”算法，但由于其算法粗糙，不適合處理一些復雜問題。本文在權衡迭代穩定性和計算效率的基礎上提出了一種新的松弛因子搜索算法，并驗證了其有效性。

1 一步成形基本理論

一步成形有限元方法的基本思想是假設變形過程是比例加載的，通過比較初始毛胚板和變形終了構型兩個狀態，確定零件中的節點在毛胚平板中的對應位置，從而得到毛胚平板的形狀和尺寸，進而得到零件中的應變、應力分布和厚度分配。在零件變形終了構型上建立虛功方程[5]：

選取松弛因子的方法有很多[6]，通常的辦法是定義一個目標函數φ，選取合適的松弛因子使該函數的值隨著迭代的進行逐漸減小，下面是常用的一種目標函數形式：

對于某些問題，松弛因子經過一到兩次減半就可以使目標函數減小，“折半搜索”在解決這類問題時效果較好，但是作為一種初級的線性搜索算法，“折半搜索”由于搜索模式單一，只有減半這一種方式，造成其在處理一些問題時太粗糙。因為松弛因子只能取1、0.5、0.25、0.125、0.0625、0.05這六個值中的某一個，這很有可能錯過一些較大的松弛因子（如0.8，0.9），例如在某個迭代步中，最佳的松弛因子取值為0.75，顯然“折半搜索”會錯過0.75，而得到一個較小的松弛因子，影響計算效率。

2 新的松弛因子搜索算法

不難看出，最理想的情況是可以找到使目標函數在區間上取最小值的松弛因子[10]，但一般情況下目標函數非常復雜，搜索最小值點的過程需要大量試算，這會增加搜索松弛因子的時間消耗，降低程序的計算效率。本文在權衡迭代收斂性和計算效率的基礎上，提出了一種基于局部極小值搜索的新松弛因子搜索算法，算法的思想是用局部極小值點代替最小值點，并通過一個簡單的線性搜索得到目標函數的局部極小值點。算法的具體操作步驟如下：

1）首先令松弛因子為1，計算目標函數的值，如果相比上一個迭代步有所下降，則取松弛因子為1，否則開啟目標函數的局部極小值搜索算法。

2）目標函數的局部極小值搜索算法：

（1）沿0向1的方向進行搜索，確定a和c的值，滿足條件；

（2）繼續搜索，確定b的值，滿足條件，同時保持的關系不變。

算法的思想是，用一個長度為0.1的區間沿著0向1的方向移動，每次移動步長，直到找到滿足條件的區間或者達到預先設置的最大值1。

算法的流程如圖1所示。

圖1 新松弛因子搜索算法流程圖

3 數值算例

為了驗證本文提出的新算法的有效性，作者編制計算機程序實現了新算法，并將其加入到KMAS/One-Step程序中。本文選取了U型梁、方盒等四個算例進行測試分析，并與“折半搜索”算法結果進行比較。測試所用硬件配置為Intel Core i5/2.53GHz，8GB內存。

各測試零件均采用三角形單元劃分有限元網格，具體劃分情況如表1所示。出于篇幅考慮，這里只給出了方盒和測試零件四的有限元模型圖如圖2所示。

表1 測試零件的有限元網格劃分情況

圖2 方盒（左）和測試零件四（右）的有限元模型

測試零件所用材料為鋼，具體的材料和幾何參數如表2所示。

表2 測試零件的材料和幾何參數

由于篇幅限制，這里只列出部分計算結果圖。在四個算例中，方盒最具有代表性，下面給出了采用本文提出的新松弛因子搜索算法計算所得的方盒的厚度分布圖和初始毛胚平板形狀如圖3所示。

圖3 方盒的厚度分布圖（左）和毛胚形狀（右）

表3為采用不同算法對四個算例進行測試的結果，包括零件的最大和最小厚度（mm），迭代收斂所需的步數和迭代消耗的CPU時間（s）。

根據表3，我們可以看出，對于選取的四個算例，采用兩種松弛因子選取算法迭代都可以收斂，而且計算所得的零件最大（最小）厚度相差無幾，這表明本文提出的新松弛因子搜索算法不影響Newton-Raphson迭代的收斂性和計算精度。

新松弛因子搜索算法可以有效減少收斂所需的迭代步數，尤其是在處理一些深拉延問題時。例如圓盤的迭代步數從11減少到5。這主要是因為新算法搜索目標函數的局部極小值點，而“折半搜索”只要求目標函數下降。但是新算法要比“折半搜索”復雜，搜索松弛因子的時間也長，因此在解決某些問題時，時間略有上升，這是正常的。為了減少不必要的附加計算，本算法并不是每次都開啟局部極小值搜索，而是分情況處理，在松弛因子取1目標函數無法下降時，才進行搜索。

表3 不同松弛因子算法的測試結果

對于規模較大的問題，解方程組的時間可以占到一步成形模擬程序總時間的90%以上[9]。本文提出的新算法可以有效減小收斂所需的迭代步數，雖然目標函數的局部極小值搜索會帶來一些附加計算，但是解線性方程組的時間復雜度為o(n3)，新松弛因子搜索算法的時間復雜度為o(n2)，因此可以預見新算法在解決大型問題時，計算效率的提高將更加明顯。新算法可以使目標函數隨著迭代的進行逐漸減小，因此可以放松對初始解的要求，提高N-R迭代的穩定性。

總體來說，本文提出的新的松弛因子搜索算法，比傳統的“折半搜索”算法有更好的迭代穩定性和計算效率，尤其是在處理規模較大的問題時。

4 結論

本文研究了一步成形模擬中非線性方程組的N-R迭代求解，探討了松弛因子取值對N-R迭代收斂性和計算效率的影響，并在權衡迭代穩定性和計算效率的基礎上，提出了一種新的基于局部極小值搜索的松弛因子搜索算法。選取了U型梁和方盒等四個實例，驗證了新算法的有效性。通過與經典的“折半搜索”算法的結果進行比較，本文發現新算法可以有效減少迭代收斂所需的步數，提高計算效率，在處理大型問題時，優勢將更加明顯。本文提出的新松弛因子搜索算法已經集成到Siemens PLM NX/One-Step和KMAS/One-Step中。

[1]那景新,高華,張麗,胡平.一步成形模擬方法中松弛因子選取算法[J].吉林大學學報(工學版),2005,3:292-296.

[2]Kobayashi S, Kim J H, Mechanics of Sheet Metal Forming,New York: Plenum Press,1978:341-365.

[3]Yang Shiyong and Nezu Kikuo,“Application of an inverse FE approach in the concurrent design of sheet stamping”,Journal of Material Processing Technology,1998:79,86-93.

[4]C. H. Lee and H. Huh,“Blank design and strain estimates for sheet metal forming processes by a finite element inverse approach with initial guess of linear deformation”, Journal of Materials Processing Technology,1998:82,145-155.

[5]Y. Q. Guo and J. L. Batoz,“Finite element procedures for strain estimations of sheet metal forming parts”,International Journal for Numerical Methods in Engineering,1990:30,1385-1401.

[6]C. T. Kelley, Solving Nonlinear Equations with Newton’s Method, Philadelphia:Society for Industrial and Applied Mathematics, 2003:1-24.

[7]Sven K. Esche, Gary L. Kinzel and Taylan Altan,“Issues convergence improvement for non-linear finite element programs”,Int. J.Numer.Meth.Engng.1967:40,4577-4594 .[8]L.Armijo,“Minimization of functions having Lipschitz-continuous first partial derivatives”,Pacific J.Math.1966:16,1-3.

[9]陸善彬.一步模擬中沖壓方向、單元模型及松弛因子算法研究[D].吉林大學,2004.

[10]Y. Q. Guo,J.L.Batoz, H. Naceur, S. Bouabdallah,F.Mercier, O. Barlet, “Recent developments on the analysis and optimum design of sheet metal forming parts using a simplified inverse approach”,Computers and Structures,2000:78,133-148.

A new relaxation factor search algorithm in one-step forming simulation

ZHANG Ya-qi1,2, WANG Shu2, LIU Wei-jie1,2, ZHANG Xiang-kui1,2, HU Ping1,3, ZHENG Guo-jun3

本文著重研究了一步成形模擬中松弛因子取值對迭代收斂性和計算效率的影響，在權衡迭代穩定性和計算效率的基礎上，提出了一種新的松弛因子搜索算法，然后用計算機程序實現了算法，并將其加入到KMAS/One-Step中，通過U型梁、方盒等四個計算實例，與經典的“折半搜索”算法結果進行了比較，驗證了新算法的有效性。

板料沖壓成形；一步成形；松弛因子

張亞奇（1989 -），男，碩士研究生，主要從事汽車覆蓋件沖壓成形仿真研究。

TP391

A

1009-0134(2013)06(上)-0014-04

10.3969/j.issn.1009-0134.2013.06(上).05

2013-04-09

863計劃項目（2013AA040501）