翼柱形藥柱燃面退移過(guò)程的變量化設(shè)計(jì)方法①

蔡 強(qiáng)，鮑福廷，丁 林，劉 旸

(西北工業(yè)大學(xué) 航天學(xué)院，西安 710072)

0 引言

固體火箭發(fā)動(dòng)機(jī)翼柱形藥柱具有三維復(fù)雜構(gòu)型，其燃面面積很難用解析公式表達(dá)。實(shí)體造型方法［1］是目前采用的主要方法，通常對(duì)商用CAD軟件進(jìn)行二次開(kāi)發(fā)，避免了復(fù)雜的數(shù)理方程推導(dǎo)，卻容易因特征定義不嚴(yán)格，而出現(xiàn)特征消失或自相交，導(dǎo)致幾何模型拓?fù)浠兌煨褪?，且商用CAD軟件計(jì)算大型藥柱的時(shí)間較長(zhǎng)，不利于方案的驗(yàn)證與藥柱的優(yōu)化設(shè)計(jì)。文獻(xiàn)［2］中的方法僅能處理具有1組翼的藥柱，其采用基于造型歷史追溯的燃面算法，容易發(fā)生統(tǒng)計(jì)疏漏和重復(fù)。

為了解決上述問(wèn)題，本文結(jié)合翼柱形藥柱構(gòu)型的特點(diǎn)，引入變量化設(shè)計(jì)方法，在幾何造型內(nèi)核上自主開(kāi)發(fā)。通過(guò)定義統(tǒng)一的翼形體特征函數(shù)，給出翼形體的關(guān)鍵設(shè)計(jì)變量，嚴(yán)格控制模型參數(shù)的變化，避免了造型失敗;調(diào)節(jié)三維圖形變換矩陣中的比例變換因子，縮小幾何模型的尺寸，顯著縮短了大型翼柱形藥柱的仿真時(shí)間;本文還研究新的燃面面積算法，避免統(tǒng)計(jì)燃面出錯(cuò)，并針對(duì)后翼、前后翼和大小翼等各種類(lèi)型的翼柱形藥柱，總結(jié)了統(tǒng)一的設(shè)計(jì)流程。

1 翼柱形藥柱變量化仿真建模

1.1 假設(shè)條件

(1)整個(gè)翼柱形藥柱表面光滑均一，沒(méi)有裂紋、尖點(diǎn)等缺陷。

(2)翼柱形藥柱的燃面同時(shí)被點(diǎn)燃，燃面各點(diǎn)的條件相同。

(3)燃面退移滿足平行層燃燒定律，即燃面各點(diǎn)都以相同燃速向藥柱內(nèi)推進(jìn)。

(4)翼柱形藥柱的各組翼形體周向均勻分布，滿足嚴(yán)格的循環(huán)對(duì)稱(chēng)條件。

1.2 建模仿真過(guò)程

本文將翼柱形藥柱分為外輪廓、回轉(zhuǎn)體和翼形體3種特征形體，其變量化設(shè)計(jì)過(guò)程可分為變量化約束草圖繪制、特征形體定義與生成和燃面退移仿真3步，如圖1所示。首先，交互式繪制外輪廓、回轉(zhuǎn)體和翼形體的二維截面草圖，并對(duì)草圖施加約束，定義設(shè)計(jì)變量與變量方程。然后，根據(jù)藥柱相應(yīng)特征的二維草圖，采用幾何造型算法生成外輪廓、回轉(zhuǎn)體和翼形體的三維幾何實(shí)體，外輪廓依次“減去”回轉(zhuǎn)體和翼形體，就得到藥柱的初始形狀，并調(diào)用統(tǒng)計(jì)函數(shù)獲得初始燃面面積及質(zhì)量特性信息。最后，仿真燃面退移時(shí)，按照平行層燃燒定律假設(shè)，每次退移一定的肉厚，通過(guò)改變回轉(zhuǎn)體和翼形體的特征變量，使得變化后的各形體表面與初始表面等距，外輪廓“減去”變化后的回轉(zhuǎn)體和翼形體，得到一系列退移后的藥柱形狀，計(jì)算燃燒面積和質(zhì)量特性，得到燃面面積、質(zhì)心、慣量、形狀的變化規(guī)律。

