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

北極冰海耦合模式對兩種不同大氣再分析資料響應的分析

2015-01-05 06:45:04牟龍江趙進平
海洋學報 2015年11期
關鍵詞:海冰實驗

牟龍江,趙進平

(1.中國海洋大學 物理海洋教育部重點實驗室,山東 青島 266100)

北極冰海耦合模式對兩種不同大氣再分析資料響應的分析

牟龍江1,趙進平1

(1.中國海洋大學 物理海洋教育部重點實驗室,山東 青島 266100)

本文中我們比較了Climate Forecast System Reanalysis(CFSR)高分辨率的再分析數據集和低分辨率的Japanese 25-year Reanalysis Project(JRA25)再分析數據集在向下短波輻射、向下長波輻射、10 m風場、近地面氣溫、降水、濕度上的不同,發現二者差異最大的為降水數據,其次為向下短波輻射數據、向下長波輻射數據。用這兩個數據集驅動同一冰海耦合模式,CFSR強迫的海冰、北冰洋中層水和加拿大海盆溫鹽結構與實測相比有很大差距,等密度面上的地轉流速在加拿大海盆和歐亞海盆比JRA25強迫的結果高20%,同時等密度面的深度偏深、位溫偏高,在弗拉姆海峽的流通量也比海洋再分析數據Simple Ocean Data Assimilation(SODA)偏多。CFSR的向下輻射數據更加接近實測,采用此數據的敏感性實驗模擬結果與實測符合的更好。對于海冰的模擬,云量起著至關重要的作用,降水帶來的淡水通量通過影響大西洋入流水攜帶的熱量進而影響到冰區。此外,CFSR過量的降水也是二者對于北冰洋溫鹽結構、弗拉姆海峽流通量以及地轉流強度模擬產生偏差的主要原因。盡管風場的分辨率不同,在海盆尺度上對于海冰和海水溫鹽結構的影響并不大。

北冰洋;CFSR;JRA25;冰海耦合模式;加拿大海盆;弗拉姆海峽

1 引言

全球氣候變暖背景下,北極的氣候系統發生了巨大變化[1—3],北極近地面氣溫的增幅是全球平均水平的2.5倍[4—6],被稱為北極放大現象。近年來不少研究針對北極放大的原因進行了探討,對強迫來自大氣還是海洋,學界有不同的觀點[5]。云量的增多、大氣水含量的增加、來自低緯度更多的大氣熱含量以及海冰的減退被視為產生北極放大現象的主要原因。有研究指出海冰的持續減少在其中扮演了最重要的角色[4],最近有研究分析Coupled Model Intercomparison Project Phase 5(CMIP5)數據時在沒有冰雪覆蓋范圍變化的模式中同樣發現了該現象,反而指出這不是主要原因[7]。總之,北冰洋作為大氣-海洋-海冰系統的重要組成部分,在北極放大中起怎樣的作用一直是研究的熱點。

由于觀測資料的稀缺,冰海耦合模式是研究北極海冰氣系統不可或缺的工具。冰海耦合模式中海冰動力學上的一些參數,如冰間應力、冰厚分布函數等采用的是半經驗性的參數化,盡管如此大多數模式對海冰尤其是海冰密集度的模擬已相當可信。然而,冰海耦合模式中海洋模塊并不能完整重現整個海洋過程,主要原因是一些物理過程的缺失,如格陵蘭冰架淡水通量的參數化、格陵蘭海深對流的參數化、結冰析鹽參數化等。

北冰洋中海洋模式的突出問題在于對中層水的模擬存在偏差。環極邊界流占據了北冰洋150~1 000 m的中層[8],其流向由輸入位渦的正負決定[9—10]。West Spitsbergen Current(WSC)攜帶了位渦和高溫高鹽水進入北冰洋洋盆,是影響整個中層水體的主要部分。在加拿大海盆北冰洋中層水的平均溫度一般維持在0.5℃左右,近年來由于中層水的增暖[11],平均溫度可能達到0.8℃左右。而很多海洋模式由于自身的漂移,在30年積分后整體平均溫度的漂移可能都達到0.5℃,更長時間的積分會產生更大的漂移。很多模式為此在計算過程中引入了溫鹽回復邊界條件[12—13]。

除了海洋模式中物理過程的參數化之外,表面數據的輸入是模擬誤差的另一個主要來源。冰海耦合模式不僅依賴于精確的大氣強迫場資料、徑流數據,區域模式同時也依靠精確的開邊界數據。海表面由于非絕熱過程和風應力摩擦形成的垂向位渦均會影響到內區的海洋結構[14]。這兩種過程都受到大氣強迫數據的影響,因此考察大氣強迫場的數據對于海洋模式的影響相當重要。

在本文中,我們對比了日本氣象廳的Japanese 25-year Reanalysis Project(JRA25)再分析數據集和美國國家環境預報中心的Climate Forecast System Reanalysis(CFSR)高精度的再分析數據集,并分別利用這兩個數據集強迫冰海耦合模式,來研究兩種不同數據集對于冰海耦合模式的影響。本文的第二部分簡要介紹了采用的區域冰海耦合模式及再分析數據集;第三部分比較了兩個數據集的不同之處;第四部分是實驗配置;第五部分是模式結果的比較以及可能的機制探討;第六部分為結論。

2 數值模式及數據集介紹

本文中我們使用的模式來自于ECCO2全球海洋環流和氣候評估計劃。模式中海洋模塊采用正方體網格的正交曲線坐標系,可避免極點附近空間步長太小造成的計算不穩定,取包含北冰洋的面建立水平分辨率為18 km的區域模式。垂向上有50層,其中300 m以上有20層,能很好地刻畫出鹽躍層的結構,往深層垂向網格逐漸增加到450 m。差分格式采用C網格、有限體積法。水平邊界采用滑動邊界條件,底摩擦為平均流平方的函數,系數取為常數。除此之外,KPP方案用來參數化上下邊界層和對流。動量方程水平黏性項的計算利用改進的雙調和Leith方案。為增強模式的穩定性,垂向的黏性項和擴散項均采用了隱格式處理。對溫鹽方程平流項的處理使用了7階的保單調差分格式。考慮到北冰洋內部較弱的水平擴散,顯式的水平擴散系數取為0。海冰模塊同樣建立在C網格之上,在熱力學上采用了零層模式,海冰動力學上采用了黏塑性(VP)參數化方案。模式積分過程中并沒有引入任何海表面溫鹽回復以及數據約束[15],因此模式結果可用來診斷北冰洋中淡水通量及熱通量的變化。地形數據取自ETOPO2數據集,初始溫鹽場由World Ocean Atlas 2005(WOA05)數據插值得到,海冰初始場取自泛北極海冰海洋模擬同化系統(Pan-Arctic Ice-Ocean Modeling and Assimilation System)。此外,模式中使用了氣候態的徑流數據。

