廉鑫, 張衛東, 毛玉明, 于哲峰*
(1.上海交通大學航空航天學院, 上海 200240; 2.上海宇航系統工程研究所, 上海 201109)
大型薄壁加筋結構具有結構重量輕、強度高等優點,廣泛應用于飛行器結構等領域,如火箭貯箱結構等[1-6]。其加筋特征可以顯著地提升承力薄殼或壁板的抗彎剛度,避免結構過早的發生屈曲失穩變形,從而有效增強結構在軸壓載荷下的承載能力[7-11]。加強筋的形狀、尺寸與布局直接影響著結構的質量與性能。因此,如何實現板殼結構的加強筋布局自動化設計一直是國內外學者與工業界高度重視的問題。
在過去的30年中,國內外學者針對相關結構優化問題做了許多工作。倪楊等[12]基于多學科優化軟件,利用粒子群優化算法對平板加筋結構的參數進行優化,使得結構在單軸壓載作用下的穩定性能與極限承載能力得到了明顯提升。Zhao等[13]針對貯箱的加筋筒殼,開展了多級三角形加筋構型優化,得到了高承載力的設計方案。相比于單級加筋筒殼,多級加筋殼在同樣的承載力約束下具有更大的輕量化設計潛力。Alinia[14]基于ANSYS有限元分析軟件研究了加筋板在剪切荷載作用下的優化設計,通過確定整體屈曲向局部屈曲轉化的臨界點來設計加強筋的最優幾何特性。Lam等[15]將板的厚度作為設計變量對加筋肋的布局進行了優化,在優化后板的厚度較薄的地方添加加強肋。Chung等[16]基于變密度方法的拓撲優化,研究了加筋肋的最優形狀和位置,使得結構的材料利用率得到了很大的提高。
從以上文獻中可見,現有的加筋結構優化方法多采用商業有限元軟件,直接求解大型有限元模型,往往需要消耗較多時間,對于含受局部集中載荷的加筋壁板優化問題也缺少研究。現提出一種基于解析公式的加筋結構優化方法,在誤差允許的范圍內,能夠較快地求解出優化后的結構參數。然后針對捆綁接頭與加筋壁板之間協調設計問題進行相關的參數研究,以期為工程設計提供參考。
結構優化指在滿足給定的幾何約束、材料約束和狀態約束的所有可能結構中,尋求使得給定目標最優的設計[17]。約束條件和目標函數共同組成了優化問題的數學模型,進行優化設計的數學模型為
(1)
式(1)中:f(x)為優化問題的目標函數;gj(x)為優化問題的第j個不等式約束;hk(x)為優化問題的第k個等式約束。
對于這個數學優化問題,可以使用拉格朗日乘子法進行求解。首先將式(1)中的不等式約束全部轉化為等式約束,即對模型中的每一個不等式約束引入一個松弛變量sj,即
(2)
轉換后的數學模型只包含等式約束,可構造出拉格朗日函數為

(3)
式(3)中:λk為與等式約束hk(x)有關的拉格朗日乘子;μj為與不等式約束gj(x)有關的拉格朗日乘子,對式(3)求一階偏導數可得

=0
(4)
引入Karush-Kuhn-Tucker(KKT)條件[18],該條件為約束優化問題的一階必要性最優條件。
如果規范點x*為約束優化問題[式(1)]的最優解,則存在非負常數μj使得在x*處滿足下列條件:①x*滿足式(1)中的所有約束;②式(4)成立;③所有在x*處取不等號的約束相關的拉格朗日乘子為0。
證明優化問題解的存在性后,需要通過具體的算法對優化問題進行求解,得到最終的優化解。
采用MATLAB軟件編程實現載荷擴散優化,采用的優化函數是fmincon函數。在默認情況下,fmincon函數主要采用直接法求解優化解。在直接法中,求式(3)中拉格朗日函數的散度可以得到Hessian函數,即


(5)
為了求解滿足KKT條件的優化解,構造如下方程[19]。
(6)
式(6)中:Jg為約束函數g的雅可比矩陣;Jh為約束函數h的雅可比矩陣;S為松弛變量組成的對角陣;λ為與約束g相關的拉格朗日乘子向量;Λ為拉格朗日乘子組成的對角陣;y為與h相關的拉格朗日乘子向量;e為與g大小相同的向量。
為了求解(Δx,Δs)的方程,該算法對矩陣進行LDL分解[20]。這是計算成本最高的一步。這種分解的一個結果是確定投影的H矩陣是否正定;如果不是,該算法使用共軛梯度方法進行求解。
載荷擴散結構優化目標是在保證結構穩定的情況下使其質量最小,采用解析公式建立優化目標函數。
平板加筋結構的幾何模型如圖1所示,該結構由m根筋條與n塊蒙皮組成,所有部件材料相同,密度設為ρ。整個平板加筋結構可按照載荷的施加狀態分成兩個部分。一類是沒有接頭的部分,只受到遠端施加的線載荷Nx1,另一類是含有接頭的部分,受到包含遠端線載荷Nx2和由接頭集中力F2擴散的線載荷Nx2,兩種部件的橫截面參數如圖1所示。