圖1 翼柱形藥柱變量化設(shè)計(jì)過(guò)程Fig.1 Variational design principle for finocyl grain

2 翼柱形藥柱變量化設(shè)計(jì)關(guān)鍵技術(shù)

在翼柱形藥柱的實(shí)際應(yīng)用過(guò)程中，采用在藥柱前段增加1組軸向前翼，或在后段周向增加另1組小翼的方法來(lái)增大初始燃面，從而為助推段提供更大推力，前者稱(chēng)為前后翼，后者稱(chēng)為大小翼。從藥柱二維截面生成三維特征形體以及燃面退移時(shí)，為了簡(jiǎn)便，通常只處理循環(huán)對(duì)稱(chēng)藥柱的1個(gè)對(duì)稱(chēng)單元，然而前后兩組翼通常不在同一對(duì)稱(chēng)平面內(nèi)，大小2組翼是周向陣列而成，因此需采取措施，確保所有翼形體都處在一份對(duì)稱(chēng)單元中，從而得到合理的布爾運(yùn)算結(jié)果。另外，還要解決翼形體的通用描述方法、穩(wěn)健準(zhǔn)確的燃面面積算法和縮短仿真時(shí)間等關(guān)鍵技術(shù)。

2.1 圖形比例變換

本文采用的是支持邊界表達(dá)方法的精確幾何造型內(nèi)核［3］，藥柱的尺寸規(guī)模越大，耗費(fèi)機(jī)時(shí)越長(zhǎng)。因此，一個(gè)直觀的想法就是在保證造型精度的條件下，對(duì)外輪廓、回轉(zhuǎn)體和翼形體的相應(yīng)二維截面進(jìn)行圖形比例變換，以縮小尺寸［4］。這樣生成的三維特征形體尺寸較小，能顯著縮短后續(xù)的布爾運(yùn)算和等距偏移運(yùn)算所耗費(fèi)的時(shí)間。

式(1)為三維圖形比例變換矩陣:

式中 (x，y，z，1)為原始截面各點(diǎn)的齊次坐標(biāo);(x*，y*，z*，1)為縮放后的特征截面相應(yīng)點(diǎn)的齊次坐標(biāo);(xf，yf，zf)為比例變換的參考點(diǎn)坐標(biāo);Sx、Sy和 Sz分別為x、y和z坐標(biāo)分量的比例變換因子，比例變換因子的取值以不損失造型精度的前提下越小越好。

翼柱形藥柱特征形體原始截面的各點(diǎn)依次乘以比例變換矩陣，就得到特征形體縮放后截面的各點(diǎn)坐標(biāo)。

2.2 翼形體特征函數(shù)定義

考慮到通用性，本文定義翼形體特征函數(shù)如下:

式中 Section為翼形體二維特征截面;w為翼形體厚度;r為翼形體頂端倒角半徑(見(jiàn)圖2);angle為該組翼形體起始角度;symN為藥柱的循環(huán)對(duì)稱(chēng)數(shù)目;npc為每個(gè)外輪廓對(duì)稱(chēng)單元中該組翼的個(gè)數(shù)。

圖2 翼形體示意圖Fig.2 Parameters of fin entity

藥柱的循環(huán)對(duì)稱(chēng)數(shù)目symN是前后2組翼個(gè)數(shù)(nf，na)或大小2組翼個(gè)數(shù)(nb，ns)的最大公約數(shù)，則

翼在外輪廓對(duì)稱(chēng)單元中的周向排列如圖3所示。翼形體起始角度angle表示第一個(gè)翼距外輪廓對(duì)稱(chēng)單元中心線的角度，其余翼的角度為

