1. 整流罩地面分離過程仿真

1.1 仿真動機

對于整流罩地面分離過程仿真的必要性，從以下幾個方面考慮：

1）為了準(zhǔn)確模擬整流罩在高空的分離過程，需要在大型真空罐內(nèi)進(jìn)行分離試驗。然而，由于整流罩真空分離試驗所需成本和對設(shè)備要求均較高，目前國內(nèi)尚不具備實驗條件。

2）隨著我國航天水平的進(jìn)步，對大推力火箭的需求日益迫切。為能承載更多有效載荷，大推力火箭要求配備大型整流罩。整流罩的尺寸增大使其結(jié)構(gòu)剛度相對減小，地面環(huán)境下氣動載荷帶來的影響相對增大，可能導(dǎo)致出現(xiàn)如彎曲、扭轉(zhuǎn)、呼吸等復(fù)雜的變形情況，不能視其為剛體。

1.2 仿真難點

整流罩地面分離過程仿真中存在的難點包括：

1）為提升仿真結(jié)果的精確性，要求整流罩模型具有一定精細(xì)度。省略一些不必要的結(jié)構(gòu)或作適當(dāng)簡化，加強關(guān)鍵結(jié)構(gòu)的準(zhǔn)確性和精確度，如何掌握這兩者的平衡是難點之一。

2）為實現(xiàn)仿真過程的高效性，必須考慮已有結(jié)構(gòu)的網(wǎng)格劃分。如何選取網(wǎng)格類型和網(wǎng)格密度，使之在不影響結(jié)果精確性的前提下實現(xiàn)較高的計算速度，又是一大難點。

3）整流罩地面分離是一典型的流固耦合問題。在此過程中，罩外大氣附加整流罩氣動力，影響其運動和變形；整流罩的運動與變形反過來又會影響罩外大氣的流動。如何選擇具有解決流固耦合問題能力的大型有限元軟件以及如何選取合理高效的算法，成為整流罩地面分離過程仿真的首要難點。

1.3 仿真路線

以上，在仿真動機驅(qū)動下，重點關(guān)注仿真難點，同時考慮甲方對于大型有限元軟件的要求，選擇基于Abaqus/Explicit求解器的耦合歐拉-拉格朗日（CEL，coupled Eulerian-Lagrangian）算法分析整流罩在氣動力作用下的分離特性，從而實現(xiàn)地面分離過程仿真。

1.4 耦合歐拉-拉格朗日算法優(yōu)勢分析

1.4.1 拉格朗日算法和歐拉算法

拉格朗日算法常用于固體力學(xué)中的受力與變形分析。其重要特征是有限元網(wǎng)格固連于材料區(qū)域且兩者共享邊界，所形成的拉格朗日單元內(nèi)充滿材料。因此，結(jié)構(gòu)變形一致反映于有限元網(wǎng)格的變化，可跟蹤節(jié)點的運動，從而簡化控制方程的求解過程。但在大變形情況下，網(wǎng)格發(fā)生嚴(yán)重畸變，此時拉格朗日算法喪失了其準(zhǔn)確性。

與拉格朗日算法相比，歐拉算法可有效應(yīng)用于大變形問題，如液體晃動、氣體流動以及滲流等。有限元網(wǎng)格固定于空間，其形狀、大小、位置不隨結(jié)構(gòu)變形而變化。一般地，歐拉算法所形成的歐拉單元難以被同種材料填滿或者無任何材料。因此，難以準(zhǔn)確描述結(jié)構(gòu)的材料邊界。圖1對比了拉格朗日算法和歐拉算法的單元特性。

a）拉格朗日算法 b）歐拉算法

圖1 拉格朗日算法和歐拉算法的單元特性

1.4.2 耦合歐拉-拉格朗日算法

1.4.2.1 概述

耦合歐拉-拉格朗日算法由學(xué)者Noh提出，最初應(yīng)用于帶有移動邊界的二維流體動力學(xué)問題。CEL算法吸收了拉格朗日算法和歐拉算法的優(yōu)點并克服了兩者的缺點。對固體建立拉格朗日模型，劃分拉格朗日網(wǎng)格；對流體建立歐拉模型，劃分歐拉網(wǎng)格。兩類網(wǎng)格重疊處是耦合區(qū)，能夠高效傳遞計算中的信息。

1.4.2.2 理論基礎(chǔ)

