飛機整體結構件加工變形的初始殘余應力- 初始幾何誤差耦合影響與控制

林鋒， 姚婉， 秦國華， 葉海潮， 陶江

(南昌航空大學 航空制造工程學院， 江西 南昌 330063)

0 引言

航空航天領域的發展水平在當前已成為衡量一個國家綜合國力的重要指標之一。人們對于航空材料性能的要求越來越高，而航空材料是一類特殊且要求比較嚴苛的材料。鋁合金由于較高的比強度、低密度、良好的加工性能以及耐腐蝕等優點，廣泛應用于航空領域[1]。此外，為了提高飛行機動性，許多骨架零件大量采用整體結構件[2]。在加工過程中，毛坯材料的去除率高達90%，使得零件的剛性變差，易發生加工變形且加工精度難以得到保證，每年因為加工變形而導致報廢的零件數不勝數，造成巨大的經濟損失，也限制了我國航空業的進一步發展。

零件在高速加工時變形的主要致因包括初始殘余應力和加工殘余應力[3-4]。Li等[5]介紹兩種殘余應力的產生原因及影響因素，并分析了這兩種殘余應力的預測、控制方法及其適用范圍。據統計，航空零件加工過程中發生變形的約90%可以歸結于毛坯內部的初始殘余應力[6]。Sun等[7]根據板變形理論進行了不同初始殘余應力作用下的銑削變形實驗，發現初始壓應力和加工殘余應力耦合下板的變形會增加，而初始拉應力和加工殘余應力作用下板的變形會減小。Huang等[8]分析毛坯初始殘余應力和銑削應力耦合作用對零件變形的影響，對比化學銑削實驗結果后認為，在高速銑削情況下，殘余應力對零件變形的影響遠遠大于銑削應力。Wang等[9]和Lu等[10]基于航空整體薄壁件加工變形提出一種僅在殘余應力作用下的加工變形預測模型，該模型可以較準確地預測出殘余應力對加工變形的影響以及殘余應力分布中行測的偏移對變形的影響。Li等[11]利用零件的簡化模型進行有限元建模，分析了雙向殘余應力對變形的影響，發現初始軸向殘余應力對變形的影響大于徑向初始殘余應力對變形的影響。Cerutti等[12]研究初始殘余應力的作用下加工順序對加工質量的影響，提出一種計算工具，利用該工具根據初始工件和加工順序來優化加工質量。Li等[13]依據兩個典型航空零件在加工后72 h內的變形規律提出三步加工工藝，以控制薄壁零件的加工變形，實驗結果表明在第1次加工后將零件自由放置72 h可以有效降低約48%變形值。秦國華等[14]依據彎曲變形理論推導出零件加工變形的力學模型，將加工要求作為加工變形的約束條件，結合夾逼準則以及引入單位變量，將加工變形轉化為齊次線性不等式，通過改變約束條件來求解方程組，以控制加工變形。

國內外學者對于整體結構件加工變形方面的研究較為深入，但目前對于毛坯的初始幾何誤差這一方面的研究較少。在實際生產加工過程中，由于材料的非均勻性，毛坯表面不可避免地存在初始幾何誤差。

為了進一步揭示毛坯不同數值大小初始幾何誤差對零件變形的影響，本文針對不同的毛坯初始幾何誤差建立相應的有限元模型，并選取零件在毛坯內的不同位置進行研究，得到不同初始幾何誤差以及零件位于毛坯內不同位置對零件加工變形的影響規律。再以零件加工后的變形最小為目標，利用變向迭代優化方法搜索零件在毛坯內的最佳位置，即當零件處于毛坯內的該位置時，零件的變形為最小，進而可以通過優化零件的加工位置來控制變形。

1 零件變形的分析方法與驗證

在毛坯制造過程中，由于材料力學性能的非均勻性，使得鋁合金厚板內不僅存在初始殘余應力，而且存在初始幾何誤差(即預拉伸板的平面度)。而在高速銑削毛坯成型為飛機整體結構件的工藝過程中，毛坯內部的初始殘余應力得以釋放，破壞了飛機整體結構件的靜力平衡狀態，造成飛機整體結構件發生變形。

