基于理想形變理論的有限元逆算法在板料成形中的應用

朱光成，吳建軍

（西北工業大學 機電學院，陜西 西安 710072）

基于理想形變理論的有限元逆算法在板料成形中的應用

朱光成，吳建軍

（西北工業大學 機電學院，陜西 西安 710072）

基于理想形變理論，研究了金屬板料成形的有限元逆算法，并開發了計算程序。采用線性三角形膜單元和厚向異性的剛塑性材料模型，計算了一個帶凸緣的方盒形件，并與Dynaform的一步逆向法和增量法的計算結果進行了比較。實例分析結果表明，在工程允許精度范圍內，本文方法能夠有效分析零件的成形性能。由于該方法計算速度快，所以可用于零件的早期設計。

機械制造；板料成形；理想形變理論；有限元

1 引言

隨著計算機技術的迅速發展和有限元理論的日趨成熟，有限元方法在坯料計算及板料成形中的應用越來越廣泛。目前，常用的板料成形數值模擬方法主要有兩種：基于增量理論（又稱流動理論）的有限元增量法和基于全量理論（又稱形變理論）的有限元逆算法。增量法可以全面考慮各種影響因素，是目前比較精確的方法。但在產品及工藝設計的早期階段，通常只有產品形狀而無模具信息，因而，使用增量法進行數值模擬比較困難。對設計人員來講，最重要的是如何根據產品或已經完成工藝補充的沖壓件幾何形狀來快速預測它的毛料展開形狀和零件的厚度分布等信息，以此來預測零件的可成形性和工藝的可加工性。因此，K.Chung和O.Richmond將形變理論與極值功路徑相結合，提出了理想形變的概念[1]。本文將其應用于板料成形過程，開發了有限元逆算法程序。采用線性三角形膜單元和厚向異性的剛塑性材料模型，計算了一個帶凸緣的方盒形件，驗證了本文方法的有效性。

2 理想形變理論[2]

理想形變的基本假設為：變形沿著最小塑性功路徑進行，材料為剛塑性，滿足R.Hill’79屈服準則。整個變形過程達到最小塑性功路徑的條件有兩個：塑性變形的主伸長對應于固定的主物質線；塑性變形的主自然應變比保持不變。

理想形變理論認為變形體在整體塑性功取得相對極值的條件下變形較為均勻。變形取得相對極值的條件為：在成形的最后一刻，變形體在邊界上僅受法向力作用而處于平衡狀態，且變形分布較為均勻，并以變形過程中消耗的整體塑性功最小作為應變優化的基礎。文獻[1]中已經對該假設做了證明。將變形過程看作一個勻速的過程，考慮到塑性變形體積不可壓縮的假設，包括不均勻變形的塑性功可表示為：

式中：σ——Cauchy應力；

D——應變速率；

σˉ——等效應力；

εˉ——等效應變；

V0——初始體積。

當物質微元的變形路徑確定之后，塑性功就只由位移確定，可看作是位移的函數W=W（U）。

3 理想形變理論的有限元表達及實現

3.1 變形分析[3][4]

在Lagrange描述法中，物質坐標系在初始構形中與整體坐標系重合。變形前，物質坐標系為直角坐標系，各坐標軸之間的夾角為直角，且基矢量為單位矢量。變形后，物質坐標系不再是直角坐標系。根據Kirchhoff假設，變形前的中面法線在變形后仍然垂直于中面。用P0和Pt分別表示物質點P在變形前后的空間位置，我們可以得到變形后物質坐標系在P點的基矢量t→1、t→2和t→3的表達式為：

由于已知毛坯形狀為平板，所以單元的右Cauchy－Green變形張量[C]可表示為：

[C]的特征值為單元面內兩個主伸長λ1、λ2的平方，則有：

根據塑性變形中體積不可壓縮的假設，可得板料上各點的厚向伸長量λ3為：

右Cauchy－Green變形張量[C]也可由坐標旋轉矩陣[M]與主伸長λi（i=1，2，3）表示為：

其中，坐標旋轉矩陣[M]為：

θ是λi與最終構形的局部坐標系x軸的夾角，求解（6）式得：

由右Cauchy－Green變形張量[C]容易得到對數主應變表達式為：

那么單元沿著厚度方向任意一點的大變形對數應變表達式為：

3.2 塑性功及其極值[5]

將理想形變應用于板料成形設計時，考慮到零件的最終構形和坯料的初始輪廓（一般是平面）是已知的，成形過程的邊界條件未指定，最終構形上的物質點必須滿足幾何約束條件x3=x3（x1，x2）,再根據U= x-X（X，x分別是物質點在初始構形和最終構形上的節點坐標），取得整體塑性功極值的靜力平衡方程：

采用Newton－Raphson迭代法，構造方程組，求{R（X）}=0，有：

定義切線剛度矩陣為：

求解方程組式（13），就可得到初始構形上節點坐標的調整量{d Xi}，調整后的節點坐標為：

式中：β——減速因子，取值一般為0～1，本文取為1。

以上各式的具體表達形式可以非常容易的求得，故這里不再詳細介紹。

3.3 程序流程圖及毛料初始解的獲取

在Visual C++6.0環境下開發有限元逆算法程序，程序的計算流程圖如圖1所示。

有限元逆算法需要一個毛料的初始猜測值。一個合理的毛料初始解將直接影響到求解的收斂性和準確性。本文采用幾何映射法[6]獲得該初始猜測值。

