復合材料編織帶平面承力層結構有限元模型

劉 敏,襲建人,房光強,劉 鈺

(1.山東大學土建與水利學院,濟南250061；2.山東大學材料科學與工程學院,濟南250061；3.上海宇航系統工程研究所,上海201108)

1 引言

2016年4月美國太空探索技術公司將Bigelow公司的BEAM充氣太空艙成功與國際空間站完成對接,并順利完成太空艙充氣膨脹,這種可充氣膨脹展開太空艙重量輕、發射成本低,并可提供較大的居住空間和更好的空間防護[1]。可充氣太空艙通常采用多層復合材料制造,滿足不同的功能性和結構性要求,如圖1所示[1]。

可充氣太空艙結構的核心部分是其承力層,要滿足整個太空艙的結構性能要求,同時要求具有折疊和膨脹展開功能。為達到這一目的,承力層通常設計為由高強度纖維帶交叉編織在一起,如圖2所示[2]。

圖1 BEAM太空艙結構示意圖[1]Fig.1 Schematic diagram of BEAM space cabin structure[1]

圖2 承力層結構示意圖[2]Fig.2 Schematic diagram of load bearing layer[2]

纖維織物本身由纖維束按照編織模式由經紗和緯紗編織而成,在承受力的作用時,由于經紗和緯紗相互作用,而使其力學性能具有復雜的幾何非線性和接觸非線性特征,導致承力層的結構計算難度大大增加。將平面編織纖維織物等效為單層均質膜可以實現圓柱面的折疊和展開的模擬[3],而將平面編織纖維織物簡化為薄板或殼,按照正交各向異性的均質連續體進行計算是前期該領域研究的普遍做法[4-6],基于編織物的細觀結構特點建立的織物膜理論則進一步體現了纖維織物的力學特性[7]。從纖維束的細觀結構入手,Parsons等[8]通過將每束纖維視為連續體,考慮各纖維束間的接觸,并允許纖維束之間具有滑動特性,更真實地反映了平面編織纖維織物特性,更好地實現了力學性能的模擬[8]。然而將文獻[8]中的方法應用于實際的工程結構,將會由于計算規模太大而難于實際應用。

鑒于現有的計算能力,本文嘗試將平面的編織纖維帶視作連續的殼體,考慮各纖維帶之間的相互接觸作用和承受載荷作用下的相互之間的滑動,并根據編織纖維帶的力學特性,探討更適于解決編織纖維帶構成的承力層結構力學計算模擬方法。

2 理論與建模

2.1 平面殼單元的有限元理論

根據平面殼單元理論[9],空間殼單元在其局部坐標系中,每個節點有5個自由度,單元任一點的位移向量f寫為式(1)。

式中,u,v,w分別為x、y和z坐標軸方向的平動位移分量,θx和θy分別為繞x、y坐標軸的轉動位移分量。若采用4節點矩形的殼單元,則節點i的位移向量fi及其對應的位移分量為式(2)。

殼單元中位移分量ui和vi對應的薄膜效應的單元剛度矩陣為式(3)。

殼單元中位移分量wi、θi和θi對應板單元的用于表示彎曲效應的單元剛度矩陣為式(4)。

式(3)和(4)中上標m、b分別表示對應薄膜效應和彎曲效應的單元剛度矩陣。將式(3)和(4)分別擴充為20×20的矩陣,式(3)和(4)中的分塊矩陣和按照式(1)右端項的位移分量順序對應放置擴充矩陣中,在沒有對應分塊矩陣的位置填充O并將兩式相加得到式(5)。

