天宮飛行器低軌控空氣動(dòng)力特性一體化建模與計(jì)算研究

李志輝，吳俊林，彭傲平，唐歌實(shí)（.中國(guó)空氣動(dòng)力研究與發(fā)展中心超高速所，綿陽(yáng)62000；2.國(guó)家計(jì)算流體力學(xué)實(shí)驗(yàn)室，北京009；.北京航天飛行控制中心航天飛行動(dòng)力學(xué)技術(shù)重點(diǎn)實(shí)驗(yàn)室，北京00094）

·工程技術(shù)·

天宮飛行器低軌控空氣動(dòng)力特性一體化建模與計(jì)算研究

李志輝1，2，吳俊林1，彭傲平1，2，唐歌實(shí)3
（1.中國(guó)空氣動(dòng)力研究與發(fā)展中心超高速所，綿陽(yáng)621000；2.國(guó)家計(jì)算流體力學(xué)實(shí)驗(yàn)室，北京100191；3.北京航天飛行控制中心航天飛行動(dòng)力學(xué)技術(shù)重點(diǎn)實(shí)驗(yàn)室，北京100094）

對(duì)非規(guī)則板艙組合體天宮飛行器300～200 km低軌道飛行過(guò)程空氣動(dòng)力特性一體化計(jì)算建模，提出考慮復(fù)雜構(gòu)型物面遮蓋效應(yīng)面元解析法與經(jīng)修正的Boettcher/Legge非對(duì)稱橋函數(shù)，發(fā)展基于三角形面元逼近復(fù)雜外形通用處理方法，建立適于天宮飛行器復(fù)雜物形處理與面元?dú)鈩?dòng)力系數(shù)計(jì)算規(guī)則；將DSMC方法與求解Boltzmann模型方程氣體運(yùn)動(dòng)論統(tǒng)一算法應(yīng)用于天宮飛行器簡(jiǎn)化外形，進(jìn)行氣動(dòng)力當(dāng)?shù)鼗P(guān)聯(lián)參數(shù)計(jì)算修正，建立針對(duì)大型復(fù)雜結(jié)構(gòu)天宮飛行器低軌道飛行控制過(guò)程空氣動(dòng)力特性一體化快速算法與程序軟件。對(duì)大尺度圓柱體外形與天宮飛行器300～200 km不同高度變軌飛行過(guò)程不同迎角與側(cè)滑角及帆板平面與本體主軸不同夾角復(fù)雜構(gòu)型氣動(dòng)力特性計(jì)算分析驗(yàn)證，表明天宮飛行器在200 km以上低軌道飛行控制過(guò)程中所受空氣動(dòng)力系數(shù)隨飛行高度發(fā)生顯著變化（8%～50%），證實(shí)長(zhǎng)期在軌運(yùn)行的大型航天器若采用統(tǒng)一固定的氣動(dòng)力系數(shù)，誤差累積巨大，需要采取防護(hù)措施，低軌道飛控大氣阻力仍是制約航天器定軌預(yù)報(bào)精度最關(guān)鍵因素。

天宮飛行器；低軌控空氣動(dòng)力特性；當(dāng)?shù)鼗瘶蚝瘮?shù)；關(guān)聯(lián)參數(shù)確定；工程計(jì)算；DSMC方法；統(tǒng)一算法

1　引言

航天器從離地面數(shù)百公里外層空間往地球低軌道再入飛行過(guò)程，是一個(gè)自數(shù)百到數(shù)十量級(jí)Knudsen（克努森）數(shù)高稀薄自由分子流動(dòng)狀態(tài)多物理場(chǎng)復(fù)雜構(gòu)型極高超聲速流動(dòng)問(wèn)題［1］。低地球軌道環(huán)境，是對(duì)地觀測(cè)衛(wèi)星、氣象衛(wèi)星、通訊遙感衛(wèi)星、載人飛船、空間站等大型復(fù)雜結(jié)構(gòu)航天器運(yùn)行區(qū)域。距離地面200～1000 km的低地球軌道環(huán)境大氣壓力很低、氣體非常稀薄，分子間斷粒子效應(yīng)與稀疏氣體分子分布變化通過(guò)對(duì)航天器產(chǎn)生力/熱沖擊作用，直接影響航天器運(yùn)行軌道、姿態(tài)和壽命。為了研究和控制航天器低地球軌道環(huán)境飛行過(guò)程，在過(guò)去計(jì)算機(jī)條件不夠先進(jìn)情況下，傳統(tǒng)做法是近似處理，不同流動(dòng)區(qū)域發(fā)展相應(yīng)的計(jì)算方法。對(duì)于Knudsen數(shù)很大的稀薄氣體自由分子流區(qū)，氣體流動(dòng)由無(wú)碰撞或近于無(wú)碰撞的自由分子流理論控制。為了進(jìn)行航天器從340～200 km低軌道飛行試驗(yàn)，通常認(rèn)為飛行器在200 km以上，氣體分子平均自由程達(dá)兩百多米以上，其Knudsen數(shù)數(shù)十以上，飛行繞流已屬于完全自由分子流，阻力、升力等空氣動(dòng)力不再隨飛行高度發(fā)生變化，而采用固定的氣動(dòng)力系數(shù)參與飛行力學(xué)軌道計(jì)算［2］。于是大型航天器開(kāi)展低軌道飛行試驗(yàn)采用固定阻力系數(shù)2.2，發(fā)現(xiàn)飛行偏差有時(shí)越來(lái)越大，以致使用其上RCS姿控系統(tǒng)強(qiáng)制控回標(biāo)稱軌道，造成燃料消耗，甚至有時(shí)難以準(zhǔn)確有效控制［3，4］。為此，提出這樣的問(wèn)題：我國(guó)天宮飛行器低軌道飛行所受大氣阻尼系數(shù)是否會(huì)變化？在軌天宮飛行器本體直徑數(shù)米、長(zhǎng)度達(dá)十余米、其攜帶的太陽(yáng)電池翼帆板展開(kāi)寬度數(shù)十米，對(duì)如此復(fù)雜巨大不規(guī)則板艙組合體飛行器實(shí)施低軌道飛行控制、軌道參數(shù)計(jì)算時(shí)，能否采用統(tǒng)一、固定的空氣動(dòng)力系數(shù)？如何開(kāi)展天宮飛行器低軌道飛行過(guò)程空氣動(dòng)力特性一體化建模與計(jì)算研究［4］，剖析天宮飛行器氣動(dòng)特性隨飛行高度變化規(guī)律？為此，本文在研究建立航天空氣動(dòng)力學(xué)統(tǒng)一計(jì)算理論與各流動(dòng)區(qū)域多物理場(chǎng)復(fù)雜流動(dòng)機(jī)理模擬方法研究基礎(chǔ)上，嘗試開(kāi)展上述問(wèn)題計(jì)算研究。

