基于協(xié)同仿真方法的U型橡膠結(jié)構(gòu)渦致振動研究

熊小慧, 唐明贊,*, 汪海燕, 朱 亮, 李小白

(1. 中南大學(xué) 軌道交通安全教育部重點實驗室, 長沙 410075; 2. 中南大學(xué) 軌道交通安全關(guān)鍵技術(shù)國際合作聯(lián)合實驗室, 長沙 410075; 3. 中南大學(xué) 軌道交通列車安全保障技術(shù)國家地方聯(lián)合工程研究中心, 長沙 410075; 4. 中車青島四方機車車輛股份有限公司, 青島 266111)

0 引 言

U型橡膠結(jié)構(gòu)外風(fēng)擋被廣泛應(yīng)用于高速列車兩車廂端部連接處，如圖1所示，使車體表面光順平滑過渡，以此來減小列車運行時的氣動阻力與氣動噪聲[1-3]。近年來，隨著旅客列車向高速化、輕量化方向發(fā)展，結(jié)構(gòu)設(shè)計中被忽略的氣動彈性問題逐漸突顯出來。當(dāng)列車高速運行時，列車表面周圍空氣的流速急劇加快，運動的空氣以列車風(fēng)的形式作用于外風(fēng)擋上，對U型橡膠結(jié)構(gòu)產(chǎn)生沿流向的阻力(拉力)和垂直于流向的升力，風(fēng)產(chǎn)生的壓力或U型橡膠外風(fēng)擋上旋渦的脫落均可能誘發(fā)結(jié)構(gòu)產(chǎn)生以流固耦合作用為特點的變形或振動現(xiàn)象[4-5]。而外風(fēng)擋的變形振動反過來又影響周圍流場運動，從而改變氣動載荷的分布和大小，將導(dǎo)致車廂端部連接處的氣動力分布不均勻，嚴(yán)重影響列車行車穩(wěn)定性。同時，外風(fēng)擋的變形振動對于結(jié)構(gòu)的疲勞壽命同樣產(chǎn)生危害。因此對于高速列車U型橡膠結(jié)構(gòu)外風(fēng)擋的氣動彈性研究顯得尤為重要。

圖1 高速列車U型橡膠結(jié)構(gòu)外風(fēng)擋

氣動彈性其本質(zhì)是氣動力、彈性力以及慣性力之間相互耦合作用發(fā)生在兩交界面上的流固耦合問題[6-7]。在工程領(lǐng)域，除了航空航天工程中的飛行器結(jié)構(gòu)外，民用工程領(lǐng)域的橋梁、煙囪、高塔、高樓等高聳結(jié)構(gòu)，機械工程領(lǐng)域的渦輪機械以及電力工程領(lǐng)域的輸電線等，都存在流固耦合問題[8-11]。對于此類較復(fù)雜結(jié)構(gòu)流固耦合的研究，主要采用風(fēng)洞試驗與數(shù)值計算方法。隨著計算機軟硬件技術(shù)的迅速提高，數(shù)值模擬方法已發(fā)展成為分析流固耦合作用的重要工具[12-13]。周岱等采用了弱耦合分區(qū)法模擬了不同雷諾數(shù)下單圓柱橫向和兩向流致振動、串列雙圓柱兩向流致振動，分區(qū)弱耦合算法計算效率較高，但難以克服強附加質(zhì)量效應(yīng)，易導(dǎo)致數(shù)值失穩(wěn)[14-15]。強耦合分區(qū)算法在每時間步內(nèi)依次迭代計算各單場方程直至收斂，具備良好的物理守恒性和程序模塊性[16]。對于低質(zhì)量比、大位移或柔性體等復(fù)雜情形可采用分區(qū)強耦合算法進(jìn)行求解。并且基于分區(qū)強耦合的模塊化優(yōu)勢，運用商業(yè)軟件計算流固耦合問題具有廣闊的工程應(yīng)用潛力[17-18]。