按照式(5)的單元剛度矩陣得到的有限元方程,可以求解包含薄膜效應和彎曲效應的力學問題,而在式(5)中去掉彎曲效應對應的剛度矩陣,即將帶有字母b的分塊矩陣置為O,則可用于求解只考慮薄膜效應的力學問題。但在進行如圖2所示編織帶結構的力學分析中,因可充氣式太空艙的承力層,需要折疊和充氣膨脹展開,其彎曲剛度一定很小。對這樣的編織帶構造的結構進行力學模擬,傳統做法常只考慮膜應力的殼單元,從而低估了編織帶的抗彎曲性能。而采用式(5)的殼單元進行力學模擬,由于將編織帶的厚度完全計入其彎曲剛度,則明顯過高地計入了編織帶的彎曲剛度,與實際情況相差較大。因此本文提出基于薄膜效應殼單元附加一定彎曲剛度的混合單元來解決可充氣膨脹展開式太空艙的承力層結構的力學計算。

2.2 平面編織織物殼單元的有限元理論

為求解編織帶結構的力學問題,在節點自由度和單元剛度矩陣元素維數不變的條件下,將式(5)修正為式(6)。

上式中右邊的第1項為按照編織帶幾何尺寸和材料參數建立的只考慮薄膜效應的單元剛度矩陣,第2項則為修正項,使原有單元增加一定的彎曲剛度。根據殼單元理論,增加彎曲剛度的同時,也帶來一定的薄膜效應。因此為便于對比分析彎曲剛度的影響規律,彎曲剛度的修正是基于編織帶的厚度,再縮減一定百分比來實現,在此基礎上探討彎曲剛度的修正量對編織帶結構力學性能的影響。對于像編織帶這樣的結構,由于彎曲剛度很小,其給出厚度的縮減百分比較大。

為了減少計算量,原本由經紗和緯紗交織而成的具有離散結構特點的編織帶,將其宏觀性能簡化為用表觀寬度和厚度表示的平板。編織帶實際上具有一定的各向異性的宏觀特性,但其各向異性參數的測試十分困難,考慮到編織帶在實際承載中以承受單向拉伸為主,因此將編織帶簡化為各向同性材料,其彈性模量按照實際的測試數據選取。

2.3 平面編織承力層結構有限元模型

圖3為平面平紋編織帶結構,依此方式建立有限元模型,模型中橫豎方向各采用7條編織帶,編織帶的橫截面尺寸為25 mm寬,1 mm厚,長度為265 mm,材料的彈性模量為30 GPa。編織帶的4個外邊緣,全部平動位移約束,轉動自由度自由。編織帶的相互接觸面不考慮摩擦效應。承力層結構主要承受內部壓力,因此只在中心的3×3共9個上表面施加0.1 MPa的均布壓力。由于只考慮薄膜效應的有限元計算涉及到幾何非線性,因此按照軟件的要求,將單元網格劃分為三角形單元,即將單元中的4節點中的2個節點合并(圖3(b))。計算采用ANSYS軟件進行,彎曲剛度的調整取不同編織帶厚度的百分比進行,混合單元的形成是用同一位置的一組節點同時生成只考慮薄膜效應的單元和具有較小彎曲剛度的殼單元。有限元模型單元參數的計算分組見表1。

3 結果與討論

3.1 計算結果

3.1.1 薄膜單元計算結果

圖4為只采用薄膜單元的垂直位移和應力分布圖。圖5為中心編織帶的x方向和z方向的應力分布圖。

從圖4(a)中可以看出,盡管編織帶結構是交叉編織而成,但垂直方向的變形就像是一張帶有孔的平板變形一樣,具有很好的整體性。最大的垂直位移為-5.999 mm,位于編織帶的中心。從圖4(b)編織帶中x方向的應力分布可以看出薄膜應力呈現出較為均勻的分布,最大數值為36.31 MPa,從中心向邊緣方向,應力數值逐漸減小。圖5結果進一步顯示,x方向的應力沿x方向和z方向,應力數值均呈現小幅的變化,最小33.73 MPa,最大36.31 MPa。編織帶中z方向的應力數值均很低,在±0.5 MPa之內,這說明編織帶相互之間的接觸作用并沒有改變其單向應力狀態的性質,這樣的結果與傳統編織帶單向承載能力的設計一致。

