負高斯曲率曲面薄壁管吸能特性研究*

馬梓鴻，張慧樂，孫澤玉，陳慧敏，岳曉麗

（1. 東華大學機械工程學院，上海 201620；2. 東華大學上海市輕質(zhì)結(jié)構(gòu)復合材料重點實驗室，上海, 201620；3. 電子科技大學長三角研究院（湖州），浙江 湖州 313000）

隨著我國汽車產(chǎn)銷量的逐年增長，交通安全問題已經(jīng)成為我國汽車工業(yè)發(fā)展所面臨的重要問題之一。薄壁吸能結(jié)構(gòu)在受到?jīng)_擊載荷作用時，能夠吸收耗散碰撞動能，因此薄壁結(jié)構(gòu)在耐撞性領(lǐng)域得到了廣泛應用。

目前，薄壁結(jié)構(gòu)表面形狀以零高斯曲率曲面為主，但自然界中還存在負高斯和正高斯曲率曲面。數(shù)學家高斯將其總結(jié)為“高斯絕妙定理”[1]：對于任意曲面，當該曲面內(nèi)任意點P的兩條主曲線曲率k的乘積為負時，則該曲面為負高斯曲率曲面，同理可得零高斯及正高斯曲率曲面。負高斯曲率曲面在建筑領(lǐng)域受關(guān)注較多[2-3]，但在薄壁結(jié)構(gòu)耐撞性領(lǐng)域鮮有研究。

在針對薄壁結(jié)構(gòu)耐撞性和吸能特性的研究中，研究者們較多關(guān)注橫截面為方形[4-5]、六邊形[6-7]、八邊形[8-11]和圓形[12-13]的薄壁管。張宗華等[14]對上述4 類薄壁管進行動態(tài)軸向沖擊模擬后發(fā)現(xiàn)從方形到六邊形，結(jié)構(gòu)壓潰力和吸能量顯著提高，但八邊形后趨于穩(wěn)定，圓形截面優(yōu)勢并不顯著。Guillow 等[15]研究了不同長徑比和徑厚比下薄壁圓管的壓潰響應，發(fā)現(xiàn)其變形模式主要包括軸對稱變形、金剛石變形、混合變形（軸對稱變形和金剛石變形共存）以及歐拉屈曲等4 類。為進一步探究薄壁圓管吸能特性，Zarei 等[16]、Hou 等[17]和Marzbanrad 等[18]基于響應面法和遺傳算法對薄壁圓管進行了多目標優(yōu)化設計，提高了吸能量的同時降低了結(jié)構(gòu)質(zhì)量。但上述研究對象表面形狀按照高斯曲率類型分類仍屬零高斯曲率曲面，對于負高斯曲率情況下的薄壁管吸能特性研究較少。

為設計出吸能特性優(yōu)良的薄壁管，本文中將負高斯曲率曲面形態(tài)引入薄壁管設計中，提出一種新型負高斯曲率曲面薄壁圓管結(jié)構(gòu)（negative Gaussian curvature surface circular tube，NGC-C）。利用經(jīng)驗證的有限元分析方法對其進行軸向動態(tài)沖擊模擬，提取各項性能指標。借助復雜比例評估（complex proportion assessment，COPRAS）法將其與傳統(tǒng)薄壁管進行性能比對。采用拉丁超立方采樣法構(gòu)建樣本空間并獲取各樣本點對應的性能響應值，進一步推導出其代理模型。基于此代理模型，借助改進非支配排序遺傳算法（non-dominated sorting genetic algorithm，NSGA-Ⅱ）進行多目標優(yōu)化設計。以期為薄壁吸能結(jié)構(gòu)的設計提供新思路。

1 有限元模型及驗證

1.1 有限元模型建立

薄壁管設計參數(shù)主要有：口徑R、厚度T、曲率參數(shù)G和長度L。固定長度L=180 mm，如圖1(a)所示。當曲率參數(shù)G分別取負數(shù)、零、正數(shù)時，薄壁管形狀分別呈負高斯曲率曲面、零高斯曲率曲面以及正高斯曲率曲面，傳統(tǒng)薄壁結(jié)構(gòu)一般為零高斯曲面。