連接微觀分子動(dòng)力學(xué)與宏觀流體力學(xué)的介觀Boltzmann（玻爾茲曼）速度分布函數(shù)方程［5］本身可描述各個(gè)流域氣體分子輸運(yùn)現(xiàn)象，該方程自1872年提出以來(lái)，一直是國(guó)際學(xué)術(shù)界追蹤研究的問(wèn)題，它作為一個(gè)高度復(fù)雜高維積分、微分、多相空間多尺度非線性剛性問(wèn)題［5，6］，精確求解描述各流域氣體流動(dòng)特征的Boltzmann方程至今未成現(xiàn)實(shí)。為此，近五十年來(lái)世界航天發(fā)展了基于微觀分子動(dòng)力學(xué)將分子運(yùn)動(dòng)與碰撞解耦進(jìn)行隨機(jī)統(tǒng)計(jì)模擬的DSMC（直接模擬Monte Carlo）方法，該方法自1963年Bird［7］將其發(fā)展用于稀薄氣體流動(dòng)模擬至今，已在稀薄氣體動(dòng)力學(xué)領(lǐng)域獲得了廣泛的應(yīng)用與檢驗(yàn)，國(guó)內(nèi)外學(xué)者將DSMC方法用于模擬包括內(nèi)能、化學(xué)反應(yīng)、電離、輻射等多種作用物理模型的稀薄流區(qū)氣體流動(dòng)問(wèn)題研究，已取得較好的應(yīng)用發(fā)展［8-17］，并利用巨型計(jì)算機(jī)研究并行算法，引入中樞網(wǎng)格系統(tǒng)開(kāi)展并行計(jì)算［18，19］。在DSMC方法對(duì)較高Knudsen數(shù)稀薄氣體流動(dòng)仿真取得巨大成功同時(shí)，該方法因受網(wǎng)格劃分、時(shí)間步長(zhǎng)和模擬分子數(shù)等方法本身模擬準(zhǔn)則所限，難以對(duì)低Knudsen數(shù)復(fù)雜大尺度飛行器氣動(dòng)問(wèn)題數(shù)值仿真［14，20，21］。

為了探索跨流域氣體流動(dòng)問(wèn)題一體化模擬方法，通過(guò)跟蹤國(guó)際上關(guān)于Boltzmann方程研究現(xiàn)狀與發(fā)展趨勢(shì)［22-27］，結(jié)合從事DSMC方法與計(jì)算流體力學(xué)有限差分法研究基礎(chǔ)，過(guò)去十余年本文第一作者等從研究Boltzmann方程碰撞松弛間隔理論出發(fā)，確立描述各流域微觀分子輸運(yùn)現(xiàn)象統(tǒng)一的Boltzmann模型方程，提出并應(yīng)用離散速度坐標(biāo)法對(duì)分子速度分布函數(shù)數(shù)值離散，發(fā)展可用于速度空間宏觀流動(dòng)取矩的離散速度數(shù)值積分法。將計(jì)算流體力學(xué)有限差分方法推廣拓展到基于時(shí)間、位置空間與速度空間的Boltzmann模型方程數(shù)值求解，先后建立起模擬各流域一維、二維、三維繞流問(wèn)題氣體運(yùn)動(dòng)論統(tǒng)一算法（GKUA）理論與系列計(jì)算技術(shù)［28-33］，并開(kāi)展了飛船返回艙、近空間飛行器、箔條云等跨流域氣動(dòng)問(wèn)題研究［33-37］，解決了系列工程研制需求問(wèn)題。由于上述求解Boltzmann模型方程統(tǒng)一算法將氣體分子速度分布函數(shù)基于位置空間與速度空間的信息均保存起來(lái)數(shù)值求解，需要在位置空間與速度空間組成的六維多相空間進(jìn)行計(jì)算，由此占用大量計(jì)算機(jī)內(nèi)存資源。如果進(jìn)行三維高超聲速氣體流動(dòng)問(wèn)題數(shù)值計(jì)算，如天宮飛行器跨流域繞流問(wèn)題研究，必須開(kāi)展大規(guī)模并行計(jì)算。

另一方面，航天器軌道計(jì)算可靠性直接依賴于軌道動(dòng)力學(xué)方程空氣動(dòng)力項(xiàng)求解的準(zhǔn)確性，需要與飛行力學(xué)數(shù)值推進(jìn)方法耦合，計(jì)算確定航天器任意時(shí)刻空間位置與速度，需要發(fā)展跨流域空氣動(dòng)力特性快速計(jì)算方法。為了能在有限時(shí)間，建立天宮飛行器低軌道飛行過(guò)程空氣動(dòng)力一體化計(jì)算平臺(tái)，擬以稀薄氣體高超聲速繞流當(dāng)?shù)鼗瘶蚝瘮?shù)理論為基礎(chǔ)，開(kāi)展天宮飛行器低軌道飛行控制過(guò)程300～200 km空氣動(dòng)力特性當(dāng)?shù)鼗瘶蚝瘮?shù)工程計(jì)算建模，研究使用DSMC方法、求解Boltzmann模型方程的氣體運(yùn)動(dòng)論統(tǒng)一算法對(duì)天宮飛行器簡(jiǎn)化外形典型繞流狀態(tài)計(jì)算，修正適于大尺度復(fù)雜結(jié)構(gòu)航天器當(dāng)?shù)鼗P(guān)聯(lián)參數(shù)，建立適于天宮飛行器低地球軌道環(huán)境不同高度、馬赫數(shù)與不同迎角、側(cè)滑角空氣動(dòng)力特性快速計(jì)算分析應(yīng)用研究框架。

2　適于天宮飛行器低軌控不同高度空氣動(dòng)力特性當(dāng)?shù)鼗焖偎惴?/h2>