模式大氣強迫場的數據分別來自JRA25和CFSR。海表面通量利用塊體公式進行計算。我們選取這兩個再分析數據集中的向下短波輻射數據、向下長波輻射數據、近地面氣溫數據、10 m風速數據、降水數據、比濕數據強迫冰海耦合模式。

JRA25是日本氣象廳發布的大氣再分析數據集,最早的版本于2007年公開,其中使用了ECMWF、NCEP、NCAR、NCDC和一些其他研究機構的數據[16]。本文中我們選取采用T106截斷的高斯網格下320×160(約1.125°)的版本,它提供了每天4個6 h一次的數據。值得注意的是JRA25的同化模式中并沒有海洋模式,海洋只是給大氣模式提供了一個下墊面。

CFSR數據集是由美國國家環境預報中心發布和維護的,采用的同化模式是一個大氣、海洋、陸表、海冰全耦合系統。我們使用的是采用T382截斷的高斯網格下1 152×576(約0.312 5°)的版本[17],其同樣提供了每天4個6 h一次的數據。

3 CFSR和JRA25數據集比較

我們分別比較兩數據集在模式海洋格點上的場平均值在時間序列上的差異、氣候態平均場在空間分布的差異。另外,由于高精度的CFSR數據集包含了更多的中尺度過程,場平均值時間序列在各個譜段能量分布的差異能體現出兩個數據集刻畫不同時間尺度過程的能力。

我們定義CFSR數據集相對JRA25的偏差為:

(1)

式中,Δ為相對偏差,DCFSR為CFSR數據集的時間序列,DJRA25為JRA25數據集的時間序列,用以比較所有物理量之間的偏差。

CFSR數據集中的向下短波輻射數據與JRA25中的相對偏差為-0.165 3(圖1),這意味著JRA25中向下短波輻射數據更強。在每年夏季二者差距達到最大,相對偏差的幅度為0.068 6,體現出明顯的季節性特征。在1980-2008年間,二者每年的極大差值都超過40 W/m2,從整個時間序列來看最大的差值達到93.47 W/m2。多年平均的結果顯示(圖2a),白令海、北歐海、拉布拉多海、北大西洋東側是二者差別最大之處,JRA25要比CFSR高20 W/m2之多。由于高緯度的原因,在北冰洋中兩個數據集的差比低緯度要小,中心區二者的差值在0~10 W/m2之間,沿著陸架附近的海區二者的差在10~20 W/m2之間。

與向下短波輻射數據在JRA25數據集中較高不同,CFSR數據集中的向下長波輻射數據始終強于JRA25,二者相對偏差為0.092(圖1)。實際偏差平均值達到了21.71 W/m2,且在每年冬季相差最大,超過30 W/m2,而在每年夏季兩者的差距最小,超過5 W/m2,這與向下短波輻射數據在夏季差距大、冬季差距小正相反。多年平均的結果(圖2b)顯示CFSR中向下的長波輻射在模式整個海區均強于JRA25,在格陵蘭海、巴倫支海、拉布拉多海、白令海峽入口二者差距較大,超過25 W/m2;而在大西洋西側、北冰洋中心區差別較小,在0~20 W/m2之間。事實上,在二者差距最大的海區與海冰邊緣線的位置接近,這是由于再分析過程中模式同化造成的。在不同的同化資料中,云的參數化方案是影響向下短波、長波輻射最顯著的因素[16—18]。CFSR數據中的云量比JRA25數據多(圖3),這在向下的長短波輻射數據比較中也體現出來。

圖1 CFSR數據集相對于JRA25數據集各分量的相對偏差Fig.1 Relative deviation of different variables between CFSR and JRA25

近地面氣溫影響到感熱、潛熱通量的交換。整個平均場的時間序列反映出CFSR和JRA25數據集中的近地面氣溫整體一致,每年近地面氣溫上升、下降速度不同導致二者的相對偏差為-0.066(圖1)。冬季CFSR中的近地面氣溫高于JRA25,超過0.75℃;夏季JRA25反而高于CFSR,超過0.87℃。多年平均的近地面氣溫空間分布場(圖2c)顯示CFSR在巴倫支海北部、弗拉姆海峽、白令海峽入口、拉布拉多海近地面溫度比JRA25高1℃之多,北冰洋中JRA25均高于CFSR,靠近加拿大群島附近JRA25近地面氣溫要高2℃之多。

比濕數據影響到模式中潛熱通量的釋放。CFSR和JRA25數據集的比濕數據相對偏差為-0.032 0,差距較小,JRA25中的數據稍大于CFSR(圖1)。在弗拉姆海峽、巴倫支海北部CFSR數據中的比濕比JRA25數據大,在其他海區弱于JRA25(圖2d)。

降水和陸地徑流影響到海表面鹽度,進而影響到整個環流[12]。但降水和陸地徑流一直沒有連續可靠的數據,再分析資料中的降水數據是一個有效的參考,但并不一定是真實可靠的。為彌補降水和陸地徑流數據的缺失,解決模式長期積分造成的鹽度漂移,一般處理的方式是采用海表面鹽度回復。對模擬真實海洋的情形來說,一些實驗結果表明,海表面鹽度回復的尺度需要經驗性的調試,這種人為引入的誤差帶來的影響是需要仔細研究評估的。理想模型的研究表明,降水引起的虛鹽度通量在平衡中層水運動中起了很大的作用[19]。通過比較這兩數據集的降水數據發現,CFSR中的降水量明顯比JRA25中的降水量大(圖1),且差距有逐年增大的趨勢。二者的相對偏差達到0.574 4,且二者相對偏差的標準差也達到了0.279 0。其中,相差較大的區域集中在北大西洋、挪威海暖水流經的海區以及白令海(圖2e)。