1.1 變形理論

由于鋁合金厚板沿厚度方向的殘余應力值非常小，完全可以忽略不計[15]。因此在自然狀態下，鋁合金厚板內部的初始殘余應力在任意截面都處于自平衡狀態，即內部的力和力矩處于自平衡(見圖1)，可以表示如下：

圖1 毛坯初始殘余應力分布Fig.1 Initial residual stress distribution of blank

(1)

(2)

式中：X、Y、Z分別為軋制方向、橫向方向、厚度方向；x、y、z分別為沿X軸方向、Y軸方向、Z軸方向的位置；FX(x)、FY(y)分別為X軸方向、Y軸方向的力；MX(y)、MY(x)分別為繞X軸方向、Y軸方向的力矩；σX(z)、σY(z)分別為沿X軸方向、Y軸方向的殘余應力；ΩOYZ、ΩOZX分別為x、y處的截面面積。

梁類結構件變形的特點為：軋制方向Y軸上彎矩產生的彎曲變形遠大于其他方向上的變形，且扭轉變形極小。根據彈性變形理論，可知軋制方向上的應變為

(3)

式中：ν、E分別為零件材料的泊松比、楊氏模量。

因此，為了更加清晰地分析軋制方向上的變形，定義零件軋制方向上的等效應力如下：

σ(z)=σX(z)-νσY(z).

(4)

在殘余應力分布σ(z)的作用下，梁類結構件發生了彎曲變形(即撓度)。由于一般整體結構件不同位置的截面是不相同的，在Y軸的不同位置，梁的每個彎矩也是不同的。因此不同截面繞X軸方向的力矩可以表示為

(5)

式中：Ωw(y)表示零件在任意位置y處的截面面積；z0(y)表示零件任意位置y處的截面中性軸高度，如圖2所示。圖2中，h為零件在毛坯中的加工位置，X′、Y′、Z′為經過零件中性層的3個坐標軸，UY(yA)、UZ(yA)、UZ(yB)分別為A點Y軸、A點Z軸B點Z軸方向的位移，UZ(y)為任意位置y處Z′軸方向的位移。

圖2 受力狀態Fig.2 Stress state of parts

根據材料力學可知，所選取截面處彎曲變形向下凹時，截面上的彎矩為正、反之為負。因此，對于梁結構，任意位置y處上彎矩引起的撓度UZ(y)可以表示為

(6)

式中：IX(y)為任意位置y對應的截面慣性矩。

1.2 三框件變形計算

圖3所示為飛機整體三框結構件，其材料選用7075-T7451鋁合金板材，彈性模量E=71 700 MPa，泊松比ν=0.33.圖3中，l、w、hf為三框結構件的長、寬、高，l1、l2、l3分別為3個框的內長，hb為腹板的厚度，tf為框的厚度。

圖3 三框結構件示意圖Fig.3 Schematic diagram of three-frame structure

選用的毛坯外廓尺寸為1 200 mm×120 mm×60 mm，毛坯初始幾何誤差Δ=0.6 mm. 在采用裂紋柔度法對殘余應力進行測量時，每隔1.5 mm測量一次，共測出40組應力值，最后得到毛坯內初始殘余應力分布曲線，如圖4所示。另外，零件在毛坯中的厚向加工位置為正向h=16.5 mm，如圖2(a)所示。

圖4 殘余應力分布曲線圖Fig.4 Distribution curves of residual stress

如圖5(a)所示，將毛坯沿厚度方向劃分為40層，由于零件的加工位置為16.5 mm，故處于毛坯的第12層～第31層，共計20層。

圖5 有限元仿真結果Fig.5 Finite element simulated results

在有限元仿真時，按照毛坯和零件尺寸進行建模，建模過程如下：

1)在無初始幾何誤差的毛坯中分割出零件和去除材料兩個部分。

2)將毛坯均勻劃分為40層。

3)選取單元類型C3D20R進行網格劃分，設置網格密度15 mm.

4)給毛坯底部施加均布載荷，使其產生一定量的初始幾何誤差。

