国产日韩欧美一区二区三区三州_亚洲少妇熟女av_久久久久亚洲av国产精品_波多野结衣网站一区二区_亚洲欧美色片在线91_国产亚洲精品精品国产优播av_日本一区二区三区波多野结衣 _久久国产av不卡

?

固體火箭發(fā)動(dòng)機(jī)噴管分離流動(dòng)流固耦合數(shù)值仿真①

2012-07-09 09:12:06吳朋朋楊月誠(chéng)高雙武張昆鵬
固體火箭技術(shù) 2012年3期
關(guān)鍵詞:氣流流場(chǎng)流動(dòng)

吳朋朋,楊月誠(chéng),高雙武,張昆鵬

(第二炮兵工程大學(xué)二系201室,西安 710025)

固體火箭發(fā)動(dòng)機(jī)噴管分離流動(dòng)流固耦合數(shù)值仿真①

吳朋朋,楊月誠(chéng),高雙武,張昆鵬

(第二炮兵工程大學(xué)二系201室,西安 710025)

針對(duì)固體火箭發(fā)動(dòng)機(jī)大膨脹比噴管出現(xiàn)的分離流動(dòng),采用數(shù)值仿真方法進(jìn)行分析,并與試驗(yàn)進(jìn)行對(duì)比。通過集成軟件平臺(tái)MpCCI,連接計(jì)算流體動(dòng)力學(xué)軟件FLUENT和有限元軟件ABAQUS,對(duì)燃?xì)饬鲃?dòng)與噴管結(jié)構(gòu)運(yùn)動(dòng)變形進(jìn)行了耦合計(jì)算。耦合計(jì)算結(jié)果發(fā)現(xiàn),此大膨脹比噴管發(fā)生氣流分離,且分離處斜激波后的氣流溫度與壓力變化較大,采用流固耦合數(shù)值方法能體現(xiàn)噴管的結(jié)構(gòu)變形,從而更準(zhǔn)確地反映噴管與燃?xì)饬飨嗷ビ绊懙恼鎸?shí)環(huán)境。耦合計(jì)算結(jié)果與試驗(yàn)進(jìn)行對(duì)比得出,耦合計(jì)算得到的分離位置能很好地?cái)M合實(shí)驗(yàn)測(cè)得的氣流分離位置,說明了流固耦合數(shù)值方法的有效性,為更深入研究大膨脹比噴管分離流動(dòng)現(xiàn)象提供了支撐。

固體火箭發(fā)動(dòng)機(jī);大膨脹比噴管;分離流動(dòng);流固耦合分析

0 引言

隨著航天科學(xué)技術(shù)的飛速發(fā)展,固體火箭發(fā)動(dòng)機(jī)的研究工作開始突飛猛進(jìn)。如美國(guó)于2010年8月31日對(duì)DM-2固體火箭發(fā)動(dòng)機(jī)進(jìn)行了測(cè)試,DM-2發(fā)動(dòng)機(jī)可產(chǎn)生360萬磅(約1 630 t)推進(jìn)力,是迄今為止用于飛行的體積最大且功率最大的固體火箭發(fā)動(dòng)機(jī),可能被用在未來的重型運(yùn)載火箭上。之所以有這么大的進(jìn)步,非常重要的一點(diǎn)就是對(duì)噴管的改進(jìn)。噴管是固體火箭發(fā)動(dòng)機(jī)的基本部件之一,在許多情況下,它決定了發(fā)動(dòng)機(jī)的外形和能量質(zhì)量完善程度。改進(jìn)噴管是提高火箭發(fā)動(dòng)機(jī)性能重要途徑之一。氣動(dòng)性能設(shè)計(jì)、結(jié)構(gòu)強(qiáng)度設(shè)計(jì)都是噴管設(shè)計(jì)的重要內(nèi)容。

隨著推進(jìn)技術(shù)的發(fā)展,運(yùn)載火箭的助推級(jí)或第一級(jí)發(fā)動(dòng)機(jī)正采用越來越大面積比的噴管,以提高高空性能。但大面積比噴管在地面試車以及發(fā)動(dòng)機(jī)的啟動(dòng)和關(guān)機(jī)過程中,都會(huì)產(chǎn)生分離流動(dòng)現(xiàn)象。使用大面積比噴管來加速燃?xì)獾陌l(fā)動(dòng)機(jī),像美國(guó)的J2發(fā)動(dòng)機(jī)、航天飛機(jī)主發(fā)動(dòng)機(jī)(SSME)、俄羅斯的RD-0120發(fā)動(dòng)機(jī)、歐洲的火神發(fā)動(dòng)機(jī)和日本的LE-7A[1-2]等發(fā)動(dòng)機(jī)的研制過程中,均遇到了分離流動(dòng)現(xiàn)象。國(guó)內(nèi)外學(xué)者對(duì)分離流已開展了許多研究,如文獻(xiàn)[3-5]。

