999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

GRAPES模式中不同陸面方案對新疆一次強降水事件的模擬

2014-02-14 01:39:54琚陳相李淑娟于曉晶李曼
沙漠與綠洲氣象 2014年6期
關鍵詞:區域

琚陳相,李淑娟,于曉晶,李曼

(中國氣象局烏魯木齊沙漠氣象研究所,新疆烏魯木齊 830002)

GRAPES模式中不同陸面方案對新疆一次強降水事件的模擬

琚陳相,李淑娟,于曉晶,李曼

(中國氣象局烏魯木齊沙漠氣象研究所,新疆烏魯木齊 830002)

用中尺度數值預報模式GRAPES_Meso V3.3.2.5版本以及NCEP的GFS資料,分別選用模式中不同陸面參數化方案(SLAB、LSM、NOAH)對2013年9月14—17日新疆強降水過程進行數值模擬試驗,結果表明:(1)陸面方案對主要雨帶的落區和大致走向影響并不大,但對降水強度的預報還是敏感的,SLAB陸面方案預報的降水場與實況場最相似;(2)不同陸面方案在選取的強降水區域預報的降水量均較實況偏小,比較各試驗后發現NOAH方案較其他方案的預報結果顯得更加穩定與合理。

GRAPES_Meso模式;暴雨預報;陸面方案

近些年數值模擬在暴雨預報中的作用越來越重要,而地形、地貌和下墊面復雜的特征,土地利用和土地覆蓋的不均勻,陸氣相互作用的多樣性,都對暴雨過程產生作用和影響。陸面作為氣候系統一個非常重要的因素,近幾十年來陸氣相互作用對氣候的影響受到了高度關注[1],并發展了許多描述陸面與大氣之間動量、熱量、水汽及其他微量氣體交換的陸面方案。前人的研究表明[2],在數值模式描述強降水的過程中,模式描述云的形成,大氣邊界層的溫度、濕度、風廓線等信息的準確性都主要依賴于下墊面的土壤、植被等溫度、濕度物體特征。因此,在中尺度模式對降水的模擬過程中,選用合適的陸面方案,對提高模式的預報準確性也非常重要。

Fei Chen等[3]在NCARmM5模式中耦合改進了的陸面模式,指出在MM5模式中耦合陸面過程對表層十分敏感,尤其是在計算熱量和水汽含量方面,同時分子子層的作用也可以提高對表面熱通量的模擬。Jimy Dudhia[4]在MM5模式對多層土壤溫度的預報中,將土壤分成兩層來預報土壤溫度,將上層土壤作為計算垂直尺度的晝夜溫度波的一個深度,將底層土壤作為晝夜平均溫度;而在能量計算中,考慮了感熱通量、潛熱通量和輻射熱通量,上層溫度的變化主要取決于它的有效深度和熱容;文中還指出有效深度和熱容均可以計算。Hsin-I Chang等[5]在WRF模式中討論陸面過程在印度孟買強降水事件(2005年7月20日)的作用,文章采用SLAB、NOAH和NOAH-GEM三種陸面過程,作者指出WRF模式對孟買強降水事件中的各種氣象要素均有較好的模擬,而陸面過程對強降水的模擬有著潛在的影響。馬紅云等[2]利用WRF模式,分別耦合不同的陸面方案(SLAB、RUC、NOAH、UCM)對2007年7月7—8日江淮地區強降水過程進行數值模擬試驗,對比模擬結果表明:陸面方案對強降水的基本位置和大致走向影響不大,但對降水強度的模擬還是比較敏感,不同的陸面方案對降水中心落點、雨量值、降水日變化、降水類型以及降水條件的模擬各有所長。李安泰等[6]對2005年7月1—2日發生在甘肅東南部的一次暴雨天氣進行高分辨率數值模擬,研究NOAH陸面過程中土壤最大容積含水量(MAXSMC)的初始擾動對此次暴雨的敏感性。模擬的結果表明,此次暴雨對陸面參數MAXSMC擾動比較敏感,MAXSMC減少20%,模擬的降水與實況更接近。

本文采用中國氣象局自主研發的中尺度氣象模式GRAPES_Meso,對2013年9月14—17日發生在新疆的一次強降水過程進行數值模擬,以期進一步了解和認識不同陸面方案對強降水預報的準確程度,為合理選擇和使用模式中陸面參數化方案提供一些依據,也為GRAPES_Meso在新疆區域開展研究工作提供一定參考。