5)重新導入具有初始幾何誤差的毛坯(不導入由均布載荷帶來的附加內部應力)，利用網格編輯功能調整局部節點位置，確保零件在毛坯內處于指定位置。

6)給毛坯施加初始殘余應力。

采用“單元生死”技術模擬銑削過程的材料去除效應。將實驗測量得到的初始殘余應力值賦到毛坯中，得到毛坯存在0.6 mm(方向為上凸)初始幾何誤差，有限元仿真中原始模型和簡化模型的結果如圖5所示。

對比圖5(b)和圖5(c)可以看出，簡化結構與原始結構的有限元仿真結果非常相近，最大的相對誤差僅有3.5%. 因此在解析法中可以利用簡化結構來進行計算。

假定毛坯具有初始幾何誤差Δ，則在極坐標系統中，根據圖4的應力值，結合(4)式計算出等效殘余應力的等高曲線為

σ(ρ)=5.445 3+2.498 9ρ-2.211 7ρ2+
0.374 2ρ3-0.027 2ρ4+0.001 02ρ5-
0.000 020 7ρ6+0.000 000 215ρ7-
9×10-10ρ8，

(7)

式中：ρ為位置z處由Δ形成的圓弧半徑；σ(ρ)為半徑ρ處的應力值。

圖6所示為應力等高線。圖6中，R為毛坯底面處由Δ形成的圓弧半徑，θ為半徑ρ與Y″軸的夾角。根據圖6的幾何關系，由于應力等高曲線為以O″為圓心的圓弧，則有

圖6 應力等高線Fig.6 Stress contour

(8)

(9)

(10)

尺寸l=1 100 mm代入(10)式后，再代入(8)式，經整理可得

(11)

簡化后零件在任意y處都具有相同的截面，因此三框結構件的中性軸高度位置和截面慣性矩分別為

(12)

(13)

三框結構件在不同位置y處，橫截面不同。因此可以計算任意位置y處截面繞X軸方向的力矩為

(14)

(14)式代入(6)式，二次積分即可得到撓度曲線。但由于在進行積分計算時無法求解出解析解，沿Y軸方向以50 mm為間距，依次從y=0 mm到y=1 100 mm進行力矩的計算，求解得到一系列力矩值，如表1所示。

由(6)式和(14)式可知，MX(y)是一個9次多項式。根據表1數據，采用多項式進行擬合即可求出彎矩MX(y)表達式，即

MX(y)=-3.996×10-35y9+2.31×10-21y8-1.016×10-17y7+1.79×10-14y6-1.6×10-11y5+7.665×10-9y4-1.91×10-8y3-1.0×10-7y2+0.013 6y-4 284.9. (15)表1 不同位置處的力矩Tab.1 Torques at different locations

由于高次項系數太小，不予考慮，進一步可得力矩表達式為

MX(y)=-1.0×10-7y2+0.013 6y-4 284.9.

(16)

(16)式代入(6)式，可以求得厚度方向上撓曲線方程為

UZ(y)=-1.66×10-6y2+αy+β，

(17)

式中：α、β為常數。

根據圖2(b)中零件首尾兩端的初始變形值條件，可求得零件在軋制方向上的撓曲線方程為

UZ(y)=-1.66×10-6y2+1.826×10-3y.

(18)

同樣，原始結構的每一段撓度曲線也可用同樣方法求得：

(19)

1.3 試驗驗證與分析

為了驗證解析法和有限元仿真方法的正確性，下面進行試驗驗證。試驗在意大利菲迪亞公司生產的高速銑床K211A上進行零件加工，在青島前哨精密機械有限責任公司生產的橋式三坐標測量機ADVANTAGE 15.30.10進行變形測量。整個試驗過程如下：

步驟1對毛坯進行無應力裝夾。

無應力裝夾是指該裝夾不對毛坯產生附加的應力，一般采用壓板沿四周頂住毛坯側面，如圖7所示。圖7中，Q、P分別為毛坯的正面和反面。

圖7 毛坯無應力裝夾Fig.7 Blank clamping without stress

步驟2對毛坯P、Q表面反復粗銑。

高速銑削毛坯Q面和P面，獲得零件在毛坯中的加工位置，具體工序如表2所示。

