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

?

不飽和聚酯包覆層流動及澆注的仿真模擬

2022-11-17 08:39張豫魯劉奔奔陳國輝曹貝貝何吉宇李向梅楊榮杰
含能材料 2022年11期
關鍵詞:入口流體速率

張豫魯,劉奔奔,陳國輝,曹貝貝,何吉宇,李向梅,楊榮杰

(1. 北京理工大學材料學院,北京 100081;2. 北京航空航天大學前沿科學技術創(chuàng)新研究院,北京 100191;3. 西安近代化學研究所,陜西 西安 710065;4. 北京理工大學化學與化工學院,北京 100081)

0 引言

包覆層用于對固體火箭推進劑的藥柱進行包覆,可起到耐燒蝕、限燃、保證固體火箭推進劑正常工作的作用。近年來,隨著固體火箭推進劑技術的發(fā)展,對裝藥包覆層的要求也越來越高,為了改善和提高包覆層的各項性能,較多文獻報道了改性環(huán)氧樹脂、硅橡膠、三元乙丙橡膠(EPDM)、聚磷腈包覆層的制備及性能研究[1-5],不飽和聚酯(UPR)由于其良好的力學性能、優(yōu)異的加工性能和良好的抗遷移性而在包覆層應用中備受青睞。目前對UPR 包覆層的研究主要集中于UPR 包覆層的性能和應用方面[6-8],而關于UPR 包覆工藝的研究較少,實際澆注包覆過程中工藝參數(shù)的選擇關乎到固體推進劑包覆過程的工藝安全性和包覆質量。具體而言,包覆澆注過程中樹脂的粘度、流動性、溫度、壓力、澆注速度等都會影響包覆層的質量控制,從而影響到推進劑的質量和工藝安全性,因此有必要對UPR 包覆層澆注過程進行深入研究。

對澆注過程進行模擬仿真,可以大大減少重復試驗的次數(shù)、優(yōu)化包覆工藝、提高包覆效率和質量。Mitani等[9]建立模型模擬環(huán)氧樹脂的流動,使用有限元方法模擬了環(huán)氧樹脂在重力下的非等溫和反應流動行為,且該模型不僅可以用于環(huán)氧樹脂的流動模擬,也可以用于模擬其他熱固性樹脂的成型過程;王慶濤等[10]通過PAM-RTM 軟件模擬比較了樹脂轉移模塑成型(RTM)和真空輔助樹脂轉移模塑成型(VARTM)2 種工藝,研究了不飽和聚酯樹脂粘度、壓力對成型灌注過程的影響;楊進水等[11]對UPR 固化過程中的流變特性進行研究,在修正雙Arrhenius 模型和工程粘度模型的基礎上建立了流變模型,并預測了低粘度工藝操作窗口;Kiehl 等[12]研究了碳酸鈣填料含量不同時的UPR 的流變性能,得到了不同碳酸鈣含 量UPR 的 粘 度-靜 置 時 間 關 系;Aktas 等[13]模 擬 了UPR 在真空澆注過程中的粘度變化,將實驗在不同溫度條件下重復,得到了粘度-時間-溫度曲線;張曼曼等[14]研究了UPR 在動態(tài)升溫及恒溫條件下的粘度變化,建立了工程粘度模型,并使用建立的模型預測了適用于真空輔助成型工藝的低粘度工藝窗口;劉奔 奔 等[15]采 用Carreau 模 型 建 立EPDM 的 流 變 特 性方程,并使用Moldflow 軟件模擬EPDM 包覆層的注射成型過程。已有的研究成果表明,熱固性樹脂的流動過程是較為復雜的過程[16-17],而目前的粘度模型各有局限性,如雙Arrhenius 模型未考慮實際反應的 級 數(shù)[11,18],Castro-Macosko 模 型 忽 略 了 凝 膠 點 附近較為明顯的熱效應[19-21],工程粘度模型不適用于凝 膠 速 度 慢 的 樹 脂[11,19,22]等 問 題。此 外,目 前 的 澆注工藝研究多以環(huán)氧樹脂為研究對象,對UPR 體系的澆注工藝研究研究較少。