稀薄氣體流動(dòng)領(lǐng)域氣動(dòng)力系數(shù)計(jì)算所用的當(dāng)?shù)鼗こ逃?jì)算方法，是基于半經(jīng)驗(yàn)橋函數(shù)理論的當(dāng)?shù)貥蚬椒ǎ?8，39］。自由分子流和牛頓連續(xù)流兩個(gè)極限流態(tài)都可以用當(dāng)?shù)貧鈩?dòng)力的解析來(lái)表述，當(dāng)?shù)乇砻婷嬖系臍鈩?dòng)力系數(shù)只依賴于來(lái)流和當(dāng)?shù)匦再|(zhì)，如當(dāng)?shù)赜恰⒈砻孀饔玫?。?dāng)?shù)剡B續(xù)流和自由分子流系數(shù)之間的間隙可以用半經(jīng)驗(yàn)的橋公式光滑搭接起來(lái)。這是一種加權(quán)函數(shù)，依賴于獨(dú)立的關(guān)聯(lián)或相似參數(shù)（如克努森數(shù)），通常由實(shí)驗(yàn)結(jié)果與數(shù)值計(jì)算分析確定。對(duì)連續(xù)流區(qū)域的壓力系數(shù)計(jì)算，采用修正的牛頓無(wú)粘流理論以及背風(fēng)真空效應(yīng)法進(jìn)行估算；對(duì)自由分子流區(qū)，采用基于不同模型材料修正的Nocilla壁面反射模型進(jìn)行壓力系數(shù)和摩擦力系數(shù)計(jì)算；對(duì)過(guò)渡流區(qū)域，基于修正Boettcher/Legge非對(duì)稱橋函數(shù)理論［40］，使用可分段描述的非對(duì)稱壓力與摩擦力系數(shù)關(guān)聯(lián)橋函數(shù)。通常針對(duì)軸對(duì)稱鈍體飛行器設(shè)計(jì)，能預(yù)測(cè)估算鈍柱形高超聲速飛行器的稀薄氣動(dòng)力特性，但本文研究大型復(fù)雜結(jié)構(gòu)航天器如天宮飛行器，其空氣動(dòng)力特性計(jì)算不能沿用常規(guī)旋成體當(dāng)?shù)鼗?jì)算方法［41］?？裳芯渴褂肈SMC方法、求解Boltzmann模型方程從稀薄流到連續(xù)流統(tǒng)一算法［28-37］，基于準(zhǔn)確計(jì)算進(jìn)行當(dāng)?shù)鼗瘶蜿P(guān)聯(lián)系數(shù)的選擇調(diào)試與檢驗(yàn)修正，確定各關(guān)聯(lián)參數(shù)的最佳組合［35，42］，發(fā)展針對(duì)縱橫數(shù)十米大尺度復(fù)雜結(jié)構(gòu)航天器，跨越高稀薄自由分子流、過(guò)渡流到大氣層內(nèi)低高度各流區(qū)、高低不同馬赫數(shù)、不同迎角與側(cè)滑角氣動(dòng)力系數(shù)快速計(jì)算方法。

2.1高超聲速空氣動(dòng)力特性當(dāng)?shù)鼗?jì)算方法

當(dāng)?shù)鼗瘶蚬椒ㄊ紫瓤紤]的是在兩個(gè)邊界流域連續(xù)流和自由分子流區(qū)域，氣動(dòng)力特性可用兩個(gè)不同的數(shù)學(xué)函數(shù)稱為連續(xù)流函數(shù)FCont和自由分子流函數(shù)FFM來(lái)描述，而在過(guò)渡流區(qū)域需要給出一個(gè)新的函數(shù)Fb，亦即橋函數(shù)［4，35，38，39］，使其在連續(xù)流或自由分子流分別趨近FCont和FFM。對(duì)給定入射角θ，當(dāng)?shù)孛嬖辖?jīng)歸一化壓力和摩擦力系數(shù)表達(dá)為式（1）：

式中，F(xiàn)b，p和Fb，τ分別是壓力和摩擦力橋函數(shù)，依賴于獨(dú)立參數(shù)Kn、TW/T0和θ。

對(duì)連續(xù)流區(qū)壓力系數(shù)計(jì)算，可采用諸如修正牛頓無(wú)粘流理論及背風(fēng)真空效應(yīng)法計(jì)算；對(duì)自由分子流區(qū)，可使用基于不同模型材料修正的Nocilla壁面反射模型進(jìn)行壓力與摩擦力系數(shù)計(jì)算；對(duì)過(guò)渡流區(qū)，研究修正的Boettcher/Legge非對(duì)稱橋函數(shù)理論［39］，發(fā)展可分段描述的非對(duì)稱壓力與摩擦力系數(shù)關(guān)聯(lián)橋函數(shù)，其中，非對(duì)稱壓力橋函數(shù)［38，39，42］可表示為式（2）：

并且Knm，τ、ΔKnτ1、ΔKnτ2為可調(diào)參數(shù)［40］。

當(dāng)?shù)鼗椒P(guān)于飛行器物形處理［38，43，44］有兩種方法，一種是用解析式分別描述飛行器的各部分，這種方法適用于外形簡(jiǎn)單且易于解析表達(dá)的軸對(duì)稱旋成體飛行器［38，41，44］，但大多數(shù)飛行器外形復(fù)雜且很不規(guī)則，難于進(jìn)行解析表示。于是當(dāng)?shù)鼗?jì)算方法所依賴的物形處理常常使用另一種較為通用而近似的面元處理方法［4，39，43］，該方法將飛行器物面劃分為若干個(gè)面元，利用四邊形或三角形來(lái)逼近復(fù)雜外形，計(jì)算出每個(gè)面元上中心點(diǎn)坐標(biāo)、中心點(diǎn)處的法向及面元面積。

飛行器外形曲面分成若干塊小的曲面，對(duì)于每一個(gè)小曲面，選用一個(gè)小的平面四邊形來(lái)代替。這樣，計(jì)算小曲面上的氣動(dòng)力就轉(zhuǎn)換成計(jì)算面元上的氣動(dòng)力，將這些面元上的氣動(dòng)力加起來(lái)就可以得到整個(gè)飛行器的氣動(dòng)力。當(dāng)面元分得很細(xì)小時(shí)，由此引起的誤差就更小［35，44］。

2.2天宮飛行器計(jì)算模型構(gòu)建及表面非結(jié)構(gòu)網(wǎng)格物形表征技術(shù)

由于天宮飛行器與其左右寬度達(dá)數(shù)十米的太陽(yáng)電池翼帆板及連接裝置構(gòu)成一個(gè)相當(dāng)復(fù)雜而極不規(guī)則巨大板艙組合體。為了解決天宮飛行器及附件復(fù)雜物形表征的困難，研究引入分區(qū)網(wǎng)格與非結(jié)構(gòu)網(wǎng)格生成技術(shù)，以此獲得能準(zhǔn)確表征復(fù)雜物面邊界的非結(jié)構(gòu)網(wǎng)格，以便進(jìn)行復(fù)雜物形表征處理。模型表面采用非結(jié)構(gòu)貼體網(wǎng)格，通過(guò)網(wǎng)格處理得到計(jì)算需要的表面網(wǎng)格，具體過(guò)程見(jiàn)文獻(xiàn)［4］，其中基于三角非結(jié)構(gòu)網(wǎng)格表征天宮飛行器物形，圖1繪出天宮飛行器計(jì)算模型表面三角非結(jié)構(gòu)網(wǎng)格布局。