表2 毛坯表面的粗銑工藝路線Tab.2 Rough milling process of blank surface

步驟3對粗銑后的零件進行有應力裝夾。

有應力裝夾是指對零件施加夾緊力，將阻止材料去除過程中應力釋放造成的零件變形。故采用壓板在零件頂面壓緊的裝夾方式。

步驟4對零件進行精加工。

在意大利菲迪亞公司生產的K211A 3500×1500型銑床上加工成形為三框結構件，如圖8所示。加工3個框的內腔時，采用直徑為16 mm、刀頭數為4、刃長為55 mm的硬質合金銑刀；加工零件外側緣條時，則采用直徑為12 mm、刀頭數為4的硬質合金銑刀。主軸轉速均為15 000 r/min，進給量均為0.083 mm/r. 為了盡量避免銑削對加工變形的影響，在銑削外側緣條時，軸向切深僅為0.5 mm.

圖8 加工后的三框結構件Fig.8 Three-frame structure after machining

步驟5對三框結構件進行變形測量。

零件變形的測量在青島前哨精密機械有限責任公司生產的橋式三坐標測量機ADVANTAGE 15.30.10上進行。測量位置為零件腹板底面中線的數據，每隔50 mm進行一組數據的測量。

圖9所示為零件在毛坯內的加工位置分別為正向(零件開口朝向Q面)16.5 mm和反向(零件開口朝向P面)16.5 mm時的試驗值、計算值和仿真值的比較。

由圖9可見：加工位置為正向16.5 mm時，最大加工變形的試驗測量值為0.585 mm，毛坯存在0.6 mm上凸初始幾何誤差、無初始幾何誤差的計算值分別為0.527 mm、0.51 mm，其最大相對誤差分別為9.54%和12.8%，對應的有限元仿真值分別為0.535 mm、0.525 mm，最大相對誤差分別為9.4%、10.42%；當零件的加工位置為反向16.5 mm時，最大加工變形的試驗測量值為-0.323 2 mm，毛坯存在0.6 mm上凸初始幾何誤差、毛坯無初始幾何誤差的零件變形計算值分別為-0.320 6 mm、-0.301 8 mm，其相對誤差分別為14.5%和17.988%，相應的有限元仿真值分別為-0.325 9 mm、-0.301 mm，相對誤差分別為11.95%和17.6%. 由此可知，不論是有限元仿真還是計算值，毛坯存在初始幾何誤差時加工變形數據會更為接近試驗值。

圖9 零件變形試驗值、仿真值與計算值Fig.9 Experimental simulated and calculated values of part deformation

誤差存在的主要原因有以下兩個方面：其一，計算值中所使用的殘余應力值是人為試驗測量得到的，設備條件和試驗方法以及測量精度的限制不可避免地帶來誤差，因此計算彎矩時，所使用的殘余應力值是通過多項式擬合得到的曲線，會帶來一定的偏差；其二，對于鋁合金厚板而言殘余應力在同一厚度方向上的數值并不相同，在有限元分析過程中假定每一層的殘余應力值是中心平面的數值，使得仿真結果存在誤差。

2 零件變形演變機制

由第1節可知，零件變形的影響因素主要有兩個方面：一是毛坯的初始殘余應力分布；二是零件的結構形狀。

顯然，如果毛坯具有初始幾何誤差，則必將對零件變形產生影響。其本質就是毛坯初始幾何誤差導致零件受到的殘余應力發生了變化，如圖10所示。

圖10 初始幾何誤差對受力狀態的影響Fig.10 Influence of initial geometric error on stress state

這里，分析毛坯幾何誤差Δ取值為0～4 mm(方向為下凹)對零件變形的影響規律。事實上，零件在毛坯中的位置不同，也會影響零件的變形，故選取4個不同的位置進行分析，如圖11所示，分別為正向16.5 mm、正向9.0 mm、反向16.5 mm和反向9.0 mm.

圖11 零件在毛坯內位置示意圖Fig.11 Schematic diagram of the positions of parts in the blank

由于受到應力的腹板是零件變形的關鍵部位，表3列出了腹板所受到的應力。

