客車側(cè)翻一步碰撞算法中的接觸修正方法研究

王童那景新閆亞坤張?zhí)O蘋李婷婷

（1.長(zhǎng)安大學(xué)；2.吉林大學(xué) 汽車仿真與控制國(guó)家重點(diǎn)實(shí)驗(yàn)室）

客車側(cè)翻一步碰撞算法中的接觸修正方法研究

王童1那景新2閆亞坤2張?zhí)O蘋2李婷婷2

（1.長(zhǎng)安大學(xué)；2.吉林大學(xué) 汽車仿真與控制國(guó)家重點(diǎn)實(shí)驗(yàn)室）

提出一種彈塑性耦合迭代的接觸修正方法，將該方法用于客車側(cè)翻一步碰撞算法中，通過對(duì)假定的最大變形構(gòu)形進(jìn)行碰撞接觸判斷、彈性回彈計(jì)算及塑性能量平衡迭代計(jì)算，可獲得計(jì)算車體結(jié)構(gòu)最終變形所需的初始解變形構(gòu)形。以某12 m公路客運(yùn)車輛的典型車身段為研究對(duì)象，利用該方法進(jìn)行了客車側(cè)翻碰撞模擬試驗(yàn)，結(jié)果表明，與Ls-dyna等增量法及側(cè)翻試驗(yàn)結(jié)果對(duì)比，該算法計(jì)算效率高，試驗(yàn)次數(shù)少，模擬時(shí)間是Ls-dyna等增量法的1/10，具有實(shí)際應(yīng)用價(jià)值。

1 前言

客車側(cè)翻一步碰撞算法是對(duì)客車側(cè)翻碰撞結(jié)構(gòu)變形進(jìn)行快速模擬并快速評(píng)價(jià)側(cè)翻碰撞安全性能[1～4]的新方法,與目前的Ls-dyna等增量法軟件相比，該方法在保證工程所需計(jì)算精度的情況下可大幅提升計(jì)算效率。

由于側(cè)翻一步碰撞算法采用全量理論，在應(yīng)用Newton-Raphson方法對(duì)初始解構(gòu)形的廣義失衡力進(jìn)行平衡迭代獲得側(cè)翻碰撞結(jié)構(gòu)最終變形前，該算法需在結(jié)構(gòu)最大變形狀態(tài)下獲得可進(jìn)行Newton-Raph?son迭代的滿足能量轉(zhuǎn)換關(guān)系與結(jié)構(gòu)變形條件的初始解構(gòu)形[5]。由于側(cè)翻碰撞過程是一個(gè)自由旋轉(zhuǎn)落體過程，理論上變形形態(tài)無限多，無法準(zhǔn)確確定初始解構(gòu)形，因此，首先需在最大變形狀態(tài)下，假定1個(gè)基本滿足側(cè)翻碰撞過程能量平衡條件的符合實(shí)際變形規(guī)律的最大變形構(gòu)形，但此構(gòu)形與地面發(fā)生碰撞接觸的部分可能已侵入地面，因此需對(duì)已侵入地面的部分節(jié)點(diǎn)進(jìn)行接觸判斷與修正，以符合工程實(shí)際情況，從而得到算法所需的初始解構(gòu)形。

為此，本文提出一種彈塑性耦合迭代的碰撞接觸修正方法，對(duì)所假定的最大變形構(gòu)形進(jìn)行碰撞接觸判斷與修正，以得到可進(jìn)行Newton-Raphson迭代的初始解構(gòu)形，從而獲得側(cè)翻碰撞車體結(jié)構(gòu)的最終變形。

2 客車側(cè)翻一步碰撞算法計(jì)算過程及碰撞接觸判斷條件

2.1 計(jì)算過程

將碰撞開始狀態(tài)的車體結(jié)構(gòu)作為原始構(gòu)形｛X｝，此時(shí)車體未發(fā)生變形。碰撞開始后車體重心下降微小，可忽略不計(jì)，此時(shí)車體結(jié)構(gòu)動(dòng)能Ed最大，其計(jì)算式為：

式中，J為車體繞固定轉(zhuǎn)軸的轉(zhuǎn)動(dòng)慣量；ω為車體角速度；M為車體質(zhì)量；g為重力加速度；Δh為車體重心下降高度。

碰撞開始狀態(tài)車體結(jié)構(gòu)各節(jié)點(diǎn)的速度{v0}計(jì)算式為：

式中，ri為各節(jié)點(diǎn)到側(cè)翻固定轉(zhuǎn)軸距離；n為節(jié)點(diǎn)數(shù)。