1 模式簡介

GRAPES是英文全稱“全球/區域同化預報系統”的縮寫(GRAPES:Global/Regionalassimilationand PrEdiction System),是我國自主研究發展的新一代數值預報系統。GRAPES以多尺度通用動力模式為核心,其框架的主要特點是標準化、模塊化、多尺度通用。GRAPES的核心技術包括三維變分,并可向四維變分拓展的新的資料同化技術;半隱式半拉格朗日時空分離技術和全可壓非靜力平衡動力模式;可自由組合的、優化的物理過程參數化方案;全球、區域一體化的同化與預報系統。

2 個例選取

受西北氣流和低值系統影響,2013年9月14—18日新疆出現強降水過程,北疆各地、天山山區、克州、阿克蘇、哈密和喀什、巴州北部、吐鄯托盆地等地小到中雨,并伴有局地強對流天氣,山區局部雨轉雪。此次降水過程前期,存在弱的徑向環流,歐洲大陸和新疆分別受淺脊控制,里咸海為低值活動區,隨著歐洲脊的不斷發展,脊前北風帶不斷引導冷空氣南下,里咸海到中亞低槽向南加深東移,受到槽前西南氣流的影響,南北疆偏西地區出現明顯降水;過程中后期,里咸海長脊,中亞低槽在東移過程中與西西伯利亞低槽攜帶的冷空氣結合,造成了天山山區及其兩側的明顯降水。

3 試驗方案設計

采用由美國國家環境預報中心(NCEP,National Centers for Environmental Prediction)提供的GFS資料作為初始場。GRAPES_Meso模式的水平網格距設為0.1°×0.1°,垂直層次為26層。模式的模擬區域為(30°~55°N,50°~105°E),中心緯度為42.5°N,初始時間為2013年9月14日20時(北京時,下同),積分時間為72 h,模式積分的時間步長為90 s,每1 h輸出一次預報結果。本文在以下分析中,將15日08時—17日20時作為一次完整的降水過程分析。

微物理過程方案選用WSM6方案,長波輻射選用RRTM方案,短波輻射采用Dudhia方案,近地面層采用Monin-Obukhov方案,邊界層選用MRF方案,積云參數化方案采用Betts-Miller-Janjic方案。采用GRAPES_Meso模式所提供的3套陸面方案:(1)SLAB方案,即5層土壤熱擴散方案,該方案可以描述包括輻射、感熱、潛熱等在內的能量平衡,但是未考慮植被作用及一些物理量隨時間變化的影響;(2)LSM方案[7]是NCAR的陸面模式,包括了能量、動量、水、大氣、土地、不同植被類型以及大氣和土地之間CO2交換的一維模型,該方案主要考慮植被、土壤、水、雪和冰等對太陽輻射的吸收、反射和透射,植物冠層的湍流輸送以及地表徑流和下滲等;(3)NOAH方案[8]是基于俄勒岡州立大學的OSU陸面模式發展而來,包括了一個四層土壤模塊和一層植被冠層的植被模塊,不僅能夠預報土壤溫度,還可以預報土壤濕度、地表徑流等。同時,該方案向邊界層方案提供感熱和潛熱通量,改進了城市冠層,并考慮了地表輻射系數。

4 預報的檢驗方法

為了定量描寫不同陸面方案預報的降水場與觀測場的相似性,擬采用相似系數[9]衡量,計算式如下:

其中ai、bi分別為從站點降水量的模式預報值和觀測值通過Cressman插值法按同一分辨率插值到同一區域(33°~50°N,72°~98°E)均勻網格上的值,n為格點數。相似系數等于1.00(-1.00)為模擬與實況完全相同(相反),等于0.00時表示完全不相似,正值越大圖形越相似,負值越大圖形越相反。

為了定量描寫不同陸面方案預報的強降水區域降水量預報值與實況值變化趨勢的吻合程度,擬采用預報區域模式輸出的逐3 h平均降水量(x1,x2,…,xi)和所在區域站點平均降水量(y1,y2,…,yi)的相關來衡量[10]:

其中降水區域預報值(x1,x2,…,xi)是用雙線性插值法插值到站點后求取區域站點平均而得。式2中,n為樣本容量(20),x、y分別為區域降水量預報與實況的樣本平均值。各檢驗區域所代表的站點見表1。

表1 相關性檢驗的區域和站點

因樣本量小于30時,故采用無偏相關系數r*加以校正。計算式如下:

用t檢驗進行顯著性檢驗。

遵從自由度v=n-2的t分布。給定顯著性水平α,若t>ta,則拒絕原假設,認為相關系數是顯著的。

5 試驗與分析比較

5.1 預報降水場與實況場的相似性比較

圖1給出了全疆105個站的60 h(9月15日 08時—17日20時)累計總降水的實況和預報結果。實際降水帶呈西南、東北向分布,主要降水區位于北疆大部、天山山區、阿克蘇、喀什和和田等地,強降水中心位于伊犁河谷及塔城附近(44°N,81°E和45° N,84.5°E),雨量值達到24mm以上;在其東側和南側也分別存在兩個雨量較小的中心(圖1a)。

圖1b~1d給出3種方案下預報105個站點的60 h總降水量分布,整體來看,不同陸面過程參數化方案均能預報出此次強降水過程的雨帶分布,預報的強降水區域均在天山南北側。然而3種陸面方案對陸面暴雨中心落點和降水強度的預報還存在一定差異,模擬的強降水范圍均較實況偏小,NOAH方案的極值降水范圍最大,其次是SLAB和LSM方案。

根據公式(1)計算相似系數,3種陸面方案(SLAB、LSM和NOAH)的降水預報場與實況場的相似系數分別為0.62、0.58和0.52。SLAB陸面方案預報的降水場與實況場最相似,其次是LSM方案和NOAH方案。

5.2 不同陸面參數化方案對降水趨勢的預報對比

圖1 新疆2013年9月15日08時—17日20時的60 h累計降水量實況和預報結果比較(單位:mm)

圖2 區域平均的觀測降水和預報降水的趨勢演變(單位:mm)

降水趨勢變化作為降水變化的最基本形式,已經成為檢驗一個數值模式正確性的關鍵,正確地模擬降水趨勢變化的振幅和位相可以為數值模式參數化提供理想的試驗基礎。圖2為區域平均降水和預報降水的趨勢演變過程比較,從區域A的實況降水的演變過程看出,此次降水過程主要呈現雙峰型,15日11時降水達到一天中的最低值,然后逐漸增強,降水峰值出現在15日20時左右,隨后降水整體呈增加趨勢,在16日08時左右過程降水達到最大值,而后降水迅速減弱,在16日21時左右,整個降水過程趨于結束。3種陸面方案預報的降水日變化整體與實況降水較接近,但從位相上看均比實況提前2~ 3 h左右。各方案預報的降水振幅與實況也存在一定偏差,3個方案預報的降水量在過程前期均較實況偏大,而后又趨于實況,在過程后均較實況偏小。對比發現,NOAH方案預報的降水演變振幅比LSM方案和SLAB方案更接近實況。而在區域B,實況降水的演變主要呈現單峰型,16日08時降水過程逐漸開始,降水峰值出現在17日00時左右,隨后降水迅速減弱,在17日21時左右,整個降水過程趨于結束。3種陸面方案預報的降水日變化整體與實況降水相差較大,在降水前期均比實況偏大,而在降水過程中又均偏??;位相上也均比實況提前2~3 h左右。各方案預報的降水振幅與實況也存在一定偏差,NOAH方案與實況最為接近,其次為LSM方案和SLAB方案??偟膩砜?,NOAH陸面參數化方案對降水演變趨勢的預報較其他2種方案更為合理和穩定。

為了客觀分析各試驗的預報效果,在此次強降水過程的中心區域伊犁河谷附近(區域A:42°~46° N,80°~85°E)和烏魯木齊附近(區域B:42°~46°N,85°~89°E)分別將各試驗預報的逐3 h平均降水量與實況作相關分析,并對主要降水帶內的60 h總降水量做誤差分析。相關系數可以反映出預報值與實況降水的相關性,絕對誤差則能估量預報值可能的誤差范圍。分析結果見表2,根據計算無偏相關系數,求得區域A的SLAB、LSM和NOAH的t分別為12.81、27.18、12.72,在區域A內,3種陸面方案的預報降水量與實況所求的相關系數均通過了0.001的顯著性檢驗,其中LSM的相關性最好。在區域B,求得SLAB、LSM和NOAH的t分別為0.92、0.43、1.78,只有耦合NOAH陸面方案的試驗通過了0.1的顯著性檢驗。整體來看,NOAH方案在區域A和區域B均通過顯著性檢驗,這說明耦合NOAH方案的試驗預報的降水量比其他方案更合理。從降水量誤差看,3種陸面方案預報的總降水量平均誤差范圍<4mm。降水量的模擬一直是比較困難的工作,特別是極值降水,能否準確的模擬出實況降水量,有待對模式進行更深入的了解與探討。