Lx為加筋板的高度,mm; Lyn為加筋板的寬度,mm; Lyf為接頭部分的寬度,mm; Nx1為遠端施加的線載荷,N/mm; F2為接頭傳遞的集中力,N; Nx2為集中力等效的線載荷,N/mm; tsk為蒙皮厚度,mm; Ly1為第一類結構部件的蒙皮寬度,mm; ly1為第一類結構部件的桁條寬度,mm; lz1為第一類結構部件的桁條高度,mm; ty1與tz1為第一類結構部件的桁條厚度,mm; Ly2為第二類結構部件的蒙皮寬度,mm; ly2為第二類結構部件的桁條寬度,mm; lz2為第二類結構部件的桁條高度,mm; ty2與tz2為第二類結構部件的桁條厚,mm
當結構承受壓縮線載荷時,失穩形式可能包括筋條局部失穩、蒙皮局部失穩、周期結構部分的失穩與結構整體的寬柱失穩。
(1)筋條局部失穩。在筋條局部失穩的情況下,加強筋可視作為三邊簡支、一邊自由的長板,此時對于兩類部件,筋條的臨界失穩應力計算公式為
(7)
式(7)中:E為材料彈性模量;ν為材料泊松比。
(2)蒙皮局部失穩。當蒙皮發生局部失穩時,可將其視作4邊簡支的板,臨界失穩應力為
(8)
(3)周期結構部分的失穩。當蒙皮與筋條組成的周期結構發生失穩時,失穩形態類似柱的失穩,可求出周期結構的總體慣性矩Ic、橫截面面積Ac與非加載邊的長度Lx,按兩端為簡支的假設,此時失穩應力為
(9)
(4)結構整體的寬柱失穩。當結構發生整體失穩時,整個結構被看作一個寬柱發生失穩,此時需要求出整個結構的總體慣性矩Ict。非加載邊的長度與上文一致,均為Lx,此時臨界失穩載荷可通過式(10)計算得到,即
(10)
對于結構中第i根筋條,失穩應力應該取式(7)中局部失穩應力與式(9)中周期結構失穩應力的最小值。
[σcr]i=min{[σcr,st]i,[σcr,c]i}。
對于結構中第j塊蒙皮,失穩應力應該取式(8)中局部失穩應力與式(9)中周期結構失穩應力的小值。
[σcr]j=min{[σcr,sk]j,[σcr,c]j}。
平板加筋結構的失穩應力集合可表示為
{[σcr]}={[σcr]i,[σcr]j},i=1,2,…,m;j=1,2,…,n。
對于第一類結構部件,可求得筋條的橫截面面積Ap1,該類部件的寬度為Ly1,蒙皮的厚度為tsk,可根據式(11)求得部件的等效厚度為
(11)
部件的應力計算公式為
(12)
對于第二類結構部件,用相同的方法求得筋條橫截面積Ap2,部件加載端的寬度Ly2,蒙皮厚度仍為tsk,部件的等效厚度為
(13)
部件的應力計算公式為
(14)
根據式(12)與式(14)可求得每個部件的應力,對于模型中的所有部件,計算局部失穩偏差相關的變量V1為

(15)
整個平板加筋結構的加載邊寬度為Lyn,接頭載荷擴散的區域寬度為Lyf,總載荷為
F=Nx1Lyn+Nx2Lyf
(16)
對于整個平板加筋模型,計算總體失穩偏差相關的變量V2為

(17)
求解加筋結構的總質量,對于結構中的第i根筋條,其橫截面積為
Ai=ly,itz,i+(lz,i-tz,i)ty,i
(18)
所以第i根筋條的質量可表示為
mst,i=ρAiLx
(19)
對于結構中的第j塊蒙皮,其質量為
msk,j=ρLy2,jtskLx
(20)
整個結構的總質量為

(21)
優化問題的目標函數為
G=(V1+wV2)M
(22)
式(22)中:w為部件局部失穩偏差和總體失穩偏差之間的權重比例。采用優化算法,求解G的最小值,即得到各個筋條橫截面的幾何參數、筋條間距和蒙皮的厚度。
以一個平板加筋結構的優化為算例,整個結構部件由鋁合金材料組成,彈性模量E=73 GPa,泊松比v=0.3,密度取ρ=2 700 kg/m3。在遠端施加均勻線載荷Nx1=15 kN/m。結構加載邊的長度Lyn=0.5 m,非加載邊的長度Lyn=1 m。結構中心位置有均勻線載荷Nx2=10 kN/m,其初始寬度LyfS=0.25 m,用來模擬由接頭處擴散的集中載荷。所以結構承擔的總載荷F=10 000 N。如圖2所示,9根筋條將整個結構劃分為8個部分。這8個部分又可分為兩個區域,一個是沒有接頭的區域,由1、2、7、8部分組成;另一個是包含接頭的區域,由3、4、5、6部分組成。其中筋條3和7屬于包含接頭的區域。

1~9代表桁條編號;①~⑧代表蒙皮編號