為此,本研究以Fontana-Kiuna 模型為基礎建立化學流變模型,得到UPR 固化時的溫度-時間-粘度關系,依據(jù)工藝要求得到了澆注工藝適合的溫度,并在此溫度下采用有限元法模擬了UPR 包覆層的澆注過程,為UPR 包覆層的生產(chǎn)工藝提供參考。

1 模型建立

1.1 UPR 流變模型建立

UPR 在固化過程中存在多種化學反應,經(jīng)歷自由基引發(fā)、微粒凝膠、過度凝膠和大凝膠等階段[23]。由于UPR 的粘度和固化反應相關,而固化反應又和溫度相關,因此粘度本質上是熱歷史的函數(shù)。Fontana-Kiuna模型是由Fontana 首先提出[24]、經(jīng)由Kiuna 修正[25]的一種流變模型,該模型可以在樹脂的固化動力學未知的情況下,在低固化度區(qū)域預測樹脂體系粘度。該模型不能精確預測某時刻的粘度,但能較好模擬整個固化過程,與工程粘度模型、雙Arrhenius 模型相比,F(xiàn)ontana-Kiuna 模型更適合描述UPR 包覆層的流動規(guī)律[25-26]。Fontana-Kiuna 模型的基本形式為[25]:

根據(jù)α與η(T,t)的關系,使用e×η0的值求出t1,即可求出k(T),由τ=k(T) ×t可求出對應的τ值。

將α和τ代入還原,得到式(6)和式(7),整理得到UPR 流變方程為式(8),式中A,B,C均為待定系數(shù)。

1.2 包覆層幾何模型建立

圖1 為澆注包覆層所使用的模具,模具中流體的流域包括圓孔入口部分、底部頂部不規(guī)則部分以及藥柱包覆層。模具整體高214 mm,其中藥柱外層的包覆層高174 mm,厚度為1.5 mm。使用Pro/Engineer軟件建立整體包覆層的幾何模型,使用Hyper Mesh軟件進行網(wǎng)格劃分[28],共49916 個網(wǎng)格,除入口部分采用四面體網(wǎng)格外,其余部位均為六面體網(wǎng)格,幾何模型如圖1 所示。

圖1 模具的幾何模型:(a)整體結構,(b)不規(guī)則部分XY 平面,(c)不規(guī)則部分ZY 平面Fig.1 Geometric model of mould:(a)General structure,(b)XY plane of irregular part,(c)ZY plane of irregular part

計算邊界條件設置:流體入口采用體積流率設定;流體出口不受阻力,可自由流出;其余邊界為壁面邊界。流體澆注過程的模擬計算選用有限元瞬態(tài)體積分數(shù)VOF 模型,流體分數(shù)傳輸任務由系統(tǒng)自動生成[29]。流體入口溫度為298.15 K,模具外表面設定為絕熱邊界。

1.3 本構方程建立

根據(jù)粘度的定義,粘度的通用格式可表示為[30]:

式中,γ為剪切速率,s-1,F(xiàn)(γ)和H(T)分別表示粘度的剪切速率依賴性和溫度依賴性的函數(shù)。

未固化的包覆層為假塑性流體,為了擬合假塑性流體的粘度-剪切速率依賴關系,常用的模型包括冪律模 型、Herschel-Bukely 模 型、Bird-Carreau 冪 律 方 程等[27]。冪律模型僅在有限的剪切速率范圍內才成立,而實際澆注中,開始時的剪切速率很大,最后趨近于零,冪律模型無法滿足要求[31]。Bird-Carreau 冪律方程適用于整個充填過程中的剪切速率分析,僅考慮物理流變(剪切速率-粘度)。熱固性樹脂流變包含化學(時間-粘度)和物理(剪切速率-粘度)2 種流變過程,但UPR 在未固化時粘度變化很小,物理流變對粘度貢獻遠遠大于化學流變,因此可簡化為僅考慮物理流變。Bird-Carreau 冪律方程[32]可以表示為:

式中,η∞為極限剪切速率時的粘度,mPa·s;η0為零剪切速率時的粘度,mPa·s;λ為材料的時間因子;n為冪律指數(shù)。

為了擬合粘度-溫度依賴關系,選用近似Arrhenius 冪律方程:

2 實驗部分