風場是模式需要輸入的重要外強迫場,在熱力學上影響到海洋、海冰感熱、潛熱的釋放,動力學上驅動整個海洋環流場、海冰漂移場,導致海冰形變。從兩個數據集氣候態的風場矢量來看,二者均體現出波弗

圖2 CFSR與JRA25數據集中6個變量在1980-2008年間平均場的空間分布差異Fig.2 Deviations of mean states averaged from 1980 to 2008 between CFSR and JRA25

圖3 1980-2008多年月平均的云量Fig.3 Monthly mean cloud data averaged from 1980 to 2008

特高壓控制的加拿大海盆反氣旋的風場結構,位于北極濤動核心區[20]內的南下經向風以及冰島低壓控制的氣旋式的風場結構。CFSR的風速于每年夏季均強于JRA25中的風速,冬季小于JRA25,兩個數據集的相對偏差為-0.015 8(圖1),相對較小。多年平均的結果顯示(圖2f)在北冰洋除歐亞海盆部分海區外的其他區域以及拉布拉多海,JRA25的風速普遍比CFSR高,而在北歐海和北大西洋的海域,CFSR的風速更高。雖然兩個數據集中風速大小的差別不大,但CFSR高精度的資料多年平均意義上在加拿大海盆上空引入了一個氣旋式的風場,這事實上削弱了波弗特高壓。而在巴倫支海,CFSR的資料引入了反氣旋式風場(圖4)。采取常通量層雷諾應力的計算方法[21]計算出風應力和風應力旋度。對比1996年2月6日兩數據集的風應力旋度,CFSR的數據在空間上具有JRA25數據的基本特征,但更多的是體現了在JRA25無法表現的精細中尺度過程。這些不同之處對于驅動海洋模式會造成什么影響是一個值得深入探討的問題。Marshall在早期的一些研究中指出,大氣的一些過程會影響到北大西洋亞極地海區底層水的形成,進而影響全球大洋環流[22]。中尺度的天氣過程同樣會加強北歐海的深層對流的深度、頻率和范圍,使得進入該海區的熱量增多,從丹麥海峽流出更多的底層水[23]。

對CFSR和JRA25數據集的時間序列進行譜分析后,可以發現二者在時間尺度上所包含的周期性過程基本一致。為體現出各分量在相同時間尺度上的能量差異,我們得到CFSR與JRA25各分量時間序列譜能量的差(圖5)。可以看到,CFSR數據集的各個分量在半年周期左右的時間尺度上均比JRA25具有更高的強度。而在1年周期左右的時間尺度上,JRA25向下的長短波輻射數據、近地面氣溫數據和比濕數據強度更大,在風速數據上CFSR具有更強的年周期變化,降水數據反而差別不大。

總之,兩個數據集中偏差較大的幾個變量分別為降水數據、向下短波輻射數據、向下長波輻射數據。從實測數據驗證的角度來看,CFSR數據集中的向下短波、長波輻射與實測值的偏差較小,更符合實際[24],而JRA25數據集中的降水數據與實測相比差距最小[24]。

圖4 CFSR與JRA25氣候態風場矢量偏差(a)和CFSR與JRA25氣候態風場風應力旋度差(b)Fig.4 Vector of the deviation of climatological wind filed between CFSR and JRA25(a) and the deviation of climatological wind stress curl between CFSR and JRA25(b)

圖5 CFSR數據集與JRA25數據集各分量時間序列的功率譜差(0.95的置信度)Fig.5 Biases of spectrum power of different variables between CFSR and JRA25

4 數值實驗配置

我們設計了6個數值實驗,有如下兩個目的:一是為研究兩個不同再分析數據集驅動同一區域冰海耦合模式會對模擬結果產生怎樣的影響,為發展冰海耦合模式選擇更好的大氣強迫場;二是為研究兩個數

據集中相對偏差最大的幾個變量:降水、向下的短波輻射、向下長波輻射、風場數據分別對海冰模式造成怎樣的效應。

表1展示了6個數值實驗的配置方式。Exp1采用JRA25數據集的大氣強迫場數據和CS510全球模擬實驗提供的開邊界;Exp2采用CFSR數據集的大氣強迫場數據,開邊界的配置和Exp1相同;Exp3采用JRA25數據集的大氣強迫場數據和SODA 2.1.6提供的開邊界;Exp4除采用CFSR數據中高精度的10 m風場外,其他強迫場以及開邊界的配置與Exp3相同。Exp5、Exp6均采用與Exp1相同的開邊界,分別使用CFSR數據中的輻射數據和降水數據,其他強迫場來自JRA25數據集。除此之外,在上述6個實驗中,模式中所有的海洋海冰參數的設置是完全一致的。

各個實驗均從1980年開始啟動積分,到2008年底結束。前10年的計算結果為模式啟動過程,將1990年之后的數據用來分析。

表1 模式配置Tab.1 Configuration of the experiments

續表1

5 實驗結果討論

5.1 海冰模擬結果比較

2007年是模式模擬時間段內夏季海冰顯著減少的一年。Zhang等[25]指出造成該現象的主要原因:一方面由于全球變暖的大背景,另一方面源于異常的風場以及海冰與表面反照率的反饋過程。因此,2007年夏季海冰密集度的結果可以用來對比6個實驗的模擬能力。

盡管使用了不同的開邊界條件,使用JRA25強迫場的Exp1和Exp3模擬的海冰密集度相差不大,與AMSR-E衛星反演的海冰密集度產品一致性較好(圖6)。Exp2使用CFSR高分辨率大氣強迫場模擬的海冰偏多,海冰外緣線與衛星數據有較大出入,海冰密集度的分布與其他結果也有很大的偏差。這充分說明了冰海耦合模式的模擬精度與大氣再分析資料關系很大,而開邊界條件對結果的影響并不顯著。很多研究也指出了不同再分析數據集驅動下的冰海耦合模式對于海冰模擬產生的巨大差異[24,26]。

為了解兩套再分析數據中哪個物理量的影響起主要作用,我們以模擬較好的JRA25強迫場為基礎,逐一更換氣象參數場。Exp4的配置與Exp3相同,唯一改變的是更換了CFSR的風場。模擬結果表明,Exp4中高精度風場模擬的結果在海冰邊緣區與實測數據吻合得更好,在向密集冰區過渡的海區模擬的密集度偏大。Exp5更換了輻射數據,模擬的海冰范圍與Exp1具有一致的形態,海冰外緣線的分布與衛星數據出入很小。Exp6更換了降水數據,模擬的海冰范圍與實測數據相比略有偏差,但仍好于Exp2的模擬結果。雖然具有不同的向下輻射通量,但模式均反映出合理的海冰生消年變化特征,在此情況下海冰的范圍由風場決定。對比圖6中Exp1、5、6不難發現,向下輻射數據影響到密集冰區海冰密集度的模擬,而降水數據顯然影響到對于海冰范圍的模擬。