1）控制方程由質(zhì)量守恒、動量守恒、能量守恒及連續(xù)性方程組成。本次仿真的歐拉材料(即罩外大氣）由本構(gòu)方程與連續(xù)性方程描述，分別為：

(1)

其中：—應(yīng)變張量

—大氣壓強

—剪切粘度

—應(yīng)變變化率

(2)

其中：—大氣壓強

—標(biāo)準(zhǔn)大氣壓

—大氣密度

—氣體常數(shù)

—大氣溫度

—絕對零度

2）有限差分形式的時間積分。當(dāng)前解由前一步獲得而不需迭代。其數(shù)值穩(wěn)定性由臨界時間步判斷

(3)

其中：—單元特征尺寸

—波速

3）罰函數(shù)耦合系數(shù)追蹤拉格朗日節(jié)點和歐拉材料界面上錨定點的相對位移，檢查每一個節(jié)點對歐拉物質(zhì)界面的貫穿，如果節(jié)點不出現(xiàn)貫穿，則無需進(jìn)行處理；如果發(fā)生節(jié)點對物質(zhì)界面的貫穿，界面力就會分布到歐拉材料的錨定點上，如圖2所示。界面力的大小與發(fā)生的貫穿相對位移成正比，即

(4)

其中：—罰剛度系數(shù)，僅取決于材料特性。

圖2 罰函數(shù)耦合算法

1.4.2.3 研究現(xiàn)狀

CEL算法在解決流固耦合問題時具有強大的優(yōu)勢，目前已有廣泛的應(yīng)用。左輝等人對比了簡單氣體壓縮過程的仿真結(jié)果和解析結(jié)果，成功證明了Abaqus中理想氣體狀態(tài)方程的可靠性與正確性，并應(yīng)用于產(chǎn)品分析。姚小虎等人采用CEL算法對水陸兩棲飛機水上降落的流固耦合問題進(jìn)行了結(jié)構(gòu)響應(yīng)分析，得到的高精度結(jié)果為結(jié)構(gòu)設(shè)計提供了十分有效的參考。徐文杰采用CEL算法模擬分析了滑坡涌浪這一復(fù)雜的流固耦合問題，模擬結(jié)果與試驗結(jié)果的高吻合度證明了CEL算法的巨大優(yōu)越性。宋慶華利用ABAQUS的CEL功能實現(xiàn)燃油箱晃動的流固耦合分析，為產(chǎn)品設(shè)計提供了寶貴的參考。

1.4.2.4 Abaqus應(yīng)用

Abaqus/Explicit求解器包含了CEL算法，且具有如下特性與優(yōu)勢：

1）高速計算。Abaqus/Explicit求解器對CEL算法進(jìn)行優(yōu)化，即使計算資源有限，仍能以較快速度運行，保證了計算的穩(wěn)定性。

2）高度集成。將流體前處理與典型的固體前處理一同集成于Abaqus/CAE模塊中，消除調(diào)用其它模塊甚至其他軟件的麻煩。

3）運動流體網(wǎng)格。于Abaqus 6.10引入的運動流體網(wǎng)格僅需用戶在初始分析步內(nèi)包含固體結(jié)構(gòu)，而非整個運動軌跡，極大減少了計算成本。

4）歐拉體積分?jǐn)?shù)（EVF，Eulerian volume fraction）。基于流體體積法，用戶可定義每個歐拉單元內(nèi)的歐拉體積分?jǐn)?shù)，若材料充滿某一單元，則體積分?jǐn)?shù)為1；若單元內(nèi)無材料，則體積分?jǐn)?shù)為0。歐拉體積分?jǐn)?shù)的定義極大方便了復(fù)雜歐拉材料的建模。

5）通用接觸設(shè)置。用戶只需在通用接觸中選擇接觸類型，即可定義流體材料與結(jié)構(gòu)幾何邊界的歐拉-拉格朗日接觸。

1.4.3 光滑粒子流體動力學(xué)（SPH，smoothed particle hydrodynamic）算法

光滑粒子流體動力學(xué)是一種無網(wǎng)格數(shù)值方法，用大量離散的光滑粒子的集合定義流體材料，從而通過插值離散化連續(xù)方程組。其理論核心為核函數(shù)，實質(zhì)是一定光滑長度范圍內(nèi)其他臨近粒子對研究粒子影響程度的權(quán)函數(shù)，如圖3所示。

圖3 核函數(shù)

SPH可有效解決流體流動和自由表面等問題，但存在以下不足：

1）屬于傳統(tǒng)的拉格朗日方法，在處理大變形問題時結(jié)果準(zhǔn)確性不如CEL。