利用非線性有限元軟件LS-DYNA 進行動態(tài)軸向壓潰數(shù)值模擬分析，并利用Ls-Prepost 4.8 軟件進行前后處理。薄壁管單元類型為Belytschko-Tsay 殼單元，采用全積分方法，厚度方向積分點設置為5。為防止管件自身在變形過程中出現(xiàn)穿透現(xiàn)象，利用Automatic single surface 單面自動接觸類型來設置管件自身接觸，利用Automatic surface to surface 定義管件與沖擊墻之間的接觸[18]，摩擦因數(shù)設置為0.25[16]。圖1(b)展示了3 種網(wǎng)格尺寸（4、3 和2 mm）的載荷-位移曲線，當網(wǎng)格尺寸為3 mm 時，其載荷波動趨勢與網(wǎng)格為2 mm 的吻合度較高，且初始峰值力、平均壓潰力以及總吸能等指標相差均在3%以內(nèi)；而網(wǎng)格尺寸為4 mm 時曲線整體偏低，尤其第2 個峰值力相較于2 mm 網(wǎng)格下降了16%。當網(wǎng)格進一步細化會大幅增加計算時間，3 mm 網(wǎng)格計算結(jié)果已經(jīng)能夠滿足分析需求，因此采用3 mm 網(wǎng)格。設置薄壁管底部邊緣為固定支撐。上端有300 kg 剛性沖擊墻，以6.0 m/s 初速度沖擊薄壁管，并設置其只存在軸向位移自由度。薄壁管材料選用A6060-T5 鋁合金[12]，使用LS-DYNA 中的MAT103 號材料模型，材料參數(shù)見表1。由此建立G=0 時的有限元模型（finite element model, FEM），如圖1(c)所示。

圖1 有限元模型以及網(wǎng)格選擇Fig. 1 Finite element model and element size selection

表1 材料參數(shù)Table 1 Material parameters

1.2 有限元模型驗證

如圖2，有限元模型求解后得出的沙漏能（hourglass energy）對總能量（total energy）的占比小于5%，求解結(jié)果收斂性良好。為驗證有限元模型的可靠性，借鑒文獻[18]的驗證方法，利用相關(guān)文獻實驗結(jié)果驗證有限元模型可靠性。首先分別構(gòu)建出文獻[6]中討論的方形、六邊形、圓形和文獻[11]中探究的八邊形截面薄壁管，這4 個模型受動態(tài)軸向壓潰后實驗與數(shù)值模擬的變形模式結(jié)果對比情況可見圖3，各模型實驗與數(shù)值模擬的變形模式相似度較高。進一步地，構(gòu)建出文獻[16]中具有代表性的S-1、S-3、S-5 和S-6 薄壁圓管模型，其參數(shù)見表2。由于文獻中只給出了S-1 模型的變形模式和力-位移曲線，因此將S-1 實驗與數(shù)值模擬結(jié)果的對比情況分別繪制成圖4 和圖5，將各模型實驗與數(shù)值模擬結(jié)果列于表3。比較后發(fā)現(xiàn)：S-1圓管受壓潰后呈現(xiàn)出混合變形模式，與實驗結(jié)果相似性較強；壓潰過程的力-位移曲線與實驗所得曲線變化趨勢吻合度較好，數(shù)值上存在差異，但偏差在可接受范圍內(nèi)。綜上，所建立的有限元模型能夠較好地反映結(jié)構(gòu)變形模式受壓潰時的載荷以及吸能情況，故可靠性較高。

圖2 能量變化情況Fig. 2 Energy change

圖3 薄壁管變形模式的實驗與模擬對比情況Fig. 3 Comparison of experiment and simulation of four kinds of thin-walled tube deformation modes

表3 實驗與數(shù)值模擬結(jié)果對比Table 3 Comparison of experimental and simulation numerical results