圖1　天宮飛行器計(jì)算模型表面非結(jié)構(gòu)網(wǎng)格布局Fig.1 Com putational m odel and sur face grid layout of the TG spacecraft

2.3天宮飛行器跨流域氣動(dòng)特性快速算法

根據(jù)稀薄過(guò)渡流流動(dòng)特征，當(dāng)?shù)鼗焖偎惴▽⑾”∵^(guò)渡流區(qū)壓力橋函數(shù)（2）式與摩擦力橋函數(shù)（3）式用非對(duì)稱誤差函數(shù)形式關(guān)聯(lián)，兩式中Knm，p、ΔKnp1、ΔKnp2及Knm，τ、ΔKnτ1、ΔKnτ2應(yīng)根據(jù)風(fēng)洞試驗(yàn)數(shù)據(jù)或數(shù)值計(jì)算結(jié)果進(jìn)行調(diào)試確定，由于天宮飛行器的大尺度復(fù)雜結(jié)構(gòu)與在軌飛行高稀薄自由分子流區(qū)繞流特點(diǎn)，風(fēng)洞實(shí)驗(yàn)條件并不現(xiàn)實(shí)，可依靠高性能計(jì)算機(jī)開(kāi)展DSMC方法［4，13，17］或求解Boltzmann模型方程統(tǒng)一算法［30-34］典型飛行狀態(tài)數(shù)值計(jì)算，調(diào)試選擇與計(jì)算修正上述六個(gè)關(guān)聯(lián)參數(shù)，使當(dāng)?shù)鼗こ逃?jì)算結(jié)果最大程度地滿足實(shí)際計(jì)算精度。由于本文旨在有限時(shí)間搭建天宮飛行器低軌道飛行過(guò)程300～200 km空氣動(dòng)力特性初步計(jì)算平臺(tái)，作為算法驗(yàn)證，橋函數(shù)中當(dāng)?shù)鼗烧{(diào)關(guān)聯(lián)參數(shù)的確定采用250 km高度與天宮飛行器具有同樣尺寸的圓柱體繞流數(shù)值計(jì)算結(jié)果進(jìn)行修正確定。首先利用數(shù)值計(jì)算得到該圓柱體繞流典型狀態(tài)所受到的氣動(dòng)阻力系數(shù)為參考對(duì)象；改變上述六個(gè)關(guān)聯(lián)參數(shù)，使用當(dāng)?shù)鼗焖偎惴ㄟM(jìn)行調(diào)試計(jì)算修正；最終獲得一套針對(duì)天宮飛行器計(jì)算適用的六個(gè)關(guān)聯(lián)參數(shù)的最佳組合，使得經(jīng)關(guān)聯(lián)參數(shù)修正后的氣動(dòng)特性當(dāng)?shù)鼗?jì)算方法能滿足本文關(guān)于天宮飛行器在不同高度、不同馬赫數(shù)、不同迎角與側(cè)滑角的氣動(dòng)力系數(shù)計(jì)算要求。圖2給出了250 km高度0°迎角上述大尺度圓柱體繞流狀態(tài)壓力等值線云圖分布，可看出圓柱體前端表面壓力值遠(yuǎn)大于后表面，氣流完全附著物面強(qiáng)擾動(dòng)，流場(chǎng)呈現(xiàn)出高稀薄自由分子流動(dòng)特征。

圖2　 H=250 km，α=0°大尺度圓柱體繞流壓力等值線分布Fig.2 Pressure distribution contours around largescale cylinder body w ith H=250 km，α=0°

數(shù)值計(jì)算得到0°迎角時(shí)圓柱體模型阻力系數(shù)為3.169025，以此為基礎(chǔ)得到當(dāng)?shù)鼗焖偎惴ㄖ袠蚝瘮?shù)的可調(diào)關(guān)聯(lián)參數(shù)如式（5）～（6）：

以上述六個(gè)當(dāng)?shù)鼗P(guān)聯(lián)參數(shù)為基礎(chǔ)，發(fā)展可用于確定大尺度飛行器在不同高度、馬赫數(shù)、迎角與側(cè)滑角的氣動(dòng)力系數(shù)快速算法。由此計(jì)算得到上述圓柱體繞流阻力系數(shù)3.169013，與數(shù)值模擬結(jié)果相當(dāng)一致，彼此偏差僅0.000379%。表1列出了該狀態(tài)數(shù)值計(jì)算與當(dāng)?shù)鼗焖偎惴ǚ謩e得到4°迎角上述圓柱體模型氣動(dòng)力系數(shù)比較情況，可看出二者吻合很好，偏差范圍0.063%～2.85%，驗(yàn)證了本文提出發(fā)展低軌道飛行過(guò)程大尺度飛行器300 km到200 km空氣動(dòng)力特性快速計(jì)算平臺(tái)可行性。

表1　當(dāng)?shù)鼗焖偎惴ㄓ?jì)算H=250 km，α=4°圓柱體繞流氣動(dòng)力系數(shù)驗(yàn)證比較Table 1 Validation and comparison of aerodynam ic coefficients around the large-scale cylinder w ith H=250 km，α=4°

3　天宮飛行器簡(jiǎn)化外形200 km至300 km空氣動(dòng)力特性計(jì)算分析

3.1算法驗(yàn)證

為了進(jìn)一步驗(yàn)證本文發(fā)展的當(dāng)?shù)鼗焖偎惴Ｐ驼_性以及200～300 km高度飛行的氣動(dòng)力系數(shù)變化規(guī)律，分別采用數(shù)值模擬方法和當(dāng)?shù)鼗焖偎惴ㄓ?jì)算了200 km和300 km上述大尺度圓柱體0°迎角氣動(dòng)力系數(shù)，表2給出了兩類方法計(jì)算0°迎角200 km和300 km阻力系數(shù)Cd比較情況?？煽闯?，兩種計(jì)算方法得到的氣動(dòng)力系數(shù)完全吻合，最大偏差不超過(guò)0.48%，證實(shí)本文發(fā)展的當(dāng)?shù)鼗焖偎惴ㄕ_可靠；在300 km飛行高度阻力系數(shù)比200 km阻力系數(shù)均高出至少7.43%，說(shuō)明近地軌道飛行過(guò)程200～300 km不同高度飛行器阻力系數(shù)是隨飛行高度增加而增大，這反映出傳統(tǒng)工程處理手段認(rèn)為200～300 km氣體分子平均自由程λ∞達(dá)235.2～2595.37m，遠(yuǎn)大于飛行器特征尺寸，來(lái)流氣體處于自由分子流態(tài)，飛行器繞流使用固定不變阻力系數(shù)統(tǒng)一處理方式對(duì)類天宮飛行器兩艙結(jié)構(gòu)大尺度圓柱體繞流會(huì)產(chǎn)生至少7.43%偏差。