高速列車U型結(jié)構(gòu)外風(fēng)擋采用橡膠材料，相同的工作環(huán)境下，其強度與剛度均低于金屬材料，在氣動載荷作用下U型橡膠結(jié)構(gòu)更容易出現(xiàn)大變形及振動問題。關(guān)于橡膠材料的有限元分析，軟件Abaqus具有較完備的材料模型，而軟件STAR-CCM+對于流體的建模仿真具有較好的分析能力。因此，需要構(gòu)建軟件之間的協(xié)同仿真環(huán)境，以充分發(fā)揮各仿真平臺的優(yōu)勢，并消除在單一平臺中進(jìn)行仿真的局限性[19]。

目前，尚未有研究采用STAR-CCM+和Abaqus相結(jié)合的方法，對高速列車U型橡膠結(jié)構(gòu)外風(fēng)擋以及車體其他部件的氣動彈性問題進(jìn)行流固耦合分析。為了驗證該協(xié)同仿真方法的可行性，本文利用高速列車上拆卸的U型橡膠結(jié)構(gòu)試件進(jìn)行風(fēng)洞試驗，并根據(jù)風(fēng)洞試驗構(gòu)建仿真模型，其目的為了說明該流固耦合協(xié)同仿真方法能否較好的分析U型橡膠結(jié)構(gòu)在氣動載荷作用下的變形振動問題。在本次U型橡膠結(jié)構(gòu)試件協(xié)同仿真中，研究了不同空氣流速條件下U型橡膠結(jié)構(gòu)變形及振動情況，探討了U型橡膠結(jié)構(gòu)振動頻率與空氣旋渦脫離頻率的相互關(guān)系。關(guān)于高速列車整車外風(fēng)擋氣動彈性問題尚未有研究，此基于流固耦合協(xié)同仿真的U型橡膠結(jié)構(gòu)試件的試驗研究，是為下一步整車外風(fēng)擋結(jié)構(gòu)氣動彈性問題研究做的基礎(chǔ)性工作，可利用該協(xié)同仿真方法去分析列車不同運行工況下外風(fēng)擋的氣動響應(yīng)狀態(tài)，為高速列車外風(fēng)擋設(shè)計及應(yīng)用評估提供分析方法與技術(shù)支撐。

1 流固耦合理論基礎(chǔ)

在流固耦合分析中，常用耦合方法有直接耦合法和分區(qū)耦合法[20]。直接耦合法可同步直接求解耦合方程，但由于流固耦合問題中高度的非線性特征，很難在一種網(wǎng)格體系中用一種格式同步推進(jìn)求解，導(dǎo)致在實際應(yīng)用中存在局限性。分區(qū)耦合法最大限度地保持了流體與固體之間的獨立性，能夠充分利用現(xiàn)有的計算流體力學(xué)(CFD)和計算固體力學(xué)(CSD)的方法和程序，加上數(shù)據(jù)交換模塊，減少了計算復(fù)雜度、簡化了隱式/顯示的處理、有利于子循環(huán)求解、程序模塊化等優(yōu)點，使得分區(qū)耦合法在流固耦合問題中有了廣泛的應(yīng)用。圖2所示是一個典型的流固耦合問題，將計算域分為流體區(qū)域和固體區(qū)域，其中流體模型和固體模型的材料屬性、邊界條件等分別定義，并且在兩個模型的交界面發(fā)生相互作用，通過耦合求解獲得流體和結(jié)構(gòu)的響應(yīng)。

圖2 流固耦合模型

1.1 動力學(xué)與運動學(xué)條件

固體模型基于拉格朗日坐標(biāo)系，流體模型通常采用歐拉坐標(biāo)系。而對于流固耦合問題，由于流固界面是可變形的，因此流體模型必須建立在ALE(任意拉格朗日-歐拉)坐標(biāo)系上。并且流體模型與固體模型分區(qū)進(jìn)行計算，其相互作用在交界面上需滿足運動學(xué)和動力學(xué)條件[21]。

其流固耦合交界面的運動學(xué)條件(或位移相容性)為：

其流固耦合交界面的動力學(xué)條件(或牽引力平衡)為：