采用全部CFSR大氣場的Exp2模擬的冰厚在海冰邊緣區比其他實驗厚(0.5±0.2)m,在一些海冰密集度較大的海區如加拿大群島和格陵蘭島北部比其他實驗薄(0.6±0.2)m,空間分布上與其他實驗出入較大。只采用CFSR風場的Exp4與其他實驗的結果相對接近。Exp5模擬的加拿大群島及格陵蘭島以北厚度超過4 m的厚冰區范圍比采用JRA25數據集的Exp1的模擬的少近1/2,更加接近ICESat衛星反演的2007年10月到11月的海冰冰厚分布場。采用CFSR降水數據的Exp6模擬的厚冰區范圍介于Exp1和Exp5之間。由此可以看出,盡管除Exp2外的各個實驗對海冰范圍的模擬相差不大,但在海冰厚度的模擬上體現出各氣象參數對模式結果產生的差異。

通過圖7a所示,1990-2008年多年平均的海冰體積可以更清楚看到各個實驗中對海冰總量的模擬狀況。Exp1、Exp2、Exp3、Exp4、Exp5、Exp6中模擬的平均海冰體積分別為(1.889±0.707)×1013m3,(1.553±0.621)×1013m3,(1.891± 0.708)×1013m3,(1.905±0.710)×1013m3,(1.624±0.741)×1013m3,(1.694±0.654)×1013m3。為更清晰顯示各實驗模擬結果的區別,我們計算了各實驗相對于6個實驗平均值的距平(圖7b)。圖7b顯示出Exp1、Exp3、Exp4的模擬結果在整個年周期上比較一致,Exp5和Exp6在冬季模擬結果接近,在夏季有很大偏差。Exp2在所有實驗中模擬的海冰體積最少。盡管Exp2中模擬的海冰范圍比其他各個實驗更廣,由于模擬的冰厚較薄,海冰體積相對偏差近17.8%。

在第3節中指出,與JRA25數據集相比,CFSR向下短波輻射通量的平均值偏低,而向下長波輻射通量平均值偏高。采用CFSR輻射通量的Exp5模擬的海冰體積比采用JRA25輻射通量的Exp1、Exp3、Exp4都少,表明向下的長波輻射起的作用更大。兩種數據輻射的差異主要是云量引起的。在70°N以北的海區,CFSR中的平均云量比JRA25中高出0.24(圖3),導致了CFSR數據中較低的向下短波輻射和較高的向下長波輻射。正是由于二者云量的差異,Exp5比Exp1模擬的海冰體積偏少2.65×1012m3,表明云量對于海冰的模擬起重要的作用。

圖6 各實驗模擬的2007年8月27日到9月6日平均海冰密集度Fig.6 Simulated sea ice concentration of each experiment averaged from August 27,2007 to September 6,2007紫線為由Bremen大學反演AMSR-E衛星海冰密集度產品所得到的海冰外緣線Sea ice edge retrieved from AMSR-E by University of Bremen is plotted in purple

圖7 1990-2008多年平均的海冰體積(a)及其距平(b) Fig.7 Annual mean sea ice volume averaged from 1990 to 2008 in each experiment (a) and anomaly (b)

如第3節所述,CFSR和JRA25數據集的降水數據存在巨大的差異。降水的增多并沒有使得Exp6模擬的海冰體積增多,反而采用CFSR降水數據的Exp6模擬的海冰體積比Exp1的模擬結果偏少1.95×1012m3。之所以出現該差異,是由于過多降水使得低緯海區的熱量釋放受到強層結的約束,由大西洋水攜帶更多的熱量進入北冰洋后釋放,減弱了海冰的生長。雖然夏季降水的增多使得更多的熱量用來融化降雪,不利于海冰的融化,而從長期積分來看降水引入的淡水通量造成的海洋升溫對海冰形成的間接影響遠超過淡水直接的影響。

Exp4與Exp1、Exp3中模擬的海冰體積、分布范圍相對一致,從另一個角度說明兩數據集中的風場差異對于海冰凍結或融化影響不大。在海冰動力學上,海冰漂流速度不僅取決于風場施加于冰面的應力,冰間應力是一個很大的量。黏塑性流變學方案(VP)中冰間應力的參數化取決于海冰密集度和冰厚。這種參數化方案決定了海冰邊緣區的冰漂流顯著受制于風場的分辨率,而在密集度較大且厚度較大的冰區,風場分辨率造成的影響并不顯著,6個實驗比較下不難得出該結論。我們比較了2007年3月和9月平均海冰漂流場,除Exp2之外的5個實驗中海冰漂流的整體模態是一致的,且與美國冰雪數據中心(NSIDC)發布的海冰漂流場數據可比,都反映出波弗特渦和穿極流控制的流場漂移形態,僅在海冰密集度差別較大的海區漂流速度有一定偏差。Exp2由于模擬海冰密集度與其他3個實驗有所出入,3月份海冰漂流場在加拿大海盆靠近阿拉斯加一側與其他實驗不一致,9月份同樣在密集度較小的海冰邊緣區有很大偏差。上述分析說明風場分辨率的提高僅影響到對海冰邊緣區冰漂流以及海冰外緣線的模擬,并沒有顯著影響到整體的海冰漂流場以及海冰體積的模擬。

5.2 海洋模擬結果比較

5.2.1 北冰洋中層水的模擬

第一斜壓Rossby變形半徑在北冰洋陸架海區大約為1~7 km,在南森海盆約為5 km,至加拿大海盆逐步增加到約15 km[27],因此我們使用的分辨率為約18 km的海洋模式并不具備分辨所有渦的能力。然而,模式中仍然可以模擬出較大尺度溫鹽和流場特征[15]。北冰洋上層環流一方面受制于太平洋入流水較高的海表面高度而形成穿極流,另一方面受制于常年存在的波弗特高壓形成典型的反氣旋式風生環流。中層水的模擬是北極冰海耦合模式中比較棘手的問題,AOMIP的各個模式均有很大的偏差[28]。海洋再分析資料也不令人滿意,同化模式中由于觀測資料的引入會產生中層水的突然升溫,不采用同化模式的再分析資料中與觀測值具有很大偏差[29]。不僅如此,北冰洋模式長期積分實驗中溫度、鹽度會有較大漂移,而中層水對此更加敏感。我們對比各個實驗在2006年12月到2008年11月模擬的450 m水深平均位溫來驗證模式對中層水的模擬效果。