利用第1節建立的有限元方法，可以獲得各個位置處毛坯初始幾何誤差對零件變形的影響，如圖12所示。

由圖12(a)和表3可知，當零件正向位于毛坯內16.5 mm時，腹板受到的等效應力隨著毛坯初始幾何誤差的增加呈現逐漸減小的趨勢，而且均為拉應力。因此，零件變形向上凸，且隨著初始幾何誤差的增加而減小。

由圖12(b)和表3可知，當零件反向位于毛坯內16.5 mm時，腹板受到的等效應力值隨著毛坯初始幾何誤差的增加呈現逐漸減小的趨勢，且為壓應力。因此，零件變形向上凸，且隨著初始幾何誤差的增加而減小。

表3 腹板受到的應力值Tab.3 Stress values of web plate

由圖12(c)和表3可知：當零件正向位于毛坯內9 mm，在初始幾何誤差為0～1.0 mm時，腹板受到的等效應力為壓應力，且呈現逐漸減小的趨勢，正向的零件變形為下凹，且變形逐漸減小；在初始幾何誤差為1.0～4.0 mm時，腹板受到的等效應力值為拉應力，且逐漸增加，則正向的零件變形也逐漸增加，變形向上凸。

由圖12(d)和表3可知：當零件反向位于毛坯內9 mm，腹板受到的等效應力均為壓應力，初始幾何誤差為0～1.0 mm時，其值呈增加趨勢，故反向的零件變形上凸，且隨著初始幾何誤差的增加而增加；初始幾何誤差為1.0～4.0 mm時，其值呈減小趨勢，故反向的零件變形上凸，且隨著初始幾何誤差的增加而減小。

圖12 不同初始幾何誤差下零件變形規律Fig.12 Deformation law of parts with different initial geometrir errors

由此可知，毛坯初始幾何誤差對零件變形的影響規律主要有以下3個方面：

1)零件在毛坯中的加工位置不同時，無論是否存在初始幾何誤差，零件變形也不相同；

2)零件在毛坯中的加工位置相同時，在相同的毛坯初始幾何誤差情況下，若零件結構出現變化，則零件變形曲線也會出現波動或拐點；

3)毛坯初始幾何誤差導致零件腹板受到的應力，若為拉應力，則變形向下凹；若為壓應力，則變形向上凸。若受到的應力隨毛坯初始幾何誤差的增大而增大，則變形也增大；若受到的應力隨毛坯初始幾何誤差的增大而減小，則變形也減小。

3 加工位置的優化方法

由變形規律可知，在毛坯初始幾何誤差一定時，零件在毛坯內的不同位置具有不同的變形量。因此，可通過優化加工位置來實現零件變形的控制[16-17]。

3.1 優化模型

在高速切削加工過程中，由于加工過程中產生的附加切削應力對零件變形的影響很小，優化切削工藝參數意義不大。因此，可以通過搜索零件在毛坯內的最佳加工位置，實現零件在加工過程中的變形控制。加工位置的優化模型建立如下：

findz，
minUX(y)，
s.t.zmin≤z≤zmax，

(20)

式中；zmin和zmax分別為零件在毛坯內的最低厚向位置和最高厚向位置。

3.2 步長遞減算法

為了求解出(16)式中z的解，可按照下列步長遞減算法來進行。

步驟1確定加工位置的初始近似值。

(21)

步驟2確定加工位置的第1個近似值。

(22)

根據加工位置的第1個近似值與初始近似值對應的變形方向，定義變形方向變化標識函數為

k1=λ1-λ0.

(23)

若k1=0，則表示變形方向未發生變化；若k1≠0，則表示變形方向發生了改變。

步驟3確定加工位置的下一個近似值。

根據相鄰兩個近似位置處變形標識的變化，確定當前近似位置的步長sn為

sn=(-1)kn-1η|kn-1|sn-1，

(24)

式中：η為0～1之間的任意數，也稱遞減系數，η∈[0，1]；kn-1=λn-1-λn-2，k0=0，n為加工位置近似值的次數，n≥1；λn-1、λn-2分別為對應于第n-1個近似值zn-1、第n-2個近似值zn-2的變形方向。

