基于拉格朗日擬序結構的黑潮區域海流運輸結構提取和分析

田豐林，任一丹，何遒，陳戈*

( 1. 中國海洋大學 信息科學與工程學院，山東 青島 266100；2. 青島海洋科學與技術試點國家實驗室 區域海洋動力學與數值模擬功能實驗室，山東 青島 266237)

1 引言

黑潮是西太平洋海域最大的暖流系統，它起源于呂宋島以東的洋面，其主干流沿臺灣以東進入東海，順東海大陸坡向東北流去，流經日本列島南部，于35°N，141°E附近海域結束，后續成為黑潮延伸體。黑潮是一種典型的西邊界流，具有高溫、高鹽等特征，對其周邊海區的氣象、水文產生巨大影響。黑潮流量以及流幅存在明顯的季節變化和年際變化，黑潮流量夏季大，秋季小，流幅遵循春夏寬，秋冬窄的變化規律[1]。目前對于黑潮的研究工作中，研究范圍多為整個流系中的某一局部區域，如選取黑潮延伸體[2-4]或者臺灣東北部地區[5]，有研究考察了流軸的形態變化[6-7]，黑潮流軸多采用局部最強速度場等屬性等值線方法獲取[8]，也有利用模型驗證模擬黑潮的動力學過程[9]，對黑潮入侵南海路徑和機理做了研究[10]。黑潮不斷與周邊水域發生強烈的物質交換和能量輸送，目前關于黑潮區域拉格朗日擬序結構（LCSs）的研究較少，LCSs作為研究海洋流場結構的有力工具，能夠獲取通過流線和軌跡信息等無法顯現的物質運輸特征，對于研究黑潮環流結構、水體性質等有著重要意義。

LCSs最早在2000年時被提出[11]，它可以精確地定義和識別具有任意時間依賴性的動力系統中的障礙[12]，是流場表面存在的一種骨架結構。利用拉格朗日的方法提取相干結構方法眾多，定義LCSs的方法可以有所不同，但均應滿足4個條件：客觀性、有限時間性、拉格朗日不變性、空間收斂性[13-14]。有研究對這些方法做了對比，探討了不同方法得到的結果[15-16]，其中有限時間里雅普諾夫指數（FTLE）和有限尺寸里雅普諾夫指數（FSLE）是應用較為廣泛的兩種LCSs方法[17]，FTLE描述了初始時刻相鄰的流體質點在有限時間范圍里的平均分離率，FSLE表示在有限距離里的平均分離率。d'Ovidio等[18]利用FSLE研究地中海區域的混合結構識別，黃高龍等[19]利用FTLE分析呂宋海峽附近環流模式和水體輸運，解釋浮標運動軌跡的變化。更進一步的研究中Haller[20]及Farazmand和Haller[21]提出了變分理論來計算LCSs的方法，LCSs表現為一些具有不同收縮和拉伸特性的物質線，可以用來區分流場中具有不同動力學特征的區域，根據流體內部物質受力情況不同，分為橢圓型、雙曲型和拋物線型LCSs，橢圓型LCSs是拉格朗日渦旋的邊界，雙曲型LCSs代表了物質收縮或拉伸線，拋物線型LCSs可以作為射流核心[22-23]，雙曲型LCSs和拋物型LCSs在理論流場中的表現形式如圖1所示。近些年LCSs在各種海洋現象的研究上取得了重大進展，LCSs在渦旋識別方面已有較為成熟的理論，在渦旋對水體輸運效果評估方面也得到應用[24-25]，在海流方面多是對緯向射流進行分析[26-28]。Olascoaga等[29]將LCSs應用于墨西哥灣地區，表明了LCSs對該地區內部署的浮標的運動施加了重要的約束，發現吸引型LCSs近似于浮標運動形成的中尺度圖案的中心，而這種物質運輸模式對于傳統的歐拉分析來說是非常困難的。