圖3 平面平紋編織帶結構及其有限元網格圖Fig.3 The geometry and FEM meshing diagram ofwoven fabrics structure

表1 單元厚度情況分組Table 1 Thickness of the membrane and hybrid elements

圖4 平面編織結構在0.1 MPa壓力作用下的位移和應力分布圖(薄膜單元)Fig.4 Distribution contour of displacement and stress in the plain woven fabrics structure at 0.1 MPa(only membrane element)

3.1.2 混合單元計算結果

圖6為采用薄膜單元加0.1 mm厚度的殼單元,即采用混合單元模型得到的x方向應力分布圖,圖中看出編織帶的正面和背面的應力數值不同,最大應力值在背面,數值為50.23 MPa,而同一點的正面,數值為45.63 MPa。圖7(a)為中心編織帶的x方向應力分布,最大應力值是位于中心部位的臨近區,而圖5(a)應力最大值則在中心部位。

圖7(b)為中心編織帶z方向應力分布,其應力數值相對于只采用薄膜效應單元的計算值要大得多,中心部位的 2個邊緣數值最高,達到20.63 MPa。

不同彎曲剛度的混合單元的應力分布與圖6所示的規律相似,只是在數值有一定的差異(表2)。由于對稱關系,另一個方向的編織帶,亦即沿z方向布置的編織帶的z方向的應力分布與圖6和圖7中的x方向應力分布相同。

圖5 中心部位編織帶在0.1 MPa壓力作用下應力分布圖(薄膜單元)Fig.5 Distribution contour of stress in the center of plain woven fabrics at 0.1 MPa(only membrane element)

圖6 平面編織結構在0.1 MPa壓力作用下的應力分布圖(混合單元)Fig.6 Distribution contour of stress in the plain woven fabrics structure at 0.1 MPa(hybrid element)

為分析殼單元的彎曲剛度對混合單元計算結果的影響,選取殼單元的厚度值從0.01～0.3 mm之間變化,即殼單元厚度的變化是編織帶厚度的1%～30%。從表2可以看出,殼單元的彎曲剛度越大,其垂直位移越小,從薄膜效應單元的-5.999 mm,減少到混合單元的-5.209 mm,垂直位移數值的絕對值減少了13.2%。相對于垂直位移影響來說,混合單元的采用對x方向和z方向應力的影響較大。混合單元中,當彎曲剛度較小時,同等條件下產生數值為52.62 MPa的x方向應力增加幅度高達45%。同樣z方向應力盡管數值在20 MPa左右,但增加幅度較大。因此只采用薄膜效應單元會高估結構的靜態承載能力。

表2中的數據也表明,殼單元的厚度值從0.01～0.3 mm變化,以0.01 mm的計算結果作為參考,則垂直位移的變化幅度的絕對值為8.412%,x方向應力變化幅度的絕對值為11.78%,z方向應力變化幅度的絕對值為15.92%。

3.2 討論分析

1)對于像編織帶這樣的平面織物來講,傳統上為計算簡化常忽略其彎曲剛度,利用只考慮薄膜效應的薄膜單元進行結構的計算分析,這樣做是可表征其主要承載特性。從計算結果中可以發現(圖5(a)),平面織物結構承受局部壓力載荷會造成編織帶的局部薄膜應力增加,這一點與將編織結構進行均勻簡化的計算結果不同。同時可以注意到,局部的承壓會導致織物在受壓處產生較大的局部下凹,且這種局部下凹不僅沿編織帶長度方向存在,同時沿編織帶平面內的橫向也有一定程度的存在,從而造成其薄膜應力不均勻分布,但應力狀態基本處于單向應力狀態。