原材料:不飽和聚酯(西安204 所提供),密度:1064.11 kg·m-3,導熱率:0.30 W·(m·K)-1,比熱容:1.60 J·(g·K)-1;環(huán)烷酸鈷(北京市通廣精細化工公司,鈷含量7.8%~8.2%,溶劑為40%~80%的石腦油);過氧化甲乙酮(上海邁瑞爾化學技術有限公司,濃度50%,溶劑為鄰苯二甲酸二辛酯);

流變性能測試:HAAKE MARS 40 應力流變儀(德國Haake 公司)。為了不破壞不飽和聚酯的固化過程,流變測試采用振蕩模式,振蕩頻率1 Hz。取決于不同溫度下的固化快慢,測量持續(xù)時間為1300~4600 s。

3 結果與討論

3.1 粘度-溫度-時間關系的研究

不飽和聚酯包覆層的配方為UPR(100 phr),過氧化環(huán)己酮(3 phr),環(huán)烷酸鈷(0.3 phr),常溫下混合均勻后,使用應力流變儀測試其粘度。測得的不同溫度下,粘度隨時間變化的曲線見圖2,由圖2 可知,隨著時間增加,粘度增大;且溫度越高,初始粘度越小,但是粘度增加更快。

圖2 在不同溫度下,UPR 包覆層粘度隨時間變化的曲線Fig.2 Viscosity vs. time curves of UPR coating layer at different temperatures

將不同溫度下不同時間的粘度值代入式(2)~式(5),計算出相應的α和τ,繪制出α-τ的函數(shù)曲線。如圖3 所示,不同溫度下α-τ的曲線均呈指數(shù)曲線形式,說明在不同溫度下該體系固化規(guī)律沒有改變。

圖3 不同溫度下UPR 包覆層α-τ 的曲線Fig.3 α vs. τ curves of UPR coating layer at different temperatures

將式(8)輸入MATLAB 軟件,使用不同溫度下不同時間的粘度值進行待定系數(shù)擬合,得到A,B,C 分別見表1。

表1 不飽和聚酯包覆層體系參數(shù)表Table 1 Parameters of UPR coating layer system

再將表1 中的A、B、C代入式(8),通過Origin 進行擬合,可分別得到式(13)~式(15),即A、B、C與溫度T的關系式。將A、B、C 帶入式(8)得到式(16):粘度、溫度和時間的關系式。不同溫度下,粘度隨時間變化的曲線,ln(η)-t曲線見圖4,由圖4 可以看出,擬合結果與實驗結果吻合性較好,說明式(8)可以在25~40 ℃內擬合UPR 包覆層的流變性。

圖4 不同溫度下包覆層的ln(η)-t 曲線Fig.4 ln(η)vs. t curve of coating layer at different temperatures

在實際生產(chǎn)過程中,溫度會隨著摩擦生熱、固化反應放熱等有所變化。為了得到25~40 ℃下任意溫度時的粘度-時間關系,使用MATLAB 軟件繪制了粘度-溫度-時間三維曲線(圖5),從而得到一定范圍內任意溫度下的UPR 包覆層粘度特性。如圖5 所示,其中粘度低于1000 mPa·s 可認為是低粘度區(qū),粘度在1000 mPa·s和3000 mPa·s之間可認為是中粘度區(qū),粘度大于3000 mPa·s可認為是高粘度區(qū)。UPR 體系在凝膠點之前粘度變化很小,在達到凝膠點后急劇增大。隨溫度升高,粘度維持在1000 mPa·s 以下的時間先增大后減小,主要原因是低溫下體系初始粘度較大,隨著溫度升高,初始粘度逐漸減小,而繼續(xù)升高溫度,固化反應加速造成粘度升高。根據(jù)包覆工藝要求,UPR 體系粘度低于3000 mPa·s 的時間至少為30 min,而35 ℃時,UPR 體系在低粘度區(qū)保持時間為18 min,粘度低于2000 mPa·s 的時間為30 min。因此,包覆層澆注仿真模擬時設定的溫度應在35 ℃以下,以防提前凝膠化。

圖5 UPR 包覆層溫度-時間-粘度曲線Fig.5 Temperature-time-viscosity curves of UPR coating layer

3.2 包覆層澆注模擬仿真

根據(jù)本構方程(式12),使用所測得的不同溫度下不同時間的粘度值進行擬合,得到:

將式(19)導入Polyflow 軟件進行澆注模擬,同時假設如下:流體為不可壓縮純粘性非牛頓流體;流體與壁面之間無滑移;根據(jù)Polyflow 軟件VOF 計算任務設置,計算時考慮慣性力和重力;澆注過程視為等溫過程[33]。

由圖5 可知,澆注模擬時的溫度應在35 ℃以下,因此,選擇澆注溫度為室溫25 ℃,在此溫度下分別分析UPR 包覆層的澆注速率、澆注壓力對建立的幾何模型充填完整性的影響,以及UPR 充填過程的動態(tài)模擬和充填完成時熔接線的位置等。

3.2.1 澆注速率的影響

受到設備最大壓力限制,澆注速率不能過快,根據(jù)經(jīng)驗,一般最大速率在200 mm3·s-1左右,澆注速率過慢則會降低生產(chǎn)效率,一般不低于50 mm3·s-1。因此,設定澆注體系溫度為25 ℃,澆注速率分別恒定為50 mm3·s-1,100 mm3·s-1,150 mm3·s-1,175 mm3·s-1和200 mm3·s-1,研究不同澆注速率下的包覆層的充填完整性和整體壓力分布。

1) 澆注速率對充填完整性影響

不同恒定澆注速率條件下澆注完成后包覆層的體積分數(shù)分布見圖6,紅色部分表示澆注完整部分,非紅色部分表示澆注不完整(流體體積分數(shù)小于1)。由圖6a 可見,在澆注速率為50 mm3·s-1時,頂部扇形區(qū)域的邊角處有未澆注完全的地方。而圖6b~6e 中流體均填滿模具型腔,沒有孔洞生成,因此,入口處的流率≥100 mm3·s-1才能保證充填完整。

不同澆注速率時流體澆注完成所需的時間見圖6,由圖6 可見,澆注速率為100 mm3·s-1時,澆注時間為400 s;增 加 速 率 到150 mm3·s-1時,澆 注 時 間 縮 短 了133.9 s,僅為266.9 s;而進一步加大速率到175 mm3·s-1,澆注時間變化不大,僅減少了37 s,再增加澆注速率,澆注時間幾乎不變??紤]加工效率,澆注速率應在150 mm3·s-1以上比較合適。

圖6 不同澆注速率時的澆注時間及包覆層體積分數(shù)分布圖Fig.6 Casting time and volume fraction distribution images of coating layer at different casting rates

2) 澆注速率對壓力分布影響

不同澆注速率條件下形成的包覆層最大壓力及壓力分布見圖7,包覆層底部壓力較大,頂部壓力較小,且隨著澆注速率的增大,壓力明顯增加。包覆層壓力的最大點位于澆注入口,當澆注速率為200 mm3·s-1時,入口處最大壓力為3.52×106Pa。

圖7 不同流速時包覆層最大壓力及壓力分布Fig.7 The max pressure and pressure distribution of coating layer at different flow rates

不同流速下的流體流動阻力不同,因而導致入口處的壓力有一定差異,且達到相同填充程度時所需的入口最大壓力不同。圖8 為恒定入口澆注速率時,入口最大壓力與填充程度之間的關系。隨著澆注速率的增加,流體流動速度變快,流動阻力相應增大,達到同樣填充程度時所需的入口最大壓力也隨之增加。由圖8 可見,在充填達到20%之前,入口處最大壓力變化相對較小。填充程度在20%~80%之間時,入口處最大壓力隨填充程度的增加呈線性增大趨勢。填充程度達到80%以后,入口處最大壓力的變化相對來說也比較小。

圖8 不同澆注速率下入口最大壓力隨填充程度變化Fig.8 The Inlet maximum pressure vs. filling degree at different casting rates

在填充過程中,入口處的壓力主要來源于重力和流體的粘性,因此入口處壓力在填充后期比初期大很多。根據(jù)實際操作工藝的要求,設備可以提供的最大壓力為3 MPa。當澆注速率為150 mm3·s-1和175 mm3·s-1時,入口處最大壓力為2.78 MPa和3.15 MPa。因此,對于這種恒速澆注模式,澆注的入口流速應小于175 mm3·s-1。

3.2.2 澆注壓力對充填完整性影響

設定澆注體系溫度為25 ℃,澆注壓力分別恒定為0.5、1、2 MPa 和3 MPa,研究澆注過程中壓力對充填完整性、入口速率分布等影響。