流體速度條件是由運動學(xué)條件決定，其方程為：

在流體域耦合界面上的節(jié)點位置由運動學(xué)條件決定，其他流體域節(jié)點的位移可由程序自動確定，以保持計算網(wǎng)格質(zhì)量，然后求解其ALE公式中的流體流動控制方程。

根據(jù)動力學(xué)條件，流體牽引力沿流固耦合界面整合成流體力，施加于結(jié)構(gòu)節(jié)點上，其方程為：

其中hd為固體位移的虛量。

1.2 雙向耦合的迭代計算

耦合分區(qū)法中的雙向耦合迭代計算不斷利用耦合系統(tǒng)中流體(或固體)場求解變量提供給固體(或流體)場進(jìn)行連續(xù)迭代求解[22]。

①從流體方程中求解流體節(jié)點解向量

此解向量利用固體模型中的節(jié)點位移在流體流動分析中得到，固體位移松弛因子λd(0<λd<1)。壓力向量解的殘差滿足判別條件。

②從固體方程中求解固體節(jié)點解向量

流體壓力松弛因子λτ(0<λτ<1)。

③在指定的耦合邊界條件上計算流體節(jié)點的位移為：

④壓力和位移迭代解還沒有收斂，則程序返回到步驟①并繼續(xù)進(jìn)行下一個迭代。在這種求解方法中，時間步長和求解時間由流體模型控制，并確定了控制耦合系統(tǒng)收斂的參數(shù)。這些參數(shù)包括壓力和位移公差、松弛因子、收斂準(zhǔn)則等。

1.3 STAR-CCM+與Abaqus協(xié)同仿真方法

軟件STAR-CCM+和Abaqus協(xié)同仿真方法建立在雙向耦合迭代計算的基礎(chǔ)上，利用SIMULIA協(xié)同仿真引擎在兩個軟件之間以耦合步間隔自動交換數(shù)據(jù)，這種協(xié)同仿真技術(shù)相對于文件交換方法的耦合更強。在協(xié)同仿真耦合步中，Abaqus優(yōu)先計算，即圖3協(xié)同仿真流程圖中灰色部分為一個耦合步。

圖3 STAR-CCM+和Abaqus協(xié)同仿真流程

2 U型橡膠結(jié)構(gòu)協(xié)同仿真建模及驗證

為驗證STAR-CCM+與Abaqus協(xié)同仿真方法的可行性，利用U型橡膠結(jié)構(gòu)外風(fēng)擋試件在風(fēng)洞進(jìn)行壓力位移測試，再根據(jù)風(fēng)洞試驗進(jìn)行流固耦合協(xié)同仿真分析，比較風(fēng)洞試驗和數(shù)值計算的測試點壓力、位移以及變化頻率，并基于流固耦合協(xié)同仿真方法分析了不同速度等級下U型橡膠結(jié)構(gòu)變形及振動情況。

2.1 U型橡膠結(jié)構(gòu)風(fēng)洞試驗

風(fēng)洞試驗中，在U型橡膠結(jié)構(gòu)表面布置壓力傳感器測量表面壓力變化，并使用紅外線激光位移計對橡膠變形振動進(jìn)行測量，風(fēng)洞試驗傳感器布置如圖4(a)所示。U橡膠結(jié)構(gòu)的固有動態(tài)特性由模態(tài)試驗所得，試驗如圖4(b)所示。力錘法模態(tài)試驗測量的一階固有頻率為7.29 Hz，二階模態(tài)16.59 Hz，三階模態(tài)23.54 Hz。

圖4 風(fēng)洞試驗測點布置

2.2 建立協(xié)同仿真計算模型

根據(jù)風(fēng)洞試驗中的模型及工況建立流固耦合協(xié)同仿真模型。該外風(fēng)擋試件以EPDM(三元乙丙橡膠)為基材，通過一系列工藝處理得到現(xiàn)U型結(jié)構(gòu)，其材料試驗所得該試件楊氏模量6 MPa，泊松比0.4995，密度1300 kg/m3。