氣流分離往往會(huì)對(duì)噴管造成危害,因?yàn)闅饬鞣蛛x的不對(duì)稱性,可能造成較嚴(yán)重的噴管側(cè)向載荷問題,導(dǎo)致噴管氣動(dòng)性能下降。因此,氣流分離的精確預(yù)測(cè)對(duì)于合理設(shè)計(jì)噴管是非常關(guān)鍵的。在很大的空氣彈性變形情況下,由于流動(dòng)和結(jié)構(gòu)相互作用,將會(huì)引起側(cè)向載荷巨大的增長(zhǎng)[6]。以往的仿真計(jì)算大多單獨(dú)采用計(jì)算流體動(dòng)力學(xué)軟件FLUENT,但沒有反映噴管結(jié)構(gòu)與燃?xì)饬飨嗷ビ绊懙恼鎸?shí)環(huán)境,因此缺乏足夠的準(zhǔn)確性。隨著數(shù)值算法的發(fā)展,為得到高質(zhì)量的數(shù)值仿真結(jié)果,采用流固耦合方法對(duì)其進(jìn)行耦合求解逐漸成為趨勢(shì)。該方法可揭示內(nèi)流場(chǎng)、結(jié)構(gòu)相互影響及規(guī)律,更真實(shí)反映發(fā)動(dòng)機(jī)工作狀態(tài),提高固體火箭發(fā)動(dòng)機(jī)設(shè)計(jì)水平[7]。

文中運(yùn)用MpCCI(Mesh-based parallel Code Coupling Interface)耦合器作為計(jì)算流體動(dòng)力學(xué)軟件FLUENT和有限元分析軟件ABAQUS的數(shù)據(jù)交換平臺(tái),對(duì)燃?xì)饬鲃?dòng)與噴管結(jié)構(gòu)運(yùn)動(dòng)變形進(jìn)行了耦合計(jì)算。通過與某固體火箭發(fā)動(dòng)機(jī)噴管實(shí)驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比,證明了流固耦合數(shù)值方法的有效性。

1 流固耦合數(shù)理模型

1.1 物理模型

基于文中主要研究大面積比噴管的氣流分離現(xiàn)象,因此把大面積比噴管和外場(chǎng)作為研究對(duì)象。文中以某固體火箭發(fā)動(dòng)機(jī)噴管為研究背景[3]??紤]該模型的對(duì)稱性,采用二維軸對(duì)稱模型。模型如圖1所示。

圖1 噴管結(jié)構(gòu)和流場(chǎng)區(qū)域Fig.1 The nozzle and flow field

該噴管模型總長(zhǎng)0.162 m,噴管面積比 ε=55.2。噴管后段為外場(chǎng)。為簡(jiǎn)化計(jì)算,文中結(jié)構(gòu)計(jì)算僅考慮燃?xì)饬鲃?dòng)對(duì)噴管結(jié)構(gòu)變形的影響。

1.2 計(jì)算模型

如圖1所示,流場(chǎng)取壓力入口邊界ab;噴管內(nèi)壁bd,外壁de;壓力遠(yuǎn)場(chǎng)邊界ef、fg、gh、hi;壓力出口邊界ij和對(duì)稱軸jk、ka圍成的封閉區(qū)域作為控制體。由于噴管結(jié)構(gòu)變形,需采用動(dòng)網(wǎng)格技術(shù)來模擬,所以接近噴管壁面bd、de的流場(chǎng)局部區(qū)域采用了三角形非結(jié)構(gòu)化網(wǎng)格,流場(chǎng)其他區(qū)域均采用四邊形網(wǎng)格。

1.2.1 假設(shè)條件

文中分析重點(diǎn)是噴管結(jié)構(gòu)與燃?xì)庵g的流固耦合作用。為便于分析,引入假設(shè):(1)燃?xì)饬鳛槔硐霘怏w;(2)不考慮噴管結(jié)構(gòu)溫度變化所造成的熱應(yīng)力。

