• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    邊界層高度和動量粗糙度的不確定性對臺風(fēng)“天鴿”(1713)模擬的影響

    2023-12-05 10:42:52譚曉偉
    海洋氣象學(xué)報 2023年1期
    關(guān)鍵詞:參數(shù)化邊界層數(shù)值模擬

    譚曉偉

    摘要:基于中國氣象局區(qū)域臺風(fēng)模式CMA-TYM模式,通過一系列敏感性試驗,對臺風(fēng)“天鴿”(1713)生命過程中不同時間階段的移動路徑及強(qiáng)度的數(shù)值模擬對模式參數(shù)化方案中邊界層高度(h)和動量粗糙度(z0)的敏感性進(jìn)行了研究。試驗結(jié)果表明,使用不同參數(shù)化方案計算的h在“天鴿”(1713)初期的熱帶風(fēng)暴階段對熱帶氣旋的移動路徑有較明顯的影響,在臺風(fēng)成熟后對熱帶氣旋的移動路徑影響不顯著;臺風(fēng)中心附近最大10 m風(fēng)速對h的變化不敏感,而最低海平面氣壓對h的變化卻非常敏感。同時在臺風(fēng)發(fā)展初期階段,邊界層過薄或過厚都不利于臺風(fēng)強(qiáng)度的發(fā)展加強(qiáng),這表明邊界層高度h在熱帶氣旋數(shù)值模擬和預(yù)報中是非常重要的,尤其是在臺風(fēng)發(fā)展的初期階段。動量粗糙度z0的變化對臺風(fēng)“天鴿”的影響主要體現(xiàn)在臺風(fēng)增強(qiáng)階段,臺風(fēng)中心附近10 m風(fēng)速最大值在臺風(fēng)增強(qiáng)階段對z0敏感,尤其是在臺風(fēng)發(fā)展的初期階段。

    關(guān)鍵詞:臺風(fēng);數(shù)值模擬;邊界層;參數(shù)化

    中圖分類號:P444;P456? ?文獻(xiàn)標(biāo)志碼:A? ? ?文章編號:2096-3599(2023)01-0001-00

    DOI:10.19513/j.cnki.issn2096-3599.2023.01.003

    Influence of the uncertainty of boundary layer height and momentum roughness on Typhoon Hato (1713) simulation

    TAN Xiaowei1,2

    (1. CMA Earth System Modeling and Prediction Centre (CEMC), Beijing 100081, China; 2. State Key Laboratory of Severe Weather, Chinese Academy of Meteorological Sciences, Beijing 100081, China)

    Abstract: Based on the CMA-TYM model, through sets of sensitivity experiments, the sensitivities of the simulated moving track and intensity of Typhoon Hato (1713) to the boundary layer height h and momentum roughness z0 in different life stages are studied. The experiment results show that the boundary layer height h from different parameterization schemes has significant influence on the moving track during the tropical storm stage, while there is no obvious influence on the moving track after Hato gets mature. The maximum 10-m wind speed near the typhoon center is not sensitive to the change of h, while the minimum sea-level pressure is very sensitive to the change of h. The study also finds that, in the early stage, it is not beneficial for the development of typhoon intensity whether h is too thin or thick, which shows that h is very important in the numerical simulation and prediction of tropical cyclones, especially in the early stage of typhoon. The influence of the change of momentum roughness z0 on Typhoon Hato is mainly reflected in the typhoon intensification stage. The maximum 10-m wind speed near the typhoon center is sensitive to z0 in the typhoon intensification stage, especially in the early stage of typhoon.

    Keywords: typhoon; numerical simulation; boundary layer; parameterization

    引言

    熱帶氣旋(tropical cyclone,TC)會對生命和財產(chǎn)造成巨大危害,為了減少這類損失,相關(guān)機(jī)構(gòu)和個人都需要獲得TC的位置和強(qiáng)度的準(zhǔn)確預(yù)報。眾所周知,邊界層是熱帶氣旋的重要組成部分,邊界層過程是大氣的主要熱量和水汽的源,是動量的匯,調(diào)節(jié)著與TC發(fā)展及增強(qiáng)密切相關(guān)的動量和焓的徑向和垂直分布,在TC的能量傳輸過程中起著重要的作用。然而TC邊界層與大部分大氣邊界層不同,TC的旋轉(zhuǎn)對其動力學(xué)有強(qiáng)烈的影響[1-5],這也進(jìn)一步增加了TC數(shù)值模擬中邊界層參數(shù)化的復(fù)雜性[6]。

    在數(shù)值模式中,邊界層高度h和動量粗糙度z0都是需要進(jìn)行參數(shù)化且具有重要意義的變量。例如,在被業(yè)務(wù)數(shù)值模式廣泛使用的一階K閉合方案(如Medium Range Forecast Model(MRF)方案和Yonsei University(YSU)方案)中,混合層的湍流擴(kuò)散系數(shù)需要使用混合層特征速度和邊界層高度h進(jìn)行參數(shù)化,而混合層特征速度的參數(shù)化也依賴于邊界層高度h。YSU方案的一個特點是包含了邊界層頂夾卷的參數(shù)化計算,而該夾卷項的參數(shù)化也依賴于邊界層高度h。因此,在數(shù)值模式中邊界層高度h影響著邊界層內(nèi)湍流熱量通量的垂直輸送,進(jìn)而影響TC的發(fā)展。另外,在大多數(shù)參數(shù)化方案中,描述地(海)氣交界面的動量、感熱和潛熱交換作用的表面通量要依賴于動量粗糙度z0,部分參數(shù)化方案會使用與z0有關(guān)的表面動量交換系數(shù)CD(也稱為曳力系數(shù))代替z0。因此,在數(shù)值模式中,邊界層高度h和動量粗糙度z0(或表面動量交換系數(shù)CD)的確定在邊界層參數(shù)化過程中十分關(guān)鍵。

    但是,由于數(shù)值計算過程中參數(shù)化方法的設(shè)計及對TC缺少足夠的科學(xué)認(rèn)識,都會為邊界層高度h和動量粗糙度z0的參數(shù)化計算帶來不確定性。利用總體理查森數(shù)(Rib)確定邊界層高度是邊界層參數(shù)化方案中廣泛使用的一種方案[7-10]。例如,MRF方案[11]和 YSU方案[12]中都使用了總體理查森數(shù)方法。從最低模式層向上進(jìn)行尋找,取總體理查森數(shù)首次等于臨界總體理查森數(shù)Ribcr(一般取Ribcr=0.5)時所在的高度作為大氣邊界層高度h,即認(rèn)為是首次達(dá)到中性層結(jié)的高度。如果Rib=Ribcr所在的高度層出現(xiàn)在兩個模式層之間,那么還需要通過線性插值的方法進(jìn)一步確定h的值。因此,在大氣邊界層高度h的這一參數(shù)化過程中,參數(shù)化方法、線性插值等數(shù)值計算方法都會導(dǎo)致h相對于真實大氣的邊界層高度具有不確定性。假設(shè)地面風(fēng)速為0,Hong等[12]用公式(1)定義了總體理查森數(shù),表示浮力與切變力的比率,它們分別負(fù)責(zé)減少和產(chǎn)生湍流。

    ,? ?(1)

    其中為地表面層與高度z之間的總體理查森數(shù),為重力加速度,為z高度處的虛位溫,表示地表面的虛位溫,且,為最低模式層上的虛位溫,為地表面附近的額外虛溫,U(z)為z高度處的水平風(fēng)速。然而,Zhang等[13-14]的研究指出,用傳統(tǒng)的臨界總體理查森數(shù)法確定的邊界層高度可能不能準(zhǔn)確地再現(xiàn)颶風(fēng)邊界層的高度尺度。在進(jìn)行TC數(shù)值模擬過程中,如何對邊界層高度h進(jìn)行合適的參數(shù)化表達(dá)尚沒有讓人十分滿意的方案,仍在探索中。

    在動量粗糙度z0的參數(shù)化過程中,一般先通過地形靜態(tài)資料識別模式格點下方所處位置是屬于陸面還是水面,并將陸面上的動量粗糙度z0設(shè)置為0.01,水面上使用比較復(fù)雜的參數(shù)化公式計算,不同的參數(shù)化方案在水面上常使用不同的參數(shù)化公式。在對水面上方的模式格點進(jìn)行參數(shù)化時,通常需要使用摩擦速度、水平風(fēng)速、重力加速度。其中,摩擦速度是在模式積分的每一時步通過參數(shù)化計算獲得,水平風(fēng)速是通過模式每一時步的積分計算獲得,重力加速度取一個近似的常數(shù)。因此,在動量粗糙度z0的參數(shù)化過程中,不同參數(shù)化設(shè)計方法會導(dǎo)致數(shù)值模式中的z0具有不確定性。但是目前水面上的動量粗糙度z0(和曳力系數(shù)CD)的參數(shù)化尚沒有讓人十分滿意的方案,仍然在進(jìn)行探尋和摸索[15-19]。

    以往的研究多集中在不同邊界層參數(shù)化方案對臺風(fēng)模擬的影響上,然而,熱帶氣旋的模擬對邊界層參數(shù)化的敏感性除對方案變化敏感以外,還可能對參數(shù)化表達(dá)的微小變化也敏感[20-23]。其影響對不同臺風(fēng)并不完全一致,因此針對更多臺風(fēng)開展數(shù)值試驗來加深邊界層中變量的參數(shù)化對臺風(fēng)模擬影響的認(rèn)識至關(guān)重要。深入了解邊界層參數(shù)化變量表達(dá)差異在臺風(fēng)過程不同階段中的影響十分必要,可為有針對性地改進(jìn)數(shù)值預(yù)報模式的邊界層參數(shù)化方案提供參考依據(jù)。因此,分析臺風(fēng)過程中不同階段的數(shù)值模擬對邊界層參數(shù)化中變量要素的敏感性,認(rèn)清這些敏感要素是如何通過影響臺風(fēng)結(jié)構(gòu)來影響臺風(fēng)的路徑和強(qiáng)度具有重要意義。

    對1713號臺風(fēng)“天鴿”開展敏感性試驗,分別對邊界層高度h和動量粗糙度z0采用不同的參數(shù)化表達(dá),考察“天鴿”生命過程中不同階段的移動路徑和強(qiáng)度的數(shù)值模擬結(jié)果對邊界層高度h和動量粗糙度z0的敏感性,為理解邊界層高度h和動量粗糙度z0參數(shù)化的差異對臺風(fēng)不同階段數(shù)值模擬的影響提供參考。

    1 臺風(fēng)“天鴿”簡介

    1713號臺風(fēng)“天鴿”(Hato)于2017年8月20日14時于西北太平洋洋面上生成,此時為熱帶風(fēng)暴級別。8月22日08時加強(qiáng)為強(qiáng)熱帶風(fēng)暴,22日15時加強(qiáng)為臺風(fēng),23日07時加強(qiáng)為強(qiáng)臺風(fēng),23日15時前后以強(qiáng)臺風(fēng)級別在廣東省珠海市登陸,23日18時減弱為強(qiáng)熱帶風(fēng)暴,23日22時開始減弱為熱帶風(fēng)暴,24日14時減弱為熱帶低壓。該臺風(fēng)最強(qiáng)時臺風(fēng)中心最小海平面氣壓為935 hPa,風(fēng)力達(dá)15級(48 m/s),強(qiáng)度變化大是該臺風(fēng)的一個主要特點。該臺風(fēng)為珠海、香港、澳門等地區(qū)帶來重大破壞,造成24人死亡和68.2億美元的經(jīng)濟(jì)損失,于2018年的第59屆臺風(fēng)委員會年度會議上被除名。

    2 數(shù)值模式介紹

    本文的研究工作均是使用CMA-TYM模式系統(tǒng)開展的。CMA-TYM模式是國家氣象中心于2010年基于中尺度數(shù)值預(yù)報系統(tǒng)GRAPES-MESO(Global/Regional Assimilation and Prediction System MESOscale model)[24]開始研發(fā)的區(qū)域模式臺風(fēng)數(shù)值預(yù)報系統(tǒng),該系統(tǒng)于2012年7月開始投入業(yè)務(wù)運行并得到持續(xù)改進(jìn),CMA-TYM熱帶氣旋路徑及強(qiáng)度預(yù)報能力不斷提升。2011年實時運行結(jié)果和2012年業(yè)務(wù)運行結(jié)果表明:CMA-TYM對臺風(fēng)強(qiáng)度具有一定的預(yù)報能力,強(qiáng)度預(yù)報平均誤差小于歐洲中期天氣預(yù)報中心(European Centre for Medium-Range Weather Forecasts,ECMWF)全球模式;CMA-TYM的路徑及強(qiáng)度預(yù)報存在明顯的系統(tǒng)偏差:臺風(fēng)強(qiáng)度預(yù)報系統(tǒng)性偏強(qiáng),尤其是在預(yù)報前期存在明顯的快速增長(spin-up)現(xiàn)象,同時對于南海臺風(fēng)的路徑預(yù)報存在明顯的系統(tǒng)性北偏[25]。

    使用的CMA-TYM模式水平分辨率為0.09°×0.09°,共901×567個格點,模擬區(qū)域覆蓋西北太平洋和中國南海(90°~171°E,0°~50.94°N),垂直方向采用50層地形追隨高度坐標(biāo),使用的主要物理方案包括:WSM6(Single-Moment 6-class)微物理方案、RRTM(Rapid and Accurate Radiative Transfer Model)長波輻射方案、Dudhia短波輻射方案、Noah陸面方案、Kain-Fritsch積云參數(shù)化方案、邊界層參數(shù)化使用YSU方案,并用莫寧-奧布霍夫(Monin-Obukhov)相似理論對表面層進(jìn)行參數(shù)化。使用美國國家環(huán)境預(yù)報中心(National Centers for Environmental Prediction,NCEP)提供的每3 h一次0.5°×0.5°的零場和預(yù)報場進(jìn)行模式驅(qū)動,并根據(jù)中央氣象臺的臺風(fēng)報文,使用渦旋重定位技術(shù)和渦旋強(qiáng)度調(diào)整技術(shù)對臺風(fēng)渦旋進(jìn)行初始化調(diào)整[25-28]。

    3試驗設(shè)計

    采用中央氣象臺提供的熱帶氣旋最佳路徑數(shù)據(jù)作為實況。研究使用CMA-TYM模式對臺風(fēng)“天鴿”進(jìn)行數(shù)值模擬的過程中,選取代表“天鴿”(1713)生命過程中不同階段的四個不同時次作為起報時間,分別進(jìn)行120 h的數(shù)值模擬。

    起報時間分別選取“天鴿”(1713)生命過程中不同階段的四個時次:2017年8月20日12:00 (世界時,下同)是臺風(fēng)初始階段,為熱帶風(fēng)暴級別,中心附近最小海平面氣壓為1 000 hPa,風(fēng)力8級(18 m/s);2017年8月22日00:00剛剛升級為強(qiáng)熱帶風(fēng)暴級別,中心附近最小海平面氣壓為985 hPa,風(fēng)力10級(25 m/s);2017年8月23日00:00剛升級為強(qiáng)臺風(fēng)不久,尚未達(dá)到最大強(qiáng)度,中心附近最小海平面氣壓為950 hPa,風(fēng)力14級(42 m/s);2017年8月23日06:00即將登陸,雖然仍為強(qiáng)臺風(fēng)級別,但已經(jīng)開始減弱一段時間,中心附近最小海平面氣壓為955 hPa,風(fēng)力14級(42 m/s)。

    基于CMA-TYM模式分別對混合層參數(shù)化中的邊界層高度h和近地面層參數(shù)化中的動量粗糙度z0使用不同的參數(shù)化表達(dá)形式,考察h和z0的變化對臺風(fēng)“天鴿”不同生命階段中路徑和強(qiáng)度模擬的影響情況。

    3.1 CMA-TYM模式的參數(shù)化方案介紹

    CMA-TYM業(yè)務(wù)模式的邊界層參數(shù)化方案使用YSU方案,采用如下的方法對邊界層高度h和動量粗糙度z0進(jìn)行參數(shù)化。

    3.1.1 邊界層高度的參數(shù)化方案介紹

    YSU 方案是一階非局地K閉合方案。在YSU方案中,邊界層高度通過以下公式來推算出:

    ,? ? ?(2)

    其中 為重力加速度,U為風(fēng)速,h為邊界層頂高度,為邊界層頂?shù)奶撐粶?,為最低模式層上的虛位溫,為地表面附近的額外虛溫,為臨界總體理查森數(shù)。該定義中的邊界層為熱力學(xué)邊界層,邊界層頂部是識別的通量最小值所在的層[12]。

    3.1.2 動量粗糙度z0的參數(shù)化方案介紹

    與HWRF模式一樣,目前CMA-TYM業(yè)務(wù)模式中使用公式(3)對動量粗糙度z0進(jìn)行參數(shù)化,該方案是對應(yīng)于曳力系數(shù)CD在颶風(fēng)級風(fēng)力下保持恒定的觀點[15]開發(fā)的,并在WRF-ARW3.4.1版本中采用[29]。

    ,? (3a)

    ,? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?(3b)

    ,? ? ? ? ? ? ? ? ? ? ? ? ? ? (3c)

    ,? ? ? ? ? ? ? ? ? ? ? ? ? (3d)

    其中是摩擦速度。

    3.2 對邊界層高度的敏感性試驗方案

    公式(2)中的是一個經(jīng)驗項,具有一定的不確定性[30],這導(dǎo)致YSU 方案中定義邊界層高度h的不確定性。進(jìn)一步,YSU方案中用于計算擴(kuò)散系數(shù)的變量項kpbl(表示邊界層頂在模式面中所在的位置:第k個模式層)也是會根據(jù)調(diào)制的邊界層高度h進(jìn)行相應(yīng)的改變。因此,針對YSU方案中定義邊界層高度的這種不確定性問題,本文中通過改變一個α系數(shù),來調(diào)節(jié)邊界層高度,進(jìn)行針對邊界層高度h的一系列敏感性試驗(表1)。

    敏感性試驗方案包括1組參考試驗(CTRL)和5組敏感性試驗(表1)。這里,用H表示CTRL中用YSU 方案(即公式(2))診斷出的邊界層高度h。將H乘以一個比率系數(shù)α,通過合理地改變α的取值構(gòu)造其他5組敏感性試驗。因此,在CTRL中,α=1.00。根據(jù)Ma等[31]的研究,這里將H的比率系數(shù)取為α=0.50,得出的邊界層高度h(即h=0.50×H)視為實際邊界層高度的下限值,設(shè)計了試驗H0.50;取α=0.75,得出的邊界層高度h(即h=0.75×H)對應(yīng)邊界層高度在一個合適的范圍里,設(shè)計了試驗H0.75;取α=1.50,得出的邊界層高度h(即h=1.50×H)對應(yīng)于一個大約1 km的典型深度,設(shè)計試驗H1.50;分別取α=2.00和α=3.00,得出的邊界層高度h(即h=2.00×H和h=3.00×H)對應(yīng)某些條件下邊界層厚度能夠比1 km深得多的情況,設(shè)計試驗H2.00和試驗H3.00。Ma等[31]在文中詳細(xì)闡述了α取值的依據(jù)。

    3.3 對動量粗糙度z0的敏感性試驗方案

    在計算z0的時候,采用與公式(3)不同的其他幾種方法,進(jìn)行針對動量粗糙度z0的一系列敏感性試驗(表2)。敏感性試驗方案包括1組參考試驗(CTRL)和4組敏感性試驗(表2),4組敏感性試驗分別記為試驗z0old1.5、試驗z0old、試驗z0new和試驗z0CRTL0.5。

    (4)

    公式(4)產(chǎn)生的曳力系數(shù)CD隨風(fēng)速單調(diào)增加(對應(yīng)Green等[29]中的公式(5)),是WRF-ARW3.0版本之前就已經(jīng)采用的一種方案。試驗z0old使用公式(4)來計算粗糙度z0。

    試驗z0new采用公式(5)的方法計算z0。該方案與Green等[29]中的第一個isftcflx選項方案一致,即粗糙度z0在低風(fēng)速下是使用本文公式(3)進(jìn)行計算,但允許風(fēng)速超過33 m/s以上時曳力系數(shù)持續(xù)增加。

    ,? ? ? ? (5)

    其中z0是由公式(3)給出的。

    為了使z0在保持隨風(fēng)速不變的情況下變化幅度有更大的差別,另外設(shè)計了兩個方案z0old1.5和z0CRTL0.5。試驗z0old1.5是將試驗z0old方案得到的z0乘以系數(shù)1.5,即將試驗z0old得到的z0放大到原來的1.5倍。試驗z0CRTL0.5是將CTRL方案得到的z0乘以系數(shù)0.5,即將試驗CTRL的z0減小到原來的0.5倍。

    對于表1中的5組敏感性試驗方案,當(dāng)使用公式(6)的摩擦速度和10 m風(fēng)速的近似關(guān)系時[ 注意,實際的GRAPES_TYM模式中使用的摩擦速度與10 m風(fēng)速之間的關(guān)系比公式(6)復(fù)雜得多,是使用公式進(jìn)行參數(shù)化計算的。忽略穩(wěn)定度函數(shù)項,并取高度z=10 m,得到簡化的公式(6)。],可得到如圖1中所示的動量粗糙度z0隨10 m風(fēng)速的變化關(guān)系。

    ,? ? ? ? ? ? ? ? ? ? ? ?(6)

    其中取,。

    4 數(shù)值試驗結(jié)果分析

    4.1 邊界層高度h對臺風(fēng)路徑和強(qiáng)度模擬的影響分析

    4.1.1 對臺風(fēng)移動路徑模擬的影響分析

    圖2a、c、e、g給出了邊界層高度h的變化對臺風(fēng)“天鴿”四個不同階段起報的移動路徑模擬的影響,圖2b、d、f、h給出了對應(yīng)的移動路徑預(yù)報誤差。在初期為熱帶風(fēng)暴的階段(2017年8月20日12:00)起報的試驗中,積分約24 h之后,各試驗中的移動路徑誤差有明顯差異,邊界層高度h較薄的H0.50、H0.75和CTRL試驗中移動路徑誤差較小,在積分84 h之前基本小于150 km,但邊界層高度h較厚的H2.00和H3.00試驗中的移動路徑誤差較大,一般超過150 km(圖2a、b)。

    而在臺風(fēng)成熟之后(強(qiáng)熱帶風(fēng)暴及后期,即2017年8月22日00:00、2017年8月23日00:00和2017年8月23日06:00)起報的試驗中,邊界層高度h對臺風(fēng)移動路徑的模擬影響差異較小(圖2c—h),沒有臺風(fēng)生命史初期(2017年8月20日12:00起報)顯著。

    4.1.2 對臺風(fēng)強(qiáng)度模擬的影響分析

    圖3和圖4分別給出了邊界層高度h的變化對臺風(fēng)“天鴿”四個不同階段起報的臺風(fēng)中心附近最大10 m風(fēng)速和最小海平面氣壓模擬的影響。由圖3a和圖4a可見,在初期為熱帶風(fēng)暴的階段(2017年8月20日12:00)起報的試驗中,h在合適的典型深度范圍(H0.75、CTRL和H1.50試驗)有利于臺風(fēng)加強(qiáng)發(fā)展,h過?。℉0.50試驗)或太厚(H2.00和H3.00試驗)都不利于臺風(fēng)的發(fā)展加強(qiáng)。

    在臺風(fēng)成熟之后(強(qiáng)熱帶風(fēng)暴及后期,即2017年8月22日00:00、2017年8月23日00:00和2017年8月23日06:00)起報的試驗中,臺風(fēng)中心附近10 m風(fēng)速最大值的模擬對h的變化不敏感(圖3b—d),但臺風(fēng)中心附近最小海平面氣壓的模擬對h的變化非常敏感(圖4b—d)。圖4b表明,在臺風(fēng)達(dá)到最大強(qiáng)度之前,較淺薄的h(H0.50、H0.75和CTRL試驗)導(dǎo)致臺風(fēng)加強(qiáng)偏快、強(qiáng)度發(fā)展過強(qiáng)。

    同時發(fā)現(xiàn),對于成熟的臺風(fēng),臺風(fēng)中心附近最小海平面氣壓對h變化的反應(yīng)約在6 h以后才顯現(xiàn)(圖4b—d);而初期階段的臺風(fēng)需要更長的響應(yīng)時間(約24 h)(圖4a)。在臺風(fēng)達(dá)到最大強(qiáng)度以后的減弱階段,臺風(fēng)中心附近最小海平面氣壓的減弱對h的變化非常敏感(圖4a—d),圖4d顯示H3.00試驗中臺風(fēng)減弱最快。

    4.2 動量粗糙度z0對臺風(fēng)路徑和強(qiáng)度模擬的影響分析

    4.2.1 對臺風(fēng)移動路徑模擬的影響分析

    圖5給出了粗糙度z0的變化對臺風(fēng)“天鴿”四個不同階段起報的移動路徑模擬的影響及其對應(yīng)的預(yù)報誤差。容易看出,移動路徑的模擬對粗糙度z0不太敏感,僅z0old1.5方案在初期為熱帶風(fēng)暴的階段(2017年8月20日12:00)起報的試驗有可分辨的移動路徑預(yù)報誤差(圖5b)。

    4.2.2 對臺風(fēng)強(qiáng)度模擬的影響分析

    圖6和圖7分別給出了粗糙度z0的變化對臺風(fēng)“天鴿”4個不同階段起報的臺風(fēng)中心附近最大10 m風(fēng)速和最小海平面氣壓模擬的影響。在初期為熱帶風(fēng)暴的階段(2017年8月20日12:00)起報的試驗中,臺風(fēng)在達(dá)到最大強(qiáng)度之前的加強(qiáng)過程中,臺風(fēng)強(qiáng)度(臺風(fēng)中心附近最大10 m風(fēng)速和最小海平面氣壓)模擬對z0敏感,在臺風(fēng)減弱過程中,臺風(fēng)強(qiáng)度對z0不敏感(圖6a、7a)。由圖7a可見,z0old1.5方案在臺風(fēng)增強(qiáng)過程中產(chǎn)生的誤差最大,使臺風(fēng)增強(qiáng)加深偏快的最多。

    在臺風(fēng)成熟之后(強(qiáng)熱帶風(fēng)暴及后期,即2017年8月22日00:00、2017年8月23日00:00和2017年8月23日06:00)起報的試驗中,臺風(fēng)中心附近最小海平面氣壓對z0的變化不敏感(圖7b—d);減弱階段的臺風(fēng)中心附近最大10 m風(fēng)速對z0的變化也不敏感(圖6b—d)。加強(qiáng)階段的臺風(fēng)中心附近最大10 m風(fēng)速對z0的變化略有差異,z0old試驗的臺風(fēng)發(fā)展明顯偏慢偏弱(圖6b)。

    5 結(jié)論與討論

    使用CMA-TYM模式系統(tǒng),選擇1713號臺風(fēng)“天鴿”生命過程中4個典型的不同階段,進(jìn)行一系列敏感性試驗,通過分別改變邊界層高度h和動量粗糙度z0,研究邊界層參數(shù)化過程中的兩個關(guān)鍵要素對該臺風(fēng)移動路徑和強(qiáng)度模擬的影響。結(jié)果顯示:

    (1)邊界層高度h在初期的熱帶風(fēng)暴階段對熱帶氣旋的移動路徑有較明顯影響,在臺風(fēng)成熟之后對熱帶氣旋的移動路徑影響不顯著。反映臺風(fēng)強(qiáng)度的兩個要素,即臺風(fēng)中心附近最大10 m風(fēng)速和最小海平面氣壓,對邊界層高度h的敏感性有很大差別。臺風(fēng)中心附近最大10 m風(fēng)速對邊界層高度h的變化不敏感,而最小海平面氣壓對邊界層高度h的變化卻非常敏感。

    (2)在臺風(fēng)發(fā)展初期階段,邊界層過薄或過厚(即邊界層高度h過高或過低)都不利于臺風(fēng)強(qiáng)度的發(fā)展加強(qiáng)。這表明邊界層高度h在熱帶氣旋數(shù)值模擬中的重要性,尤其是在臺風(fēng)發(fā)展初期階段。

    (3)動量粗糙度z0的變化對臺風(fēng)“天鴿”移動路徑模擬的影響不顯著,對臺風(fēng)減弱過程中的強(qiáng)度影響也不顯著;動量粗糙度z0的變化對臺風(fēng)“天鴿”的影響主要體現(xiàn)在臺風(fēng)增強(qiáng)階段。在臺風(fēng)增強(qiáng)階段,臺風(fēng)強(qiáng)度(特別是臺風(fēng)中心附近最大10 m風(fēng)速)的加強(qiáng)過程對z0很敏感,尤其是在臺風(fēng)發(fā)展的初期階段。

    反映臺風(fēng)強(qiáng)度的臺風(fēng)中心附近最大10 m風(fēng)速和最小海平面氣壓的模擬對邊界層高度h的敏感性非常不同,這是CMA-TYM模式改進(jìn)的一個方向。未來還需要通過討論和分析大量臺風(fēng)個例來進(jìn)一步理解臺風(fēng)生命過程不同階段的模擬對邊界層高度h和動量粗糙度z0的敏感性問題,例如是否個例依賴、參數(shù)化方案依賴或模式依賴,希望未來能總結(jié)出普遍適用的結(jié)論來幫助數(shù)值模式中邊界層參數(shù)化改進(jìn)。邊界層高度對模擬的臺風(fēng)增強(qiáng)的影響具有非線性特點,在臺風(fēng)發(fā)展初期階段,邊界層高度h過薄或過厚都不利于臺風(fēng)強(qiáng)度的發(fā)展加強(qiáng),不平衡過程的作用可能是一個主要原因。基于K廓線的YSU邊界層方案本質(zhì)上是用熱力學(xué)定義邊界層的一種方案,而動力學(xué)的重要作用也是臺風(fēng)邊界層區(qū)別于其他邊界層問題的主要特點之一,臺風(fēng)的動力學(xué)特點如何在邊界層參數(shù)化模式中合理地體現(xiàn)也是未來的模式邊界層改進(jìn)研究中需要著重考慮的一個方向。

    參考文獻(xiàn):

    [1] ROSENTHAL S L. A theoretical analysis of the field of motion in the hurricane boundary layer[R]. Nat Hurr Res Proj Rep,1962,NO.56:12.

    [2] ELIASSEN A. On the Ekman layer in a circular vortex[J]. J Meteor Soc Japan,1971,49:784-789.

    [3] ELIASSEN A, LYSTAD M. The Ekman layer of a circular vortex: A numerical and theoretical study[J]. Geophys Norv,1977,31(7):1-16.

    [4]KEPERT J. The dynamics of boundary layer jets within the tropical cyclone core. Part I: Linear theory[J]. J Atmos Sci,2001,58(17): 2469-2484.

    [5] KEPERT J, WANG Y. The dynamics of boundary layer jets within the tropical cyclone core. Part II: Nonlinear enhancement[J]. J Atmos Sci,2001,58(17):2485-2501.

    [6] KEPERT J D. Choosing a boundary layer parameterization for tropical cyclone modeling[J]. Mon Wea Rev,2012,140(5):1427-1445.

    [7] TROEN I B, MAHRT L. A simple model of the atmospheric boundary layer: sensitivity to surface evaporation[J]. Bound-Layer Meteor,1986,37(1):129-148.

    [8] VOGELEZANG D H P, HOLTSLAG A A M. Evaluation and model impacts of alternative boundary-layer height formulations[J]. Bound-Layer Meteor,1996,81(3):245-269.

    [9] NOH Y, CHEON W G, HONG S Y, et al. Improvement of the K-profile model for the planetary boundary layer based on large eddy simulation data[J]. Bound-Layer Meteor,2003,107(2):401-427.

    [10] BENDER M A, GINIS I, TULEYA R, et al. The operational GFDL coupled hurricane–ocean prediction system and a summary of its performance[J]. Mon Wea Rev,2007,135(12):3965-3989.

    [11] HONG S Y, PAN H L. Nonlocal boundary layer vertical diffusion in a medium-range forecast model[J]. Mon Wea Rev,1996,124(10):2322-2339.

    [12] HONG S Y, NOH Y, DUDHIA J. A new vertical diffusion package with an explicit treatment of entrainment processes[J]. Mon Wea Rev,2006,134(9):2318-2341.

    [13] ZHANG J A, MARKS F D, MONTGOMERY M T, et al. An estimation of turbulent characteristics in the low-level region of intense hurricanes Allen (1980) and Hugo (1989)[J]. Mon Wea Rev,2010,139(5):1447-1462.

    [14] ZHANG J A, ROGERS R F, NOLAN D S, et al. On the characteristic height scales of the hurricane boundary layer[J]. Mon Wea Rev,2011,139(8),2523-2535.

    [15] POWELL M D, VICKERY P J, REINHOLD T A. Reduced drag coefficient for high wind speeds in tropical cyclones[J]. Nature,2003,422(6929):279-283.

    [16] MAKIN V K. A note on the drag of the sea surface at hurricane winds[J]. Bound-Layer Meteor,2005,115(1):169-176.

    [17] KUDRYAVTSEV V N, MAKIN V K. Aerodynamic roughness of the sea surface at high winds[J]. Bound-Layer Meteor,2007,125(2):289-303.

    [18] GOLBRAIKH E, SHTEMLER Y M. Foam input into the drag coefficient in hurricane conditions[J]. Dyn Atmos Oceans,2016,73:1–9.

    [19] DONELAN M A. On the decrease of the oceanic drag coefficient in high winds[J]. Journal of Geophysical Res Oceans,2018,123(2):1485–1501.

    [20] BRAUN S A, TAO W K. Sensitivity of high-resolution simulations of Hurricane Bob (1991) to planetary boundary layer parameterizations[J]. Mon Wea Rev,2000,128(12):3941–3961.

    [21] NOLAN D S, ZHANG J A, STERN D P. Evaluation of planetary boundary layer parameterizations in tropical cyclones by comparison of in situ observations and high-resolution simulations of Hurricane Isabel (2003). Part I: Initialization, maximum winds, and the outer-core boundary layer[J]. Mon Wea Rev,2009,137(11):3651–3674.

    [22] NOLAN D S, STERN D P, ZHANG J A. Evaluation of planetary boundary layer parameterizations in tropical cyclones by comparison of in situ observations and high-resolution simulations of Hurricane Isabel (2003). Part II: Inner-core boundary layer and eyewall structure[J]. Mon Wea Rev,2009,137(11):3675–3698.

    [23] ZHANG J A, NOLAN D S, ROGERS R F, et al. Evaluating the impact of improvements in the boundary layer parameterizations on hurricane intensity and structure forecasts in HWRF[J]. Mon Wea Rev,2015,143(8):3136–3155.

    [24] 陳德輝,沈?qū)W順.新一代數(shù)值預(yù)報系統(tǒng)GRAPES研究進(jìn)展[J].應(yīng)用氣象學(xué)報,2006,17(6):773-777.

    [25] 張進(jìn),麻素紅,陳德輝,等. GRAPES_TYM改進(jìn)及其在2013年西北太平洋和南海臺風(fēng)預(yù)報的表現(xiàn)[J].熱帶氣象學(xué)報,2017,33(1):64-73.

    [26] KURIHARA Y, BENDER M A, TULEYA R E, et al. Improvements in the GFDL hurricane prediction system[J]. Mon Wea Rev,1995,123(9):2791-2801.

    [27] KURIHARA Y, BENDER M A, ROSS R J. An initialization scheme of hurricane models by vortex specification[J]. Mon Wea Rev,1993,121(7):2030-2045.

    [28] WANG Y Q. On the bogusing of tropical cyclones in numerical models: The influence of vertical structure[J]. Meteor Atmos Phys,1998,65:153-170.

    [29] GREEN B W, ZHANG F Q. Impacts of air-sea flux parameterizations on the intensity and structure of tropical cyclones[J]. Mon Wea Rev,2013,141(7):2308–2324.

    [30] SHIN H H, HONG S Y, DUDHIA J. Impacts of the lowest model level height on the performance of planetary boundary layer parameterizations[J]. Mon Wea Rev,2012,140(2):664–682.

    [31] MA Z, FEI J, HUANG X, et al. Sensitivity of the Simulated Tropical Cyclone Intensi?cation to the Boundary-Layer Height Based on a K-pro?le boundary-Layer parameterization scheme[J]. JAMES,2018,10(11):2912–2932.

    猜你喜歡
    參數(shù)化邊界層數(shù)值模擬
    基于HIFiRE-2超燃發(fā)動機(jī)內(nèi)流道的激波邊界層干擾分析
    Pro/E的三維往復(fù)壓縮機(jī)參數(shù)化模型庫的建立
    一種懸架運動仿真快速建模方法研究
    汽車科技(2016年6期)2016-12-19 20:32:56
    跨音速飛行中機(jī)翼水汽凝結(jié)的數(shù)值模擬研究
    科技視界(2016年18期)2016-11-03 20:38:17
    姚橋煤礦采空區(qū)CO2防滅火的數(shù)值模擬分析
    雙螺桿膨脹機(jī)的流場數(shù)值模擬研究
    科技視界(2016年22期)2016-10-18 14:53:19
    一種基于液壓緩沖的減震管卡設(shè)計與性能分析
    科技視界(2016年20期)2016-09-29 11:08:27
    基于Pro/ENGINEER的圓柱齒輪參數(shù)化虛擬設(shè)計
    一類具有邊界層性質(zhì)的二次奇攝動邊值問題
    非特征邊界的MHD方程的邊界層
    国产在线男女| 一进一出好大好爽视频| av福利片在线观看| 精品久久久久久久久久久久久| 国产精品伦人一区二区| 免费看日本二区| 欧美成人免费av一区二区三区| 老司机影院成人| 一本一本综合久久| 卡戴珊不雅视频在线播放| 床上黄色一级片| 国产精品三级大全| 国产黄a三级三级三级人| 九色成人免费人妻av| 男女边吃奶边做爰视频| 国产精品免费一区二区三区在线| 香蕉av资源在线| 亚洲美女视频黄频| 国产男人的电影天堂91| 国产男人的电影天堂91| 精品午夜福利在线看| 亚洲欧美成人综合另类久久久 | 成人三级黄色视频| 波多野结衣高清无吗| 2021天堂中文幕一二区在线观| 2021天堂中文幕一二区在线观| 91av网一区二区| 此物有八面人人有两片| 久久欧美精品欧美久久欧美| 国产 一区精品| 美女免费视频网站| 国产精品av视频在线免费观看| 午夜a级毛片| 国产精品国产三级国产av玫瑰| www日本黄色视频网| 国内揄拍国产精品人妻在线| 日韩制服骚丝袜av| 亚洲精品色激情综合| a级毛片免费高清观看在线播放| 免费观看人在逋| 一本一本综合久久| 亚洲人与动物交配视频| 美女高潮的动态| 国产aⅴ精品一区二区三区波| 成人欧美大片| 国产av麻豆久久久久久久| 麻豆精品久久久久久蜜桃| 亚洲欧美日韩高清专用| 国内精品宾馆在线| 狂野欧美白嫩少妇大欣赏| 美女高潮的动态| 国产在视频线在精品| 日韩欧美国产在线观看| 精品久久久久久久久久免费视频| 黄色欧美视频在线观看| 国产精品野战在线观看| 久久亚洲国产成人精品v| 欧美一级a爱片免费观看看| 精品少妇黑人巨大在线播放 | 亚洲高清免费不卡视频| 高清毛片免费观看视频网站| 1024手机看黄色片| 久久久久性生活片| 秋霞在线观看毛片| 欧美国产日韩亚洲一区| 黄色欧美视频在线观看| 不卡视频在线观看欧美| 小说图片视频综合网站| 麻豆久久精品国产亚洲av| 国产爱豆传媒在线观看| 精品人妻视频免费看| 一a级毛片在线观看| 国产私拍福利视频在线观看| 久久久精品94久久精品| 国产精品日韩av在线免费观看| 观看美女的网站| 久久韩国三级中文字幕| 久久久久九九精品影院| 乱系列少妇在线播放| 亚洲精品456在线播放app| 最近最新中文字幕大全电影3| 18+在线观看网站| 婷婷精品国产亚洲av| av.在线天堂| 精华霜和精华液先用哪个| 五月玫瑰六月丁香| 久久久久久大精品| 一个人免费在线观看电影| 日日摸夜夜添夜夜爱| 男人的好看免费观看在线视频| 精品人妻熟女av久视频| 欧美高清性xxxxhd video| 国产在线精品亚洲第一网站| 中国美白少妇内射xxxbb| 婷婷亚洲欧美| 一级黄色大片毛片| 亚洲激情五月婷婷啪啪| 三级经典国产精品| 亚洲欧美精品自产自拍| 国产不卡一卡二| 国产欧美日韩精品亚洲av| av在线播放精品| 能在线免费观看的黄片| 成人亚洲精品av一区二区| 成人一区二区视频在线观看| 白带黄色成豆腐渣| 国产精品一区二区性色av| 五月伊人婷婷丁香| 免费电影在线观看免费观看| 久久精品综合一区二区三区| 国模一区二区三区四区视频| 亚洲国产色片| 乱系列少妇在线播放| 成人二区视频| 亚洲三级黄色毛片| 国产精品国产高清国产av| 日韩一本色道免费dvd| 在线免费观看的www视频| 99九九线精品视频在线观看视频| 国产片特级美女逼逼视频| 亚洲人成网站在线播| 在线播放无遮挡| 99久久久亚洲精品蜜臀av| 亚洲美女视频黄频| 99九九线精品视频在线观看视频| 99久国产av精品国产电影| 91精品国产九色| 人妻丰满熟妇av一区二区三区| 午夜免费男女啪啪视频观看 | 狠狠狠狠99中文字幕| 在线观看午夜福利视频| 晚上一个人看的免费电影| 亚洲在线自拍视频| 欧美日韩国产亚洲二区| 精品人妻偷拍中文字幕| 菩萨蛮人人尽说江南好唐韦庄 | 少妇的逼水好多| 国产高清不卡午夜福利| 亚洲美女黄片视频| 国产av不卡久久| 一区二区三区四区激情视频 | 成人性生交大片免费视频hd| 久久草成人影院| 久久久精品大字幕| 亚洲最大成人中文| 国内少妇人妻偷人精品xxx网站| 国产私拍福利视频在线观看| avwww免费| 亚洲精品日韩av片在线观看| 欧美日韩在线观看h| 欧美高清成人免费视频www| 51国产日韩欧美| av女优亚洲男人天堂| 成人漫画全彩无遮挡| 床上黄色一级片| 我的女老师完整版在线观看| 日本 av在线| av国产免费在线观看| 婷婷色综合大香蕉| 久久精品人妻少妇| 亚洲av二区三区四区| 亚洲欧美成人精品一区二区| av专区在线播放| 国产亚洲精品久久久久久毛片| 免费看a级黄色片| 久久久久久久午夜电影| 久久九九热精品免费| 亚洲欧美精品综合久久99| 尾随美女入室| 午夜福利在线在线| 亚洲自偷自拍三级| 欧美成人a在线观看| 晚上一个人看的免费电影| 在线播放无遮挡| 麻豆久久精品国产亚洲av| 天天一区二区日本电影三级| 男女之事视频高清在线观看| 最近最新中文字幕大全电影3| 国产高清激情床上av| 99久国产av精品| 99riav亚洲国产免费| 国产精品三级大全| 18禁在线无遮挡免费观看视频 | 22中文网久久字幕| 亚洲av.av天堂| 夜夜爽天天搞| 久久久国产成人精品二区| 免费一级毛片在线播放高清视频| 99视频精品全部免费 在线| avwww免费| 狠狠狠狠99中文字幕| 亚洲欧美精品综合久久99| 久久久欧美国产精品| 一区二区三区高清视频在线| 亚洲综合色惰| 91久久精品国产一区二区成人| 热99re8久久精品国产| 真实男女啪啪啪动态图| 亚洲人成网站在线播放欧美日韩| 少妇人妻一区二区三区视频| 黑人高潮一二区| 少妇猛男粗大的猛烈进出视频 | 我要看日韩黄色一级片| 国产成人aa在线观看| 啦啦啦韩国在线观看视频| 熟女电影av网| 男女下面进入的视频免费午夜| 夜夜爽天天搞| 午夜精品国产一区二区电影 | 午夜激情福利司机影院| 日日摸夜夜添夜夜添小说| 国产熟女欧美一区二区| 国产蜜桃级精品一区二区三区| 欧美成人a在线观看| 亚洲精品乱码久久久v下载方式| 黄色日韩在线| 日韩制服骚丝袜av| 免费人成视频x8x8入口观看| 日本黄色片子视频| 18+在线观看网站| 22中文网久久字幕| 九色成人免费人妻av| 99在线视频只有这里精品首页| 天堂网av新在线| 三级国产精品欧美在线观看| 天堂动漫精品| 午夜精品在线福利| 亚洲在线观看片| 亚洲国产高清在线一区二区三| 亚洲av.av天堂| 久久国内精品自在自线图片| 女的被弄到高潮叫床怎么办| 搞女人的毛片| 少妇丰满av| av视频在线观看入口| 亚洲精品456在线播放app| 亚洲国产精品合色在线| www日本黄色视频网| 狂野欧美激情性xxxx在线观看| 99热6这里只有精品| 中文资源天堂在线| 午夜福利18| 欧美一区二区亚洲| 国产91av在线免费观看| 亚洲人成网站在线观看播放| 国产三级在线视频| 美女xxoo啪啪120秒动态图| 免费看日本二区| 天堂av国产一区二区熟女人妻| 三级国产精品欧美在线观看| 日韩欧美一区二区三区在线观看| 色播亚洲综合网| 亚洲欧美精品综合久久99| 亚洲国产精品合色在线| 久久精品久久久久久噜噜老黄 | 日韩国内少妇激情av| av中文乱码字幕在线| 久久久精品大字幕| av在线亚洲专区| 午夜福利高清视频| 波多野结衣高清作品| 亚洲欧美中文字幕日韩二区| 午夜精品一区二区三区免费看| 欧美三级亚洲精品| 亚洲综合色惰| 日日干狠狠操夜夜爽| 91精品国产九色| 精品少妇黑人巨大在线播放 | 成熟少妇高潮喷水视频| 国产精品99久久久久久久久| 一区二区三区高清视频在线| 日韩欧美一区二区三区在线观看| 成人午夜高清在线视频| 美女高潮的动态| 内射极品少妇av片p| 婷婷色综合大香蕉| 欧美性猛交╳xxx乱大交人| 午夜激情欧美在线| 国产精品人妻久久久久久| av天堂在线播放| 亚洲性夜色夜夜综合| 中出人妻视频一区二区| 搞女人的毛片| 成人三级黄色视频| 在线观看美女被高潮喷水网站| 色视频www国产| 日韩精品中文字幕看吧| 日韩中字成人| 日本爱情动作片www.在线观看 | 国产高清视频在线播放一区| 看非洲黑人一级黄片| 嫩草影院新地址| 女的被弄到高潮叫床怎么办| 看免费成人av毛片| 搡女人真爽免费视频火全软件 | 亚洲第一电影网av| 亚洲电影在线观看av| 亚洲av.av天堂| 最近手机中文字幕大全| 精品99又大又爽又粗少妇毛片| 久久久久久久久久成人| 国产三级在线视频| a级毛片免费高清观看在线播放| 午夜精品在线福利| 亚洲国产精品久久男人天堂| 乱码一卡2卡4卡精品| 亚洲在线观看片| 欧美日本视频| 亚洲国产精品合色在线| 变态另类成人亚洲欧美熟女| 99久久久亚洲精品蜜臀av| 无遮挡黄片免费观看| 日韩欧美一区二区三区在线观看| 亚洲av美国av| 嫩草影院精品99| 国产精品久久久久久av不卡| 毛片女人毛片| 在线免费十八禁| 最新在线观看一区二区三区| 蜜桃久久精品国产亚洲av| 中出人妻视频一区二区| av在线老鸭窝| 日日摸夜夜添夜夜添小说| 少妇高潮的动态图| 精品不卡国产一区二区三区| 国产精品无大码| 日日干狠狠操夜夜爽| 国产精品永久免费网站| 在线观看66精品国产| 久久久久久国产a免费观看| 毛片女人毛片| 国产91av在线免费观看| 国产精品一区二区三区四区免费观看 | 国产成人一区二区在线| 国产真实伦视频高清在线观看| 男女之事视频高清在线观看| 美女高潮的动态| 国产亚洲精品综合一区在线观看| 我的老师免费观看完整版| 美女黄网站色视频| 亚洲婷婷狠狠爱综合网| 桃色一区二区三区在线观看| 欧美日本视频| 午夜福利高清视频| 国产av不卡久久| 午夜激情福利司机影院| 成人精品一区二区免费| 成人av一区二区三区在线看| 日韩欧美三级三区| 国产精品久久视频播放| 亚洲天堂国产精品一区在线| 亚洲经典国产精华液单| 看非洲黑人一级黄片| 69人妻影院| 极品教师在线视频| 看免费成人av毛片| 精品人妻偷拍中文字幕| 亚洲av美国av| 免费搜索国产男女视频| 如何舔出高潮| 久久鲁丝午夜福利片| 精品久久久久久久久亚洲| 18禁在线播放成人免费| 色视频www国产| 成人精品一区二区免费| 大型黄色视频在线免费观看| 亚洲内射少妇av| 最近2019中文字幕mv第一页| 99久久精品一区二区三区| 免费在线观看影片大全网站| 国产精品国产三级国产av玫瑰| 日韩高清综合在线| 亚洲国产精品久久男人天堂| 亚洲性久久影院| 久久精品国产亚洲av香蕉五月| 日韩av在线大香蕉| 亚洲av电影不卡..在线观看| 国产熟女欧美一区二区| 成人亚洲精品av一区二区| 日本在线视频免费播放| 十八禁国产超污无遮挡网站| 少妇的逼水好多| 又黄又爽又免费观看的视频| 国产精品久久久久久久电影| 亚洲精品日韩在线中文字幕 | 亚洲欧美成人精品一区二区| 三级经典国产精品| 人妻夜夜爽99麻豆av| 亚洲精品成人久久久久久| 国产一区二区在线av高清观看| 美女xxoo啪啪120秒动态图| 亚洲国产精品成人综合色| 日韩亚洲欧美综合| 色av中文字幕| 天天一区二区日本电影三级| 日韩精品青青久久久久久| 国产一区二区在线观看日韩| 免费大片18禁| 秋霞在线观看毛片| 看黄色毛片网站| 精品久久久久久久久亚洲| 日韩亚洲欧美综合| 国内精品一区二区在线观看| 黑人高潮一二区| 国产亚洲欧美98| 高清午夜精品一区二区三区 | 国产黄a三级三级三级人| 欧美丝袜亚洲另类| 在线国产一区二区在线| 日韩精品有码人妻一区| 久久午夜亚洲精品久久| 人人妻人人看人人澡| 国产高清激情床上av| 在线免费十八禁| 欧美绝顶高潮抽搐喷水| 美女黄网站色视频| 美女 人体艺术 gogo| 婷婷精品国产亚洲av| 国产一区二区在线av高清观看| 午夜福利成人在线免费观看| 色哟哟哟哟哟哟| 欧美在线一区亚洲| 啦啦啦韩国在线观看视频| 日本色播在线视频| 久久午夜亚洲精品久久| 久久久久免费精品人妻一区二区| 国产精品一及| 韩国av在线不卡| av福利片在线观看| 97在线视频观看| 亚洲欧美精品综合久久99| 色综合色国产| 免费高清视频大片| 干丝袜人妻中文字幕| 寂寞人妻少妇视频99o| 欧美一区二区亚洲| 国产色爽女视频免费观看| 又粗又爽又猛毛片免费看| 卡戴珊不雅视频在线播放| 日韩大尺度精品在线看网址| 村上凉子中文字幕在线| 亚洲精品影视一区二区三区av| 国产伦精品一区二区三区四那| 色噜噜av男人的天堂激情| 最近在线观看免费完整版| 人妻制服诱惑在线中文字幕| 一级毛片电影观看 | 国产精品久久久久久亚洲av鲁大| 一个人观看的视频www高清免费观看| 看片在线看免费视频| 国产午夜精品久久久久久一区二区三区 | 亚洲自偷自拍三级| 成人二区视频| 亚洲国产高清在线一区二区三| 91在线精品国自产拍蜜月| 国产欧美日韩精品一区二区| 国产精品久久久久久精品电影| 色视频www国产| 久久热精品热| 午夜福利视频1000在线观看| 又黄又爽又免费观看的视频| 亚洲精品国产成人久久av| 精品久久久久久久久久免费视频| av国产免费在线观看| 午夜爱爱视频在线播放| 国内精品美女久久久久久| 亚洲国产精品成人综合色| 久久久欧美国产精品| 黑人高潮一二区| 久久草成人影院| 欧美一区二区精品小视频在线| 午夜福利在线观看吧| 中文资源天堂在线| 国产老妇女一区| 大香蕉久久网| 一级a爱片免费观看的视频| 久久久久久伊人网av| 午夜免费激情av| 又黄又爽又刺激的免费视频.| 色综合色国产| 看非洲黑人一级黄片| 91在线精品国自产拍蜜月| 最新在线观看一区二区三区| 一个人看视频在线观看www免费| 日韩中字成人| 一区二区三区四区激情视频 | 国产中年淑女户外野战色| 变态另类成人亚洲欧美熟女| 听说在线观看完整版免费高清| 国产aⅴ精品一区二区三区波| 国产一区二区三区av在线 | 国产精品一区二区三区四区久久| 97热精品久久久久久| 日日摸夜夜添夜夜添av毛片| 日韩欧美三级三区| 欧美成人免费av一区二区三区| 99久久久亚洲精品蜜臀av| 国产精品久久久久久av不卡| 搡老熟女国产l中国老女人| 国产精品免费一区二区三区在线| 成熟少妇高潮喷水视频| 一进一出抽搐gif免费好疼| 日韩av在线大香蕉| 国产在线男女| 亚洲一级一片aⅴ在线观看| 精品久久久久久成人av| 欧美+亚洲+日韩+国产| 国产欧美日韩一区二区精品| 国产高清不卡午夜福利| 最近中文字幕高清免费大全6| 成人欧美大片| 久久久成人免费电影| 欧美成人免费av一区二区三区| 亚洲成人久久性| 别揉我奶头~嗯~啊~动态视频| 国产熟女欧美一区二区| 两个人视频免费观看高清| 中文资源天堂在线| 国国产精品蜜臀av免费| 欧美一级a爱片免费观看看| 免费av毛片视频| 中出人妻视频一区二区| 午夜福利18| 成年免费大片在线观看| 国产三级中文精品| 国产欧美日韩一区二区精品| 国产精品亚洲一级av第二区| 日日摸夜夜添夜夜爱| 亚洲av二区三区四区| 日本a在线网址| 欧美日本视频| a级毛片免费高清观看在线播放| 联通29元200g的流量卡| 天堂√8在线中文| 精品久久久噜噜| 欧美精品国产亚洲| 在线播放国产精品三级| 国产精品无大码| 成人鲁丝片一二三区免费| 精品人妻偷拍中文字幕| 日本撒尿小便嘘嘘汇集6| 日本与韩国留学比较| 亚洲一区高清亚洲精品| 亚洲自拍偷在线| 欧美中文日本在线观看视频| 精华霜和精华液先用哪个| 成年女人毛片免费观看观看9| 欧美+亚洲+日韩+国产| 人妻少妇偷人精品九色| 国产视频一区二区在线看| 午夜福利视频1000在线观看| 两性午夜刺激爽爽歪歪视频在线观看| 综合色av麻豆| 亚洲国产日韩欧美精品在线观看| 天堂动漫精品| 尤物成人国产欧美一区二区三区| 久久99热这里只有精品18| 色尼玛亚洲综合影院| 美女黄网站色视频| 老师上课跳d突然被开到最大视频| 热99在线观看视频| 久久久久精品国产欧美久久久| 亚洲精品色激情综合| 中文字幕久久专区| 亚洲熟妇熟女久久| 日韩欧美在线乱码| 激情 狠狠 欧美| 亚洲无线观看免费| 欧美日韩乱码在线| 一个人看的www免费观看视频| 亚洲婷婷狠狠爱综合网| 欧美日韩一区二区视频在线观看视频在线 | 日韩 亚洲 欧美在线| 亚洲激情五月婷婷啪啪| 老司机影院成人| 久久久午夜欧美精品| 日本一本二区三区精品| 亚洲国产欧美人成| 成人欧美大片| 久久久久久大精品| 欧美绝顶高潮抽搐喷水| 久久精品影院6| 亚洲人成网站高清观看| 高清日韩中文字幕在线| 黄色欧美视频在线观看| 亚洲精品色激情综合| 成人综合一区亚洲| 哪里可以看免费的av片| 日韩成人伦理影院| 亚洲av第一区精品v没综合| 久久精品国产亚洲av涩爱 | 午夜福利在线观看吧| 婷婷亚洲欧美| 一本久久中文字幕| 不卡一级毛片| 久久久国产成人免费| 一进一出抽搐gif免费好疼| 久久精品夜色国产| 国产午夜精品久久久久久一区二区三区 | 一级黄片播放器| 国产三级中文精品| 国产亚洲av嫩草精品影院| 两个人的视频大全免费| 非洲黑人性xxxx精品又粗又长| 日产精品乱码卡一卡2卡三| 成年免费大片在线观看| 一级a爱片免费观看的视频| 国产精品永久免费网站| 国产私拍福利视频在线观看| 欧美一级a爱片免费观看看| 露出奶头的视频| 国产探花极品一区二区| 毛片一级片免费看久久久久| 午夜福利在线观看免费完整高清在 | АⅤ资源中文在线天堂|