1) 計算模型及測點布置。圖5所示為U型橡膠試件計算模型及測點位置，為了便于數(shù)值仿真與風(fēng)洞試驗結(jié)果的比較，計算模型的尺寸參數(shù)以及工藝孔位置與風(fēng)洞試驗試件一致，其中壓力測點p1’和位移測點d1’的位置與風(fēng)洞試驗測點布置一致，空間壓力測點p2’用于分析U型結(jié)構(gòu)上空氣旋渦脫離的頻率。

圖5 計算模型及測點布置

2) 模型離散化。STAR-CCM+和Abaqus協(xié)同仿真分析中，需要準(zhǔn)備流體和固體兩個計算模型，幾何一致且相對于坐標(biāo)軸處于同一坐標(biāo)位置。圖6(a)為Abaqus軟件中固體模型的網(wǎng)格，采用C3D8RH單元；圖6(b)為STAR-CCM+軟件中流體計算區(qū)域部分壁面網(wǎng)格示意，在近U型橡膠結(jié)構(gòu)區(qū)域進(jìn)行網(wǎng)格加密并在結(jié)構(gòu)表面和風(fēng)洞壁面都設(shè)置附面層。

(a) 固體模型離散

3) 計算區(qū)域及邊界條件。計算區(qū)域的尺寸與風(fēng)洞試驗的試驗區(qū)段尺寸一致，斷面寬1 m，高0.8 m。計算區(qū)域速度入口為均勻風(fēng)速來流，四面設(shè)為壁面，為避免壓力出口對所關(guān)心區(qū)域流場結(jié)構(gòu)的影響，并保證區(qū)域內(nèi)流場的充分發(fā)展，壓力出口設(shè)置較遠(yuǎn)，速度入口到U型橡膠結(jié)構(gòu)壁面為3 m，U型橡膠結(jié)構(gòu)壁面到壓力出口為6 m。U型橡膠結(jié)構(gòu)的固定邊界與風(fēng)洞安裝支架的固定方式相同。計算區(qū)域及邊界條件的設(shè)置如圖7所示。

圖7 流場計算區(qū)域及邊界條件

4) 模型參數(shù)設(shè)置。在STAR-CCM+中設(shè)置流場計算模型，采用三維非定常黏性不可壓縮k-ω湍流模型，網(wǎng)格變形基于RBF(徑向基函數(shù))方法，二階精度時間離散，耦合時間步長0.001 s，并設(shè)置內(nèi)迭代10次，速度入口空氣流速分別為20 m/s、30 m/s、40 m/s。

在Abaqus設(shè)置中，楊氏模量6 MPa，泊松比0.4995，密度1300 kg/m3，瑞利阻尼alpha=0.125，beta=0.00159，采用動隱式求解，耦合時間步長0.001 s，并設(shè)置內(nèi)迭代10次。并將Abaqus中固體模型的配置參數(shù)，輸出“.inp”文件，在STAR-CCM+中的聯(lián)合仿真設(shè)置中指定該文件路徑。

2.3 數(shù)值計算方法驗證

數(shù)值計算中壓力測點p1’在不同風(fēng)速下的平均值與風(fēng)洞試驗該測點位置的壓力平均值如表1所示，位移測點結(jié)果如表2所示。由于U型橡膠結(jié)構(gòu)在不同風(fēng)速的作用下產(chǎn)生變形振動，為了比較測點的響應(yīng)過程，選取30 m/s風(fēng)速條件下，位移測點風(fēng)洞試驗與數(shù)值計算中的彎曲變形振動時程曲線進(jìn)行比較，并進(jìn)行了功率譜分析，如圖8所示。

(a) 位移測點振動時程曲線

表1 U型橡膠結(jié)構(gòu)表面壓力測點結(jié)果

表2 U型橡膠結(jié)構(gòu)表面位移測點結(jié)果