圖4 S-1 模型實驗與數(shù)值模擬變形模式對比Fig. 4 Comparison of experiment and simulation deformation modes of S-1 model

圖5 S-1 模型實驗與數(shù)值模擬力-位移曲線對比Fig. 5 Comparison of force-displacement curves between experiment and simulation of S-1 model

表2 有限元驗證模型尺寸參數(shù)[16]Table 2 Finite element verification model size parameters[16]

2 研究結(jié)果及討論

吸能結(jié)構(gòu)一般依靠其材料的變形來對外部沖擊能進行吸收耗散，其性能一般由總吸能E、比吸能 δ 、初始峰值壓潰力F0、平均壓潰力F、壓潰力效率 η 等指標進行衡量，同時對于類似落錘沖擊的動態(tài)工況，還應考慮有效壓潰長度l。著重關(guān)注比吸能以及有效壓潰長度，并且期望薄壁管兼具高比吸能和較低的有效壓潰長度，暨能夠在一定質(zhì)量下吸收更多沖擊能，同時具備較強的軸向抗沖擊能力，能夠在更短的進程中使外部沖擊物減速制動。各指標示意可見圖6，計算式為：

圖6 沖擊過程中力-位移曲線Fig. 6 Force-displacement curve during impact

式中：F為壓潰過程結(jié)構(gòu)受力，m為結(jié)構(gòu)質(zhì)量。

2.1 薄壁管變形模式分析

構(gòu)建等質(zhì)量條件下正高斯曲率曲面（positive Gaussian curvature, PGC），零高斯曲率曲面（zero Gaussian curvature, ZGC）和負高斯曲率曲面（negative Gaussian curvature, NGC）等3 類曲率情況下的正方形（square, S）、六邊形（hexagon, H）、八邊形（octagon, O）和圓形（circular, C）橫截面薄壁管，共計12 個。由于材料密度相同，使模型體積相同即可保證質(zhì)量相同。

將其中某一模型沿縱向剖切，所得剖面的一半會產(chǎn)生兩段曲線，將兩段曲線進行數(shù)學表達及處理，之后將相應結(jié)構(gòu)橫截面面積S沿任意一條曲線方向進行積分即可求得其體積V。以NGC-S 模型為例，如圖7 所示：沿縱向剖切得剖切面(紅線所圍區(qū)域)，其上兩曲線分別為y1,y2，則NGC-S 體積VNGC-S等于其橫截面積SS沿任一曲線積分為：

圖7 NGC-S 剖切面示意圖Fig. 7 Schematic diagram of section plane of NGC-S

因此，求得各薄壁管橫截面面積與邊長a或半徑r之間關(guān)系式關(guān)系式，即可進一步推得其體積：

為確保結(jié)構(gòu)變形從初始受沖擊端發(fā)生，進而獲得完整變形形態(tài)，在薄壁管頂端添加長1.5 mm、與水平面呈37°角的斜坡觸發(fā)結(jié)構(gòu)[19]，見圖8。

圖8 觸發(fā)結(jié)構(gòu)Fig. 8 Trigger structure

各模型沖擊壓潰后變形模式如圖9 所示。首先從截面類型角度分析可知：正方形截面下的各結(jié)構(gòu)呈軸對稱模式，且產(chǎn)生的褶皺間距較大；六邊形與八邊形結(jié)構(gòu)則表現(xiàn)出軸對稱模式和混合模式；對于圓截面，只有NGC-C 結(jié)構(gòu)表現(xiàn)為軸對稱模式，其余均為混合模式。其次從高斯曲率類型角度分析可知：負高斯曲率情況下各結(jié)構(gòu)均呈現(xiàn)軸對稱模式；正高斯曲率結(jié)構(gòu)中除PGC-S外均為混合模式；零高斯曲率情況下除ZGC-C外均為軸對稱模式。由上述現(xiàn)象可知，當前工況下方形截面結(jié)構(gòu)變形模式并不會隨著高斯曲率類型的改變而改變。在負高斯曲率條件下這4 類不同截面類型薄壁管更易產(chǎn)生軸對稱變形。