對比圖8a全球海洋數據集(World Ocean Database,WOD)實測結果,長期積分實驗中,使用更接近實測輻射數據的Exp5模擬的450 m深中層水與實測符合的最好。大西洋水流經弗拉姆海峽的一支(FSB)是中層水的重要來源,Exp5模擬出該支流在歐亞海盆暖水范圍和溫度與實測數據相同的大致分布。在拉普捷夫海、東西伯利亞海、楚科奇海臺外側由于觀測數據的稀少和缺失,模式的結果無法得到有效驗證。在加拿大海盆以及加拿大群島、格陵蘭島以北,Exp5的模擬結果很好的展示出實測溫度的分布。這體現出云量通過影響向下輻射進而對北冰洋中層水模擬產生的顯著影響。

Exp1中由于接受了更多的向下短波輻射,模擬的中層水位溫在加拿大海盆的環極邊界流路徑上比Exp5高出0.3℃,在歐亞海盆暖水具有更大的范圍。直接比較Exp1和Exp2的模擬結果可以發現兩個數據集不僅對海冰的模擬有很大的誤差,在450 m深的中層,CFSR數據集模擬的位于歐亞海盆的暖水已擴散到加拿大海盆的大部分海區,使得加拿大海盆大部分位溫超過1.3℃。由北冰洋流出的東格陵蘭流水溫因此上升,影響到格陵蘭海的低溫環境,格陵蘭海氣旋式環流中心區在該深度的位溫提高了將近1℃。

盡管Exp1和Exp3模擬的海冰體積幾乎一致,但不同開邊界條件輸入熱通量的不一致導致采用SODA開邊界的模擬結果與Exp1相比超過1℃的暖水范圍更廣。Exp3模擬的暖水在歐亞海盆與Exp1的分布一致,但在楚科奇海臺外沿以及阿拉斯加沿岸海區高出0.3℃。

在較強降水的Exp6模擬結果中,整個北冰洋被暖水覆蓋,與實測產生了巨大的偏差。產生上述結果的原因是由于強降水影響到模式中的環流結構,CFSR強降水在海表面易形成障礙層,減弱了垂向對流,使得海洋向大氣的熱釋放受到抑制,過高的熱量使由弗拉姆海峽進入北冰洋的這一支水溫度過高、厚度偏大。我們的另一個對比實驗中,保持其他強迫場不變,加倍了JRA25數據集中的降水過程,結果與Exp6一致,流經弗拉姆海峽的大西洋入流水的溫度和厚度顯著升高和增大,從而影響到北冰洋內部的溫鹽性質,這支持了我們的上述解釋,充分說明北冰洋冰海耦合模式對于表層淡水通量的敏感性。

圖8 2006年11月至2008年11月WOD實測450 m深位溫平均結果(a)以及各個實驗模擬450 m深位溫結果(b)Fig.8 Mean potential temperature in 450 m averaged from November 2006 to November 2008 from WOD and each experiment

圖9 加拿大海盆中心區位溫垂直剖面(a)以及T-S圖(b)Fig.9 Vertical profiles of potential temperature (a) and T-S diagrams (b) in Canada basinCHIN(黑色)為中國第三次北極科學考察在2008年8月14日位于75.992°N,147.989°W的實測數據;WOD(淺灰色)為WOD中ITP在2007年5月17日位于75.979°N,148.553°W的實測數據;各實驗模擬的結果是2006年12月至2008年11月的平均值,位于75.972°N,148.369°WThe black curve indicated by CHIN shows CTD data got during the third Chinese National Arctic Research Expedition on August 14,2008 at 75.992°N,147.989°W,the dashed grey curve indicated by WOD is observed at 75.979°N,148.553°W on May 17,2007,other curves at 75.972°N,148.369°W averaged between December 2006 and November 2008 are derived from each experiment

5.2.2 加拿大海盆溫鹽結構模擬

加拿大海盆受波弗特高壓直接控制,由于海冰融化凍結、陸地徑流、太平洋入流水和大西洋水互相混合形成了復雜的溫鹽結構,可檢驗模式模擬能力。圖9中,我們選取了各個模式中位于加拿大海盆75.972°N,148.369°W的網格點,再選取與之相鄰指示夏季、冬季加拿大海盆的實際觀測數據,分別來自中國第三次北極科學考察(CHIN)和WOD中的冰基錨系浮標(ITP)提供的數據。CHIN垂直剖面和T-S圖顯示出夏季加拿大海盆海冰融化形成的淡水和陸地徑流在表層混合形成的高溫低鹽水,深度在100 m左右的夏季太平洋水(鹽度31~33,位溫-0.55℃),水深180 m左右的冬季太平洋水以及中層高溫高鹽的大西洋水(鹽度高于34.5)。WOD數據具有同樣的垂直結構,與夏季不同的是表層水的低溫特征以及觀測到溫度較高的夏季太平洋水。

從圖9位溫的垂直剖面上來看,Exp1、Exp3、Exp4、Exp5均模擬出實測數據中觀測到的夏季太平洋水、冬季太平洋水以及大西洋水的垂直結構特征。相比之下,Exp5模擬的結果與實測數據符合的較好,T-S圖上顯示出模擬的偏差主要在于冬季太平洋水溫偏高0.3℃,另外表現在夏季太平洋水的鹽度比實測高0.1~0.2。Exp1的模擬的大西洋水的溫度偏高0.3℃,冬季太平洋水偏高0.5℃,夏季太平洋水在T-S圖上的表現與Exp5相差不大,基本體現出應有的垂直特征。Exp3、Exp4雖然模擬出了垂向位溫結構,T-S圖中也體現出相應的躍層,但太平洋水的鹽度溫度與實測相比均有較大的偏差,中層大西洋水的溫度也比實測數據偏高0.5℃。