2）其計算效果依賴于對于光滑粒子的建模，若粒子布置不均勻，可能會造成流體材料質(zhì)量分布不均。

3）適合單線程計算，當(dāng)粒子數(shù)量較多時對多線程計算有限制。

4）Abaqus使用SPH的步驟較為繁瑣，首先在Abaqus/CAE中創(chuàng)建質(zhì)量單元，再寫入輸入文件。

1.5 Abaqus/Explicit優(yōu)勢分析

有限元仿真領(lǐng)域?qū)α鞴恬詈蠁栴}一般采用計算流體動力學(xué)（CFD，Computational Fluid Dynamics），目前較為主流的求解器包括CFX、Fluent、COMSOL。這些傳統(tǒng)求解器相對于Abaqus/Explicit的劣勢包括：

1）以流體為核心，流體材料產(chǎn)生運動或變形后再將信息傳遞給固體結(jié)構(gòu)，與整流罩地面分離過程截然相反，難以滿足初始時刻的情況。

2）不同程度簡化了氣動載荷，對結(jié)果精度產(chǎn)生較大影響。落?壽曾利用Fluent進(jìn)行整流罩地面分離試驗仿真；張小偉等分析了整流罩分離過程的模態(tài)。兩學(xué)者基于相似的簡化原則，未考慮氣動力的時變影響和空間分布不均勻性。

3）對固體材料仿真的求解能力較差。由于流固耦合分析對邊界條件、模型尺寸、單元大小等因素極為敏感，要求同其他求解器聯(lián)合求解，不利于仿真的二次開發(fā)。

1.6 技術(shù)方案

1.6.1 概述

圖4為整流罩地面分離仿真的基本步驟。其中，固體結(jié)構(gòu)的建模如前所述，在此不贅述。下文將重點闡述Abaqus CAE中的流體材料建模、裝配部件、設(shè)置相互作用及施加載荷環(huán)節(jié)，為二次開發(fā)作鋪墊。

圖4 整流罩地面分離仿真的基本步驟

1.6.2 關(guān)鍵步驟

1） 流體材料建模。

流體材料為罩外大氣，在創(chuàng)建部件時命名為“atmosphere”，類型為三維空間下的歐拉區(qū)域。由于整流罩半罩為半圓柱，宜采用旋轉(zhuǎn)方式建模（見圖5）。考慮到二次開發(fā)的實現(xiàn)，應(yīng)簡化建模操作，不妨設(shè)為半圓柱并進(jìn)行partition操作（見圖6）。為創(chuàng)建重疊網(wǎng)格區(qū)域并設(shè)置耦合歐拉-拉格朗日接觸，兩模型應(yīng)有重合區(qū)域，故取分割面位置由內(nèi)靠近整流罩蒙皮內(nèi)表面。

圖5 圖6

        采用理想氣體狀態(tài)方程（2）式描述流體材料的屬性

(2)

選擇Abaqus材料庫EOS模型中的Ideal Gas類型（見圖7）。其中：Gas Constant為氣體常數(shù)，取為；Ambient Pressure為環(huán)境壓強，取為。此外還需設(shè)置氣體密度為（見圖8），動力粘度系數(shù)為（見圖9），和絕對零度為-273（見圖10）。

圖7

圖8

圖9

圖10

    選擇EC3D8R型歐拉單元為其劃分網(wǎng)格（見圖11）。

圖11

2） 裝配部件

如前所述，調(diào)整部件相對位置，使固體結(jié)構(gòu)與流體材料有重合部分，重合部分為整流罩蒙皮。

3） 設(shè)置相互作用

對接觸性質(zhì)進(jìn)行編輯（見圖12），選擇罰函數(shù)定義切向行為，其中Friction Coeff為罰剛度系數(shù)，僅取決于材料特性。在通用接觸中定義相互作用，賦予接觸性質(zhì)，關(guān)鍵是定義相互接觸面（見圖13）。對于固體結(jié)構(gòu)，接觸面為整流罩半罩蒙皮外表面，類型為幾何；對于流體材料，接觸面為整體，類型為網(wǎng)格（見圖14）。

圖12

圖13

圖14

4） 施加載荷

僅為整流罩固體結(jié)構(gòu)施加重力場。創(chuàng)建預(yù)定義場（見圖15），為已進(jìn)行partition操作的罩外大氣進(jìn)行材料設(shè)置，罩外部分設(shè)為1，罩內(nèi)部分使用缺省值0（見圖16）。

圖15 圖16