忽略整個(gè)過程能量損失，隨著車體結(jié)構(gòu)與地面發(fā)生碰撞，車體動(dòng)能逐漸轉(zhuǎn)化為結(jié)構(gòu)形變能W，即滿足Ed=W。對(duì)于最大變形時(shí)刻的車體結(jié)構(gòu)變形，由于理論上變形形態(tài)無限多，無法準(zhǔn)確確定初始解構(gòu)形，因此需假定1個(gè)基本滿足能量平衡條件及實(shí)際變形規(guī)律的最大變形構(gòu)形{x′}，然后對(duì)該構(gòu)形進(jìn)行碰撞接觸判斷與修正，重新構(gòu)造出滿足碰撞接觸條件的最大變形構(gòu)形，將滿足能量平衡、變形條件及碰撞接觸條件的最大變形構(gòu)形{x}作為Newton-Raphson迭代初始解構(gòu)形，此時(shí)車體結(jié)構(gòu)各節(jié)點(diǎn)的位移{u}為：

由于車體結(jié)構(gòu)在空間內(nèi)變形過程除接觸碰撞力及旋轉(zhuǎn)軸支反力外，無其它外力作用[6]，此時(shí)節(jié)點(diǎn)失衡力已處于不平衡狀態(tài)，即

應(yīng)用Newton-Raphson法,按照式（5）和式（6）對(duì)此時(shí)初始解構(gòu)形下的各節(jié)點(diǎn)失衡力進(jìn)行平衡迭代求解，得到車體結(jié)構(gòu)的最終變形，使式（4）達(dá)到平衡：

2.2 碰撞接觸判斷條件

為對(duì)車體結(jié)構(gòu)進(jìn)行接觸判斷，需對(duì)側(cè)翻碰撞過程進(jìn)行簡(jiǎn)化。根據(jù)ECE R66法規(guī)，客車從側(cè)翻臨界位置開始，在重力作用下繞翻轉(zhuǎn)臺(tái)的旋轉(zhuǎn)軸自由旋轉(zhuǎn)下落，至車體結(jié)構(gòu)最大變形狀態(tài)時(shí)刻，輪胎與旋轉(zhuǎn)軸的相對(duì)位置幾乎未發(fā)生變化。為方便計(jì)算，忽略其極微小的滑動(dòng)，將整個(gè)側(cè)翻碰撞過程簡(jiǎn)化為定軸轉(zhuǎn)動(dòng)，如圖1所示。

如圖1所示，旋轉(zhuǎn)軸與X軸平行，距地面的高度為800 mm。車體與地面發(fā)生碰撞后，車體結(jié)構(gòu)各節(jié)點(diǎn)在X方向的位移幾乎不發(fā)生變化，Y方向的位移與地面侵入量沒有直接聯(lián)系，而Z方向的坐標(biāo)值直接決定車體結(jié)構(gòu)是否侵入地面。故將Z方向坐標(biāo)值作為側(cè)翻碰撞變形結(jié)構(gòu)與地面的接觸判斷條件， 即各節(jié)點(diǎn)的Z向坐標(biāo)值均不能低于地板平面的Z向坐標(biāo)值。假設(shè)A點(diǎn)（圖1）為旋轉(zhuǎn)軸上的一點(diǎn)，則地板平面的Z向坐標(biāo)值z(mì)地為：

式中，zA為A點(diǎn)的Z向坐標(biāo)值。

將最大變形構(gòu)形各節(jié)點(diǎn)的Z向坐標(biāo)值與z地對(duì)比，判斷各節(jié)點(diǎn)是否已經(jīng)侵入地面，并對(duì)相應(yīng)節(jié)點(diǎn)進(jìn)行接觸修正。接觸判斷條件為zi＜z地，zi為最大變形構(gòu)形任意節(jié)點(diǎn)i的Z向坐標(biāo)。

3 彈塑性耦合迭代的接觸修正方法

變形結(jié)構(gòu)經(jīng)過碰撞接觸判斷后，需對(duì)已經(jīng)侵入地面的節(jié)點(diǎn)進(jìn)行接觸回彈修正。參考結(jié)構(gòu)碰撞接觸修正的罰函數(shù)法[7，8]，提出一種彈塑性耦合迭代的側(cè)翻碰撞接觸修正方法。為應(yīng)用有限元結(jié)構(gòu)平衡方程計(jì)算回彈修正位移，將車體結(jié)構(gòu)碰撞后的力學(xué)模型簡(jiǎn)化為懸臂結(jié)構(gòu)模型（圖1），以側(cè)翻固定旋轉(zhuǎn)軸線的節(jié)點(diǎn)作為固定支點(diǎn)進(jìn)行碰撞接觸修正。