3.2天宮飛行器200～300 km氣動(dòng)力特性計(jì)算分析

為使本文當(dāng)?shù)鼗焖偎惴ㄓ?jì)算300～200 km期間天宮飛行器氣動(dòng)特性結(jié)果能盡可能真實(shí)可靠，本文對(duì)任一電池翼帆板平面與天宮飛行器本體主軸夾角0°、22.5°、45°、67.5°、90°確定的計(jì)算模型，進(jìn)行300～200 km不同高度300 km、280 km、250 km、220 km、200 km，不同迎角α= -13.5°，0°，13.5°，22°，38°與側(cè)滑角β=5.5°，2.6°，0°，-2.6°，-4.1°的天宮飛行器氣動(dòng)特性計(jì)算，提供了大量可供軌道控制彈道飛行力學(xué)所用空氣動(dòng)力系數(shù)計(jì)算數(shù)據(jù)。下面僅列兩組典型狀態(tài)。

表2　 200 km和300 km飛行高度0°迎角圓柱體繞流阻力系數(shù)驗(yàn)證比較Table 2 Com parison of drag coefficients around the large-scale cylinder w ith H=200、300 km，α=0°

3.2.1帆板平面與天宮飛行器主軸夾角為0°氣動(dòng)特性計(jì)算

擬定天宮飛行器250 km、高度13.5°迎角、0°側(cè)滑角外形，根據(jù)數(shù)值計(jì)算調(diào)試確定的當(dāng)?shù)鼗瘶蚝瘮?shù)關(guān)聯(lián)參數(shù)為Knm，p=0.9、ΔKnp1=3.5、ΔKnp2=3.5，Knm，τ=0.85、ΔKnτ1=3.5、ΔKnτ2=3.0。圖3繪出天宮飛行器簡(jiǎn)化外形展向水平面與縱向豎直平面內(nèi)壓力等值線分布，可看出繞流流場(chǎng)附著物面流動(dòng)，僅在物面附近出現(xiàn)強(qiáng)擾動(dòng)，逐漸擴(kuò)展到遠(yuǎn)前方，表現(xiàn)出高稀薄自由分子流動(dòng)特征，等值線分布范圍相當(dāng)廣，達(dá)天宮飛行器本體尺寸2～3倍。

圖3　天宮飛行器簡(jiǎn)化外形帆板與主軸夾角0°高度250 km迎角13.5°、側(cè)滑角0°軸對(duì)稱面壓力等值線分布Fig.3 Pressure distribution contours around the sim p lified TG-spacecraft w ithθ=0°，H=250 km，α= 13.5°，β=0°

表3列出了帆板平面與本體主軸夾角θ=0°時(shí)，考慮帆板與天宮飛行器兩艙本體三腳架連接（見(jiàn)圖1）天宮一號(hào)目標(biāo)飛行器計(jì)算模型，采用本文算法計(jì)算得到α=-13.5°、β=5.2°條件下，不同飛行高度氣動(dòng)力系數(shù)，包括軸向力Ca、法向力Cy、側(cè)向力Cz、阻力Cd、升力Cl、側(cè)力Czl，其中“最大差異”指300 km的氣動(dòng)力系數(shù)與200 km相應(yīng)氣動(dòng)力系數(shù)值之間的相對(duì)偏差?？煽闯?，軸向力系數(shù)與法向力系數(shù)最大差異在5%～10%，而側(cè)向力系數(shù)最大差異在7%～13%。但是，由軸向力、法向力和側(cè)向力系數(shù)經(jīng)坐標(biāo)變換得到的風(fēng)軸系下阻力、升力和側(cè)力系數(shù)的最大差異分別達(dá)到了10%、30%和80%。由此可見(jiàn)，在200～300 km的大型復(fù)雜結(jié)構(gòu)航天器飛行過(guò)程，若采用固定統(tǒng)一的氣動(dòng)力系數(shù)進(jìn)行飛行控制會(huì)產(chǎn)生較大誤差。

表3　天宮飛行器帆板與本體0°夾角、-13.5°迎角、5.2°側(cè)滑角200～300 km各高度氣動(dòng)力系數(shù)Table 3 Aerodynam ic coefficients around the TG-spacecraft from 200～300 km withθ=0°，α=-13.5°，β=5.2°

3.2.2帆板平面與天宮飛行器主軸夾角為45°氣動(dòng)特性計(jì)算

對(duì)θ=45°，天宮飛行器與帆板形成計(jì)算模型見(jiàn)圖4（a），根據(jù)數(shù)值計(jì)算確定Knm，p=0.9、ΔKnp1=3.2、ΔKnp2=3.5，Knm，τ=0.65，ΔKnτ1=3.1，ΔKnτ2=3.0。圖4（b）繪出天宮飛行器簡(jiǎn)化外形繞流展向水平面與縱向豎直平面內(nèi)壓力等值線分布云圖，同樣看出繞流流場(chǎng)附著物面流動(dòng)，僅在物面附近出現(xiàn)強(qiáng)擾動(dòng)，并逐漸擴(kuò)展到遠(yuǎn)前方，表現(xiàn)出高稀薄自由分子流動(dòng)特征。

表4介紹了本文方法計(jì)算天宮飛行器θ= 45°、α=0°、β=0°條件下，200～300 km間隔20 km各飛行高度氣動(dòng)力系數(shù)?？煽闯?，300 km與200 km的氣動(dòng)力系數(shù)差別很明顯，隨飛行迎角狀態(tài)不同，這種差異也不同。對(duì)于某些飛行狀態(tài)，出現(xiàn)側(cè)力系數(shù)正負(fù)符號(hào)變化，法向力系數(shù)“最大差異”達(dá)25%以上。

圖4　天宮飛行器簡(jiǎn)化外形帆板與主軸夾角θ= 45°，H=250 km，α=0°，β=0°模型與對(duì)稱面壓力分布Fig.4 Pressu re distribution contours aroundthe simp lified TG-spacecraft withθ= 45°，H=250 km，α=0°，β=0°