通過表1中數(shù)據(jù)可以看出：風(fēng)洞試驗壓力監(jiān)測點p1的結(jié)果重復(fù)性較好，將三次風(fēng)洞試驗測量數(shù)據(jù)求取平均值并與仿真計算結(jié)果進(jìn)行比較，在風(fēng)速為40 m/s條件下數(shù)值計算與試驗的壓力最大相對誤差為2.99%。因此，對壓力場的計算該仿真方法具有一定的可靠性。

在表2中，通過對比位移測點的仿真計算與風(fēng)洞試驗結(jié)果可以發(fā)現(xiàn)，數(shù)值計算結(jié)果相對偏小，并且隨著風(fēng)速的增大，相對誤差隨之增大，當(dāng)風(fēng)速達(dá)到40 m/s時，相對誤差為8.92%。從位移比較結(jié)果來看，該流固耦合協(xié)同仿真方法用來研究U型橡膠結(jié)構(gòu)振動變形具有一定的可靠性。

通過表2與圖8可以得出，30 m/s風(fēng)速條件下，風(fēng)洞試驗位移測點的平均位移為16.5 mm，振動主頻為20.1 Hz；數(shù)值計算的平均位移為15.3 mm，振動主頻為20.5 Hz。風(fēng)洞試驗與數(shù)值計算的位移相對誤差為7.27%，屬于工程應(yīng)用可接受范圍。

3 計算結(jié)果及討論分析

3.1 不同速度等級下結(jié)構(gòu)位移及頻率

不同風(fēng)速條件下結(jié)構(gòu)由靜止到振動平穩(wěn)，測點d1’的計算全過程位移時程曲線如下圖9所示；取0～2 s時間段位移做功率譜分析如圖10所示；取4～6 s時間段位移做功率譜分析如圖11所示。

不同風(fēng)速條件下位移測點d1’在0～2 s和4～6 s內(nèi)兩個階段的振動主頻以及平均位移見下表3所示。

由圖9中可以看出，U型橡膠結(jié)構(gòu)受到風(fēng)載荷作用由靜止產(chǎn)生彎曲振動。在0～2 s時間內(nèi)，結(jié)構(gòu)受氣動力、彈性力、慣性力以及阻尼的作用，彎曲振動的振幅逐漸變小。并且隨著風(fēng)速的增加U型橡膠結(jié)構(gòu)產(chǎn)生的變形位移量越大。在4～6 s時間內(nèi)，U型橡膠結(jié)構(gòu)維持在彎曲變形狀態(tài)下等幅振動，且風(fēng)速越大，振動幅值越大。

圖9 不同速度等級下測點d1’的位移時程曲線

通過圖10和圖11中及表3中得出：1) 在0～2 s時間內(nèi)，風(fēng)速為20 m/s、30 m/s及40 m/s條件下U型橡膠結(jié)構(gòu)振動主頻分別為7.0 Hz、7.2 Hz和7.5 Hz。而該U型結(jié)構(gòu)的一階彎曲固有頻率為7.29 Hz，因此，在不同速度等級下，結(jié)構(gòu)受到氣動力作用由靜止產(chǎn)生變形振動，其振動特點為該結(jié)構(gòu)產(chǎn)生以接近自身一階彎曲固有頻率的彎曲振動。(2)在4～6 s時間內(nèi)，U型橡膠結(jié)構(gòu)在氣動力、彈性力、慣性力以及結(jié)構(gòu)阻尼的作用下變形彎曲振動幅值變小，其結(jié)構(gòu)處于一定彎曲變形狀態(tài)下的振動，風(fēng)速為20 m/s、30 m/s及40 m/s條件下U型橡膠結(jié)構(gòu)振動主頻分別為7.0 Hz、20.5 Hz和26.5 Hz。

表3 測點d1’的平均位移及振動主頻

圖10 不同速度等級下測點d1’的位移功率譜分析(0～2 s)

圖11 不同速度等級下測點d1’的位移功率譜分析(4～6 s)

3.2 渦脫頻率與結(jié)構(gòu)振動關(guān)系

圖12為不同速度等級下空間監(jiān)測點p2’的壓力時程曲線和壓力功率譜。表4中為測點p2’的壓力功率譜變化主頻和4～6s時間內(nèi)U型橡膠結(jié)構(gòu)位移測點d1’的振動位移主頻。