表2 3種試驗下60 h總降水量預報降水量與實況降水量的相關分析和誤差分析結果

5.3 不同陸面方案對云微物理降水和對流降水的預報比較

模式中描繪的降水由云微物理降水和對流降水組成,分別由云微物理方案和積云對流參數化方案計算得出。圖3給出了各方案預報的對流降水與云微物理降水之間的比較。對流降水的分布較云微物理降水分布廣;而云微物理降水的分布主要集中在伊犁河谷及沿天山一帶的強降水區域。對比3種方案,NOAH方案模擬的對流降水在南疆分布比SLAB和LSM方案分布更廣,云微物理降水的預報也更向東延伸??梢钥闯觯琋OAH方案較其他方案預報出更多、更廣的對流性降水信息,這與其模擬最大降水值高于其他方案的結果較為一致。

圖3 云微物理降水(等值線)和對流性降水(陰影)的比較(單位:mm)

為了定量分析幾種陸面過程參數化方案預報的對流性降水貢獻情況,依次對3種方案在主要降水區伊犁河谷附近(區域A;42°~46°N,80°~85°E)求算對流性降水占降水總量比例的平均值。由計算結果可知,NOAH方案預報的對流性降水占總降水量的比例最高,達到7.3%,其他依次為:SLAB方案6.5%、LSM方案6.4%。雖然強降水的發生與對流系統密切相關,然而在模擬結果中,對流性降水不超過8%。值得指出,模式的分辨率為0.1°×0.1°,對流過程可能有部分被網格尺度所描述,因此,各試驗中對流性降水占總降水的比例均不高。

5.4 對影響系統以及動力和水汽條件的預報比較

在15日08時500 hPa的高空圖上,歐亞中高緯環流維持“兩脊一槽型”,高壓脊之間的西伯利亞一帶為一低壓槽區,槽下為分裂短波,槽線位于70° ~80°E之間。受槽前西南氣流影響,南北疆偏西地區出現了明顯降水,隨后中亞低槽在東移的過程中與西西伯利亞低槽攜帶的冷空氣結合,又造成了天山山區及其兩側的明顯降水。通過試驗結果分析,不同陸面方案均能較好地模擬出這一低值系統,并且模擬的高低空水平環流形勢均與實況較吻合,模擬出15日08時的500 hPa高度場與實況較一致,850hPa在伊犁河谷、北疆沿天山一帶有強西南和西北氣流輻合,伴隨著中亞低槽南伸東移,造成了此次強降水過程。

從圖2可以看出,伊犁河谷附近(區域A;42°~ 46°N,80°~85°E)在15日18時左右開始了新一輪強降水。圖4為16日14時沿44°N的垂直剖面圖,比較各方案模擬暴雨發生時動力和水汽條件可以發現:3種陸面方案模擬的水汽通量較好,在400 hPa左右的水汽通量均能達到0.001 g/s;在對垂直速度的模擬方面,NOAH陸面方案在區域A模擬的垂直速度,不僅在低層存在強烈上升運動,在中高層也有很好的體現。強烈的上升運動和充分的水汽條件對該區域的暴雨發生至關重要。

6 結論

目前,陸面過程對暴雨等短時天氣過程影響的研究尚處于發展階段,盡管大氣本身的動力過程對降水事件的產生至關重要,但從模擬試驗的結果發現,陸面過程的影響也同樣有著重要的作用。本文通過耦合GRAPES_Meso模式中不同陸面方案模擬降水過程進行比較,結果發現:

(1)總體上,60 h累計降水的預報對不同陸面方案是敏感的,不同陸面方案對降水場的預報存在偏差,其中SLAB陸面方案預報的降水場與實況場最相似。

圖4 模擬2013年9月16日14時沿44°N垂直速度(陰影為垂直速度>0.1m/s的上升區)和水汽通量(等值線,單位:g/s)

(2)3種陸面方案模擬的高低空水平環流特征與實況均有較好的一致性,但對動力和水汽條件的模擬還存在一定差異,這將對降水預報的準確性產生影響。