在車體結(jié)構(gòu)與地面碰撞過程中，按照接觸判斷條件（zi＜z地）進(jìn)行判斷后，以侵入地面各節(jié)點(diǎn)的侵入量作為參考指標(biāo)，通過對(duì)已侵入地面的節(jié)點(diǎn)施加彈性回彈修正載荷，應(yīng)用有限元結(jié)構(gòu)平衡方程進(jìn)行彈性回彈修正計(jì)算。若經(jīng)過彈性回彈修正后的車體結(jié)構(gòu)形變能還不滿足側(cè)翻碰撞過程中的能量轉(zhuǎn)換關(guān)系（即W≠Ed），則需重新計(jì)算形變能并修正車體變形，使其滿足能量平衡條件，最終得到不侵入地面且滿足能量轉(zhuǎn)換關(guān)系的初始解構(gòu)形。

3.1 接觸修正方法的彈性回彈計(jì)算

假設(shè)最大變形構(gòu)形中的節(jié)點(diǎn)j已侵入地面，以該節(jié)點(diǎn)為例進(jìn)行碰撞接觸修正，如圖2所示。將Z向坐標(biāo)變化量qjz作為節(jié)點(diǎn)與地板平面的侵入量指標(biāo)，其計(jì)算式為：

式中，zj為最大變形構(gòu)形任意節(jié)點(diǎn)j的Z向坐標(biāo)。

對(duì)比所有侵入地面節(jié)點(diǎn)的侵入量，可得到最大侵入量qzmax為：

式中，j，k，…，m為侵入地面的節(jié)點(diǎn)號(hào)。

為保證彈性回彈修正的載荷值基本在合理范圍內(nèi)，提出一種基于值域歸一化的方法，以便精確獲得合理載荷，減少計(jì)算次數(shù)，提高計(jì)算效率。

根據(jù)式（9）所得結(jié)果可找到一個(gè)β值，使得10β≤qzmax≤10β+1。令侵入地面的各節(jié)點(diǎn)法向回彈修正載荷（Z向）與侵入量成正比，并引入懲罰因子α，定義節(jié)點(diǎn)j的侵入量qjz與法向回彈修正載荷fjn之間的函數(shù)關(guān)系為：

通過調(diào)整懲罰因子α大小，對(duì)彈性修正的法向回彈載荷進(jìn)行調(diào)節(jié)，確定比較合理的節(jié)點(diǎn)法向載荷值。

懲罰因子α的大小可由式（10）估算和調(diào)節(jié)：

式中，μ為步長(zhǎng)因子。

μ的取值范圍為-3～3，當(dāng)μ取3時(shí)，步長(zhǎng)較小，節(jié)點(diǎn)j的法向回彈修正載荷fjn較小，彈性回彈修正次數(shù)增加，算法模擬速度降低；當(dāng)μ取-3時(shí)，步長(zhǎng)較大，節(jié)點(diǎn)j的法向回彈修正載荷fjn較大，彈性回彈修正次數(shù)減少，計(jì)算效率提高，但容易導(dǎo)致計(jì)算結(jié)果奇異，迭代計(jì)算不收斂。

根據(jù)車體結(jié)構(gòu)碰撞接觸修正的罰函數(shù)法，對(duì)已侵入地面的節(jié)點(diǎn)反向施加由式（10）計(jì)算的法向回彈修正載荷，應(yīng)用彈性修正有限元結(jié)構(gòu)平衡方程對(duì)最大變形各節(jié)點(diǎn)位移進(jìn)行彈性回彈修正。

彈性修正有限元結(jié)構(gòu)平衡方程為：

式中，[K]為當(dāng)前構(gòu)形的結(jié)構(gòu)整體剛度矩陣；{u′}為回彈修正位移；{fn}為侵入地面各節(jié)點(diǎn)的法向回彈修正載荷。

利用計(jì)算獲得的回彈修正位移{u′}進(jìn)行回彈修正，計(jì)算修正后的結(jié)構(gòu)構(gòu)形{X′}：

將修正后的結(jié)構(gòu)構(gòu)形進(jìn)行碰撞接觸判斷，若不滿足接觸判斷條件（zi＜z地），則繼續(xù)進(jìn)行回彈修正計(jì)算過程，直至車體變形結(jié)構(gòu)不侵入地面。