圖12 不同速度等級下空間點的壓力時程曲線及功率譜

表4 渦脫頻率與結(jié)構(gòu)振動主頻

由圖12中可以看出，測點p2’處的壓力變化隨風(fēng)速變化明顯，壓力變化幅值隨風(fēng)速增加而增加。壓力變化頻率明顯隨風(fēng)速的增加而增大。

在4～6 s時間內(nèi)，U型橡膠結(jié)構(gòu)的振動受旋渦脫落的影響產(chǎn)生高頻振動，當(dāng)風(fēng)速較低時，旋渦脫落的湍動能較小，不足以強迫結(jié)構(gòu)產(chǎn)生與渦脫頻率相同的振動，因此20 m/s風(fēng)速下空間點壓力變化頻率為14.2 Hz，而U型結(jié)構(gòu)振動頻率為7.0 Hz；而當(dāng)風(fēng)速較大時，渦脫頻率強迫結(jié)構(gòu)做相近頻率的振動，如30 m/s風(fēng)速下空間點壓力變化頻率為20.2 Hz，而U型結(jié)構(gòu)振動頻率為20.5 Hz。

3.3 流場壓力與速度矢量云圖

圖13、圖14所示分別為速度入口30 m/s風(fēng)速條件下，U型橡膠結(jié)構(gòu)處于4～6 s時間內(nèi)的流場壓力云圖及矢量云圖。

(a) t=4.0 s

壓力云圖及矢量云圖展示旋渦的產(chǎn)生、發(fā)展及脫離。由圖13中可以看出，在U型橡膠結(jié)構(gòu)的圓弧部位可以觀察到渦的產(chǎn)生，隨著時間的推移，逐漸從圓弧頂部脫離，并沿著流場方向向后發(fā)展，從而導(dǎo)致圓弧部位的擺動。圖14矢量云圖中通過速度矢量顯示了渦的旋轉(zhuǎn)方向及流動方向，旋渦的產(chǎn)生和脫離導(dǎo)致弧頂部位產(chǎn)生擺動，形成U型橡膠結(jié)構(gòu)的渦致振動現(xiàn)象。

(a) t=4.0 s

4 結(jié) 論

本文基于STAR-CCM+和Abaqus協(xié)同仿真方法，研究了不同空氣流速條件下U型橡膠結(jié)構(gòu)變形及振動情況，探討了U型橡膠結(jié)構(gòu)振動頻率與空氣旋渦脫離頻率的相互關(guān)系。并利用U型橡膠結(jié)構(gòu)試件風(fēng)洞試驗，驗證了該協(xié)同仿真方法用于研究此類問題的可行性。得到的主要結(jié)論有：

1) STAR-CCM+和Abaqus協(xié)同仿真方法對于橡膠結(jié)構(gòu)大變形振動問題的研究滿足工程應(yīng)用要求，具有一定的可靠性。

2) U型橡膠結(jié)構(gòu)受到氣動載荷作用由靜止產(chǎn)生彎曲變形，風(fēng)速越大其結(jié)構(gòu)產(chǎn)生的彎曲變形量越大。氣動力作用下結(jié)構(gòu)彎曲變形同時產(chǎn)生彈性力、慣性力以及阻尼力，各作用力相互耦合導(dǎo)致結(jié)構(gòu)產(chǎn)生振動，其振動特點為接近該結(jié)構(gòu)自身一階彎曲固有頻率的彎曲振動。

3) 當(dāng)結(jié)構(gòu)受到氣動力、彈性力、慣性力趨于平衡，其結(jié)構(gòu)處于一定彎曲變形狀態(tài)下的振動。風(fēng)速較低，旋渦脫離的能量不足以強迫U型橡膠結(jié)構(gòu)的振動與渦脫頻率相一致。隨著風(fēng)速的增大，U型橡膠結(jié)構(gòu)的振動是由旋渦脫離導(dǎo)致的受迫振動。