3.4 收斂準則及收斂性的提高[7][8]

采用Newton－Raphson迭代法必然存在兩次迭代之間滿足什么準則可以終止本次迭代的問題。根據有限元逆算法的實際計算情況，本文選擇位移準則來進行收斂性判斷，計算節點位移向量的范數：

式中：n——節點總數；

ρ——給定的收斂判斷因子，對位移收斂準則ρ一般取10－3。

當網格中出現畸變單元或單元大小突變比較嚴重時，可能使切線剛度矩陣呈現病態，程序算法無法收斂。為了提高逆算法的網格處理能力，添加阻尼矩陣來改善切線剛度矩陣的病態性。

由式（13）和式（14）可得：

圖1 程序計算流程圖

在切線剛度矩陣中加入阻尼矩陣δi[I]，[I]為單位矩陣，δi[I]為對角矩陣，則方程組的解為：

在本文程序算法中取：

式中：Kjj——第i迭代步中切線剛度矩陣第j行第j列的對角線元素值。

4 有限元逆算法應用實例

圖2 劃分有限元網格的產品最終構形

以尺寸為50mm×50mm×20mm、各處圓角半徑均為6mm的帶凸緣的方盒形件為例進行計算。圖2是劃分三角形單元后產品的最終構形。單元尺寸為2.5mm，盒形件劃分為2384個單元和1237個節點。材料的應力應變曲線σˉ=520.4εˉ0.232MPa，厚向異性指數為r= 1.65，毛料的初始厚度為t=1.0mm。采用本文方法、Dynaform的一步逆向法和增量法的計算結果見3所示，三種方法計算的零件厚度的云圖分布規律基本上是一致的，只存在個別區域的差異。本文方法與Dynaform一步逆向法相比：增厚最大誤差為2.60%，減薄最大誤差為 0.82%；本文方法與Dynaform增量法相比：增厚最大誤差為12.84%，減薄最大誤差為6.77%。采用Dynaform增量法用時1976s，而本文程序迭代6次就收斂，用時僅493s，本文方法明顯比增量法快很多。本例計算是在Windows XP系統上完成的，CPU是AMD Processor 5050e/2.61GHz，內存為2GB。

在板料成形中，零件的厚度變化是一重要的因素，因為變薄嚴重會有破裂的可能，增厚嚴重會有起皺的可能，所以算例中主要對零件的厚度進行了分析和比較。

5 結論

根據理想形變理論結合板料成形的具體情況，給出了有限元逆算法的數學公式和有限元表達，開發了工程分析模塊。同Dynaform的一步逆向法和增量法的計算結果進行對比，本文方法計算的零件厚度分布情況同實際情況更加吻合。由于本文方法在計算過程中簡化了工藝參數對板料成形的影響，這樣的結果在工程分析中是可以接受的。通過實例可以證明，有限元逆算法能夠較好的滿足設計的需要，是一種高效的數值模擬工具。雖然逆算法模擬精度相對于增量法低了一些，但隨著算法的不斷完善，其模擬精度還可以得到進一步的提高。

圖3 零件的厚度分布

[1]K.Chung，O.Richmond.Idealforming-Ⅱ.Sheetforming withoptimum deformation[J].Int.J.Mech.Sci.，1992，34（8）:617-633.

[2] 蘭 箭，董湘懷，李志剛.用有限元逆算法計算板料成形毛坯形狀和應變分布[J].塑性工程學報，2001，8（2）:60-62.

[3] 吳建軍，陳衛彬，李順平.復雜形狀拉深件快速展開與成形模擬[J].機械科學與技術，2004，23（3）:370-375.

[4] 陳衛彬.基于反向模擬法的拉深成形性分析技術研究[D].西安:西北工業大學出版社，2003.

[5] K.Chung，S.Y.Lee，F.Barlat.Finite elementsimulation ofsheetforming based on a planaranisotropic strain-rate potential[J].InternationalJournalofPlasticity，1996，12（1）:93-115.

[6] 吳建軍，楊漢平.考慮外法向的鈑金零件展開方法研究[J].中國機械工程，2006，17（15）:1546-1549.

[7] Y.Q.Guo，J.L.Batoz，H.Naceur.Recent developments on the analysis and optimum design of sheet metal forming parts using a simplified inverse approach[J].Computers and Structures，2000，78:133-148.

[8] 徐國艷，高 峰，杜發榮，張立玲.基于快速有限元分析的沖壓件毛料展開[J].機械科學與技術，2006，25（4）:386-389.

Application of finite element inverse approach in sheet metal forming process based on ideal forming theory

ZHU Guangcheng，WU Jianjun
（Northwestern Polytechnical University, School of Mechatronics, Xi'an 710072，Shanxi China）

Based on ideal forming theory, the finite element inverse approach for sheet metal forming has been explored, and the computer program has been implemented. By use of the linear triangular membrane element and the thickness anisotropy rigid-plasticity material model，a flanged square box has been analyzed and its result has been compared with the results obtained from MSTEP and incremental approach of Dynaform software. The results show that the finite element inverse approach is an effective approach for sheet metal forming only if its error is allowable during the engineering calculation. Because of its rapid computing speed, this method can be used at the early design stage of sheet metal forming parts.

Sheet metal forming；Ideal forming theory；Finite element

TG386

A

1672－0121（2011）03－0049－03

國家863計劃資助項目（2008AA04Z120）

2011－02－22

朱光成（1985－），男，碩士在讀，主攻板料成形數值模擬技術