禮者,以財(cái)物為用,以貴賤為文,以多少為異,以隆殺為要。文理繁,情用省,是禮之隆也;文理省,情用繁,是禮之殺也;文理、情用相為內(nèi)外表里,并行而襍,是禮之中流也。故君子上致其隆,下盡其殺,而中處其中。

1.2.2 結(jié)構(gòu)參數(shù)

噴管結(jié)構(gòu)的主要性能參數(shù):密度ρ=7 850 kg/m3,泊松比 γ =0.33,彈性模量E=2.0 ×1011Pa。

1.2.3 初始條件

燃燒室壓強(qiáng)為6.0×106Pa。噴管流場(chǎng)取海平面的壓強(qiáng)與溫度作為流固耦合初始條件。噴管頭部壁面bc為絕熱固壁邊界;噴管結(jié)構(gòu)表面bd、de為耦合邊界。壓力遠(yuǎn)場(chǎng)邊界條件和壓力出口邊界條件取海平面的壓強(qiáng)和溫度。

1.2.4 模型驗(yàn)證

數(shù)值模擬有自己的學(xué)科哲學(xué),由于分離流很難算準(zhǔn),這是當(dāng)前學(xué)術(shù)界公認(rèn)的問題,分離流對(duì)網(wǎng)格密度有很強(qiáng)的依賴性。因此,模型驗(yàn)證工作較重要。為選取合適的計(jì)算模型網(wǎng)格,文中對(duì)4種不同計(jì)算模型網(wǎng)格數(shù),通過單獨(dú)使用計(jì)算流體動(dòng)力學(xué)軟件FLUENT進(jìn)行網(wǎng)格無關(guān)性驗(yàn)證,范圍從網(wǎng)格6 826(網(wǎng)格A)到58 581個(gè)網(wǎng)格(網(wǎng)格D)。由于文中主要研究噴管分離流動(dòng),故通過對(duì)計(jì)算模型網(wǎng)格是否能準(zhǔn)確捕捉分離點(diǎn)及馬赫盤位置(x軸坐標(biāo))進(jìn)行對(duì)比。計(jì)算模型網(wǎng)格無關(guān)性驗(yàn)證在表1中進(jìn)行了顯示。

表1 燃燒室壓強(qiáng)為6 MPa時(shí)計(jì)算模型網(wǎng)格數(shù)Table 1 Grid distributions under 6 MPa combustion pressure

分析發(fā)現(xiàn),網(wǎng)格C與網(wǎng)格D在分離點(diǎn)位置的差異為4%,在馬赫盤位置的差異為0.8%,基于上述分析和考慮到計(jì)算成本,網(wǎng)格數(shù)為38 105(網(wǎng)格C)的計(jì)算模型網(wǎng)格已可較好地用于文中研究?jī)?nèi)容,故采用網(wǎng)格C作為發(fā)動(dòng)機(jī)內(nèi)流場(chǎng)計(jì)算模型網(wǎng)格。

1.2.5 耦合方法

在整個(gè)耦合過程中,采取時(shí)間步為10-6s。發(fā)動(dòng)機(jī)內(nèi)流場(chǎng)由FLUENT軟件計(jì)算,采用耦合隱式求解方法,湍流模型采用RNGk-ε模型,對(duì)流項(xiàng)和粘性項(xiàng)的離散采用二階迎風(fēng)格式,整個(gè)計(jì)算具有二階精度。噴管結(jié)構(gòu)運(yùn)動(dòng)變形計(jì)算采用ABAQUS隱式求解模塊——ABAQUS/Standard模塊。雙方在耦合區(qū)域部分的網(wǎng)格可不匹配,而網(wǎng)格數(shù)據(jù)之間的轉(zhuǎn)換是通過MPCCI的插值來實(shí)現(xiàn),從而將FLUENT軟件和ABAQUS軟件每一個(gè)迭代步的計(jì)算結(jié)果進(jìn)行數(shù)據(jù)交換。

2 結(jié)果分析

2.1 流場(chǎng)分析

在耦合過程中,通過計(jì)算流體動(dòng)力學(xué)軟件FLUENT計(jì)算結(jié)果,顯示了噴管不同時(shí)刻的流場(chǎng)變化。