表4　θ=45°，α=0°，β=0°飛行過(guò)程200～300 km各高度氣動(dòng)力系數(shù)Table 4 Aerodynam ic coefficients around the TG-spacecraft from 200～300 km w ithθ=45°，α=0°，β =0°

4　結(jié)論

針對(duì)在軌服役天宮飛行器板艙組合體大型復(fù)雜構(gòu)型特點(diǎn)，研究修正的Boettcher/Legge非對(duì)稱橋函數(shù)，發(fā)展了可靠模擬不同高度、不同迎角與側(cè)滑角及不同復(fù)雜外形高超聲速空氣動(dòng)力特性快速計(jì)算方法；利用三角形逼近復(fù)雜物形，研制適于天宮飛行器復(fù)雜物面表征的非結(jié)構(gòu)網(wǎng)格生成技術(shù)；研究通用面元處理法，發(fā)展適于天宮飛行器復(fù)雜物形處理與面元?dú)鈩?dòng)力系數(shù)計(jì)算規(guī)則；建立了適于天宮飛行器近地軌控200～300 km空氣動(dòng)力特性一體化快速算法應(yīng)用研究平臺(tái)，揭示大尺度天宮飛行器空氣動(dòng)力變化規(guī)律，并得到如下結(jié)論：

1）天宮飛行器在300 km到200 km變軌飛行過(guò)程中，氣動(dòng)力系數(shù)隨飛行高度變化幅度較大：體軸系中軸向力系數(shù)與風(fēng)軸系中阻力系數(shù)變化幅度5%～9%；法向力、側(cè)向力等其它氣動(dòng)力系數(shù)變化幅度更大，有的甚至超過(guò)50%，某些狀態(tài)（如帆板平面與本體主軸夾角67.5°，以0°迎角飛行）法向力系數(shù)“最大差異”達(dá)48%，而一些飛行狀態(tài)升力和側(cè)力系數(shù)出現(xiàn)了正負(fù)符號(hào)與方向差別。對(duì)長(zhǎng)期在軌飛行大型航天器，這種誤差累積后果難以控制，須采取防護(hù)措施；

2）對(duì)天宮飛行器這類由太陽(yáng)電池翼與主體構(gòu)成大型復(fù)雜板艙結(jié)構(gòu)航天器低軌道飛行300～200 km氣動(dòng)力系數(shù)隨飛行高度、姿態(tài)與太陽(yáng)電池翼繞連接軸轉(zhuǎn)動(dòng)外形變化而改變，變化幅度與飛行狀態(tài)相關(guān)；不同迎角、側(cè)滑角氣動(dòng)力系數(shù)存在較大差異；同時(shí)，太陽(yáng)能電池翼圍繞主體連接軸旋轉(zhuǎn)對(duì)天宮整器氣動(dòng)力特性也產(chǎn)生較大影響；

3）天宮飛行器近地軌道飛行過(guò)程，氣動(dòng)力系數(shù)不僅發(fā)生變化，且變化幅度與方向不同。如軸向力和阻力系數(shù)隨飛行高度300～200 km變化的“最大差異”在5%～10%，該類氣動(dòng)力系數(shù)本身基數(shù)大，隨高度變化量絕對(duì)值較大（一般在0.18左右）；與之相應(yīng)法向力系數(shù)“最大差異”有時(shí)達(dá)48%，但其隨高度變化量絕對(duì)值為0.06；而某些狀態(tài)，如太陽(yáng)電池翼平面與本體主軸夾角45°，以0°迎角、2.6°側(cè)滑角所受側(cè)力系數(shù)200 km為-0.0001，而300 km變?yōu)?.0032，出現(xiàn)正負(fù)側(cè)力方向的改變，它們之間“最大差異”就會(huì)出現(xiàn)異常嚴(yán)重，而該側(cè)力系數(shù)隨高度變化量絕對(duì)值僅0.0033。這種由太陽(yáng)電池翼旋轉(zhuǎn)所致氣動(dòng)力方向上的變化會(huì)直接危及航天器精細(xì)軌控。

（致謝：本文寫(xiě)作過(guò)程中和李革非研究員、李勰工程師及課題組蔣新宇助理工程師等進(jìn)行了交流討論，部分計(jì)算在總參三部北方計(jì)算中心完成，特此感謝！）

［1］PilinskiM D，Argrow BM，Palo SE.Drag coefficientsofsatellites with concave geometries：comparingmodels and observations［J］.Journal of Spacecraft and Rockets，2011，48（2）：312-325.

［2］Moe K，Bowman B R.The effects of surface composition and treatment on drag coefficients of spherical satellites［C］// AAS/AIAA Astrodynamics Specialist Conference，Americal Astronautical Soc.，Paper 05-258，Springfield，VA，2005.

［3］Moe K，Moe M M，Wallace SD.Improved satellite drag coefficient calculations from orbitalmeasurements of energy accommodation［J］.Journal of spacecraft and rockets，1998，35（3）：266-272.

［4］李志輝，吳俊林，彭傲平.天宮飛行器低軌道飛行控制過(guò)程空氣動(dòng)力建模與計(jì)算研究［R］.中國(guó)國(guó)防科技報(bào)告，綿陽(yáng)：中國(guó)空氣動(dòng)力研究與發(fā)展中心超高速研究所，2013. Li Z H，Wu JL，Peng A P.Aerodynamicsmodeling and calculation during low-orbit flying control of TG-1 target spacecraft［R］.Technical Report，Hypersonic Aerodynamics Institute，China Aerodynamic Research and Development Centre，2013.（in Chinese）

［5］Cercignani C.The Boltzmann Equation and its Applications［M］.Springer Verlag，Berlin，1988：55-58.

［6］應(yīng)純同.氣體輸運(yùn)理論及應(yīng)用［M］.北京：清華大學(xué)出版社，1990：123-136. Ying C T.Gas transport theory and application［M］.Tsinghua Univ.Press，Beijing，1990：123-136.（in Chinese）

［7］Bird G A.Approach to translational equilibrium in a rigid sphere gas［J］.Physics of Fluids（1958-1988），1963，6（10）：1518-1519.

［8］Pham-Van-Diep G，Erwin D，Muntz E P.Nonequilibrium molecular motion in a hypersonic shock wave［J］.Science，1989，245（4918）：624-626.

［9］Carlson A B，Hassan H A.Direct simulation of reentry flows with ionization［R］.AIAA 90-144，1990.