9根筋條的初始幾何參數在表1中列出,蒙皮的初始厚度為tsk=1 mm。整個結構的初始質量為m0=3.513 kg。考慮到裝配的需要,所有筋條橫截面上ly的長度固定為30 mm。優化過程所有的設計參數為:接頭寬度Lyf、9根筋條各自的ty、tz、lz和蒙皮厚度tsk。邊界條件為:施加線載荷Nx1的邊界y與z方向位移被約束,其對邊的x、y、z方向位移被約束,其兩個鄰邊的y、z方向位移被約束。表2列出

表1 結構初始幾何參數

表2 初始結構各部分的失穩應力和估算的應力
了初始結構各個部分的失穩應力和估算的應力,可見估算應力小于兩類失穩應力。該參數下結構的總體失穩臨界載荷Fcr=119 668 N。
對該算例進行優化,優化后接頭寬度Lyf=0.230 m;所有筋條的幾何參數在表3中列出;蒙皮的厚度tsk=0.350 mm,結構總質量為m1=2.537 kg。表4列出了優化后各個部分的失穩應力和估算的應力,可見部件應力與失穩應力中的小值更為接近,而是結構失穩形態是總體的寬柱失穩載荷,該參數下結構的總體失穩臨界載荷Fcr=10 282 N,與結構總載荷F=10 000 N很接近。

表3 優化后結構的幾何參數

表4 優化后各部分的失穩應力和估算的應力
使用有限元來分析優化前后結構的穩定性,在有限元建模軟件ABAQUS中按照優化前后的幾何參數建立有限元模型,考慮到結構的對稱性,只取結構的一半進行建模分析。采用S4R殼單元建立網格模型,網格種子的全局尺寸參數為3.0,網格單元個數為22 400,網格結點總數為22 715。實際的火箭艙體結構承受的載荷來自遠端的集中力F1和接頭的集中力F2,而這會在艙體結構處擴散為均勻的線載荷Nx1與Nx2,為了在有限元模型中更好地模擬這個過程,在結構外部的一個節點上施加等效集中力F1和F2,并將此節點與艙體結構擴散端的節點以剛體的接觸形式相連接,這樣能夠保證承受線載荷的端部的所有節點能夠有相同的y方向位移,接觸方式的施加形式如圖3所示,整個模型的邊界條件為四邊簡支,如圖4所示。

圖3 模型載荷施加方式

圖4 對稱模型的邊界條件
分別對優化前后的結構有限元模型進行線性屈曲分析后,所得結果如圖5和圖6所示。結構第一階失穩特征值由最初的3.261降低至1.418,失穩形式是中心部分的寬柱失穩,而有限元模型一邊是簡支的,所以不完全符合柱模型,這也是造成特征值較大的原因。

圖5 初始加筋壁板結構的屈曲分析結果

圖6 優化后的平板加筋結構的屈曲分析結果
從以上分析結果可知,初始結構的一階失穩特征值較高,結構重量有較大的冗余;而優化后的結構在保證一階失穩特征值大于1的前提下減少了27.78%的質量,使得材料利用率得到了提升。


圖7 優化結構質量與接頭載荷所占比例的關系曲線
處理圖7中的數據可以發現,優化結構的質量方差數值為0.12,可見接頭載荷所占比例變化時,優化結構質量變化并不明顯。這個結果表明在總載荷一定的情況下,可以較為自由地設計接頭載荷所占比例而不必擔心結構質量會因此大幅波動。
在核心艙段的參數化設計中,接頭寬度也是經常關注的重要參數,在優化設計中,都是將其作為優化參數,使用優化程序來獲得寬度值。而不同寬度接頭的質量也有所不同,確定優化結果時應該考慮不同寬度接頭的質量,所以進一步研究接頭寬度對優化結構質量的影響,計算其在一定范圍內變化時優化得到的結構質量。取接頭寬度的變化范圍為0.175~0.325 m,求出各個接頭寬度下優化結構的質量,二者的關系曲線如圖8所示。

圖8 優化結構質量與接頭寬度的關系
從圖8中可以發現,隨著接頭寬度的增加,優化結構質量呈現出下降的趨勢,這符合越均勻的結構性能更優的特點,但如上所述,捆綁接頭寬度增加會造成自身重量的增加,所以應同時考慮接頭重量的優化結果來確定接頭寬度。
研究受集中載荷的加筋壁板結構優化方法。基于板和柱的穩定性公式建立結構模型,提出了同時考慮穩定性和結構質量的目標函數,優化算例和有限元模擬結果驗證了該方法的有效性。
針對實際應用中主要關注的集中載荷所占比例和捆綁接頭寬度兩個參數的影響進行了參數分析。結果表明,在考慮結構輕量化的設計中,不必過于強調接頭載荷占比,可以較為自由地進行設計;接頭寬度的增加會使得優化結構的質量呈現出下降的趨勢,在結構設計時可以在綜合考慮接頭質量情況下,通過增加接頭寬度來進行減重。
在根據穩定性解析公式進行建模時,沒有考慮不同部件之間的應力傳遞,這是導致應力計算誤差和優化后結構特征值仍大于1的主要原因。在未來的研究工作中,將根據結構載荷擴散的相關公式對模型進行修正,進一步提高計算結果的準確度。