完全采用CFSR數據集的Exp2以及采用CFSR降水數據的Exp6并沒有模擬出合理的躍層結構,上層太平洋水的結構沒有反應出來,中層模擬的大西洋水超過1℃,與實測值有接近1℃的偏差。北冰洋中層水在近百年來雖然有增暖的趨勢[11],但在2006年11月至2008年11月加拿大海盆450 m深平均位溫并未超過1℃。T-S圖上均反映出Exp2、Exp6中大量的大西洋水進入北冰洋海盆,導致了中層鹽度、溫度的升高,同時影響到表層的鹽度。

5.2.3 弗拉姆海峽大西洋入流水的模擬

定義一個與等密度面相差很小的虛擬等密度面[30],該面上的地轉流可由Montgomery流函數梯度得到,我們采用與Aksenov等[31]類似的地轉流函數[32]:

(2)

式中,r定義為一虛擬等密度面,r=ρ[θ,S,p(z)]-ρ[θ0,S0,p(z)]+σ0(θ0,S0),σ0(θ0,S0)=ρ(θ0,S0,0)為位勢密度,θ0=0℃,S0=35,ρ0=1 027.5 kg/m3,p(z)可由靜力近似得到。這樣即有:

(3)

虛擬等密度面r與等位密面相差不大,在陸架海區的差別小于0.05 kg/m3,內區小于0.01 kg/m3。由弗拉姆海峽進入北冰洋的大西洋水(FSB)可用r=27.9 kg/m3來指示[31]。

我們考察各實驗模擬的FSB在北冰洋的分布。圖10是模式1990-2004年模擬結果的平均。根據地轉關系式(3),由Montgomery流函數的梯度可以判斷出環極邊界流流動的方向及強度。圖10g~k顯示出在27.9 kg/m3等密度面上的大西洋水由弗拉姆海峽進入北冰洋后,沿陸架外緣逆時針流動,在流動的過程中深度逐漸加深(圖10a~e),并與周圍水體不斷混合,熱量逐步損耗,溫度逐漸降低(圖10m~q)。這說明前5個實驗的模擬結果均大致符合目前對環極邊界流的一般認識。圖11a為各實驗模擬的弗拉姆海峽流通量結果與SODA再分析數據的計算結果。對比SODA數據,各實驗的結果均模擬出海峽通量與再分析海洋數據一致的季節性和年際變化特征。

采用CFSR數據集的Exp2模擬的弗拉姆海峽大西洋入流水在歐亞海盆和加拿大海盆比Exp1深120 m左右。在27.9 kg/m3等密度面上Exp2模擬的溫度在靠近太平洋一側也比Exp1高將近0.5℃。Montgomery流函數的梯度顯示出Exp2模擬的地轉流速比Exp1高約20%,在弗拉姆海峽的流通量也高出0.66×106m3/s(圖11b)。

Exp3和Exp1盡管使用了不同的開邊界條件,在27.9 kg/m3等密度面上Montgomery流函數、位溫和等密度面的深度相差不大。Exp3和Exp4的模擬結果整體上具有一致的形態,差別體現在Exp3模擬的等密度面的深度,在加拿大海盆達350 m以上的范圍比Exp4更廣(圖10c、d)。事實上,在第三部分我們比較兩個數據集的風應力旋度時發現,CFSR數據集在多年平均上比JRA25引入了一個正的渦度(圖4a、b),在Ekman抽吸下必然造成海洋中等密度面的抬升。除此之外,在歐亞海盆弗拉姆海峽入口處Exp3等密度面上的溫度比Exp4高0.2℃左右(圖10o、p)。多年月平均的結果顯示二者在弗拉姆海峽的流通量模擬結果相差并不明顯,只有0.036×106m3/s(圖11b)。這說明兩數據集中的風場數據對FSB模擬的影響不大,采用完整JRA25數據集的Exp1和CFSR數據集的Exp2流通量偏差的原因并不是由于風場造成的。

考察向下輻射通量對于FSB模擬結果的影響。對比Exp5和Exp1的模擬結果,Exp5模擬的北冰洋中FSB位溫在歐亞海盆、加拿大海盆明顯低于Exp1(圖10m、q),同時模擬的地轉流速也相對較弱。Exp5在弗拉姆海峽的流通量比Exp1低0.26×106m3/s(圖11b),與SODA數據的結果相對一致。

圖10 各實驗模擬27.9 kg/m3等密度面所在深度(a~f),其上的Montgomery流函數(g~l),位溫(m~r)結果在1990-2004年間的平均Fig.10 Simulated mean result of each experiment from 1990 to 2004,depth (a~f),montgomery streamfunction (g~l),potential temperature (m~r) on pseudo isopycnal surface r=27.9 kg/m3

圖11 各實驗模擬的弗拉姆海峽位于78.75°N斷面的體積通量(a)以及1990-2004年多年月平均的體積通量(b)Fig.11 Time series for volume flux (a) and monthly mean values averaged (b) between 1990 and 2004 of each experiment

與其他實驗的結果相比,Exp6采用了CFSR較高的降水,27.9 kg/m3等密度面的深度因此相對其他實驗中偏深(圖10f),Montgomery流函數也相對較小(圖10l)。Exp6與Exp2模擬的弗拉姆海峽流通量相對一致,比其他實驗以及SODA數據均偏高(圖11b)。這充分表明降水顯著影響到北大西洋水在弗拉姆海峽的凈通量,進而會影響到大西洋水在北冰洋中層的分布。

6 結論

通過比較兩個不同分辨率大氣再分析數據集,我們發現JRA25和CFSR數據集中用來驅動冰海耦合模式的6個變量在時間和空間尺度上有很大差別。時間序列上來看,CFSR的向下長波輻射、降水數據普遍高于JRA25;向下的短波輻射、比濕數據弱于JRA25數據集;而近地面氣溫在冬季高于JRA25,夏季相對較小;風速數據在冬季弱于JRA25,夏季強于JRA25。各變量在空間分布上的差異均不相同,值得注意的是,CFSR高精度的風場資料多年平均意義上在加拿大海盆上空加入了一個氣旋式的風場,削弱了波弗特高壓的影響強度;除此之外,還在巴倫支海引入了反氣旋的風場。在能量分配上,CFSR數據集中的各個變量在半年周期的尺度上強度均強于JRA25;在年周期尺度上除風場強于JRA25外,其他均相對較弱。總之,兩個數據集中降水偏差最大,其次為向下短波輻射數據、向下長波數據。