圖9 變形模式對比Fig. 9 Comparison of deformation modes

各結(jié)構(gòu)力-位移曲線如圖10 所示，其中方形截面結(jié)構(gòu)曲線在3 類曲率情況中均位于最下方，且其有效壓潰長度最大。負高斯曲率情況下NGC-C 結(jié)構(gòu)載荷波動程度較低且其有效壓潰長度最小，如圖10(a)所示。正高斯曲率結(jié)構(gòu)中，PGC-C 結(jié)構(gòu)載荷曲線位于PGC-H 以及PGC-O 之下，且初期載荷波動較大，如圖10(b)所示。零高斯曲率結(jié)構(gòu)中，ZGC-H 和ZGC-O 載荷變化趨勢相似性較強且曲線明顯位于ZGC-S 和ZGC-C 之上，如圖10(c)所示。由上述現(xiàn)象可知，方形截面結(jié)構(gòu)的吸能以及軸向抗沖擊性較差，圓形截面結(jié)構(gòu)在負高斯曲率下表現(xiàn)出了較小的有效壓潰長度，能夠在更短進程中使沖擊墻制動，抗沖擊性較好。

圖10 3 類高斯曲率情況下的模型力位移曲線Fig. 10 Force-displacement curves of model under three kinds of Gaussian curvature

2.2 薄壁管吸能情況分析

為進一步探究各結(jié)構(gòu)吸能特性，從力-位移曲線中提取各性能指標，將結(jié)果繪制成圖11。分析可得，對于零高斯曲率情況下的各模型，其截面類型從四邊形變化至六邊形后，總吸能、比吸能、初始峰值力以及壓潰力效率均大幅增加，有效壓潰長度則顯著減小，這是由于隨著截面邊數(shù)增加結(jié)構(gòu)抗變形能力增強，變形區(qū)域增多，變形模式發(fā)生改變。但隨著邊數(shù)進一步增加至圓截面，性能指標差異逐漸減小，上述指標變化情況與文獻[14]所得一致，其中ZGC-O 結(jié)構(gòu)的比吸能和有效壓潰長度最優(yōu)，如圖11(a)所示。上述由于截面邊數(shù)增加所產(chǎn)生的指標變化情況同樣在正高斯曲率和負高斯曲率情況下出現(xiàn)，但在正高斯曲率組中PGC-H 結(jié)構(gòu)的比吸能和有效壓潰長度最優(yōu)，如圖11(b)所示。負高斯曲率組中NGC-O 結(jié)構(gòu)的比吸能最優(yōu)；NGC-C 的有效壓潰長度指標最優(yōu)，如圖11(c)所示。

圖11 3 類高斯曲率情況下的模型性能指標比較Fig. 11 Comparison of model performance indexes under three kinds of Gaussian curvature

為在3 組模型中實現(xiàn)結(jié)構(gòu)優(yōu)選，鑒于文獻[20]的成功經(jīng)驗，本文中借助復雜比例評估法（complex proportion assessment，COPRAS）首先對比吸能、有效壓潰長度、初始峰值力和壓潰力效率四項性能指標進行了權(quán)重分配，其次對各結(jié)構(gòu)進行了綜合性能評分，從而實現(xiàn)針對各結(jié)構(gòu)的客觀評價以及準確擇優(yōu)，COPRAS 法流程如圖12 所示，圖中：wj代表性能指標的權(quán)重因子；S+、S-分別代表有利值和不利值；Qi、Ui分別代表綜合性能評分和相對分值，Ui越大說明該結(jié)構(gòu)綜合性能越好。

圖12 COPRAS 法流程圖Fig. 12 Process flow chart of COPRAS method