3.2 接觸修正方法的塑性能量平衡修正

經(jīng)過彈性回彈修正后的結(jié)構(gòu)構(gòu)形的形變能會(huì)產(chǎn)生變化，不再滿足W=Ed的能量平衡關(guān)系，因此需對(duì)各節(jié)點(diǎn)位移再次進(jìn)行修正，使能量重新達(dá)到平衡。

以碰撞開始狀態(tài)的原始構(gòu)形{X}作為位移計(jì)算基準(zhǔn)，計(jì)算經(jīng)過彈性回彈修正后的各節(jié)點(diǎn)位移{u″}：

對(duì)各節(jié)點(diǎn)位移{u″}進(jìn)行等比例修正：

式中，γ為等比例修正系數(shù)。

調(diào)整等比例修正系數(shù)γ，按照修正后的節(jié)點(diǎn)位移{u}重新計(jì)算當(dāng)前構(gòu)形下的結(jié)構(gòu)形變能W（式（16）），使得Ed=W。

式中，{σS}為塑性應(yīng)力；{εS}為塑性應(yīng)變。

對(duì)于重新達(dá)到能量平衡的結(jié)構(gòu)可能再次侵入地面情況，還需不斷重復(fù)進(jìn)行接觸判斷與修正，即重復(fù)式（8）～式（16）的的接觸修正過程，直至獲得滿足能量平衡且不侵入地面的初始解構(gòu)形，最后對(duì)初始解構(gòu)形的廣義失衡力進(jìn)行Newton-Raphson迭代，得到側(cè)翻碰撞車體最終變形。

4 實(shí)例分析

為檢驗(yàn)所提出的修正方法在客車側(cè)翻一步碰撞算法中的應(yīng)用效果，選擇某12 m公路客運(yùn)車輛的典型車身段作為研究對(duì)象（圖3）,對(duì)車體結(jié)構(gòu)的變形量與計(jì)算效率進(jìn)行分析。

利用CATIA軟件建立該車身段幾何模型，在Hy?permesh軟件中選擇shell單元對(duì)車體結(jié)構(gòu)進(jìn)行離散，單元大小為10 mm。模型共離散四邊形單元259 976個(gè)，節(jié)點(diǎn)258 368個(gè)，如圖4所示。應(yīng)用側(cè)翻一步碰撞算法進(jìn)行模擬，模擬結(jié)果如圖5所示。同時(shí)也進(jìn)行了Ls-dyna等增量法的分析和實(shí)際側(cè)翻試驗(yàn)，結(jié)果如圖6和圖7所示。

另外，對(duì)該車體結(jié)構(gòu)最終變形進(jìn)行了定量分析，通過側(cè)翻一步碰撞算法、Ls-dyna及側(cè)翻試驗(yàn)3種方式得到的該車身段中封閉環(huán)A和封閉環(huán)B（圖3）兩側(cè)立柱變形量如表1所列，兩側(cè)立柱變形量對(duì)比結(jié)果如圖8所示。

表1 封閉環(huán)A和封閉環(huán)B兩側(cè)立柱變形量統(tǒng)計(jì)結(jié)果

通過上述對(duì)比可知，側(cè)翻一步碰撞算法與Ls-dy?na等增量法及側(cè)翻試驗(yàn)的車體結(jié)構(gòu)變形趨勢(shì)非常吻合，變形量誤差均小于10%，符合工程計(jì)算誤差要求，保證了工程所需計(jì)算精度。

為驗(yàn)證該方法的計(jì)算效率，將圖4所示的典型車身段有限元模型作為原始模型進(jìn)行了側(cè)翻碰撞模擬和Ls-dyna等增量法模擬，結(jié)果表明，側(cè)翻一步碰撞算法模擬時(shí)間為20 min，Ls-dyna等增量法模擬時(shí)間為180 min，側(cè)翻一步碰撞算法計(jì)算時(shí)間約為L(zhǎng)s-dyna等增量法的1/10，大大提高了側(cè)翻碰撞分析的效率。

為進(jìn)一步探討車身段有限元模型單元尺寸對(duì)側(cè)翻一步碰撞算法模擬效率的影響，應(yīng)用5 mm單元重新構(gòu)造圖4中的典型車身段有限元模型，共離散四邊形單元519 846個(gè)，節(jié)點(diǎn)518 932個(gè)。利用側(cè)翻一步碰撞算法與Ls-dyna等增量法對(duì)新構(gòu)造的車身段有限元模型進(jìn)行模擬，結(jié)果表明，側(cè)翻一步碰撞算法模擬時(shí)間為30 min，Ls-dyna等增量法模擬時(shí)間為280 min。雖然單元尺寸改為5 mm后，側(cè)翻一步碰撞算法的計(jì)算時(shí)間比單元為10 mm時(shí)有所延長(zhǎng)，但仍明顯快于Ls-dyna等增量法的分析時(shí)間。

