• 
    

    
    

      99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

      復雜荷載環(huán)境下海上風力機的建模及動力學特性分析

      2016-12-21 03:31:30王青占趙建中郭興明
      上海大學學報(自然科學版) 2016年5期
      關鍵詞:三階風力機固有頻率

      王青占,趙建中,郭興明

      (上海大學上海市應用數(shù)學和力學研究所,上?!?00072)

      復雜荷載環(huán)境下海上風力機的建模及動力學特性分析

      王青占,趙建中,郭興明

      (上海大學上海市應用數(shù)學和力學研究所,上海200072)

      通過分析處于海風、海浪、海流、土壤力等復雜荷載作用下的海上風力機支撐結(jié)構(gòu),采用赫維賽德階躍函數(shù)和狄拉克δ函數(shù)建立了連續(xù)統(tǒng)一的、頂端帶有集中質(zhì)量塊的懸垂梁風力機動力學模型.基于對控制方程的Galerkin截斷,得到離散化的常微分方程組,使用四階Runge-Kutta方法求解,得到了模型的動力學行為特性云圖曲線.通過懸垂梁風力機模型的時程曲線、龐加萊映射對風力機模型進行動力學分析,給出了位移和速度的幅值隨激振力頻率變化的幅頻特性曲線,并研究了垂向激振、自重、變剛度參數(shù)對風力機結(jié)構(gòu)振動特性及其穩(wěn)定性的影響.

      海上風力機;變軸力;復雜荷載;變剛度;動力響應

      近年來,各國對風能、太陽能等可替代的潔凈新能源的開發(fā)利用力度加大,尤其是對海上風力機的研究越發(fā)引起各國科研工作者的重視[1].以英國、丹麥和德國為代表的歐洲走在了海上風力發(fā)電研究的前端,近兩年的裝機總量不斷增加,發(fā)電量所占的比重相應提高.圖1為海上風力機系統(tǒng)的示意圖.海上風力機系統(tǒng)由風力機、支撐結(jié)構(gòu)和地基基礎三部分組成,其中風力機由葉片、輪轂和機艙構(gòu)成;支撐結(jié)構(gòu)包括塔筒和下部結(jié)構(gòu),下部結(jié)構(gòu)分為固定式和漂浮式兩種形式.無論支撐結(jié)構(gòu)采用固定式還是漂浮式,海上風力機與陸地風力機、海上油氣平臺工程都有較大區(qū)別[2],在系統(tǒng)結(jié)構(gòu)、環(huán)境條件、荷載特征等方面都具有特殊性.首先,海上風力機與塔架要經(jīng)受臺風的嚴峻考驗,其下部結(jié)構(gòu)還要受到波浪、海流以及與土體之間的耦合作用[3],且風、浪、流是有一定耦合的.國內(nèi)外的相關學者和研究機構(gòu)也開始對海上風電工程和其他浮體項目進行研究.南京航空航天大學胡文瑞等[4]進行了大型風力機的空氣動力學基礎研究,中國科學院力學所周濟福等[5]針對海上風電工程和地基的關鍵力學問題進行了相關研究.國內(nèi)一些學者通過有限元建模仿真計算對海上風力機進行氣動力學、結(jié)構(gòu)靜力學和動力學的研究[6-9].國外學者主要對風力機結(jié)構(gòu)進行建模分析[10-11],包括對不同基礎、浮動風機以及海床結(jié)構(gòu)的參數(shù)化模擬分析.但是考慮復雜海洋環(huán)境的綜合作用,以及風力機結(jié)構(gòu)剛度變化的報道較少.本工作建立的帶有集中質(zhì)量塊的懸垂梁風力機動力學模型較好地考慮了這些因素,研究的結(jié)果也可為后續(xù)風力機的優(yōu)化設計和安全評估提供重要參考.

      圖1 海上風力機系統(tǒng)示意圖Fig.1 Sketch of offshore wind turbines

      1 海上風力機的結(jié)構(gòu)分析與簡化模型

      對于海上風力機模型的研究,采用較多的是有限元仿真計算.例如,用有限元進行非線性結(jié)構(gòu)動力學、流固耦合或者土體與結(jié)構(gòu)耦合[12]、葉輪氣動力學研究等.王磊等[13]通過計算機仿真,結(jié)合水動力學模型和風輪空氣動力學模型,建立了“風輪-機艙-塔筒-系泊系統(tǒng)”組成的多柔體系統(tǒng)動力學模型,并對漂浮式風力機系統(tǒng)和近海的定樁式風力機系統(tǒng)進行動力學分析對比,研究結(jié)果表明氣動載荷與水動力相互耦合對整機結(jié)構(gòu)動力響應及功率波動有著明顯影響.空間太陽能電板、飛機機翼、海上風力機塔筒等大尺寸細長物體,在某種程度下都可以簡化為梁模型.Piana[14]和Virgin[15]通過建立類似的彈性結(jié)構(gòu)模型研究了上述結(jié)構(gòu)的動態(tài)穩(wěn)定性和振動之間的相互作用,進一步分析了結(jié)構(gòu)承受類似于軸向壓力的載荷(盡管可能不會直接導致結(jié)構(gòu)失穩(wěn))對結(jié)構(gòu)的固有頻率產(chǎn)生的影響,并以此結(jié)果為參照對結(jié)構(gòu)施加控制[16].

      Lajimi等[17]用Frobenius冪級數(shù)法研究了質(zhì)量塊的質(zhì)量、旋轉(zhuǎn)慣性以及偏心距對頂端承受集中質(zhì)量塊固支懸垂梁的固有頻率的影響,研究結(jié)果表明該結(jié)構(gòu)的固有頻率隨著集中質(zhì)量塊質(zhì)量的增加而降低.Auciello[18]研究了水對風力機模型固有頻率的影響,結(jié)果表明水只對較高的固有頻率產(chǎn)生顯著影響,而對于通??紤]的前三階固有頻率的影響不大.U′sci′lwska等[19]建立了半浸水頂端承受集中質(zhì)量的懸垂梁模型,并以水面為界限建立了兩段梁的運動方程,將集中質(zhì)量塊歸類于邊界條件,求解六階奇異矩陣得到了控制其固有頻率的特征方程.Andersen等[20]通過彈簧模擬樁體與土壤之間的p-y曲線,研究了海床對風力機支撐機構(gòu)的影響.

      工程設計中海上風力機塔架的固有頻率必須要避開風力機轉(zhuǎn)子的自激頻率.通過對海上風載以及海浪沖擊荷載的觀測和統(tǒng)計,找出風力機設計自然頻率的安全頻域,通過本工作的參數(shù)化計算和模擬,可以選擇合適的風力機設計參數(shù),提高風力機的效率和使用壽命.通過本工作所考慮的機艙重量和風力機工作動荷載進行調(diào)節(jié)和振動控制,也可以結(jié)合海上環(huán)境氣候統(tǒng)計[21]分析來調(diào)節(jié)海上復雜荷載的組合參數(shù),進而更加真實模擬復雜的海洋作用環(huán)境.綜上所述,本工作所建立的模型和計算方法有很好的適應性,可用于不同海洋環(huán)境和風力機參數(shù)的普適計算.

      2 懸垂梁風力機模型的動力學方程

      2.1動力學一般方程的建立

      海上風力機支撐結(jié)構(gòu)模型如圖2所示.整個艙重和葉輪用一集中質(zhì)量塊M表示,葉片承受的風載等效為作用在輪轂處的集中力f1(x,t),塔體所受的風載為f2(x,t),海浪潮汐的作用為f3(x,t),海流的作用為f4(x,t),海床巖土對塔體的作用為f5(x,t).同時引入隨徑向和時間變化的軸力,以及由填充鋼筋混凝土引起的變剛度高度,依據(jù)能量法建立Euler-Bernoulli懸垂梁模型的運動方程[22]:移,M為塔筒頂端(機艙和葉片)質(zhì)量,ε和σ分別為風力機振動的動質(zhì)量參數(shù)和額定工作頻率, l1為塔筒內(nèi)部混凝土填充高度,m為單位長度空心塔筒的質(zhì)量,mc為單位長度填充鋼筋混凝土的質(zhì)量,dc為塔筒內(nèi)徑,H(x)表示赫維賽德階躍函數(shù),δ(x)為狄拉克脈沖函數(shù).

      圖2 海上風力機支撐結(jié)構(gòu)模型Fig.2 Structure model of offshore wind turbines

      2.2控制方程及邊界條件的無量綱化處理

      為方便計算,對懸垂梁風力機模型進行無量綱化處理.引入無量綱參數(shù),令

      式中,Υ為位移標定參量,v(ξ,τ)為橫向振動位移w(x,t)的無量綱表示,ξ1為鋼筋混凝土填充變量,ξ2為反映海域深淺的量,ξ3為反映海床環(huán)境的量,κ為混凝土剛度調(diào)節(jié)參數(shù),κ1為混凝土與塔筒旋轉(zhuǎn)慣性比參數(shù),η1為塔筒與艙體的質(zhì)量比(又稱固有頻率影響因子),η2為旋轉(zhuǎn)慣性參數(shù),κ3為混凝土質(zhì)量填充參數(shù),α為集中質(zhì)量塊慣性影響參數(shù),τ為無量綱時間尺度,σ為風機激振頻率,ζ為結(jié)構(gòu)阻尼比.

      風力機模型的無量綱控制方程為

      邊界條件的無量綱化為

      2.3復雜荷載及其無量綱化處理

      風力發(fā)電機組運行時,其葉片上的風荷載和風力機偏航引起的荷載通過結(jié)構(gòu)和傳動機構(gòu)作用在塔架頂端,所以相關規(guī)范規(guī)定海上風電機組基礎結(jié)構(gòu)設計應考慮風電機組的荷載.這部分荷載包括風輪上的靜風壓引起的荷載、湍流和尾流引起的荷載、風力發(fā)電機偏航引起的荷載和風力發(fā)電機組的重力荷載等.動荷載包括風力機運轉(zhuǎn)荷載,以及風、浪、流、冰和地震荷載.為了便于分析計算,分別將風輪荷載簡化為輪轂集中力和風力機垂向激振荷載,將塔筒承受的風荷載簡化為均勻脈動的荷載,將海面處承受的海浪等沖擊荷載用脈沖荷載表示,土壤與塔基之間的作用力[12]用近似簡化的p-y曲線表示.具體化簡后的荷載表示如下.

      海平面處的波浪、海冰船舶沖擊荷載(集中):

      式中,p5(ξ,τ)=?5H(ξ3?ξ)tanh(2πξ sin(2ξ))或?5H(ξ3?ξ)tanh(2πξ?1(2ξ)),A=1(與循環(huán)荷載有關),B與土壤質(zhì)地(密度)有關,Pu為土體極限抗力.

      2.4Galerkin截斷選取的模態(tài)函數(shù)

      模態(tài)函數(shù)為

      2.5無量綱控制方程的Galerkin截斷

      利用選取的模態(tài)函數(shù)可求得方程(3)的三階Galerkin離散方程:

      式中,φi,j,k(i,j=1,2,3,k=1,2,…,8)與梁的剛度有關,是反映剛度變化的量.

      2.6控制方程的矩陣化表示

      為了從整體上更加清晰地了解該懸垂梁的計算模型,對控制方程進行如下的矩陣化歸一表示

      ξ1,ξ2,ξ3為反映風力機海洋環(huán)境和風力機剛度控制的設計參量,Λ矩陣清晰反映了風力機所承受的海洋環(huán)境荷載.如能獲得更加詳盡和真實的海洋環(huán)境和風力機的資料,理論上可以通過這些參數(shù)組合實現(xiàn)對風力機設計的操控和優(yōu)化,為設計更加合理的風力機提供重要的分析依據(jù).

      3 特定條件下求解結(jié)果的對比

      3.1Galerkin截斷誤差及有效性驗證

      為了驗證Galerkin截斷的收斂性和有效性,取自由振動的均勻懸臂梁進行對比驗證.

      圖3為懸臂梁自由振動控制方程近似解析解與三階Glerkin截斷數(shù)值解的時空振動曲線.可以看出,二者整體對比的結(jié)果基本一致.圖4為懸臂梁自由振動控制方程近似解析解與三階Galerkin截斷數(shù)值解的絕對誤差和相對誤差.可以看出,位移最大的頂端處的相對誤差<3%.因此,該模態(tài)下的Galerkin截斷完全可以滿足工程中的精度要求.

      3.2數(shù)值結(jié)果與有特定解的解析解對比

      圖3 方程近似解析解與三階Galerkin截斷數(shù)值解的時空振動曲線Fig.3 Vibration curves of approximate analytical and numerical solution of 3-order Galerkin truncation

      圖4 方程近似解析解與三階Galerkin截斷數(shù)值解的絕對誤差和最大相對誤差Fig.4 Absolute and maximum relative errors of approximate analytical and numerical solution of 3-order Galerkin truncation

      由此求得無量綱控制方程的解析解為w(ξ,τ)=?cosτ(1+ξ?0.679 7e-0.9934ξ?0.3463e0.9909ξ+e-0.00125ξ(0.0209cos(1.0079ξ)?1.3216sin(1.0079ξ))).此處w(ξ,τ)與v(ξ,τ)類似,均為相同單位標定的無量綱量.

      圖5 方程數(shù)值結(jié)果與解析結(jié)果的對比Fig.5 Comparisons between numerical and analytical solution for the equation

      圖5為懸臂梁自由振動控制方程數(shù)值結(jié)果與解析結(jié)果的對比.可以看出,兩種方法得出的結(jié)果吻合較好.產(chǎn)生誤差的原因可能是由于模型頂端集中質(zhì)量塊處邊界條件的近似處理,這也從側(cè)面反映了機艙的重量對其振動特性有顯著的影響.優(yōu)化機艙重量或者在機艙中設計減震與控制裝置[16],可以提高風力機的使用年限和安全系數(shù).

      圖6為特殊荷載下懸垂梁風力機模型的幅頻特性曲線.根據(jù)文獻[24]的海上風力機頻域設計區(qū)域理論,可以看出,在強有力的海上荷載作用下,該參數(shù)下的風力機處于Soft-Stiff區(qū)域,是比較安全的設計.計算結(jié)果表明,風力機在頻率fo=0.37 Hz(轉(zhuǎn)化為無量綱參數(shù)λo=2.375 6)的慣常荷載作用下是安全的.

      圖6 特殊荷載下風力機模型的幅頻特性曲線Fig.6 Amplitude-frequency curve of the wind turbine model under specific loads

      4 數(shù)值仿真算例

      本工作以文獻[24]中給出的海上風力機支撐結(jié)構(gòu)為算例,取結(jié)構(gòu)額定功率為3 MW,塔筒主體材料為Q345鋼且塔筒內(nèi)部填充了鋼筋混凝土.

      4.1原始物理參數(shù)

      將以上物理量代入與其對應的無量綱表達式(3)中,并將所得的結(jié)果代入式(7)~(9)中,再通過式(10)計算得到控制方程的系數(shù)矩陣(剛度矩陣)和荷載矩陣:

      4.2動力學數(shù)值分析

      通過四階Runge-Kutta數(shù)值方法求解矩陣微分方程組,可得到穩(wěn)態(tài)時前三階動力學時程曲線、時空響應曲線(見圖7),可以看出u2與u3重合.這些基本曲線反映了梁的動力學特性[25].

      圖7 穩(wěn)態(tài)下海上風力機頂端處的時程曲線和幅頻曲線Fig.7 Time-history and amplitude-frequency curves in stationary state for the top end of offshore wind turbines

      從圖7(b)可以看出,危險頻段為1.8~2.2,折算為實際的頻率為0.28~0.34Hz.如果荷載頻率在此頻段內(nèi),說明設計不合理,就需要調(diào)節(jié)其他影響固有頻率的參數(shù).從梁受迫振動的前三階時程曲線(見圖7(a))、相圖(見圖8)和龐加萊映射(見圖9)可以看出,一階截斷存在較大的誤差,而二、三階截斷誤差要小得多.因此用Galerkin截斷來計算梁的動力曲線是收斂的.從圖9中可以看出,龐加萊映射是一極限環(huán)狀的,可見在此參數(shù)下梁的運動是準周期的、穩(wěn)定的響應.

      圖8 頂端處的穩(wěn)態(tài)前三階Galerkin截斷相圖Fig.8 Steady-state first 3-order Galerkin truncation phase diagram at top end

      圖9 頂端處的穩(wěn)態(tài)前三階Galerkin截斷龐加萊映射Fig.9 Steady-state first 3-order Galerkin truncation Poincar′e map at top end

      4.3變剛度參數(shù)影響分析

      圖10為載荷條件相同時無填充和半填充鋼筋混凝土塔筒情況下風力機結(jié)構(gòu)的動力特性.

      圖10 荷載條件相同時無填充和半填充鋼筋混凝土塔筒的動力特性Fig.10 Dynamic characteristics of no-filled and half-filled wind turbines under the same load

      從圖10可以看出,填充混凝土后風力機結(jié)構(gòu)的最大振幅(0.165個標準)要比無填充時(0.465個標準)小得多.因此,塔筒填充鋼筋混凝土對風力機結(jié)構(gòu)的振動特性有著顯著的影響,通過調(diào)節(jié)海上風力機塔架設計的固有頻率,可以避開海上荷載的常見頻段組合.

      5 結(jié)束語

      本工作建立了海上風力機支撐結(jié)構(gòu)的懸垂梁動力學模型.研究結(jié)果表明,該模型可以較好地模擬風力機支撐結(jié)構(gòu)的工作現(xiàn)狀.通過分析計算復雜載荷下的結(jié)構(gòu)動力特性,研究各種荷載、風力機安裝環(huán)境參數(shù)、混凝土的填充高度、適宜海域的水深和海床條件等對風力機安全性的影響,不僅可以驗證風力機在特定海洋環(huán)境的安全性能,同時也可以對風力機的參數(shù)進行優(yōu)化設計.相對于有限元計算,這種數(shù)值計算模型大大縮減了優(yōu)化的時間成本;相對于模擬試驗,較好地降低了試驗成本;而與解析條件下的情況對比,則大大提高了適應求解的范疇.

      [1]林鶴云,郭玉敬,孫蓓蓓.海上風電的若干關鍵技術(shù)綜述[J].東南大學學報(自然科學版),2011,41(4):882-888.

      [2]尚景宏.海上風力機基礎結(jié)構(gòu)設計選型研究[D].哈爾濱:哈爾濱工程大學,2010:11-19.

      [3]李益,張宏戰(zhàn),馬震岳.基于土-結(jié)構(gòu)相互作用的海上風機模態(tài)分析[J].水電與新能源,2013(1):70-73.

      [4]胡文瑞,王同光,李曄.“大型(兆瓦級)風力機的空氣動力學問題”專刊前言[J].應用數(shù)學和力學, 2013,34(10):1-2.

      [5]周濟福,林毅峰.海上風電工程結(jié)構(gòu)與地基的關鍵力學問題[J].中國科學:物理學力學天文學, 2013,43(12):1590-1593.

      [6]陳前,付世曉,鄒早建.海上風力發(fā)電機組支撐結(jié)構(gòu)動力特性分析[J].振動與沖擊,2012,31(2):86-90.

      [7]顧岳飛.大兆瓦風電機組塔架的有限元分析與優(yōu)化設計[D].上海:上海交通大學,2012:18-45.

      [8]李立崗.風力發(fā)電機塔架結(jié)構(gòu)動力特性分析[D].哈爾濱:哈爾濱工業(yè)大學,2011:21-29.

      [9]嚴磊.風力發(fā)電機支撐體系結(jié)構(gòu)設計研究[D].天津:天津大學,2008:37-41.

      [10]LEBLANc C.Design of offshore wind turbine support structures[D].Aalborg:Aalborg University, 2009:16-27.

      [11]AGARwAL P.Structural reliability of offshore wind turbines[D].Texas:The University of Texas at Austin,2008.

      [12]LOMBARDI D,BHATTAcHARyA S,WOOD D M.Dynamic soil-structure interaction of monopile supported wind turbines in cohesive soil[J].Soil Dynamics and Earthquake Engineering,2013, 49:165-180.

      [13]王磊,何玉林,金鑫,等.漂浮式海上風電機組動力學仿真分析[J].中南大學學報(自然科學版), 2012,43(4):1309-1314.

      [14]PIANA G.Vibrations and stability of axially and transversely loaded structures[D].Turin:Polytechic University Turin,2013.

      [15]VIRGIN L N.Vibration of axially-loaded structures[M].Cambridge:Cambridge University Press, 2007.

      [16]HE W,GE S S.Vibration control of a nonuniform wind turbine tower via disturbance observer[J].IEEE/ASME Transactions on Mechatronics,2014,20(1):237-244.

      [17]LAjIMI S A M,HEppLER G R.Free vibration and buckling of cantilever beams under linearly varying axial load carrying an eccentric end rigid body[J].Transactions of the Canadian Society for Mechanical Engineering,2013,37(1):89-109.

      [18]AUcIELLO N M.Dynamic analysis of immersed tapered column carrying an eccentric mass at the free[J].Advanced Research in Scientific Areas,2013(1):412-417.

      [19]U′ScI'LOwSKA A,KO'LODzIEj J A.Free vibration of immersed column carrying a tip mass[J]. Journal of Sound and Vibration,1998,216(1):147-157.

      [20]ANDERSEN L V,VAHDATIRAD M J,SIcHANI M T,et al.Natural frequencies of wind turbines on monopile foundations in clayey soils—a probabilistic approach[J].Computers and Geotechnics, 2012,43(11):1-11.

      [21]BERAN J,CALvERI L,LANGE B,et al.Offshore wind modelling and forecast[C]//Proceedings of the 6th WRF/15th MM5 Users’Workshop.2005:1-9.

      [22]ADHIKARI S,BHATTAcHARyA S.Vibrations of wind-turbines considering soil-structure interaction[J].Wind and Structures,2011,14(2):85-112.

      [23]劉延柱,陳文良,陳立群.振動力學[M].北京:高等教育出版社,1998:130-135.

      [24]ADHIKARI S,BHATTAcHARyA S.Dynamic analysis of wind turbine towers on flexible foundations[J].Shock and Vibration,2012,19(1):37-56.

      [25]丁虎,嚴巧赟,陳立群.軸向加速運動黏彈性梁受迫振動中的混沌動力學[J].物理學報,2013, 62(20):200502.

      本文彩色版可登陸本刊網(wǎng)站查詢:http://www.journal.shu.edu.cn

      Dynamics modeling and analysis of offshore wind turbines under complicated loads

      WANG Qingzhan,ZHAO Jianzhong,GUO Xingming
      (Shanghai Institute of Applied Mathematics and Mechanics,Shanghai University, Shanghai 200072,China)

      By analyzing the structure and complexity of offshore wind turbines(OWTs), a beam model with free end carrying a concentrated mass body was built to analyze dynamical characteristics of OWTs subject to adverse working environment loads. Factors including vertical excitation,self-weight,rotary inertia and variety of stiffness were considered to study the effect on offshore wind turbines.A general and uniform governing equation was built by taking advantages of the Heaviside step function and Dirac delta function.Based on the Galerkin truncation and the Runge-Kutta time discretization, numerical solutions of the kinematic governing equation were obtained.By comparing with the analytical solution under specified conditions,validity of the method was checked.The time history of the beam’s free end was chosen to represent motion of the beam.Based on the steady time history of the beam,a Poincar′e map was constructed to study its periodic motion.Furthermore,an amplitude-frequency curve was given to find the dangerous frequency range where OWTs exist.

      offshore wind turbines;varying axial force;complicated load;varying rigidity;dynamic response

      TH 43

      A

      1007-2861(2016)05-0573-13

      10.3969/j.issn.1007-2861.2015.02.016

      2015-06-16

      國家重點基礎研究發(fā)展計劃(973計劃)資助項目(2014CB0462003)

      郭興明(1964—),男,教授,博士生導師,博士,研究方向為連續(xù)介質(zhì)力學與力學中的數(shù)學方法. E-mail:xmguo@shu.edu.cn

      猜你喜歡
      三階風力機固有頻率
      三階非線性微分方程周期解的非退化和存在唯一性
      現(xiàn)場測定大型水輪發(fā)電機組軸系的固有頻率
      基于UIOs的風力機傳動系統(tǒng)多故障診斷
      三類可降階的三階非線性微分方程
      大型風力機整機氣動彈性響應計算
      小型風力機葉片快速建模方法
      太陽能(2015年6期)2015-02-28 17:09:35
      總溫總壓測頭模態(tài)振型變化規(guī)律研究
      三階微分方程理論
      A novel functional electrical stimulation-control system for restoring motor function of post-stroke hemiplegic patients
      轉(zhuǎn)向系統(tǒng)固有頻率設計研究
      桐柏县| 卢龙县| 金平| 石家庄市| 泽州县| 响水县| 延寿县| 东丰县| 五大连池市| 钟祥市| 兴文县| 蚌埠市| 东乌珠穆沁旗| 永和县| 无极县| 哈尔滨市| 新邵县| 安阳市| 游戏| 平度市| 安新县| 怀安县| 乌拉特中旗| 济南市| 桐梓县| 南昌县| 平安县| 沙田区| 西宁市| 皋兰县| 温州市| 莱阳市| 高台县| 甘泉县| 中宁县| 靖边县| 旌德县| 万年县| 高雄市| 古浪县| 南投县|