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

    Nudging資料同化對(duì)北極海冰密集度預(yù)報(bào)的改進(jìn)

    2016-06-01 06:57:52趙杰臣楊清華李明李群李春花田忠翔張林
    海洋學(xué)報(bào) 2016年5期

    趙杰臣,楊清華,李明,李群,李春花,田忠翔,張林

    (1.中國海洋大學(xué)海洋與大氣學(xué)院,山東青島266100;2.國家海洋環(huán)境預(yù)報(bào)中心國家海洋局海洋災(zāi)害預(yù)報(bào)技術(shù)研究重點(diǎn)實(shí)驗(yàn)室,北京100081;3.中國極地研究中心,上海200136)

    ?

    Nudging資料同化對(duì)北極海冰密集度預(yù)報(bào)的改進(jìn)

    趙杰臣1,2,楊清華2,李明2,李群3,李春花2,田忠翔2,張林2

    (1.中國海洋大學(xué)海洋與大氣學(xué)院,山東青島266100;2.國家海洋環(huán)境預(yù)報(bào)中心國家海洋局海洋災(zāi)害預(yù)報(bào)技術(shù)研究重點(diǎn)實(shí)驗(yàn)室,北京100081;3.中國極地研究中心,上海200136)

    摘要:北極夏季海冰的快速減少使得北極航道提前開通成為可能。為了給北極冰區(qū)船運(yùn)活動(dòng)提供及時(shí)可靠有效的海冰預(yù)報(bào)保障,急需提高海冰預(yù)報(bào)水平。本文基于麻省理工大學(xué)通用環(huán)流模式(MIT-gcm),使用牛頓松弛逼近(Nudging)資料同化方法將德國不萊梅大學(xué)的第二代先進(jìn)微波輻射成像儀(A M SR2)海冰密集度資料同化到模式中,建立了北極海冰數(shù)值預(yù)報(bào)系統(tǒng)。設(shè)計(jì)試驗(yàn)對(duì)比3種不同Nudging系數(shù)計(jì)算方案的改進(jìn)效果,結(jié)果表明選擇合適參數(shù)后,不同方案均能顯著改進(jìn)海冰密集度初始場(chǎng)。通過設(shè)計(jì)有無Nudging同化的兩組預(yù)報(bào)試驗(yàn),結(jié)合衛(wèi)星遙感海冰密集度及中國第五次北極科學(xué)考察期間“雪龍”船的走航海冰密集度觀測(cè)數(shù)據(jù),定量分析了Nudging同化方案對(duì)北極海冰密集度的24~120 h預(yù)報(bào)結(jié)果的改進(jìn)效果。結(jié)果表明,Nudging同化對(duì)120 h內(nèi)全北極海冰密集度的空間分布和移動(dòng)單點(diǎn)目標(biāo)的海冰密集度預(yù)報(bào)結(jié)果均有顯著改善;但在海冰變化很小的情況下,Nudging同化試驗(yàn)的24~120 h預(yù)報(bào)結(jié)果均劣于慣性預(yù)報(bào)結(jié)果,說明基于Nudging同化的數(shù)值預(yù)報(bào)系統(tǒng)還需進(jìn)一步提高預(yù)報(bào)技巧。

    關(guān)鍵詞:北極海冰;密集度預(yù)報(bào);Nudging;資料同化

    1 引言

    自1979年有連續(xù)海冰衛(wèi)星觀測(cè)資料以來,北極海冰的范圍、厚度和多年冰的比例都在顯著減?。?—2]。特別是北半球夏季,北極海冰范圍減少的趨勢(shì)更為明顯,以9月為例,海冰范圍的減小速率在1979-1998年間為每年約3.2×104k m2,而在1999-2010年間急劇增加到每年約15.4×104k m2[3—4]。利用ICESat衛(wèi)星資料反演得到的海冰厚度數(shù)據(jù)顯示,冬季北極海盆尺度上海冰厚度的減小速率在2003-2008年間達(dá)到每年0.17 m[5]。而融冰季節(jié)的時(shí)間長(zhǎng)度1979-2007年間增加了20 d[6]。北極海冰的這些快速變化,特別是夏季海冰范圍、密集度、厚度的持續(xù)減小,加速了北極地區(qū)科學(xué)考察和商業(yè)航運(yùn)的發(fā)展。為了更好地把握航運(yùn)機(jī)會(huì),控制風(fēng)險(xiǎn),如何及時(shí)準(zhǔn)確地為船舶提供北極海冰預(yù)報(bào)信息成為當(dāng)前亟需解決的問題。

    在此背景下,國家海洋環(huán)境預(yù)報(bào)中心近年來初步發(fā)展了北極海冰數(shù)值預(yù)報(bào)系統(tǒng)[7—8]。該系統(tǒng)基于MITgcm冰-海耦合模式,目前采用衛(wèi)星觀測(cè)海冰密集度數(shù)據(jù)直接替代為數(shù)值預(yù)報(bào)的海冰密集度初始場(chǎng)的方法。這種初始化方法可以融合所有的衛(wèi)星海冰密集度信息,較好地提高初始場(chǎng)的準(zhǔn)確度,但也導(dǎo)致出現(xiàn)初始化后計(jì)算不穩(wěn)定的問題[7]。因此,需要引入一種省時(shí)、高效、穩(wěn)定的海冰密集度資料同化方法來同化衛(wèi)星觀測(cè)的海冰密集度,進(jìn)而改進(jìn)海冰預(yù)報(bào)的密集度初始場(chǎng)。

    國際上對(duì)于通過同化海冰密集度來改進(jìn)冰-海耦合模式的預(yù)報(bào)水平已有較多研究,包括利用(局地)集合卡曼濾波方法、最優(yōu)插值方法、Nudging方法等,將SS M/I、SS MIS和歐盟海洋海冰衛(wèi)星應(yīng)用中心(OSISA F)等衛(wèi)星海冰密集度數(shù)據(jù)同化到H Y C O M、R O M S、MITgcm、Hadley氣候模式和南森中心T OP A Z等模式中[9—15]。前人的工作中,使用(局地)集合卡曼濾波方法時(shí)涉及到的程序較多,數(shù)據(jù)量較大,當(dāng)模式在每一個(gè)時(shí)間步長(zhǎng)上都進(jìn)行大數(shù)據(jù)量的計(jì)算和交換時(shí),就需要占用較多的計(jì)算資源或者消耗較長(zhǎng)的計(jì)算時(shí)間;同化使用SS M/I和SS MIS等衛(wèi)星數(shù)據(jù)時(shí),因其較低的分辨率(25 k m)會(huì)造成某些區(qū)域的信息失真。考慮到天氣尺度的海冰預(yù)報(bào)系統(tǒng)不僅需要準(zhǔn)確的預(yù)報(bào)結(jié)果,更注重快速及時(shí)的預(yù)報(bào)時(shí)效性,因此通過比較前人的同化方法,并結(jié)合預(yù)報(bào)中心實(shí)際計(jì)算資源,本文建立了基于Nudging同化方法和MIT-gcm冰-海耦合模式的北極海冰預(yù)報(bào)系統(tǒng)。選擇將分辨率更高(6.25 k m)的準(zhǔn)實(shí)時(shí)衛(wèi)星海冰密集度產(chǎn)品A M SR2同化到預(yù)報(bào)系統(tǒng)中,并通過設(shè)計(jì)對(duì)比試驗(yàn)來驗(yàn)證,Nudging同化對(duì)天氣尺度上北極海冰密集度預(yù)報(bào)的改進(jìn)作用。

    2 模式和數(shù)據(jù)

    本文所用的MITgcm冰-海耦合模式[16]是CheckPoint63k版本,海洋和海冰模式均采用Arakawa C網(wǎng)格,采用有限體積離散方法。海冰動(dòng)力學(xué)使用黏-塑迭代算法[17],海冰熱力學(xué)使用含雪蓋的“零層”方案[18]。本試驗(yàn)中設(shè)置420×384個(gè)格點(diǎn)覆蓋全北極區(qū)域,網(wǎng)格局地正交,水平網(wǎng)格較均勻,平均間距18 k m,考慮了主要的北極海峽。開邊界取在大西洋和太平洋的55°N附近[19—20],邊界條件采用月平均的ECC O2全球模擬結(jié)果,主要包括位溫、鹽度、流場(chǎng)等。海洋模式垂直分為50層,垂直分辨率為10~450 m不等,上層海洋分層更多。海冰和海洋模式采用相同的水平網(wǎng)格。地形數(shù)據(jù)采用美國國家地球物理數(shù)據(jù)中心(N G D C)的2′全球地貌數(shù)據(jù)E T O P O2。模式的積分時(shí)間步長(zhǎng)為1 200 s。

    MITgcm模式需要的大氣強(qiáng)迫場(chǎng)主要包括10 m風(fēng)速、2 m氣溫和比濕、向下長(zhǎng)波和短波輻射、降水等。本文使用的大氣強(qiáng)迫數(shù)據(jù)有兩種,控制試驗(yàn)使用的是美國國家環(huán)境預(yù)報(bào)中心(N CEP)發(fā)布的時(shí)間間隔為6 h,水平分辨率為2.5°的再分析資料,預(yù)報(bào)試驗(yàn)使用的是N CEP發(fā)布的時(shí)間間隔為6 h,水平分辨率為0.5°的全球預(yù)報(bào)系統(tǒng)(G FS)資料。

    用于Nudging同化和預(yù)報(bào)結(jié)果對(duì)比的衛(wèi)星數(shù)據(jù)是德國不萊梅大學(xué)提供的日平均A M SR2海冰密集度數(shù)據(jù),該數(shù)據(jù)基于全球變化觀測(cè)任務(wù)衛(wèi)星(G C O MW 1)上搭載的第二代先進(jìn)微波輻射成像儀(A M SR2,頻率89 G Hz),采用ASI算法[21],具有比SS M/I更高的空間分辨率(6.25 k m),同時(shí)在北半球夏季海冰表面有融池存在的情況下,A M SR2比SS M/I具有更高的準(zhǔn)確度。另外本文利用2012年中國第五次北極科學(xué)考察期間(7月25-30日和8月28-31日,圖1)“雪龍”船的走航海冰密集度觀測(cè)數(shù)據(jù)作為獨(dú)立數(shù)據(jù)來驗(yàn)證模式預(yù)報(bào)結(jié)果,走航觀測(cè)參考ASPeCt(Antarctic Sea ice Processes and Climate)的船基海冰觀測(cè)標(biāo)準(zhǔn)進(jìn)行[22],人工海冰密集度觀測(cè)方式是觀測(cè)員站在駕駛臺(tái)高處目測(cè)視線范圍1 k m內(nèi)海冰覆蓋區(qū)域所占的比例。冰區(qū)航行時(shí)“雪龍”船船速大約18 k m/h,走航觀測(cè)為半小時(shí)一次,對(duì)走航觀測(cè)的數(shù)據(jù)進(jìn)行三點(diǎn)平滑處理得到水平距離約18 k m范圍的平均走航觀測(cè)值。A M SR2產(chǎn)品空間分辨率為6.25 k m,取“雪龍”船位置處最接近的9個(gè)網(wǎng)格求平均值,得到約18 k m×18 k m范圍的平均海冰密集度值。模式的水平網(wǎng)格約為18 k m,這樣3種數(shù)據(jù)可以在相同的水平尺度上進(jìn)行比較。

    3 Nudging資料同化

    Nudging同化方法的基本的原理是,數(shù)值模式讀入日平均的衛(wèi)星觀測(cè)海冰密集度數(shù)據(jù),然后線性插值到每一個(gè)時(shí)間步長(zhǎng)上,在進(jìn)行計(jì)算時(shí),每一積分步長(zhǎng)下的海冰密集度由以下公式進(jìn)行同化:

    Nudging系數(shù)G在以往的研究中有不同的取法,Lindsay和Zhang[10]將該系數(shù)表達(dá)為海冰密集度模擬值和觀測(cè)值之差的非線性函數(shù)(表1,Nadging-A),其中σm是模擬的海冰密集度標(biāo)準(zhǔn)偏差,σo是衛(wèi)星觀測(cè)的海冰密集度標(biāo)準(zhǔn)偏差。W ang等[14]則引入海冰密集度模擬值和觀測(cè)值的標(biāo)準(zhǔn)偏差來估計(jì)該系數(shù)(表1,Nudging-B)。MITgcm模式CheckPoint63k版本中Nudging系數(shù)G由表1,Nudging-C公式計(jì)算得到。本文中τo取值為海冰熱力學(xué)積分的時(shí)間尺度與海冰密集度Nudging的時(shí)間尺度的比值,指數(shù)部分的作用是使得G在同化過程中隨著時(shí)間和空間不斷變化,以往的數(shù)值試驗(yàn)表明,變化的Nudging系數(shù)通常比常值更有效[14]。本文選取表1中3種不同的Nudging方案分別進(jìn)行一系列的敏感性試驗(yàn),確定各方案中經(jīng)驗(yàn)性參數(shù)的最佳取值,然后探討不同方案對(duì)海冰初始場(chǎng)同化結(jié)果的影響。

    圖1 第五次中國北極科學(xué)考察期間“雪龍”船航跡圖[23]。本文將7月25-30日低緯區(qū)域的航線定義為第一航段,將8月26-31日高緯區(qū)域的航線定義為第二航段,藍(lán)色和綠色分別是是7月22日和8月24日衛(wèi)星觀測(cè)到的海冰覆蓋區(qū)域Fig.1 The R/V Xuelong route during the 5th Chinese National Arctic Research Expedition(C HIN A R E).The route along the low latitude during 25th to 30th July was named Route I and the route along the high latitude during 26th to 31th August was named Route II.The blue area showed the seaice extent on 22nd July and the green area showed the seaice extent on 24th August

    表1 不同Nudging同化方案中系數(shù)G的計(jì)算方法Tab.1 The calculation formula of parameter G in different Nudging schemes

    本試驗(yàn)建立的基于的MITgcm模式和Nudging資料同化方法的北極海冰數(shù)值預(yù)報(bào)系統(tǒng)的運(yùn)行流程包括兩部分:同化模塊和預(yù)報(bào)模塊。在同化模塊,MITgcm模式使用Nudging資料同化方法將A M SR2衛(wèi)星海冰密集度數(shù)據(jù)同化到模擬結(jié)果中,得到預(yù)報(bào)部分需要的海冰密集度初始場(chǎng),然后再進(jìn)行24 h、72 h、120 h的海冰密集度預(yù)報(bào)(圖2)。

    圖2 基于MITgcm模式和Nudging資料同化的北極海冰數(shù)值預(yù)報(bào)系統(tǒng)流程圖Fig.2 The flow chart of the operational sea ice forecasting system based on MITgcm model and Nudging assimilation

    4 試驗(yàn)設(shè)計(jì)

    控制試驗(yàn)中,MITgcm模式從初始狀態(tài)開始,使用1992年的N CEP再分析資料對(duì)模式進(jìn)行10 a的循環(huán)強(qiáng)迫,使模擬的海冰和海洋狀態(tài)達(dá)到基本穩(wěn)定,再以此模擬結(jié)果為初始場(chǎng),以1992-2012年的N CEP再分析資料為大氣強(qiáng)迫場(chǎng),從1992年1月1日積分到2012年12月31日。

    預(yù)報(bào)試驗(yàn)中分別設(shè)置兩組對(duì)比試驗(yàn),來比較Nudging同化對(duì)模式海冰密集度初始場(chǎng)的改進(jìn),如表2所示。為利用“雪龍”船走航觀測(cè)數(shù)據(jù)驗(yàn)證預(yù)報(bào)結(jié)果,模式分別從2012年7月25日00時(shí)和8月26日00時(shí)刻起報(bào),進(jìn)行120 h的海冰密集度預(yù)報(bào)。以第一航段7月25日00時(shí)刻起報(bào)的情況為例,試驗(yàn)I:無Nudging同化方案的預(yù)報(bào)試驗(yàn)中,模式以控制試驗(yàn)輸出的2012年7月23日的模擬結(jié)果作為初始場(chǎng),以7 月24日G FS資料00時(shí)刻的數(shù)據(jù)驅(qū)動(dòng)模式,得到的7 月24日的模擬結(jié)果再作為初始場(chǎng),以7月24日的G FS資料中00~120 h時(shí)刻的預(yù)報(bào)數(shù)據(jù)作為大氣強(qiáng)迫場(chǎng),驅(qū)動(dòng)模式得到7月25-29日的海冰密集度預(yù)報(bào)結(jié)果。試驗(yàn)Ⅱ:有Nudging同化方案的預(yù)報(bào)試驗(yàn)中,同化的目的是得到更真實(shí)、更接近觀測(cè)的模式初始場(chǎng)。因此試驗(yàn)中運(yùn)行Nudging同化模塊,以控制試驗(yàn)輸出的2012年7月23日的模擬結(jié)果作為初始場(chǎng),以7月24日G FS資料中00時(shí)刻的大氣數(shù)據(jù)作為大氣強(qiáng)迫場(chǎng),運(yùn)行24 h,將7月24日的A M SR2海冰密集度數(shù)據(jù)同化到模式中,得到新的7月24日的模擬結(jié)果作為預(yù)報(bào)試驗(yàn)中的初始場(chǎng);然后再運(yùn)行預(yù)報(bào)模塊,使用和無Nudging方案的預(yù)報(bào)試驗(yàn)相同的大氣強(qiáng)迫場(chǎng),得到7月25-29日的海冰密集度預(yù)報(bào)結(jié)果。而第二航段8月26日00時(shí)刻起報(bào)的試驗(yàn)設(shè)置情況和第一航段一致。

    為更好地比較Nudging同化對(duì)120 h內(nèi)海冰密集度預(yù)報(bào)的改進(jìn),本文引入慣性預(yù)報(bào)來輔助評(píng)判改進(jìn)效果。對(duì)于變化緩慢的事物,可以利用其初始場(chǎng)外推的方法進(jìn)行短期預(yù)報(bào),稱為慣性預(yù)報(bào)。本文利用7月25日和8月26日的海冰密集度初始場(chǎng)外推24~120 h,分別得到兩個(gè)時(shí)段的24~120 h慣性預(yù)報(bào)結(jié)果,并與Nudging同化的預(yù)報(bào)結(jié)果進(jìn)行比較。

    表2 預(yù)報(bào)試驗(yàn)的設(shè)置情況Tab.2 The configurations of forecasting experiments with and without Nudging

    5 結(jié)果分析

    5.1Nudging資料同化對(duì)模式初始場(chǎng)的改進(jìn)

    在預(yù)報(bào)試驗(yàn)中,Nudging同化的作用是優(yōu)化MITgcm模式預(yù)報(bào)時(shí)的初始場(chǎng),通過比較無Nudging試驗(yàn)和有Nudging試驗(yàn),可得到Nudging同化對(duì)初始場(chǎng)的改進(jìn),而不同Nudging同化方案對(duì)模式初始場(chǎng)的優(yōu)化效果也不同。圖3是不同試驗(yàn)得到的2012年7 月25日的預(yù)報(bào)初始場(chǎng)。圖3a是無Nudging試驗(yàn)得到2012年7月25日的海冰密集度初始場(chǎng),與當(dāng)天的A M SR2衛(wèi)星海冰密集度數(shù)據(jù)(圖3b)比較,較好地反映了北極高緯中心區(qū)密集度高,低緯周邊區(qū)域密集度低的空間分布特征,但在85°N以北模擬的密集度為0.6~0.8,比觀測(cè)的0.8~1.0偏低。在波弗特海和巴倫支海部分邊緣區(qū)域觀測(cè)基本無冰,但模擬結(jié)果顯示有0.6~0.8海冰覆蓋;在東北航道的楚科奇海和東西伯利亞海觀測(cè)顯示有0.6~0.8海冰,但模擬結(jié)果只有0.2~0.4,明顯偏低;圖3c-e分別為采用了3種不同Nudging同化方案后得到的模式初始場(chǎng),結(jié)果顯示Nudging同化對(duì)海冰密集度的空間分布有顯著改進(jìn),85°N以北高緯區(qū)域的海冰密集度模擬值明顯增加,范圍幾乎接近于觀測(cè)結(jié)果;海冰邊緣區(qū)的范圍和密集度也和觀測(cè)符合較好;特別是本文使用的Nudging-C方案模擬結(jié)果和觀測(cè)最為符合,且在波弗特海和巴倫支海的海冰邊緣仍然保留有0.2~0.4的低密集度區(qū)域,該Nudging方案是漸進(jìn)式的優(yōu)化模擬結(jié)果,這有助于模式自身其他要素的調(diào)整,也有助于計(jì)算的穩(wěn)定性,更為合理。因此表1所列的3種不同Nudging同化方案對(duì)海冰密集度初始場(chǎng)的改進(jìn)效果均比較顯著,但Nudging-C方案更合理。

    圖3 2012年7月25日海冰密集度初始場(chǎng)Fig.3 The initial sea ice concentration field for 25th July,2012

    對(duì)于8月26日00時(shí)刻起報(bào)的個(gè)例,N udging同化同樣發(fā)揮重要作用;無N uding情況下,模式在高緯區(qū)域模擬的海冰密集度偏低,海冰范圍偏小,N udging同化有效地改進(jìn)了海冰在高緯區(qū)域的分布和密集度,其中N udging-C方案較另外兩方案更合理(圖略),在第二航段的預(yù)報(bào)試驗(yàn)中依然采用N udging-C方案。本文的預(yù)報(bào)試驗(yàn)中將采用N udging-C同化方案。

    為進(jìn)一步定量分析不同實(shí)驗(yàn)方案對(duì)海冰密集度初始場(chǎng)的改進(jìn)效果,本文以7月25日初始場(chǎng)為例計(jì)算了不同試驗(yàn)方案得到的初始場(chǎng)與A M SR2衛(wèi)星海冰密集度的偏差沿不同緯度的分布(圖4)。對(duì)于無Nudging同化方案的情況,在80°N以南區(qū)域,平均偏差小于-0.2,而均方根偏差位于0.3~0.4之間,這表明在此區(qū)域模擬的初始場(chǎng)與觀測(cè)場(chǎng)相比,偏大和偏小情況都存在,但偏小的情況多。而在80°N以北區(qū)域,平均偏差為-0.2~-0.3,均方根偏差為0.2~0.3,這表明在此區(qū)域以模擬的初始場(chǎng)以比觀測(cè)偏小的情況為主。3種同化方案得到的海冰密集度初始場(chǎng)與觀測(cè)的偏差和均方根偏差均明顯減小,說明Nudging同化使得模擬的初始場(chǎng)在各個(gè)緯度均更接近于觀測(cè)。同時(shí)Nudging-C方案的平均偏差和均方根誤差在各個(gè)緯度均小于其他兩種同化方案,特別是在80°N以北的高緯區(qū)域尤為明顯。

    圖4 不同試驗(yàn)方案得到的7月25日海冰密集度初始場(chǎng)與A M SR2數(shù)據(jù)的平均偏差(a)和均方根偏差(b)隨緯度的分布Fig.4 The average deviation(a)and root mean square error(b)calculated along latitudes for different Nudging experiments on 25th July

    5.2Nudging同化對(duì)24~120 h預(yù)報(bào)結(jié)果的改進(jìn)

    為了評(píng)估Nudging同化對(duì)MITgcm模式預(yù)報(bào)能力的改進(jìn),設(shè)計(jì)兩組試驗(yàn),如表2所示。兩組對(duì)比試驗(yàn)結(jié)果的差異均由模式中有無Nudging資料同化所引起,因此,分析兩組試驗(yàn)結(jié)果就能得到Nudging同化對(duì)北極海冰數(shù)值預(yù)報(bào)系統(tǒng)的影響。

    5.2.12012年7月25-29日

    由于沒有Nudging資料同化,試驗(yàn)I的海冰密集度初始場(chǎng)在波弗特海、東西伯利亞海、巴倫支海和阿蒙森海盆附近和觀測(cè)相比,存在較大誤差,這導(dǎo)致了120 h預(yù)報(bào)結(jié)果偏差偏大;而試驗(yàn)Ⅱ在有Nudging同化的情況下,海冰密集度在上述區(qū)域的偏差減小,120 h內(nèi)的全北極海冰密集度預(yù)報(bào)結(jié)果的平均偏差由試驗(yàn)I 的-0.14減小至-0.03,均方根誤差由0.42減小到0.26(表3)。這組試驗(yàn)表明,初始場(chǎng)準(zhǔn)確對(duì)于天氣尺度預(yù)報(bào)的重要性,而Nudging同化則有效地改善了海冰密集度初始場(chǎng),進(jìn)而提高了120 h短期預(yù)報(bào)結(jié)果的準(zhǔn)確性。

    圖5 2012年7月25日00時(shí)刻起報(bào)的試驗(yàn)I(a、b、c)和Ⅱ(d、e、f)對(duì)應(yīng)的預(yù)報(bào)結(jié)果與A M SR2衛(wèi)星觀測(cè)數(shù)據(jù)(g、h、i)的比較Fig.5 The results of forecasting experiments I(a,b,c)andⅡ(d,e,f)initialized at 00 h on 25th July,2012,compared with A M SR2 data(g,h,i)

    5.2.22012年8月26-30日

    無Nudging同化時(shí),試驗(yàn)I海冰密集度初始場(chǎng)在歐亞大陸扇區(qū)的海冰覆蓋范圍誤差較大,海冰密集度偏小,這同樣導(dǎo)致了120 h預(yù)報(bào)結(jié)果的較大偏差;而加入Nudging同化后,試驗(yàn)Ⅱ的海冰密集度和范圍在上述區(qū)域的偏差得到顯著改善,全北極海冰密集度預(yù)報(bào)結(jié)果的平均偏差由試驗(yàn)I的-0.15減小至0.03,均方根誤差由0.47減小到0.27(表3),120 h預(yù)報(bào)結(jié)果更為接近觀測(cè)。這組試驗(yàn)進(jìn)一步表明Nudging同化有助于改善海冰密集度初始場(chǎng),進(jìn)而改善24~120 h預(yù)報(bào)結(jié)果。

    圖6是全北極海冰面積的120 h預(yù)報(bào)結(jié)果和衛(wèi)星觀測(cè)數(shù)據(jù)的比較。無N udging同化時(shí),海冰的分布情況和觀測(cè)比較有較大差別(圖5,圖6)。海冰密集度大范圍的低估也導(dǎo)致24~120 h的預(yù)報(bào)海冰面積均偏低(圖7),而加入N udging同化后,預(yù)報(bào)結(jié)果的海冰分布和海冰密集度都得到很大改善,因此計(jì)算得到的海冰面積和A M SR2衛(wèi)星觀測(cè)資料亦符合較好。

    圖6 2012年8月26日00時(shí)刻起報(bào)的試驗(yàn)I(a、b、c)和Ⅱ(d、e、f)對(duì)應(yīng)的預(yù)報(bào)結(jié)果與A M SR2衛(wèi)星觀測(cè)數(shù)據(jù)(g、h、i)的比較Fig.6 The results of forecasting experiments I(a,b,c)andⅡ(d,e,f)initialized at 00h on 26th August,2012,compared with A M SR2 data(g,h,i)

    海冰的分布和密集度在短期天氣尺度上通常變化較小,因此慣性預(yù)報(bào)也有一定可行性。表3分別計(jì)算了基于7月25日00時(shí)刻和8月26日00時(shí)刻初始場(chǎng)的120 h慣性預(yù)報(bào)結(jié)果和A M SR2的平均偏差和均方根偏差。結(jié)果表明,7月25日的預(yù)報(bào)試驗(yàn)中,無Nudging同化的試驗(yàn)I劣于慣性預(yù)報(bào),但有Nudging同化的試驗(yàn)Ⅱ優(yōu)于慣性預(yù)報(bào);而8月26日的預(yù)報(bào)試驗(yàn)中,有、無Nudging同化的試驗(yàn)都劣于慣性預(yù)報(bào)。分析圖5、圖6中A M SR2海冰密集度在120 h內(nèi)的變化可以發(fā)現(xiàn),7月25-29日太平洋扇區(qū)的楚科奇海、東西伯利亞海區(qū)域的海冰密集度發(fā)生較顯著變化,因此考慮Nudging同化的預(yù)報(bào)結(jié)果優(yōu)于慣性預(yù)報(bào)結(jié)果;而8月26-30日的A M SR2海冰密集度數(shù)據(jù)顯示,海冰外緣線收縮至80°N以北,海冰密集度沒有顯著變化,因此慣性預(yù)報(bào)的結(jié)果較好。

    圖7 試驗(yàn)I(a)和Ⅱ(b)得到的120 h海冰面積結(jié)果與衛(wèi)星觀測(cè)結(jié)果的比較Fig.7 Comparisons between 120 h sea ice area forecasting results of the experiment I(a)andⅡ(b)with A M SR2 data

    為進(jìn)一步評(píng)估海冰密集度預(yù)報(bào)的準(zhǔn)確性,本文利用2012年7月25-30日和8月26-31日期間中國第五次北極科學(xué)考察期間“雪龍”船走航觀測(cè)的海冰密集度和模式結(jié)果進(jìn)行比較。

    表3 試驗(yàn)I和Ⅱ預(yù)報(bào)結(jié)果與A M SR2的比較Tab.3 Comparisons between forecasting results of experiment I andⅡand A M SR2 data

    第一航段“雪龍”船沿北極東北航道經(jīng)過楚科奇海、東西伯利亞海、拉普捷夫海和喀拉海,航線處于大陸邊緣和海冰邊緣區(qū)域(圖1)。預(yù)報(bào)試驗(yàn)結(jié)果顯示,試驗(yàn)Ⅱ中7月25日的Nudging同化對(duì)海冰密集度預(yù)報(bào)有顯著改進(jìn),密集度的平均偏差/均方根誤差從0.21/0.40提高到0.12/0.26(表4);在無Nudging的試驗(yàn)I中,7月25日模擬的航跡上的海冰密集度基本均為0,無海冰覆蓋,而走航觀測(cè)和A M SR2均顯示該段大部分區(qū)域海冰密集度為0.5以上。加入Nudging同化后,試驗(yàn)Ⅱ的25日海冰密集度有較大改善,24 h、 48 h、72 h結(jié)果與走航觀測(cè)及A M SR2衛(wèi)星資料均吻合較好。

    圖8 第一航段“雪龍”船不同時(shí)間段所處的緯度(a)和2012年7月25日00時(shí)刻起報(bào)的試驗(yàn)I和Ⅱ預(yù)報(bào)結(jié)果與“雪龍”船航線上走航觀測(cè)資料及衛(wèi)星數(shù)據(jù)的比較(b)Fig.8 The latitude of R/V Xuelong during her Route I(a)and comparisons between forecasting results of experiment I andⅡinitialized at 00h on 25 July,2012 with ship-based observations and A M SR2 data(b)

    表4 A M SR2、試驗(yàn)I和Ⅱ與走航觀測(cè)結(jié)果的比較Tab.4 Comparisons between A M R2 data,experiment I andⅡ,with ship-based observations

    第二航段由于北極海冰面積急劇減小,“雪龍”船沿高緯80°N以北航線由西向東航行,8月28日之前“雪龍”船航線處于海冰邊緣區(qū)域(圖1),28日開始“雪龍”船調(diào)整航線,直線向北進(jìn)入密集冰區(qū),至30日抵達(dá)本航次最北點(diǎn)87.6°N附近。

    預(yù)報(bào)試驗(yàn)結(jié)果顯示,試驗(yàn)Ⅱ中的8月26日的Nudging同化對(duì)海冰密集度有很大改進(jìn),密集度的平均偏差/均方根誤差從0.38/0.55減小到-0.03/ 0.21(表4);在無Nudging的試驗(yàn)I中,8月26日模擬的航跡上的海冰密集度偏高,達(dá)0.8左右,而走航觀測(cè)和A M SR2均顯示該段大部分區(qū)域海冰密集度為0.2以上;加入Nudging同化后,試驗(yàn)Ⅱ得到的26日海冰密集度預(yù)報(bào)結(jié)果有較大改善。初始場(chǎng)在海冰范圍上的較大偏差導(dǎo)致試驗(yàn)I模擬的24~120 h海冰密集度基本為0,而在試驗(yàn)Ⅱ的24~120 h預(yù)報(bào)結(jié)果與走航觀測(cè)和A M SR2衛(wèi)星資料均符合較好,并且成功模擬出28-29日期間的低密集度水道區(qū)域(圖8)。這說明加入Nudging同化的MITgcm模式可較好地預(yù)報(bào)海冰密集度的空間分布,也提高了對(duì)移動(dòng)單點(diǎn)的海冰密集度預(yù)報(bào)結(jié)果。但無Nudging同化的預(yù)報(bào)結(jié)果則與觀測(cè)結(jié)果偏差較大。

    圖9 第二航段“雪龍”船不同時(shí)間段所處的緯度(a)和2012年8月26日00時(shí)刻起報(bào)的試驗(yàn)I和Ⅱ預(yù)報(bào)結(jié)果與“雪龍”船航線上走航觀測(cè)資料及衛(wèi)星數(shù)據(jù)的比較(b)Fig.9 The latitude of R/V Xuelong during her RouteⅡ(a)and comparisons between forecasting results of experiment I andⅡinitialized at 00h on 26 August,2012 with ship-based observations and A M SR2 data(b)

    6 結(jié)論和討論

    本文預(yù)報(bào)試驗(yàn)表明,Nudging資料同化方法簡(jiǎn)單、實(shí)用、節(jié)省計(jì)算資源且時(shí)效性高,在預(yù)報(bào)中心IB M刀片機(jī)系統(tǒng)上使用32個(gè)CP U可以在2 h內(nèi)完成整個(gè)預(yù)報(bào)系統(tǒng)120 h的計(jì)算,并通過網(wǎng)站公布相關(guān)預(yù)報(bào)結(jié)果。而通過本文設(shè)計(jì)的敏感性試驗(yàn)發(fā)現(xiàn),選擇合適的參數(shù),不同的Nudging公式均可以獲得較理想的同化效果,顯著優(yōu)化預(yù)報(bào)試驗(yàn)的海冰密集度初始場(chǎng),使其更接近于衛(wèi)星觀測(cè)數(shù)據(jù)。本文選取的3種不同Nudging同化方案的對(duì)比試驗(yàn)表明,MITgcm模式自帶的Nudging-C方案在合適的調(diào)整相關(guān)參數(shù)后得到最好的海冰密集度初始場(chǎng);與觀測(cè)相比,其在低緯低海冰密集度區(qū)域和高緯高海冰密集度區(qū)域的平均偏差和均方根偏差均最小,因此Nudging-C方案應(yīng)該更適合應(yīng)用于目前的海冰數(shù)值預(yù)報(bào)系統(tǒng)。Nudging同化的改進(jìn)不僅體現(xiàn)在海冰密集度的初始場(chǎng)上,在24~120 h的預(yù)報(bào)試驗(yàn)中,海冰面積、密集度和空間分布的預(yù)報(bào)結(jié)果均得到顯著改善,這說明預(yù)報(bào)初始場(chǎng)的改進(jìn)有效的提高了短期預(yù)報(bào)的效果,因此對(duì)觀測(cè)數(shù)據(jù)的同化在短期海冰密集度的預(yù)報(bào)上有非常重要的作用。而通過和“雪龍”船走航觀測(cè)海冰密集度比較可知,針對(duì)某一移動(dòng)目標(biāo)物的單點(diǎn)(小區(qū)域)的預(yù)報(bào)準(zhǔn)確度也在有Nudging同化的試驗(yàn)中得到顯著提高。全北極空間分布和雪龍船單點(diǎn)數(shù)據(jù)對(duì)比結(jié)果,都表明Nudging資料同化顯著減小了預(yù)報(bào)誤差。針對(duì)海冰在短時(shí)間內(nèi)變化較小的特點(diǎn),本文也引入慣性預(yù)報(bào)來輔助評(píng)判同化方案的改善效果,結(jié)果表明:在8月26-30日海冰變化極小的情況下,24~120 h慣性預(yù)報(bào)優(yōu)于有Nudging同化的試驗(yàn)Ⅱ,而在7月25-29日太平洋扇區(qū)海冰發(fā)生顯著變化的情況下,24~120 h慣性預(yù)報(bào)劣于試驗(yàn)Ⅱ。這表明,基于Nudging同化的預(yù)報(bào)系統(tǒng)還需要進(jìn)一步優(yōu)化改進(jìn),以提高其在海冰變化較小情況下的預(yù)報(bào)技巧。

    現(xiàn)階段,準(zhǔn)確的海冰分布和密集度預(yù)報(bào)對(duì)于我國北極科學(xué)考察和商業(yè)航行均具有重要的現(xiàn)實(shí)意義。對(duì)于“雪龍”船航行來說,提前獲知3~5 d后的北極大范圍海冰分布有助于其規(guī)劃航線。以圖1中的去程穿越東北航道為例,雪龍船的目的是穿過順利穿過東北航道到達(dá)冰島,因此通過衛(wèi)星遙感手段獲得當(dāng)天的北極的海冰信息后,結(jié)合3~5 d的全北極海冰外緣線變化趨勢(shì),“雪龍”船可以根據(jù)其航次目的,選擇盡可能靠近海冰邊緣又可以避開密集海冰區(qū)的航線,這樣既可以節(jié)省路程,少量的浮冰也可以減少涌浪對(duì)船的晃動(dòng)。以“雪龍”船回程嘗試抵達(dá)北極點(diǎn)為例,3 ~5 d的海冰密集度預(yù)報(bào)信息可以幫助雪龍船提前獲知較低密集度區(qū)域,從而選取合適的路線向北航行,圖1中8月28日沿125°E附近的北進(jìn)路線即為“雪龍”船根據(jù)海冰預(yù)報(bào)信息選擇的位于密集冰區(qū)的冰間水道,“雪龍”船沿此水道順利航行至87°N附近。對(duì)于抗冰能力很低的商業(yè)船只,準(zhǔn)確的海冰預(yù)報(bào)更為重要。2013年8-9月中遠(yuǎn)集團(tuán)“永盛”輪從國內(nèi)出發(fā)通過北極東北航道到達(dá)歐洲,意味著我國的北極商業(yè)航行已經(jīng)拉開序幕。以“永盛”輪為例,為了船舶安全其必須在遠(yuǎn)離浮冰區(qū)的清水區(qū)航行,因此海冰外緣線的準(zhǔn)確預(yù)報(bào)至關(guān)重要。利用3~5 d的海冰密集度預(yù)報(bào)信息提前獲知海冰范圍,規(guī)劃大致的航行方向,保證船舶盡量在遠(yuǎn)離浮冰的清水區(qū)航行。由此可見,不同類型的船舶對(duì)海冰預(yù)報(bào)信息的利用角度不同,但準(zhǔn)確的預(yù)報(bào)結(jié)果是共同期待的。本試驗(yàn)中結(jié)合Nudging同化和MITgcm冰-海耦合模式的北極海冰預(yù)報(bào)系統(tǒng),可以提供北極夏季東北航道和高緯度航道的準(zhǔn)確海冰分布和密集度預(yù)報(bào),幫助北極航行船舶找到最佳的航行路線,提高航行的安全性。

    同時(shí)也要指出的是,目前的北極海冰預(yù)報(bào)系統(tǒng)還不能完全滿足需求,如較低的18 k m水平分辨率,不能獲取足夠精細(xì)的單點(diǎn)目標(biāo)的海冰信息;僅考慮了海冰密集度的資料同化無法準(zhǔn)確的預(yù)報(bào)海冰厚度等重要信息等。目前,除海冰密集度外,我們也能獲取準(zhǔn)實(shí)時(shí)的海表面溫度(如AV H R R,網(wǎng)址http://noaasis.noaa.gov/N O A ASIS/ml/avhrr.html),海冰漂移產(chǎn)品(如OSISA F,http://w w w.osi-saf.org),和冰厚產(chǎn)品(如S M OS,http://icdc.zmaw.de和Cryosat-2,http://w w w.cpom.ucl.ac.uk/csopr/seaice.html)。今后的工作中,我們將嘗試把上述海洋-海冰要素同化進(jìn)MITgcm模式中,并將水平分辨率提高到4 k m,以進(jìn)一步增強(qiáng)數(shù)值預(yù)報(bào)能力。

    致謝:感謝中國第五次北極科學(xué)考察隊(duì)走航海冰觀測(cè)隊(duì)員辛勤細(xì)致的工作及“雪龍”船全體隊(duì)員對(duì)海冰觀測(cè)工作的支持。

    參考文獻(xiàn):

    [1]Stocker T F,Qin D,Plattner G K.Climate Change 2013:The Physical Science Basis.Contribution of W orking Group I to the Fifth Assessment Report of the Intergovern mental Panel on Climate Change[M].Cambridge:Cambridge U niversity Press,2013.

    [2]Shu Q,Song Z Y,Qiao F L.Assessment of sea ice simulations in the C MIP5 models[J].The Cryosphere,2015(9):399-409.

    [3]Cavalieri D J and Parkinson C L.Arctic sea ice variability and trends,1979-2010[J].The Cryosphere,2012,6(4):881-889.

    [4]Stroeve J C,Serreze M C,H olland M M,et al.The Arctic's rapidly shrinking seaice cover:a research synthesis[J].Climatic Change,2012,110 (3/4):1005-1027.

    [5]K wok R,Cunningham G F,W ensnahan M,et al.Thinning and volu me loss of the Arctic Ocean sea ice cover:2003-2008[J].Journal of Geophysical Research,2009,114(C7):C07005.

    [6]M arkus T,Stroeve J C and Miller J.Recent changesin Arctic seaice melt onset,freezeup,and melt season length[J].Journal of Geophysical Research,2009,114(C12):C12024.

    [7]楊清華,劉驥平,張占海,等.北極海冰數(shù)值預(yù)報(bào)的初步研究——基于海冰-海洋耦合模式MITgcm的模擬試驗(yàn)[J].大氣科學(xué),2011,35(3):473-482.Yang Qinghua,Liu Jiping,Zhang Zhanhai,et al.A preliminary study of the Arctic sea ice nu mericalforecasting:Coupled sea ice-ocean modeling experiments based on MITgcm[J].Chinese Journal of Atmospheric Science,2011,35(3):473-482.

    [8]楊清華,李春花,邢建勇,等.2010年夏季北極海冰數(shù)值預(yù)報(bào)試驗(yàn)[J].極地研究,2012,24(1):87-94.Yang Qinghua,Li Chunhua,Xing Jianyong,et al.Arctic sea ice forecasting experimentsin the su m mer of 2010[J].Chinese Journal of Polar Research,2012,24(1):87-94.

    [9]Lis?ter K,Rosanova J,Evensen G.Assimilation ofice concentration in a coupled ice-ocean model,using the Ensemble Kalman filter[J].Ocean Dynamics,2003,53(4):368-388.

    [10]Lindsay R W,Zhang J.Assimilation ofice concentration in an ice-ocean model[J].Journal of Atmospheric and Oceanic Technology,2006,23 (5):742-749.

    [11]Stark J D,Ridley J,M artin M,et al.Seaice concentration and motion assimilation in a seaice-ocean model[J].Journal of Geophysical Research,2008,113(C5):C05S91.

    [12]Bertino L,Lis?ter K A.The T O P A Z monitoring and prediction system for the Atlantic and Arctic Oceans[J].Journal of Operational Oceanogra-phy,2008,1(2):15-18.

    [13]Li Qun,Zhang Zhanhai,W u H uiding.Ice concentration assimilation in a regionalice-ocean coupled model and its application in seaice forecasting [J].Advances in Polar Science,2013,24(4):258-264.

    [14]W ang K,Debernard J,Sperrevik A K et al.A combined optimalinterpolation and nudging scheme to assimilate OSISA F seaice concentration into R O M S[J].Annals of Glaciology,2013,64:8-12.

    [15]Yang Q H,Losa N S,Losch M,et al.Assimilating su m mer seaice concentration into a coupled ice-ocean model using a localized SEIK filter[J].Annals of Glaciology,2014,56(69):38-44.

    [16]M arshall J,Adcroft A,Hill C,et al.A finite-volu me,incompressible Navier Stokes model for studies of the ocean on parallel computers[J].Journal of Geophysical Research,1997,102(C3):5753-5766.

    [17]Zhang Jinlun,Hibler III W D.On an efficient nu merical method for modeling sea ice dynamics[J].Journal of Geophysical Research,1997,102 (C4):8691-8702.

    [18]Zhang J,Hibler W D,Steele III M,et al.Arctic ice-ocean modeling with and without climate restoring[J].Journal of Physical Oceanography,1998,28(2):191-217.

    [19]Losch M,M enemenlis D,Campin J M,et al.On theformulation of seaice models.Part 1:Effects of different solverimplementations and parameterizations[J].Ocean M odelling,2010,33(1):129-144.

    [20]Nguyen A T,M enemenlis D,K wok R.Arctic ice-ocean simulation with optimized model parameters:Approach and assessment[J].Journal of Geophysical Research,2011,116:C04025.

    [21]Spreen G,Kaleschke L,Heygster G.Sea ice remote sensing using A M SR-E 89 G Hz channels[J].Journal of Geophysical Research,2008,113:C02S03.

    [22]W orby A,Ian A,Vito D.Technique for making ship-based observations of antarctic sea ice thickness and characteristics[S].H obart:Australia Antarctic Division,1999.

    [23]馬德毅.中國第五次北極科學(xué)考察報(bào)告[M].北京:海洋出版社,2013.M a Deyi.The Fifth Chinese National Arctic Expedition Team Report[M].Beijing:China Ocean Press,2013.

    趙杰臣,楊清華,李明,等.Nudging資料同化對(duì)北極海冰密集度預(yù)報(bào)的改進(jìn)[J].海洋學(xué)報(bào),2016,38(5):70-82,doi:10.3969/j.issn.0253-4193.2016.05.007

    Zhao Jiechen,Yang Qinghua,Li Ming,et al.Improving Arctic seaice concentration forecasts with a Nudging data assimilation method [J].Haiyang Xuebao,2016,38(5):70-82,doi:10.3969/j.issn.0253-4193.2016.05.007

    Improving Arctic sea ice concentration forecasts with a Nudging data assimilation method

    Zhao Jiechen1,2,Yang Qinghua2,Li Ming2,Li Qun3,Li Chunhua2,Tian Zhongxiang2,Zhang Lin2
    (1.Collegeof Oceanic and Atmospheric Sciences,Ocean University of China,Qingdao 266100,China;2.Key Laboratory of Research on Marine Hazards Forecasting of State Oceanic Administration,National Marine Environmental Forecasting Center,Beijing 100081,China;3.Polar Research Institute of China,Shanghai 200136,China)

    Abstract:The rapid decrease of Arctic seaicein su m mer makes shipping in the Arctic possible.The accurate seaice forecasts are urgently required to well service the Arctic shipping activities.A nu merical Arctic forecasting system was built based on MIT generalcirculation model(MITgcm)ice-ocean coupled modeland the Nudging data assimilation method was applied into this model and assimilate the Advanced Microwave Scanning Radiometer 2 (A M SR2)sea ice concentration data.Three different kinds of Nudging assimilation schemes were firstly accessed and the results showed that all three nudging schemes can largely improve the initial sea ice concentration fields.For comparison,two forecasting experiments with and without Nudging assimilation but with the same forcing were designed to evaluate the role of nudging data assimilation.By comparing with the assimilated satellite-derived data and the ship-based insitu sea ice concentration observations,it was shown that the nudging assimilation significantly improved the 24-120 h sea ice concentration forecasts.The results showed thatimprovements occurred not only in the whole Arctic sea ice concentration forecasts,but also in the single point forecasts.The persistence forecasts performed better in 24-120 h forecast than Nudging experiments when sea ice chance little in August.

    Key words:Arctic sea ice;concentration forecast;Nudging;data assimilation

    作者簡(jiǎn)介:趙杰臣(1984—),男,山東省乳山市人,博士研究生,主要從事極地海冰觀測(cè)和預(yù)報(bào)研究。E-mail:zhaojc@n mefc.gov.cn

    基金項(xiàng)目:國家極地考察專項(xiàng)(C HIN A R E-03-01);國家自然科學(xué)基金(41376188,41206184,41406218);海洋公益性行業(yè)科研專項(xiàng)(201205007)。

    收稿日期:2015-07-14;

    修訂日期:2015-10-25。

    中圖分類號(hào):P731.32

    文獻(xiàn)標(biāo)志碼:A

    文章編號(hào):0253-4193(2016)05-0070-13

    日本wwww免费看| 精品99又大又爽又粗少妇毛片| 亚洲av福利一区| 亚洲第一区二区三区不卡| 久久久欧美国产精品| 国产综合懂色| 直男gayav资源| 在现免费观看毛片| 国产免费视频播放在线视频 | 成人av在线播放网站| 大片免费播放器 马上看| 日本av手机在线免费观看| 美女黄网站色视频| 又大又黄又爽视频免费| 纵有疾风起免费观看全集完整版 | 精品国产露脸久久av麻豆 | 国产精品久久久久久久久免| 在现免费观看毛片| 一区二区三区四区激情视频| 夫妻午夜视频| xxx大片免费视频| 国产 亚洲一区二区三区 | 神马国产精品三级电影在线观看| 最近中文字幕2019免费版| 99热这里只有是精品50| 超碰97精品在线观看| 国产精品美女特级片免费视频播放器| 春色校园在线视频观看| 高清欧美精品videossex| 大片免费播放器 马上看| 国产精品久久久久久av不卡| av在线播放精品| 国产综合精华液| 听说在线观看完整版免费高清| 狂野欧美激情性xxxx在线观看| 国产精品一区二区三区四区久久| 免费av毛片视频| 亚洲真实伦在线观看| 久久99蜜桃精品久久| 永久免费av网站大全| 色网站视频免费| 看免费成人av毛片| 久久精品国产亚洲av涩爱| 亚洲欧美成人精品一区二区| 亚洲三级黄色毛片| 精品久久久噜噜| 精品人妻一区二区三区麻豆| 国产美女午夜福利| 亚洲最大成人手机在线| 亚洲欧美日韩卡通动漫| 卡戴珊不雅视频在线播放| 性插视频无遮挡在线免费观看| 久久6这里有精品| 日日摸夜夜添夜夜爱| 日韩av在线大香蕉| 欧美zozozo另类| 欧美性感艳星| 91久久精品国产一区二区三区| 精品久久久久久电影网| 久久久成人免费电影| 国产午夜精品久久久久久一区二区三区| 午夜日本视频在线| 男人舔奶头视频| freevideosex欧美| 欧美日韩在线观看h| 亚洲欧美精品专区久久| 亚洲成人精品中文字幕电影| 免费观看a级毛片全部| 哪个播放器可以免费观看大片| 亚洲乱码一区二区免费版| 97人妻精品一区二区三区麻豆| 国产精品久久久久久久电影| 国产精品.久久久| 女人久久www免费人成看片| av福利片在线观看| 九九在线视频观看精品| 日本熟妇午夜| 亚洲精品一二三| 欧美高清成人免费视频www| 久久精品熟女亚洲av麻豆精品 | 麻豆久久精品国产亚洲av| 亚洲精品视频女| 永久免费av网站大全| 汤姆久久久久久久影院中文字幕 | 综合色av麻豆| 乱人视频在线观看| 深爱激情五月婷婷| av在线天堂中文字幕| 床上黄色一级片| 搡老妇女老女人老熟妇| 亚洲综合精品二区| 国产亚洲av嫩草精品影院| 好男人在线观看高清免费视频| 亚洲人与动物交配视频| 久久久久性生活片| 成人二区视频| 亚洲精品中文字幕在线视频 | 亚洲国产精品专区欧美| 26uuu在线亚洲综合色| 99久久精品热视频| 亚洲国产av新网站| 亚洲自偷自拍三级| 国产精品日韩av在线免费观看| 亚洲自偷自拍三级| 免费av毛片视频| 色吧在线观看| 床上黄色一级片| 亚洲av一区综合| 亚洲欧美日韩卡通动漫| 国产一级毛片在线| 日日摸夜夜添夜夜爱| 免费大片18禁| 校园人妻丝袜中文字幕| 激情 狠狠 欧美| 久久99热这里只频精品6学生| 丝袜美腿在线中文| 99久国产av精品国产电影| 日本一二三区视频观看| 男人和女人高潮做爰伦理| 亚洲av免费高清在线观看| 国产黄色视频一区二区在线观看| 禁无遮挡网站| 简卡轻食公司| 嘟嘟电影网在线观看| 久久久久久久国产电影| 蜜桃亚洲精品一区二区三区| 免费看av在线观看网站| 高清毛片免费看| 成人美女网站在线观看视频| 久久久精品免费免费高清| 日韩成人伦理影院| 免费大片黄手机在线观看| 久久久久久久久久人人人人人人| 精品酒店卫生间| 九草在线视频观看| 免费观看在线日韩| 91午夜精品亚洲一区二区三区| 22中文网久久字幕| 黄片无遮挡物在线观看| eeuss影院久久| 最近最新中文字幕免费大全7| 亚洲内射少妇av| 一夜夜www| 插阴视频在线观看视频| 乱系列少妇在线播放| 久久99热6这里只有精品| 日本黄色片子视频| 亚洲精品国产av蜜桃| 日日啪夜夜爽| 国产淫语在线视频| 精品少妇黑人巨大在线播放| 免费av不卡在线播放| 亚洲精品aⅴ在线观看| 国产精品福利在线免费观看| 亚洲成人中文字幕在线播放| av播播在线观看一区| 一本—道久久a久久精品蜜桃钙片 精品乱码久久久久久99久播 | 久久人人爽人人片av| 亚洲一区高清亚洲精品| 免费黄网站久久成人精品| 亚洲丝袜综合中文字幕| 免费av观看视频| 久久久久久久久久人人人人人人| 高清午夜精品一区二区三区| 午夜福利网站1000一区二区三区| 街头女战士在线观看网站| 亚洲av成人精品一二三区| 午夜激情久久久久久久| 色播亚洲综合网| 在现免费观看毛片| 熟女人妻精品中文字幕| 久久精品国产亚洲网站| 中文在线观看免费www的网站| kizo精华| 国产精品.久久久| 精品人妻偷拍中文字幕| eeuss影院久久| 综合色av麻豆| 2021少妇久久久久久久久久久| 如何舔出高潮| 日韩电影二区| 国产伦理片在线播放av一区| 国产免费福利视频在线观看| 国产三级在线视频| 亚洲精品,欧美精品| 日韩中字成人| 亚洲欧美日韩无卡精品| 在线免费观看的www视频| 一个人观看的视频www高清免费观看| 日韩 亚洲 欧美在线| 国产精品久久久久久久久免| 国产午夜精品一二区理论片| av福利片在线观看| 国产 一区精品| 色哟哟·www| 一级二级三级毛片免费看| 成人午夜精彩视频在线观看| 久久久久国产网址| 国产黄色免费在线视频| 久久久久久伊人网av| 最近视频中文字幕2019在线8| 男插女下体视频免费在线播放| 一级av片app| 赤兔流量卡办理| 亚洲经典国产精华液单| 日产精品乱码卡一卡2卡三| 日本猛色少妇xxxxx猛交久久| 色播亚洲综合网| 干丝袜人妻中文字幕| 日日摸夜夜添夜夜爱| 亚洲va在线va天堂va国产| 少妇熟女欧美另类| 国产综合懂色| 波野结衣二区三区在线| 国产精品无大码| 久久久精品94久久精品| 日日干狠狠操夜夜爽| 婷婷色综合www| 久久草成人影院| 亚洲人成网站在线观看播放| 亚洲色图av天堂| 国产av在哪里看| 黑人高潮一二区| 国产黄片视频在线免费观看| 亚洲欧美一区二区三区国产| 国产一区有黄有色的免费视频 | 亚洲经典国产精华液单| 91午夜精品亚洲一区二区三区| 亚洲美女视频黄频| 国产亚洲最大av| 99久久精品国产国产毛片| 99久国产av精品| 欧美97在线视频| 精品不卡国产一区二区三区| 亚洲国产精品sss在线观看| 97在线视频观看| 亚洲国产最新在线播放| 亚洲欧美一区二区三区国产| 免费黄频网站在线观看国产| 亚洲欧美日韩东京热| 亚洲av中文av极速乱| 99久久精品热视频| 亚洲欧美日韩卡通动漫| 免费黄色在线免费观看| 国产乱人偷精品视频| 亚洲av中文字字幕乱码综合| 大话2 男鬼变身卡| 亚洲一级一片aⅴ在线观看| 国产av国产精品国产| 亚洲一区高清亚洲精品| 又黄又爽又刺激的免费视频.| 在线播放无遮挡| 亚洲电影在线观看av| 免费播放大片免费观看视频在线观看| 99九九线精品视频在线观看视频| 国内精品宾馆在线| 尤物成人国产欧美一区二区三区| h日本视频在线播放| 美女脱内裤让男人舔精品视频| 久久99精品国语久久久| 亚洲国产欧美人成| 日本一二三区视频观看| 如何舔出高潮| 搞女人的毛片| 日韩av不卡免费在线播放| 一级毛片 在线播放| 久久6这里有精品| 又爽又黄无遮挡网站| 免费看不卡的av| 国产免费一级a男人的天堂| 99热网站在线观看| 免费观看无遮挡的男女| 三级国产精品片| 欧美高清成人免费视频www| 亚洲精品亚洲一区二区| 国产亚洲91精品色在线| 精品久久久精品久久久| 乱码一卡2卡4卡精品| 亚洲av电影不卡..在线观看| 伦精品一区二区三区| 国产高潮美女av| 日本免费a在线| 欧美高清性xxxxhd video| 成年版毛片免费区| 国产一区亚洲一区在线观看| 久久久色成人| 中国美白少妇内射xxxbb| 五月伊人婷婷丁香| 欧美日韩国产mv在线观看视频 | 国产极品天堂在线| 日韩av免费高清视频| 男插女下体视频免费在线播放| 乱码一卡2卡4卡精品| 亚洲经典国产精华液单| av免费观看日本| 成人午夜高清在线视频| 亚洲精品国产av蜜桃| 一个人看视频在线观看www免费| 国产亚洲5aaaaa淫片| 国产亚洲午夜精品一区二区久久 | 在现免费观看毛片| 中文字幕人妻熟人妻熟丝袜美| 看免费成人av毛片| 亚洲真实伦在线观看| 欧美一区二区亚洲| 久久精品国产亚洲网站| 伊人久久精品亚洲午夜| 看非洲黑人一级黄片| 一个人免费在线观看电影| 99热这里只有是精品在线观看| 欧美高清成人免费视频www| 男女边摸边吃奶| 美女xxoo啪啪120秒动态图| 国产精品美女特级片免费视频播放器| 国产激情偷乱视频一区二区| 国产女主播在线喷水免费视频网站 | 国产真实伦视频高清在线观看| 综合色丁香网| 如何舔出高潮| 夫妻性生交免费视频一级片| 国产精品一及| 午夜老司机福利剧场| 日本三级黄在线观看| 夜夜看夜夜爽夜夜摸| 精品午夜福利在线看| 插阴视频在线观看视频| 久久久a久久爽久久v久久| 人体艺术视频欧美日本| 国产欧美另类精品又又久久亚洲欧美| av天堂中文字幕网| 精品久久久久久久久亚洲| 免费观看无遮挡的男女| 嘟嘟电影网在线观看| 日韩,欧美,国产一区二区三区| 男人舔女人下体高潮全视频| videossex国产| 日本一本二区三区精品| 午夜福利网站1000一区二区三区| 国产国拍精品亚洲av在线观看| 成年人午夜在线观看视频 | 少妇高潮的动态图| 黄片wwwwww| 欧美xxxx黑人xx丫x性爽| 国产成人a区在线观看| 日韩av不卡免费在线播放| 少妇裸体淫交视频免费看高清| 亚洲图色成人| 高清午夜精品一区二区三区| 亚洲精品乱久久久久久| 国产精品嫩草影院av在线观看| 国产淫语在线视频| freevideosex欧美| 久久精品久久精品一区二区三区| 如何舔出高潮| 国产有黄有色有爽视频| 日本猛色少妇xxxxx猛交久久| 建设人人有责人人尽责人人享有的 | 可以在线观看毛片的网站| 亚洲激情五月婷婷啪啪| 少妇裸体淫交视频免费看高清| 亚洲美女搞黄在线观看| 国产精品.久久久| 亚洲精品色激情综合| 久久久久精品性色| 国产精品人妻久久久久久| 少妇人妻一区二区三区视频| 搡老妇女老女人老熟妇| 高清毛片免费看| 欧美xxⅹ黑人| 久久久精品欧美日韩精品| 中文字幕免费在线视频6| 你懂的网址亚洲精品在线观看| 精品久久国产蜜桃| 晚上一个人看的免费电影| 亚洲成人av在线免费| 国产欧美另类精品又又久久亚洲欧美| 水蜜桃什么品种好| 伊人久久精品亚洲午夜| 久久精品夜夜夜夜夜久久蜜豆| 在线免费十八禁| 国产探花极品一区二区| 丰满人妻一区二区三区视频av| 三级毛片av免费| 99久久九九国产精品国产免费| 国产成人精品婷婷| 久久久色成人| 尾随美女入室| 久久久久久伊人网av| 国产成人aa在线观看| 亚洲在久久综合| 日韩大片免费观看网站| 亚洲人与动物交配视频| 久久久久久久久久久免费av| 一级a做视频免费观看| 蜜桃亚洲精品一区二区三区| 欧美日韩综合久久久久久| 2018国产大陆天天弄谢| 男人狂女人下面高潮的视频| 日本-黄色视频高清免费观看| 欧美性猛交╳xxx乱大交人| 久久久久久久久久久丰满| 成人亚洲欧美一区二区av| 久久国产乱子免费精品| 久久精品综合一区二区三区| 99热网站在线观看| 国产免费一级a男人的天堂| 久久久精品94久久精品| 好男人在线观看高清免费视频| 国产精品麻豆人妻色哟哟久久 | 国产精品一区二区性色av| 波多野结衣巨乳人妻| 国产高潮美女av| 嫩草影院精品99| 精品久久久久久成人av| 精品人妻偷拍中文字幕| 一二三四中文在线观看免费高清| 欧美一级a爱片免费观看看| 欧美日韩在线观看h| 亚洲最大成人手机在线| 不卡视频在线观看欧美| 三级男女做爰猛烈吃奶摸视频| 亚洲av福利一区| 免费黄网站久久成人精品| 日韩av在线免费看完整版不卡| 69av精品久久久久久| 成人一区二区视频在线观看| av国产免费在线观看| 亚洲欧美一区二区三区国产| 九草在线视频观看| 国产国拍精品亚洲av在线观看| 国产91av在线免费观看| 日韩av在线免费看完整版不卡| 成人无遮挡网站| 国产伦一二天堂av在线观看| 亚洲熟妇中文字幕五十中出| 国产精品久久久久久久久免| 永久免费av网站大全| 亚洲aⅴ乱码一区二区在线播放| 成人欧美大片| 国产成人91sexporn| 久久精品久久久久久久性| 一级毛片黄色毛片免费观看视频| 深爱激情五月婷婷| av在线老鸭窝| 黄片无遮挡物在线观看| 国产片特级美女逼逼视频| 99久久九九国产精品国产免费| 色综合站精品国产| 国产亚洲av片在线观看秒播厂 | 偷拍熟女少妇极品色| 国产成人精品一,二区| 深夜a级毛片| 97精品久久久久久久久久精品| 国产老妇女一区| 在线观看美女被高潮喷水网站| 欧美潮喷喷水| 国内揄拍国产精品人妻在线| 国产精品国产三级国产专区5o| 国产有黄有色有爽视频| 18禁在线播放成人免费| 大又大粗又爽又黄少妇毛片口| 亚洲av成人精品一二三区| 国产精品一区www在线观看| 美女脱内裤让男人舔精品视频| 免费大片18禁| 可以在线观看毛片的网站| 九九久久精品国产亚洲av麻豆| .国产精品久久| 日本午夜av视频| 午夜爱爱视频在线播放| 精品国产一区二区三区久久久樱花 | 五月玫瑰六月丁香| av线在线观看网站| 黄片无遮挡物在线观看| 成人性生交大片免费视频hd| 全区人妻精品视频| 久久97久久精品| 亚洲精华国产精华液的使用体验| 免费电影在线观看免费观看| 蜜桃亚洲精品一区二区三区| 中文字幕av成人在线电影| 禁无遮挡网站| 国产精品一区二区在线观看99 | 欧美丝袜亚洲另类| 精品国产一区二区三区久久久樱花 | 日日干狠狠操夜夜爽| 看黄色毛片网站| 亚洲自偷自拍三级| 精品久久久久久久久久久久久| 亚洲精品日韩av片在线观看| 在线免费观看不下载黄p国产| 永久免费av网站大全| 男插女下体视频免费在线播放| 白带黄色成豆腐渣| 国产综合精华液| 亚洲丝袜综合中文字幕| 99久久中文字幕三级久久日本| 中文字幕亚洲精品专区| 91午夜精品亚洲一区二区三区| 美女cb高潮喷水在线观看| 国产白丝娇喘喷水9色精品| 视频中文字幕在线观看| 91精品一卡2卡3卡4卡| 国国产精品蜜臀av免费| videos熟女内射| 亚洲欧美成人综合另类久久久| 久热久热在线精品观看| 久久99精品国语久久久| 亚洲精品第二区| 中文天堂在线官网| 亚洲欧美清纯卡通| 亚洲成人av在线免费| 欧美zozozo另类| 亚洲欧美日韩东京热| 我要看日韩黄色一级片| 国产激情偷乱视频一区二区| 久久久国产一区二区| 老女人水多毛片| 亚洲欧洲国产日韩| 女的被弄到高潮叫床怎么办| 亚洲成人精品中文字幕电影| 日韩一区二区三区影片| 日韩av不卡免费在线播放| 国产老妇伦熟女老妇高清| 永久网站在线| 午夜精品在线福利| 欧美三级亚洲精品| 国产亚洲精品久久久com| 日本免费在线观看一区| 最近手机中文字幕大全| 插逼视频在线观看| 最近手机中文字幕大全| 国产男女超爽视频在线观看| av卡一久久| 亚洲经典国产精华液单| 一二三四中文在线观看免费高清| 伊人久久精品亚洲午夜| 深爱激情五月婷婷| 看黄色毛片网站| 久久久精品94久久精品| 午夜免费男女啪啪视频观看| 国产一区二区亚洲精品在线观看| 国产淫语在线视频| 亚洲国产欧美在线一区| 亚洲在线自拍视频| 国产精品一及| 高清在线视频一区二区三区| 色综合亚洲欧美另类图片| 少妇人妻精品综合一区二区| 国产毛片a区久久久久| 18禁在线播放成人免费| 免费看不卡的av| 激情五月婷婷亚洲| 久久精品熟女亚洲av麻豆精品 | 国产精品熟女久久久久浪| 日本免费在线观看一区| 国产精品熟女久久久久浪| 免费观看在线日韩| 五月天丁香电影| 22中文网久久字幕| 国产美女午夜福利| 日日啪夜夜撸| 欧美3d第一页| 五月伊人婷婷丁香| 99热这里只有精品一区| 99久久精品热视频| 高清毛片免费看| 99久久中文字幕三级久久日本| 99久久人妻综合| 免费av毛片视频| 中文字幕人妻熟人妻熟丝袜美| 最近最新中文字幕免费大全7| 观看免费一级毛片| 免费观看a级毛片全部| 国产精品久久久久久精品电影小说 | 人人妻人人澡人人爽人人夜夜 | 99久久人妻综合| av免费在线看不卡| 国产一区有黄有色的免费视频 | 国产老妇伦熟女老妇高清| av网站免费在线观看视频 | 精品99又大又爽又粗少妇毛片| 欧美bdsm另类| 婷婷色av中文字幕| 男女下面进入的视频免费午夜| 久久精品国产鲁丝片午夜精品| 色吧在线观看| 啦啦啦韩国在线观看视频| 亚洲av成人av| 亚洲av中文字字幕乱码综合| 亚州av有码| 欧美潮喷喷水| 成人漫画全彩无遮挡| www.av在线官网国产| 好男人视频免费观看在线| 高清午夜精品一区二区三区| 99久国产av精品国产电影| 精品一区二区三区人妻视频| 亚洲欧美一区二区三区国产| 国产真实伦视频高清在线观看| 亚洲精品国产av蜜桃| 亚洲国产精品专区欧美| 成人综合一区亚洲| 我要看日韩黄色一级片| 黄色配什么色好看| av卡一久久| 寂寞人妻少妇视频99o| 最近的中文字幕免费完整| 2021天堂中文幕一二区在线观| 久热久热在线精品观看| 精品久久久久久电影网| 精品久久久噜噜| 亚洲最大成人av| 人妻制服诱惑在线中文字幕| 欧美成人一区二区免费高清观看|