隨著全球海洋觀測技術日趨成熟，能夠獲取到更高質量的網格化海洋數據集，這些數據具有高時空分辨率、更長時間跨度等特性，將這些數據用于黑潮的研究，可以更加精確和直觀的表達出海流信息。本文創新性地利用變分方法計算了黑潮區域的LCSs，并且將得到的LCSs簡化合并為LCSs束，基于長時間序列的海流數據得到了清晰的氣候態黑潮運輸模式。

2 數據和方法

2.1 數據介紹

本次采用的數據是絕對動力地勢（Absolute Dynamic Topography & Absolute Geostrophic Velocities，ADT）[30]，是利用在給定時間內可用的所有衛星計算出來的一種網格狀海面高度和速度數據，計算方法為ADT=SLA+MDT，其中，SLA為海平面高度異常，MDT為平均動力地勢。ADT數據質量高，能夠突出海洋環流細節，反映出典型的海流特征，對于研究黑潮系統的結構十分有利。該數據的空間分辨率為0.25°×0.25°，時間分辨率為1 d。本研究采用了25年的數據，時間范圍為1993-2017年。

我們獲取了用于驗證工作的浮標軌跡數據，數據來自NOAA的全球漂流浮標計劃（Global Drifter Program）[31]，選取了全球范圍內1993-2017年共計25年的浮標數據，按月份提取了浮標軌跡信息。

用于驗證的溫鹽數據是來自CMEMS發布的ARMOR3D數據[32]，是一種結合了衛星觀測數據和現場數據而得到的數據產品，數據空間分辨率為0.25°×0.25°，采用的時間范圍也是1993-2017年。

2.2 研究方法

2.2.1 計算LCSs的變分方法

在二維不穩定流場中，速度場可以表示如下：

柯西格林應變張量場可以表示為[21]

圖1 雙曲線型 LCSs（a）和拋物線型 LCSs（b）Fig. 1 Hyperbolic LCSs（a）and parabolic LCSs（b）

應變線處處與單位應變特征向量場n1相切，可以根據下式中的常微分方程計算得到平滑的軌跡

為使此方程的數值解為一個全局平滑的應變線集，依據一般張量線的尺度變換方式，將式（5）進一步表示為[21]

利用變分方法從速度場中計算得到LCSs是一些特殊的物質線，一條物質線處處切于時為收縮線，處處切于時為拉伸線[23,30]。對于，當后向積分（T＜0）時，吸引型 LCSs定義為局部最大收縮線，排斥型LCSs定義為局部最大拉伸線，前向積分（T＞0）時，吸引型LCSs定義為局部最大拉伸線，排斥型LCSs定義為局部最大收縮線。圖1所示的雙曲型LCSs，中心橫線代表的是拉伸線，中心豎線代表的是收縮線，虛線表示周圍的流形，它們將流場劃分為4個具有不同運輸特征的區域，可見收縮線作為流場中的運輸屏障，阻止了其兩邊的水體和物質交換[23,33]。

2.2.2 計算氣候態LCSs

在氣候態下研究LCSs具有更廣泛的意義，能夠反映當前海流區域的綜合特征和集中趨勢，為不同月份和不同季節提供對比分析。我們利用1993-2017年總共25年的速度場數據，通過年際求平均得到氣候態海流速度場，用計算得到的氣候態海流速度場做進一步的研究。為方便計算，取每年前360 d的數據，每個月按照30 d計算，1年共計12個月[34]。

利用密集流體粒子在流場中運動一段時間后計算得到柯西格林應變張量場，張量場描述了粒子在當前位置下的運動趨勢和形變特性。我們在整個區域內均勻選取N×N個網格點，在每一個網格點周圍放置4個距離較近的輔助點，如圖2所示。本文中輔助點與主網格點之間的距離δx1、δx2設定為兩個網格點之間距離的0.05倍，使用龍格-庫塔四階積分計算輔助點處的粒子運動后的位置，使用有限差分方法計算出每個網格點處的形變梯度矩陣，得到柯西格林應變張量，進而計算特征值和特征向量，積分得到LCSs。