圖2為在6.0×106Pa的燃燒室壓強(qiáng)條件下,噴管內(nèi)不同時(shí)刻的速度云圖分布。由圖2可看出,氣流向下游傳播,在0~1 ms隨著流動(dòng)的發(fā)展,馬赫盤尺寸不斷增長(zhǎng),產(chǎn)生了不利于噴管流動(dòng)的壓強(qiáng)梯度。由0~9 ms的流動(dòng)分析可知,在2 ms時(shí)噴管流動(dòng)達(dá)到穩(wěn)定狀態(tài),噴管內(nèi)部有激波,但激波未達(dá)到喉部,說明工作在6.0×106Pa的燃燒室壓強(qiáng)下,該噴管處在超臨界流動(dòng)狀態(tài),噴管流場(chǎng)呈現(xiàn)馬赫盤激波模態(tài)[8](Mach disk shock pattern)。受噴管內(nèi)馬赫盤的強(qiáng)烈阻擋,氣流有繞開其流動(dòng)的趨勢(shì),從而使馬赫盤下游形成了低速區(qū),低速區(qū)內(nèi)的氣流溫度較高,這在噴管溫度云圖上反映很明顯。

圖3為噴管壁面9 ms時(shí)刻壓強(qiáng)曲線圖。此時(shí),噴管內(nèi)流場(chǎng)處于自由激波分離模態(tài)。壁面壓強(qiáng)首先隨著氣流膨脹逐漸降低,在分離點(diǎn)由于分離斜激波的作用,壁面壓強(qiáng)迅速升高到一個(gè)相對(duì)穩(wěn)定的平臺(tái)壓強(qiáng)。分離點(diǎn)以后,外界氣流進(jìn)入分離區(qū)形成回流,壁面壓強(qiáng)逐漸升高到接近環(huán)境壓強(qiáng)。氣流分離后沒有再附到壁面,區(qū)域產(chǎn)生不穩(wěn)定的壓強(qiáng)差,勢(shì)必造成噴管流場(chǎng)的壓強(qiáng)震蕩,導(dǎo)致噴管產(chǎn)生側(cè)向載荷。

圖2 不同時(shí)刻噴管內(nèi)流場(chǎng)速度云圖Fig.2 Velocity magnitude contour of the nozzle at different times

圖4為噴管流場(chǎng)在9 ms時(shí)的溫度云圖分布,圖5為分離點(diǎn)周圍區(qū)域壁面在9 ms時(shí)的溫度變化曲線。由圖4和圖5可看出,流場(chǎng)溫度沿軸線逐漸減小,在分離點(diǎn)附近突然升高,然后又快速下降。在氣流分離處存在明顯的溫度變化梯度,這是由于氣流分離處形成斜激波,氣流經(jīng)過該激波減速增壓,導(dǎo)致溫度升高。但在分離斜激波后為環(huán)境空氣,空氣溫度較低,導(dǎo)致溫度梯度較大。在馬赫盤附近區(qū)域,可看出由于馬赫盤對(duì)流動(dòng)的阻擋作用,上下游存在很大的溫度差。

圖3 9 ms時(shí)噴管分離點(diǎn)壁面區(qū)域壓力曲線圖Fig.3 Pressure graph of the nozzle on the wall of seperation point at 9 ms

圖4 9 ms時(shí)噴管的溫度云圖分布Fig.4 Temperature contour of the nozzle at 9 ms

圖5 9 ms時(shí)分離點(diǎn)壁面區(qū)域溫度曲線圖Fig.5 Temperature graph of the nozzle on the wall of seperation point at 9 ms

2.2 結(jié)構(gòu)分析

在耦合過程中,通過有限元軟件ABAQUS計(jì)算結(jié)果顯示了噴管結(jié)構(gòu)在不同時(shí)刻的應(yīng)力、應(yīng)變分布圖。圖6為9 ms時(shí)噴管結(jié)構(gòu)應(yīng)力圖,圖7為噴管外壁de中點(diǎn)(見圖1)的位移-時(shí)間變化曲線。