這樣，加工位置的下一個近似值可表示為

zn=zn-1+sn-1.

(25)

步驟4判斷當前近似值的范圍。

根據(16)式可知，加工位置的約束條件為[zmin，zmax]，則當前近似值zn應該滿足下列關系：

zmin≤zn≤zmax.

(26)

若znzmax，則zn=zmax+(-1)ηsn-1.

步驟5判斷終止條件。

由(20)式可知，當前步長可以進一步描述為

(27)

則當且僅當

|sn|≤ε

(28)

時迭代過程終止。

此時，最佳加工位置zo即為加工位置第n次的近似值，即

zo=zn.

(29)

4 應用與分析

以典型飛機整體三框結構件的高速銑削為例，說明其在無初始幾何誤差的毛坯和有初始幾何誤差的毛坯中加工變形的控制策略。

零件結構尺寸如圖3所示，所使用的毛坯結構尺寸和初始殘余應力曲線如圖4所示。毛坯和零件的厚度分別為60 mm和30 mm.若毛坯無初始幾何誤差，則零件在毛坯內的變動范圍是zmin=0 mm，zmax=30 mm.給定遞減系數η=0.5，閾值ε=0.000 1，初始步長為s=30，根據步長遞減算法，最佳加工位置的搜索過程如表4所示。由表4可知：當毛坯不存在初始幾何誤差時，步長迭代算法在n=27時，當前步長為sn=0.000 1.由于|sn|<ε，優化過程達到收斂。從理論上講，零件在毛坯內的最佳位置z=23.265 0 mm.因此，從工程實際考慮，最終選取最佳位置z=23.265 0 mm.表5所示為有初始幾何誤差時零件最佳位置搜索過程。

表4 無初始幾何誤差時零件最佳位置搜索過程Tab.4 Searching process of optimal positions of partswithout initial geometric error

由表5可知：毛坯存在初始幾何誤差Δ=0.6 mm時，零件在毛坯內的變動范圍是zmin=0.6 mm，zmax=29.4 mm. 同樣地，給定遞減系數η=0.5，閾值ε=0.000 1，初始步長s=30；步長迭代算法在n=30時，當前步長sn=-0.000 1，由于|sn|<ε，優化過程達到收斂。因此，由算法可知零件在毛坯內的最佳位置z=22.754 3 mm，故實際應用中選取最佳位置z=22.754 3 mm.

表5 有初始幾何誤差時零件最佳位置搜索過程Tab.5 Searching process of optimal positions ofparts with initial geometric error

通過有限元仿真得到毛坯存在0.6 mm下的凹初始幾何誤差，零件位于毛坯內正向23.265 0 mm時零件的變形為-0.144 mm；當毛坯不存在初始幾何誤差，零件位于毛坯內正向22.754 3 mm時，零件的最終變形為0.043 mm. 這兩組數據進一步表明毛坯的初始幾何誤差對零件加工變形具有重要影響。

5 結論

1)利用變形力學模型的計算值、實驗測量結果與有限元仿真值進行比較分析，結果顯示三者無論在變形趨勢還是變形幅值上均非常吻合。而且，考慮毛坯初始幾何誤差的零件變形比不考慮初始幾何誤差的結果更加接近實際結果。

2)零件位于毛坯內不同位置，零件加工變形大小不同。主要原因是，零件的變形與零件腹板位置所受的應力有關，若受到的應力隨初始幾何誤差增大而增大則變形也增大，反之則減小。

3)零件處于毛坯內一定位置時，零件變形曲線存在拐點或波動。原因有兩個：對于毛坯不同數值初始幾何誤差大小，拐點的存在是由于腹板位置受力的方向發生改變；對于毛坯同一數值初始幾何誤差大小，拐點存在與零件結構相關。

4)毛坯存在初始幾何誤差，會影響零件在毛坯內變形最小的位置。當毛坯存在0.6 mm的下凹初始幾何誤差時，其最佳加工位置為22.754 3 mm，對應的變形為-0.007 mm；在相同的位置時，若毛坯不存在初始幾何誤差，則零件的加工變形值為0.043 mm.