文章分別用CFSR、JRA25兩套再分析數據集強迫同一冰海耦合模式。對于海冰的模擬,二者在海冰外緣線、海冰冰量、冰厚上有巨大的偏差。產生這些偏差的首要因素是兩數據集中由于云量計算差別導致的向下長短波輻射不同,其次是降水通過熱力學過程和淡水通量到影響海冰融化凍結,而風場對于海盆尺度的海冰模擬影響較小。

此外,對于海洋的模擬,由于CFSR的向下長短波輻射與實測更接近,采用該數據的Exp5模擬的北冰洋中層水、加拿大海盆溫鹽結構與觀測符合的更好,對弗拉姆流通量的模擬更加接近海洋再分析數據SODA的結果。Exp2和Exp6模擬的溫鹽結構嚴重失真,模擬的FSB在北冰洋中的分布、等密度面上的地轉流以及弗拉姆海峽流通量與實測有很大偏差,文中發現CFSR中過量降水是其最主要的原因。因此更合理的表面輻射通量數據以及淡水通量數據對北冰洋冰海耦合模式的模擬十分必要。

[1] Dickson B. Oceanography: All change in the Arctic[J]. Nature,1999,397(6718): 389-391.

[2] Serreze M C,Holland M M,Stroeve J. Perspectives on the Arctic’s shrinking sea-ice cover[J]. Science,2007,315(5818): 1533-1536.

[3] Parkinson C L,Cavalieri D J,Gloersen P,et al. Arctic sea ice extents,areas,and trends,1978—1996[J]. Journal of Geophysical Research: Oceans(1978—2012),1999,104(C9): 20837-20856.

[4] Screen J A,Simmonds I. The central role of diminishing sea ice in recent Arctic temperature amplification[J]. Nature,2010,464(7293): 1334-1337.

[5] Serreze M C,Francis J A. The Arctic amplification debate[J]. Climatic Change,2006,76(3/4): 241-264.

[6] Holland M M,Bitz C M. Polar amplification of climate change in coupled models[J]. Climate Dynamics,2003,21(3/4): 221-232.

[7] Pithan F,Mauritsen T. Arctic amplification dominated by temperature feedbacks in contemporary climate models[J]. Nature Geoscience,2014,7:181-184.

[8] Rudels B,Friedrich H J,Quadfasel D. The Arctic circumpolar boundary current[J]. Deep-Sea Research Part Ⅱ: Topical Studies in Oceanography,1999,46(6): 1023-1062.

[9] Yang J. The Arctic and Subarctic Ocean flux of potential vorticity and the Arctic Ocean circulation[J]. Journal of Physical Oceanography,2005,35(12): 2387-2407.

[10] Karcher M,Kauker F,Gerdes R,et al. On the dynamics of Atlantic Water circulation in the Arctic Ocean[J]. Journal of Geophysical Research: Oceans (1978-2012),2007,112(C4).

[11] Polyakov I V,Alekseev G V,Timokhov L A,et al. Variability of the intermediate Atlantic water of the Arctic Ocean over the last 100 years[J]. Journal of Climate,2004,17(23): 4485-4497.

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

[13] Aksenov Y,Bacon S,Coward A C,et al. Polar outflow from the Arctic Ocean: A high resolution model study[J]. Journal of Marine Systems,2010,83(1): 14-37.

[14] Haynes P H,McIntyre M E. On the evolution of vorticity and potential vorticity in the presence of diabatic heating and frictional or other forces[J]. Journal of the Atmospheric Sciences,1987,44(5): 828-841.

[15] Nguyen A T,Menemenlis D,Kwok R. Arctic ice-ocean simulation with optimized model parameters: Approach and assessment[J]. Journal of Geophysical Research: Oceans (1978-2012),2011,116(C4).

[16] Onogi K,Tsutsui J,Koide H,et al. The JRA-25 reanalysis[J]. 気象集誌. 第2輯,2007,85(3): 369-432.

[17] Saha S,Moorthi S,Pan H L,et al. The NCEP climate forecast system reanalysis[J]. Bulletin of the American Meteorological Society,2010,91(8): 1015-1057.

[18] Iziomon M G,Mayer H,Matzarakis A. Downward atmospheric longwave irradiance under clear and cloudy skies: Measurement and parameterization[J]. Journal of Atmospheric and Solar-Terrestrial Physics,2003,65(10): 1107-1116.

[19] Spall M A. On the circulation of Atlantic Water in the Arctic Ocean[J]. Journal of Physical Oceanography,2013,43(11): 2352-2371.

[20] Zhao J,Cao Y,Shi J. Core region of Arctic Oscillation and the main atmospheric events impact on the Arctic[J]. Geophysical Research Letters,2006,33(22): L22708.

[21] Large W G,Pond S. Sensible and latent heat flux measurements over the ocean[J]. Journal of Physical Oceanography,1982,12(5): 464-482.

[22] Marshall D. Subduction of water masses in an eddying ocean[J]. Journal of Marine Research,1997,55(2): 201-222.

[23] Condron A,Renfrew I A. The impact of polar mesoscale storms on northeast Atlantic Ocean circulation[J]. Nature Geoscience,2013,6(1): 34-37.

[24] Lindsay R,Wensnahan M,Schweiger A,et al. Evaluation of seven different atmospheric reanalysis products in the Arctic[J]. Journal of Climate,2014,27(7): 2588-2606.

[25] Zhang J,Lindsay R,Steele M,et al. What drove the dramatic retreat of arctic sea ice during summer 2007?[J]. Geophysical Research Letters,2008,35(11):L11505.

[26] Hunke E C,Holland M M. Global atmospheric forcing data for Arctic ice-ocean modeling[J]. Journal of Geophysical Research: Oceans (1978-2012),2007,112(C4):C04S14.

[27] Nurser A J G,Bacon S. The Rossby radius in the Arctic Ocean[J]. Ocean Science,2014,10(6): 967-975.

[28] Lique C,Steele M. Where can we find a seasonal cycle of the Atlantic water temperature within the Arctic Basin?[J]. Journal of Geophysical Research: Oceans (1978-2012),2012,117(C3):C03026.

[29] Li X,Su J,Zhao J. An evaluation of the simulations of the Arctic Intermediate Water in climate models and reanalyses[J]. Acta Oceanologica Sinica,2014,33(12): 1-14.