由于本文中期望性能優(yōu)異的吸能結(jié)構(gòu)應兼具高比吸能和低有效壓潰長度，因此考慮指標 δ 、l同時享有最高優(yōu)先級，其次為 η 、F0。根據(jù)COPRAS 法，高優(yōu)先級指標與低優(yōu)先級指標對比時，前者賦3 分，后者賦1 分；相同優(yōu)先級對比時，二者各賦2 分，最終將各指標得分數(shù)據(jù)統(tǒng)計運算后給出權(quán)重因子wj。上述4 項性能指標的權(quán)衡賦分過程以及被分配的權(quán)重因子wj列于表4，各結(jié)構(gòu)通過COPRAS 法所獲取的相關(guān)計算值以及最終綜合性能排名情況列于表5。

表4 指標 δ 、l、F0、η 的權(quán)衡賦分過程及其權(quán)重因子w jTable 4 Weighing and scoring process of four indicators ( δ , l, F0, η ) and their weighting factors ( w j )

表5 各結(jié)構(gòu)COPRAS 法相關(guān)計算值Table 5 Relevant calculated values of COPRAS method for each structure

COPRAS 法分析結(jié)果表明，NGC-C 結(jié)構(gòu)在3 組模型中綜合性能最優(yōu)，其次為ZGC-O 和PGC-H。

為進一步分析NGC-C 結(jié)構(gòu)各項性能指標，將NGC-C、PGC-H 和ZGC-O 三者的力-位移曲線繪制成圖13(a)，發(fā)現(xiàn)NGC-C 有效壓潰長度最小，載荷波動程度最小。從曲線中提取出各項性能指標，繪制成圖13(b)，發(fā)現(xiàn)NGC-C 結(jié)構(gòu)的比吸能比ZGC-O 高0.1 %，比PGC-H 低0.6 %。但有效壓潰長度比ZGC-O低12.0%，比PGC-H 低15.0%。可知在相同工況下若壓潰相同距離，NGC-C 結(jié)構(gòu)具備較大能量吸收潛力。進一步地，利用應變云圖來表現(xiàn)NGC-C 受軸向沖擊后所發(fā)生的變形情況，并取ZGC-O 為對照，如圖15 所示。其中褶皺波長定義及計算方法參照文獻[21]，褶皺波長越小說明結(jié)構(gòu)對于沖擊能吸收越充分，吸能效率越高。其次引入平均應變值來衡量結(jié)構(gòu)綜合應變程度，相同工況下該值越大說明結(jié)構(gòu)變形程度越大，能夠吸收更多沖擊能，平均應變值根據(jù)總應變與單元總數(shù)相除得出:

圖13 各組優(yōu)選模型結(jié)果比較Fig. 13 Comparison of optimal model results of each group

式中：b為單元平均應變值，εi為第i個單元所產(chǎn)生的應變值，n為單元總數(shù)。

比較后發(fā)現(xiàn)，NGC-C 結(jié)構(gòu)變形區(qū)域形狀為條帶狀，對應軸對稱變形模式下所產(chǎn)生的均勻褶皺，褶皺波長均值為25 mm，b=0.118，如圖14(a)所示；ZGC-O 結(jié)構(gòu)變形區(qū)域主要集中在8 條棱邊上且軸向分布不均勻，同樣存在條帶狀變形區(qū)域，褶皺波長均值為33 mm，吸能情況很大程度受棱邊數(shù)量的影響，同時b=0.103，如圖14(b)所示。因此，NGC-C 結(jié)構(gòu)單元平均應變值更大，褶皺波長更小。綜合以上分析可知，NGC-C 結(jié)構(gòu)具有吸能效率高，軸向抗沖擊性好等優(yōu)點，更適合作為能量吸收結(jié)構(gòu)。

圖14 NGC-C 與ZGC-O 結(jié)構(gòu)的應變情況Fig. 14 Strain of NGC-C and ZGC-O structure

3 薄壁圓管優(yōu)化設計

3.1 代理模型建立

步驟1：優(yōu)化問題。以比吸能最大化和有效壓潰長度最小化為優(yōu)化目標，參數(shù)T、G、R為優(yōu)化變量，參考薄壁管定義（外徑與壁厚之比大于20），規(guī)定參數(shù)T、G、R取值范圍分別為1～5 mm、-10～10 mm、60～100 mm。綜上，優(yōu)化問題為：

