上節所述之 TLF 理論係以to時刻為變形前之座標Xi為基準所建 立之虛功原理方程式。本節則討論每一變形後時刻t，以變形後之座 標xi為基準而加以考慮之情形。因此前節所述之座標Xj以xj取代，Vo 及So分別以V 及S置換，而 f•oi與b•oi則以 f•i與b•i取代，即可得 Updated Lagrangian Formulation（ULF）虛功原理方程式如下：
dV
間之應變增量求其與應力之關係之應變增量理論。
換算為對降伏狀態或應變硬化等效的單向應力狀態，F,G,H,L,M,N
σ
及
圖2-1 連續體在t =t0與t=t時之變化狀況
圖 2-2 定義應力之連體座標
圖 2-3 平衡時之連體體積與曲面定義
三、數值分析
3.1簡介
近年來，有限元素法應用於解決工程的問題上，已發展的相當成 熟，而金屬板金是由無限個質點所組成的連續體，因此將其分割成有 限個元素亦稱離散化（dispersed），並假設這些元素間僅在節點上相 互連接與傳遞力，建立起節點上作用力與位移之間的矩陣關係式，並 進行計算，以求出各節點處變形和應力的近似解。
應 用 前 述 之 大 變 形 應 力 及 應 力 率 關 係 、 有 限 變 形 之 update
Lagrangian formulation和材料之構成關係式，可逐步推導得到適用於
各種金屬加工製程之統制方程式。
本研究中，係以 Cauchy應力之 Jaumann率作為構成關係式的應
力率，則updated Lagrangian形式的虛速度原理，省略體積力之作用
可以下式表示之
∫
∫
VE (σo ij−2σikε•kj)δε•ijdV +∫
VEσjkLikδLijdV = Sf •f iδνidS （3.1） 有限元素法常用之關係式為} ]{
[ }
{ν = N d• （3.2） }
d [B]{
}
{ε• = • （3.3） }
d [E]{
{L}= • （3.4） 其中，[N]為形狀函數（shape function），{d•}為節點速度（nodal
velocity），[B]為應變率-速度關係矩陣（strain rate-velocity matrix），[E] 為速度梯度-速度關係矩陣（velocity gradient-velocity matrix）。
由於虛速度原理之方程式與材料之構成關係式之變化率為線性
本節主要推導元素剛性矩陣[K]E過程，虛速度原理的離散化，以
其中k 為應力係數。
對於連體的離散化而言，標準的形狀函數為Nα，並透過此函數 將元素內部的速度νi與虛速度δνi表示如下：
∑
α ανα=
νi N i （3.18）
∑
α αδνα=
δνi N i （3.19）
αi
ν 及δν 的附標部分表示，在元素第 α 個節點上其節點速度與節點αi 虛速度的第i 個分量。
3.3 退化殼元素（Degenerated Shell Element）之分析
本研究以退化殼元素為例，推導其形狀函數，由圖3-1 可知，對 殼而言，ζ 定義為薄層（layers of lamina），而ξ及η兩常數線定義為纖 維（fibers），纖維一般而言未必與層相垂直。殼的運動是藉由使纖維 不可伸張的磚形元素（brick element）運動而獲得。每一個殼節點相 對兩個磚形節點透過不可伸張的纖維來連接。由於一個位移拘束條件 加諸於其上，使得兩個磚形元素的六個平移自由度不再各自獨立，而 僅需要三個平移與兩個旋轉座標，如圖3-2 所示。假設在殼變形的過 程中，節點的纖維保持直線，則變形幾何可藉由中間曲面與節點纖維 來描述。在此分析中，是選擇由上頂與下底連體節點的平均節點位置 所描述的中間曲面做為參考薄層曲面如圖3-1。
節點纖維向量定義如下：
bot
5
其中
由（3.33）式組合成有限元素法（3.5）之基本型式。
其中µ 為庫侖摩擦係數，而（3.37）式稱為降伏曲面或滑動曲面方程
j
陣、標準構成剛性矩陣與初始應力剛性矩陣，vαi與 f•αi的附標部份表 示在元素第α 節點之速度與表面力率的第i個分量。
處理摩擦的邊界條件使用到接觸彈簧（spring）元素之觀念如圖 3-3 所示，接觸元素的構成關係描述如（3.43）式。將（3.43）式與（3.45）
式結合，並針對所有節點的彈簧元素積分，可表示如下：
) } { } { ](
[ } { ]
[K ∆u =− Kfr ∆u − ∆ut （3.46）
其中{ u∆ 與} {∆ut}分別為為節點位移向量與模具位移向量，將（3.46）
式重新整理可得有限元素方程式之最後形式如下所示：
} ]{
[ } { ) ] [ ] [
( K + K fr ∆u = Kfr ∆ut （3.47）
而[K fr]摩擦剛性矩陣(其矩陣大小為3× ）僅存在於與模具接觸時之3 節點自由度。
3.6 廣義 rmin法之增量步驟的計算
每一次時間增量的過程中，板金成形是以其增量開始之座標為參 考基準，其採用ULF 理論並結合 Yamada 等所提之廣義之 rmin方法來 決定每一次計算步驟的增量值的大小，下列六種情況下 r 值中之最小 者定為rmin。
(1) 元素最大容許應變增量 (2) 元素最大容許旋轉增量 (3) 元素最大容許等效應力增量
(4) 材料降伏之判定 (5) 自由節點與模具接觸 (6) 接觸節點與模具分離
應用r 方法的程序及增量大小的計算包含下列步驟: min
1. 在增量開始前，先給予一個假設模具位移增量∆a；
2. 相對於前述模具的位移增量∆a，利用剛性方程式解出其虛擬 解∆u；
3. 經由前述(1)~(6)的r值，找出最小的值即為rmin的值；
4. 將所得的解 u∆ 乘上r 係數(0＜min r ＜1)作為其加權，故此增min 量的大小將變為rmin∗ u∆ ，此將用來改善成形的變形狀態、總 位移、總應力及每一個元素降伏極限的數值。
由於本研究採增量理論，每一步驟與材料之幾何形狀及元素所處 狀態有關，故以rmin決定該步驟之實際增量。在第一次的虛擬增量計 算後，經由上述各種增量限制之判定，每一步驟僅允許一元素或節點 改變其原有狀態，並選擇最小之rmin去決定該步驟之實際增量。
圖3-1 退化殼元素之座標
圖 3-2 殼元素之自由度
圖3-3 接觸彈簧元素
四、結果與討論
本計畫為整合型計畫之子計畫三，主要以開發金屬板材成形三維 解析 CAE 軟體有關應力與應變的部份為主。本計畫已如期將增量型 彈塑性大變形三維有限元素分析應力與應變之理論與程式，與其他各 子計畫所建構之程式進行整合，且完成測試。此外，亦完成 IDEAS 軟體之升級與本計畫分析程式之整合應用。且為使本程式之可靠度提 高，於去年亦完成實驗之驗證，因此本年度的計畫，將程式針對實際 加工所發生的狀況，做大幅修改及改善執行效率，已如期完成泛用型 CAE 分析軟體之目標。
今年為本程式完成階段，本計畫最後亦嘗試加入判斷材料破裂之 理論，藉由變形履歷與板材厚度變化，可判斷成形中材料最早產生破 裂之位置，應該位於料片與沖頭或料片與沖模接觸之隅角處。將本計 畫所推導之應力與應變程式，與其他子計畫之前後處理程式進行整合 後，以金屬板材彎曲成形加工之解析為例進行分析。在模具方面有限 元素網格分割是採用三角形元素，其目的是解決模具之幾何形狀問 題，以避免網格分割時在模具表面同時產生三角形和四邊形元素，增 加分析程式的複雜程度。而工件之有限元素網格分割則是採用四邊形 元素。上述模具與工件之幾何圖形皆由套裝軟體I-DEAS 建立，並透 過 Simulation 模組產生有限元素網格資料，再轉出 Universal 檔案，
以作為本整合計畫程式之前處理輸入資料，經分析後再經由 I-DEAS 軟體作後處理。以下為金屬板材彎曲成形加工為例，進行分析：
4.1 實驗與材料參數
為了瞭解金屬板材彎曲成形時沖頭負荷與沖程關係、變形歷程、
厚度分佈、應力與應變分佈關係及驗證有限元素法之適用性，本節實 施了一系列 V 型彎曲馬鞍成形實驗，實驗設備與板金材料由中國鋼 鐵公司提供。圖4-1 表示模具之模具幾何尺寸。至於，實驗之板金材 料，經使用標準的拉伸試驗方法及最小平方法以得到應力-應變特性 曲線方程式。拉伸試驗所得的材料參數與板材尺寸如下所示：
板材原始厚度：t = 2.0 mm
降伏應力：σy = 291 MPa
應力應變關係式：σ =550.87(0.0135+εp)0.1483 Mpa 楊氏係數：E = 218100 MPa
柏以松比：ν = 0.3 板金長度(L)：L = 60 mm
在一般文獻，指出板金寬度(W)於彎曲成形時，影響馬鞍外形甚 為明顯，故將實驗依板金寬度(W)：分為 12、16、24mm 三組。實驗 中以鐵氟龍作為板材之潤滑處理。照片4-1 表示 V 型彎曲馬鞍成形在 不同成形階段之成品。
4.2 數值分析
本研究是採用四節點四邊形殼元素推導剛性矩陣，並以 CAD 套 裝軟體進行前後處理，將所建立的模具與板材經過網格分割後轉換成 數據檔，輸入增量型彈塑性大變形三維有限元素分析程式中進行數值 解析，將所模擬出之沖頭負載與沖頭沖程之結果與實驗結果做比較，
並將分析結果輸出至 CAD 套裝軟體顯示變形圖以及應力與應變分佈 圖。
V 型彎曲馬鞍成形模擬時所建構沖頭與下模之網格分割如圖 4-2、圖 4-3 所示。在模具方面有限元素網格分割是採用三角形元素。
圖4-4 顯示將板材(Blank 3)進行網格分割，產生四邊形元素，而在彎 曲變形激烈的沖頭與下模間之部份作較細的元素分割，以求得精確的 計算結果。在板材模擬時，需給予適當之邊界條件，邊界條件須設在 節點上，即於節點不能發生旋轉或位移的方向上給予拘束。在板材X 軸上的節點為Y 方向的位移拘束，而 Y 軸上的節點為 X 方向的位移 拘束，此外中心節點為X 與 Y 方向的位移拘束，Z 方向為旋轉拘束。
表4-1 顯示模具與板材有限元素網格分割之相關數據。
V 形彎曲馬鞍成形過程中，板材會與沖頭及下模接觸，因此需瞭 解成形過程中節點是否與模具接觸或分離，故可將節點分為接觸節點 與自由節點。圖 4-1 顯示板材與模具系統組合時之邊界條件，與模具
不接觸之節點為自由節點採用空間固定座(X， Y， Z)，而與模具接 觸之節點則採用埋入座標(ξ，η，ζ)，座標方向採用右手定則，故 Y 與η方向朝紙內。
由於各個節點之接觸狀況，會隨成形過程中板材的變形而有所改 變，故在每一位移增量步驟計算時，均需對接觸之節點檢查其節點力 之法線分量是否小於或等於零，若小於或等於零時，則下一位移增量 步驟需將此節點之邊界條件更改為自由節點之邊界件。而原來為自由 節點亦需檢查其幾何位置是否與模具接觸，若接觸則在下一位移增量 步驟計算時，將此節點之邊界條件更改為接觸節點之邊界條件。以上 所述之檢查方式是以 extended r-minimum 來處理。
由於各個節點之接觸狀況，會隨成形過程中板材的變形而有所改 變，故在每一位移增量步驟計算時，均需對接觸之節點檢查其節點力 之法線分量是否小於或等於零，若小於或等於零時，則下一位移增量 步驟需將此節點之邊界條件更改為自由節點之邊界件。而原來為自由 節點亦需檢查其幾何位置是否與模具接觸，若接觸則在下一位移增量 步驟計算時，將此節點之邊界條件更改為接觸節點之邊界條件。以上 所述之檢查方式是以 extended r-minimum 來處理。