(3)不同陸面方案(SLAB、LSM、NOAH)預報的降水量整體偏差在1.2~4mm之間,NOAH方案在本次暴雨過程的預報較其他2個物理方案更穩定,在選取的兩個區域中均通過了顯著性水平檢驗,較其他方案描繪出更多、更廣的對流性降水信息,對強降水區域動力和水汽條件的模擬也更為合理。

首先GRAPES_Meso模式對區域降水的預報能力主要取決于模式基本的動力框架和物理過程參數化的選擇,陸面過程對降水的影響也是比較有限的。其次,本文討論陸面參數化方案對降水的敏感性僅僅是一個個例,還需更多的個例進行研究和驗證。另外,本文的工作也將為今后GRAPES_Meso在新疆區域開展研究提供一定參考。

[1]曾新民,吳志皇,宋帥,等.WRF模式不同陸面方案對一次暴雨事件模擬的影響[J].地球物理學報,2012,55(1):16-28.

[2]馬紅云,郭品文,宋潔.耦合不同陸面方案的WRF模式對2007年7月江淮強降水過程的模擬過程[J].大氣科學,2009,33(3):557-567.

[3]FeiChen,JimyDuhia.CouplinganAdvancedLand Surface/Hydrologymodel with the Penn State/NCARmM5modeling System.Part II:Preliminarymodel validation[J].monthly Weather Review,2001,129(4):587-604.

[4]Jimy Dudhia.Amulti-layer soil temperaturemodel formM5.TheSixthPSU/NCARMesoscaleModelUsers’Workshop,22-24 July 1996,Boulder,Colorado,49-50.

[5]Chang H I,Kumara,Niyogid,etal.The role of land surface processes on themesoscale simulation of the July 26,2005 heavy rain event overmumbai,India[J].Globaland Planetary Change,2009,67(1):87-103.

[6]李安泰,何宏讓,張云.WRF模式陸面參數擾動對一次西北暴雨影響的數值模擬[J].高原氣象,2012,31(1):65-75.

[7]Bonan.Alandsurfacemodel(LSMversion1.0)for ecological,hydrological,andatmosphericstudies:technical descriptionanduser’sguide.NCARTechnicalNote NCAR/TN-417+STR.National Center foratmospheric Research 1-150.

[8]Chen F,Dudhia J.Couplinganadvanced land surfacehydrology model with the PennState-NCARMM5modelingsystem.PartI:Modelimplementationand sensitivity[J].Monthly Weather Review,2001,129(4):569-585.

[9]曾新民,張強.一次暴雨天氣對陸面參數擾動的敏感性數值影響[J].解放軍理工大學學報(自然科學版),2009,8(10):384-390.

[10]王穎,施能,顧駿強,等.中國雨日的氣候變化[J].大氣科學,2006,30(1):160-170.

[11]魏風英.現代氣候統計診斷與預測技術[M].北京,氣象出版社,2009:18-31.

[12]陳潛,趙鳴.地形對降水影響的數值模擬試驗[J].氣象科學,2006,26(5):484-491.

[13]王薇,張瑛.陸面過程模式的研究進展簡介[J].氣象與減災研究,2010,33(3):1-6.

[14]王秋云,嚴明良,包云軒,等.基于不同陸面參數化方案的高溫天氣數值模擬[J].氣象科技,2011,39(5):537-544.

[15]辛渝,湯劍平,趙逸舟,等.模式不同分辨率對新疆達坂城一小草湖風區地面風場模擬結果的分析[J].高原氣象,2010,29(4):884-893.

[16]李夢婕,申雙和,李雨鴻,等.北京一次下擊暴流的三維數值模擬分析[J].沙漠與綠洲氣象,2013,7(6):22-29.

[17]馬玉芬,趙玲,趙勇.天山地形對新疆強降水天氣影響的數值模擬研究[J].沙漠與綠洲氣象,2012,6(5):40-45.

Simulation ofa Heavy Rainfall Case in Xinjiang Using the GRAPESmodel with Different Land Surface Schemes

JU Chenxiang,LI Shujuan,YU Xiaojing,LIman
(Institute of Desertmeteorology,CMA,Urumqi 830002,China)