步驟2：樣本點選取。為確保優(yōu)化變量取值范圍的全覆蓋，采用拉丁超立方抽樣法進行樣本點選取，并獲取各樣本點的響應值，如表6 所示。

表6 樣本點及其響應Table 6 Sample points and responses

表6（續(xù)）Table 6 (Continued)

步驟3：響應面代理模型建立。鑒于文獻[22]中的成功案例，采用多項式響應面代理模型表達出設計變量與響應之間的關(guān)系。為驗證代理模型的準確性，采用調(diào)整后的擬合優(yōu)度系數(shù)R2Adj來衡量，該系數(shù)越接近1 說明模型準確性越好，具體模型表達式如下：

3.2 優(yōu)化結(jié)果及討論

基于上述代理模型，借助改進非支配排序遺傳算法（NSGA-Ⅱ）進行多目標優(yōu)化。優(yōu)化后所得出的非劣解組成了帕累托解集，該解集在空間上形成了帕累托前沿（Pareto front）。為從解集中選取最合適的解，使用最小距離法得到最優(yōu)點B，如圖15 所示。該點對應最優(yōu)解為T=2.56 mm、G=-9.71 mm 和R=75.48 mm，對應的響應值 δ =26.24 J/g、l=69.29 mm。根據(jù)最優(yōu)解參數(shù)建立有限元模型，將所得數(shù)值模擬結(jié)果與上述代理模型預測結(jié)果進行對比，對比結(jié)果見表7。比吸能和有效壓潰長度的數(shù)值模擬與預測結(jié)果相對誤差均在5%以內(nèi)，優(yōu)化結(jié)果準確度較高。相較于2.2 節(jié)中通過COPRAS 法評定的最優(yōu)NGC-C 結(jié)構(gòu)，優(yōu)化后所得出的NGC-C 結(jié)構(gòu)的比吸能提高了16.47%，有效壓潰長度降低了12.4%，同時結(jié)構(gòu)質(zhì)量降低了20.18%。設計人員可通過更改 α 角來改變比吸能與有效壓潰長度的權(quán)重，進而選取能夠滿足工藝性能需求的最優(yōu)解。

圖15 Pareto 解集Fig. 15 Pareto solution set

表7 數(shù)值模擬與代理模型預測結(jié)果對比Table 7 Comparison of simulation and agent model prediction results

4 結(jié) 論

本文中提出了一種負高斯曲率薄壁吸能結(jié)構(gòu)，通過文獻驗證有限元模型的可靠性，對其進行軸向動態(tài)沖擊模擬，分析了該結(jié)構(gòu)的吸能特性。借助COPRAS 法將其與各類非負高斯曲率薄壁結(jié)構(gòu)進行了性能對比。最終利用拉丁超立方抽樣法構(gòu)建樣本空間，借助響應面法建立代理模型，通過NSGA-Ⅱ遺傳算法求解出最優(yōu)解集，對設計參數(shù)進行了優(yōu)選，主要結(jié)論如下。

(1)通過數(shù)值模擬對薄壁圓管進行動態(tài)軸向壓潰所得出的力-位移曲線和變形模式均與文獻中相應結(jié)果一致性較好，所建立的有限元模型的可靠性較高。

(2)在等質(zhì)量條件下與各類非負高斯曲率薄壁結(jié)構(gòu)比對后發(fā)現(xiàn)，負高斯曲率圓管綜合性能最優(yōu)。相同動態(tài)沖擊工況下，該結(jié)構(gòu)能夠在更短的進程中使沖擊墻制動，同時能夠吸收更多的沖擊能，故其較適合作為能量吸收結(jié)構(gòu)。

(3)優(yōu)化后負高斯曲率圓管結(jié)構(gòu)的比吸能提高了16.47%，有效壓潰長度降低了12.40%，同時結(jié)構(gòu)質(zhì)量降低了20.18%。設計人員能夠根據(jù)不同工藝性能要求，從最優(yōu)解集中選取最適合的設計參數(shù)組合。