［10］Haas B L，Boyd ID.Models for direct Monte Carlo simulation of coupled vibration-dissociation［J］.Physics of Fluids A：Fluid Dynamics，1993，5（2）：478-489.

［11］Bird G A.Molecular Gas Dynamics and Direct Simulation of Gas Flows［M］.Oxford University Press，1994：78-189

［12］樊菁，沈青.過(guò)渡領(lǐng)域高超聲速圓柱繞流的直接模擬［J］.空氣動(dòng)力學(xué)學(xué)報(bào)，1995，13（4）：405-413. Fan J，Shen Q.Direct simulation of hypersonic transitional flow around circular cylinder，Acta Aerodynamica Sinica，1995，13（4）：405-413.（in Chinese）

［13］李志輝，吳振宇.阿波羅指令艙稀薄氣體動(dòng)力學(xué)特征的蒙特卡羅數(shù)值模擬［J］.空氣動(dòng)力學(xué)學(xué)報(bào)，1996，14（2）：230-233. Li Z H，Wu Z Y.DSMC simulation of hypersonic rarefied flow past Apollo-CM，Acta Aerodynamica Sinica，1996，14（2）：230-233.（in Chinese）

［14］Ivanov M S，Gimelshein S F.Computational hypersonic rarefied flows［J］.Annual Review of Fluid Mechanics，1998，30（1）：469-505.

［15］吳其芬，陳偉芳.高溫稀薄氣體熱化學(xué)非平衡流動(dòng)的DSMC方法［M］.長(zhǎng)沙：國(guó)防科技大學(xué)出版社，1999：165-176. Wu Q F.Chen W.F.DSMC method for high-temperature thermochemical nonequilibrium flow of rarefied gas［M］.National University of Defense Technology Press，ChangSha，1999：165-176.（in Chinese）

［16］Sun Q，F(xiàn)an J，Boyd ID.Improved sampling techniques for the direct simulation Monte Carlomethod［J］.Computers& Fluids，2009，38（2）：475-479.

［17］梁杰，閻超，李志輝，等.稀薄過(guò)渡流區(qū)橫向噴流干擾效應(yīng)數(shù)值模擬研究［J］.空氣動(dòng)力學(xué)學(xué)報(bào)，2013，31（1）：27-33. Liang J，Yan C.，Li Z H，et al.Numerical investigation of lateral jet interaction effects in rarefied transition flow regime，Acta Aerodynamica Sinica，2013，31（1）：27-33.（in Chinese）

［18］Wong B C，Long L N.Direct simulation Monte Carlo（DSMC）on the connection machine［R］.AIAA 92-0564，1992.

［19］Scanlon T J，RoohiE，White C，etal.An open source，parallel DSMC code for rarefied gas flows in arbitrary geometries［J］.Computers&Fluids，2010，39（10）：2078-2089.

［20］李志輝，楊彥廣，梁杰，等.飛船返回艙高空稀薄氣動(dòng)特性研究［C］.全國(guó)載人航天工程氣動(dòng)工作總結(jié)暨關(guān)鍵問(wèn)題研討會(huì)，6月，北京，2004. Li Z H.Yang Y G.，Liang J，etal.Study of rarefied aerodynamics for spacecraft capsule in high-altitude，Symposium of front key problems for national manned spaceflight，Beijing，June，2004.（in Chinese）

［21］Yoon S，Gnoffo P A，White JA，et al.Computational challenges in hypersonic flow simulations［R］.AIAA 2007-4265，2007.

［22］Bhatnagar P L，Gross E P，Krook M.A model for collision processes in gases.I.Small amplitude processes in charged and neutral one-component systems［J］.Physical review，1954，94（3）：511-525.

［23］Holway Jr L H.New statistical models for kinetic theory：methods of construction［J］.Physics of Fluids（1958-1988），1966，9（9）：1658-1673.

［24］Shakhov E M.Generalization of the Krook kinetic relaxation equation［J］.Fluid Dynamics，1968，3（5）：95-96.

［25］Yang JY，Huang JC.Rarefied flow computations using nonlinearmodel Boltzmann equations［J］.Journal of Computational Physics，1995，120（2）：323-339.

［26］Titarev V A，Shakhov E M.Heat transfer and evaporation from a p lane surface into a half-space upon a sudden increase in body temperature［J］.Fluid dynamics，2002，37（1）：126-137.

［27］Mieussens L.Discrete-velocitymodels and numerical schemes for the Boltzmann-BGK equation in plane and axisymmetric geometries［J］.Journal of Computational Physics，2000，162（2）：429-466.

［28］李志輝，張涵信.稀薄流到連續(xù)流的氣體運(yùn)動(dòng)論統(tǒng)一數(shù)值算法初步研究［J］.空氣動(dòng)力學(xué)學(xué)報(bào)，2000，18（3）：255-263. Li ZH.，Zhang H X.Study on gas-kinetic algorithm for flows from rarefied transition to continuum，Acta Aerodynamica Sinica，2000，18（3）：255-263.（in Chinese）

［29］Li Z，Zhang H.Study on gas kinetic algorithm for flows from rarefied transition to continuum［C］//Proc.of 22ndInternational Symposium on Rarefied Gas Dynamics，Sydney，Australia，Jul.9-14，2000：628-636.

［30］李志輝.從稀薄流到連續(xù)流的氣體運(yùn)動(dòng)論統(tǒng)一數(shù)值算法研究［D］.綿陽(yáng)：中國(guó)空氣動(dòng)力研究與發(fā)展中心研究生部，2001. Li Z H.Study on gas kinetic unified algorithm for flows from rarefied transition to continuum，Ph.D thesis，China Aerodynamics Research and Development Center，Jan.，2001.（in Chinese）

［31］Li ZH，Zhang H X.Study on gas kinetic unified algorithm for flows from rarefied transition to continuum［J］.Journal of Computational Physics，2004，193（2）：708-738.

［32］李志輝，張涵信.基于Boltzmann模型方程的氣體運(yùn)動(dòng)論統(tǒng)一算法研究［J］.力學(xué)進(jìn)展，2005，35（4）：559-576. Li Z H.，Zhang H X.Study on gas kinetic numerical algorithm using Boltzmann model equation.Advances in Mechanics 2005；35（4）：559-576.（in Chinese）

［33］Li Z，Zhang H.Gas-kinetic numericalmethod for solvingmesoscopic velocity distribution function equation［J］.Acta Mechanica Sinica，2007，23（2）：121-132.