圖3 外輪廓對(duì)稱(chēng)單元中多個(gè)翼的角度分布Fig.3 Angles of multiple sets of fins defined in one cyclic symmetry cell

每組翼的第1個(gè)翼作為幾何造型的樣本，其余通過(guò)將樣本的副本繞X軸進(jìn)行三維圖形旋轉(zhuǎn)變換得到［4］，即

式中 (x，y，z，1)為樣本點(diǎn)的齊次坐標(biāo);(x*，y*，z*，1)為副本點(diǎn)的齊次坐標(biāo)。

經(jīng)過(guò)如上定義以后，即可確保所有翼形體均處在外輪廓實(shí)體對(duì)稱(chēng)單元中，從而得到嚴(yán)格旋轉(zhuǎn)對(duì)稱(chēng)的1/symN份藥柱實(shí)體，然后周向陣列即可得到完整的藥柱信息。

2.3 顏色標(biāo)識(shí)的燃面面積計(jì)算方法

本文采用邊界表達(dá)的幾何模型，其數(shù)據(jù)結(jié)構(gòu)完整記錄了模型的拓?fù)?、幾何和屬性信息?］，顏色作為基本屬性，附加在藥柱外輪廓實(shí)體和內(nèi)腔特征實(shí)體的幾何模型上，并隨著實(shí)體模型之間的布爾運(yùn)算及圖形變換，保持其合理性，具體體現(xiàn)在:

式中 藥柱實(shí)體Sgrain的表面來(lái)源于2部分:由外輪廓實(shí)體Sout的幾何模型表面分裂而成，記為SFa;由內(nèi)腔特征實(shí)體Score的幾何表面分裂而成，記為SFb。

本文用不同的顏色標(biāo)記燃面和非燃面:

(1)外輪廓實(shí)體所有表面標(biāo)識(shí)為藍(lán)色，內(nèi)腔特征實(shí)體表面標(biāo)識(shí)為紅色。因此，所有的SFa為藍(lán)色，所有的SFb為紅色。

(2)按照式(6)進(jìn)行布爾差運(yùn)算，得到的藥柱實(shí)體的表面中，源自外輪廓實(shí)體的為藍(lán)色，源自?xún)?nèi)腔特征實(shí)體的為紅色。顯然，紅色的面即為翼柱形藥柱的燃面，統(tǒng)計(jì)所有紅色面的面積，即為藥柱的初始燃面面積。

(3)內(nèi)腔特征實(shí)體擴(kuò)張以后，仍保持所有表面為紅色。同理，與不發(fā)生變化的外輪廓實(shí)體進(jìn)行布爾差運(yùn)算之后，統(tǒng)計(jì)得到的擴(kuò)張后藥柱實(shí)體的紅色表面面積，即為退移后的燃面面積。

這樣，只需在初始時(shí)進(jìn)行1次顏色標(biāo)識(shí)，就可保證燃面統(tǒng)計(jì)不發(fā)生疏漏，得到精確的計(jì)算結(jié)果。

3 翼柱形藥柱的變量化設(shè)計(jì)流程

在解決以上關(guān)鍵技術(shù)后，整個(gè)翼柱形藥柱的變量化設(shè)計(jì)流程總結(jié)見(jiàn)圖4。

圖4 翼柱形藥柱的變量化設(shè)計(jì)流程Fig.4 Variational design flow chart of finocyl grain

(1)在變量化設(shè)計(jì)環(huán)境中，交互式繪制外輪廓、回轉(zhuǎn)形內(nèi)孔、前后翼或大小翼的約束草圖，得到各個(gè)特征實(shí)體的原始二維截面。

(2)按照2.1節(jié)的方法對(duì)特征實(shí)體的截面執(zhí)行比例變換，縮小圖形的尺寸規(guī)模，以節(jié)省求解時(shí)間。