2)編織帶由于采用經紗和緯紗交叉編織在一起,具有較高的單向承載拉伸的能力,同時也具有較小的彎曲剛度,以便折疊彎曲之用。這個較小的彎曲剛度盡管對編織帶結構的剛度影響較小(13.2%以內),但對工作應力會帶來較大的影響。如果殼單元的厚度只取原薄膜單元厚度的10%,即殼單元的彎曲剛度只取原厚度剛度的千分之一,則沿編織帶縱向的最大應力計算值會增加45%,而編織帶橫向的最大應力計算值會從很小數值增加為20 MPa左右(表2)。同時應力分布規律也與只采用薄膜效應單元的完全不同,縱向應力數值較高的區域并不位于中心編織帶的中心區域,而是在相互垂直編織帶交叉重疊區域,且最大縱向應力位于長方形的交叉重疊區的4個角點處,而中心區域的應力值低于4個角點處的應力值(圖7(a))。橫向應力較大的區域仍位于中心編織帶的中心區域,但最大的橫向應力出現在相互垂直編織帶交叉重疊區域的2個邊緣上(圖7(b))。

圖7 中心部位編織帶在0.1 MPa壓力作用下應力分布圖(混合單元)Fig.7 Distribution contour of stress in the center of plain woven fabrics structure at 0.1 MPa(hybrid element)

表2 不同單元計算結果比較Table 2 Comparison of calculation results of different elements

3)對比圖5和圖7,可以看出編織帶的整體變形規律基本相同,即壓力作用區出現明顯的下凹,中心區域較為平坦,垂直位移數值較大,但其局部曲率較小；壓力作用區的邊緣垂直位移數值較小,但其局部曲率較大,這也是具有彎曲剛度的編織帶其應力的高水平點位于中心編織帶重疊區的邊緣的原因,而只采用薄膜效應單元則無法模擬這種力學特性。

以上分析表明,對于編織帶構成的平面織物結構,采用不同的單元類型對結構靜態應力變形模擬計算有較大的影響,采用混合單元更能反映編織帶的受力特點。從混合單元的計算結果可以看出(表2),垂直位移和應力計算值對殼單元厚度值的敏感性較低,其變化幅度在16%以內。在編織帶的彎曲剛度無法準確獲知的條件下,利用基于混合單元的有限元模型進行編織帶結構的靜態力學分析可更準確地反映出編織帶固有的力學特性。由于殼單元厚度值不僅增加彎曲剛度,同時增加薄膜剛度,會對編織帶中的薄膜應力計算帶來較大的影響,所以取相對較小的殼單元的厚度值更接近編織帶的實際承載狀態。

因此對于像編織帶這樣的結構件而言,目前用于結構分析的彎曲剛度的實驗方法測定仍屬空白,在此條件下,選取具有較小彎曲剛度的混合單元有限元模型進行其結構力學計算分析是值得推薦的可行方法,可以更準確地模擬編織帶結構的力學響應。

4 結論

1)采用薄膜效應單元并考慮編織帶之間相互接觸的有限元模型,在非均勻承載條件下,每根編織帶內會呈現不均勻的縱向應力分布,更準確體現膜的受力特征。但由于薄膜效應單元不能精確地反映編織帶的彎曲剛度,導致編織帶基本呈現單向應力狀態,無法真實地反映局部接觸區域的應力狀態,從而高估編織帶的靜態承載能力。

2)采用以薄膜效應單元為基礎,疊加具有較小彎曲剛度的殼單元形成混合單元來建立有限元模型,更能準確地反映編織帶的實際承載效應,更真實地評定其靜態受力狀況,在目前缺少編織帶彎曲剛度實驗方法的條件下,是一種更好的評估編織帶織物結構的靜態力學性能的有限元模型。

3)基于混合單元的有限元模型,對殼單元彎曲剛度取值的敏感性較低,可采用較小數值彎曲剛度的混合單元有限元模型進行編織帶結構力學分析,以更準確地模擬編織帶結構的靜態力學響應。