5 結(jié)束語(yǔ)

針對(duì)客車側(cè)翻一步碰撞算法，在車身結(jié)構(gòu)與地面發(fā)生碰撞的最大變形位置，提出一種彈塑性耦合迭代的客車側(cè)翻碰撞接觸修正方法，將滿足接觸判斷條件、結(jié)構(gòu)變形規(guī)律及能量轉(zhuǎn)換關(guān)系的結(jié)構(gòu)最大變形，應(yīng)用New?ton-Raphson迭代得到結(jié)構(gòu)的最終變形。以某12 m公路運(yùn)輸車輛的典型車身段為研究對(duì)象，將側(cè)翻一步碰撞算法模擬結(jié)果與Ls-dyna等增量法及側(cè)翻試驗(yàn)結(jié)果進(jìn)行對(duì)比，驗(yàn)證了所提出的接觸修正方法的有效性。

1 Guo Yingqiao,Batoz J L,EL Mouatassim M,et al.On the Es?timation of Thickness Strains in Thin Car Panels by The In?verse Approach.NUMIFORM’92，1992,1403～1408.

2 Naceur H,Guo Y Q,Batoz J L,et al.Design of Process Pa?rameters in Deep Drawing of Thin Sheets Using The Simpli?fied Inverse Approach.NUMISHEET’99,vol.1,France 1999,517～522.

3 Guo Y Q,Naceur H,Debray K，et al.Initial Solution Esti?mation to Speed Up Inverse Approach in Stamping Model?ing.Engineering Computations，2003，20（7）:810～834.

4 Na Jingxin,Wang Tong,Xu Ziwen.Research on a One-step Fast Simulation Algorithm for Bus Rollover Collision Based on Total Strain Theory.International Journal of Crashworthi?ness,2014,19（3）:275～287.

5 Lee C H,Huh H.Blank Design and Strain Estimates for Sheet Metal Forming Processes by a Finite Element Inverse Approach With Initial Guess of Linear Deformation,Journal of Materials Processing Technology,1998,145～155.

6 Tang B T,Zhao Z,Hagenah H,et al.Energy Based Algo?rithms to Solve Initial Solution in One-step Finite Element Method of Sheet Metal Stamping,Computer Methods in Ap?plied Mechanics and Engineering,2007,196,2187～2196.

7 胡美燕，姜獻(xiàn)峰.有限元分析法在接觸現(xiàn)象中的應(yīng)用研究.機(jī)電工程，2003，20（5）：160～162.

8 韓青，張毅剛，趙凱紅.結(jié)構(gòu)工程中接觸問題的數(shù)值計(jì)算方法.北京工業(yè)大學(xué)學(xué)報(bào)，2006,32（4）：321～326.

9 ECE R66,關(guān)于核準(zhǔn)大型客運(yùn)車輛上部結(jié)構(gòu)強(qiáng)度的統(tǒng)一技術(shù)法規(guī).

（責(zé)任編輯文 楫）

修改稿收到日期為2015年1月16日。

Research on the Contact Modification Method for One-step Algorithm for Bus Rollover Collision

Wang Tong,Na Jingxin,Yan Yakun,Zhang Pingping,Li Tingting
（1.Chang’an University；2.State Key Laboratory of Automotive Simulation and Control,Jilin University）

An elastic-plastic coupling iteration method is proposed,which is used in one-step algorithm for bus rollover collision.By collision contact estimation of the assumed maximum deformation contour,elastic resilience calculation and plastic energy balance iteration calculation,the original solution deformation configuration needed to calculate body structure’s final deformation can be obtained.The typical bus body section of a 12-meter bus is selected as the research object,and this method is applied in bus rollover collision simulation test,which shows that compared with Ls-dyna and rollover test,the one-step algorithm features high calculation efficiency,less test cycles,and simulation duration is only one tenths of Ls-dyna method,which has practical application value.

Bus;One-step algorithm for rollover collision;Elastic-plastic coupling iteration; Contact modification

客車 側(cè)翻一步碰撞算法 彈塑性耦合迭代 接觸修正

U462.3

A

1000-3703（2015）06-0031-04