(3)由特征截面分別生成外輪廓的對(duì)稱(chēng)單元、回轉(zhuǎn)形內(nèi)孔的對(duì)稱(chēng)單元、前后翼或大小翼的三維實(shí)體，并按照2.2節(jié)中的處理方法，確保所有翼形體在外輪廓對(duì)稱(chēng)單元中。

(4)按照2.3節(jié)中介紹的顏色標(biāo)識(shí)的燃面面積計(jì)算方法，計(jì)算初始燃面面積(肉厚Web=0)和退移一定距離(肉厚Web＞0)以后的燃面面積，每退移一步依次計(jì)算燃面面積及質(zhì)量特性。

(5)判斷藥柱是否完全燃燒，即判斷當(dāng)前藥柱質(zhì)量是否為0，如果完全燃燒，就輸出燃面面積及質(zhì)量特性隨燒去肉厚的變化關(guān)系;否則，繼續(xù)執(zhí)行退移仿真。

4 實(shí)例驗(yàn)證

本文以LEDAS約束求解器［6］，構(gòu)建變量化設(shè)計(jì)環(huán)境，實(shí)現(xiàn)了上述針對(duì)翼柱形藥柱的處理方法。

4.1 可行性實(shí)例驗(yàn)證及結(jié)果分析

為了驗(yàn)證該方法的可行性，以1個(gè)有3組翼形體的藥柱作為示例，截面尺寸如圖5所示。3組翼分為前翼、后大翼、后小翼，其參數(shù)見(jiàn)表1。

圖5 大小翼和前后翼藥柱二維截面草圖Fig.5 2-D geometric constraint sketch of rear big fin and rear small fin

表1 3組翼的參數(shù)Table 1 Parameters of three sets of fins

藥柱燃面退移過(guò)程及燃面面積-肉厚關(guān)系曲線如圖6和圖7所示。

圖6 大小翼和前后翼藥柱燃面退移過(guò)程仿真及燃面面積-肉厚曲線Fig.6 Grain burnback simulation and profile of burning surface area vs web

這個(gè)燃面退移過(guò)程分為3段，初始階段燃面不斷增大，直至肉厚Web=30 mm時(shí)，后大翼開(kāi)始燒穿;Web=30～60 mm為減面燃燒，前翼和后小翼逐漸擴(kuò)張，后大翼逐漸消失;Web=60 mm時(shí)，后大翼完全消失，前翼和后小翼開(kāi)始燒穿，直至Web=85 mm藥柱完全消失，這一階段也為減面燃燒。仿真時(shí)，比例變換因子Sx=Sy=Sz=1，肉厚推進(jìn)步長(zhǎng)為1 mm，85步計(jì)算共耗時(shí)221.63 s，如果令 Sx=Sy=Sz=0.02，即縮小 50 倍，計(jì)算時(shí)間縮短至151.16 s(CPU主頻2.8 GHz)。

圖7 參考文獻(xiàn)實(shí)例的燃面退移仿真過(guò)程及燃面面積-肉厚曲線Fig.7 Grain burnback simulation and profile of burning surface area vs web

4.2 計(jì)算精度實(shí)例驗(yàn)證及結(jié)果分析

為了驗(yàn)證本文計(jì)算方法的精度，對(duì)文獻(xiàn)［7］中17#藥柱進(jìn)行燃面退移仿真，并結(jié)合推進(jìn)劑參數(shù)和噴管參數(shù)，進(jìn)行了內(nèi)彈道計(jì)算。整個(gè)藥柱肉厚為21.35 mm，燃燒初始階段(Web=0～4 mm)內(nèi)孔和翼形體同時(shí)擴(kuò)張，為增面燃燒;在Web=4 mm時(shí)，后端翼形體開(kāi)始燒穿，直至Web=12 mm翼形體完全燒完，這一階段為減面燃燒;燃燒最后階段(Web=12～21.35 mm)為近似等面燃燒過(guò)程。該例本文設(shè)置的比例變換因子為Sx=Sy=Sz=2/3，仿真肉厚步長(zhǎng)與文獻(xiàn)［7］相同，都為1 mm，22步仿真時(shí)間僅為7.344 s，而文獻(xiàn)［7］中仿真耗費(fèi)了將近30 min。