圖2 網格點與輔助點的位置Fig. 2 Locations of grid points and auxiliary points

本文采用后向積分的方法，將吸引型LCSs作為研究對象，從每個月里選取第30天代表計算LCSs的初始時刻，積分時間為30 d，計算張量場，基于此張量場計算得到的12個不同動力系統的LCSs，分別代表了各月的氣候態LCSs。具體步驟如下：

（1）生成均勻密集網格點，求每點處的柯西格林應變張量，并計算特征值和特征向量；

（3）將所有極值點按特征值大小遞減排序；

（4）以特征值極值點為種子點，在對應的特征矢量場中積分，得到吸引型LCSs。

通過調整尋找極值點的范圍和積分線的密集程度等參數可以保留部分吸引型LCSs，如圖3所示，得到了代表12個月份的氣候態LCSs，圖中的線代表了求取的吸引型LCSs，表示其法向排斥率，我們將線的顏色用ln值劃分。利用FTLE計算得到的LCSs體現為標量場的脊，在有限時間[]內，FTLE定義如下：

我們將獲取吸引型LCSs與FTLE場進行疊加，可以發現抽取到的吸引型LCSs與FTLE的脊基本相吻合，在這種條件下該吸引型LCSs可以作為黑潮兩側邊界。

2.2.3 提取 LCSs束

由于過少的LCSs無法全面地展示海流結構，過多的LCSs容易引起視覺雜亂，我們獲取整個研究區域內的全部吸引型LCSs，如圖4a所示，當LCSs足夠密集時，把積分得到的所有吸引型LCSs進行簡化合并。保留前20%和40%的吸引型LCSs（圖4b，圖4c），觀察其覆蓋的范圍，發現黑潮區域的特征值要明顯大于周邊區域，選取的吸引型LCSs不僅要表現出黑潮區域的整體結構，還摒棄大部分周圍海流結構的干擾。

圖3 氣候態下的吸引型LCSs和FTLE場Fig. 3 Monthly climatological attracting LCSs and FTLE fields

圖4 保留 100%（a），40%（b）和 20%（c）的吸引型 LCSsFig. 4 Reserved attracting LCSs for 100% (a), 40% (b) and 20% (c)

選取前40%的吸引型LCSs進行研究，經過簡化合并后的LCSs我們稱之為LCSs束，具體合并流程以及示例如圖5所示，對于每一條吸引型LCSs，計算它與區域內其他線的相似度，將相似度最高的兩條線合并為一條線。在對任意的兩條線計算相似度時，由于線是由有序排列的頂點構成，對于一條線中的每個頂點坐標，去尋找另一條線中距離最近的點進行配對并計算距離，最終得到以該條線為基準的平均距離，用同樣的方法得到以另一條線為基準的平均距離，對這兩個距離求取平均，得到的平均值作為這兩條線相似度的數值描述[35]。通過試驗設定容許兩線合并的最大相似度S=1°，將兩條線按照逐點對取中間點或者取較長線的方法進行合并，隨著保留的線越來越少，計算復雜度逐漸降低，最終能夠得到反映流場突出特征的更為簡潔的LCSs束。

圖5 LCSs簡化合并流程圖（a）和相似度表示例（b，c）Fig. 5 Flow chart of simplifing and merging LCSs (a) and an example of similarity table (b, c)

將吸引型LCSs進行簡化合并后，得到的效果如圖6所示，合并的線的條數越多，LCSs束在圖中表現得越粗，原來密集的LCSs經過合并后能夠清晰地表現出特征較強的流場區域，反映出黑潮的主軸結構，該方法能夠保持大范圍流場中突出的海洋特性，也有利于海流結構的二維可視化研究。

3 結果與討論