1) 澆注壓力對充填完整性影響

不同澆注壓力條件下,流體填充的體積分數(shù)分布如圖9 所示。澆注壓力為0.5 MPa時,澆注11 min 后包覆層仍沒有填充完全。澆注壓力為1 MPa 時,澆注時間達到500 s 后,除頂部不規(guī)則區(qū)域的扇形區(qū)域有小部分未填充完全外,其余部分均能較好填充;當澆注壓力為2 MPa 和3 MPa 時填充完全,沒有空隙產(chǎn)生。因此,當入口處采用壓力控制條件時壓力應該在1 MPa以上。由圖9 可見,在澆注壓力為3 MPa 時,受到重力和壓力的雙重影響,藥柱底部和頂部的壓力差較大,在固化后可能引起內應力不均勻,影響藥柱的性能。

圖9 不同澆注壓力對體積分數(shù)分布影響Fig.9 Effect of different casting pressure on volume fraction distribution

2) 澆注壓力對入口流速分布影響

圖10 為恒壓條件下澆注開始時,不同澆注壓力條件下的入口流速分布。入口處流速呈現(xiàn)中心大邊緣小的特點。隨著澆注壓力的增加,流體在入口處流速明顯增大,中心位置流動速度最快。當壓力為0.5 MPa時,橫截面整體的流速都很低,中心的流速僅為0.001877 m·s-1;當壓力達到1 MPa 時,中心區(qū)域的流速達到1.211 m·s-1,基本滿足澆注的需求。澆注壓力也不宜過大,若壓力過大,在包覆層的同一高度上,橫截面各部分速度差異大,可能引起內應力不均勻,進而導致固化后收縮不均勻,產(chǎn)生質量問題。

圖10 不同澆注壓力下初始時入口流速分布圖Fig.10 The distribution images of inlet flow rate at different casting pressures

3.2.3 澆注各階段模擬仿真

將包覆層按結構分為三部分,分別為底部流體入口部分、中間圓柱部分和頂部流體出口部分。以澆注速率為150 mm3·s-1的計算結果為例,分析流體的澆注過程。

1) 底部流體入口的澆注模擬

圖11 為底部澆注過程的示意圖和矢量圖,矢量圖根據(jù)速度方向分析流體澆注的軌跡,能夠更加詳細地分析流體澆注過程中的流動情況。由圖11 可知,隨著流體不斷流入型腔,入口部分的流速最大,在澆注入口底部時,包覆層中間圓柱部分已有部分流體流入(圖11a 和圖11d)。在澆注底部流體入口部分時,會形成兩股流體相向匯流的情況(圖11b 和圖11e)。流體在向另外2 個凹槽處流動時,底部還未填充滿,頸部已經(jīng)開始匯流,從而導致另兩個凹槽底部還存在部分空隙(圖11c 和圖11f),由于上方流體已匯流,導致此處空氣無法排出,形成缺陷。

圖11 包覆層底部澆注過程(a~c:示意圖;d~f:矢量圖)Fig.11 Casting process at the bottom of coating layer(a-c:schematic diagram;d-f:vector diagram)

2) 中間圓柱部分的澆注過程模擬

包覆層中間部分是規(guī)則的圓柱,流動過程較為簡單。由圖12a~12c 可知,流體剛開始澆注到中間部分時會分為兩股,并在對側融合,隨后開始向上均勻澆注。由于體系粘度小,這部分會很快融合均一,再加上尺寸比較規(guī)整,所以流動較為穩(wěn)定,不會出現(xiàn)缺陷。

圖12 包覆層中部澆注過程(a~c:示意圖;d~f:矢量圖)Fig. 12 Casting process at the middle of coating layer(a-c:schematic diagram;d-f:vector diagram)

3) 頂部流體出口的澆注過程模擬

流體充滿中部圓柱后,繼續(xù)向頂部凹槽部分流動。頂部凹槽部分的流體流動情況相對來說比較復雜,如圖13所示。先從中間圓柱部位向頸部區(qū)域填充(圖13a 和圖13d),由于流體在中部的包覆層流動時較為平穩(wěn),沒有明顯的高度落差,流體的上升也呈現(xiàn)出整體均勻穩(wěn)定的特點。在經(jīng)過頸部位置后,再次繼續(xù)向上流動逐漸填充凹槽頂部位置,流動狀態(tài)如圖13c 和圖13f。