本文內(nèi)彈道計(jì)算結(jié)果與文獻(xiàn)［7］中壓力-時(shí)間進(jìn)行對(duì)比，如圖8所示。由圖8可看出，由于本文內(nèi)彈道計(jì)算未考慮點(diǎn)火段和拖尾段的非定常效應(yīng)，點(diǎn)火段和拖尾段與試驗(yàn)曲線差距較大，其余吻合相當(dāng)好(最大誤差僅3.1%)。

圖8 壓力-時(shí)間曲線計(jì)算結(jié)果與試驗(yàn)曲線對(duì)比Fig.8 Results of pressure vs time compared with experiment data

5 結(jié)論

(1)引入翼的起始角度，合理布置大小翼、前后翼的位置，使其在1個(gè)外輪廓實(shí)體對(duì)稱(chēng)單元中。調(diào)節(jié)三維比例變換矩陣的比例變換因子，以縮小藥柱的尺寸規(guī)模，顯著縮短求解時(shí)間。

(2)在基于造型歷史追溯的燃面面積統(tǒng)計(jì)方法的基礎(chǔ)上，提出改進(jìn)的顏色標(biāo)識(shí)燃面面積計(jì)算方法，該方法具有執(zhí)行效率高、操作簡(jiǎn)單、不會(huì)統(tǒng)計(jì)疏漏等優(yōu)點(diǎn)。

(3)歸納總結(jié)了翼柱形藥柱的變量化設(shè)計(jì)流程，包括特征形體變量化約束草圖的定義、比例變換、燃面退移過(guò)程仿真、統(tǒng)計(jì)燃面面積等步驟。并通過(guò)2個(gè)實(shí)例，對(duì)其有效性及計(jì)算精度進(jìn)行了驗(yàn)證。

(4)文中研究工作為翼柱形藥柱的設(shè)計(jì)、參數(shù)優(yōu)化和燃面退移仿真提供一套實(shí)用的解決方法。同時(shí)，本文的研究思路也適用于其他復(fù)雜的任意三維藥型。

［1］董新剛，陳林泉，侯曉.基于Pro/E平臺(tái)下的固發(fā)裝藥CAD軟件［C］//昆明:2002年中國(guó)宇航學(xué)會(huì)固體推進(jìn)專(zhuān)業(yè)委員會(huì)年會(huì)論文集(上)，2002:109-114.

［2］劉旸，鮑福廷，蔡強(qiáng)，等.變量化約束草圖驅(qū)動(dòng)的藥柱三維特征設(shè)計(jì)研究［J］.固體火箭技術(shù)，2010，33(1).

［3］蔡強(qiáng)，鮑福廷.基于ACIS幾何造型平臺(tái)的固體火箭發(fā)動(dòng)機(jī)裝藥設(shè)計(jì)［J］.固體火箭技術(shù)，2008，31(3).

［4］孫家廣，等.計(jì)算機(jī)圖形學(xué)(第三版)［M］.北京:清華大學(xué)出版社，1998:358-373.

［5］詹海生，李廣鑫，馬志欣.基于ACIS的幾何造型技術(shù)與系統(tǒng)開(kāi)發(fā)［M］.北京:清華大學(xué)出版社，2002:30-49.

［6］LEDAS Ltd.LEDAS geometric solver 2D overview［R］.http://ledas.com/products/lgs2d.

［7］Püskülcü G，Ulas A.3-D grain burnback analysis of solid propellant rocket motors:Part 2—modeling and simulations［J］.Aerospace Science and Technology，2008(12):585-591.