由圖6可發(fā)現(xiàn),噴管結(jié)構(gòu)在分離點(diǎn)附近的應(yīng)力變化很大,將引起結(jié)構(gòu)震蕩。由圖7可發(fā)現(xiàn),曲線也較好反映了流場(chǎng)參數(shù)變化對(duì)結(jié)構(gòu)造成的影響。由于縮比例模型和全尺寸噴管之間的兩個(gè)主要差別是氣體特性和噴管尺寸。對(duì)于此噴管,結(jié)構(gòu)尺寸較小,而噴管厚度大,結(jié)構(gòu)變形相應(yīng)較小。噴管結(jié)構(gòu)的微小變形對(duì)噴管壁壓也相當(dāng)敏感[9]。對(duì)于大型噴管,如美國(guó)的J2S發(fā)動(dòng)機(jī)噴管,以及軍用固體火箭發(fā)動(dòng)機(jī)噴管需追求輕質(zhì)殼體,均會(huì)導(dǎo)致結(jié)構(gòu)變形較大。結(jié)構(gòu)發(fā)生變化又必定對(duì)流場(chǎng)造成影響。這也是造成分離流不穩(wěn)定性的一個(gè)重要因素,可作為解釋在燃燒室壓強(qiáng)不變時(shí),分離位置不斷地移動(dòng)的一個(gè)原因。

2.3 仿真結(jié)果和試驗(yàn)結(jié)果對(duì)比

文中對(duì)試驗(yàn)數(shù)據(jù)[3]與流固耦合數(shù)值計(jì)算結(jié)果(t=9 ms時(shí))進(jìn)行了比較,圖8顯示了燃燒室壓強(qiáng)為6.0×106Pa時(shí),試驗(yàn)和數(shù)值計(jì)算得到的分離點(diǎn)位置基本吻合,驗(yàn)證了流固耦合數(shù)值計(jì)算的準(zhǔn)確性。

圖6 9 ms時(shí)噴管結(jié)構(gòu)應(yīng)力圖Fig.6 Stress contour of the nozzle at 9 ms

圖7 噴管外壁中點(diǎn)的位移-時(shí)間變化曲線Fig.7 Displacement of the center of the outer wall

圖8 仿真結(jié)果和試驗(yàn)結(jié)果對(duì)比Fig.8 Simulation result compared with the experiment data

3 結(jié)論

(1)對(duì)大膨脹比噴管流動(dòng),采用了二維軸對(duì)稱流固耦合數(shù)值模擬,展現(xiàn)了噴管流場(chǎng)與結(jié)構(gòu)相互影響,并對(duì)噴管結(jié)構(gòu)的應(yīng)力、位移進(jìn)行分析。經(jīng)分析表明,由于分離點(diǎn)附近應(yīng)力、速度、溫度變化都較大,噴管在這些復(fù)雜載荷作用下,導(dǎo)致在分離點(diǎn)附近各項(xiàng)參數(shù)變化劇烈。因此,對(duì)大膨脹比的噴管設(shè)計(jì)時(shí),應(yīng)考慮噴管擴(kuò)張段的結(jié)構(gòu)防護(hù)。

(2)經(jīng)分析顯示,MpCCI軟件能很好地將計(jì)算流體動(dòng)力學(xué)軟件FLUENT和有限元軟件ABUQUS聯(lián)合進(jìn)行流固耦合計(jì)算,同步對(duì)噴管流場(chǎng)和噴管結(jié)構(gòu)進(jìn)行分析。由于其功能較好,這種方法可被用來構(gòu)造多種工況下發(fā)動(dòng)機(jī)的失效機(jī)理,尤其是準(zhǔn)確呈現(xiàn)發(fā)動(dòng)機(jī)的工作過程,為合理設(shè)計(jì)發(fā)動(dòng)機(jī)并延壽提供準(zhǔn)確的技術(shù)支持,節(jié)約大量研制經(jīng)費(fèi)。

文中流固耦合數(shù)值模擬方法可在下階段用于對(duì)實(shí)際應(yīng)用進(jìn)行三維數(shù)值仿真。需指出的是由分離流產(chǎn)生的側(cè)向載荷引起噴管破壞,當(dāng)前國(guó)際上還沒有實(shí)質(zhì)性的解決方法,這種不足導(dǎo)致了發(fā)動(dòng)機(jī)壽命的減少和質(zhì)量的增加。分離流狀態(tài)下的流固耦合分析是非常復(fù)雜的問題,文中方法僅重點(diǎn)考慮了在一定條件下噴管流場(chǎng)與結(jié)構(gòu)的相互影響,后續(xù)還應(yīng)通過與試驗(yàn)對(duì)比不斷修改完善分析模型,為工程實(shí)際應(yīng)用打下基礎(chǔ)。