［34］Li Z H，Zhang H X.Gas-kinetic numerical studies of threedimensional complex flows on spacecraft re-entry［J］.Journal of Computational Physics，2009，228（4）：1116-1138.

［35］李志輝，梁杰，李四新，等.箔條云跨流域整體氣動(dòng)特性計(jì)算研究［J］.空氣動(dòng)力學(xué)學(xué)報(bào)，2011，29（01）：59-67. Li Z H.，Liang J，Li SX，et al.Computing study on holistic aerodynamics of chaff cloud covering various flow regimes，Acta Aerodynamica Sinica，2011，29（1）：59-67.（in Chinese）

［36］Li Z H.Gas-kinetic unified algorithm for re-entering complex flows covering various flow regimes by solving boltzmann model equation［M］//Hall J.Advances in Spacecraft Technologies.Croatia：InTech，2011：273-332.

［37］李志輝，彭傲平，吳俊林，等.稀薄氣體自由分子流到連續(xù)流跨流域氣動(dòng)力熱繞流統(tǒng)一算法研究［J］.載人航天，2013，19（2）：81-91. Li ZH.，Peng A P，Wu JL，etal.Gas-kinetic unified algorithm for aerothermodynamics from free molecular regime to continuum，Manned Spaceflight，2013，19（2）：81-91.（in Chinese）

［38］Gentry A E，Smyth D N，OliverW R.The Mark IV supersonic-hypersonic arbitrary-body program，Vol.11：Program Formulation［R］.AFFDL TR-73-159，1973.

［39］Kashkovskyt A V.Computational tools for rarefied aerodynamics［C］//Rarefied Gas Dynamics：Technical Papers from the Proceedings of the Eighteenth International Symposium on Rarefied Gas Dynamics，University of British Columbia，Vancouver，British Columbia，Canada，1992：115-126.

［40］Mosanov SV，F(xiàn)reedlender O G，Nikiforov A P，et al.Experimental determination ofmomentum transfer coefficients in hypersonic freemolecular flow and distribution function recovery of reflected molecules［C］//Proceedings of the XIII International Symposium Rarefied Gas Dynamics，Vol.1，Plenum，New York，1985：669-676.

［41］謝硯儒，李志輝.尖鈍錐誘餌再入阻力特性工程計(jì)算研究［R］.綿陽(yáng)：中國(guó)空氣動(dòng)力研究與發(fā)展中心，中國(guó)國(guó)防科技報(bào)告，1992. Xie Y R.，Li ZH.，Study ofengineering calculation of re-entry drag characteristics for sharp-cone body［R］.Technical Report，Hypersonic Aerodynamics Institute，China Aerodynamic Research and DevelopmentCentre，1992.（in Chinese）

［42］李志輝.箔條及箔條云整體氣動(dòng)特性建模方法研究［R］.綿陽(yáng)：中國(guó)空氣動(dòng)力研究與發(fā)展中心，中國(guó)國(guó)防科技報(bào)告，2005. Li Z H.，Modeling study on holistic aerodynamic characteristic of foil and foil-cloud cluster［R］.Technical Report，Hypersonic Aerodynamics Institute，China Aerodynamic Research and Development Centre，2005.（in Chinese）

［43］謝硯儒，劉濤.彈頭稀薄氣動(dòng)特性的工程計(jì)算方法研究［R］.綿陽(yáng)：中國(guó)空氣動(dòng)力研究與發(fā)展中心，中國(guó)國(guó)防科技報(bào)告，1990.Xie Y R.，Liu T.，Study ofengineering calculation method of rarefied aerodynamics for bullet［R］.Technical Report，Hypersonic Aerodynamics Institute，China Aerodynamic Research and Development Centre，1990.（in Chinese）

［44］粱杰，李志輝.目標(biāo)稀薄氣體特性工程計(jì)算［R］.綿陽(yáng)：中國(guó)空氣動(dòng)力研究與發(fā)展中心，國(guó)防科技報(bào)告，1999. Li Z H.，Liang JEngeering calculation of target rarefied aerodynamic characteristics［R］.Technical Report，Hypersonic Aerodynamics Institute，China Aerodynamic Research and Development Centre，1999.（in Chinese）

Unified Modeling and Calculation of Aerodynam ics Characteristics during Low-O rbit Flying Control of the TG Vehicle

LIZhihui1，2，WU Junlin1，PENG Aoping1，TANG Geshi3
（1.Hypervelocity Aerodynamics Institute，China Aerodynamics Research and Development Center，Mianyang 621000，China；2.National Laboratory for Computational Fluid Dynamics，Beijing 100191，China；3.Science&Technology on Aerospace Flight Dynamics Laboratory，Beijing Aerospace Control Center，Beijing100094，China）

The surface analyticalmethod，considering surface shielding effect of the complex structure and the modified Boettcher/Legge nonsymmetric bridge correction function，was proposed to computationally model aerodynamic characteristics of the large-scale Tiangong spacecraft of irregular plate-capsule assembly.The complex configuration processing and computing rules of surface element aerodynamic coefficientswere setup by developing a general triangle elementapproximation for complex shapes.A unified fast algorithm and computing software for aerodynamic characteristics of large-scale complex spacecraft structures were developed，by computing correction of local correlation parameters during Tiangong spacecraft’s low-earth orbit flight control，in which the DSMC method and the gas-kinetic unified algorithm solving the Boltzmann model equation were applied.Itis indicated that，by computing the aerodynamics of the Tiangong spacecraft in 300～200 km altitude-orbit flight process，with various flying heights，various angles of attack and various angles between the panel and the principal axis，the aerodynamic coefficients varied remarkably，with a range of 8%～50%，with the altitude change，and the atmospheric drag in low-earth orbit flight control was the key factor of the spacecraft orbit prediction accuracy.It is validated that，if the fixed and constant unity aerodynamic coefficient is used for long-term orbit flight of a large spacecraft，the error accumulation will be huge，and protectivemeasures are needed.

TG-1 target spacecraft；aerodynamic characteristic during low-orbit control；local bridge function；correlation parameter evaluation；engineering calculation；DSMCmethod；Gas-Kinetic U-nified Algorithm（GKUA）

V411.4

A

1674-5825（2015）02-0106-09

2014-09-12；

2015-02-26

國(guó)家973計(jì)劃（2014CB744100）；國(guó)家自然科學(xué)基金（11325212、91016027）；國(guó)防基礎(chǔ)科研基金（51313030104）作者簡(jiǎn)介：李志輝（1968-），男，博士，研究員，研究方向?yàn)榭缌饔蚩諝鈩?dòng)力學(xué)。E-mail：zhli0097@x263.net