[30] McDougall T J,Klocker A. An approximate geostrophic streamfunction for use in density surfaces[J]. Ocean Modelling,2010,32(3): 105-117.

[31] Aksenov Y,Ivanov V V,Nurser A J,et al. The Arctic circumpolar boundary current[J]. Journal of Geophysical Research: Oceans (1978-2012),2011,116(C9):C09017.

[32] Marshall J,Jamous D,Nilsson J. Entry,flux,and exit of potential vorticity in ocean circulation[J]. Journal of physical oceanography,2001,31(3): 777-789.

Analysis on the response of an Arctic ice-ocean coupled model to two different atmospheric reanalysis datasets

Mu Longjiang1,Zhao Jinping1

(1.KeyLaboratoryofPhysicalOceanography,MinistryofEducation,OceanUniversityofChina,Qingdao266100,China)

The downward radiative fluxes, wind speed, near surface temperature, precipitation, humidity got from Climate Forecast System Reanalysis (CFSR) and the Japanese 25-year Reanalysis Project (JRA25) are compared in this article. We find that most significant biases between the two datasets are precipitation, downward fluxes for both shortwave and longwave radiation. Driving by these two datasets, model results forced by CFSR shows big differences on sea ice, Atlantic water and thermohaline structure in Canada basin compared to in situ observations, with the simulated geostrophic current on isopycnal surface 20% higher than that forced by JRA25 and a larger volume fluxes than that derived from SODA data. Sensitivity experiment forced by downward radiative fluxes from CFSR, which have been evaluated to be close to observed values, demonstrates comparable results to observational results. The cloud data plays a key role in modeling sea ice while freshwater flux brought by precipitation can change the heat transport of Atlantic inflow prominently and carry a further effect on sea ice in the Arctic. The overestimated precipitation in CFSR is the major for large biases of volume flux through Fram Strait, geostropic current on isopycnal surface and thermohaline structure in central Arctic. Although reanalysis wind have different resolution for the two datasets, our results indicates that it carries an ignorable effect on modeling sea ice and thermohaline structure on basin scale.

Arctic Ocean; CFSR; JRA25; ice ocean coupled model; Canada Basin; Fram strait

2015-04-16;

2015-06-18。

國家重大科學研究計劃項目(2015CB953903);國家自然科學委員會重點基金項目(41330960)。

牟龍江(1989—),男,山東省日照市人,博士研究生,從事冰海耦合數值模型研究。E-mail:479689286@qq.com

10.3969/j.issn.0253-4193.2015.11.008

P732.6

A

0253-4193(2015)11-0079-13

牟龍江,趙進平. 北極冰海耦合模式對兩種不同大氣再分析資料響應的分析[J]. 海洋學報,2015,37(11): 79-91,

Mu Longjiang,Zhao Jinping. Analysis on the response of an Arctic ice-ocean coupled model to two different atmospheric reanalysis datasets[J]. Haiyang Xuebao,2015,37(11): 79-91,doi:10.3969/j.issn.0253-4193.2015.11.008

猜你喜歡
海冰實驗
記一次有趣的實驗
微型實驗里看“燃燒”
末次盛冰期以來巴倫支海-喀拉海古海洋環境及海冰研究進展
海洋通報(2021年3期)2021-08-14 02:20:38
近三十年以來熱帶大西洋增溫對南極西部冬季海冰變化的影響
做個怪怪長實驗
基于SIFT-SVM的北冰洋海冰識別研究
NO與NO2相互轉化實驗的改進
實踐十號上的19項實驗
太空探索(2016年5期)2016-07-12 15:17:55
累積海冰密集度及其在認識北極海冰快速變化的作用
應用MODIS數據監測河北省近海海域海冰
河北遙感(2014年4期)2014-07-10 13:54:59
主站蜘蛛池模板: 在线观看亚洲人成网站| 午夜精品久久久久久久2023| 成人午夜精品一级毛片| 成年女人18毛片毛片免费| 91久久天天躁狠狠躁夜夜| 福利在线一区| 亚洲国产看片基地久久1024| 丁香五月激情图片| 中文字幕在线免费看| 丁香五月激情图片| 国产高清在线精品一区二区三区 | 成人综合久久综合| 少妇精品在线| 草草影院国产第一页| 亚洲va在线∨a天堂va欧美va| 日韩精品亚洲精品第一页| 亚洲精品手机在线| 婷婷久久综合九色综合88| 在线看片中文字幕| 97亚洲色综久久精品| 一本大道香蕉久中文在线播放| 91精品伊人久久大香线蕉| 天天综合网色| 日韩毛片免费| 在线看AV天堂| 国产不卡在线看| 在线免费不卡视频| 在线日本国产成人免费的| 日韩在线视频网| 无码'专区第一页| 欧美怡红院视频一区二区三区| 国产成人91精品| av在线无码浏览| 国产高清无码麻豆精品| 成人久久18免费网站| 亚洲精品动漫| 色成人亚洲| 婷婷丁香在线观看| 免费人成又黄又爽的视频网站| 日韩精品资源| 亚洲三级网站| 91福利免费| 久久99国产乱子伦精品免| 女人一级毛片| 国产在线无码一区二区三区| 久久久久久高潮白浆| 国产在线观看91精品| 天堂成人在线| 97在线免费视频| 午夜小视频在线| 伊人成人在线| 狠狠色狠狠综合久久| 又污又黄又无遮挡网站| 国产精品亚洲天堂| 91在线一9|永久视频在线| 欧美日韩亚洲国产主播第一区| 国产美女视频黄a视频全免费网站| 伊人天堂网| 国产精品刺激对白在线| 国产9191精品免费观看| 激情乱人伦| 亚洲中文字幕av无码区| 色综合天天综合| 在线观看国产黄色| 欧美在线视频不卡第一页| 亚洲欧美在线精品一区二区| 成人午夜网址| 久久青草免费91线频观看不卡| 国产18页| 亚洲日韩精品无码专区97| 国产女人18水真多毛片18精品 | 日韩国产欧美精品在线| 亚洲中文字幕手机在线第一页| 91久久偷偷做嫩草影院精品| 欧美中日韩在线| 永久免费av网站可以直接看的 | 国产精品99久久久久久董美香| 波多野结衣的av一区二区三区| 99热这里只有精品免费| 精品無碼一區在線觀看 | 凹凸精品免费精品视频| 精品国产成人av免费|