圖13 包覆層頂部澆注過程(a~c:示意圖;d~f:矢量圖)Fig.13 Casting process at the top of coating layer(a-c:schematic diagram;d-f:vector diagram)

為了進一步觀察中間薄層部分和頂部凹槽部分交界處的流動狀態(tài),每隔2 mm 建立頂部凹槽包覆層的截面。圖14 為所取各截面有流體流入時的體積分數(shù)分布圖,以及流體開始流入各截面對應的時間,圖中194~204 mm 依次表示包覆層頂部凹槽的底面到頂面的位置。由圖14 可知,流體首先填充包覆層的外圍部分,之后向內側流動,再從內側向上流動。這種情況下就會出現(xiàn)多股流體融合的情況,尤其在204 mm 的頂部截面處進行填充時,包覆部分的邊界處會出現(xiàn)明顯的多股流體交融情況(圖14f 中紅色與橙色交接部分),這時可能會出現(xiàn)熔接痕。

圖14 包覆層頂部各截面填充體積分數(shù)分布(a~f:194、196、198、200、202、204 mm 處截面)Fig.14 Distribution of filling volume fraction of the top section of coating layer(a-f:194,196,198,200,202,and 204 mm section)

3.2.4 熔接線模擬預測

根據(jù)上面對包覆層澆注過程中的體積分數(shù)、流速、壓力分布等情況的模擬分析,推測可能會出現(xiàn)缺陷的部位在頂部和底部凹槽部分(圖15)。底部凹槽包覆部分在澆注時有兩股流體相融的情況,在兩股流體相匯部位有可能會出現(xiàn)熔接線。頂部凹槽包覆部分澆注過程中流體流動情況更復雜,會出現(xiàn)多股流體相融合的情況,亦有可能產(chǎn)生熔接線。因此,在受到外力作用時這些部位產(chǎn)生缺陷的可能性較大。

圖15 包覆層的缺陷位置Fig.15 Location of defects of the coating layer

4 結論

本研究以UPR 為基體研究其澆注包覆時的生產(chǎn)工藝,可對包覆層連續(xù)自動化生產(chǎn)起到指導作用。

(1)研究了UPR 熱固性包覆層的化學流變性,得到了粘度-溫度-時間方程,預測了UPR 包覆層的操作工藝窗口,模型結果與實驗能夠很好的吻合。滿足操作工藝要求時,澆注體系的最佳溫度在35 ℃以下,為澆筑過程的模擬提供了合適的溫度參考。

(2)通過Bird-Carreau 冪律方程建立了UPR 包覆層基體的本構方程,將其引入到Polyflow 仿真軟件中模擬包覆層澆注過程中的流動規(guī)律。得到了澆注溫度25 ℃時,恒定不同流率和壓力模式下包覆層澆注的完整性及壓力、速率分布,預測了澆筑完成時熔接線的形成位置。

(3)通過仿真模擬得到了UPR 包覆層澆注工藝參數(shù)為壓力大于1 MPa,入口流速大于150 mm3·s-1且小于175 mm3·s-1。同時包覆層澆注過程中可能出現(xiàn)缺陷的部位在其頂部和底部凹槽部分,但可通過澆注前期降低澆注速率來避免。

猜你喜歡
入口流體速率
高速公路入口疏堵解決方案及應用
納米流體研究進展
山雨欲來風滿樓之流體壓強與流速
秘密入口
第九道 靈化閣入口保衛(wèi)戰(zhàn)
找準入口,打開思路的閘門
蓮心超微粉碎提高有效成分的溶出速率
2011年盈江5.8級地震前近場流體異常初探
乌兰浩特市| 辽中县| 石屏县| 仁寿县| 温州市| 仪陇县| 青海省| 玉山县| 米易县| 新宾| 浑源县| 凌海市| 鱼台县| 尼木县| 白银市| 珲春市| 普陀区| 绥中县| 平顺县| 贡山| 昌吉市| 林芝县| 新化县| 宁海县| 邢台县| 耒阳市| 土默特右旗| 台南县| 咸宁市| 灵寿县| 丹江口市| 康平县| 荆州市| 九寨沟县| 深泽县| 岳西县| 文登市| 平度市| 永修县| 安图县| 昭苏县|