3.1 結果分析

用氣候態速度場計算得到的吸引型LCSs，在黑潮主干區域能夠得到較連貫黑潮流結構，在黑潮延伸體區域，雖然不同時間的流場變化幅度大，復雜且不穩定，但是通過氣候態數據的計算也可以反映出典型的運輸特性。圖7是提取到的黑潮系統的流場結構圖，揭示了12個月份的拉格朗日運輸模式。

目前，黑潮流軸和路徑的獲取方法多種多樣，最為常用的有以下3種方法：其一是基于流線的識別方法，黑潮表面流速最大的部分看作是表面黑潮流軸；其二是基于等溫線的識別方法，黑潮海水的溫度顯著高于周圍海流，日本南部海域處的主軸采用 16℃等溫線進行識別，而黑潮延伸體處的主軸則用12℃等溫線作為識別依據；其三是基于高度數據的識別方法，例如對于ADT數據，ADT=240 cm處的等值線可以作為黑潮主軸[36]。圖8展示了利用流速和海表面高度識別的黑潮流軸。

圖6 合并前的LCSs (a)和合并后的LCSs束(b)示意圖Fig. 6 LCSs before merging (a) and LCSs bundles after merging (b)

圖7 12個月的氣候態流場結構Fig. 7 Monthly climatological flow field structures

圖8 基于多年平均流速提取的黑潮主軸[8]（a）和 170 cm海表面高度輪廓線定義的黑潮延伸體路徑[4]（b）Fig. 8 Kuroshio axis extracted based on average velocity[8] (a) and paths of the Kuroshio Extension defined by the 170 cm contours in sea surface height fields[4] (b)

本文提出了一種用LCSs束識別黑潮流路的方法，LCSs束本身代表了兩側海流具有較強的排斥作用，LCSs束可以作為提取黑潮主軸的有力工具，將其應用于海洋流場中具有重要價值。從圖7中可以發現黑潮主干結構穩定，但在局部地區存在月際變化。在流經呂宋海峽時，產生不連續的隔離現象，這可能是因為黑潮入侵南海所導致。在黑潮主干區域，冬季提取到的流結構較為緊密和穩健，可能是由于冬季黑潮流幅窄而平穩的緣故。黑潮在經過臺灣東北部時由于海底地形的作用產生向東的轉折，分出一支臺灣暖流，從圖7中發現春季在這一地區的轉折相較于其他時期更為平緩一些，而后主干沿著大陸架的約束繼續向東北方向前進，直到達到日本九州島西南部時又出現轉折，在此之前一直存在著穩定和強健的吸引型LCSs。根據圖7中線的顏色可以看出在轉折區域多會發現特征強度較弱的信號，可能是因為產生流系分支的原因，這些區域的流軸不是很健壯。由于黑潮延伸體區域流場的復雜多變性，提取到的LCSs束的連續性較弱，但是從圖7中可以發現延伸體在143°～150°E之間，不同月份有著相似彎曲形態，均在146°E處形成谷狀，且冬季彎曲幅度更大。黑潮延伸體區域存在中尺度及亞中尺度渦旋等復雜流場結構，海流在流動過程中存在較強的物質交換，使得海流邊界更為曲折復雜，難以提取到完整的結構。

3.2 驗證

為了評估以上結果的可靠性，我們對結果進行了驗證，驗證方法包括拉格朗日粒子運輸驗證、浮標軌跡驗證以及溫度鹽度驗證。

3.2.1 拉格朗日粒子的運輸

采用1月份的數據進行拉格朗日粒子運輸的驗證，通過示蹤粒子在不穩定流場中運動一段時間后觀察其是否與LCSs束吻合。選定6個不同位置且間距合理的區域，在直徑為1°的圓形區域內均勻地放置示蹤粒子，它們的初始位置和形態如圖9a所示，讓粒子在1月的不穩定流場中分別運動10 d、20 d和30 d，得到的結果如圖9b至圖9d所示?？梢园l現粒子輸送總是沿著LCSs的方向拉伸，初始時刻的圓形組合隨著時間的增加逐漸變為回旋鏢狀，最后拉伸成絲狀，證明LCSs束可以代表黑潮的主軸結構，也可預測物質運輸的形態。