[1]Hadjadj A,Onofri M.Nozzle flow separation[J].Shock Waves,2009,19(3):163-169.

[2]Yasuhide W,Norio S.LE-7A engine nozzle flowseparation phenomenon and the possibility of RSS suppression by the step inside the nozzle[R].AIAA 2004-4014.

[3]胡海峰,鮑福廷,王藝杰,等.噴管分離流流動(dòng)-熱結(jié)構(gòu)順序耦合數(shù)值模擬及試驗(yàn)研究[J].宇航學(xué)報(bào),2011,32(7):1534-1541.

[4]Francesco Nasuti,Marcello Onofri.Numerical analysis of flow separation structures in rocket nozzles[R].AIAA 2007-5473.

[5]Ralf H Stark,Chloe Genin.Experimental study on rocket nozzle side load reduction[R].AIAA 2011-389.

[6]Zhang S J,F(xiàn)uchiwaki T.Aeroelastic coupling and side loads in rocket nozzles[R].AIAA 2008-4064.

[7]Dennis C W,Belfield J.Towards identifying rocket motor failure modes using coumputational fluid dynamic and finite element codes[R].AIAA 2008-5209.

[8]Frey M,Hagemann G.Status of flow separation prediction in rocket nozzles[C]//34th AIAA/ASME/SAE/ASEE Joint Propulsion Conference Exhibit.Cleveland,Oh,USA,July 13-15,1998.

[9]Nave L H,Coffey G A,孫國(guó)慶.大面積比火箭發(fā)動(dòng)機(jī)噴管的側(cè)向載荷[J].導(dǎo)彈與航天運(yùn)載技術(shù),1980(6).

Fluid-structure coupled numerical simulation of flow separation in SRM nozzle

WU Peng-peng,YANG Yue-cheng,GAO Shuang-wu,ZHANG Kun-peng
(Second Artillery Engineering University No.201 Staff Room,Xi'an 710025,China)

The flow separation in the overexpanded nozzle was studied by using numerical simulation method and compared with experiment result.The MpCCI software is used to link the FLUENT CFD code and the ABUQUS FE code to analyze gas flow and the nozzle deformation.The results show that the gas in the nozzle is separated and the dramatic temperature and pressure alteration occur in the separation zone.The location of the separation point of fluid-structure coupled simulation was similar to the experimental data and can display the nozzle deformation,which support the accuracy of the method of numerical simulation.The simulation provides the base for further study.

solid rocket motor;over-expanded rocket nozzle;flow separation;fluid-structure coupled analysis

V435

A

1006-2793(2012)03-0344-04

2011-12-05;

2012-01-09。

吳朋朋(1982—),男,碩士生,主要研究航空宇航推進(jìn)理論與工程。E-mail:wupeng_ch@163.com

(編輯:崔賢彬)

猜你喜歡
氣流流場(chǎng)流動(dòng)
氣流的威力
大型空冷汽輪發(fā)電機(jī)轉(zhuǎn)子三維流場(chǎng)計(jì)算
流動(dòng)的光
流動(dòng)的畫
轉(zhuǎn)杯紡排雜區(qū)流場(chǎng)與排雜性能
基于HYCOM的斯里蘭卡南部海域溫、鹽、流場(chǎng)統(tǒng)計(jì)分析
為什么海水會(huì)流動(dòng)
固體運(yùn)載火箭變軌發(fā)動(dòng)機(jī)噴管氣流分離研究
飛片下的空氣形成的“超強(qiáng)高速氣流刀”
基于瞬態(tài)流場(chǎng)計(jì)算的滑動(dòng)軸承靜平衡位置求解
湟中县| 大同县| 金湖县| 肥东县| 卓尼县| 安西县| 攀枝花市| 夏河县| 尼玛县| 封丘县| 都昌县| 乌恰县| 平阳县| 余江县| 石台县| 泰兴市| 揭东县| 海丰县| 青神县| 绵阳市| 尼玛县| 大竹县| 乐安县| 北京市| 六盘水市| 大渡口区| 罗江县| 辽宁省| 富源县| 教育| 仁怀市| 将乐县| 墨玉县| 肥东县| 黔西县| 抚远县| 遵义县| 泸州市| 秭归县| 东莞市| 遂平县|