Themesoscale numerical predictionmodel GRAEPSmeso V3.3.2.5 versionand National Centers for Environmental Prediction(NCEP)GFS dataare used to simulatea heavy rainfall event in Xinjiang during 14-17 September,2013,and the precipitation sensitivity to land surface parameterization schemes is tested.The results show that:(1)The forecast precipitation patterns from GRAPESmodel generallyagree with the observations,and the sensitivity of the precipitation forecast to the land surface physical process is distinct.The precipitation field forecasted by the SLAB scheme is themost similar to the observation field.(2)The precipitation forecasted by different land surface schemes in selected heavy rainfallarea is smaller than the observation. Comparingall results of precipitation forecast with different land surface parameterization schemes, NOAH scheme givesmore reasonableand stable forecast effect than the others.

GRAPESmesomodel;heavy rainfall forecast;land surface parameterization scheme

P456.7

B

1002-0799(2014)06-0016-07

10.3969/j.issn.1002-0799.2014.06.003

2014-03-31;

2014-09-15

新疆氣象局科學技術研究項目:天山準靜止鋒鋒面三維結構的數值模擬(Q201402)和新疆快速更新循環同化分析數值預報系統技術研究(ZD201404)共同資助。

琚陳相(1988-),男,實習研究員,現從事數值天氣預報與模擬工作。E-mail:juchenxiang@126.com

琚陳相,李淑娟,于曉晶,等.GRAPES模式中不同陸面方案對新疆一次強降水事件的模擬[J].沙漠與綠洲氣象,2014,8(6):16-22.

猜你喜歡
區域
分割區域
探尋區域創新的密碼
科學(2020年5期)2020-11-26 08:19:22
基于BM3D的復雜紋理區域圖像去噪
軟件(2020年3期)2020-04-20 01:45:18
小區域、大發展
商周刊(2018年15期)2018-07-27 01:41:20
論“戎”的活動區域
敦煌學輯刊(2018年1期)2018-07-09 05:46:42
區域發展篇
區域經濟
關于四色猜想
分區域
公司治理與技術創新:分區域比較
主站蜘蛛池模板: 国产女人18水真多毛片18精品| 91av国产在线| 亚州AV秘 一区二区三区 | 国产精品女熟高潮视频| 色男人的天堂久久综合| 欧美一级99在线观看国产| 一本久道久综合久久鬼色| 国产丝袜第一页| 美臀人妻中出中文字幕在线| 亚洲精品无码在线播放网站| 污视频日本| 亚洲性影院| 亚洲伊人电影| 国产免费怡红院视频| 亚洲天堂日本| 亚洲天堂免费在线视频| 欧美在线网| 综合色区亚洲熟妇在线| 欧美亚洲激情| 国产精品 欧美激情 在线播放| 成年人福利视频| 成人在线视频一区| 伊人91视频| 中文字幕在线永久在线视频2020| 久久无码av三级| 中文字幕日韩视频欧美一区| 日韩无码黄色| 中文字幕乱码二三区免费| a级毛片毛片免费观看久潮| 国产一区二区三区在线观看免费| 亚洲国产日韩视频观看| 久久人搡人人玩人妻精品| 亚洲女同欧美在线| 国产aⅴ无码专区亚洲av综合网| 在线99视频| 成人国产精品视频频| 亚洲无码精品在线播放| 亚洲一区二区三区在线视频| 在线观看免费AV网| 在线观看国产精品一区| 亚洲一区二区视频在线观看| 四虎永久在线| 91网在线| 欧美亚洲综合免费精品高清在线观看 | 欧美α片免费观看| 国产肉感大码AV无码| 国产综合欧美| 久久亚洲国产一区二区| 国产麻豆精品久久一二三| 18黑白丝水手服自慰喷水网站| 国产在线观看一区二区三区| 免费看美女自慰的网站| 久久精品无码一区二区日韩免费| 午夜一区二区三区| 久久大香伊蕉在人线观看热2| 日韩精品欧美国产在线| 国产精品一区在线麻豆| 国产免费羞羞视频| 亚洲妓女综合网995久久| 亚洲看片网| 精品久久久久久久久久久| 色综合色国产热无码一| 亚洲日韩图片专区第1页| 69视频国产| 国产菊爆视频在线观看| 经典三级久久| 69视频国产| 久久夜夜视频| 亚洲日本韩在线观看| 欧洲欧美人成免费全部视频| 久久婷婷五月综合97色| 亚洲国产高清精品线久久| 国产主播在线一区| 午夜精品国产自在| 她的性爱视频| 黄色一级视频欧美| 国产簧片免费在线播放| 一本久道热中字伊人| 国产精品久久久久久久伊一| 亚洲综合香蕉| 青草精品视频| 美女毛片在线|