3.2.2 浮標軌跡

我們通過1993-2017年的浮標軌跡信息對結果進行驗證，提取了25年中1月和6月的全部浮標軌跡，與相應月份的LCSs束進行比較，結果如圖10所示，浮標的運輸模式符合LCSs束的預測情況。雖然單獨的浮標運動由于受局部因素影響不確定性大，軌跡較為曲折和雜亂，但從整體上看它們的軌跡分布基本覆蓋了整個黑潮流域，浮標沿著黑潮蜿蜒前行，在臺灣島東北部、九州島南部等區域發生轉折，在黑潮主干區域與LCSs束的提取結果保持了高度的一致性，在黑潮延伸體區域也存在部分相似性。

3.2.3 溫鹽分布

此前相關研究多采用溫鹽等值線來確定黑潮的位置[37]，基于這種方法我們分析了溫鹽分布與提取到的流場結構的相關性。采用海表面和50 m深的溫度和鹽度數據，得到了25年的1月份溫鹽數據的平均值，計算了溫鹽異常，將其與1月份提取到的LCSs束疊加，結果如圖11所示。從溫度異常圖中可以看出，黑潮北側的溫度異常值變化明顯，LCSs束有效區分了不同溫度的區域。鹽度異常圖中，黑潮主干大部分地區的LCSs束存在于鹽度異常變化相對較大的區域?？傮w上看，黑潮的北側溫鹽異常梯度大，邊界明顯，南側溫鹽邊界較為模糊，與提取的LCSs束具有一致性。而在黑潮延伸體區域，溫鹽異常分布與LCSs束的相似性較小，這部分區域是海洋渦旋活躍地帶，物質交換劇烈，難以獲得精確結果。

圖9 示蹤粒子初始位置（a）和運動 10 d（b）、20 d（c）、30 d（d）后的 LCSs輸運圖Fig. 9 Initial position (a) and advected images for 10 days (b), 20 days (c) and 30 days (d) of particles

圖10 1月（a）和 6 月（b）的浮標軌跡Fig. 10 Trajectories of drifters for January (a) and June (b)

圖11 氣候態溫鹽異常分布Fig. 11 Distributions of climatological temperature and salt anomalies

4 總結

本文研究了黑潮區域的拉格朗日運輸模式，面對復雜多變的流場環境，提出了利用氣候態LCSs束來計算流場結構的方法。整合了25年的速度場數據，通過變分方法計算得到了氣候態下的吸引型LCSs，并根據特征值的強度選取了黑潮附近的LCSs，將其簡化合并為LCSs束。LCSs束可以作為描述黑潮形態和流場結構的工具，用于研究黑潮的時間和空間變化規律等。我們提取并分析了12個月份的氣候態流場結構圖，討論了不同月份和季節中物質運輸的一致性和區別，最后用拉格朗日虛擬粒子、真實浮標軌跡以及溫度鹽度數據作為驗證，表明了其合理性和準確性。

同時本文也存在諸多不足，首先我們只關注了黑潮及其延伸體的主干部分，未來可以考慮黑潮分支的提取，如臺灣暖流、對馬暖流等，將研究拓展到完整的黑潮流系結構。其次，在黑潮延伸體區域，可以結合中尺度渦旋的位置、特性以及其他海洋現象展開分析，進而獲得更為連續和完整的海流結構。同時，可以將本方法用于其他大范圍的海域，利用并行方案來解決大數據問題，對全球海流的給定精確描述，為物質運輸狀態做指導，應用于海上溢油、海上垃圾堆積和擴散區域的識別和預測。